A micromechanically-based model for predicting dynamic damage evolution in ductile polymers

A micromechanically-based model for predicting dynamic damage evolution in ductile polymers

Mechanics of Materials 33 (2001) 177±184 www.elsevier.com/locate/mechmat A micromechanically-based model for predicting dynamic damage evolution in ...

236KB Sizes 3 Downloads 104 Views

Mechanics of Materials 33 (2001) 177±184


A micromechanically-based model for predicting dynamic damage evolution in ductile polymers David H. Allen *, Chad R. Searcy Department of Aerospace Engineering, Texas A&M University, College Station, TX 77843-3141, USA Received 7 July 2000; received in revised form 6 November 2000

Abstract A new micromechanically-based viscoelastic cohesive zone model is proposed that extends the features of a quasistatic model previously developed by the authors. A statistical distribution of ®bril radii is assumed, and a critical ®bril radius is introduced as a fracture criterion. Results are presented for various loading rates, and it is found that this dynamic damage model is both physically motivated and capable of modeling complex phenomena in cohesive zones. Ó 2001 Elsevier Science Ltd. All rights reserved. Keywords: Dynamic fracture; Interfaces; Cohesive elements; Viscoelasticity

1. Introduction Solid media, such as metals, ceramics, polymers, and polymeric composites, exhibit numerous damage modes that can degrade structural performance. These damage modes have been observed for both monotonic and fatigue loading, and under both quasi-static and impact conditions. Generally, these damage modes are detrimental to design criteria by way of decreases in sti€ness, residual strength, and life. Therefore, it is attractive to produce models that can predict the e€ects of damage on these performance measures, thereby decreasing the necessity for expensive experiments. These models can also provide valuable insight for both material and structural designers in the production of more cost-e€ective products. * Corresponding author. Tel.: +1-979-845-7541; fax: +1-979845-6051. E-mail address: [email protected] (D.H. Allen).

Models that have been developed to predict these damage phenomena generally fall into one of two categories: continuum damage models (Kachanov, 1958; Rabotnov, 1969); and discrete fracture models (Grith, 1920). Continuum damage models are based on the assumption that the modeling of every discrete crack is unfeasible. Instead, a volumetric approach to damage is employed. In the case of laminated composites, the damage mechanisms are typically matrix cracks, ®ber/matrix debonds, delaminations, and ®ber fracture. Damage within the plies is generally quite small compared to the scale of the structure, so that in this case continuum damage models have been utilized with some success. These models may be formulated phenomenologically (Talreja, 1985, 1994; Allen et al., 1987; Allen, 1994), or by homogenizing a micromechanical solution (Boyd et al., 1993; Costanzo et al., 1996; Allen and Yoon, 1998). Pragmatic modeling methods for predicting the evolution of

0167-6636/01/$ - see front matter Ó 2001 Elsevier Science Ltd. All rights reserved. PII: S 0 1 6 7 - 6 6 3 6 ( 0 0 ) 0 0 0 6 9 - 7


D.H. Allen, C.R. Searcy / Mechanics of Materials 33 (2001) 177±184

several of these damage modes up to failure have also been developed (Ladeveze, 1990). In contrast, fracture models contain discrete internal boundaries. These models utilize some form of fracture criterion to predict the evolution of these internal cracks in the continuum. One form of the discrete fracture approach that has been receiving increasing attention is the socalled cohesive zone model (Dugdale, 1960; Barenblatt, 1962; Needleman, 1987; Tvergaard, 1990; Tvergaard and Hutchinson, 1993). In this approach, the zone ahead of the crack tip is viewed as a surface of dissipation having zero volume in the undeformed state. Recently, cohesive zone models have emerged as ecient methods for predicting the evolution of multiple damage modes simultaneously. Historically, the primary detraction against using cohesive zone models has been the lack of an underlying physical interpretation; however, the authors have developed a micromechanically based cohesive model with material inelasticity (Allen and Searcy, 2000, 2001). A unique feature of the approach employed by the authors is that in their previously developed micromechanical cohesive zone model, rate-dependent behavior is predicted in addition to the well known nonconvexity that is necessary in order to model crack extension. In this paper, the methodology for constructing a micromechanically based constitutive model for cohesive zones is reviewed. A brief description of the quasi-static model is presented. This is then followed by the solution of a micromechanics problem representing the dynamic damage for a single ®bril. The solution to this problem leads to an explicit form of the damage evolution law which permits a critical ®bril radius fracture criterion. An example problem will be shown where a statistical distribution of ®bril radii is assumed. 2. Cohesive zone models In viscoelastic media the nature of microscale damage zones near macroscale crack-tips can be quite complex. Therefore, it is not mathematically tractable to account for the geometry of every

microscale internal boundary in the medium of interest when attempting to predict the propagation of an existing macroscale crack. However, some simplifying assumptions can be made provided that the damage zone meets certain dimensional requirements. If the microscale damage that exists near the tip of a much larger crack has a thickness that is very small as compared its planform dimensions, then it is mathematically accurate to model this damaged volume as a cohesive zone. A physically idealized planar cohesive zone is produced by ®rst, extracting a representative volume element (RVE) from the damage zone; and then, homogenizing the thermomechanical constitutive properties of that RVE. For this model, an RVE is de®ned as a volume of the damage zone with a statistical distribution of ®bril radii where the standard deviation is suciently less than the mean ®bril radius. This homogenization process results in a nonlinear and history-dependent traction±displacement relation that governs the separation of the plane ahead of the macroscale crack tip. The aim is to construct a constitutive model for this problem. In general, the traction vector for the cohesive zone, Ti , will be represented by a nonmonotonic functional of the crack opening displacements, di , that is of the general form Ti …xk ; t† ˆ f fdi …xk ; t†g

on oXcz ;


where t is the time of interest and the braces signify the functional mapping of the entire displacement history onto the current value of the cohesive traction. Furthermore, oXcz indicates that part of the boundary upon which the cohesive zone model is applied. The mapping must be nonmonotonic so that eventual unloading occurs, thereby converting this part of the boundary of the body from a Dirichlet-type to a Neumann-type boundary, which signi®es crack extension. Functional mappings of the type of Eq. (1) are especially dicult to construct for the case where the response is nonmonotonic. Therefore, a pragmatic, although not as mathematically attractive, alternative that may be used is to propose a slight variation of the above Ti …xk ; t† ˆ f …fdi …xk ; t†g; a…xk ; t††

on oXcz ;


D.H. Allen, C.R. Searcy / Mechanics of Materials 33 (2001) 177±184

where a is an internal variable representing the microstructural damage state in the cohesive zone. Although the above form is mathematically redundant, it is a physically attractive and convenient descriptor, both experimentally and numerically. In order to complete this alternative description, it is necessary to supply an additional constitutive equation describing the evolution of the damage variable. This is usually assumed to be of the form a_ …xk ; t† ˆ g…di …xk ; t†; a…xk ; t†† on oXcz :


In this paper the alternate formulation provided by Eqs. (2) and (3) will be utilized. Deploying the above model, in conjunction with the ®eld equations used to describe the bulk medium will result in a well posed initial boundary value problem whose solution is capable of predicting crack extension (Foulk et al., 2000). Whether this procedure bears any relationship to reality depends on the accuracy of the model for the structural component at hand. However, some a priori observations can be made regarding a physically acceptable form of Eqs. (2) and (3). 3. Quasi-static model This model utilizes a micromechanical description of the damaged zone ahead of a crack tip to construct the cohesive zone model. Solid media


often dissipate substantial energy even without macroscale crack growth. Of course, a signi®cant portion of this energy dissipation can be attributed to the viscoelastic nature in the bulk material. Additional energy is usually dissipated due to the development of highly damaged microscale regions ahead of crack tips. These regions are often hackled or ®brillated in nature. In viscoelastic media the nature of microscale damage zones near macroscale crack-tips can be quite complex. Therefore, it is not mathematically tractable to account for the geometry of every microscale internal boundary in the medium of interest when attempting to predict the propagation of an existing macroscale crack. However, if the aforementioned dimensional requirements are met, then the damage zone can be idealized as a cohesive zone. For the quasi-static model, a homogenization technique is utilized that results in a nonlinear and history-dependent traction±displacement relation. This relation governs the separation of the plane ahead of the macroscale crack tip. For example, the homogenized solution for the physical representation described by a RVE in Fig. 1 would be obtained by taking boundary averages of the tractions and displacements on the upper and lower surfaces. For the case where the material in the representative volume is linear viscoelastic (except for damage) and the response is assumed to be quasi-static, the homogenization of the

Fig. 1. Damaged region ahead of a crack tip and RVE within the damaged region.


D.H. Allen, C.R. Searcy / Mechanics of Materials 33 (2001) 177±184

micromechanics problem will result in a single integral cohesive zone traction±displacement law of the general form  Z t 1 di ok cz Ti …t† ˆ ds ;  ‰1 a…t†Š  E …t s† k di os 0 …4† where Ti are the components of the traction vector in the local coordinates of the macroscale crack, a is a damage parameter resulting from the geometry of hackling and/or crazing, Ecz …t† is the viscoelastic relaxation modulus of the undamaged material, di are the opening displacements of the cohesive zone also in the local coordinates of the macroscale crack, and di are the components of a material constant representing the fracture mode mixity. Finally, k is the Euclidean norm of the opening displacements which is given by "   2  2 #1=2 2 dn dr ds k… t † ˆ ‡ ‡ : …5† dn dr ds The geometric details of the solution to the micromechanics problem will dictate the precise nature of the parameters above. For the case where the ®brils shown in Fig. 1 are assumed to be prismatic rods, PN p A pˆ1 A …t † a…t† ˆ ; …6† A where A is the total planform area of the RVE, N is the number of ®brils in the RVE, and Ap is the planform area of the pth ®bril. It is assumed that the pth ®bril can fracture, so that the range of summation in Eq. (7) gradually decreases as the ®brils break. Thus, the values for a range between zero and one. When a is zero there is no damage in the RVE and when a attains its upper limit of unity, Eq. (6) produces zero tractions ahead of the macroscale crack tip. At this point the macroscale crack tip necessarily advances. It is this nonmonotonicity in the traction±displacement relation that replaces the Grith criterion (1920) in the macroscale fracture problem, so that no energy release rate criterion need be applied at the macroscale, although the concept of a critical energy

release rate is still applicable (Costanzo and Allen, 1993, 1996). It should be pointed out that in some cases it might be dicult to measure the critical cohesive traction accurately. In this case, it may be preferable to utilize a phenomenological damage evolution law. Previous research has indicated that a phenomenological evolution law of the following general form may be applicable (Yoon and Allen, 1999): a_ ˆ c…1 a_ ˆ 0;


a† kn ;

k_ > 0 and a < 1;

k_ 6 0 or a ˆ 1;

…7† …8†

where c, m, and n are microscale phenomenological material constants and k is the aforementioned Euclidean norm of the crack opening displacement vector. 4. Dynamic damage model In forming the dynamic damage model, the general physical description proposed in the quasistatic model is retained. Suppose we return to the physical depiction shown in Fig. 1, where it is assumed that the cohesive zone is represented by a set of parallel ®brils of varying cross-section. In order to present a more accurate physical depiction of a typical ®bril, it is now assumed that the loading is applied at the end of each ®bril at such a rate that inertial e€ects cannot be neglected. Consider a single ®bril, as shown in Fig. 2. Assuming a one-dimensional response, the conservation of momentum for in®nitesimal motions is now given by or o2 u ˆq 2; ox ot


where r is the longitudinal stress in the ®bril, u is the longitudinal displacement, and q is the mass density of the ®bril. The undamaged constitutive behavior of the ®brils is assumed to be linearly viscoelastic and characterized by a convolutiontype stress±strain relationship, i.e., Z t oe rˆ …10† E…t s† ds; ot 0

D.H. Allen, C.R. Searcy / Mechanics of Materials 33 (2001) 177±184

o2 u q ˆ su; ox2 E



where the overbars represent the Laplace transform of those variables, and s is the Laplace transform variable. This results in the following solution in the Laplace transform domain:    r 1 ti sq u…x; s† ˆ u0i  exp x : …14† s2 s E

Fig. 2. A single ®bril in a cohesive zone subjected to end displacement.

where e is the axial strain and E…t† is the relaxation modulus of the ®bril. The solution to this initial boundary value problem is characterized by the propagation of a discontinuity in both the stress and strain, but not the displacement. This discontinuity propagates with wave speed, c0 , determined by satisfying jump conditions at the wave front (Achenbach, 1994) s E … 0† c0 ˆ : …11† q For coordinate locations x P c0 …t ti † the longitudinal displacement will thus necessarily be zero. A solution can be obtained to the above initialboundary value problem for certain sets of initial and boundary conditions. For example, consider the following imposed boundary condition u…0; t† ˆ u0i  …t

ti †  H …t

ti †;


where u0i is a constant, and H …t† is the Heaviside step function. This particular Dirichlet-type condition was chosen because it is amenable to the discretization of numerical algorithms. A steadystate solution can be found for this case by transforming this problem into Laplace space. Laplace transforming (9) gives

In general, performing the inverse Laplace transform of the above can be quite dicult due to the nature of the relaxation modulus, E…t†. For the case where E…t† is not homogeneous in time, the inversion may not be obtainable. In many realistic materials, the relaxation modulus has been shown to be realistically ®t with a Prony series   n X t Ej exp E…t† ˆ E1 ‡ ; …15† sj jˆ1 where E1 , Ej , and sj are material constants obtained by ®tting to experimental relaxation data. The value of n is essentially equivalent to the number of decades of time over which the material exhibits time dependence in a relaxation test. Since there is no exact Laplace transform inverse to the above equation, an approximate solution may be inferred by using the quasi-elastic method for cases where the material exhibits a small loss tangent. This method yields the following approximate solution:    Z t  x x 0 0 u… x; t† ˆ ui H s ds ui ti H t : c0 c0 0 …16† The relative response of this quasi-elastic solution is compared to the response to its elastic counterpart in Fig. 3. For an arbitrary location x, the viscoelastic response is delayed as compared to the elastic response. Furthermore, the viscoelastic response diverges due to creep. Returning to the larger problem, it is assumed that the displacement condition described in Eq. (16) applies to the entire top surface of the RVE. Thus, assuming that the upper and lower surfaces are parallel to one another, and that the


D.H. Allen, C.R. Searcy / Mechanics of Materials 33 (2001) 177±184

the pth ®bril at the critical location x ˆ L is approximated by   u…L; t† ; …18† rpf …L; t†  rpf …L; 0†  1 ‡ m  2L

Fig. 3. Comparison of elastic and viscoelastic responses for dynamic loading.

viscoelastic constitution of each ®bril is identical, the stresses in the ®brils will be identical due to conservation of momentum and kinematic compatibility. The average traction on the upper surface can then be obtained by using Cauchy's formula and boundary averaging this traction to an expression that is equivalent to Eq. (4). Suppose now that a similar de®nition of the damage parameter (6) is utilized, together with a ®bril failure model to construct a damage evolution law. A precise description of this phenomenon will require microscale experiments on the material to determine details regarding the ®bril strength and geometry in the RVE. Lacking precise details at this point, we seek only to construct a general form of the damage evolution law. Therefore, certain assumptions will be made towards this end. First, assume that the pth ®bril fractures when at some point along its length its radius reaches a f critical value, rcr , i.e., f rpf …t† 6 rcr ;

f where rcr ˆ constant:

where m is the Poisson's ratio. Note that the plus sign enters the equation due to the fact that u is negative for tensile conditions. Finally, it is assumed that the distribution of ®bril radii, rpf , is governed by a Gaussian distribution function, f …l; f†, such that ! 1 … y l†2 f …l; f† ˆ p exp ; …19† 2f2 f 2p where l is the initial mean ®bril radius for the RVE at the coordinate location x ˆ L, and f is the standard deviation of the ®bril radii. A graphical depiction of the normalized ®bril distribution used in this study is o€ered in Fig. 4. Note that a critical ®bril radius has been established. The total area under the curve is unity. The area under the curve and left of the critical ®bril radius represents the measure of ®bril breakage or damage. This leads to following expression for the damage parameter a: Z 1 a…t† ˆ 1 f …l; f† dr: …20† rcr

The time-evolution of damage of the ®brils is now introduced through Poisson e€ects. Although both


Furthermore, it is assumed that the Poisson's ratio for the ®brils is time independent. This is a fairly accurate approximation for some viscoelastic materials (Schapery, 1962). The coordinate location x where each ®bril fractures is a variable that complicates the model development. In the absence of experimental research, it will be assumed that failure occurs at the midpoint of the ®bril (x ˆ L, where 2L is the ®bril length). Thus, the radius of

Fig. 4. Normalized Gaussian distribution of ®bril radii.

D.H. Allen, C.R. Searcy / Mechanics of Materials 33 (2001) 177±184


the mean and standard deviation could potentially vary as a result of Poisson e€ects, for this study only the mean ®bril radius is a€ected. Thus,  u  l…t† ˆ l0 1 ‡ m  ; 2L


where l0 is the mean ®bril radius at t ˆ 0 and u is the displacement as described by Eq. (16). With this ®nal step, the dynamic damage model is complete. Numerical simulations of this model were performed for a single ramp displacement. The damage accumulation for three rates of loading is shown in Fig. 5. The ®gure demonstrates the direct relationship between loading rate and damage evolution rate: as the rate of loading increases, so does the rate of damage. Also of note, slower rates of loading naturally e€ect a longer delay in damage initiation; the delay in damage evolution for the displacement rate u0i ˆ 0:1 is greater than that for u0i ˆ 0:3. The dynamic damage also has a strong e€ect on the cohesive traction±displacement relationship. Retaining the same traction±displacement framework as prescribed for the quasi-static model (4), numerical simulations were performed for the same three displacement rates. The e€ect of the dynamic damage is shown in Fig. 6. In each of the three cases, the cohesive tractions are able to grow unimpeded by damage. Any initially strain-soft-

Fig. 5. Accumulated damage vs. time for dynamic cohesive zone model.

Fig. 6. Cohesive traction vs. time for delayed damage model.

ening is the result of viscoelastic relaxation. However, once the stress wave reaches the critical location x ˆ L, the onset of damage occurs, and the cohesive traction quickly decays.

5. Conclusions In this paper the authors have proposed a method for building a viscoelastic cohesive zone model that incorporates dynamic e€ects. It has been demonstrated that this method accounts for viscoelastic response exhibited in ductile media. It has also been shown that this method accounts for the delayed response in damage accumulation typically seen in dynamic load experiments. It is also interesting to note that the micromechanical model reported herein provides a sucient framework for identifying and characterizing the material properties necessary to predict crack growth using a cohesive zone model. Note that these test, while dicult, indeed can obviate the necessity to perform fracture toughness tests. This is indeed attractive, as macroscale fracture toughness test performed on ductile media almost always include energy dissipation in the bulk material and away from the crack tip, thus making these tests subject to experimental error. The utilization of this model in the prediction of crack growth in solid media remains a future effort.


D.H. Allen, C.R. Searcy / Mechanics of Materials 33 (2001) 177±184

Acknowledgements The authors are grateful for the support provided by the US Army ± Army Research Oce under research grant no. DAAG55-98-1-0119.

References Achenbach, J.D., 1994. Wave Propagation in Elastic Solids. North-Holland, Amsterdam. Allen, D.H., 1994. Damage evolution in laminates. In: Talreja, R. (Ed.), Damage Mechanics of Composite Materials. Elsevier, Berlin, pp. 79±114. Allen, D.H., Harris, C.E., Groves, S.E., 1987. A thermomechanical constitutive theory for elastic composites with distributed damage ± Part I: Theoretical development. International Journal of Solids and Structures 23 (9), 1301± 1318. Allen, D.H., Searcy, C.R., 2000. Numerical aspects of a micromechanical model of a cohesive zone. Journal of Reinforced Plastics and Composites 19, 240±248. Allen, D.H., Searcy, C.R., 2001. A micromechanical model for a viscoelastic cohesive zone. International Journal of Fracture (to appear). Allen, D.H., Yoon, C., 1998. Homogenization techniques for thermoviscoelastic solids containing cracks. International Journal of Solids and Structures 35, 4035±4054. Barenblatt, G.I., 1962. The mathematical theory of equilibrium cracks in brittle fracture. Advances in Applied Mechanics 7, 55±129. Boyd, J.G., Costanzo, F., Allen, D.H., 1993. A micromechanics approach for constructing locally averaged damage dependent constitutive equations in inelastic composites. International Journal of Damage Mechanics 2, 209±228. Costanzo, F., Allen, D.H., 1993. A continuum mechanics approach to some problems in subcritical crack propagation. International Journal of Fracture 63 (1), 27±57. Costanzo, F., Allen, D.H., 1996. A continuum thermodynamic analysis of cohesive zone models. International Journal of Engineering Science 33 (15), 2197±2219. Costanzo, F., Boyd, J.G., Allen, D.H., 1996. Micromechanics and homogenization of inelastic composite materials with

growing cracks. Journal of the Mechanics and Physics of Solids 44 (3), 333±370. Dugdale, D.S., 1960. Yielding of steel sheets containing slits. Journal of the Mechanics and Physics of Solids 8, 100±104. Foulk, J.W., Allen, D.H., Helms, K.L.E., 2000. Formulation of a three dimensional cohesive zone model for application to a ®nite element algorithm. Computational Methods in Applied Mechanics and Engineering 183, 51±66. Grith, A.A., 1920. The phenomena of rupture and ¯ow in solids. Philosophical Transactions of the Royal Society of London A 221, 163±197. Kachanov, L.M., 1958. On creep rupture time. Izvestiya Akademii Nauk SSSR, Otd. Techn. Nauk, No. 8, pp. 26± 31. Ladeveze, P., 1990. Towards a fracture theory. In: Owen, D.R.J., Onate, E., Hinton, E. (Eds.), Proceedings of the Third International Conference on Computational Plasticity, Part 2. Pineridge Press, Cambridge, UK, pp. 1369± 1400. Needleman, A., 1987. A continuum model for void nucleation by inclusion debonding. Journal of Applied Mechanics 54, 525±531. Rabotnov, Y.N., 1969. Creep Problems in Structural Members. North-Holland, Amsterdam. Schapery, R.A., 1962. Approximation methods of transform inversion for viscoelastic stress analysis. In: Proceedings of the Fourth US National Congress on Applied Mechanics. p. 1075. Talreja, R., 1985. Transverse cracking and sti€ness reduction in composite laminates. Journal of Composite Materials 19, 355±375. Talreja, R., 1994. Damage characterization by internal variables. In: Talreja, R. (Ed.), Damage Mechanics of Composite Materials. Elsevier, Berlin, pp. 53±78. Tvergaard, V., 1990. Micromechanical modelling of ®bre debonding in a metal reinforced by short ®bres. In: Dvorak, G.J. (Ed.), Proceedings of the IUTAM Symposium on Inelastic Deformation of Composite Materials. Springer, New York, pp. 99±111. Tvergaard, V., Hutchinson, J.W., 1993. The In¯uence of plasticity on mixed mode interface toughness. Journal of the Mechanics and Physics of Solids 41, 1119±1135. Yoon, C., Allen, D.H., 1999. Damage dependent constitutive behavior and energy release rate for a cohesive zone in a thermoviscoelastic solid. International Journal of Fracture 96, 56±74.