Spacetime dimensional analysis and self-similar solutions of linear elastodynamics and cohesive dynamic fracture

International Journal of Solids and Structures 48 (2011) 2076–2087

Contents lists available at ScienceDirect

International Journal of Solids and Structures journal homepage: www.elsevier.com/locate/ijsolstr

Spacetime dimensional analysis and self-similar solutions of linear elastodynamics and cohesive dynamic fracture R. Abedi, R.B. Haber ⇑ Department of Mechanical Science and Engineering, University of Illinois at Urbana-Champaign, 1206 West Green Street, Urbana, IL 61801, USA

a r t i c l e

i n f o

Article history: Received 17 September 2010 Received in revised form 6 March 2011 Available online 30 March 2011 Keywords: Dimensional analysis Similitude Elastodynamics Fracture mechanics Cohesive model Traction–separation relation Ductile-to-brittle transition

a b s t r a c t We present a dimensional analysis and self-similar solutions for linear elastodynamics with extensions to dynamic fracture models based on cohesive traction–separation relations. We formulate the problem using differential forms in spacetime and show that the scaling rules expressed in terms of forms are simpler and more uniform than those obtained for tensor representations of the solution. In the extension to cohesive elastodynamic fracture, we identify and study the influence of certain intrinsic cohesive scales on dynamic fracture behavior and describe a fundamental set of nondimensional groups that uniquely identifies families of self-similar solutions. We present numerical studies of the influence of selected nondimensional parameters on dynamic fracture response to verify the dimensional analysis, including the identification of the fundamental set for cohesive fracture mechanics. We show that distinct values of a widely-used nondimensional quantity can produce self-similar solutions. Therefore, this quantity is not fundamental, and it cannot parameterize dynamic, cohesive-fracture response. Published by Elsevier Ltd.

1. Introduction There are numerous applications of dimensional analysis1 and similarity methods2 in the fields of fluid mechanics, thermomechanics, electromagnetics and astronomy by Langhaar (1951), Sedov (1959), Huntley (1967), Isaacson and Isaacson (1975) and Szirtes (1998). Although applications of these methods to solid mechanics exist, they are less common and tend to be more limited in scope. For example, the analyses of elastodynamics in Miles (1960) and Norwood (1973) yield both the similarity variables and the complete similarity solutions, but only for specific planar configurations. Historically, the application of dimensional analysis to modeling of fatigue and fracture of materials was limited by inadequate knowledge of the significant variables that govern these phenomena (Langhaar, 1951; Wagner, 1984). Nonetheless, several applications to Linear Elastic Fracture Mechanics (LEFM) can be ⇑ Corresponding author. Tel.: +1 217 333 3826; fax: +1 217 244 5707. E-mail addresses: [email protected] (R. Abedi), [email protected] (R.B. Haber). 1 Dimensional analysis is a method by which we can deduce information about a given phenomenon by assuming only that the phenomenon can be described by dimensionally-consistent relations between a selected set of variables (Langhaar, 1951). It can generate a partial solution to nearly any problem with relatively little effort, even when a complete mathematical formulation of the problem is not available (Wagner, 1984). 2 Similarity methods attempt to represent families of solutions that share a common form when expressed in terms of certain nondimensional similarity variables. Techniques used to identify the similarity variables include dimensional analysis (Birkhoff, 1948), group theory (Morgan, 1952), universal graph methods (Szata, 2001) and integral transforms (Norwood, 1973). 0020-7683/$ - see front matter Published by Elsevier Ltd. doi:10.1016/j.ijsolstr.2011.03.013

found in the literature. Carpinteri (1982) and Wagner (1984) derive complete sets of nondimensional parameters using the Buckingham theorem (Buckingham, 1914), including the familiar Griffith’s form in the latter work. Setien and Varona (1996) discuss the computation of stress intensity factors in the context of dimensional analysis, and Szata (2001) derives fatigue crack growth rates in isotropic bodies via the universal graph method, which differs slightly from the Buckingham method. Progress in understanding the microscopic mechanisms of material failure enable new applications of dimensional analysis. For example, dimensional analysis has been used to determine size effects and the dominant failure modes in fracture. Kysar (2003) considers dislocation-induced deformations to obtain a set of nondimensional parameters that control crack-tip energy dissipation and to identify the dominant failure mode, ductile or brittle, at the onset of crack propagation. Cohesive models are among the most effective, and currently the most popular, class of continuum numerical models for dynamic fracture. They developed from the cohesive zone models first introduced by Dugdale (1960) and Barenblatt (1962). Cohesive models simulate crack initiation and extension by modeling the macroscopic effects of various nonlinear damage processes in the neighborhood of the crack tip. A constitutive relation, called a traction–separation relation (TSR), describes the tractions acting across a cohesive interface as nonlinear, bounded functions of the interface separation. A limited literature on dimensional analysis of cohesive fracture models exists. Carpinteri (1989, 1991) demonstrated that a nondimensional brittleness factor, obtained from cohesive scales and a

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

domain dimension, determines the transition from ductile to brittle failure. Tvergaard and Hutchinson (1992) examine crack growth in an elastic–plastic solid with an idealized TSR and express the critical stress intensity factor required to advance the crack in terms of other nondimensional parameters. Xu and Needleman (1994) emphasize the importance of nondimensional parameters and propose a ‘‘key dimensionless group’’ for the cohesive fracture problem. Camacho and Ortiz (1996) derive an intrinsic cohesive ~, for specific loads and a specific geometric configuratime scale s tion, in a study of dynamic spall strength. Pandolfi et al. (1999) ~ influences the minimum time step required to observe that s ensure a convergent time-stepping algorithm for dynamic fracture. Rahulkumar et al. (2000) express the macroscopic fracture energy in viscoelastic solids in terms of bulk and cohesive nondimensional parameters. Carpinteri et al. (2003) use fractal topology to further investigate the brittleness factor and its relation to scale-invariant cohesive models. In this work (see also Abedi, 2010), we present a systematic dimensional analysis for linearized elastodynamics, without restrictions on the spatial dimension and configuration of the spacetime analysis domain. We seek conditions for scaling the fields of an elastodynamic system so that the resulting system also satisfies balance of momentum, kinematic compatibility and the governing constitutive relation. We extend our dimensional analysis of linearized elastodynamics to include cohesive models of dynamic fracture, and to the authors’ knowledge, a similarly complete treatment has not been reported before. Our analysis reveals the fundamental set of nondimensional parameters required to represent all families of self-similar cohesive fracture solutions. Using this information, the solution to any problem with specific material properties, cohesive parameters and loading can be easily generalized to a family of self-similar solutions. We use the results of our dimensional analysis to derive intrinsic scales that govern cohesive fracture problems. In particular, we ~. In the case obtain intrinsic cohesive length and time scales, e L and s of self-similar solutions, the sizes of the spatial domain and the cohesive process zone scale with e L. We show that the cohesive ~ is a function of the cohesive scales for separation time scale s ~ , respectively. In fact, the value of s ~ we derive and strength, ~ d and r coincides with the one reported for the spall-strength example in Camacho and Ortiz (1996). We investigate changes in the self-similar solutions due to ~; ~ ~ and changes in various intrinsic cohesive scales, including r d; s ~ Moreover, we introduce two a scale for work of separation /. groups of nondimensional cohesive parameters that describe the problem data. Members of the first group relate the spatial and ~, while members temporal scales of the analysis domain to e L and s of the second group relate measures of the loading components to corresponding cohesive scales. We use the spacetime discontinuous Galerkin (SDG) method (Abedi et al., 2006a,b, 2009) in numerical studies of cohesive fracture mechanics to verify our dimensional analysis. These studies confirm the predicted oneto-one correspondence between specific selections of the fundamental nondimensional set and families of self-similar solutions. They also show that a widely-used nondimensional parameter in the cohesive fracture literature is not fundamental. That is, distinct values of this parameter can generate self-similar solutions, so taken alone, it does not provide a useful parameterization of cohesive-fracture response. We use differential forms and the exterior calculus on manifolds to formulate the initial and boundary-value problem of elastodynamics directly on spacetime manifolds; cf. (Abedi et al., 2006a, 2009; Miller et al., 2009). This approach has several advantages in the spacetime setting relative to conventional tensorial representations, as described in the following section. However, both for the convenience of the reader unfamiliar with

2077

this branch of mathematics and to demonstrate one of the advantages of differential forms notation, we present the main results of our dimensionless analysis in Section 3 in both conventional tensorial and forms notations. Our presentation of numerical results in Section 4 uses conventional notation and does not require knowledge of exterior calculus. 2. Formulation In this section, we formulate the initial and boundary-value problem for linearized elastodynamics, and then extend it to problems that include cohesive fracture interfaces modeled with traction–separation relations (TSRs). Our formulation systematically combines space and time quantities, and, following the development in Abedi et al. (2006a, 2009) and Miller et al. (2009), it uses the notation of differential forms on spacetime manifolds. This approach provides a direct, coordinate-free notation that can be used to express fluxes across spacetime interfaces with arbitrary orientation. For example, we combine the stress and linear momentum density fields in a single form that delivers momentum flux density across any spacetime d-manifold. This leads to concise representations of the governing equations that emphasize the notion of conservation on spacetime control volumes. Although differential forms are not widely used in solid mechanics, their use is well justified in the spacetime setting. In contrast to traditional tensor notation, for example, we can use differential forms to express the Rankine–Hugoniot jump conditions without referring to unit normal vectors on spacetime manifolds. This is a significant advantage, since no objective metric is available in classical mechanics to define magnitude and the orthogonality property for spacetime vectors. Similarly, the spacetime Stokes theorem has a simple and elegant structure when written for differential forms, but its expression for spacetime tensor fields is problematic due to the absence of an objective definition for spacetime normal vectors. Perhaps the most important advantage of the differential forms notation is that certain intrinsic relations between the spacetime mechanics fields, while obscured by tensorial notation, become clearly evident under differential forms notation. For example, in the dimensional analysis below, only four independent scalings appear when forms notation is used, while eight scalings are required with tensorial representations. 2.1. Spacetime analysis domain and differential forms notation Let D be the reference spacetime domain, an open (d + 1)-manifold in Ed  R, where d is the spatial dimension. The coordinates (x1, . . . , xd, t) = (x, t) in D are defined with respect to the ordered basis (e1, . . . , ed, et) and are understood to be material coordinates associated with the undeformed configuration of a body followed by the time coordinate. The dual basis is denoted (e1, . . . , ed, et). From here on, we adopt the standard summation convention with Latin indices ranging from 1 to d. We employ differential forms with scalar, vector, covector, tensor and cotensor coefficients and follow the convention that symbols displayed in italic bold fonts denote differential forms, while symbols in upright bold fonts denote vector, covector, tensor or cotensor fields. Thus, s and s denote, respectively, the differential form for stress and the stress tensor field, as explained below. The standard basis for spacetime 1-forms is {dx1, . . . , dxd, dt}. The spacetime volume element is the (d + 1)-form given by X :¼ dx1 ^  ^ dxd ^ dt, where ‘‘^’’ is the exterior product operator on forms; cf. (Spivak, 1965; Fleming, 1964; Arnold, 1989). We have the standard basis for d-forms, {wdxj, wdt}, in which w is the Hodge star operator and the indices of wdxj shall, from here on, be treated

2078

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

as subindices for purposes of the summation convention. These satisfy dxi ^ Hdxj ¼ dij X; dt ^ Hdxj ¼ 0; dt ^ Hdt ¼ X and dxi ^ wdt = 0. For example, in the case d = 2, we have X = dx1 ^ dx2 ^ dt, wdx1 = dx2 ^ dt, wdx2 = dx1 ^ dt and wdt = dx1 ^ dx2. The temporal insertion is defined in terms of the standard insert is a vector field tion (contraction) operator as, i :¼ iet , in which e on D with uniform value et. For d = 2, for example, we have iwdx1 = dx2, iwdx2 = dx1 and iwdt = 0. Let a and b be r- and s-forms on D, respectively, let a and b be m- and n-tensor fields on D with m P n, and let w be a scalar field on D. We write aa and bb to denote an r-form and an s-form with tensor coefficients of order m and n, respectively. The exterior product of aa and bb is the (r + s)-form with tensor coefficients of order m  n given by

aa ^ bb :¼ aðbÞða ^ bÞ;

ð1Þ

in which a(b) is the standard tensor mapping of b by a into a (m  n)-tensor field. We introduce a special 1-form with vector coefficients, dx :¼ ei dxi, and a corresponding d-form with covector coefficients, wdx :¼ eiwdxi. Let a be an r-form defined on a t-manifold Q, and let C be a smanifold, such that r 6 s 6 t and C  Q. We use ajC to denote the restriction of a to C, a r-form on C. The interpretation of the restriction operation is straightforward when s = t; it only involves restriction to the submanifold C with no alteration of the form. It has more subtle implications when s < t. First, we must sometimes interpret the restriction to C in the sense of the trace operator. This is the case, for example, when Q  D is open and C  @Q. Second, we must account for the fact that the cotangent spaces, T P C and T P Q, are distinct at any P 2 C. The cases where t  s = 1 and r = 1, s are of particular interest in this work. In these cases, the dimension of T P C is one less than the dimension of T P Q, and we must use a suitably reduced basis for covectors (s-covectors) to express the restriction as a 1-form (s-form) on C. See Section 2.3 for a more detailed discussion of restrictions of spacetime forms on D to vertical d-manifolds. We work with forms whose coefficients might suffer jumps, so we must interpret the exterior derivative operator weakly, in the sense of distribution theory (Arnold, 1989). Thus, the exterior derivative of a form, indicated by the exterior derivative operator d, generally contains a diffuse part and a jump part. Following the convention used in Miller et al. (2009), we use da to denote the diffuse part of the exterior derivative of any form a and write the jump part separately and explicitly. As shown below, the Rankine–Hugoniot jump conditions for elastodynamics arise naturally as the jump parts of exterior derivatives that appear in the spacetime governing equations for kinematic compatibility and momentum balance. Spacetime d-manifolds play an important role in our formulation, and certain special configurations of d-manifolds are of particular interest. A d-manifold with a uniform time coordinate is horizontal, and a d-manifold on which the time coordinate is independent of the spatial coordinate is vertical. The intermediate case, i.e., d-manifolds where the spatial coordinates locally parameterize the temporal coordinate, are inclined. The jump part of the exterior derivative of a across any d-manifold embedded in D is a function of the restriction of a from opposing sides of the manifold. For any open Q  D, we use undecorated symbols and symbols decorated with a superscript ‘‘+’’ to denote, respectively, the interior and exterior restrictions of a differential form to @Q. Then, for any Q  D, we define the jump in f across @Q as,

sf t@Q :¼ f

þ

j@Q  f j@Q :

ð2Þ

2.2. Mechanics fields Let the ordered set PðDÞ ¼ fQa gNa¼1 be a partition of the spacetime domain D into N open subdomains with regular boundaries S such that D ¼ a Qa . Let L2 ðDÞ and H1 ðQÞ be the Hilbertian Sobolev spaces of order 0 on D and order 1 on Q, respectively. We define a broken Sobolev space on P; V :¼ fw 2 L2 ðDÞ : wjQa 2 H1 ðQa Þ; a ¼ 1; . . . ; Ng, in which w is a covector field on D, (i.e., a 0-form with covector coefficients), and note that V admits covector fields with jumps between adjacent subdomains. 2.2.1. Kinematic quantities Let u denote the displacement covector field on D; u ¼ ui ei . The velocity v and the linearized strain E are 1-forms on D with covector coefficients given by v :¼ vdt = vieidt and E :¼ E ^ dx = E(ek) dxk, in which v = viei is the velocity covector field and the linearized strain cotensor field, E = Eijei  ej: Eij = Eji, is defined such that the value of E(ek) is a covector field. The velocity–strain is the 1-form with covector coefficients defined by,

e :¼ v þ E:

ð3Þ

2.2.2. Force-like quantities The force-like quantities include two d-forms with vector coefficients: the linear momentum density, p = pwdt, and the stress, s = s^wdx = s(ek)wdxk. The vector field p and the second-order tensor field s (under the assumption of balance of angular momentum) have the Cartesian component expansions, p = piei and s = sijei  ej:sij = sji, on D. We combine the linear momentum density and the stress in a single d-form,

M :¼ p  s;

ð4Þ

called the spacetime momentum flux. The momentum flux M acts on any oriented, spacetime d-manifold embedded in D to deliver the flux of linear momentum across the manifold. We also introduce the body force as a (d + 1)-form with vector coefficients given by b = bX, in which b = biei is the vector field on D for body force per unit mass. The corresponding form for body force per unit vol :¼ qb, where q is the mass density per unit volume is given by b ume in the reference configuration. 2.3. Governing equations 2.3.1. Kinematic compatibility The kinematic compatibility relations couple the independent displacement, velocity and strain fields. The displacement–velocity relation requires that for all Q  D,

½ðdu  v Þ ^ Hdt jQnCJ ¼ 0; u  sut@Q Hdt @Q ¼ 0;

ð5aÞ ð5bÞ

where CuJ is the jump set of u.3 The diffuse part of the displacement-velocity relation (5a) enforces the characteristic relation, u_  v ¼ 0. In lieu of (5b), we enforce the stronger condition that, for all Q  D,

   u@Q  uj@Q Hdtj@Q ¼ 0;

ð6Þ

u@Q

where is the restriction to @Q of a target displacement field that is uniquely defined on every non-vertical d-manifold embedded in D. We note that (5b) is trivially satisfied on any vertical d-manifold (because Hdtj@Q ¼ 0 on vertical manifolds), and observe that summing Eq. (6) for adjacent subdomains Q on opposing sides of any non-vertical manifold on the interior of D implies (5b). That is, in 3 The jump set of a form is the set of all points where the form’s coefficient field is discontinuous. This definition also applies to tensor fields when viewed as 0-forms.

2079

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

addition to enforcing the jump condition across any interior nonvertical manifold, we also require the jumps from both sides to vanish independently with respect to a common target field on the manifold. Consistent with the principle of causality, the target value u@Q on non-vertical manifolds is taken as the trace of u from the earlier side of the manifold or is computed directly from initial data on @D (Abedi et al., 2006a; Miller et al., 2009). The value of u@Q is immaterial, and need not be specified, on vertical parts of @Q. The tensorial representation of the velocity–strain relation is,

r v  E_ ¼ 0;

ð7Þ



in which r is the symmetric part of the spatial gradient operator. This equation holds wherever e is continuous. However, we need an associated jump condition on e to complete the exterior derivative. It is convenient to express the complete velocity–strain compatibility relation, including the jump part, in terms of the weak exterior derivative of e. For all open regions Q  D and for all symmetric, second-order tensor fields T on D, we have

ðde ^ TÞjQnCJ ¼ 0;  e set ^ T  ¼ 0; @Q

@Q

ð8aÞ ð8bÞ

in which CJe is the jump set of e, and T :¼ T ^ iwdx. It is easily shown that (8a) is equivalent to (7). Let the target velocity-strain, e@Q , be a single-valued 1-form on @Q, similar to u@Q , but computed from prescribed initial/boundary data for e or from the solution to a local Riemann problem, as explained later. Then, parallel to our treatment of (5b), we replace (8b) with the stronger condition, for all Q  D,





e@Q  ej@Q ^ Tj@Q ¼ 0:

ð9Þ

In addition to enforcing initial/boundary conditions and (8b) across interior d-manifolds (9) requires the solution to preserve the characteristic structure of the governing system across all interior interfaces between adjacent spacetime subdomains. In other words, any jumps that arise in the solution for e must be consistent with the principle of causality. The tensorial strain–displacement relation is, ðE  r uÞjQnCJ ¼ 0. u If the initial data satisfy this condition, then (5) and (8) are sufficient to enforce the strain–displacement relation everywhere. From here on, we assume that the initial data satisfy this condition, and do not include explicitly the strain–displacement relation in our governing system of equations. 2.3.2. Constitutive relation We introduce a linear transformation, denoted by C, that maps 1-forms with covector coefficients into d-forms with vector coefficients such that

M ¼ CðeðE; v ÞÞ :¼ q1ðv ÞHdt  CðEÞ ^ Hdx on D

et al., 2006a). Let CJM be the jump set of M on D. Then for all Q  D, the system

   J ¼ 0; ðdM  bÞ QnC

ð12aÞ

sMtj@Q ¼ 0;

ð12bÞ

M

enforces (11) via the Stokes theorem while accounting for possible   pX _ jumps in M. The tensorial form of (12a) is ½r  s þ b ¼ 0 on J Q n CM . Thus, (12a) is the forms representation of the equation of motion. It can also be shown that (12b) is the forms representation of the Rankine-Hugoniot condition that governs shocks in p and s. That is, the Rankine–Hugoniot condition is simply the jump part of the equation of motion. Once again, we replace the basic jump condition (12b) with a stronger condition,

M @Q  Mj@Q ¼ 0 8Q  D;

in which the target momentum flux M⁄ is defined uniquely on every d-manifold embedded in D through the initial/boundary conditions or the solution to a local Riemann problem (Abedi et al., 2006a). In addition to enforcing the boundary and initial conditions, writing the jump condition from each side with respect to M @Q requires that the solution for M preserves the characteristic structure of the elastodynamic system across all interior d-manifolds embedded in D. In other words, the solution must satisfy simultaneously the Rankine– Hugoniot condition and the principle of causality. 2.4. Initial and boundary data  as well as the initial and Given suitable prescribed data for b boundary data, and subject to the symmetry constraint on the stress tensor, the governing equations, (5), (8), (10) and (12), fully define the elastodynamic initial boundary value problem. In general, the spacetime domain boundary @D might include inclined segments that model moving boundaries. For simplicity, we focus here on a more typical problem, depicted in Fig. 1, where @D is comprised of a horizontal d-manifold at the initial time, denoted by @Di , a horizontal d-manifold at the terminal (final) time, denoted by @Df , and a collection of vertical d-manifolds that connect @Df and @Di . The latter collection, denoted by @Db , represents the spacetime extension of the boundary of a fixed reference spatial domain. In this subsection, we consider the specification of initial data on @Di , boundary data on @Db and free, unconstrained conditions on @Df for the simplified spacetime domain geometry. We stress, however, that the dimensional analysis to follow is valid for more general domain configurations, including those whose boundary includes inclined segments.

ð10Þ

in which 1 = dijei  ej, the scalar q > 0 is the mass density field, and the positive fourth-order elasticity tensor field, C = Cijklei  ej  ek  el, exhibits the standard major and minor symmetries. Eqs. (3), (4) and (10) imply the familiar tensorial component relations, pi = qdijvj and sij = CijklEkl. 2.3.3. Balance of momentum Balance of linear momentum requires that for all Q  D,

Z @Q

M

Z

 ¼ 0: b

ð11Þ

Q

Eq. (11) also implies balance of angular momentum under our a priori assumption that the stress tensor field s is symmetric (Abedi

ð13Þ

Fig. 1. Cylindrical spacetime domain and boundary partition for d = 2.

2080

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

The restrictions of the velocity–strain and the momentum flux to the horizontal, initial boundary @Di simplify according to, ej@Di ¼ Ej@Di and Mj@Di ¼ pj@Di . Thus, we express the initial conditions for displacement, strain and linear momentum density on @Di through the jump conditions (6), (9) and (13), by setting  ; e i ¼ E and M @Di ¼ p  , in which an overbar denotes preu@Di ¼ u @D  and E are kinescribed initial data. Recall that, by assumption, u matically compatible. The solution is unconstrained on the terminal boundary @Df , another horizontal d-manifold where the simplifications, ej@Df ¼ Ej@Df and Mj@Df ¼ pj@Df , hold. We model the free conditions by setting u@Df ¼ uj@Df ; e@Df ¼ Ej@Df and M @Df ¼ pj@Df . Thus, (6), (9) and (13) are trivially satisfied on @Df . The remainder of @D consists of the collection of vertical dmanifolds that comprise @Db . Since (6) is trivially satisfied on any vertical manifold, there is no need to specify u@Db . In preparation for enforcing the boundary conditions, we choose a disjoint partition of @Db into a Dirichlet part @De and a Neumann part @DM . On any vertical d-manifold C, we have (e ^ T)jC = (v ^ T)jC (for any T defined as above) and MjC = sjC. Thus, we enforce the  in (9) and Dirichlet boundary conditions by setting e@De ¼ v  is the prescribed boundary velocity. M @De ¼ sj@De in (13), where v We similarly enforce the Neumann boundary conditions via (13) and (9) by setting M @DM ¼ s and e@DM ¼ v j@DM , in which s denotes the prescribed boundary traction.

s :¼ s ^ Hdx ¼ s ^ ek Hdxk ¼ sðek ÞHdxk :

ð14Þ

However, the restriction, rjC, is a top form on C that involves only a single d-form:

sjC ¼ s ^ e1 Hdx1 ¼ sðe1 ÞHdx1 ;

ð15Þ

k

since wdx jC = wdtjC = 0 for k – 1. Recalling the Cauchy relation and that the unit covector e1 is everywhere normal to the cotangent space ðT  CÞjEd , it is clear that s(e1) and wdx1 are, respectively, the surface traction field and the singleton basis for top d-forms, both on the d-manifold C. In other words, sjC is the surface traction dform on C. Note that the restriction operation on s maps a d-form with vector coefficients on the (d + 1)-manifold D into a top d-form with vector coefficients on the d-manifold C; i.e., it involves more than a simple trace operation. e denote the collection of all the cohesive interfaces in a giLet C e – ;.4 ven elastodynamic model, and consider all Q  D 3 @Q \ C e Since C is entirely vertical, we have, uHdtje ¼ 0 for arbitrary C u; eje ¼ v je and Mje ¼ sje . For each subdomain Q, we choose C

ee

C \@Q

C

C

¼ v je

C \@Q

C

ð16Þ

:

Thus, the kinematic jump conditions, (6) and (9), are trivially satisfied, allowing independent motion of the material on each side of the cohesive interface. It only remains, then, to specify the function e according to the TSR. that uniquely determines s on C e C In general, we write,

2.5. Extension to incorporate cohesive models

se ¼ ~sðu; v ; sÞ;

Cohesive interfaces are often added to the elastodynamics problem to model crack nucleation and growth. In its most basic form, a cohesive interface is a material surface embedded in the interior of the analysis domain across which jumps in the kinematic fields are permitted and the momentum flux is described by a cohesive traction–separation relation. In the spacetime setting, as is the case for all material surfaces, cohesive interfaces are modeled as vertical dmanifolds. It is useful, then, to consider the relation between the spacetime momentum flux M and the tensorial representation of surface traction on vertical d-manifolds. Let C be a vertical d-manifold embedded in Ed  R, as illustrated in Fig. 2 for the case d = 2. To facilitate a description of the restricted cotangent space T⁄C, we define on C a local frame, fei ; et gdi¼1 : e1 ? T  CjEd , with local coordinates fxi ; tgdi¼1 , in which underlined symbols denote items referred to the local frame. In contrast to the (d + 1)-manifold, D  Ed  R, where the standard basis for d-forms in local coordinates is fHdxk ; Hdtgdk¼1 , the basis for d-forms on the d-manifold C is the singleton set, {wdx1}. Thus, the expansion of the stress form s with respect to the local frame involves a linear combination of d-forms:

e whose coefin which ~s ¼ ~sHdx1 is a top form on the d-manifold C ficient ~s is the surface traction field generated by the cohesive TSR. Thus, (13) enforces the cohesive traction condition. From here on, but without affecting the generality of our dimensional analysis, we focus on the simplest and most common constitutive relation e as, for ~s, a TSR written on C

ð17Þ

C

  ~s sute ; r ~ ; ~d; oe ¼ r ~ fð~d1 sute ; oC Þ; C

C

C

ð18Þ

in which f is a normalized TSR function that relates normalized traction to normalized separation, the covector field sute is the cohesive e and o is the local orientation on C e . CThe parameters separation on C e C r~ and ~d are reference scales for the cohesive traction and the ~ is sometimes taken to be the cohesive separation. For example, r cohesive strength, with ~ d the cohesive separation that corresponds ~ in the TSR. to r We shall use the notation, (k); k = 1, . . ., d, where (1) :¼ I, (2) :¼ II and (3) :¼ III, to map the local, coordinate directions xk; k = 1,. . .,d e into the normal (I) and tangential (II, III) modes of cohesive on C separation. The cohesive modal works of separation are given by,

~ ðkÞ :¼ /

Z

1

~sðnek Þ  ek dn ðno sum on kÞ;

ð19Þ

0

in which we have suppressed the secondary parameters in ~s, cf. (18). Eq. (18) does not imply equal cohesive strengths, equal critical separations or equal works of separation across the separation ~ and ~ modes. The cohesive parameters, r d, only provide dimensional scaling for the TSR. However, (18) does imply fixed ratios between the modal cohesive strengths and between the modal critical separations in the dimensional analysis presented in Section 3. Cohesive models can also be parameterized by the work of ~ or ~ separation and either of r d. Thus, alternative forms of (18) are given by

Fig. 2. Local coordinate system on a vertical 2-manifold C embedded in E2  R.

4 From here on, we use a superposed tilde ‘ ’ to denote a quantity defined on, restricted to or associated with a cohesive interface.

2081

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

~ ~ 1 sut ; oC Þ ¼ / hð~d1 sut ; oC Þ; ~sðsute ; oC Þ ¼ r ~ gðr ~/ e e C C C ~d

ð20Þ

~ is a reference scale for work of separation, and the values of where / e. g and h are vector fields on C 3. Dimensional analysis 3.1. Elastodynamic processes We define an elastodynamic process as an ordered set  u  ; e; M; ðx; tÞg that satisfies the governing equafu; e; M; C; b; tions of momentum balance and kinematic compatibility subject to the the corresponding jump conditions, including the initial and boundary conditions as well as the constitutive relations. We similarly define an elastodynamic cohesive process as an ordered  u  ; e; M; ðx; tÞg that additionally satisfies set fu; e; M; C; ~s; b; the cohesive TSR (18). We categorize the members of elastodynamic processes as follows. 1. Kinematic fields: The kinematic solution fields are the displacements u and the velocity–strain e = v + E. The prescribed initial  and e ¼ v  þ E comprise the kinematic and boundary data in u loading. 2. Force-like fields: The spacetime momentum flux M = p  s, is the force-like solution field. The force-like loading includes the  as well as the initial and boundary body force per unit volume b   data in M ¼ p  s. 3. Constitutive parameters: These include the linear transformation C, defined by the elasticity tensor C and mass density q. When a cohesive fracture model is included, we also have the cohesive traction ~s, whose coefficient ~s can be expressed in terms of the ~ and ~ normalized TSR function f and the cohesive scales, r d; cf. (18). 4. Spacetime coordinates: The spacetime coordinates (x,t) parameterize the spacetime analysis domain D. We consider below independent scalings of the members of an elastodynamic process and investigate the conditions under which the scaled system is itself an elastodynamic process. We also describe a similar dimensional analysis of elastodynamic cohesive fracture processes. These analyses identify key nondimensional groups that govern the response of elastodynamic systems in general and dynamic cohesive fracture models in particular, and we investigate the influence of each nondimensional group on the system response. We carry out our dimensional analyses using differential forms notation, but consider scalings of both the forms overall and of just the form coefficients, as would be done in analyses using conventional tensor notation. We show that to obtain a scaled elastodynamic process, the velocity scalings implied by the spatial and temporal coordinate scalings and by the scalings of the material properties must be identical. We also find that only four distinct scalings are required to define an elastodynamic process when we work directly with forms, while eight are required when tensor notation is used. This simplification provides new evidence of the elegance and insight afforded by the use of differential forms and the exterior calculus in spacetime mechanics analysis. 3.2. Scaled elastodynamic processes

x ¼ kx x0

t ¼ kt t 0

ð22Þ

and the resulting scalings between the bases for 1-forms, d-forms and (d + 1)-forms in the two systems are,

dxi ¼ kx dx0 i

Hdx ¼

i

0

dx ¼ kx dx

0i kd1 x kt Hdx

dt ¼ kt dt0 ;

Hdx ¼

0 kd1 x kt Hdx

ð23aÞ Hdt ¼

kdx Hdt0 ;

X ¼ kdx kt X0 :

ð23bÞ ð23cÞ

Since

@ui @ðku u0 i Þ ku @u0i ¼ ¼ ; @xj @ðkx x0 j Þ kx @x0 j

@ui @ðku u0 i Þ ku @u0 i ; ¼ ¼ @ðkt t0 Þ @t kt @t0

ð24Þ

and according to (21) and (23), if the scaled system satisfies the kinematic compatibility relations, then the kinematic forms and their tensor coefficients must transform as,

u ¼ ku u0

u ¼ k u u0 ; k v ¼ ku v 0 v ¼ u v 0 ; kt ku 0 0 E ¼ ku E E ¼ E ; kx ) e ¼ ku e0 ;

ð25aÞ ð25bÞ ð25cÞ ð25dÞ

so that ku = kv = kE = ke, while kv = ku/kt and kE = ku/kx. Thus, all of the kinematic forms share a common scaling, ku, while the three scalings for their vector and tensor coefficients are distinct. To ensure that the kinematic jump conditions are satisfied everywhere in the scaled system, it is necessary and sufficient that the prescribed kinematic initial and boundary data satisfy

 ¼ ku u 0 u

 ¼ ku u 0; u k v ¼ ku v 0 v ¼ u v 0 ; kt ku 0 0 E ¼ ku E E ¼ E ; kx ) e ¼ ku e0 :

ð26aÞ ð26bÞ ð26cÞ ð26dÞ

Next we consider conditions which guarantee that the scaled system satisfies the constitutive relations. Combining (10) with results from (23) and (25), we obtain scaling relations that govern the body force and the momentum-flux components,

 ¼ k  kd kt b 0 b b x p¼

kq ku kdx kt

p0

0 s ¼ kC ku kd2 x kt s

 ¼ k b 0 b b ; kq ku 0 p¼ p; kt kC ku 0 s¼ s: kx

ð27aÞ ð27bÞ ð27cÞ

so that kp = kqku/kt and ks = kCku/kx.

 u  ; e; M; ðx; tÞg be an elastodynamic Let fu; e; M; C; b; 0 ; u  0 ; e0 ; M 0 ; ðx0 ; t0 Þg be a scaled process, and let fu0 ; e0 ; M 0 ; C0 ; b system defined according to the relations,

a ¼ ka a0 ;

in which a and a0 are members (or member components) of the original process and the scaled system, and ka is a scaling factor specific to a. For generality, we provide independent scalings for the components of e, M and C. For each member of the process, we will use the font of the subscript a to distinguish the scaling factor for the overall form from the scaling factor for its vector or tensor coefficient. For example, kb and kb denote, respectively, the scaling fac and the vector field b.  Next, we develop tors for the (d + 1)-form b the necessary conditions on the set {ka} to ensure that the scaled system is also an elastodynamic process. According to (21) the scalings of the spacetime coordinates are5

ð21Þ

5 For simplicity, we restrict our attention to equal scalings in all spatial directions. Our results can easily be extended to accommodate distinct scales for each direction. However, the scalings of the cohesive parameters introduced in Section 3.3 would then depend on the local orientation of the spacetime manifolds that model cohesive interfaces.

2082

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

Finally, we identify conditions which ensure that the scaled system satisfies the equation of motion. We use (23), (27) and the Chain Rule to rewrite the equation of motion for the unscaled system in terms of scaled quantities:



"



  p_ X ¼ rsþb

#

kq ku ks r  s0 þ kb b 0  2 p_ 0 kdx kt X0 ¼ 0: kx kt

ð28Þ

In order for the scaled system to be an elastodynamic process, we  0  p_ 0 X0 ¼ 0. Combining this result with (27) must have ½r  s0 þ b and (28) delivers,

kC ku k2x

ð29aÞ

k2 stffiffiffiffiffi kx kC ) ¼ : kt kq

ð29bÞ

Combining these results with (27), we obtain

 ¼ k b 0 b b

p ¼ kb p0

3.3.1. Scalings  u  ; e; M; ðx; tÞg be an elastodynamic cohesive Let fu; e; M; C; ~s; b; 0 ; u  0 ; e0 ; M 0 ; ðx0 ; t 0 Þg be a scaled sysprocess, and let fu0 ; e0 ; M 0 ; C0 ; ~s0 ; b tem that is defined according to (18). We seek the conditions under which the scaled system is also an elastodynamic cohesive process. That is, in addition to the conditions for a general elastodynamic system, it must also satisfy (17) and (18). Recalling eje ¼ v je and C

s ¼ kb s0

) M ¼ kb M 0 ;

ð30Þ

C

Mje ¼ sje , and according to 9, 13 and (33), we get C

kq ku

¼ kb ¼

intrinsic scales and describe how each of these influences the elastodynamic cohesive response.

C

k~s ¼ kb ;

k~d ¼ ku :

ð34Þ

The cohesive scaling factors, k~s and k~d , are natural choices for defining the scaling of an elastodynamic cohesive process. We must select two additional scaling factors to determine all the scaling factors in (32). Here, we choose two that govern the constitutive relation for the bulk material, kc and kq.

so that kb ¼ kM ¼ kp ¼ ks ¼ ¼ That is, all the force-like forms share a common scaling factor, kb , while the scaling factors for their vector and tensor coefficients are distinct. To ensure that the momentum-flux jump conditions are satisfied everywhere in the scaled system, it is necessary and sufficient that the prescribed initial and boundary data satisfy

3.3.2. Nondimensional groups We obtain a nondimensional representation of the solution by setting

M ¼ kb M 0

in which q0 and c0 are reference scales for the wave speed and mass density of the bulk material.6 From here on, we assume that the bulk material is homogeneous, and we set q0 = q and c0 = cd, where cd is the dilatational wave speed,7

kC ku kd2 x kt

 ¼ kp p 0; p

kq ku kdx =kt .

s0 ¼ ks s0 :

ð31Þ

Let kCk be a measure of the elasticity tensor. Then the elastic pffiffiffiffiffiffiffiffiffiffiffiffiffiffi speeds are proportional to kCk=q. As a result, pffiffiffiffiffiffiffiffiffiffiffiffi kc1 :¼ kC =kq is the scaling factor for the elastic wave speeds due to the scaling of the material properties. Let us also define kc2 :¼ kx =kt as the velocity scaling implied by the spacetime coordinate mapping. For the scaled problem to be an elastodynamic process, (29b) requires kc1 ¼ kc2 . That is, the scaling of the elastic wave speeds due to scaling the material properties must match the velocity scaling implied by the scaling of the spacetime coordipffiffiffiffiffiffiffiffiffiffiffiffi nates. From here on, we write this common scaling as kc ¼ kC =kq . In summary, the scaled system is an elastodynamic process if the scalings of the spacetime coordinates and the tensorial representations of the mechanical fields satisfy wave

kx ¼ kc kt ; ku ¼ ku kb ¼

kq ku k2t

ð32aÞ ku ku kE ¼ kE ¼ ; kt kx kq ku kC ku ks ¼ ks ¼ kp ¼ kp ¼ : kx kt

kv ¼ kv ¼

ð32bÞ ð32cÞ

When working with forms, the same requirements simplify to

kx ¼ kc kt ; ke ¼ ke ¼ ku ; kM ¼ kM ¼ kb :¼

ð33aÞ ð33bÞ kq ku kdx : kt

ð33cÞ

All the scaling values in (32) and (33) are written with respect to kC, kq, kt and ku. In fact, we can choose any other set of four scalings from which to derive the rest. 3.3. Elastodynamic cohesive processes: scalings, fundamental set of nondimensional groups and intrinsic scales In this subsection, we discuss the scaling of elastodynamic cohesive processes and investigate a fundamental set of nondimensional groups that uniquely identifies families of self-similar solutions. We also identify certain dimensional parameters as

~ kr~ ¼ k~s ¼ r

()

~ kb ¼ kd1 x kt r;

k~d ¼ ~d;

kc ¼ c 0 ;

kq ¼ q0 ; ð35Þ

cd ¼ max

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi C ijkl ni nj nk nl

q

n:jnj¼1

:

ð36Þ

~0 ¼ ~ Eq. (35) implies that r d0 ¼ 1. Moreover, the wave speeds and the density are scaled by c0 and q0, respectively. In particular, c0d ¼ q0 ¼ 1. All the scales in (32) can be expressed in terms of those in (35). The nondimensional, tensorial representation of the solution is

x0 ¼

x

s~cd

u ~d  u 0  ¼ u ~d

t t0 ¼ ; s~

v v0 ¼ ~ v v 0 v ¼ ~ v

ð37aÞ

cd E; v~ cd E0 ¼ E; v~ cd s 0 0 p ¼ p s ¼ ; r~ r~ s ~  c 0 ¼ s 0 ¼ d p  s0 ¼ ; b p b qv~ r~ r~ u0 ¼

E0 ¼

ð37bÞ

ð37cÞ

~ and v~ are intrinsic time and velocity scales associated with where s the TSR, as explained below. When working with forms, the same relations simplify to

x0 ¼ u0 ¼

x

s~cd u ~d

t t0 ¼ ; s~  e e u 0 ¼ u e0 ¼ ~ e0 ¼ ~ ; d~ d d  b M M 0 ¼ b M0 ¼ ; M0 ¼ ~ ~ ~ m m m

ð38aÞ ð38bÞ ð38cÞ

6 ~ and ~d are reference scales for the cohesive traction and the cohesive Recall that r separation in the TSR, cf. Section 2.5. 7 However, we can extend our results to cover non-homogeneous materials by substituting c0 and q0 for cd and q.

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

~e ~s ~ :¼ r where m L d1 is the cohesive model’s intrinsic momentum ~ is an intrinsic length scale for the cohesive scale, in which e L ¼ cd s process-zone size (see below). Eq. (37) demonstrates that the domain dimensions, the loads and the material properties must scale in a particular way to obtain a cohesive elastodynamic process. For example, the nondimen 0 and s0 , must alsional Dirichlet and Neumann boundary data, v ~ and r ~. ways scale, respectively, with v  be measures for the load data, ; v  ; E; p ; r  and b Let u  Also, let L and s be length and time measures  ; s and b. ; v  ; E; p u for the spacetime domain D, and assume that the ratios,  v  =~  =v ~ ; L=e ~ ; E= e =p ~; r  =r ~; s ~b= ~, are fixed for a family u d; v E; p L and s=s of scalings of a given cohesive elastodynamic process. Then, according to our previous discussion, the family of scaled systems are cohesive elastodynamic processes with self-similar solutions. We define a fundamental nondimensional set as a collection of independent, nondimensional parameters having a one-to-one correspondence between families of self-similar solutions and the values of the set. That is, fixing the values of all members of the fundamental set identifies a unique family of self-similar solutions, and any two problems with self-similar solutions must have identical values for all members of the fundamental set. Moreover, the members of the set must be independent in the sense that no member of the set can be expressed as a function of the other members. Accordingly, the set of nondimensional ratios listed in the preceding paragraph constitutes a fundamental nondimensional set for TSR-based cohesive elastodynamic fracture. In general, the specific choice of non-dimensional parameters in the fundamental set is not unique. However, the above list of ratios is a natural choice because each member is the ratio of a spacetime domain or loading measure to a corresponding cohesive scale. 3.3.3. Intrinsic cohesive scales ~ 8 in ~; ~ We showed in Section 2.5 that any two members of fr d; /g, combination with the bulk material properties, q and cd, are sufficient to describe the scaling properties of a TSR; cf. (20). Here we ~ and ~ choose r d as the independent parameters. In general, ~ ¼ cr ~~ / d, in which the dimensionless factor c depends on the choice of the independent cohesive scales and on the type of the cohesive model. We simplify this relation for purposes of dimensional analy~ ¼r ~~ ~ ; ~d and sis to / d. The reference scales for the TSR are then, r

~ ¼r ~ ~d; /

qc ~d s~ ¼ ~d ; r ~d r~ v~ ¼ ~ ¼ ; s qc d v~ r~ r~ e ; E¼ ¼ 2/ cd qcd kCk r~ ~ ¼ qv~ ¼ ; p cd qc2 ~d kCk/~ e ~¼ d / : L ¼ cd s r~ r~ 2

ð39aÞ ð39bÞ ð39cÞ ð39dÞ ð39eÞ ð39fÞ

~; v ~ and e ~; e in which s E; p L are, respectively, dimensional scales for time, velocity, strain, linear momentum density, and length. The second, alternative expression for e L in (39f) arises from (39a) and the definition of cd in Section 3.2. The above equations and cohesive scales are valid for general mixed-mode conditions. While they remain valid for single-mode problems, different choices for the wave speeds and cohesive parameters that define the cohesive scales might be more appropriate in a particular setting. For example, in pure mode-II and 8

~ is a reference scale for work of separation; cf. Section 2.5. Recall that /

2083

mode-III settings, the shear-wave speed cs is the natural choice for the wave speed, and we would substitute cs for cd in all equations following (36). Moreover, the appropriate normal or tangential components of cohesive traction and separation should be used ~ and ~ for r d in each case. Next, we discuss the influence of the dimensional parameters on various aspects of the cohesive response. The cohesive parameters and the bulk material properties determine the cohesive length scale, e L. Although e L and ~ d share a common dimension, length, they have distinct physical interpretations. The cohesive length scale e L reflects the cohesive process-zone size, while the cohesive separation scale ~ d influences the amplitude of the displacement field. Both cohesive length scales are independent of the spatial length scale L of the analysis domain. There are several estimates in the literature for the quasi-static process-zone size (Rice, 1968, 1980). For isotropic materials they can be presented in the combined form,

K0 ¼ fp

l 1m

~ / ~ r2

! ¼ fe L;

ð40Þ

~ is set where l and m are the shear modulus and Poisson ratio, r equal to the cohesive strength and the constant factors, f and f, depend on the details of the TSR and on the convention that defines the process-zone size. For instance, under mode I conditions, Rice (1968) estimated f = 1/4 for the Dugdale model (Dugdale, 1960), and proposed f = 9/16 for potential-based TSRs (Rice, 1980). The process-zone size in the dynamic setting is proportional to K0. Furthermore, it depends on the crack speed, v^ , and approaches zero as the crack speed approaches the Rayleigh wave speed (Freund, 1990; Yu and Suo, 2000). There is a close connection between L=e L and the brittleness of the response. We are not aware of any discussion of this issue in the dynamic fracture literature, but we note here two related results from studies of quasi-static fracture. Carpinteri (1989), Carpinteri (1991) and Carpinteri et al. (2003) define a nondimensional ~ r ~ hÞ, in which h is the beam brittleness number for beams, sE :¼ /=ð depth. They demonstrate that the nondimensional quantity, ~ r ~ 2 lÞ, in which eu :¼ r ~ =E and l is the beam span, sE h=ðeu lÞ ¼ E/=ð determines the transition from ductile to brittle fracture; cf. Eqs. (15) and (17) in Carpinteri et al. (2003). Recalling (39f), we have sE h=ðeu lÞ / e L=L when l is set equal to the domain length scale, L. Harder (1991) investigates the ductile-to-brittle transition for an initial crack of length a in an infinite domain with tensile loading normal to the crack direction. He defines a brittleness modulus, b = a/lc, in which the characteristic cohesive length scale is given by ~ r ~2 / e lc ¼ E/= L; cf. (39f). Based on his analysis, Harder predicts a transition from ductile to brittle fracture when b exceeds a critical brittleness modulus, bc. Thus, the reciprocal of the nondimensional quantity of Carpenteri et al. and Harder’s brittleness modulus are proportional to the dimensionless parameter, L=e L (cf. Section 3.3.2), when the beam span l and the crack length a are used, respectively, as the domain length scale, L. Next, we hold fixed the material parameters, C and q, and hence the wave speeds, to investigate the influence of the cohesive parameters in (39) on the system response. Consequently, the pairs ~; e ~; v ~ Þ are interchangeable for this study. Subject to these ðs LÞ and ðr conditions, we categorize the applied loads in (37) into three groups, each with a distinct scaling. First, the initial displacement  is proportional to qv  scales with ~ ~; ~ =s u d; second, the body force b and third, all boundary conditions as well as the initial linear ~. momentum density and the initial strain are proportional to r ~ and discuss ~; r ~; ~ We next consider the cohesive scales, s d and /, the impact on the system scalings when the bulk material parameters and one cohesive scale are held fixed and the others ~ be fixed, are allowed to vary. First, let the cohesive time scale s

2084

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

and let the separation scale ~ d vary. This corresponds to uniform scaling of the TSR on both the separation and stress axes, because ~ /~ ~ implies r ~ /~ ~/~ fixed s d. We also have v d and / d2 . By construction and by (32a), the spatial and temporal dimensions of the domain are unchanged. However, all of the load values scale linearly with ~ d. ~ is fixed and ~ Now consider the case where r d is varied. This corresponds to scaling the separation axis in a graph of the TSR. In this ~ /~ ~; / case, e L; s d, and the displacement and initial-displacement fields (37b) also scale linearly with ~ d. That is, both the undeformed and deformed spacetime configurations as well as the temporal axis scale proportionally to ~ d, so that e and M do not change. This invariance is also evident in the initial–boundary conditions, which  however, scales as 1=~ are unchanged. The body force b, d. ~/r ~ ~ Finally, consider varying r with ~ d held fixed. This yields / ~ / 1=r ~ . Now the scaling is on the stress axis of the TSR and e L; s graph, with the separation axis fixed. Clearly, the displacement and the initial displacement do not change. However, all of the ~ , and we other initial–boundary conditions scale linearly with r /r ~ 2. have b ~ is fixed. This situation has Finally, we consider the case where / significant practical importance, since we are often able to measure directly a material’s specific fracture energy, from which we can ~ We must then identify the remaining cohesive paramcompute /. ~ held fixed. Let r ~ be the control eters with the computed value of / ~, ~ ~ , and ~/r for the remaining cohesive parameters. Then v d / 1=r s~; eL / r~ 2 . The displacement solution and the initial displacements ~ , while the spatial and temporal coordinates are proportional to 1=r ~ 2 . All initial–boundary conditions scale linearly with r ~, scale as 1=r /r ~ 3. and b 4. Numerical verification This section presents numerical results that verify our dimensional analysis. We demonstrate the one-to-one correspondence between families of self-similar solutions and specific values of the fundamental set. That is, all of the nondimensional values in the fundamental set are equal for two problems if and only if their solutions are self-similar. Clearly, not all nondimensional parameters pertinent to a given physics are fundamental. Typically, a dimensional analysis begins with the identification of a set of key physical parameters. Then a set of nondimensional parameters are derived by applying dimensional consistency arguments to this set. For example, the ~ H :¼ r ~ =E, is a nondimensional quannormalized cohesive strength, r tity that Xu and Needleman (1993) proposed to parameterize general CFM response and that is often used to characterize the artificial compliance of intrinsic cohesive models; see, for example (Zhang et al., 2007). However, we demonstrate below that distinct values of this quantity can produce self-similar solutions; hence it is not fundamental and cannot distinguish cohesive models with ~ H is proportional to e distinct response. Indeed, the facts that r E (see (39d)) and that the ratio, E= e E, is fundamental (cf. Section ~ H cannot be fundamental. In fact, the normalized 3.3.2) imply that r  0 :¼ r  =r ~ or, equivalently, its reciprocal, is the funload intensity, r damental parameter that characterizes cohesive strength; we com~ H in this numerical study. pare and contrast it with r Numerical modeling of dynamic CFM can be very challenging. A robust numerical model must resolve sharp wave fronts and capture very large gradients in the vicinity of a cohesive process zone whose size vanishes, in the limit, as the crack speed approaches the Rayleigh wave speed. We use the SDG finite element method, as described in Abedi et al. (2006a,b, 2009), to meet these requirements. The SDG method’s element-wise balance properties, linear computational complexity, and powerful spacetime adaptive meshing capabilities, combined with two adaptive error indicators

that limit numerical dissipation in the bulk and ensure an accurate rendering of the TSR along the cohesive interface, provide the necessary resolution. These adaptive SDG models accurately determine the extent of the process zone, and they do not suffer the nonphysical crack speeds generated by some models due to liftoff effects; see, for example (Xu and Needleman, 1994). 4.1. Problem description Fig. 3 shows a plate subjected to in-plane stress loading; it contains a half-plane crack, but is otherwise unbounded. We choose the coordinate system (x1, x2), so that the semi-infinite crack coincides with the negative part of the x2-axis and the crack tip position is at the origin at times t < 0. We prescribe traction-free conditions on the crack faces at all times and a suddenly applied, spatially uniform, mode-I far-field traction loading with intensity s1 . Sharp wavefronts generated by the far-field loading reach the crack plane from both sides at t = 0. We use a bulk material model that approximates the elastic properties of polymethyl methacrylate (PMMA): Young’s modulus, E = 3.24 GPa; Poisson’s ratio, m = 0.35; and mass density, q = 1190 kg/m3. The corresponding dilatational, shear and Rayleigh wave speeds are cd = 2090 m/s, cs = 1004 m/s, and cR = 938 m/s. We introduce a cohesive interface along an assumed crack path, {(x1, x2): x1 = 0, x2 > 0}, to model fracture. The extent of the computational domain is taken sufficiently large to prevent reflected waves from affecting the cohesive response. We exploit symmetry about the x2-axis and use a dimensionally consistent regularization of the Heaviside temporal variation of s1 , cf. (Abedi, 2010), to reduce computational cost. Our cohesive constitutive relation is a specialization of the exponential TSR developed by Xu and Needleman (1994), in which the critical separations and the works of separation are equal for the normal and tangential directions. Although more physically realistic models are available, this simple idealization suffices for purposes of verifying our dimensional analysis. The function f in (18) then takes the form,



su1 t f su2 t

 ¼

D1 expð1  D1  D22 Þ 2D2 ð1 þ D1 Þ expð1  D1  D22 Þ

! ;

ð41Þ

~ are, respectively, normalized in which D1 :¼ su1 t=~ d and D2 :¼ su2 t=d separations in the normal and tangential directions. The cohesive ~ , is attained at the critical normal strength in the normal direction, r separation, su1 t ¼ ~ d. Following Xu and Needleman (1994), we use ~ d ¼ 4:0  104 mm. The analysis domain is unbounded in space and time, and all loading and initial data vanish except the far-field traction, s1 .  0 is the only relevant fundamental nondimensional paramThus, r eter in this problem; cf. Section 3.3.2. In the absence of the crack, the normal traction acting on the crack plane at t = 0+ is 2s1 , due to additive interference between the two incident waves. Thus, r ¼ 2s1 is an appropriate choice for the load scale, and we obtain r 0 ¼ 2s1 =r~ for the normalized load intensity.

Fig. 3. Domain and load description for a plate containing a half-plane crack.

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

2085

Table 1 Normalized load intensities and cohesive strengths for analysis cases. Case I II III

r 0

r~ H

(normalized load intensity)

(normalized cohesive strength)

101 101/4 101/4

101 101 102

We analyze three cases to check the predictions of our dimensionless analysis and to investigate the suitability of the normal~ H , for parameterizing cohesive response. Each ized strength, r  0 and r ~ H , as listed in Table case is identified by specific choices for r 1. Cases I and II have distinct values of the fundamental parameter r 0 to model low- and high-amplitude loading. According to our dimensional analysis, they should generate distinct, non-self-similar solutions. On the other hand, they share the same value of r~ H used by Xu and Needleman (1994), which suggests, if our analysis is correct, that this widely-used nondimensional parameter cannot distinguish problems with fundamentally different solu 0 . Thus, our dimentions. Cases II and III share the same value of r sional analysis predicts that they produce self-similar solutions ~ H. despite having distinct values of r

Fig. 5. Normalized crack-tip position versus normalized time. Normalized results for self-similar Cases II and III coincide.

4.2. Comparison of results Figs. 4–6 present nondimensional numerical results for cracktip velocity and trajectory and for the size of the cohesive process

Fig. 6. Normalized CPZ size versus normalized crack-tip velocity. Normalized results for self-similar Cases II and III coincide.

Fig. 4. Crack-tip velocity.

zone (CPZ). The crack-tip position, denoted as a, is the current x2coordinate of the crack tip, and the crack-tip velocity is the time _ According to common convention, we define the derivative, a. CPZ size, denoted by K, as the distance between the nominal crack-tip position and the trailing edge of the CPZ, where the nominal crack-tip position is the location where the critical separation, ~ d, is attained. The Xu and Needleman TSR generates tractions that only vanish asymptotically for large separations, so we define the trailing edge of the CPZ as the location behind the crack tip where ~ . Alternathe normal traction falls to 1% of the cohesive strength, r tive definitions of K are possible, but reasonable options generally scale K by a factor that is only Oð1Þ. Fig. 4(a) presents histories of crack-tip velocity in dimensional form. In all cases, the crack velocity approaches the Rayleigh wave speed, cR = 938 m/s. We observe apparently distinct results for the three cases, with the fastest acceleration and earliest onset of crack propagation occurring in Case II, with high-amplitude loads and higher cohesive strength, and the slowest acceleration and latest onset in Case I, with low-amplitude loads and higher cohesive strength. Fig. 4(b) presents the same results in nondimensional form, ~ (cf. (37a)), and normalized crack-tip using normalized time, t 0 ¼ t=s _ R . In this format, the results for Cases II and III velocity, a_ 0 :¼ a=c are indistinguishable, with identical, immediate onsets and accelerations, indicating that these cases produce self-similar solutions.

2086

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087

The results for the low-amplitude loading Case I are clearly distinct. The delayed onset of crack growth and the lower initial crack-tip velocity of Case I make clear that there is no combination of scalings under which the solution to this problem is self-similar to the others. Fig. 5 presents crack-tip trajectories in nondimensional form, where the nondimensional crack-tip position is given by a0 :¼ a=e L, and Fig. 6 presents nondimensional results for normalized CPZ size, defined by K0 :¼ K/K0, as a function of normalized crack-tip velocity, a_ 0 .9 Both figures confirm that the solutions for Cases II and III are self-similar, even though they use distinct values ~ H , because they share the same value for the fundamental for r  0 . The solution for Case I is fundamentally different, with parameter, r lower values of K0 at lower crack-tip velocities, despite sharing the ~ H with Case II. same value of r Fig. 6 also compares our transient numerical results to the estimate for dynamic CZP size reported by Freund (1990) and Yu and Suo (2000), based on caveats of steady-state propagation and the LEFM small-scale-yielding (SSY) assumption. As expected, the estimate delivers K0 equal to unity for the quasi-static condition, a_ 0 ¼ 0. The estimate also predicts that the CPZ size vanishes as a_ 0 ! 1. Although the estimate was derived under the assumption of steady-state propagation, it is expected to remain valid for unsteady conditions when the SSY caveat is satisfied and the normalized crack speed, a_ 0 , does not change appreciably during crack extensions on the order of the CPZ-size (Freund, 1990). The results in Fig. 6 support this expectation. The SSY assumption is fairly well satisfied in all three cases considered here, but we expect increas 0 , increases.10 ing violation of the SSY caveat as the load intensity, r Thus, the high-amplitude loading in Cases II and III produces larger discrepancies with respect to the steady estimate than the lowamplitude loading in Case I. The time required for the crack tip to traverse a distance equal to the CPZ size decreases as a_ 0 increases, due to both the reduction in K with increasing crack velocity and the higher crack-tip velocity itself. Moreover, the crack-tip acceleration tends to zero at high velocities as a_ 0 approaches its steady asymptotic limit; cf. Fig. 4(b). These two effects combine to satisfy Freund’s second condition for validity of the steady estimate at high crack-tip velocities. Accordingly, all three transient solutions in Fig. 6 are in good agreement with the estimate as a_ 0 ! 1. 5. Conclusions We have presented a complete dimensional analysis for the general initial and boundary-value problem of linearized elastodynamics using differential forms notation. We obtained families of self-similar solutions and derived relations between the corresponding scalings of the mechanical fields. In particular, we showed that the velocity scalings between self-similar solutions implied by the spacetime coordinate scalings and by the material property scalings must be identical. In addition to certain advantages previously described by Abedi et al. (2006a, 2009), we showed that working with spacetime differential forms halves the number of distinct scalings between self-similar solutions that are required when working with tensorial representations. In particular, for any pair of self-similar solutions, we showed that all of the kinematic forms share a common scaling, as do all of the force-like forms. We extended our dimensional analysis to address cohesive models of dynamic fracture, for which we obtained complete sets of intrinsic cohesive scales and nondimensional parameters. ~ was previously introduced by Although the cohesive time scale s 9 We normalize K by the quasi-static CPZ size, K0, cf. (40), which differs from e L by only a constant factor, f, cf. (39). 10 Although this trend is expected, we present a rigorous analysis of the relation  0 and satisfaction of the SSY assumption in forthcoming publications. between r

Camacho and Ortiz (1996) in a study of dynamic spall strength, the systematic derivation of the complete set of intrinsic cohesive scales has not, to the authors’ knowledge, been previously reported. In addition, we obtained the complete set of fundamental nondimensional parameters that uniquely characterize any family of self-similar cohesive fracture solutions and verified this finding numerically. We emphasize that not all nondimensional parameters are fundamental. For example, we showed both analytically (see the penultimate paragraph of Section 3.3) and numerically ~ H , is (cf. Section 4) that the nondimensional cohesive strength, r not fundamental. The set of fundamental nondimensional parameters identified in Section 3.3.2 provides a useful framework for the design of experimental, numerical and analytical studies of cohesive fracture. For modeling purposes, and subject to the applicability of the linearized theory, it suffices to match the set of fundamental parameters used in an experimental or numerical study to those of the target elastodynamic system. We can avoid redundant parametric studies by varying only the fundamental parameters, rather than some larger set of model parameters. The intrinsic cohesive scales introduced in Section 3.3.3 also facilitate analysis of CFM. For example, we show in a forthcoming paper that the incubation time, the peak velocity and the diameter of the quasi-singular velocity fields observed by Abedi et al. (2009), scale with the ~; v ~ and e intrinsic parameters, s L. Several extensions of the present work are possible. First, parametric studies of other members of the fundamental nondimensional set advanced in Section 3.3.2 would be of interest. Our methodology for dimensional analysis can be applied to other classes of fracture models and physical processes. For example, in an upcoming paper we present a dimensional analysis for a fracture model that replaces the TSR with an interfacial damage model that incorporates a contact model for crack closure; see also Abedi (2010). Our analysis generates exactly the same set of cohesive scales and fundamental nondimensional parameters, except the ~ now appears as an explicit parameter in cohesive time scale s the delay-type interfacial damage model. The cohesive separation scale ~ d no longer appears as an explicit parameter in the fracture model; instead, it derives from Eq. (39b). Nonlinear elastodynamics, plasticity and thermoelasticity present alternative descriptions of the bulk response where the use of differential forms notation in the spacetime setting and dimensional analysis might prove beneficial. Acknowledgment Support from the US National Science Foundation via Grant DMR 01-21695 is gratefully acknowledged. References Abedi, R., 2010. Spacetime damage-based cohesive model for elastodynamic fracture with dynamic contact. Ph.D. thesis. Department of Theoretical and Applied Mechanics, University of Illinois at Urbana–Champaign. Abedi, R., Haber, R.B., Petracovici, B., 2006a. A spacetime discontinuous Galerkin method for elastodynamics with element-level balance of linear momentum. Computer Methods in Applied Mechanics and Engineering 195, 3247–3273. Abedi, R., Haber, R.B., Thite, S., Erickson, J., 2006b. An h-adaptive spacetime– discontinuous Galerkin method for linearized elastodynamics. Revue Européenne de Mécanique Numérique (European Journal of Computational Mechanics) 15, 619–642. Abedi, R., Hawker, M.A., Haber, R.B., Matouš, K., 2009. An adaptive spacetime discontinuous Galerkin method for cohesive models of elastodynamic fracture. International Journal for Numerical Methods in Engineering 1, 1–42. Arnold, V.I., 1989. Mathematical Methods of Classical Mechanics, second ed. Springer, New York. Barenblatt, G.I., 1962. The mathematical theory of equilibrium of cracks in brittle fracture. Advanced Applied Mechanics 7, 55–129. Birkhoff, G., 1948. Dimensional analysis of partial differential equations. Electrical Engineering 67, 1185–1188.

R. Abedi, R.B. Haber / International Journal of Solids and Structures 48 (2011) 2076–2087 Buckingham, E., 1914. On physically similar systems; illustrations of the use of dimensional equations. Physics Review 4, 345–376. Camacho, G.T., Ortiz, M., 1996. Computational modelling of impact damage in brittle materials. International Journal of Solids and Structures 33, 2899–2938. Carpinteri, A., 1982. Notch sensitivity in fracture testing of aggregative materials. Engineering Fracture Mechanics 16, 467–481. Carpinteri, A., 1989. Cusp catastrophe interpretation of fracture instability. Journal of the Mechanics and Physics of Solids 37, 567–582. Carpinteri, A., 1991. Size-scale transition from ductile to brittle failure. structural response vs. crack growth resistance curve. International Journal of Fracture 51, 175–186. Carpinteri, A., Cornetti, P., Barpi, F., Valente, S., 2003. Cohesive crack model description of ductile to brittle size-scale transition: dimensional analysis vs. renormalization group theory. Engineering Fracture Mechanics 70, 1809–1839. Dugdale, D.S., 1960. Yielding of steel sheets containing slits. Journal of the Mechanics and Physics of Solids 8, 100–104. Fleming, W.H., 1964. Functions of Several Variables. Addison-Wesley, Reading, Massachusetts. Freund, L.B., 1990. Dynamic Fracture Mechanics. Cambridge University Press, Cambridge, UK. Harder, N., 1991. Dimensional analysis for concrete in fracture. Materiaux et Constructions 24, 202–209. Huntley, H.E., 1967. Dimensional Analysis. Dover Publications, New York. Isaacson, E., Isaacson, M., 1975. Dimensional Methods in Engineering and Physics. John Wiley & Sons, New York. Kysar, J.W., 2003. Energy dissipation mechanisms in ductile fracture. Journal of the Mechanics and Physics of Solids 51, 795–824. Langhaar, H.L., 1951. Dimensional Analysis and Theory of Models. John Wiley & Sons, New York. Miles, J., 1960. Homogeneous solutions in elastic wave propagation. Quarterly of Applied Mathematics 18, 37–59. Miller, S.T., Kraczek, B., Haber, R.B., Johnson, D.D., 2009. Multi-field spacetime discontinuous Galerkin methods for linearized elastodynamics. Computer Methods in Applied Mechanics and Engineering 199, 34–47. Morgan, A.J.A., 1952. The reduction by one of the number of independent variables in some systems of partial differential equations. Quarterly Journal of Mathematics 3, 250–259. Norwood, F.R., 1973. Similarity solutions in plane elastodynamics. International Journal of Solids and Structures 9, 789–803.

2087

Pandolfi, A., Krysl, P., Ortiz, M., 1999. Finite element simulations of ring expansion and fragmentation: the capturing of length and time scales through cohesive models of fracture. International Journal of Fracture 95, 279–297. Rahulkumar, P., Jagota, A., Bennison, S., Saigal, S., 2000. Cohesive element modeling of viscoelastic fracture: application to peel testing of polymers. International Journal of Solids and Structures 37, 1873–1897. Rice, J., 1968. Mathematical analysis in the mechanics of fracture. In: Liebowitz, H. (Ed.), Fracture, An Advanced Treatise, vol. 2. Academic Press, New York, United States, pp. 191–311. Rice, J., 1980. The mechanics of earthquake rupture. In: Physics of the Earths Interior. North-Holland Publ. Co., Amsterdam. Italian Physical Society, pp. 555– 649. Sedov, L.I., 1959. Similarity and Dimensional Methods in mechanics. Academic Press, New York. Setien, J., Varona, J., 1996. On the use of dimensional analysis in fracture mechanics. In: Proceedings of the 11th Biennial European Conference on Fracture – ECF, vol. 11, pp. 125–132. Spivak, M., 1965. Calculus on Manifolds. W.A. Benjamin, New York. Szata, M., 2001. Dimensional analysis as a tool for modeling fracture-mechanics processes. Materials Science 37, 442–446. Szirtes, T., 1998. Applied Dimensional Analysis and Modeling. McGraw Hill, New York. Tvergaard, V., Hutchinson, J.W., 1992. The relation between crack growth resistance and fracture process parameters in elastic–plastic solids. Journal of the Mechanics and Physics of Solids 40, 1377–1397. Wagner, H., 1984. On the use of the dimensional analysis in material failure. International Journal of Fracture 25, R83–R86. Xu, X.P., Needleman, A., 1993. Void nucleation by inclusion debonding in a crystal matrix. Modelling and Simulation in Materials Science and Engineering 1, 111– 132. Xu, X.P., Needleman, A., 1994. Numerical simulations of fast crack growth in brittle solids. Journal of the Mechanics and Physics of Solids 42, 1397–1434. Yu, H.H., Suo, Z., 2000. Intersonic crack growth on an interface. Proceedings of the Royal Society of London, Series A (Mathematical, Physical and Engineering Sciences) 456, 223–246. Zhang, Z., Paulino, G.H., Celes, W., 2007. Extrinsic cohesive modelling of dynamic fracture and microbranching instability in brittle materials. International Journal for Numerical Methods in Engineering 72, 893–923.