# On fracture criteria for dynamic crack propagation in elastic materials with couple stresses

## Abstract

The focus of the article is on fracture criteria for dynamic crack propagation in elastic materials with microstructures. Steady-state propagation of a Mode III
semi-infinite crack subject to loading applied on the crack surfaces is considered. The micropolar behavior of the material is described by the theory of couple-stress
elasticity developed by Koiter. This constitutive model includes the characteristic lengths in bending and torsion, and thus it is able to account for the underlying
microstructures of the material. Both translational and micro-rotational inertial terms are included in the balance equations, and the behavior of the
solution near to the crack tip is investigated by means of an asymptotic analysis. The asymptotic fields are used to evaluate the dynamic
J-integral for a couple-stress material, and the energy release rate is derived by the corresponding conservation law. The propagation stability is studied according to the
energy-based Griffith criterion and the obtained results are compared to those derived by the application of the maximum total shear stress criterion.

*Keywords:* Couple-stress elasticity, Dynamic fracture, Steady-state propagation, Energy release rate, Fracture criterion.

On fracture criteria for dynamic crack propagation

in elastic materials with couple stresses

L. Morini, A. Piccolroaz, G. Mishuris and E. Radi

*Dipartimento di Scienze e Metodi dell’Ingegneria, Universitá di Modena e Reggio Emilia,*

*Via Amendola 2, 42100, Reggio Emilia, Italy.*

*Department of Civil, Environmental and Mechanical Engineering, University of Trento,*

*Via Mesiano 77, 38123, Trento, Italy.*

*Institute of Mathematical and Physical Sciences, Aberystwyth University,*

*Ceredigion SY23 3BZ, Wales, U.K.*

Corresponding author. Tel.: +390461282583, email address: lorenzo.morini@unitn.it.

## 1 Introduction

In many experimental analyses it has been shown that the mechanical behavior of brittle materials such as ceramics, composites, cellular materials, foams, masonry, bones tissues, glassy and semicrystalline polymers, is strongly affected by the the microstructure.

The influence of materials inhomogeneities and defects on the mechanical properties of the materials and their interactions with cracks have been extensively studied in the framework of classical theory of elasticity by means of homogenization theories (Hashin, 1959; Eshelby, 1976; Budiansky, 1965; Budiansky and O’ Connell, 1976), and the effective elastic moduli of bodies containing several ensembles of microstructures have been determined using for instance self-consistent methods (Kachanov, 1987, 1992; Huang et al., 1994). Furthermore, modern multi-scale simulation approaches have been developed for modelling microstructural properties of the materials (Askes et al., 2008, 2009; Silberschmidt, 2009; Andrade et al., 2011).

In Bigoni and Drugan (2007) it has been shown that classical homogenization results describe accurately elastic properties of heterogeneous materials in situations where the bodies are subjected to slowly-varying loading and the displacement and stress gradients are small. If high gradients are present, standard homogenized materials cannot represent the physical response of composite elastic media. Moreover, approaches based on classical elasticity theory cannot always predict enough accurately the size effect experimentally observed when the representative scale of the deformation field becomes comparable to the length scale of the microstructure (Lakes, 1986). For example, it has been detected that in presence of stress concentration, such as near inclusions and holes, the strength of the material is higher if the grain size is smaller, and that the bending and torsional strength of beams and wires are greater if their cross-section is thinner (Fleck et al., 1994). Since stress concentration factors derived applying standard homogenization procedures to Cauchy elastic materials depend only on the shape of the inhomogeneities and not on its size, they cannot describe accurately these effects. On the other hand, the utilization of multi-scale techniques for studying microstructural properties of the materials implies challenging numerical computations, and the validation of the results requires a critical comparison with analytical solutions and experimental data.

Generalized theories of continuum mechanics, such as micropolar elasticity (Cosserat and Cosserat, 1909), indeterminate couple stress elasticity (Koiter, 1964) and strain gradient theories (Fleck and Hutchinson, 2001; Aifantis, 2011; Dal Corso and Willis, 2011), may be considerate as an effort to include rigorously defined characteristic length scales and to study influence of microstructures on the material behavior avoiding strong numerical calculations required by multi-scale approaches. Indeed, exact analytical formulas for characteristic lengths in couple stress elastic materials have been derived via higher order homogenization of heterogeneous Cauchy elastic materials by Bigoni and Drugan (2007) and via numerical computations by Askes and Aifantis (2011). Analytical solutions derived for couple stress and strain gradient solids can also be used for validating results of numerical simulations and analyzing experimental data obtained for materials with microstructures (Lakes, 1995).

Indeterminate couple stress elasticity theory developed by Koiter (1964) provides two distinct characteristic length scales for bending and torsion. Moreover, it includes the effects of the microrotational inertia, which can be considered as an additional dynamic length scale. Therefore, in order to study crack propagation stability in couple-stress elastic materials, new fracture criteria accounting for both effects of scale lengths and microrotational inertia must be formulated (Morozov, 1984). For antiplane crack problems, in Georgiadis (2003) and Radi (2008) a critical level for the maximum total shear stress ahead of the crack tip at which the crack starts propagating has been proposed as fracture criterion. This is known as the maximum total shear stress criterion, and later in the article we will refer to it as the criterion.

The J-integral for static crack problems in couple stress elasticity has been derived by Lubarda and Markenscoff (2000) in the case of an homogeneous material, and by Piccolroaz et al. (2012) for an interfacial crack. The energy release rate can be evaluated by means of the conservation of this integral, and further in the paper we will refer to it as ERR. Nevertheless, energy-based dynamic fracture criteria for this kind of materials are still unknown in literature. For that reason, the principal aim of this paper is to generalize the static energy release rate expression to the case of dynamic steady state crack propagation, and to study the effects of the microstructure on the propagation stability by applying the energy Griffith criterion (Freund, 1998). The results are compared with those obtained in Mishuris et al. (2013), where the criterion has been adopted.

The structure of the paper is organized as follows: in Section 2 the problem of a semi-infinite Mode III steady state propagating crack in couple stress elastic materials is formulated. Both translational and micro-rotational inertial terms are included in the balance equations, and a distributed loading applied on the crack surfaces is assumed. In Section 3, the dynamic conservation laws derived by Freund and Hutchinson (1985) and Freund (1998) for classical elasticity are generalized to couple stress materials. General expressions for the dynamic J-integral and energy release rate associated to steadily propagating cracks in couple stress elastic solids are derived. Explicit forms are obtained for the case of a Mode III crack, and the path-independence of the J-integral is demonstrated in Appendix A.

An asymptotic analysis of the stress and displacement fields near to the crack tip is performed in Section 4. The contribution of the asymptotic terms to the dynamic energy release rate is analyzed in details, the leading term corresponding to finite non-zero energy is individuated and an explicit formula for the J-integral evaluated along a circular path surrounding the crack tip is derived. The obtained formula involves a constant term depending on the boundary conditions of the problem and indicating the amplitude of the leading order term of the asymptotic shear stress. This term is evaluated in closed form in Section 5 by performing an asymptotic expansion of the full-field solution derived in Mishuris et al. (2013) for the same loading configuration applied at crack faces. In this Section, the asymptotics expansion of the full-field solution is also used for deriving an alternative equivalent formula for the dynamic J-integral, calculated considering the square-shaped path around the crack tip introduced by Freund (1998) and used in Georgiadis (2003); Gourgiotis et al. (2011); Aravas and Giannakopoulos (2009). The energy release rate associated to a steady propagating Mode III crack in couple stress elastic solids is compared to the corresponding expression in classical elastic materials.

In Section 6, the obtained expression for the energy release rate is used for studying subsonic crack propagation stability. Assuming the energy-based Griffith criterion (Willis, 1971, 1967; Obrezanova et al., 2002), under the considered loading conditions the steady state propagation turns out to be unstable regardless of the values the microstructural parameters. This result appears to be in contrast with those detected in Mishuris et al. (2013) adopting the criterion, which instead shows a stabilizing effect in presence of relevant microstructures contribution. In the authors’ opinion, this discrepancy may be due to the fact that the energy release rate depends only on the leading term of the asymptotic expansion of the stresses, which dominates very close to the crack tip but provides unphysical features such as negative total shear stress ahead of the crack tip. Therefore, at a characteristic distance from the crack tip the sole leading order term may not describe the correct behavior of stresses and displacements. Differently, the total shear stress involved in the criterion is calculated by means of the full-field solution, that takes fully into account the microstructural contributions. As a consequence, in order to study the crack propagation stability in elastic materials with microstructure, fracture criteria including also higher order terms of the asymptotic stresses and involving two or more characteristic parameters should be used. Note that in classical elasticity fracture criteria including higher order terms contributions such as T-stress criterion (Hancock and Du, 1991; Smith et al., 2006) have also been proposed.

## 2 Steady-state cracks in couple stress elastic materials

In this Section the problem of the steady-state dynamic propagation of a Mode III crack in elastic materials with microstructures is formulated by means of the fully dynamical version of the couple-stress elastic model, accounting both translational and micro-rotational inertial terms into the balance equations. Reference is made to a fixed Cartesian coordinate system centred at the crack tip at the initial time . Under antiplane shear deformation, the indeterminate theory of couple stress elasticity (Koiter, 1964) adopted in the present study provides the following kinematical compatibility conditions between the out-of-plane displacement , rotation vector , strain tensor and rotation gradient tensor :

(1) |

(2) |

Therefore, the rotations are derived from displacement. The rotation gradient tensor , also known as deformations curvature tensor or torsion-flexure tensor (Koiter, 1964), is defined in the general three-dimensional case as . The vanishing of the Saint Venant tensor (or incompatibility tensor) requires to be irrotational:

(3) |

Using expressions (2), it can be immediately verified that relation (3) is satisfied. According to the indeterminate couple stress theory the non-symmetric Cauchy stress tensor can be decomposed into a symmetric part and a skew-symmetric part , namely . In addition, the couple stress tensor is introduced as the work-conjugated quantity of . The reduced tractions vector and the couple stress tractions vector are defined as

(4) |

respectively, where denotes the outward unit normal and . The conditions of dynamic equilibrium of forces and moments, taking into consideration rotational inertia, and neglecting body forces and body couples, write

(5) |

where is the mass density and is the rotational inertia.

Within the context of small deformations theory, the total strain and the deformation curvature are connected to stress and couple stress through the following isotropic constitutive relations

(6) |

where is the elastic shear modulus, and the couple stress parameters, with . Note that for antiplane deformations . Both material parameters and depend on the microstructure and can be connected to the material characteristic length in bending and in torsion, namely

(7) |

Typical experimental values of and for some classes of materials with microstructure can be found in Lakes (1986, 1995), and analytical expressions for these moduli have been derived via homogenization of heterogeneous Cauchy elastic materials by Bigoni and Drugan (2007).The limit value of corresponds to vanishing characteristic length in torsion, which is typical of polycristalline metals. Moreover, from the definitions (7) it follows that for and for . The constitutive equations of the indeterminate couple stress theory do not define the skew-symmetric part of the total stress tensor , which instead is determined by the equilibrium equations (5). Constitutive equations (6) together with the compatibility relations (1) and (2) give stresses and couple stresses in terms of the out of plane displacement :

(8) |

(9) |

The introduction of (9) into (5) yields:

(10) |

where denotes the Laplace operator. By means of (8) and (10), the equation of motion (5) becomes

(11) |

We assume that the crack propagates with a constant velocity along the -axis. In this case it is convenient to introduce a moving framework , and the out of the plane displacement can be assumed in the form:

(12) |

It follows that the time derivative of the displacement can be written in terms of the derivative with respect to , namely and thus . Therefore the equation of motion (11) under steady-state conditions becomes:

(13) |

where is the normalized crack velocity, is the shear wave speed for classical elastic materials, the characteristic length is defined as with , and (see Mishuris et al. (2013)).

According to (4), the non-vanishing components of the reduced traction and couple stress traction vectors along the crack line can be written as

(14) |

respectively. By using (8), (9), (10) and (14), the skew-symmetry of the Mode III crack problem requires ahead of the crack tip:

(15) |

On the crack surface, vanishing of the reduced traction and couple stress traction yield to the following boundary conditions for the function :

(16) |

where is the loading applied on the crack faces, which is assumed to have the following form:

(17) |

Note that although we discuss here only a specific loading condition, the main conclusions reported in this paper have been confirmed for other types of loading.

## 3 Dynamic energy release rate

In this Section, the dynamic conservation laws obtained for linear elastic media by Freund (1998) and Freund and Hutchinson (1985) are generalized to couple stress elastic materials. An explicit integral expression for the dynamic energy release rate associated to steady state cracks propagation in elastic solids with presence of couple stress is derived.

The energy release rate for a dynamic crack has been defined by Freund (1998) as the following limit:

(18) |

where is the total energy flux through the contour surrounding the crack tip. For couple stress elastic materials, the energy flux for a dynamic crack propagating along -axis is given by:

(19) |

where is an outward unit normal on , denotes the strain-energy density

(20) |

and is the kinetic energy density

(21) |

Since a Mode III steady-state crack propagating at constant velocity along the -axis is considered, the expression (12) for the out-of-plane displacement is used and then in the moving framework the energy flux (19) assumes the special form

(22) |

then the generalized J-integral for an antiplane dynamic steady-state crack in couple-stress elastic materials can be defined:

(23) | |||||

(24) |

The expressions (23) and has been proved to be path independent, the details of the demonstration are reported in Appendix A. Moreover, the equivalence of the two alternative forms of the J-integral (23) and (24), the first written in function of the tractions and the second of the reduced tractions, is demonstrated in Appendix B. The (23) is the generalization of the static expressions derived by Atkinson and Leppington (1974, 1977), and Lubarda and Markenscoff (2000) to the antiplane dynamic steady state case. The dynamic energy release rate is then defined by the limit:

(25) |

Using definitions (1) and (2), the strain energy density (20) becomes

(26) |

whereas for steady state propagation the kinetic energy density is given by

(27) |

A polar coordinates system centered at the crack tip is assumed, and a circular contour of radius around the crack tip with is considered. Then, substituting expressions (26) and (27) the J-integral (23) becomes

(28) | |||||

According to the definition (25), the energy release rate can be evaluated as the limit for of the integral (28). In order to evaluate explicitly the J-integral (28) and to investigate the variation of the energy release rate (25) in function of the crack propagation velocity and its implications on propagation stability, the behavior of the out-of-plane displacement near to the crack tip is studied by means of an asymptotic analysis in the next Section.

## 4 Asymptotic crack tip fields

The following standard asymptotic expression for out-of-plane displacement in separate variables form is considered:

(29) |

We are interested in finding the terms of the asymptotic solution (29) corresponding to finite and non-zero contributions to the J-integral (28) in the limit . Since the displacement should be bounded and symmetrical, it follows immediately that , and then no zero-order terms are present in the asymptotic expansion. Substituting the expression (29) in the generalized J-integral formula (28) and using the following derivative rules which hold for an arbitrary function :

(30) |

we get:

(31) | |||||

where the superscript denotes the total derivative respect to the variable .

Observing expression (31), in agreement with the results reported in Radi (2008) for the static case, we deduce that the finiteness of the energy release rate towards the crack tip requires that for any non integer number.

Including higher order terms in the asymptotic expansion in the form and using this expression in (28), terms of order and where are detected. These terms involve both and and can have a non vanishing impact to the value of the energy release rate. Therefore we need to consider several asymptotic terms in the form of (29) and analyze the possible correlation between them in the nonlinear functional (28). According with this discussion and in order to find the terms corresponding to finite and non-zero contributions to the energy release rate, we assume .

It can be easily demonstrated that if , only the leading order term of the governing equation (13) can be considered, while if more terms are required with exponents differing by 2 or more than 2, the full equation (13) must be considered in the analysis (Piccolroaz et al., 2012). As a consequence, assuming , we can then keep only the leading term of the evolution equation (13)

(32) |

Introducing the expression (29) in (32), the general asymptotic solution of the equation of motion has been obtained, the derivation is illustrated in details in the Section containing the supplementary material. Referring to this general solution, since we are assuming values in the range , the terms corresponding to needs to be considered and the asymptotic expression for the out of plane displacement turns out to be:

(33) | |||||

where

(34) |

and the constants and define the amplitude of each asymptotic term the sum and thus must to be evaluated according to the boundary conditions (15) and (16).

It has been verified by symbolic calculations that the terms of the (33) corresponding to and , similarly to what has been detected by Radi (2008) and Piccolroaz et al. (2012) for a stationary crack case do not contribute to the J-integral and to energy release rate. Although these terms are relevant for evaluating displacement and total shear stress at the crack tip, the only finite and non-vanishing contribution to the generalized J-integral (28) is associated to the order of the asymptotic expression (33), and then the energy release rate is given by:

(35) | |||||

where

(36) | |||||

The proposed procedure, based on the assumption and then on the analysis of the leading order term of the governing equation (32), can in principle be performed considering a different range of values of , as , corresponding to terms , but observing the expression (31) and remembering that the other contributions to the J-integral are of orders and and that is excluded because finite energy is required, it is easy to deduce that all terms associated to do not contribute to the energy release rate. Therefore we can conclude that provides effectively the only term contributing to the J-integral and that the choice of considering the leading term of the governing equation (13) in the asymptotic analysis is correct for our purpose.

In order to evaluate the energy release rate (35), the constant must be explicitly determined. In the next Section is calculated starting from the asymptotic expansion of the full-field solution derived in Mishuris et al. (2013) for the same loading conditions (15)-(17) by means of Wiener-Hopf technique.

## 5 Explicit evaluation of the energy release rate

In order to derive an explicit expression for the constant , we perform the asymptotic analysis of the Fourier transform of the full-field solution derived in Mishuris et al. (2013) for the same loading conditions (15)-(17) by means of Wiener-Hopf technique. In the limit , the Fourier transforms of stress, couple stress field and displacements assume the following behavior:

(37) | |||||

(38) | |||||

(39) |

where is a constant determined starting from the the boundary condition (15) and applying the Liouville theorem and:

(40) |

the explicit expression for the factorization function is given by Mishuris et al. (2013). Substituting the (40) into (37) and (39) we obtain:

(41) | |||||

(42) | |||||

(43) |

Further, we consider the following transformation formula (Roos, 1969):

(44) |

where is the gamma function and the symbol indicates that the quantities on the two sides of the (44) are connected by means of unilateral Fourier transform. Applying the (44) to expressions (41) and (43), we get:

(45) | |||||

(46) | |||||

(47) |

On the crack surface, for and , the term of order of the asymptotic expression in polar coordinates (33) becomes

(48) |

where the definition has been used. Equating expressions (47) and (48), we get:

(49) |

The explicit expression (49) can the be used into the (35) for studying the variation of the energy release rate in function of the crack tip speed and of the microstructural parameters and .

In order to check the validity of the obtained results, an alternative expression for the energy release rate is derived considering the rectangular-shaped contour with vanishing height along the y-direction and with reported in Fig.1 (Freund, 1998). The Cartesian components of the outward unit vector normal to are . Considering the moving framework in Fig.1 with the origin at the crack tip, for the steady-state antiplane crack problem the generalized J-integral (23) becomes:

(50) | |||||

where and are components of the total non-symmetric Cauchy stress tensor, including both symmetric and skew-symmetric part.

In order to evaluate the energy release rate, we allow the height of the rectangular path reported in Fig.1 to vanish. In this limit, the first integral of the (50) is zero and the following expression for the dynamic energy release rate is derived:

(51) |

It is important to note that anti-symmetry conditions (15) together with boundary conditions (16) provide that the reduced traction is zero along the whole crack line , where . Consequently, the energy release rate (51) becomes:

(52) | |||||

For evaluating the integral (52), solely asymptotic expressions for the traction ahead of the crack tip , the couple stress field , and the displacement are required. Then, by substituting equations (45), (46) and (47) in the general formula (52), we finally derive:

(53) | |||||

where and are distributions of the bisection type (Arfken and Weber, 2005). For any real with the exception of , this particular type of distributions is defined as follows:

Moreover, the product of distributions inside the integral in (53) is evaluated through the application of Fisher’s theorem (Fischer, 1971), that leads to the operational relation:

(54) |

where is the Dirac delta distribution. Then, by using the relation (54) into (53) and considering the fundamental property of the Dirac delta distribution (Arfken and Weber, 2005), we finally get:

(55) |

Using this alternative procedure, we have derived the explicit expression (55) for the energy release rate corresponding to a Mode III steady state propagating crack in a couple stress elastic material. Equation (55) can be compared with the energy release rate associate to an antiplane steady state crack in a classical elastic material, derived assuming the same loading configuration (17):

(56) |

the ratio between the two expressions (55) and (56) is given by:

(57) |

## 6 Results and discussion

In this Section, the variation of the energy release rate is analyzed applying the classical Griffith criterion (Willis, 1967) in order to study the propagation stability. The results are compared with those detected using the criterion, adopted by Mishuris et al. (2013).

The normalized variation of the energy release rate versus the crack tip speed is reported in Fig. 2 for the same value of
the ratio , three different values of and four different values of the rotational inertia . The range of
the normalized crack tip velocity has been chosen in such a way that the propagation is subsonic and the conditions of validity of the full-field solution obtained in Mishuris et al. (2013)
and used for evaluating the J-integral constant are satisfied. A similar behavior is observed for all different set of parameters: the energy release rate is initially constant for small values of the crack tip speed
, then increases monotonically until its limiting value. For small values of and , the limit value of the energy release rate corresponds to ,
while as the microstructural parameters increases the limit value of the crack tip speed becomes smaller than the shear waves speed , and thus .

On the basis of the Griffith criterion (Freund, 1998), crack initiation requires that the energy release rate achieves a critical value , where is the energy needed to form a unit of new material surface, and is a constant depending only by the properties of the medium. In our case, once this critical value for the crack initiation is achieved, the energy release rate always increases in function of the velocity, this means that if the applied loading provides the energy necessary for starting the fracture process, the crack has enough energy to accelerate rapidly up to the limiting values of the speed (Willis, 1971; Obrezanova et al., 2002). This implies that, if we analyze the results shown in Fig. 2 by means of Griffith criterion, the steady state propagation of a Mode III crack in couple-stress elastic material can be considered unstable for any value of the rotational inertia and of .

The ratio between the energy release rate in couple stress materials and the energy release rate in classical elastic materials (57) is plotted in Fig. 3 as a function of the normalized crack tip speed . For small values of the rotational inertia, this ratio decreases as tends to . This is due to the fact that for , the energy release rate corresponding to a classical elastic material (56) diverges, while in presence of couple stress it has a finite limiting value. Differently, for large values of , the ratio increases monotonically until a limiting value corresponding to and thus to a crack tip speed smaller than . Therefore, in order to propagate the crack at constant velocity, more energy than in a classical elastic material must be provided by the loading if the contribution of the rotational inertia is relevant. Observing Fig. 3, we can also note that in the static limit, namely as tends to zero, the ratio is close to one for , this behavior is in agreement with the results reported by Radi (2008), which illustrate that for the solution approaches the classical elastic case.

In Mishuris et al. (2013) the criterion (Georgiadis, 2003) is applied to the same crack problem. This alternative criterion states that the maximum total shear stress possesses a critical level at which the crack starts propagating. The behavior of as a function of the crack tip speed has been studied extensively for several sets of microstructural parameters. For each values of the ratio , it has been individuated a critical value of the rotational inertia , such that for the maximum shear stress increases very rapidly in function of