- Email: [email protected]
Computational aspects of disturbed state constitutive models
Computer methods in ap@bd mea#mUcs and englneatlng ELSEYIER
Comput.
Methods
Appl. Mech. Engrg.
151 (1998) 361-376
Computational aspects of disturbed state constitutive models Chandra S. Desai*, Wu Zhang Departmenr
of Civil
Engineering In honor
and Engineering
of Professor
Mechanics,
The University
J. Tinsley Oden, on the occasion
Received
7 January
of Arizona,
Tucson, AZ 85721,
USA
of his 60th birthday
1997
Abstract Disturbed state concept (DSC) provides a unified basis for constitutive modelling including elastic, plastic and creep deformations, microcracking, damage and softening, stiffening, and cyclic fatigue under thermomechanical loading. It includes intrinsically regularization, localization, characteristic dimension and avoidance of spurious mesh dependence. It also leads to a new procedure for mesh adaptivity, particularly in zones experiencing significant microcracking and softening. The idea of the critical disturbance allows identification of the initiation of fracture and its growth as a consequence of microcrack coalescence. The DSC is implemented in nonlinear finite element procedures and a number of problems have been solved, and results are evaluated and compared with solutions by other methods so as to illustrate the foregoing capabilities of the model.
1. Introduction Many problems in engineering involve nonlinear material behavior. Here, elastic, plastic and creep strains, microcracking and softening, stiffening and fatigue failure under thetmomechanical loading are among the important factors that influence the reliability of solutions from computational procedures. For instance, when microcracking leading to growth of cracks and fractures, localization and softening occur, it becomes necessary to develop consistent procedures beyond the traditional ones based on continuum theories. The disturbed state concept (DSC) is a unified constitutive modelling approach [l-9] that provides a hierarchical framework to include elastic, plastic and creep strains, microcracking and softening, stiffening and fatigue under thermomechanical loading. 1.1.
Scope
This paper includes (1) a brief description of the DSC, (2) physical aspects such as microcracking, localization and characteristic dimension, (3) spurious mesh dependence, (4) simplified procedures based on critical disturbance criterion to identify onset of fatigue failure and growth of fracture, and (5) adaptive mesh refinement by incorporation of the critical disturbance criterion. A number of problems are solved to illustrate the capabilities of the model. 2. The disturbed
state concept
Details of the DSC are given elsewhere [l-9]; a brief description is presented below. It is assumed that applied mechanical and environmental forces cause disturbances or changes in the material’s microstructure with respect to its behavior under its relative intact (RI) and fully adjusted (FA) states * Corresponding
author.
0045-7825/98/$19.00 0 1998 Elsevier Science S.A. All rights reserved PII SOO45-7825(97)00159-X
362
C.S. Desui, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
in the material. The RI state refers to the behavior of the material part that excludes the influences of factors such as friction, anisotropy, microcracking, damage, softening, and stiffening or growth. The material in the RI state is often modelled by using a continuum theory. The FA state refers to the asymptotic (‘equilibrium’) state to which the material tends to in the ultimate stages of deformation. The material in the RI state is transformed continuously into the FA state through a process of natural self-adjustment of its microstructure. At any stage during deformation, the material is treated as the mixture of the RI and FA parts, which are distributed (randomly) over the material elements, Fig. l(a). The RI and FA states of the material are called the reference states. The material is initially in a full or a partial RI state, depending upon the initial disturbance, due to factors such as initial stresses, anisotropy and manufacturing. During the deformation, the volume of the material in the FA state increases, and the volume of the RI state material decreases. This process involves continuous interaction between the material parts in the RI and FA states. The observed response of the mixture is expressed in terms of the responses of the materials in the RI and FA states by using the disturbance function, which can be considered to denote deviation of the observed behavior from those of the two reference states; a _ are the second invariants of the deviatoric stress symbolic representation is shown in Fig. 1(b). Here, J?,, and I,,] and strain tensors, respectively. The behavior of the material in the RI state can be simulated by using a theory based on continuum mechanics such as linear or nonlinear elasticity, elastoplasticity, viscoplasticity or thermoviscoplasticity. The FA state can be represented as ( 1) the ‘failed’ material which acts like a ‘void’ as in the classical damage theory and
inltlal
Intrrmrdlrte
/
” Fallurr”
0 Schomatlc Local
(I
Matlvr
Intact
Zones
Fully Adju~trd (Crltlcrl) Zones
InstabilIty
intact O
A
0 D
Ri
observed
-I-
Q R=
Q 4” c -
T
fully (critical
Q RS
dlsturbed sate)
0
Fig.
I. Representation
of disturbed
state concept
(R ti Response).
C.S. Desai, U! Zhang I Comput. Methoa!s Appl. Mech. Engrg. 151 (1998) 361-376
363
can carry no stress at all [lo], (2) it can carry hydrostatic stress but no shear stress, like a constrained liquid, or (3) it can continue to carry the shear stress for a given hydrostatic stress, reached up to that point and can deform under constant volume, as in the case of the critical state soil mechanics concept; here, the material acts like a constrained liquid-solid. As the FA material is constrained by the surrounding RI material, the two latter simulations are considered to be more realistic than the first one. At this time, the disturbance, D, is assumed to be a scalar; however, it can be expressed as a tensor [5]. The DSC includes the coupled (observed) response, as it is influenced by the collective behavior of the interacting mechanisms in the RI and FA states. Hence, the material response is represented in an integrated manner based upon the responses of the materials in the RI and FA states. Thus, it is not necessary to measure and define particle level constitutive response, as in the case of micromechanics models. Also, since the microcrack interaction is included implicitly in the model, it is not necessary to superimpose effects of forces and kinematics in individual or collected microcracks with the damage theory [ 111.
3. Constitutive
equations
Based on the equilibrium of forces in the observed stress tensor, da:, is derived as du;
= (1 - D) doT: + D da;
(mixture),
RI and FA states, the incremental
+ dD(a;; - al)
observed
(la)
in which
(lb)
a;=(l-D)afj+Da;,
and where a, i and c denote observed, RI and FA states, respectively, D is the scalar disturbance function (described later), and d denotes increment. In Eq. (la), a; = a; - a: denotes the relative stress as the difference between the stresses in the FA and RI states; Fig. 2 shows a symbolic representation in which the randomly distributed clusters in the RI and FA states are collected in a weighted sense. Eq. (la) can be written as da;
= (1 - D)C:,,
de;,
+ DC;,,
de;, + dDu;
where Cijke is the fourth order tensor related to the constitutive behavior. In the DSC, the relative stress (u;) can cause relative motions (translations, rotations). Together with the second term on the right side in Eq. (2), they incorporate the effect of microcrack interaction. As a consequence of the relative stress, a;, the strains, deij, are different in the FA and RI states. Thus, I$ and r$ are the tensors of strains in the RI and FA parts, respectively. It is usually difficult to establish a relationship between these two strains; in the finite element analysis, it is evaluated during the incremental analysis through an iterative
Fig. 2. Stresses in RI and FA
parts and relative motion.
C.S. Desai, W. Zhang
364
procedure expressed
starting from the assumption as
I Comput. Methods Appi. Mech. Engrg. 1.51 (1998) 361-376
that they are initially
equal. A simple relation between the two can be
deIi = (1 + a) de;,
(3)
where (Y is the relative motion parameter, which can be the function trajectory and the disturbance. Also, dD can be derived as
of deformation
dD = R,, de:,
history, e.g. plastic strain
(4)
Details of R;, are given elsewhere Now, Eq. (2) can be written as
[7].
do:; = [( 1 - D)C;,A, + D( 1 + ‘y)C;,,, + a;, . R,,] d&
(W
dcr;; = (Lilk( + L;;,,) de;,
(5b)
du:; = C,,,, . de;,
(5c)
or
or
where
Lijki = (1 3. I. Disturbance
D)CI,,, + DC1+ a)C:,k/
and
L:+( = u:,R,,
.function
The disturbance which is influenced by the physical properties of the material and loading is considered to denote the deviation of the observed behavior with respect to those of the material’s reference states. The disturbance function, D, can be defined based on observed stresses, volumetric response, or nondestructive (ultrasonic) velocities during loading, unloading and reloading; details are given elsewhere [5]. Here, D, based on stresses, is defined as (Fig. 3)
D=
5’ - ii”
(6)
5’ _ (TL
Fig. 3. Schematic
of stress-strain
response
and disturbance.
C.S. Desai, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
365
1
0.8 0.8 L3 0.4 0.2 0 Fig. 4. Schematic
0
0.2
of disturbance
function
0.4
0.8
as dependent
0.8
1
on parameters
1.2
1.4
such as mean pressure
and temperature.
where 5 is the equivalent stress such as dzo, the second invariant of the deviatoric stress tensor, S,, or T, the octahedral shear stress. At this time, it is assumed that the RI response is stiffer than the observed response; hence, D is positive. If it is softer than the observed response, D will be negative, indicating stiffening or strengthening in the observed response; such a behavior can be included with additional considerations [9]. Now, D is expressed as
D = DU[l - exp(-A&E)]
(7)
where A and Z are material parameters and D, is the value of D corresponding to the residual condition, Fig. 3. A schematic of D vs.5, is shown in Fig. 4; here, &, is the trajectory of deviatoric plastic strains, Ez,as
(8) Here, 5, is based on irreversible or plastic strains corresponding to the observed response. However, as a simplification, it can be approximated in terms of plastic strains corresponding to the RI response based on the hierarchical (a,,) plasticity model [ 121. Comparisons with and differences between the DSC, and the classical damage models, damage models with microcrack interaction, micromechanics models and self-organized criticality (SOC) concept are discussed in ]4,51.
4. Finite element equations Based on the virtual work principle,
the finite element equations
corresponding
to Eq. (5) are derived as [7]
or
j- ([~Tlt~l[N) dV{dq’I= {Q>- I, WITb”I dv V
(9b)
where [B] is the strain-displacement transformation matrix in {de’} = [B] {dq’}, {q} is the vector of nodal displacement, and {Q} is the external load vector. Eq. (9a) or (9b) involves two unknowns, {dq’}, {w’}. The details of solution schemes are given in [7,8]; they are briefly discussed below.
C.S. Desai, W. Zhang / Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
366
4.1.
Solution schemes
Eqs. (9) can be solved using the following schemes: ( 1) An incremental iterative procedure under strain controlled condition is used. Here, the tangent matrix [C] is not positive definite, particularly in the strain softening zone. However, with the strain controlled loading, the procedure provides convergent and accurate solutions, which is illustrated through solutions of a number of problems, which are discussed later; (2) the term related to the relative stress, a:, is taken to the right-hand side, to lead to
Here, the term related to [L’] is evaluated as the load vector at the previous increment n - 1. Then, the system matrix on the left-hand side is positive definite because 0 G D 5 D,,, and a 2 0. Such a transformation does not change the positive or negative definite character of the problem; however, it can provide an approximate and simpler procedure for the finite element analysis. The algorithm in Eq. (10) is considered similar to that in the nonlocal damage models [13]; and (3) the term {v’} can be treated as an independent (constraint) unknown, together with the displacements, and the problem is treated as coupled [5,8].
5. Physical characteristics Internal microstructural changes due to factors such as relative motions between material particles, can lead to but involves discontinuous and microcracking. As a result, the material may no longer be a continuum, nonhomogeneous conditions. For consistent solutions, it is necessary to account for the discontinuous nature of the material. Incorporation in the model of the characteristic dimension, and procedures to define strain in a weighted sense over specific zones in the neighborhood of a joint, are among the procedures to account for the discontinuous material. Nonlocal formulations [ 131, inclusion of microcrack interaction [ 111 and Cosserat and gradient enrichments of the damage nodes are among the procedures [14-161. It has been shown that with the (1) inclusion of the second term, and the term dLM:, leading to the relative motion due to the third term in Eqs. (1) (2) the expression of the disturbance (A, 2 and D,,), Eq. (7), in terms of factors such as the initial density, confining pressure and characteristic length ratio (say, LIB where L is the length and d the diameter of material specimen), and (3) a simple averaging procedure for the (plastic) strains, the DSC model allows for the nonlocal effects, regularization, localization, and avoidance of spurious mesh dependence [7,8,]. Typical examples are given subsequently.
6. Microcracking,
fatigue with critical disturbance
Theories from linear fracture mechanics and certain modifications to account for nonlinear effects are often used to identify initiation and propagation of fracturing in deforming materials. However, it is felt that the microcracking due to the microstructural changes is a process in the material structure involving coalescence of microcracks that can lead to the formation and growth of distinct or finite-sized cracks, Indeed, such a state is the consequence of the microcrack coalescence into ‘weak’ zones in the material mixture due to the accumulation of plastic or irreversible deformations. The DSC allows consideration of microcracking as a process due to interacting responses of the material parts in the RI and FA states. Thus, it is possible to define critical levels of disturbance, D,, Fig. 4, at which the microcracks culminate into states that can be treated as causing initiation and final fracture or ‘failure’ involving severe intensity of cracks, at which the engineering service function of the material ceases. Once the critical state of microcracking is identified, which can lead to ‘crack’ initiation, the cracking can grow during subsequent loading with further increase in the disturbance. When the ultimate disturbance, D,, is reached, the material can be considered to have ‘cracked fully’.
C.S. Desai, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
361
Such critical states of disturbance can also define instability in the material’s microstructure. In the case of cyclic thermal loading, such a state indicates fatigue failure [6,7,17]. Another example is liquefaction in saturated sands under cyclic (dynamic) loading [ 181. It has been shown that when the critical disturbance, D,,is reached under cyclic loading, the corresponding accumulated energy or energy/cycle experiences an ‘abrupt’ change toward stabilization or saturation. Fig. 5 shows plots of disturbance vs number of thermal cycles from the analysis of an electronic chip-substrate system involving ceramic chip, Pb-Sn solder as the joining material and plastic printed wire board [7,17]. It can be seen that the value of D,, at which its rate of change approaches stabilization or saturation, corresponds with the similar state of energy/cycle around thermal cycles of about 400; which compared well with the observed cycles to failure in laboratory tests 1191. Fig. 6 shows disturbance and accumulated energy vs number of cycles for a saturated sand tested under cyclic shear loading [20]. It was found in the laboratory that initial (instability) or liquefaction occurred at about 11 cycles, and the final (ultimate) liquefaction at about 15 cycles. The disturbance and energy plots show similar values. For example, the critical disturbance, D,,corresponds to about 11 cycles at initial liquefaction, and the ultimate disturbance, D,, to about 15 cycles at the final liquefaction. Thus, the critical disturbance criteria can provide a simplified procedure for the identification of crack initiation, and growth and fatigue failure or instability. As the disturbance is computed as the part of the numerical computation, it is relatively easy to use this procedure compared to other procedures such as energy and fracture criteria.
(a)
Fig. 5. Disturbance
D-
vs tempaum cycles
and energy density per temperature
cycle in solder joint [7].
368
w. Zhang
C.S. Desai,
/ Comput. Methods Appl. Mech. Engrg.
151 (1998) 361- 376
0.8
cc %
0.8
0.4 0.2
CrcJqN
(a) D vs N
(b) Energy vs N [20] Fig. 6. Disturbance and energy vs cycles for cyclic sand behavior.
7. Adaptive
mesh strategy
It is possible to employ the adaptive mesh refinement procedures to identify the initiation and growth of localization and cracking. It can be used in conjunction with the critical disturbance criterion. Many previous adaptivity schemes have been developed for linear problems. Hence, it is necessary and required to develop them for nonlinear problems. Here, the disturbance approach is used in a new scheme to develop mesh adaptivity procedures for nonlinear problems, particularly in zones where microcracking and softening prevail. Fig. 3 shows a schematic of the stress-strain response that is affected by microcracking leading to strain softening. The mesh adaptivity procedure proposed involves two parts. For the prepeak behavior, where microcracking and disturbance are relatively small, a modified scheme based on an available procedure based on uniform error distribution [21,22] is used. However, for the postpeak behavior, where disturbance grows rapidly,
C.S. Desai, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
a new procedure given below. 7.1. Prepeak
based on critical disturbance
augmented
by the uniform
369
error scheme is used; brief details are
response
The uniform degree of error distribution (UDED) scheme [21] is used. In this scheme, a hybrid finite element procedure is used to evaluate the exact stress {g*} for the evaluation of errors, e(i) e(i) = abs( l”*2ns:“‘)
(11)
where i denotes an element, {a} is the vector of computed is the vector ‘exact’ stresses, described below. 7.2. Hybrid FE procedure The hybrid potential lr
Pm=
;
” {dTKl{d
formulation,
and {g*}
for exact stresses [23,24]
energy is expressed
I
stresses from the displacement
dV +
as
$ Iv C(E)-[Cl -'bHTKI(k~
-
[Cl
-‘b>) dV + I,h2)‘WdV
where {u} and {E} are the displacement and strain vectors in the displacement formulation, load vector, [C] is the constitutive matrix and {a} is the assumed stress given by
(12)
{Q] is the external
{a] = uw3
(13)
Here, [P] is the matrix of interpolation function and {/3} is the vector of stress parameters. four-node isoparametric element, assumed [P] is given by
For the case of
(144 where [F,] = [l r s] and r and s are local coordinates. The resultant [F,] = [I r* s*] where r* = r - a and s* = s - b, and a and b are determined based on the constant stress patch test such that
It is shown [23] that introduction of the second term with appropriate selection of [P] provides ‘softening’ of the traditional displacement formulation and leads to improved calculation of the stresses {a*}. Also, the second term represents a least-energy fit between the strains {E} and assumed stresses {u*}. The exact stresses {a*} are derived as [25] {a*] = U’IW-‘[Glbd where [H] - ’ = J, [P]T[C]-‘[P]
(15) dV; [G] = J, [PITIB] dV and {q} is the computed
displacement
vector.
7.3. Strategy for mesh refinement The objective
is to generate
a nearly optimal
mesh such that
ese
(16)
where e is the computed error, Eq. (ll), and P is the prescribed permissible approximation error for each element in an m-element mesh is given by
error, e.g. 5%. Now, the
e, = Zllhl(m”* where I(h(l= h,,,
is the norm or maximum
(17)
size of element.
The local mesh enrichment
indicator
is given by
370
C.S. Desai, W Zhang
/ Comput. Methods Appl. Mech. Engrg. 1.51 (1998)
361-376
(18) where l/E]/, is the energy norm for element
llEl(,=
Iv{e(~)~TICl _'[WI1
Then the element
i, given by
dV
(19)
size for the optimal mesh is calculated
approximately
as
(20) where hi is the original maximum size of elements, and p denotes the order of interpolation elements near singularities of order k, the new element size is estimated from
function.
For
(21) where k represents
the strength of the singularity,
which is often adopted to be equal to 0.5.
7.4. Postpeak response In this zone, the microcracking and resulting softening can be significant, and the foregoing UDED scheme alone is found to be not sufficient. Hence, the disturbance function, D, in the DSC, Eq. (7), is used as the ‘error’ indicator for remeshing. Here, a critical value of disturbance, D * = 0.5, 0.75 or 0.90 can be specified. As the disturbance is dependent on the plastic strains, Eqs. (7) (8), it represents the intensity of strain localization. Thus, the greater the value of D, the greater is the softening and localization. As a result, use of D* as the indicator for remeshing is consistent with the physical state of the structure as microcracking, and its growth follows the growth of disturbance, Fig. 4. In the postpeak region, the error indicator is expressed as DaD*
(22)
Here, the error ratio for stress is given by
IleA
(23a)
l=ile,ll where e,, is the lowest error and e, is the error in any element. following strategy:
h, = i’k =
The mesh is refined by using the
Ile,,ll’
( ) ile,ll
(23b)
k#
h,, is the size of the element results.
with the lowest error and t is an exponent.
In general, t = 2 provides
satisfactory
8. Validations and applications In order to verify the DSC model with FE procedure, a number of problems have been solved. They include ( 1) simple (metal) tension bar with imperfection for localization, (2) plate (metal) with imperfection for localization, shear banding and mesh dependence, (3) mesh dependence in a concrete footing, (4) mesh adaptivity in a footing problem, (5) localization in electronic chip-substrate systems involving Pb-Sn solders including cyclic fatigue, and (6) dynamic soil-structure interaction including liquefaction. Here, the RI state was characterized using a plasticity model (6,) in the hierarchical single surface (HISS) models [ 121, and the FA
C.S. Desai, W Bang
I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
state was defined using the constrained below.
liquid or constrained
liquid-solid
371
[3-g]. Only typical problems are given
8.1. Tension bar with impe$ection The tension-bar problem with imperfection considered by de Borst et al. [16] using the gradient enrichment with a classical plasticity model, was solved by using the DSC. Here, the RI response was simulated by using the 6,-plasticity model [4,12] and the disturbance parameters A, 2, D, were expressed as function of the ratio t/L, where t is the characteristic length and L is the length of the bar. It was shown that the DSC yields as good as or better predictions of localization and its zones, compared to those by the gradient enrichment procedure PI. 8.2. Mesh adaptivity:
Footing problem
Zienkiewicz et al. [22] used an adaptive remeshing procedure for the analysis of localization in which the material was modelled by using a classical theory of plasticity, with softening simulated by using negative hardening modulus. The authors stated that the procedures based on the acoustic tensor for the analysis of existence and localization were not reliable. They proposed that the discontinuity and localization can be captured by appropriate formulation of the finite element method, and adaptive remeshing. Here, we consider a typical problem of footing on a strain-softening material, and show that the DSC with mesh adaptivity can provide a consistent procedure to analyze localization. Fig. 7(a) shows a test problem in which uniform displacement increment = 0.001 m applied at the top, while the right-hand side boundary is free to move. The RI material is assumed to be elastoplastic with continuous
StlWS
t
(a) Footing
,/
_,___,__.._._....._.. v .. .. ..__.. .,__.. -Intact
(h) Schematic of stress-strain curve
(4 (c) Adapted meshes at typical locations (A, B, C) Fig. 7. Mesh adaptivity
for footing problem.
C.S. Desai, W. Zhang
312
yielding [4,12], and the FA material used are Elasticity:
I Comput. Methods Appl. Mech. Engrg. I51 (1998) 361-376
is assumed
to carry only the hydrostatic stress. The material parameters
E= 2x IOXMPa 1’ = 0.49
Plasticity: y = 0.30 x lo-J, p = 0 Ultimate: Phase change: II = 2. I Yielding: = lo- ‘(I. 7, = 0.2 Disturbance: :‘= IO-‘, 2 = 2.0, n = 0.8 Critical disturbance for remeshing: D*“= 0.75 The initial mesh (point A), and the adopted meshes using the procedure proposed herein, at the peak (point B) and during softening (point C), are shown in Figs. 7(b) and (c). It can be seen that the mesh adapts in the expected direction of slip and localization. Fig. 8 shows the growth of disturbance with applied displacement at two points (a, b) as the remeshing progresses. Point (a) lies within the localization zone and experiences higher disturbance beyond D* = 0.75, at peak and beyond peak. This indicates that the mesh adapts in the zones that experience critical disturbance and beyond. The disturbance at point (b), which is away from the localization zone, does not experience higher disturbance. Fig. 9 shows nondimensionalized load-displacement curves obtained by using the final mesh Fig. 7(c), Here, the analyses were performed using (a) the plasticity model and (b) the DSC with small level of softening (D, = 0.15). Various terms in Fig. 9 are: F = computed load on the top, B = width of the loaded area, c,, = cohesive strength, U, = vertical displacement at the top, and E = elastic modulus. It can be seen that the DSC model provides softer and consistent response compared to that by the plasticity model. 8.3. Localization
und mesh drpendmt~c:
Plutr MYth imperfection
A material specimen under biaxial test with width of 60 mm and height of 180 mm, solved by de Borst et al. [ f6] using gradient enrichment with the Drucker-Prager plasticity model, is shown in Fig, IO. An imperfect element with 10% reduction in the strength is provided, as shown in the figure. The plate was subjected to incremental displacement = 0.0001 m at the top. de Borst et al. [ 161 used four-node elements with bilinear displacement field, and a Hermitian (bicubic) interpolation for the plastic multiplier.
I
D -Disturbance
1
function
-
At location B
0.75.- D* At location b
A-MeshA
(a) Typical points
B-Mesh3
C-McshC
(b) D vs displacement Fig. 8. Growth of’ disturbance at typical points.
C.S. Desai, W. Bang
I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
373
U,EIBC. Fig. 9. Comparison
of predictions
The DSC model was used with RI response material parameters used are given below: Elasticity: Plasticity: Disturbance:
from DSC (0, = 0.15) and plasticity
characterized
as elastoplastic
models.
with continuous
yielding.
The
E = 11920 MPa Y = 0.49 y = 0.001, p = 0.0, Iz = 2.10 a, = lo-lo, v, = 0.20 A = 530, Z = 1.35, D, = 0.60
In the present analysis with the DSC model, two mesh layouts were used-6 X 18 mesh and 13 X 36 mesh with 8-noded quadrilateral elements. Figs. 11(a) and (b) show deformation patterns at the peak and after the peak for the 6 X 18 mesh, while Figs. 12(a) and (b) show similar results for the 13 X 36 mesh. Although the deformation, localization and shear banding patterns in the present analysis are somewhat different from those in de Borst et al. [16], the trends are essentially similar. The differences can be due to factors such as different material model and element approximation.
Incremental displacement
3 imperfection
Fig. 10. Plate with imperfection
VW.
C.S. Desai, W. Zhang
Iii
I Comput. Methods Appl. Mech. Engrg. I51 (1998) 361-376
i i i 1
(a) At PeaL
(a) At peak
@) After+ Fig.
I I.
Deformation
Fig. 12. Deformation
patterns for 6 X 18 mesh. patterns for I3 X 36 mesh.
Fig. 13 shows the load displacement curves from the two meshes, and indicates no significant mesh dependence. Here, F = load induced at the top, B = width of the plate, 6 = observed or average stress, V,op = displacement at the top and H = height of the plate. The above results show that the DSC models provide consistent and satisfactory localization and shear banding, and also avoids the spurious mesh dependence.
8.4. Mesh dependence:
Concrete
,footing
Fig. 14(a) shows the concrete block subjected to prescribed displacements. The footing is unconfined or free at the top and on the right side. The FE analysis with the DSC involved meshes with 4-, 16-, 64- and 2568-noded isoparametric elements. The strains were averaged over the patch of elements in the part of the block, ABCD. It was found that averaging over a single element can also provide satisfactory results. Fig. 14(b) shows typical results in terms of axial stress, a,., vs. strain, E,,, with the number of eIements. It can be seen that the results are essentially independent of the mesh.
Fig. 13. Load-displacement
@
6 x 18 8-nodeclernents
a
13 x 36 8-node elements
curves for 6 X 18 and 13 X 36 meshes.
C.S. Desai, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
+
4Ebm.
1*
(I) Mesh
forcmcme
16Ebm.
-6464.
315
block
-=-256Ebm.
@) 4 “S 5
Fig. 14. Mesh dependence
for concrete
block.
9. Conclusions The disturbed state concept provides a unified approach for constitutive modelling of nonlinear materials and interfaces. It allows a hierarchical framework from which simple to sophisticated models can be adopted, depending upon the specific material and need. The DSC intrinsically provides for regularization, characteristic dimension, localization and avoidance of spurious mesh dependence. One of the important attitudes of the DSC is that it can be used to develop ‘error’ indicators for mesh adaptivity for nonlinear behavior involving microcracking and softening. Overall, the DSC is found to be a versatile modelling approach which, when implemented in numerical (finite) element procedures, can provide consistent and reliable computations.
Acknowledgments Parts of the results here were supported by Grant Nos. DDM-9102177 and DMI-9313204 from the National Science Foundation, Washington, D.C. Drs. J. Chia, C. Basaran, and T. Dishongh participated in the results presented in this paper.
References [l] C.S. Desai and Y. Ma, Modelling 16(9) (1992) 623-653.
of joints and interfaces
using the disturbed
state concept, Int. J. Numer. Analyt. Methods Geomech.
316
C.S. Desai, W. Zhang
I Comput. Methods Appl. Mech. Engrg. ISI
(1998)
361-376
[2] S. Armaleh and C.S. Desai. Modelling and testing of a cohesionless material using the disturbed state concept, Int. J. Mech. Behavior Mater. S(3) (1994) 2799295. [3] D.R. Katti and C.S. Desai. Modelhng and testing of cohesive soil using the disturbed state concept, J. Engrg. Mech., AXE 121(S) (1995) 648-658. [4] C.S. Desai. Constitutive modelling using the disturbed state concept, m: H.B. Miihlhaus, ed., Continuum Models for Materials with Microstructure, Chap. 8 (John Wiley & Sons. UK, 1995). [S] C.S. Desai and J. To& Disturbed state constitutive modelling based on stress-strain and nondestructive behavior, Int. J. Solids Struct., 33(1 I) (1996) 1619-1650. [6] C.S. Desai, J. Chia, T. Kundu and J.L. Prince, Thermomechamcal response of materials and interfaces in electronic packaging: parts I and II. J. Elect. Packaging. ASME ( 1996) in press. [7] C. Basaran, C.S. Desai and T. Kundu, Thermomechanical tinite element analysis usmg the disturbed state concept, J. Elect. Packaging. ASME (1995) submitted. ]8] C.S. Desai. C. Basaran and W. Zhang. Numerical algorithms and mesh dependence in the disturbed state concept, Int. J. Numer. Methods Engrg. ( 1996) accepted. 191 T. Dishongh and C.S. Desai. Disturbed state concept for materials and interfaces with applicattons in electronic packaging, Report to NSF, Dept. of Civil Engrg. and Engrg. Mechs., The University of Arizona, Tucson, Arizona 1996. IlO] L.M. Kachanov. Introduction to Continuum Damage Mechanics (Martinus Nijhoft Publishers. The Netherlands, 1986). [ 1I ] Z.P. Bazant. Nonlocal damage theory based on micromechanics of crack interactions, J. En&r&. Mech., ASCE 120(3) (1994) 593-617. (121 C.S. Desai. S. Somasundaram and G. Frantziskonis, A hierarchical approach for constitutive modelling of geologic materials, Int. J. Numer. Analyt. Methods Geomech. 10 (1986) 225-257. [ 131 Z.P. Bazant and F.B. Lin, Non-local yield limit degradation, Int. J. Numer. Methods Engrg. 26 ( 1988) 1805-1823. [ 14) R. de Borst and H.B. Miihlhaus, Gradient dependent plasticity: formulation and algorithmic aspects, Int. J. Numer. Met. Engrg. 35(3) (1992) 521-539. 1IS] R.H.J. Peerlings, R. de Borst, N.A.M. Brekelmans and J.H.P. de Vree, Gradient enhanced damage for quasi-brittle materials. Int. J. Numer. Methods Engrg. 39 (1996) 3391-3403. [ 161 R. de Borst, L.J. Sluys, H.B. Miihlhaus and J. Pamin, Fundamental issues in finite element analyses of localization of deformation. Engineering Comput. IO (1993) 99-121. [ 171 C.S. Desai, C. Basaran. T. Dishongh and J. Prince, Thermomechanical analysis in electronic packaging with unified constitutive model for materials and joints. Components. Packaging and Manufacturing Technology, Part B: Advanced Packaging, IEEE (1996) submitted. [IS] C.S. Desai. C. Shao and D. Rigby. Discussion to evaluation of liquefaction by energy principles, J. Geotech. Engrg., ASCE 122(3) ( 1996) 241-244. [ 191 P.M. Hall and W.M. Sherry, Materials, structures and mechanics of solder-Joints for surface-mount microelectronics technology, Proc., 3rd Int. Conf. on Techniques de Connecion en Electronique, Welding Society, Fellbach, Germany, February 1986. [20] L. Figueroa, A.S. Saada. L. Liang and N.M. Dahisaria, Evaluation of liquefaction by energy principles. J. Geotech. Engrg.. ASCE l20(9) (1994). 1211 O.C. Zienkiewicz and J.Z. Zhu. A simple error estimation and adaptive procedure for practical engineermg analysis, Int. J. Numer. Methods Engrg. 24 (1987) 337-357. (221 O.C. Zienkiewicz, M. Huang and M. Pastor, Localization problems in plasticity using finite elements with adaptive remeshing, Int. J. Numer. Analyt. Methods Geomech. I9 (1995) 127-148. (231 W. Zhang, Theoretical basis and general optimal formulations of isoparametric generalized hybrid/mixed element model for improved stress analysis, Acta Mech. Sinica 9 (1993) 277-288 (English edition). [24] W. Zhang and D.P. Chen, The patch test conditions and some multivariant finite element formulations, Int. J. Numer. Methods Engrg. ( 1996) accepted. [25] W. Zhang and C.S. Desai, Adaptive mesh refinement with the disturbed state model and application NSF, Dept. of Civil En&r&. and Engrg. Mechs., The Univ. of Arizona, Tucson, Arizona, (19%‘).
in electronic
packaging,
Report to
Comput.
Methods
Appl. Mech. Engrg.
151 (1998) 361-376
Computational aspects of disturbed state constitutive models Chandra S. Desai*, Wu Zhang Departmenr
of Civil
Engineering In honor
and Engineering
of Professor
Mechanics,
The University
J. Tinsley Oden, on the occasion
Received
7 January
of Arizona,
Tucson, AZ 85721,
USA
of his 60th birthday
1997
Abstract Disturbed state concept (DSC) provides a unified basis for constitutive modelling including elastic, plastic and creep deformations, microcracking, damage and softening, stiffening, and cyclic fatigue under thermomechanical loading. It includes intrinsically regularization, localization, characteristic dimension and avoidance of spurious mesh dependence. It also leads to a new procedure for mesh adaptivity, particularly in zones experiencing significant microcracking and softening. The idea of the critical disturbance allows identification of the initiation of fracture and its growth as a consequence of microcrack coalescence. The DSC is implemented in nonlinear finite element procedures and a number of problems have been solved, and results are evaluated and compared with solutions by other methods so as to illustrate the foregoing capabilities of the model.
1. Introduction Many problems in engineering involve nonlinear material behavior. Here, elastic, plastic and creep strains, microcracking and softening, stiffening and fatigue failure under thetmomechanical loading are among the important factors that influence the reliability of solutions from computational procedures. For instance, when microcracking leading to growth of cracks and fractures, localization and softening occur, it becomes necessary to develop consistent procedures beyond the traditional ones based on continuum theories. The disturbed state concept (DSC) is a unified constitutive modelling approach [l-9] that provides a hierarchical framework to include elastic, plastic and creep strains, microcracking and softening, stiffening and fatigue under thermomechanical loading. 1.1.
Scope
This paper includes (1) a brief description of the DSC, (2) physical aspects such as microcracking, localization and characteristic dimension, (3) spurious mesh dependence, (4) simplified procedures based on critical disturbance criterion to identify onset of fatigue failure and growth of fracture, and (5) adaptive mesh refinement by incorporation of the critical disturbance criterion. A number of problems are solved to illustrate the capabilities of the model. 2. The disturbed
state concept
Details of the DSC are given elsewhere [l-9]; a brief description is presented below. It is assumed that applied mechanical and environmental forces cause disturbances or changes in the material’s microstructure with respect to its behavior under its relative intact (RI) and fully adjusted (FA) states * Corresponding
author.
0045-7825/98/$19.00 0 1998 Elsevier Science S.A. All rights reserved PII SOO45-7825(97)00159-X
362
C.S. Desui, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
in the material. The RI state refers to the behavior of the material part that excludes the influences of factors such as friction, anisotropy, microcracking, damage, softening, and stiffening or growth. The material in the RI state is often modelled by using a continuum theory. The FA state refers to the asymptotic (‘equilibrium’) state to which the material tends to in the ultimate stages of deformation. The material in the RI state is transformed continuously into the FA state through a process of natural self-adjustment of its microstructure. At any stage during deformation, the material is treated as the mixture of the RI and FA parts, which are distributed (randomly) over the material elements, Fig. l(a). The RI and FA states of the material are called the reference states. The material is initially in a full or a partial RI state, depending upon the initial disturbance, due to factors such as initial stresses, anisotropy and manufacturing. During the deformation, the volume of the material in the FA state increases, and the volume of the RI state material decreases. This process involves continuous interaction between the material parts in the RI and FA states. The observed response of the mixture is expressed in terms of the responses of the materials in the RI and FA states by using the disturbance function, which can be considered to denote deviation of the observed behavior from those of the two reference states; a _ are the second invariants of the deviatoric stress symbolic representation is shown in Fig. 1(b). Here, J?,, and I,,] and strain tensors, respectively. The behavior of the material in the RI state can be simulated by using a theory based on continuum mechanics such as linear or nonlinear elasticity, elastoplasticity, viscoplasticity or thermoviscoplasticity. The FA state can be represented as ( 1) the ‘failed’ material which acts like a ‘void’ as in the classical damage theory and
inltlal
Intrrmrdlrte
/
” Fallurr”
0 Schomatlc Local
(I
Matlvr
Intact
Zones
Fully Adju~trd (Crltlcrl) Zones
InstabilIty
intact O
A
0 D
Ri
observed
-I-
Q R=
Q 4” c -
T
fully (critical
Q RS
dlsturbed sate)
0
Fig.
I. Representation
of disturbed
state concept
(R ti Response).
C.S. Desai, U! Zhang I Comput. Methoa!s Appl. Mech. Engrg. 151 (1998) 361-376
363
can carry no stress at all [lo], (2) it can carry hydrostatic stress but no shear stress, like a constrained liquid, or (3) it can continue to carry the shear stress for a given hydrostatic stress, reached up to that point and can deform under constant volume, as in the case of the critical state soil mechanics concept; here, the material acts like a constrained liquid-solid. As the FA material is constrained by the surrounding RI material, the two latter simulations are considered to be more realistic than the first one. At this time, the disturbance, D, is assumed to be a scalar; however, it can be expressed as a tensor [5]. The DSC includes the coupled (observed) response, as it is influenced by the collective behavior of the interacting mechanisms in the RI and FA states. Hence, the material response is represented in an integrated manner based upon the responses of the materials in the RI and FA states. Thus, it is not necessary to measure and define particle level constitutive response, as in the case of micromechanics models. Also, since the microcrack interaction is included implicitly in the model, it is not necessary to superimpose effects of forces and kinematics in individual or collected microcracks with the damage theory [ 111.
3. Constitutive
equations
Based on the equilibrium of forces in the observed stress tensor, da:, is derived as du;
= (1 - D) doT: + D da;
(mixture),
RI and FA states, the incremental
+ dD(a;; - al)
observed
(la)
in which
(lb)
a;=(l-D)afj+Da;,
and where a, i and c denote observed, RI and FA states, respectively, D is the scalar disturbance function (described later), and d denotes increment. In Eq. (la), a; = a; - a: denotes the relative stress as the difference between the stresses in the FA and RI states; Fig. 2 shows a symbolic representation in which the randomly distributed clusters in the RI and FA states are collected in a weighted sense. Eq. (la) can be written as da;
= (1 - D)C:,,
de;,
+ DC;,,
de;, + dDu;
where Cijke is the fourth order tensor related to the constitutive behavior. In the DSC, the relative stress (u;) can cause relative motions (translations, rotations). Together with the second term on the right side in Eq. (2), they incorporate the effect of microcrack interaction. As a consequence of the relative stress, a;, the strains, deij, are different in the FA and RI states. Thus, I$ and r$ are the tensors of strains in the RI and FA parts, respectively. It is usually difficult to establish a relationship between these two strains; in the finite element analysis, it is evaluated during the incremental analysis through an iterative
Fig. 2. Stresses in RI and FA
parts and relative motion.
C.S. Desai, W. Zhang
364
procedure expressed
starting from the assumption as
I Comput. Methods Appi. Mech. Engrg. 1.51 (1998) 361-376
that they are initially
equal. A simple relation between the two can be
deIi = (1 + a) de;,
(3)
where (Y is the relative motion parameter, which can be the function trajectory and the disturbance. Also, dD can be derived as
of deformation
dD = R,, de:,
history, e.g. plastic strain
(4)
Details of R;, are given elsewhere Now, Eq. (2) can be written as
[7].
do:; = [( 1 - D)C;,A, + D( 1 + ‘y)C;,,, + a;, . R,,] d&
(W
dcr;; = (Lilk( + L;;,,) de;,
(5b)
du:; = C,,,, . de;,
(5c)
or
or
where
Lijki = (1 3. I. Disturbance
D)CI,,, + DC1+ a)C:,k/
and
L:+( = u:,R,,
.function
The disturbance which is influenced by the physical properties of the material and loading is considered to denote the deviation of the observed behavior with respect to those of the material’s reference states. The disturbance function, D, can be defined based on observed stresses, volumetric response, or nondestructive (ultrasonic) velocities during loading, unloading and reloading; details are given elsewhere [5]. Here, D, based on stresses, is defined as (Fig. 3)
D=
5’ - ii”
(6)
5’ _ (TL
Fig. 3. Schematic
of stress-strain
response
and disturbance.
C.S. Desai, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
365
1
0.8 0.8 L3 0.4 0.2 0 Fig. 4. Schematic
0
0.2
of disturbance
function
0.4
0.8
as dependent
0.8
1
on parameters
1.2
1.4
such as mean pressure
and temperature.
where 5 is the equivalent stress such as dzo, the second invariant of the deviatoric stress tensor, S,, or T, the octahedral shear stress. At this time, it is assumed that the RI response is stiffer than the observed response; hence, D is positive. If it is softer than the observed response, D will be negative, indicating stiffening or strengthening in the observed response; such a behavior can be included with additional considerations [9]. Now, D is expressed as
D = DU[l - exp(-A&E)]
(7)
where A and Z are material parameters and D, is the value of D corresponding to the residual condition, Fig. 3. A schematic of D vs.5, is shown in Fig. 4; here, &, is the trajectory of deviatoric plastic strains, Ez,as
(8) Here, 5, is based on irreversible or plastic strains corresponding to the observed response. However, as a simplification, it can be approximated in terms of plastic strains corresponding to the RI response based on the hierarchical (a,,) plasticity model [ 121. Comparisons with and differences between the DSC, and the classical damage models, damage models with microcrack interaction, micromechanics models and self-organized criticality (SOC) concept are discussed in ]4,51.
4. Finite element equations Based on the virtual work principle,
the finite element equations
corresponding
to Eq. (5) are derived as [7]
or
j- ([~Tlt~l[N) dV{dq’I= {Q>- I, WITb”I dv V
(9b)
where [B] is the strain-displacement transformation matrix in {de’} = [B] {dq’}, {q} is the vector of nodal displacement, and {Q} is the external load vector. Eq. (9a) or (9b) involves two unknowns, {dq’}, {w’}. The details of solution schemes are given in [7,8]; they are briefly discussed below.
C.S. Desai, W. Zhang / Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
366
4.1.
Solution schemes
Eqs. (9) can be solved using the following schemes: ( 1) An incremental iterative procedure under strain controlled condition is used. Here, the tangent matrix [C] is not positive definite, particularly in the strain softening zone. However, with the strain controlled loading, the procedure provides convergent and accurate solutions, which is illustrated through solutions of a number of problems, which are discussed later; (2) the term related to the relative stress, a:, is taken to the right-hand side, to lead to
Here, the term related to [L’] is evaluated as the load vector at the previous increment n - 1. Then, the system matrix on the left-hand side is positive definite because 0 G D 5 D,,, and a 2 0. Such a transformation does not change the positive or negative definite character of the problem; however, it can provide an approximate and simpler procedure for the finite element analysis. The algorithm in Eq. (10) is considered similar to that in the nonlocal damage models [13]; and (3) the term {v’} can be treated as an independent (constraint) unknown, together with the displacements, and the problem is treated as coupled [5,8].
5. Physical characteristics Internal microstructural changes due to factors such as relative motions between material particles, can lead to but involves discontinuous and microcracking. As a result, the material may no longer be a continuum, nonhomogeneous conditions. For consistent solutions, it is necessary to account for the discontinuous nature of the material. Incorporation in the model of the characteristic dimension, and procedures to define strain in a weighted sense over specific zones in the neighborhood of a joint, are among the procedures to account for the discontinuous material. Nonlocal formulations [ 131, inclusion of microcrack interaction [ 111 and Cosserat and gradient enrichments of the damage nodes are among the procedures [14-161. It has been shown that with the (1) inclusion of the second term, and the term dLM:, leading to the relative motion due to the third term in Eqs. (1) (2) the expression of the disturbance (A, 2 and D,,), Eq. (7), in terms of factors such as the initial density, confining pressure and characteristic length ratio (say, LIB where L is the length and d the diameter of material specimen), and (3) a simple averaging procedure for the (plastic) strains, the DSC model allows for the nonlocal effects, regularization, localization, and avoidance of spurious mesh dependence [7,8,]. Typical examples are given subsequently.
6. Microcracking,
fatigue with critical disturbance
Theories from linear fracture mechanics and certain modifications to account for nonlinear effects are often used to identify initiation and propagation of fracturing in deforming materials. However, it is felt that the microcracking due to the microstructural changes is a process in the material structure involving coalescence of microcracks that can lead to the formation and growth of distinct or finite-sized cracks, Indeed, such a state is the consequence of the microcrack coalescence into ‘weak’ zones in the material mixture due to the accumulation of plastic or irreversible deformations. The DSC allows consideration of microcracking as a process due to interacting responses of the material parts in the RI and FA states. Thus, it is possible to define critical levels of disturbance, D,, Fig. 4, at which the microcracks culminate into states that can be treated as causing initiation and final fracture or ‘failure’ involving severe intensity of cracks, at which the engineering service function of the material ceases. Once the critical state of microcracking is identified, which can lead to ‘crack’ initiation, the cracking can grow during subsequent loading with further increase in the disturbance. When the ultimate disturbance, D,, is reached, the material can be considered to have ‘cracked fully’.
C.S. Desai, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
361
Such critical states of disturbance can also define instability in the material’s microstructure. In the case of cyclic thermal loading, such a state indicates fatigue failure [6,7,17]. Another example is liquefaction in saturated sands under cyclic (dynamic) loading [ 181. It has been shown that when the critical disturbance, D,,is reached under cyclic loading, the corresponding accumulated energy or energy/cycle experiences an ‘abrupt’ change toward stabilization or saturation. Fig. 5 shows plots of disturbance vs number of thermal cycles from the analysis of an electronic chip-substrate system involving ceramic chip, Pb-Sn solder as the joining material and plastic printed wire board [7,17]. It can be seen that the value of D,, at which its rate of change approaches stabilization or saturation, corresponds with the similar state of energy/cycle around thermal cycles of about 400; which compared well with the observed cycles to failure in laboratory tests 1191. Fig. 6 shows disturbance and accumulated energy vs number of cycles for a saturated sand tested under cyclic shear loading [20]. It was found in the laboratory that initial (instability) or liquefaction occurred at about 11 cycles, and the final (ultimate) liquefaction at about 15 cycles. The disturbance and energy plots show similar values. For example, the critical disturbance, D,,corresponds to about 11 cycles at initial liquefaction, and the ultimate disturbance, D,, to about 15 cycles at the final liquefaction. Thus, the critical disturbance criteria can provide a simplified procedure for the identification of crack initiation, and growth and fatigue failure or instability. As the disturbance is computed as the part of the numerical computation, it is relatively easy to use this procedure compared to other procedures such as energy and fracture criteria.
(a)
Fig. 5. Disturbance
D-
vs tempaum cycles
and energy density per temperature
cycle in solder joint [7].
368
w. Zhang
C.S. Desai,
/ Comput. Methods Appl. Mech. Engrg.
151 (1998) 361- 376
0.8
cc %
0.8
0.4 0.2
CrcJqN
(a) D vs N
(b) Energy vs N [20] Fig. 6. Disturbance and energy vs cycles for cyclic sand behavior.
7. Adaptive
mesh strategy
It is possible to employ the adaptive mesh refinement procedures to identify the initiation and growth of localization and cracking. It can be used in conjunction with the critical disturbance criterion. Many previous adaptivity schemes have been developed for linear problems. Hence, it is necessary and required to develop them for nonlinear problems. Here, the disturbance approach is used in a new scheme to develop mesh adaptivity procedures for nonlinear problems, particularly in zones where microcracking and softening prevail. Fig. 3 shows a schematic of the stress-strain response that is affected by microcracking leading to strain softening. The mesh adaptivity procedure proposed involves two parts. For the prepeak behavior, where microcracking and disturbance are relatively small, a modified scheme based on an available procedure based on uniform error distribution [21,22] is used. However, for the postpeak behavior, where disturbance grows rapidly,
C.S. Desai, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
a new procedure given below. 7.1. Prepeak
based on critical disturbance
augmented
by the uniform
369
error scheme is used; brief details are
response
The uniform degree of error distribution (UDED) scheme [21] is used. In this scheme, a hybrid finite element procedure is used to evaluate the exact stress {g*} for the evaluation of errors, e(i) e(i) = abs( l”*2ns:“‘)
(11)
where i denotes an element, {a} is the vector of computed is the vector ‘exact’ stresses, described below. 7.2. Hybrid FE procedure The hybrid potential lr
Pm=
;
” {dTKl{d
formulation,
and {g*}
for exact stresses [23,24]
energy is expressed
I
stresses from the displacement
dV +
as
$ Iv C(E)-[Cl -'bHTKI(k~
-
[Cl
-‘b>) dV + I,h2)‘WdV
where {u} and {E} are the displacement and strain vectors in the displacement formulation, load vector, [C] is the constitutive matrix and {a} is the assumed stress given by
(12)
{Q] is the external
{a] = uw3
(13)
Here, [P] is the matrix of interpolation function and {/3} is the vector of stress parameters. four-node isoparametric element, assumed [P] is given by
For the case of
(144 where [F,] = [l r s] and r and s are local coordinates. The resultant [F,] = [I r* s*] where r* = r - a and s* = s - b, and a and b are determined based on the constant stress patch test such that
It is shown [23] that introduction of the second term with appropriate selection of [P] provides ‘softening’ of the traditional displacement formulation and leads to improved calculation of the stresses {a*}. Also, the second term represents a least-energy fit between the strains {E} and assumed stresses {u*}. The exact stresses {a*} are derived as [25] {a*] = U’IW-‘[Glbd where [H] - ’ = J, [P]T[C]-‘[P]
(15) dV; [G] = J, [PITIB] dV and {q} is the computed
displacement
vector.
7.3. Strategy for mesh refinement The objective
is to generate
a nearly optimal
mesh such that
ese
(16)
where e is the computed error, Eq. (ll), and P is the prescribed permissible approximation error for each element in an m-element mesh is given by
error, e.g. 5%. Now, the
e, = Zllhl(m”* where I(h(l= h,,,
is the norm or maximum
(17)
size of element.
The local mesh enrichment
indicator
is given by
370
C.S. Desai, W Zhang
/ Comput. Methods Appl. Mech. Engrg. 1.51 (1998)
361-376
(18) where l/E]/, is the energy norm for element
llEl(,=
Iv{e(~)~TICl _'[WI1
Then the element
i, given by
dV
(19)
size for the optimal mesh is calculated
approximately
as
(20) where hi is the original maximum size of elements, and p denotes the order of interpolation elements near singularities of order k, the new element size is estimated from
function.
For
(21) where k represents
the strength of the singularity,
which is often adopted to be equal to 0.5.
7.4. Postpeak response In this zone, the microcracking and resulting softening can be significant, and the foregoing UDED scheme alone is found to be not sufficient. Hence, the disturbance function, D, in the DSC, Eq. (7), is used as the ‘error’ indicator for remeshing. Here, a critical value of disturbance, D * = 0.5, 0.75 or 0.90 can be specified. As the disturbance is dependent on the plastic strains, Eqs. (7) (8), it represents the intensity of strain localization. Thus, the greater the value of D, the greater is the softening and localization. As a result, use of D* as the indicator for remeshing is consistent with the physical state of the structure as microcracking, and its growth follows the growth of disturbance, Fig. 4. In the postpeak region, the error indicator is expressed as DaD*
(22)
Here, the error ratio for stress is given by
IleA
(23a)
l=ile,ll where e,, is the lowest error and e, is the error in any element. following strategy:
h, = i’k =
The mesh is refined by using the
Ile,,ll’
( ) ile,ll
(23b)
k#
h,, is the size of the element results.
with the lowest error and t is an exponent.
In general, t = 2 provides
satisfactory
8. Validations and applications In order to verify the DSC model with FE procedure, a number of problems have been solved. They include ( 1) simple (metal) tension bar with imperfection for localization, (2) plate (metal) with imperfection for localization, shear banding and mesh dependence, (3) mesh dependence in a concrete footing, (4) mesh adaptivity in a footing problem, (5) localization in electronic chip-substrate systems involving Pb-Sn solders including cyclic fatigue, and (6) dynamic soil-structure interaction including liquefaction. Here, the RI state was characterized using a plasticity model (6,) in the hierarchical single surface (HISS) models [ 121, and the FA
C.S. Desai, W Bang
I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
state was defined using the constrained below.
liquid or constrained
liquid-solid
371
[3-g]. Only typical problems are given
8.1. Tension bar with impe$ection The tension-bar problem with imperfection considered by de Borst et al. [16] using the gradient enrichment with a classical plasticity model, was solved by using the DSC. Here, the RI response was simulated by using the 6,-plasticity model [4,12] and the disturbance parameters A, 2, D, were expressed as function of the ratio t/L, where t is the characteristic length and L is the length of the bar. It was shown that the DSC yields as good as or better predictions of localization and its zones, compared to those by the gradient enrichment procedure PI. 8.2. Mesh adaptivity:
Footing problem
Zienkiewicz et al. [22] used an adaptive remeshing procedure for the analysis of localization in which the material was modelled by using a classical theory of plasticity, with softening simulated by using negative hardening modulus. The authors stated that the procedures based on the acoustic tensor for the analysis of existence and localization were not reliable. They proposed that the discontinuity and localization can be captured by appropriate formulation of the finite element method, and adaptive remeshing. Here, we consider a typical problem of footing on a strain-softening material, and show that the DSC with mesh adaptivity can provide a consistent procedure to analyze localization. Fig. 7(a) shows a test problem in which uniform displacement increment = 0.001 m applied at the top, while the right-hand side boundary is free to move. The RI material is assumed to be elastoplastic with continuous
StlWS
t
(a) Footing
,/
_,___,__.._._....._.. v .. .. ..__.. .,__.. -Intact
(h) Schematic of stress-strain curve
(4 (c) Adapted meshes at typical locations (A, B, C) Fig. 7. Mesh adaptivity
for footing problem.
C.S. Desai, W. Zhang
312
yielding [4,12], and the FA material used are Elasticity:
I Comput. Methods Appl. Mech. Engrg. I51 (1998) 361-376
is assumed
to carry only the hydrostatic stress. The material parameters
E= 2x IOXMPa 1’ = 0.49
Plasticity: y = 0.30 x lo-J, p = 0 Ultimate: Phase change: II = 2. I Yielding: = lo- ‘(I. 7, = 0.2 Disturbance: :‘= IO-‘, 2 = 2.0, n = 0.8 Critical disturbance for remeshing: D*“= 0.75 The initial mesh (point A), and the adopted meshes using the procedure proposed herein, at the peak (point B) and during softening (point C), are shown in Figs. 7(b) and (c). It can be seen that the mesh adapts in the expected direction of slip and localization. Fig. 8 shows the growth of disturbance with applied displacement at two points (a, b) as the remeshing progresses. Point (a) lies within the localization zone and experiences higher disturbance beyond D* = 0.75, at peak and beyond peak. This indicates that the mesh adapts in the zones that experience critical disturbance and beyond. The disturbance at point (b), which is away from the localization zone, does not experience higher disturbance. Fig. 9 shows nondimensionalized load-displacement curves obtained by using the final mesh Fig. 7(c), Here, the analyses were performed using (a) the plasticity model and (b) the DSC with small level of softening (D, = 0.15). Various terms in Fig. 9 are: F = computed load on the top, B = width of the loaded area, c,, = cohesive strength, U, = vertical displacement at the top, and E = elastic modulus. It can be seen that the DSC model provides softer and consistent response compared to that by the plasticity model. 8.3. Localization
und mesh drpendmt~c:
Plutr MYth imperfection
A material specimen under biaxial test with width of 60 mm and height of 180 mm, solved by de Borst et al. [ f6] using gradient enrichment with the Drucker-Prager plasticity model, is shown in Fig, IO. An imperfect element with 10% reduction in the strength is provided, as shown in the figure. The plate was subjected to incremental displacement = 0.0001 m at the top. de Borst et al. [ 161 used four-node elements with bilinear displacement field, and a Hermitian (bicubic) interpolation for the plastic multiplier.
I
D -Disturbance
1
function
-
At location B
0.75.- D* At location b
A-MeshA
(a) Typical points
B-Mesh3
C-McshC
(b) D vs displacement Fig. 8. Growth of’ disturbance at typical points.
C.S. Desai, W. Bang
I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
373
U,EIBC. Fig. 9. Comparison
of predictions
The DSC model was used with RI response material parameters used are given below: Elasticity: Plasticity: Disturbance:
from DSC (0, = 0.15) and plasticity
characterized
as elastoplastic
models.
with continuous
yielding.
The
E = 11920 MPa Y = 0.49 y = 0.001, p = 0.0, Iz = 2.10 a, = lo-lo, v, = 0.20 A = 530, Z = 1.35, D, = 0.60
In the present analysis with the DSC model, two mesh layouts were used-6 X 18 mesh and 13 X 36 mesh with 8-noded quadrilateral elements. Figs. 11(a) and (b) show deformation patterns at the peak and after the peak for the 6 X 18 mesh, while Figs. 12(a) and (b) show similar results for the 13 X 36 mesh. Although the deformation, localization and shear banding patterns in the present analysis are somewhat different from those in de Borst et al. [16], the trends are essentially similar. The differences can be due to factors such as different material model and element approximation.
Incremental displacement
3 imperfection
Fig. 10. Plate with imperfection
VW.
C.S. Desai, W. Zhang
Iii
I Comput. Methods Appl. Mech. Engrg. I51 (1998) 361-376
i i i 1
(a) At PeaL
(a) At peak
@) After+ Fig.
I I.
Deformation
Fig. 12. Deformation
patterns for 6 X 18 mesh. patterns for I3 X 36 mesh.
Fig. 13 shows the load displacement curves from the two meshes, and indicates no significant mesh dependence. Here, F = load induced at the top, B = width of the plate, 6 = observed or average stress, V,op = displacement at the top and H = height of the plate. The above results show that the DSC models provide consistent and satisfactory localization and shear banding, and also avoids the spurious mesh dependence.
8.4. Mesh dependence:
Concrete
,footing
Fig. 14(a) shows the concrete block subjected to prescribed displacements. The footing is unconfined or free at the top and on the right side. The FE analysis with the DSC involved meshes with 4-, 16-, 64- and 2568-noded isoparametric elements. The strains were averaged over the patch of elements in the part of the block, ABCD. It was found that averaging over a single element can also provide satisfactory results. Fig. 14(b) shows typical results in terms of axial stress, a,., vs. strain, E,,, with the number of eIements. It can be seen that the results are essentially independent of the mesh.
Fig. 13. Load-displacement
@
6 x 18 8-nodeclernents
a
13 x 36 8-node elements
curves for 6 X 18 and 13 X 36 meshes.
C.S. Desai, W. Zhang I Comput. Methods Appl. Mech. Engrg. 151 (1998) 361-376
+
4Ebm.
1*
(I) Mesh
forcmcme
16Ebm.
-6464.
315
block
-=-256Ebm.
@) 4 “S 5
Fig. 14. Mesh dependence
for concrete
block.
9. Conclusions The disturbed state concept provides a unified approach for constitutive modelling of nonlinear materials and interfaces. It allows a hierarchical framework from which simple to sophisticated models can be adopted, depending upon the specific material and need. The DSC intrinsically provides for regularization, characteristic dimension, localization and avoidance of spurious mesh dependence. One of the important attitudes of the DSC is that it can be used to develop ‘error’ indicators for mesh adaptivity for nonlinear behavior involving microcracking and softening. Overall, the DSC is found to be a versatile modelling approach which, when implemented in numerical (finite) element procedures, can provide consistent and reliable computations.
Acknowledgments Parts of the results here were supported by Grant Nos. DDM-9102177 and DMI-9313204 from the National Science Foundation, Washington, D.C. Drs. J. Chia, C. Basaran, and T. Dishongh participated in the results presented in this paper.
References [l] C.S. Desai and Y. Ma, Modelling 16(9) (1992) 623-653.
of joints and interfaces
using the disturbed
state concept, Int. J. Numer. Analyt. Methods Geomech.
316
C.S. Desai, W. Zhang
I Comput. Methods Appl. Mech. Engrg. ISI
(1998)
361-376
[2] S. Armaleh and C.S. Desai. Modelling and testing of a cohesionless material using the disturbed state concept, Int. J. Mech. Behavior Mater. S(3) (1994) 2799295. [3] D.R. Katti and C.S. Desai. Modelhng and testing of cohesive soil using the disturbed state concept, J. Engrg. Mech., AXE 121(S) (1995) 648-658. [4] C.S. Desai. Constitutive modelling using the disturbed state concept, m: H.B. Miihlhaus, ed., Continuum Models for Materials with Microstructure, Chap. 8 (John Wiley & Sons. UK, 1995). [S] C.S. Desai and J. To& Disturbed state constitutive modelling based on stress-strain and nondestructive behavior, Int. J. Solids Struct., 33(1 I) (1996) 1619-1650. [6] C.S. Desai, J. Chia, T. Kundu and J.L. Prince, Thermomechamcal response of materials and interfaces in electronic packaging: parts I and II. J. Elect. Packaging. ASME ( 1996) in press. [7] C. Basaran, C.S. Desai and T. Kundu, Thermomechanical tinite element analysis usmg the disturbed state concept, J. Elect. Packaging. ASME (1995) submitted. ]8] C.S. Desai. C. Basaran and W. Zhang. Numerical algorithms and mesh dependence in the disturbed state concept, Int. J. Numer. Methods Engrg. ( 1996) accepted. 191 T. Dishongh and C.S. Desai. Disturbed state concept for materials and interfaces with applicattons in electronic packaging, Report to NSF, Dept. of Civil Engrg. and Engrg. Mechs., The University of Arizona, Tucson, Arizona 1996. IlO] L.M. Kachanov. Introduction to Continuum Damage Mechanics (Martinus Nijhoft Publishers. The Netherlands, 1986). [ 1I ] Z.P. Bazant. Nonlocal damage theory based on micromechanics of crack interactions, J. En&r&. Mech., ASCE 120(3) (1994) 593-617. (121 C.S. Desai. S. Somasundaram and G. Frantziskonis, A hierarchical approach for constitutive modelling of geologic materials, Int. J. Numer. Analyt. Methods Geomech. 10 (1986) 225-257. [ 131 Z.P. Bazant and F.B. Lin, Non-local yield limit degradation, Int. J. Numer. Methods Engrg. 26 ( 1988) 1805-1823. [ 14) R. de Borst and H.B. Miihlhaus, Gradient dependent plasticity: formulation and algorithmic aspects, Int. J. Numer. Met. Engrg. 35(3) (1992) 521-539. 1IS] R.H.J. Peerlings, R. de Borst, N.A.M. Brekelmans and J.H.P. de Vree, Gradient enhanced damage for quasi-brittle materials. Int. J. Numer. Methods Engrg. 39 (1996) 3391-3403. [ 161 R. de Borst, L.J. Sluys, H.B. Miihlhaus and J. Pamin, Fundamental issues in finite element analyses of localization of deformation. Engineering Comput. IO (1993) 99-121. [ 171 C.S. Desai, C. Basaran. T. Dishongh and J. Prince, Thermomechanical analysis in electronic packaging with unified constitutive model for materials and joints. Components. Packaging and Manufacturing Technology, Part B: Advanced Packaging, IEEE (1996) submitted. [IS] C.S. Desai. C. Shao and D. Rigby. Discussion to evaluation of liquefaction by energy principles, J. Geotech. Engrg., ASCE 122(3) ( 1996) 241-244. [ 191 P.M. Hall and W.M. Sherry, Materials, structures and mechanics of solder-Joints for surface-mount microelectronics technology, Proc., 3rd Int. Conf. on Techniques de Connecion en Electronique, Welding Society, Fellbach, Germany, February 1986. [20] L. Figueroa, A.S. Saada. L. Liang and N.M. Dahisaria, Evaluation of liquefaction by energy principles. J. Geotech. Engrg.. ASCE l20(9) (1994). 1211 O.C. Zienkiewicz and J.Z. Zhu. A simple error estimation and adaptive procedure for practical engineermg analysis, Int. J. Numer. Methods Engrg. 24 (1987) 337-357. (221 O.C. Zienkiewicz, M. Huang and M. Pastor, Localization problems in plasticity using finite elements with adaptive remeshing, Int. J. Numer. Analyt. Methods Geomech. I9 (1995) 127-148. (231 W. Zhang, Theoretical basis and general optimal formulations of isoparametric generalized hybrid/mixed element model for improved stress analysis, Acta Mech. Sinica 9 (1993) 277-288 (English edition). [24] W. Zhang and D.P. Chen, The patch test conditions and some multivariant finite element formulations, Int. J. Numer. Methods Engrg. ( 1996) accepted. [25] W. Zhang and C.S. Desai, Adaptive mesh refinement with the disturbed state model and application NSF, Dept. of Civil En&r&. and Engrg. Mechs., The Univ. of Arizona, Tucson, Arizona, (19%‘).
in electronic
packaging,
Report to