5 Nonlinear Interaction Problems
COMPUTATIONAL METHOD FOR SOIL/STRUCTURE INTERACTION PROBLEMS?
Y. MARVIN i~of and RUSSIZLLH. ENGLAND8 California Research & Technology, Inc.. Woodland Hill, California, U.S.A. and
RICHARDB. NELSO~~ ljniversity of California. Los Angeles, California, U.S.A. (Receiued 9 h&s 1980) Abstract -This paper describes a ~rnpu~tio~l method for soil/structure interaction problems, with particular emphasis on reinforced concrete structural response. The Computational method involves the modeling capabilities of the CRT/NONSAP finite-element computer code which has been extended to address soil/structure interaction problems. Some of the relevant improvements include: *Refined elastic-plastic material models for concrete and geologic materials. llnttrface element to model shear transfer and sliding at soil/structure interfaces. *Element deletion (removal) procedure for representation of cracking in brittle materials. or failure in ductile material caused by large strains. The details of these improvements regarding their applicability to dynamic analysis of rdnforced concrete structures and the effects of ~il/~ru~~ int~~on are described. Atso, the nonlinear response analysis of a buried reinforced concrete cylinder to an air blast environment is presented.
Problems
INTRODUCTION involving nonlinear response of soil/structure
systems are complicated by three features: (1) nonlinearity in the structural components. e.g. ~~asti~ty in concrete and steel reinforcements of reinforced concrete structures; (2) nonlinearity in soils, e.g. compaction and (3) slip and/or gap at the interface between the soil and structure. In this paper a nonlinear fmite-element computer program developed at California Research & Technology, Inc. (CRT) is used to anaiyze problems of this type. This code, termed CRT/ NONSAP, is an extensively modified version of the code originally developed at the University of CaIifornia, Berkeley [l. 21. and has the capability of dealing with elastic-plastic large deformation dynamical response of complex structure/media systems. The nonlinear material models for soil and structure now available in this code include elasti~plastic models for complex pressure sensitive materials. These models are generalized versions of Drucker-Prager and’Mohr-Coulomb associated flow models [3,4], and can represent soil, rock or concrete types of materials, including. if necessary, tension Failures. cracking, effects of compaction and shear failure dependency on pressure. These models have successfully predicted the nonlinear response of geologic media and reinforced concrete structures to dynamical blast and shock loading. tThis work was partially supported by the Defense N&ear Agency. $Principal Engmeer. &search Engineer. Vrofessor of Engineering and Applied Science.
The combined soil/structure problem is complicated by the interface between the soil and structure. In fact, experiments [SJ have shown that sliding between concrete/soil and steel/soil interfaces occurs at interface shear levels which are significantly less than the shear limit of the soil. Thus, analyses which assume a perfectly bonded interface condition overpredict the shear transfer and depending on the specific application, either over- or underestimate structural response. Typical interface data indicates that shear associated with interface slipping resembles a Coulomb law. That is. the shear limit for bonding is dependent in some fashion on the interface normal pressure. Various plasticity ~n~tu~ve models exist which fit the general form of the interface shear/normal sliding stress data. For example, the Drucker-Prager and Mohr-Coulomb yield criterion give a good approximate representation of interface sliding data. Interfaoe behavior can be simulated numerically by explicitly representing the interface boundary with a Unite element having a small thickness and using a constitutive model that reflects the interface sh~r/no~al stress limits. Since the interface element is to simulate the shear/pressure limits on the surface between two dissimilar materials, its material constituents should reflect interface sliding data and not necessarily the constituent laws for either of the two adjoining materials. Theoretical calculations performed wherein interface behavior is simulated using either a Drucker-Prager or Mohr-Coulomb yield law and an associated flow rule generally do not compare well with experimental data. This is due to the fact that the sliding which occurs in an interface typically results in large interface shear strain. According to the associated flow rule this shear strain 157
158
Y.
MARVINho et al
causes volumetric expansion in the interface element of the order of the shear strains to occur during the process of interface slipping. This volumetric expansion can cause significant, unrealistic buildup of pressure in the structure and surrounding media. especially if the problem is confined. To avoid this difficulty, a special interface element was developed. The yield law employed is a simple Coulomb friction law, and the plastic flow rule used assumes the direction of plastic strain increment lies colinear with the interface sliding, i.e. sliding causes plastic shear strains only. This flow rule is nonassociated and prevents plastic dilatation. The selection of a nonassociated flow rule to alleviate plastic volumetric dilatation results in a loss of symmetry in the incremental constitutive law (incremental stress-strain matrix) and hence leads to nonsymmetric interface tangent stiffness matrices. Furthermore, nonassociated flow rules can lead to indefinite or negative definite systems. An incremental solution procedure is employed in which only the symmetric portion of the interface constitutive law is used to represent the system incremental tangent stiffnesses. The errors induced by neglecting the nonsymmetric terms are eliminated through iteration and the accumulation of errors is prevented by reapplication of the errors obtained from the previous cycle to the next load step. This approach is computationally efficient and permits the retention of the usual symmetric linear equation solution procedure. In the following sections the detailed mathematical features of the numerical modeling and solution techniques are described. Also, the nonlinear response analysis of a buried reinforced concrete cylinder to an air blast environment is presented.
Inelastic deformation in geologic media and concrete materials results from two seg8rate mechanisms: (1) plastic strains which accumulate when slippage occurs due to excess shear stress and (2) compaction (the cldsing of gaps and voids) which occurs under hydrostatic states of compression. Both these mechanisms can absorb significant energy, thus modeling their effects is essential if theoretical nonlinear models are to accurately predict behavior in soil/structure interaction problems. In the analyses performed herein, the geologic material is characterized using a hysteretic pressure dependent elastioplastic model. With this model, shear failure is governed by a pressure/plastic work dependent yield law and its associated flow law. The general form of the law in terms of the stress invariants J, and J; is as follows: (1)
where r=a(J,, K=K(Jl,
w,,
HYDROSTATIC
PRESSURE
Fig. 1. Typical shear falltire data for geologx materials..
The constitutive model represents the hysteretic effects of compaction using a variable loading and unloading bulk moduli formulation where the incremental bulk modulus depends both on the current volumetric strain state and the peak volumetric strain occurring during the loading history. Figure 2 shows typical compaction behavior generally obtained from either uniaxial strain or hydrostatic compressions tests which are used to fit the variable butk modulus model. The effect of compaction is most -significant in soils and may be ignored for certain types of rock and concrete. Similarly, plastic work has little signilicance on the yield/flow surfaces in soils. Thus this general material model can be specialized to represent various soils, rock and other geologic material such as concrete.
-7
GENERALIZED DRUCKER-PRAGER MATRRlAL MODEL
F=ctJ,+,/J;-K
function of the hydrostatic pressure ( -J, ) permits an accurate representation of triaxial compression test data. Typical shear failure data for geologic materials are shown in Fig. 1.
(24
W,,
J, =~a,
W, denotes plastic work, and a and J are scalar parameters. Expressing the material parameters a and K as a
DILATATIONAL
STRAIN
Fig. 2. Typical compaction data for geoloepc matenals. REINFORCED The
CONCRETE
MATERIAL
MODEL
reinforced concrete model developed for the numerical calculations presented herein includes an explicit representation of the steel reinforcement bars (one-dimensional elastic-plastic elements) and the plain concrete (elastic-plastic continuum elements). Plain concrete is idealized using an isotropic pressure dependent yield law with strain hardening and softening capabilities and an associated Row law. The yield criterion has the form given in eqn (1). Expressing the yield parameters x and K in terms of J, (hydrostatic pressure) and plastic work ( W,)allows for an analytical representation of both the pressure dependency on yielding and post-yield strain hardening and softening Formulations. The material parameters for the concrete model were evaluated to fit concrete biaxial stress data. This stress data was preferred to triaxial because of the basic
Computational
method for soil/structure interaction problems
of the structures analyzed (thin wall cylinder) where states of stress remain essentially biaxial (a,,+~,). Figure 3 shows a typical failure surface for plain concrete in states of plane stress. This biaxial failure data can be used to evaluate corresponding values of the invariants .I, and ,/J;; thus when the failure data is mapped to J,, ,/J; space (see Fig. 4), the material parameters a and K can be evaluated to fit experimental data. geometry
q
1.59
plastic flow rule from which the incremental constitutive laws are derived. Experimental data suggests that sliding in structure/ soil interfaces in governed by a Coulomb-type criterion, i.e. the critical shear stress which causes debonding is dependent on the normal interface pressure. Furthermore, experimental data shows that the critical shear stress limit is essentially a linear function of the normal pressure for a rather large range of interface pressures. A general form of the criterion fpr interface slip which defines all combinations of stress states corresponding to sliding is as follows: F(a,, r)=O
(4)
where u,=interface normal stress; ~=interface shear stress; and F(cJ,, r)=O denotes interface sliding, and F(u,, 7)
n
-J1/3 HYDROSTATIC
PRESSURE
Fig. 3. Biaxial failure data for plain concrete. -u2
Fig. 4. Plain concrete biaxial failure data mapped to J,, 4’5;; space. This concrete constitutwe model gives an accurate representation of concrete behavior until either compressive crushing and spalling or tensile cracking occurs. In this event an element deletion (removal) procedure has been developed where at any cycle of the solution procedure. failed elements are physically removed from the grid. Two distinct approaches can be exercised: (1) removal of both the mass and stiffness of the element (simulating the complete removal of the material), or (2) removal of the stiffness but retaining the mass (idealization of a crack section). This element removal technique is by no means an attempt to simulate crack propagation. It is. however, an effective means for evaluating the effects on structural response caused by changes in load path resulting from localized failures occurring in concrete (and steel reinforcing bars).
(5)
where p = interface coefficient of friction and C = interface cohesion. During interface sliding a flow rule is required which specifies the direction of the plastic flow given the interface stresses. In order to prevent plastic dilatation the following nonassociated flow law was selected: de;=0
(W
dyP=d;‘P
(6b)
where de: = incremental interface plastic normal strain and dyP= incremental interface plastic shear strain. Again, an associated flow rule will produce considerable volumetric expansion if sliding occurs on the interface. Figure 5 illustrates the difference between interface element response for associated and nonassociated flow rules. From the definition of the criterion for slip (i.e. that during sliding F ~0) and the plastic flow rule, the incremental constitutive relations are obtained. During sliding, the sliding criterion requires F(a,+du,
T+dt)=O
(7a)
or ^ zdu,+$dr=O n or (7c)
STRUCTUREiMEDlA INTERFACE MATERIAL MODEL
This section describes the basic fundamental relations used to evaluate the incremental constitutive equations for the structure/media interface element. This includes the stress criterion for slippage and the
Fig. 5. Effect of flow rule on distortion of interface element during sliding.
160
Y MARVIN ITO et al
The stress increments are evaluated from the eiastic strain increments and the elastic constitutive relations dr = G(dy -d,+Y
@a)
da, = &de, -de;).
(8))
Use of the nonassociated flow rule which stipulates that sliding causes plastic shear strain only (i.e. ds[= 0) and eqns (8) and (7) leads to the following expression for d{C ‘7 E (9) Finally, the incremental constitutive matrix is determined by substituting eqn (9) into (8),
incremental tangent stiffness and account for the errors induced by this approx~ation either through an iterative mod&d Ne~o~Raphson approach and/or by reapplication of the unbalanced loads which result from neglecting the nonsymmetric terms and from linearization of a nonlinear system to the next load step. This solution approach was used for the calculation presented herein and was found to give excellent results, i.e. this approach gave very stable behavior and generally, satisfactory solutions could be obtained without iteration. The steps outlined below illustrate the basic procedures employed to evaluate the interface stress state and the approximate constitunve laws used in evaluating the interface tangent stiffness matrix:
(10) where (11) Note that the elastic-plastic constitutive law is nonsymmetric and indefinite. Finally, it is recognized that positive normal strain (tension on the interface is an indication of a gap forming between the two material boundaries. When this occurs. a free surface exists on both material surfaces and the materials are no longer coupled through the interlace. This is simulated by zeroing the interface stress and the constitutive laws. NUMERICAL SOLUTION PROCEDURE
Perhaps the most versatile and preferred implicit solution procedure for nonlinear problems is an incremental technique wherein effects of nonlinearity are taken into account by means of a series of piece-wise linear analyses. Each analysis is based through linearization of the nonlinear constitutive law and straindisplacement relations (for large displacement theory) which gives a linear system incremental tangent stiffness for use in the-analysis. Errors induced from linearization can be eliminated by an iteration scheme in which the errors in system equilibrium are reapplied using the previous tangent stillness to reduce errors in the equilibrium equations (modified Newton-Raphson method). If iteration is not performed, then the errors in system ~~lib~urn can be reapplied to the next load increments to prevent the accumulation of error buildup. When nonlinear material behavior is idealized using plasticity models with associated flow rules, the system tangent stiffness is symmetric even for large displacement theory, provided the external loading is conservative. The symmetry of the system tangent stiffness is advantageous since it greatly reduces storage requirements and linear equation solution. If nonassociated flow rules are required linearization leads to nonsymmetric constitutive laws and nonsymmetric tangent stiffness matrices. Storage requirements for nonsymmetric matrices are nearly double the requirements for the symmetric stiffness and solution times may increase by an order of magnitude due to additional data management problems in large systems. In many appfications it is ~mpu~tion~ly more ethcient to neglect the nons~met~~ portion of the
t 1JFrom the previous interface normal strain and the current incremental normal strain. determme closure of the interface. la) If is, +A&,) O (open interface), zero out interface normal and shear stresses. If interface tangent stiffness is to be calculated. set ail terms in constitutive matrix to zero and return. (2) From previous interface strains, stresses, and the current incremental strain, calculate the updated interface stress state from elastic analysis: AQ, = EAsE, Ar=GAy. (3) If F(o,+-AC,,. r+Ar)
If interface tangent stiffness is to be calculated, select the elastic constitutive matrix and return. (a) If Ffa, + A@,,t + At) > 0, and if previous cycle was plastic, then set RATIO to zero. Otherwise, a transttion occurs between elastic and plastic response. In this case determine RATIO from the equation F(a,i- RATIO Aa,, r+ RATIO A+=0 Update the stresses to obtain the stress state prior to yielding CS,= a,, + RATIO Aa, r=s+RATIO Ar. (b) Use the remainmg portion of the strain increment associated with an elastic-plastic response to evaluate the corresponding stress increment. Since the constitutive law is dependent on the stress state, divide the strain increment subintervals and determine the updated stress increments through integration of the ~nstitutive equation, i.e.
Update the total stress and strain states: un= (T,+ A5 * r=r+As
Computational method for soil/structure interaction problems DISCUSSIOK
&,=&,+A&” ‘i=y+Ay If the interface tangent stiffness is to be updated, use the approximate symmetric elastic-plastic constitutive law:
NUMERICAL APPLICATION
This section presents the results of a numerical simulation to demonstrate theability of theanalytical models to predict buried structure response, Figure 6 shows the numerical configuration of a buried reinforced concrete model experiment under explosive air blast loading (conducted by SRI International). Figure 7 gives the details of the element representation of the test model. The results of pre-test calibration experiments were used to select appropriate concrete, reinforcement, soil and interface material parameters. Figure 8 shows the corresponding numerical results for concrete strain histories at various axial stations. The numerical results compare very well with those observed in the test. The failure in the text model corresponds to the calculated region of peak concrete strain. Additional sensitivity analyses show the effect of the soil/structure interface condition. As shown in Fig. 9, moderate changes in the coefficient of friction has a small effect on concrete response. However. a high value of cohesion gives significantly different results.
161
The
results presented in the last section give a clear indication of the ability of the CRT,/NONSAP code to simulate a complex nonlinear dynamical soil/structure response, involving elastic-plastic behavior of soil and structure constituents and the soil/structure interface behavior. This last effect was modeled using a structure/media interface element which incorporates incremental plasticity theory to predict interface sliding stress states. A nonassociated flow rule was selected to prevent volumetric expansion of the interface when sliding occurs. The nonsymmetric incremental system
I
P(T)
Fig. 7. Finite element description of the reinforced concrete structure.
FIN.6. Finite element description of surrounding media and interface of a buried reinforced concrete cylinder.
tangent stiffness which results from the selection of a nonassociated flow rule was neglected. and the errors associated with this approximation were eliminated through iteration and by reapplication of the errors to the next load step to prevent error accumulation. This approach.was very successful. In fact, in most applications. satisfactory results were obtained without any iteration. The basic idea of limiting the surface stress conditions along boundaries of dissimilar materials using a thin interface element whose material constituents reflect typical interface behavior can be extended to other applications. For example. the basic approach can be modified to account for initial gaps existing between material boundaries. That is. closing of a gap can be detected by monitoring the normal strain in a gap element. Subsequently, stress along the interface after closure can be limited using the identical procedure as employed by the interface element. This gap element will be implemented into the CRTjNONSAP finite element code for both two- and three-dimensional applications.
Y
162
MAFWR~ ITo et ui
Time,
mti Itseconds
Fig. 8. Calculated axial concrete strain histories.
Interface 0-c
condition
c = 0.00 MP, =OO?MP,p
~1 10.00 =0.18
3
Peak axial strain,
per cent
Fig. 9 REFERENCES
1. K. J. Bathe, E. L. Wilson and R. H. Iding, NONSAP, A structural analysis program for static and dynamic response of nonlinear systems. SESM Rep. 74-3, Dept. of Civil Engineering, University of California Berkeley, California (Feb. 19741 3. K. J. Bathe, H. Ozpemir and E. L. Wilson, Static and dynamic geometric and material nonlinear analysis. SESM Rep. 74-4, Dept. of Civil Engineering, University of California Berkeley, Calif. (Feb. 1974).
3. D. C. Drucker and W. Prager. SotI mechanics and plasttc analysis of limit design. Q. Appl. Mech. 10,157-163 (1952). 4. G. C. Nayak and 0. C. Zienkiewicz, Convement form of stress invariants for plasticity. Proc. AXE 98(ST4\. 949-954 (1972). 5. P Huck, T. Liber, R. L. Chiapetta, N. T. Tbomopouios. Jr and M. M. Singh, Dynamic response of soil/ooncrete interfaces at bigb pressures. AFWL-TR-73-264, Air Force Weapons Laboratory. Kirtland A.F.B.. New Mexico (Apr 1974)