Microsoft Word - numero 19 articolo 1.doc K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 5 Damage localization and rupture with gradient damage models K. Pham Université Pierre et Marie Curie, Institut Jean le Rond d'Alembert, F75005 Paris pham@lmm.jussieu.fr J.-J. Marigo Ecole Polytechnique, Laboratoire de Mécanique des Solides, F91128 Palaiseau Cedex marigo@lms.polytechnique.fr ABSTRACT. We propose a method of construction of non homogeneous solutions to the problem of traction of a bar made of an elastic-damaging material whose softening behavior is regularized by a gradient damage model. We show that, for sufficiently long bars, localization arises on sets whose length is proportional to the material internal length and with a profile which is also characteristic of the material. The rupture of the bar occurs at the center of the localization zone when the damage reaches there the critical value corresponding to the loss of rigidity of the material. The dissipated energy during all the damage process up to rupture is a quantity cG which can be expressed in terms of the material parameters. Accordingly, cG can be considered as the usual surface energy density appearing in the Griffith theory of brittle fracture. All these theoretical considerations are illustrated by numerical examples. KEYWORDS. Damage Mechanics; Gradient Damage Models; Variational Methods; Crack Initiation. INTRODUCTION t is possible to give an account of rupture of materials with damage models by the means of the localization of the damage on zones of small thickness where the strains are large and the stresses small. However the choice of the type of damage model is essential to obtain valuable results. Thus, local models of damage are suited for hardening behavior but cease to give meaningful responses for softening behavior. Indeed, in this latter case the boundary-value problem is mathematically ill-posed (Benallal et al. [1], Lasry and Belytschko, [5]) in the sense that it admits an infinite number of linearly independent solutions. In particular damage can concentrate on arbitrarily small zones and thus failure arises in the material without dissipation energy. Furthermore, numerical simulation with local models via Finite Element Method are strongly mesh sensitive. Two main regularization techniques exist to avoid these pathological localizations, namely the integral (Pijaudier-Cabot and Bažant [13]) or the gradient (Pham and Marigo[10, 11]) damage approaches, see also [6] for an overview. Both consist in introducing non local terms in the model and hence a characteristic length. We will use gradient models and derive the damage evolution problem from a variational approach based on an energetic formulation. The energetic formulations, first introduced by Nguyen [9] and then justified by Marigo [7, 8] by thermodynamical arguments for a large class of rate independent behavior, constitute a very promising way to treat in a unified framework the questions of bifurcation and stability of solutions to quasi-static evolution problems. Francfort and Marigo [4] and Bourdin, Francfort and Marigo [3] have extended this approach to Damage and Fracture Mechanics. Considering the one-dimensional problem of a bar under traction with a particular gradient damage model, Benallal and Marigo [2] apply the variational formulation and emphasize the scale effects in the bifurcation and stability analysis: the instability of the homogeneous response and the localization of damage strongly depend on the ratio between the size of I http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 6 the body and the internal length of the material. The goal of the present paper is to extend a part of the results (the questions of stability will no be investigated) of [2] for a large class of elastic-softening material. Specifically, we propose a general method to construct localized solutions of the damage evolution problem and we study the influence of the constitutive parameters on the response. Several scenarii depending on the bar length and on the material parameters enlighten the size effects induced by the non local term. The paper is structured as follows. Section Setting of the gamage problem is devoted to the statement of the damage evolution problem. In Section Non homogeneous solutions of the damage problem we describe, perform and illustrate the method of construction of localized solutions and conclude by the different scenarii of responses. The following notation are used: the prime denotes either the spatial derivative or the derivative with respect to the damage parameter, the dot the time derivative, e.g. = /  u u x , ( ) = ( ) /  E dE d , = /   t . SETTING OF THE DAMAGE PROBLEM The gradient damage model e consider a one-dimensional gradient damage model in which the damage variable  is a real number growing from 0 to 1, = 0 is the undamaged state and = 1 is the full damaged state. The behavior of the material is characterized by the state function W which gives the energy density at each point x . It depends on the local strain ( )u x (u denotes the displacement and the prime stands for the spatial derivative), the local damage value ( ) x and the local gradient ( ) x of the damage field at x . Specifically, we assume that W takes the following form 2 2 20 1 1 ( , , ) = ( ) ( ) 2 2          W u E u w E (1) where 0E represents the Young modulus of the undamaged material, ( )E the Young modulus of the material in the damage state  and ( )w can be interpreted as the density of the energy dissipated by the material during a homogeneous damage process (i.e. a process such that ( ) x = 0) where the damage variable of the material point grows from 0 to  . The last term in the right hand side of (1) is the ``non local" part of the energy which plays, as we will see later, a regularizing role by limiting the possibilities of localization of the damage field. For obvious reasons of physical dimension, it involves a material characteristic length  that will fix the size of the damage localization zone. The local model associated with the gradient model consists in setting = 0 and hence in replacing W by 0W : 2 0 1 ( , ) := ( ) ( ) 2     W u E u w (2) The qualitative properties of the (gradient or local) model, in particular its softening or hardening character, strongly depend on some properties of the stiffness function ( )  E , the dissipation function ( )  w , the compliance function ( ) = 1 / ( )   S E and their derivatives. From now on we will adopt the following hypothesis, the importance of which will appear later: Hypothesis 1 (Constitutive assumptions) ( )  E and ( )  w are non negative and continuously differentiable with (1) = 0E , (0) = 0w , ( ) < 0E and ( ) > 0w for all [0,1)  . Moreover ( ) / ( )  w E is increasing to  while ( ) / ( )  w S is decreasing to 0 when  grows from 0 to 1. Let us comment this Hypothesis before to give an example 1. The interval of definition of  can always be taken as [0,1] after a change of the damage variable; 2. The condition < 0E denotes the decrease of the material stiffness when the damage grows; 3. The condition (1) = 0E ensures the total loss of stiffness when = 1 ; 4. The positivity and the monotonicity of w is natural since ( )w represents the energy dissipated during a damage process where the damage grows homogeneously in space from 0 to  ; 5. The boundedness of w is characteristic of strongly brittle materials with softening; this condition disappears in the case of weakly brittle materials with softening or in the case of brittle material with hardening; W http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 7 6. The condition of monotonicity of / w E is introduced by sake of simplicity; it is unessential and denotes that the strain does not decrease when the damage grows; 7. The condition of monotonicity of / w S is essential; it denotes the softening property, i.e. the decreasing of the stress when the damage grows. 8. The condition 1 ( ) / ( ) = 0lim    w S ensures that the material cannot sustain any stress when its damage state is 1. Example 1 A family of models which satisfy the assumptions above is the following one, when > > 0q p : 2 0 0 ( ) = (1 ) , ( ) = (1 (1 ) ) 2       q pc q E E w pE (3) It contains five material parameters: the sound Young modulus 0 > 0E , the dimensionless parameters p and q , the critical stress > 0c and the internal length > 0 whose physical interpretation will be given in Section The homogeneous solution and the issue of uniqueness. The condition > 0q is necessary and sufficient in order that ( )  E be decreasing from 0E to 0 while the condition > 0p is necessary and sufficient in order that ( )  w be increasing from 0 to a finite value. If > 0p and > 0q , then the condition >q p is necessary and sufficient in order that ( ) / ( )    w E be increasing to  while ( ) / ( )    w S is automatically decreasing to 0 . Example 2 Another interesting family of models which satisfy the assumptions above is the following one 2 0 0 0 (1 ) ( ) = , ( ) = ( ) (1 ) 2         q p E E w p q E (4) where 1p and 1q are two constants playing the role of constitutive parameters and 0 represents the critical stress of the material. The damage problem of a bar under traction Let us consider a homogeneous bar whose natural reference configuration is the interval (0, )L and whose cross-sectional area is S . The bar is made of the nonlocal damaging material characterized by the state function W given by (1). The end = 0x of the bar is fixed, while the displacement of the end =x L is prescribed to a non negative value tU (0) = 0, ( ) = 0, 0 t t tu u L U t (5) where, in this quasi-static setting, t denotes the loading parameter or shortly the ``time", tu is the displacement field of the bar at time t . The evolution of the displacement and of the damage in the bar is obtained via a variational formulation, the main ingredients of which are recalled hereafter, see [2] for details. Let Ut  and  be respectively the kinematically admissible displacement fields at time t and the convex cone of admissible damage fields:      0= : (0) = 0, ( ) = , = : (0) = 0, ( ) = 0 , = : ( ) 0,   U tt v v v L U v v v L x x   (6) where 0 is the linear space associated with Ut  . The precise regularity of these fields is not specified here, we will simply assume that there are at least continuous and differentiable everywhere. Then with any admissible pair ( , )u at time t , we associate the total energy of the bar 2 2 2 00 0 1 1 ( , ) := ( ( ), ( ), ( )) = ( ( )) ( ) ( ( )) ( ) 2 2                 L L u W u x x x Sdx E x Su x w x S E S x dx (7) For a given initial damage field 0 , the damage evolution problem reads as: For each > 0t find ( , )t tu in Ut  such that For all ( ) ( , ) , '( , )( , ) 0        t t t tU tv u v u   (8) http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 8 with the initial condition 00 ( ) = ( ) x x . In (8), '( , )( , ) u v denotes the derivative of  at ( , )u in the direction ( , )v and is given by 2 200 1 '( , )( , ) = ( ) ( ) ( ) 2                           L u v E Suv E Su w S E S dx The set of admissible displacement rates u can be identified with ( )U t , while the set of admissible damage rates  can be identified with  because the damage can only increase for irreversibility reasons. Inserting in (8) = t and = tv u w with 0w  , we obtain the variational formulation of the equilibrium of the bar, 00 ( ( )) ( ) ( ) = 0,     L t tE x u x w x dx w  (9) From (9), we deduce that the stress along the bar is homogeneous and is only a function of time = 0, = ( ( )) ( ), (0, )     t t t tE x u x x L (10) Dividing (10) by ( )tE , integrating over (0, )L and using boundary conditions (5), we find 0 ( ( )) =  L t t tS x dx U (11) The damage problem is obtained after inserting (9)-(11) into (8). That leads to the variational inequality governing the evolution of the damage 2 200 0 0 ( ) 2 ( ) 2 0                L L L t t t tS dx w dx E dx (12) where the inequality must hold for all   and becomes an equality when = t . After an integration by parts and using classical tools of the calculus of variations, we find the strong formulation for the damage evolution problem: For (almost) all 0t , Irreversibility condition: 0 t (13) Damage criterion: 2 20( ) 2 ( ) 2 0        t t t tS w E (14) Loading/unloading condition: 2 20( ( ) 2 ( ) 2 ) = 0         t t t t tS w E (15) Remark 1 We can deduce also from the variational approach natural boundary conditions and regularity properties for the damage field. In particular, we obtain that  t must be continuous everywhere. As boundary conditions at = 0x and =x L we will simply take (0) = ( ) = 0  t t L although the more general ones induced by the variational principle correspond to a combination of inequalities and equalities like (14)-(15). These regularity properties of the damage field (and consequently the boundary conditions) hold only for the gradient model ( 0 ) and disappear for the local model ( = 0 ). As long as the regularity in time is concerned, we will only consider evolution such that  tt is at least continuous. In terms of energy, we have the following property Property 1 (Balance of energy) Let us assume that the bar is undamaged and unstretched at time 0, i.e. 0 = 0 and 0 = 0U . By definition, the work done by the external loads up to time t is given by 0 ( ) =   t e s st U ds (16) the total dissipated energy in the bar during the damage process up to time t is given by 2 20 0 0 1 ( ) = ( ) ( ( )) 2     L L d t tt E x dx w x dx (17) while the elastic energy which remains stored in the bar at time t is equal to http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 9 2 0 ( ) = ( ( )) . 2   L t e tt S x dx (18) By virtue of the conditions (13)--(15) that the fields have to satisfy, the following balance of energy holds true at each time: ( ) = ( ) ( )e e dt t t   Proof. By virtue of the equilibrium condition and the definition of the elastic energy, the work done by the external load can read as  0 0 0( ) = ( ( )) ( ( )) ( )         t L L e s s s s s st S x dx S x x dx ds 2 0 0 0 = ( ) ( ( )) ( ) 2       t t L s e s ss ds S x x dxds Using the initial condition and the consistency condition in the bulk, one gets 2 00 0 0 0 ( ) = ( ) ( )          t L t L e e s s s st t w dxds E dxds  2 20 00 0 0 0 = ( ) ( ) ' ( ( ) ( ) (0) (0))                  L t L t e t s s s s s st w dx E dxds E L L ds Using once more the initial condition and the consistency condition at the boundary, we obtain the desired equality. The homogeneous solution and the issue of uniqueness If we assume that the bar is undamaged at = 0t , i.e. if 0 ( ) = 0 x for all x , then it is easy to check that the damage evolution problem admits the so-called homogeneous solution where t depends on t but not on x . Let us construct this particular solution in the case where the prescribed displacement is monotonically increasing, i.e. when =tU tL . From (11), we get = ( ) t tE t . Inserting this relation into (14) and (15) leads to 2 2( ) ( ) , = 0 2 ( ) 2 ( )               t t t t t w wt t E E (19) Since 0 = 0 , t remains equal to 0 as long as 0 = 2 (0) / (0)   t w E . That corresponds to the elastic phase. For 0> t , since / w E is increasing by virtue of Hypothesis 1, the first relation of (19) must be an equality. Therefore t is given by 1 2 = 2             t w t E and grows from 0 to 1 when t grows from 0 to  . During this damaging phase, the stress  t is given by 2 ( ) = ( )      t t t w S Since / w S is decreasing to 0 by virtue of Hypothesis 1,  t decreases to 0 when t grows from 0 to  . This last property corresponds to the softening character of the damage model. Note that  t tends only asymptotically to 0, what means that an infinite displacement is necessary to break the bar in the case of a homogeneous response. In terms of energy, the dissipated energy during the damage process is given by ( ) = ( )d tt w L Hence, it is proportional to the length of the bar. The total energy spent to obtain a full damaged state is equal to (1)w L and hence is finite by virtue of Hypothesis 1. http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 10 The non local term has no influence on the homogeneous solution which is solution both for the gradient and the local damage models. Let us now examine the issue of the uniqueness of the response. In the case of the local damage model, it is well known that the evolution problem admits an infinite number of solution. Does the gradient term force the uniqueness? The answer to this fundamental question essentially depends on the ratio / L of the internal length with the bar length, as it is proved in [1] in the case = 2p , = 0q of Example 2. Specifically it was shown that the homogeneous solution is the unique solution of the evolution problem when 0 0  L E , i.e. when the bar is small enough, while there exists an infinite number of solutions otherwise. However, when the bar is long enough, although the number of solutions is infinite, the fundamental difference between the local and the gradient models is that the length of the damaged zone is bounded from below for the gradient model while it can be chosen arbitrarily small for the local model. The main goal of the next section is to extend these results for a large class of gradient models and to study the properties of non homogeneous solutions. Let us remark that any solution of the evolution problem contains the same elastic phase, i.e. = 0t as long as 0t . Therefore, localizations can appear only when 0> t . Example 3 For the family of models of Example 1, the homogeneous response is given by 0 0 0 2 < = 0 < = = = 1                                               c c c c p q p q p q c c c c c E if if E E and if if Since > > 0q p , the stress is a decreasing function of the strain in the damaging phase what corresponds to a property of softening. For a given > 0p , the exponent of the power law goes from  to 1 when q goes from p to  . The area under the curve, i.e. the energy dissipated during the full process of homogeneous damage, is finite what corresponds to a strongly brittle behavior, see [12]. In the limit case where =p q , the damage evolves while the strain remains constant and equal to c . That corresponds to a perfectly brittle behavior, see Fig. 1. In the limit case where = 0p or = q , the stress-strain curve is an arc of hyperbola in its softening part. The area under the curve, i.e. the energy dissipated during the full process of homogeneous damage, is infinite. This change of the boundedness of the dissipated energy marks the transition between a strongly brittle and a weakly brittle behavior, see [12]. Figure 1: Left: The stress-strain response (black curve) associated with the homogeneous evolution in the case of the models of Example 1 with > > 0q p . The limit cases of perfectly brittle material ( =p q ) and weakly brittle material ( = 0p ) are in gray. Right: Graphical interpretation of the dissipated energy at the end of a homogeneous damage process. NON HOMOGENEOUS SOLUTIONS OF THE DAMAGE PROBLEM The method of construction of non homogeneous solutions et us consider one solution of the evolution problem. Setting L http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 11 0 0 2 (0) := = (0) (0)     w E S (20) we deduce from (14), that 00   t . Indeed, 0 t by virtue of (5) and (11). Then, integrating (14) over (0, )L and using the boundary conditions (0) = ( ) = 0  t t L , we obtain 2 0 0 ( ( )) 2 ( ( ))     L L t t tS x dx w x dx (21) But, since / w S is a decreasing function of  by virtue of Hypothesis 1 and since 0 t by virtue of the irreversibility condition, we have 202 ( ( )) ( ( )), (0, )     t tw x S x x L Integrating over (0, )L and inserting the result into (21) gives 2 20 t . Therefore 0 is the maximal stress that the material can sustain. The point of departure in the construction of localized damage solutions is to seek for solutions for which the equality in (14) holds only in some parts of the bar. For a given 0> t , the localized damage field will be characterized by its set =  it ti  of localization zones it where it is an open interval of [0, ]L of the form ( , ) i t i tx D x D (the independence of its length on i will be proved). In [0, ] \ tL  , the material is supposed to be sound and therefore these parts will correspond to elastic zones where = 0t . By sake of simplicity, we will not consider localization zones centered at the boundary of (0, )L , i.e. with = 0ix or =ix L . Therefore, we force the damage to vanish at = 0x and =x L . The successive steps of the construction are as follows: 1. For a given t , assuming that  t is known, we determine the profile of the damage field in a localization zone; 2. For a given t , we obtain the relation between  t and tU ; 3. We check the irreversibility condition. Damage profile in a localization zone Since t is fixed, we omit the index t in all quantities which are time-dependent. Let 0(0, )  be the supposed known stress and = ( , ) i i ix D x D be a putative localization zone. The damage field  must satisfy 2 20( ) 2 ( ) 2 = 0 .         iS w E in  (22) Since we assume by construction that the localization zone is matched to an elastic zone and since  and  must be continuous, see Remark 1, the damage field has also to satisfy the boundary conditions ( ) = ( ) = 0.  i ix D x D (23) Multiplying (22) by  and integrating with respect to x , we obtain the first integral 2 2 2 0( ) 2 ( ) =       iS w E C in  (24) where C is a constant. Using (23) and Hypothesis 1, we get 2 0= /C E and (24) can read as 2 2( ) = ( , ( ))   ix H x in  (25) with 20 0 1 ( , ) := 2 ( ) ( ) [0,1)              E H w S for E (26) Since 20 ( , ) = 2 ( ) ( )         H E w S and since, by virtue of Hypothesis 1, ( ) > 0w and 2 ( ) 1 2 ( )       S w decreases from 2 201 / > 0  to  when  grows from 0 to 1, H is first increasing from 0, then decreasing to  . Hence, there exists a unique positive value of  , say ( )  , where H vanishes: http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 12 ( , ( )) = 0, 0 < ( ) < 1    H (27) ( )  corresponds to the maximal value of the damage (at the given time), taken at the center of the localization zone. Therefore, we have obtained the Property 2 (Profile of a localized damage field) For a given stress (0, )  c , the damage field in an inner localized damage zone ( ( ), ( ))  i ix D x D is given by (31) while the half-length ( )D of the localized damage zone is finite, proportional to the internal length  and given by (28). The damage profile is symmetric with respect to the center ix of the localized damage zone, maximal at the center, the maximal value ( )  being given by (27). The damage profile is a continuously differentiable function of x , decreasing from ( )  at the center to 0 at the boundary of the localized damage zone. The matching with the undamaged part of the bar is smooth, the damage and the gradient of damage vanishing at the boundary of the localized damage zone, see Fig. 2. Figure 2: A typical damage profile in an inner localization zone and in a boundary localization zone when 0 < < c Concerning the dependence of ( )  on  , we have: Property 3 (The dependence on the stress of the maximal value of the damage in the localization zone.) When  decreases from 0 to 0, ( )  increases from 0 to 1. Proof. Indeed, let 1 2 00 < <   . Since 1 1 2 2 1 20 = ( , ( )) = ( , ( )) < ( , ( ))        H H H , and, since 1( , ) < 0 H when 1( ) < < 1   , we have 1 2( ) > ( )    . Hence ( )   is decreasing. Since 0/ ( , ) < 0   H for > 0 , we have 0( ) = 0  . Let us prove that 0 ( ) = 1lim   . Let 0= ( )lim  m (the limit exists and is positive since ( )  is decreasing). If < 1m , passing to the limit in (27) when  goes to 0 gives 0 = (0, ) = 2 ( ) m mH w , a contradiction. Hence = 1m . The size of the localization zone is deduced from (25) by integration. It depends also on  and is given by ( ) 0 ( ) = ( , )        d D H (28) ( )D is proportional to the internal length and is finite because the integral is convergent. (Indeed, ( , ) H behaves like ( , 0)     H near = 0 and like ( , ( ))( ( ))          H near = ( )   . Since ( , 0) > 0    H and ( , ( )) < 0      H , the integral is convergent.) Provided that 2 ( )L D , it is really possible to insert a localization zone of size 2 ( )D inside the bar. Concerning the dependence of ( )D on  , we obtain the following fundamental property the proof of which is not given here (it is based on a careful study of the behavior of the integral giving ( )D ): Property 4 (Dependence on the stress of the size of the localization zone.) The size ( )D of the localization zone varies continuously with  , 0( )D and 00 = ( )lim D D are finite and given by 1 0 0 0 02 0 0 2 ( ) = , = (0) 2 (0) 2 ( )         E E D D d S w w (29) http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 13 The size function ( )  D is not necessarily decreasing. In particular 0( ) < 0 dD d if and only if the following inequality holds 2 20 0(0)( (0) 2 (0)) > (0)( (0) 2 (0))       S S w S S w (30) The position ix of the center can be chosen arbitrarily in the interval [ ( ), ( )] D L D . We finally deduce from (25) that, in the localization zone, the damage field is given by the following implicit relation between x and  : ( ) | |= ( , )        i d x x H (31) It is easy to see that the damage field is symmetric with respect to the center of the localization zone, decreasing continuously from ( )  at the center to 0 at the boundary. Remark 2 The size of the localization zone and the profile of the damage field inside depends only on  . Since  is a global quantity, all the localization zones have the same size and the same profile at a given time. The maximal number of localization zones that can exist at a given time depends on the length of the bar: the longer the bar, the greater the maximal number of localization zones. Example 4 In the case of the family of models introduced in Example 2 with 1q and 1p the size of the localization zone at 0=  or 0 are given by 0 0 2 00 2 2 ( ) = , = ( )         D D p qp q q p The necessary condition (30) of growing of the localization zone when the stress decreases reads as 2 2( ) > ( ) ( 2)   q p q p q p . It is in particular satisfied for = 2q and 4p , but it is never satisfied when 1q and 1 2 p . Figure 3: Left: The damage profile for a given  and its evolution for different  in the case of the model of Example 2 with = 2q and = 4p (the lower  , the higher  ). Right: We check numerically that ( )  D is decreasing. Rupture of the bar When = 0 our previous construction of the damage profile is no more valid. Indeed, the differential system (22)-(23) becomes 20> 0 ( ) = 0 , = = 0        i iand w E in on  Integrating the differential equation over i and using the boundary conditions leads to ( ) = 0 i w dx  what is impossible by Hypothesis 1. As this is suggested by the fact that the maximal value of the damage tends to 1 when  goes to 0, one has to search a profile such that the damage field takes the value 1 at the center of the zone. Since some http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 14 quantities like the compliance function ( )  S and its derivatives become infinite when  goes to 1, the regularity of the damage field is lost and ( ) x is no more defined at = ix x but undergoes a jump discontinuity. So the differential system reads now as 20> 0 ( ) = 0 \ , ( ) = 1, = = 0         i i i iand w E in x x on  Multiplying by  the differential equation valid on each half-zone and taking into account the boundary conditions at the ends, one still obtains a first integral 2 20 ( ) = 2 ( ( )) E x w x in \i ix . Since > 0 in i , denoting by 0D the half- length of the localization zone, one necessarily has 0 0 0 0 2 ( ) ( , ) = 2 ( ) ( , )          i i i i S w in x D x S w in x x D Since ( ) = 1 ix , the jump of  at ix is equal to 02 2 (1) / S w . By integration we obtain the damage profile and the half-length of the localization zone 1 1 0 0 0 0 | |= , = 2 ( ) 2 ( )        i d d x x D S w S w (32) One can remark that this solution can be obtained formally by taking = 0 and (0) = 1 in (28)-(31). We have proved the following Property 5 (Rupture of the bar at the center of a localization zone) At the end of the damage process, when the stress has decreased to 0, the damage takes the critical value 1 at the center of the localized damage zone. The damage profile and the half length 0D of the damage zone are then given by (32). The profile is still symmetric and continuously decreasing to 0 from the center to the boundary, but its slope is discontinuous at the center. Example 5 In the cases of the family of models of Example 1, the half-length of the damage zone and the amplitude of the damage profile when the bar breaks are given by 1 1 0 0 0 | |= , = 1 1            i p p c c p dv p dv x x D q qv v For = 1p , the profile is made of two symmetric arcs of parabola: 2 0 0 | | 2 ( ) = 1 , =        i c x x x D D q For = 2p , the profile is made of two symmetric arcs of sinusoid: 0 0 | | ( ) = 1 sin , = 2 2       i c x x x D D q The greater p , the greater the size of the damage zone and the damage field, see Fig. 4. Dissipated energy in a localization zone By virtue of Property 1 and (17) the energy dissipated in an inner localization zone when the stress is  is given by ( ) 2 2 0( ) 1 ( ) = ( ) ( ( )) 2               x Di d x Di E x w x dx By symmetry, it is twice the dissipated energy in a half-zone. Using the change of variable x and (25), we obtain 2 ( ) 0 0 2 0 0 0 4 ( ) ( ( ) ) ( ) = 2 ( ) ( ( ) )               d w S S d S w S S S  (33) http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 15 Figure 4: Damage profile in the localization zone when the bar breaks, for = 4q and different values of the parameter p ( = 1 / 2,1, 2, 4p ) in the family of brittle materials of Example 1. Figure 5: The damage profile for a given t and its evolution with t by assuming that  tt is decreasing in the case of the model of Example 1 with = 2p and = 4q . The rupture occurs when = 0 t and ( ) = 1  t . We check numerically that ( )  D is decreasing. It is easy to check that ( )  d is decreasing with ( ) = 0d c while (0)d represents the dissipated energy in a localization zone during to the process of damage up to the rupture. Let us call fracture energy and denote by cG this energy by reference to the Griffith surface energy density in Griffith theory of fracture. Since (0) = 1 , we have Property 6 (Fracture energy) The dissipated energy in an inner localization zone during the damage process up to the rupture is a material constant cG which is given by 1 00 = 8 ( ) cG E w d (34) Because of the lack of constraint on the damage at the boundary, the dissipated energy in a boundary localization zone up to the rupture is / 2cG . Example 6 In the case of the family of strongly materials of Example 1 the fracture energy is given by 1 0 = 2 , = 1  pc p c p q G J J v dv p (35) Thus cG is proportional to the product of the critical stress by the internal length, the coefficient of proportionality depending on the exponents p and q . The force-displacement relation The time is fixed and we still omit the index t . Let U be the prescribed displacement, = / U L the average strain and  the stress in the bar which contains 1n localization zones. Using (11), recalling that = 0 outside the localization zones and that all localization zones have the same size 2 ( )D and the same profile, we get http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 16 0 1 0 2 ( ) = ( ( )) = ( ( ))              L L nD U S x dx n S x dx E The integral 1 ( ( )) S x dx can be transformed into an integral over the range of  by using (25). Indeed, by symmetry it can reads as 1 ( )1 2 ( ( ))    x x D S x dx . Making the change of variables x , since = ( , )  d H dx , we obtain ( ) 01 ( ) ( ( )) = 2 ( , )          S d S x dx H Recalling (28), we finally obtain the overall stress-strain relation 1= ( ) ( )     e dn L (36) with ( ) 1 0 0 0 1 ( ) = , ( ) = 2 ( ) ( , )                  e d d S E E H (37) Remark 3 For a given n , (36) gives the average strain in term of the stress. That corresponds to a curve in the   plane, parametrized by  varying from 0 to 0 . The curve = 0n is the segment corresponding to the elastic phase. Thus  can be decomposed into two terms, one associated with the elastic part of bar, the other with the localization zones. Note that 1 ( )   d depends neither on the length of the bar nor on the internal length of the material. The properties of monotonicity of the function 1 ( )   d play an important role on the presence of snap-backs in the overall response of the bar, see the next subsection. Since 1 0( ) = 0  d and since 1 ( ) > 0  d for 0<  , 1 ( )   d is necessarily decreasing in the neighborhood of 0 . We have in particular Property 7 (Behaviour of 1 ( )   d near 0 ). 1 0( ) = 0  d and 5/ 2 2 2 1/ 2 0 01 0 2 3/ 2 0 2 (0) ( ) = ( (0) 2 (0))          d S Ed d S w (38) On the other hand, the behavior of 1 ( )  d when 0/  is small is very sensitive to the constitutive parameters as it is shown in the following example and on Fig. 6. Figure 6: Graph of the function 1 ( )   d giving the contribution of a localized zone on the overall strain in the case of the model of Example 2 with = 4p and different values of q (dashed: = 1q , thick: = 2q , thin: = 3q ). Example 7 In the case of the family of models of Example 2, we have http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 17 5/ 2 2 / 21 0 0 12 3/ 2 0 0 < 2 2 ( ) ( ) = , ( ) = 2 = 2lim (( ) ) > 2                     d d p if q d p q if q d p q q p if q (39) Consequently, when 2q , the overall strain remains finite when the stress goes to 0, contrary to the homogeneous response where the strain becomes infinite. Checking of the irreversibility It remains to check that the localized damage fields that we have constructed at different values of  leads to an evolution in time which satisfies the irreversibility condition 0  . Let us reintroduce the time and the index t in the notation. Since the center of the localization zone is fixed, the condition of irreversibility is satisfied only if ( ) = ( )t t it x   is not decreasing. Since ( )   is decreasing, it is possible only if tt  is not increasing. Since ( , ( )) = 0t i i tx x D  and since ( ) > 0t x for | |< ( )i tx x D  by construction, the condition of irreversibility is satisfied only if ( )tt D  is not decreasing. That requires that ( )D  is not increasing, condition which is not automatically satisfied by the damage model, see Example 4. When this condition is not satisfied, our construction of the localized solution is no more valid. We must consider an evolution of the damage where a part of the localization zone is unloaded and reenters in a non damaging phase, the size of the still damaging part decreasing with time. To avoid such a situation we make the following hypothesis Hypothesis 2 We assume that ( )E  and ( )w  are such that ( )D  is decreasing. Note that this hypothesis is satisfied in the class of models of Example 2 when = 2q and 4p  . Under this condition, it is possible to obtain the following property Property 8 Under Hypothesis 2, in order that tt  given by (31) in a localization zone (and equal to 0 otherwise) is not decreasing it is necessary and sufficient that tt  is not increasing. Proof. We know that it is necessary, it remains to prove that it is sufficient. Let us assume that tt  is not increasing. Then ( )tt   and ( )tt D  are not decreasing. Let 1 2 2 ( ) ( ) m E L n d nL nD w      (40) Remark 4 Under Hypotheses 1 and 2, we have really obtained a damage evolution tt  which satisfies the evolution problem (13)--(15) if the bar is long enough and if we can control the loading in such a manner that the stress is continuously decreasing. But the continuity of tt  is not automatically ensured as we show in the next subsection. Size effects and the different scenarii Let us consider a loading process where =tU tL , i.e. such that the displacement of the end =x L is monotonically increasing. Consequently t corresponds to the average strain of the bar, = t , and (36) reads now 1 0 = ( )dt tt n E L     http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 18 It remains to study under which condition this relation between t and  is invertible. In other words, we have to find when the overall curve  -- does not contain snap-backs. Specifically, in order that there is no snap-back, we must have / ( ) 0,d d     . By virtue of (36), that gives an upper bound for L : 10 (0, ]0 ( ) :=inf d M d L n E nL d             (41) Of course, this condition is never satisfied when 1 ( ) d   is not decreasing. (For example, in the case of the family of models of Example 2, it is never satisfied if < 2q , cf Example 7.) Depending on the properties of the model, the length of the bar and the number of localization zones, we can obtain different scenarii. Let us consider only the situation with one localization zone ( = 1n ), to simplify the presentation. We can distinguish two cases: 1. Case 0> 2 ( )m ML D L  . In such a case we have three situations (a) For very short bars, i.e. 02 ( )L D  , no solution with localization is possible. The homogeneous solution is the unique solution; (b) For short bars, i.e. 02 ( ) < mD L L  , a solution with localization is possible just after the elastic phase, but the localization zone will progressively cover all the bar and a snap-back is possible; (c) For long bars, i.e. > mL L , a solution with localization is possible, but it is necessarily discontinuous in time because of the presence of a snap-back in the overall stress-strain response. 2. Case ML L , a solution with localization is possible, but it is necessarily discontinuous in time because of the presence of a snap-back. Example 8 In the case of the family of models of Example 2, if < 2q , then < 0ML . It is not possible to find a non homogeneous solution without snap-back. If = 2q and = 4p , then 0 0 4 72 = < = 6 17 17 m ML L      . Therefore, for bars with an intermediate length we can find a continuous in time localized solution, cf Fig. 7 (left), while for long bars a localized solution is necessarily discontinuous in time just after the elastic phase, cf Fig. 7 (right). Figure 7: Overall stress-strain relations for a law of Example 2 with = 4p and = 2q (Thin curve=homogeneous response; Thick curve=localized response). Left: For a bar of intermediate length; Right: For a long bar. CONCLUSION AND PERSPECTIVES e have proposed a method of construction of non homogeneous solutions for the one-dimensional damage evolution problem of a bar under traction. We have shown that the properties of such a localized solution is very sensitive to the parameters of the model. This strong dependence could be very interesting from an W http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true K. Pham et alii, Frattura ed Integrità Strutturale, 19 (2012) 5-19; DOI: 10.3221/IGF-ESIS.19.01 19 experimental viewpoint to identify the right law. From a theoretical viewpoint, the presence of the gradient of damage in the model has a regularizing role as expected and limits the possibility of localization of the damage since the size of a localization zone is necessarily greater than a value fixed by the internal length of the material and the others parameters of the model. Moreover, this non local term induces size effects in the response of the bar with in particular a necessarily discontinuous response and even a brutal onset of the damage after the elastic phase for long bars. When the bar is sufficiently long, our construction gives several solutions for the evolution problem, the number of localizations zones being only restricted by the length of the bar. We could probably construct many other solutions as it was made in [2]. This drastic lack of uniqueness that the introduction of a non local term has not removed needs to add a selection criterion. A good candidate is of course the stability criterion introduced in [2] and used to study the stability of homogeneous responses in [12]. The next challenge is to find which solutions among all those we have constructed satisfy the stability criterion. REFERENCES [1] Benallal, R. Billardon, G. Geymonat, In: C.S.I.M Lecture Notes on Bifurcation and Stability of Dissipative Systems, Q.S. Nguyen, editor, Springer-Verlag, (1993). [2] Benallal, J.-J. Marigo, Modelling Simul. Mater. Sci. Eng., 15 (2007) 283. [3] Bourdin, G.A. Francfort, J.-J. Marigo, J. Elasticity, 91 (2008) 5. [4] G.A. Francfort, J.-J. Marigo, Eur. J. Mech. A/Solids, 12 (1993) 149. [5] Lasry, Belytschko, Int. J. Solids Structures, 24 (1988) 581. [6] Lorentz, S. Andrieux, Int. J. Solids Struct., 40 (2003) 2905. [7] J.-J. Marigo, Nuclear Engineering and Design, 114 (1989) 249. [8] J.-J. Marigo, In: Continuum Thermodynamics: the art and science of modelling material behaviour, volume 76 of Solids Mechanics and Its Applications: Paul Germain's anniversary, G. A. Maugin, R. Drouot, and F. Sidoroff, editors, Kluwer Acad. Publ., (2000). [9] Q. S. Nguyen, Stability and Nonlinear Solid Mechanics, Wiley & Son, London, (2000). [10] K. Pham, J.-J. Marigo, Comptes Rendus Mécanique, 338(4) (2010) 191. [11] K. Pham, J.-J. Marigo, Comptes Rendus Mécanique, 338(4) (2010) 199. [12] K. Pham, J.-J. Marigo, C. Maurini, J. Mech. Phys. Solids, 59(6) (2011) 1163. [13] G. Pijaudier-Cabot, Z. P. Bažant, J. Eng. Mech., 113 (1987) 1512. http://www.gruppofrattura.it http://dx.medra.org/10.3221/IGF-ESIS.19.01&auth=true