Employment of Fracture Mechanics Criteria for Accurate Assessment of the Full Set of Elastic Constants of Orthorhombic/Tetragonal Mono-Crystalline YBCO

1. Introduction

Elastic constants of engineering materials are crucial for understanding the deformation and failure behaviors of structural components from both macroscopic and microscopic points of view. From microscopic perspective, their importance arises from their intimate relationship to such solid state phenomena as specific heat, Debye temperature and Grunelsen parameter [1].

Discovery of superconductivity in La-Ba-Cu-O system in the1980’s by Bednorz and Muller [2] has spurred an enormous amount of activities in search for high (i.e., above the boiling point of liquid nitrogen) T_c superconductors (HTS) such as mono-crystalline YBCO (yttrium barium copper oxide), which is also called 1-2-3 superconductor (because of the presence of one yttrium, two barium and three copper atoms in its unit cell). Modern applications include Josephson junctions, which can act as a switch for magnetic fields, or alternatively, perform the function of a magnetic detector, called the superconducting quantum interference device (SQUID) [3]. Debye temperature, Θ_D, of a superconductor can be determined from the knowledge of the elastic constants, c_ij, in a manner described by Eqs. (2) and (3) of Lei et al. [1], which, in combination with the electron-phonon coupling parameter, λ^∗, can be used to compute the superconducting transition temperature, T_c [4,5].

Practical applications of such mono-crystalline superconductors are, however, limited (at cryogenic temperature) by their poor fracture toughness [6,7,8,9,10,11,12,13,14]. Asymptotic behavior of two-dimensional stress fields at the tips of cracks and wedges, has been studied extensively in the literature [15,16,17,18,19,20]; see Nejati et al [21] for an extensive literature survey. Significant progress in research on a two-dimensional cracked anisotropic solid notwithstanding, the corresponding progress in its three-dimensional counterpart appears to, till recently [22], be marked by its complete absence. More important, the above studies employ the Lekhnitskii [15] and Stroh [16] type formulations, which are all based on complex variables-based methodology. Since the three-dimensional space is too small to accommodate the next higher dimensional analog of complex variables (for which at least a four-dimensional space will be required; e.g., quaternion [23]), these complex variables-based analyses are not by themselves adequate for analysis of three-dimensional cracked anisotropic solids.

The mathematical difficulties posed by the three-dimensional crack/anticrack type problems are substantially greater than their two-dimensional counterparts (to start with, the governing PDE’s are much more complicated). In the absence of the knowledge of the strength of singularity, in regions where the elastic stresses become unbounded, the majority of weighted residual type methods, e.g., the finite elements, finite difference and boundary elements, which are generally employed to solve fracture mechanics problems, encounter overwhelming numerical difficulties, such as lack of convergence, and oscillation resulting in poor accuracy [24]. There are several classes of problems pertaining to the issue of three-dimensional stress singularity [25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49]: (i) through crack/anticrack as well as their bi- and tri-material interface counterparts, (ii) bi-material free edge, (iii) tri-material junction, (iv) penny shaped crack/anticrack and their bi-material counterparts, (iv) bi-material hole (and inclusion), (v) matrix cracking and fiber breaks in composites among others. Only the penny shaped crack/anticrack [50] (and their bi-material counterparts [51]) and the bimaterial hole [52,53] and inclusion problems [54] had earlier been adequately addressed in the literature. Earlier attempts to solve the three-dimensional through crack problem resulted in controversies that lasted for about a quarter century [25,55]. A unified three-dimensional eigenfunction approach has recently been developed by Chaudhuri and co-workers [25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,55] to address the three-dimensional stress singularity problems covering all the aspects mentioned above.

The above separation of variables approach has recently been extended to the case of cracked/anticracked transversely isotropic (smeared-out composite) [56] as well as cubic/orthorhombic/diamond cubic mono-crystalline plates under mode I/II far-field loadings [22,57,58] and cubic/orthorhombic/monoclinic/diamond cubic mono/tri-crystalline plates subjected to mode III loading [22,58,59,60,61] by making use of the separation of the thickness variable technique in combination with an affine transformation, that is similar (but not identical) in spirit to that due to Eshelby et al. [62] and Stroh [16]. This eigenfunction expansion approach has also recently employed to obtain three-dimensional asymptotic stress fields in the vicinity of the front of the kinked carbon fiber-matrix junction [63] (see also Ref. [64] for its 2D counterpart with isotropic glass fibers).

The above-mentioned importance notwithstanding, relatively fewer attempts at experimental determination of elastic constants of monocrystalline YBCO have been reported in the solid state physics literature [1,65,66,67,68,69,70,71,72,73]; these are summarized in Table 1 of Lei et al. [1]. Golding et al. [66] and Saint-Paul and co-workers [70,71] have reported experimental results on c₁₁ and c₃₃, and c₃₃, c₄₄, c₆₆ and c₁₂, respectively, by employing the ultrasound technique, while Baumgart et al. [68,69] and Zouboulis et al. [72] have resorted to Brillouin spectroscopy/scattering to determine c₁₁, c₃₃ and c₄₄, and c₄₄, c₅₅ and c₆₆, respectively. Only three investigations [1,67,73] report complete sets of elastic constants, accessible to the present author. Worse still, those reported by Ledbetter and Lei [73] are just estimates (marked ** in Table 1(a)), while their experimental counterparts due to Reichard et al. [67] are based on the assumption of tetragonal symmetry; see Table 1(b). This only leaves the experimental measurements (by resonant ultrasound spectroscopy, described in detail by Migliori et al. [74,75]) due to Lei et al. [1], marked * in Table 1. However, their $c_{12}$ value appears to be excessively high. This is because, according to these authors themselves, "no wave speed in the crystal depends only upon $c_{12}$, it is no way to estimate it directly." It also has been known for some times that while $c_{12}$ and $c_{66}$ can be measured independently by static tests [76], these constants are always coupled in vibrations-based measurements [77,78].

The above literature review reveals an absence of reliable and accurate experimentally measured complete sets of nine elastic constants needed for characterization of the deformation/fracture as well as other solid state (e.g., Debye temperature, T_c, etc.) behaviors of superconducting (orthorhombic) YBCO single crystals. This calls for a reliable criterion for assessment of the measured data that would allow us to come up with a reasonably accurate complete set of nine elastic constants, which is the primary objective of the present investigation. It is proposed here that one way to address this important issue is to analytically examine the effect of elastic constants on crack deflection in monocrystalline YBCO (YBa₂Cu₃O_7−δ) and compare with the experimental results for easy cleavage planes, reported by Cook et al. [6], Raynes [9], and Granozio and di Uccio [14] among others. In what follows, the above-mentioned modified eigenfunction expansion technique, based in part on separation of the thickness-variable and partly on the Eshelby-Stroh type affine transformation, is developed to derive three-dimensional asymptotic stress field in the vicinity of the front of a semi-infinite through-thickness crack weakening an infinite orthorhombic mono-crystalline plate, of finite thickness and subjected to far-field mode I/II loadings. Crack-face boundary conditions and those that are prescribed on the top and bottom (free or fixed) surfaces of the plate are exactly satisfied. The present investigation considers six through-crack systems ― (010)[001] with the [100] length direction, (0$\overline{1}$0)[100] with the [001] length direction, ($\overline{1}$00)[001] with the [010] length direction, (100)[010] with the [001] length direction, (001)[0$\overline{1}$0] with the [100] length direction, (001)[100] with the [010] length direction ― weakening orthorhombic mono-crystalline plates. Explicit expressions for the singular stresses in the vicinity of the front of a through-thickness crack weakening an orthorhombic mono-crystalline plate, subjected to far-field mode I/II loadings, are presented. In addition, through-thickness distribution of the stress intensity factors and energy release rates are also presented. Next, the important issue of easy or difficult cleavage plane and the related question of crack deflection criterion is discussed. The latter is based on the relative fracture energy (or the energy release rate) available for possible “fracture paths” [18]. This said, it is noteworthy that any fracture criterion derived from equilibrium theories such as the Griffith criterion can only be regarded as a necessary condition but not as sufficient [79]. This calls for establishment of a sufficient condition for determination of easy or difficult cleavage plane, and the associated question of the crack deflection criterion, which is the second objective of the present study.

The present study, although to a smaller extent a review on earlier work on this topic, is largely based on original research on this subject. The topic, which covers mathematics (e.g., asymptotic, solution to 3D mixed boundary-value problem, and necessary and sufficient condition for fracture), solid state physics (e.g., superconductivity, single crystal cleavage) and engineering (e.g., 3D fracture mechanics), has, so far, remained largely unexplored in the literature.

2. Formulation of the Problem

The Cartesian coordinate system (x, y, z) is convenient to describe the deformation behavior in the vicinity of a semi-infinite through-thickness crack, weakening an infinite orthotropic plate of thickness, 2h (Figure 1 and Figure 2). Here, the z-axis is placed along the straight crack front, while the coordinates x, y, are used to define the directions along the length of the crack and transverse to it, respectively, in the plane of the plate. u, v and w represent the components of the displacements in the x, y and z directions, respectively.

The stress-strain relations for an orthorhombic single crystal are given as follows:

$$\left\{\begin{array}{c}{\sigma }_x\\ {\sigma }_y\\ {\sigma }_z\\ {\tau }_{yz}\\ {\tau }_{xz}\\ {\tau }_{xy}\end{array}\right\}=\left[\begin{array}{cccccc}{c}_{11}& c_{12}& c_{13}& 0& 0& 0\\ c_{12}& c_{22}& c_{23}& 0& 0& 0\\ c_{13}& c_{23}& c_{33}& 0& 0& 0\\ 0& 0& 0& c_{44}& 0& 0\\ 0& 0& 0& 0& c_{55}& 0\\ 0& 0& 0& 0& 0& c_{66}\end{array}\right]\left\{\begin{array}{c}{\epsilon }_x\\ {\epsilon }_y\\ {\epsilon }_z\\ {\gamma }_{yz}\\ {\gamma }_{xz}\\ {\gamma }_{xy}\end{array}\right\}.$$

(1a)

where c_ij, i, j = 1,...,6, denotes the elastic stiffness constants of an orthorhombic mono-crystalline plate. σ_x, σ_y, σ_z represent the normal stresses, and τ_xy, τ_xz, τ_yz denote the shear stresses, while ε_x, ε_y, ε_z denote normal strains, and γ_xy, γ_xz, γ_yz represent the shear strains. For the special case of a tetragonal single crystal,

$$c_{11}=c_{22},\hspace{1em}{c}_{13}=c_{23},\hspace{1em}{c}_{44}=c_{55}.$$

(1b)

The three equilibrium equations for a linear elastic solid, made of an orthotropic/ orthorhombic material, can be expressed in terms of the displacement components, u, v, and w, as follows:

$$c_{11}\frac{{\partial }^{2}u}{\partial x^2}+c_{66}\frac{{\partial }^{2}u}{\partial y^2}+c_{55}\frac{{\partial }^{2}u}{\partial z^2}+(c_{12}+c_{66})\frac{{\partial }^{2}v}{\partial x\partial y}+(c_{13}+c_{55})\frac{{\partial }^{2}w}{\partial x\partial z}=0,$$

(2a)

$$(c_{12}+c_{66})\frac{{\partial }^{2}u}{\partial x\partial y}+c_{66}\frac{{\partial }^{2}v}{\partial x^2}+c_{22}\frac{{\partial }^{2}v}{\partial y^2}+c_{44}\frac{{\partial }^{2}v}{\partial z^2}+(c_{23}+c_{44})\frac{{\partial }^{2}w}{\partial y\partial z}=0,$$

(2b)

$$(c_{13}+c_{55})\frac{{\partial }^{2}u}{\partial x\partial z}+(c_{23}+c_{44})\frac{{\partial }^{2}v}{\partial y\partial z}+c_{55}\frac{{\partial }^{2}w}{\partial x^2}+c_{44}\frac{{\partial }^{2}w}{\partial y^2}+c_{33}\frac{{\partial }^{2}w}{\partial z^2}=0,$$

(2c)

The boundary conditions include those at the plate faces and crack-side surfaces. The boundary conditions on the plate faces, z = +h, are given by [22,25]

$${\sigma }_z={\tau }_{xz}={\tau }_{yz}=0,$$

(3a)

while those at the crack-side surfaces are more conveniently expressed in local cylindrical polar coordinates (Figure 1), which are given as follows:

$${\sigma }_{\theta }={\tau }_{r\theta }={\tau }_{\theta z}=0,\hspace{1em}\theta =\pm \pi$$

(3b)

where σ_r, σ_θ, σ_z represent the normal stresses, and τ_rθ, τ_rz, τ_θz are the shear stresses, while ε_r, ε_θ, ε_z denote the normal strains, and γ_rθ, γ_rz, γ_θz are the shear strains in the cylindrical polar coordinate system (r, θ, z). u_r and u_θ represent the components of the displacement in r and θ directions, respectively.

3. Singular Stress Fields in the Vicinity of a Crack Front Weakening an Orthotropic/Orthorhombic Lamina/Single Crystal under General Loading

The assumed displacement functions for the three-dimensional crack problem under consideration are selected on the basis of separation of z-variables. These are as given below [22,56,57,58,64]:

$$u(x,y,z)=e^{ikz}U(x,y),\hspace{1em}v(x,y,z)=e^{ikz}V(x,y),\hspace{1em}w(x,y,z)=e^{ikz}W(x,y).$$

(4a, b, c)

It may be noted that since the z-dependent term and its first partial derivative can either be bounded and integrable at most admitting ordinary discontinuities, or the first partial derivative at worst be square integrable (in the sense of Labesgue integration) in its interval z ∈ [-h, h], i.e., admitting singularities weaker than square root (i.e., z ^(-1/2+ε), ε > 0), it can be best represented by Fourier series [22,25,58]. The latter case is justified by Parseval’s theorem [80], and its physical implication is that of satisfying the criterion of finiteness of local strain energy and path independence [81]. Substitution of Eqs. (4) into Eqs. (2) yields the following system of coupled partial differential equations (PDE’s):

$$c_{11}\frac{{\partial }^{2}U}{\partial x_1^2}+c_{66}\frac{{\partial }^{2}U}{\partial y_1^2}+c_{55}U+(c_{12}+c_{66})\frac{{\partial }^{2}V}{\partial x_1\partial y_1}+(c_{13}+c_{55})\frac{\partial W}{\partial x_1}=0,$$

(5a)

$$(c_{12}+c_{66})\frac{{\partial }^{2}U}{\partial x_1\partial y_1}+c_{66}\frac{{\partial }^{2}V}{\partial x_1^2}+c_{22}\frac{{\partial }^{2}V}{\partial y_1^2}+c_{44}V+(c_{23}+c_{44})\frac{\partial W}{\partial y_1}=0,$$

(5b)

$$(c_{13}+c_{55})\frac{\partial U}{\partial x_1}+(c_{23}+c_{44})\frac{\partial V}{\partial y_1}+c_{55}\frac{{\partial }^{2}W}{\partial x_1^2}+c_{44}\frac{{\partial }^{2}W}{\partial y_1^2}+c_{33}W=0,$$

(5c)

where

$$x_1=ikx,\hspace{1em}{y}_1=iky.$$

(6a,b)

The solution to the system of coupled partial differential equations (5) subjected to the most general loading, can now be sought in the form of the following modified Frobenius type series in terms of the variable x1+py1 as follows:

$$U(x_1,y_1)={\displaystyle \sum _{n=0}^{\infty }{a}_{s+n}^{\prime }}{(x_1+p{y}_1)}^{s+2n+1}+{\displaystyle \sum _{n=0}^{\infty }{a}_{s+n}}{(x_1+p{y}_1)}^{s+2n},$$

(7a)

$$V(x_1,y_1)={\displaystyle \sum _{n=0}^{\infty }{b}_{s+n}^{\prime }}{(x_1+p{y}_1)}^{s+2n+1}+{\displaystyle \sum _{n=0}^{\infty }{b}_{s+n}}{(x_1+p{y}_1)}^{s+2n},$$

(7b)

$$W(x_1,y_1)={\displaystyle \sum _{n=0}^{\infty }{c}_{s+n}^{\prime }}{(x_1+p{y}_1)}^{s+2n}+{\displaystyle \sum _{n=0}^{\infty }{c}_{s+n}}{(x_1+p{y}_1)}^{s+2n+1}.$$

(7c)

Out of the various combinations, such as (a′, b′, c′), (a, b, c), (a′, b, c), (a, b′, c), (a, b, c′), (a′, b′, c), (a′, b, c′), and (a, b′, c′), only the first two groupings can produce meaningful solutions, for the mode I/II and mode III loading cases, respectively. This step permits separation of the mode III from the modes I/II. The first grouping is described below, while the second one has already been employed for the antiplane shear case [59,60,61].

4. Singular Stress Fields in the Vicinity of a (010)[001] Through-Crack Front Propagating under Mode I (Extension/Bending) and Mode II (Sliding Shear/Twisting) in [100] Direction

The solution to the system of coupled partial differential equations (5), subjected to the far-field mode I/II loading, can now be sought in the form of the following modified Frobenius type series in terms of the variable x1+py1 as follows [22,53,54,60], although unlike in the case of isotropic materials [25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,55], the x₁ and y₁ variables are no longer separable:

$$U(x_1,y_1)={\displaystyle \sum _{n=0}^{\infty }{a}_{s+n}}{(x_1+p{y}_1)}^{s+2n},\hspace{1em}V(x_1,y_1)={\displaystyle \sum _{n=0}^{\infty }{b}_{s+n}}{(x_1+p{y}_1)}^{s+2n},\hspace{1em}W(x_1,y_1)={\displaystyle \sum _{n=0}^{\infty }{c}_{s+n}}{(x_1+p{y}_1)}^{s+2n+1}.$$

(8a, b, c)

Here, the combined variable x1+py1 represents an affine transformation in the same spirit as that by Eshelby et al. [62], Stroh [16] and Shih et al. [17], although these authors have employed completely different techniques. On substitution of Eqs. (8) into Eqs. (5), and equating the coefficients of ${(x_1+p{y}_1)}^{s+2n-2}$, the following set of recurrent relationships can be derived:

$$\left(s+2n\right)\left(s+2n-1\right)\left(c_{11}+c_{66}{p}^2\right)a_{s+n}+c_{55}{a}_{s+n-1}+$$

$$\left(s+2n\right)\left(s+2n-1\right)\left(c_{12}+c_{66}\right)p{b}_{s+n}+\left(s+2n-1\right)\left(c_{13}+c_{55}\right)c_{s+n-1}=0,$$

(9a)

$$\left(s+2n\right)\left(s+2n-1\right)\left(c_{12}+c_{66}\right)p{a}_{s+n}+\left(s+2n\right)\left(s+2n-1\right)\left(c_{22}+c_{66}{p}^2\right)b_{s+n}$$

$$+c_{44}{b}_{s+n-1}+\left(s+2n-1\right)\left(c_{23}+c_{44}\right)c_{s+n-1}=0,$$

(9b)

which for n = 0, reduces to

$$\left[\begin{array}{cc}{c}_{11}+c_{66}{p}^2& \left(c_{12}+c_{66}\right)p\\ \left(c_{12}+c_{66}\right)p& c_{22}{p}^2+c_{66}\end{array}\right]\left\{\begin{array}{c}{a}_s\\ b_s\end{array}\right\}=\left\{\begin{array}{c}0\\ 0\end{array}\right\}, \mathrm{for} \mathrm{s}\ne 0,1.$$

(10)

$$-A_1\left(\eta +H_2\right)+A_2\left(\xi +H_1\right)=0,$$

The characteristic equation for the coupled partial differential equations (2) or (5) can now be written as follows:

$$p^4+2\chi p^2+\frac{c_{11}}{c_{22}}=0,$$

(11)

in which the normalized elastic parameter, $\chi$, is given by

$$\chi =\frac{\left(c_{11}{c}_{22}-c_{12}^2-2c_{12}{c}_{66}\right)}{2c_{22}{c}_{66}}=\frac{1}{\left(1-{\nu }_{13}{\nu }_{31}\right)}\left[\frac{E_1}{2G_{12}}-\left({\nu }_{12}+{\nu }_{32}{\nu }_{13}\right)\right],$$

(12a)

in which E₁ is Young’s modulus in the x direction, G₁₂ is the shear modulus in the x-y plane, while ν₁₂ is the major Poisson’s ratio in the x-y plane. ν₁₃ and ν₃₁ denote the major and minor Poisson’s ratios in the x-z plane, while ν₃₂ represents the minor Poisson’s ratio in the y-z plane. $\chi$ can also be expressed in terms of the inverse anisotropic ratio (in the x-y plane), $\lambda$, as follows:

$$\chi =\frac{\lambda \left(\sqrt{c_{11}{c}_{22}}+c_{12}\right)-c_{12}}{c_{22}}=\sqrt{\frac{c_{11}}{c_{22}}}\lambda +\frac{c_{12}\left(\lambda -1\right)}{c_{22}}.$$

(12b)

where $\lambda$, is defined as

$$\lambda =\frac{\sqrt{c_{11}{c}_{22}}-c_{12}}{2c_{66}}.$$

(13)

Eq. (11) has either (a) four complex or (b) four imaginary roots, depending on whether: (a) λ < 1 or equivalently,

$$\chi <\sqrt{\frac{c_{11}}{c_{22}}}=\sqrt{\frac{E_1\left(1-{\nu }_{23}{\nu }_{32}\right)}{E_2\left(1-{\nu }_{13}{\nu }_{31}\right)}},$$

(14a)

or

(b) λ > 1 or equivalently,

$$\chi >\sqrt{\frac{c_{11}}{c_{22}}}=\sqrt{\frac{E_1\left(1-{\nu }_{23}{\nu }_{32}\right)}{E_2\left(1-{\nu }_{13}{\nu }_{31}\right)}}.$$

(14b)

λ= 1 or χ = 1 represents the degenerate isotropic material case, for which the solution is available in Chaudhuri and Xie [25].

4.1. Case (a): Complex Roots

$$p_{1,2}=\xi \pm i\eta ,\hspace{1em}{p}_{3,4}=-\xi \pm i\eta ,$$

(15a, b)

where

$$\xi =\frac{1}{\sqrt{2}}{\left[{\left(\frac{c_{11}}{c_{22}}\right)}^{1/2}-\chi \right]}^{1/2},\hspace{1em}\eta =\frac{1}{\sqrt{2}}{\left[{\left(\frac{c_{11}}{c_{22}}\right)}^{1/2}+\chi \right]}^{1/2},$$

(16a, b)

valid for $\chi <\sqrt{c_{11}/c_{22}}$.

The final results that satisfy the equilibrium equations (2) can be expressed in the following form:

$$u(x,y,z)=\left({\overline{\overline{D}}}_{1}i\sin \left(kz\right)+{\overline{\overline{D}}}_2\cos \left(kz\right)\right){\left(ik\right)}^s[{\overline{A}}_1{\left(x+(\xi +i\eta )y\right)}^s+{\overline{A}}_2{\left(x+(\xi -i\eta )y\right)}^s$$

$$+{\overline{A}}_3{\left(x+(-\xi -i\eta )y\right)}^s+{\overline{A}}_4{\left(x+(-\xi +i\eta )y\right)}^s],$$

(17a)

$$v(x,y,z)=\left({\overline{\overline{D}}}_{1}i\sin \left(kz\right)+{\overline{\overline{D}}}_2\cos \left(kz\right)\right){\left(ik\right)}^s[{\overline{B}}_1{\left(x+(\xi +i\eta )y\right)}^s+{\overline{B}}_2{\left(x+(\xi -i\eta )y\right)}^s$$

$$+{\overline{B}}_3{\left(x+(-\xi -i\eta )y\right)}^s+{\overline{B}}_4{\left(x+(-\xi +i\eta )y\right)}^s],$$

(17b)

$$w(x,y,z)=\left({\overline{\overline{D}}}_1\cos \left(kz\right)+{\overline{\overline{D}}}_{2}i\sin \left(kz\right)\right){\left(ik\right)}^{s+1}[{\overline{C}}_1{\left(x+(\xi +i\eta )y\right)}^{s+1}+{\overline{C}}_2{\left(x+(\xi -i\eta )y\right)}^{s+1}$$

$$+{\overline{C}}_3{\left(x+(-\xi -i\eta )y\right)}^s+{\overline{C}}_4{\left(x+(-\xi +i\eta )y\right)}^s],$$

(17c)

wherein ξ and η are as given in Eqs. (16), and ${\overline{A}}_k, {\overline{B}}_k, {\overline{C}}_k$, k = 1,...,4, are undetermined coefficients. It may be noted that ${\overline{B}}_k$ can be expressed in terms of the corresponding ${\overline{A}}_k$, k = 1,...,4 by using Eqs. (10) and (15).

$${\overline{B}}_1=\left(H_1+i{H}_2\right){\overline{A}}_1,\hspace{1em}{\overline{B}}_2=\left(H_1-i{H}_2\right){\overline{A}}_2,$$

(18a, b)

$${\overline{B}}_3=-\left(H_1+i{H}_2\right){\overline{A}}_3,\hspace{1em}{\overline{B}}_4=\left(-H_1+i{H}_2\right){\overline{A}}_4,$$

(18c, d)

in which

$$H_1=-\frac{\xi \left(\sqrt{c_{11}{c}_{22}}+c_{66}\right)}{\left(c_{12}+c_{66}\right)},\hspace{1em}{H}_2=\frac{\eta \left(\sqrt{c_{11}{c}_{22}}-c_{66}\right)}{\left(c_{12}+c_{66}\right)}.$$

(19a, b)

The corresponding stress field can easily be obtained from Eqs. (17). It is convenient to express the components of the displacement vector and stress tensor, in terms of the cylindrical polar coordinate system (r, θ, z). Expressing

$$\rho \cos (\psi )=r\left(\cos (\theta )+\xi \sin (\theta )\right),\hspace{1em}\rho \sin (\psi )=r\left(\eta \sin (\theta )\right),$$

(20a)

$${\rho }^{\prime }\cos (\psi )=r\left(\cos (\theta )-\xi \sin (\theta )\right),\hspace{1em}{\rho }^{\prime }\sin (\psi )=r\left(-\eta \sin (\theta )\right),$$

(20b)

in which

$$\rho =r{\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{1/2},$$

(21a)

$${\rho }^{\prime }=r{\left\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{1/2},$$

(21b)

and

$$\cos \left(\psi (\theta )\right)=\frac{\cos (\theta )+\xi \sin (\theta )}{{\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{1/2}},$$

(22a)

$$\sin \left(\psi (\theta )\right)=\frac{\eta \sin (\theta )}{{\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{1/2}},$$

(22b)

$$\cos \left({\psi }^{\prime }(\theta )\right)=\frac{\cos (\theta )-\xi \sin (\theta )}{{\left\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{1/2}},$$

(22c)

$$\sin \left({\psi }^{\prime }(\theta )\right)=\frac{-\eta \sin (\theta )}{{\left\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{1/2}},$$

(22d)

the general asymptotic form for the displacement and stress fields can be written as follows:

$$u(r,\theta ,z)=r^s{D}_b\left(z\right){\left(ik\right)}^s[{\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{s/2}\left\{A_1\cos \left(s\psi \right)+A_2\sin \left(s\psi \right)\right\}$$

$$+{\left\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{s/2}\left\{A_3\cos \left(s{\psi }^{\prime }\right)+A_4\sin \left(s{\psi }^{\prime }\right)\right\}]+O\left(r^{s+2}\right),$$

(23a)

$$v(r,\theta ,z)=r^s{D}_b\left(z\right){\left(ik\right)}^s[{\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{s/2}\{(H_1{A}_1$$

$$+H_2{A}_2)\cos \left(s\psi \right)+\left(H_1{A}_2-H_2{A}_1\right)\sin \left(s\psi \right)\}+\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2$$

$${+{\eta }^2{\sin }^2(\theta )\}}^{s/2}\left\{-\left(H_1{A}_3+H_2{A}_4\right)\cos \left(s{\psi }^{\prime }\right)-\left(H_1{A}_4-H_2{A}_3\right)\sin \left(s{\psi }^{\prime }\right)\right\}$$

$$+O\left(r^{s+2}\right),$$

(23b)

$$w(r,\theta ,z)=O\left(r^{s+1}\right),$$

(23c)

and

$${\sigma }_x\left(r,\theta ,z\right)=r^{s-1}{D}_b\left(z\right){\left(ik\right)}^{s}s\langle {\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{(s-1)/2}[(A_1\{c_{11}+$$

$$+\left(\xi H_1-\eta H_2\right)c_{12}\}+A_2\left\{\left(\eta H_1+\xi H_2\right)c_{12}\right\})\cos \left(\left(s-1\right)\psi \right)+(-A_1\left\{\left(\eta H_1+\xi H_2\right)c_{12}\right\}$$

$$+A_2\{c_{11}+\left(\xi H_1-\eta H_2\right)c_{12}\})\sin \left(\left(s-1\right)\psi \right)]+\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2$$

$${+{\eta }^2{\sin }^2(\theta )\}}^{(s-1)/2}[(A_3\{c_{11}+\left(\xi H_1-\eta H_2\right)c_{12}\}+A_4\left\{\left(\eta H_1+\xi H_2\right)c_{12}\right\})\cos \left(\left(s-1\right){\psi }^{\prime }\right)$$

$$+(-A_3\left\{\left(\eta H_1+\xi H_2\right)c_{12}\right\}+A_4\{c_{11}+\left(\xi H_1-\eta H_2\right)c_{12}\})\sin \left(\left(s-1\right){\psi }^{\prime }\right)]\rangle$$

$$+O\left(r^{s+1}\right),$$

(24a)

$${\sigma }_y\left(r,\theta ,z\right)=r^{s-1}{D}_b\left(z\right){\left(ik\right)}^{s}s\langle {\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{(s-1)/2}[(A_1\{c_{12}$$

$$+\left(\xi H_1-\eta H_2\right)c_{22}\}+A_2\left\{\left(\eta H_1+\xi H_2\right)c_{22}\right\})\cos \left(\left(s-1\right)\psi \right)+(-A_1\left\{\left(\eta H_1+\xi H_2\right)c_{22}\right\}$$

$$+A_2\{c_{12}+\left(\xi H_1-\eta H_2\right)c_{22}\})\sin \left(\left(s-1\right)\psi \right)]+\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2$$

$${+{\eta }^2{\sin }^2(\theta )\}}^{(s-1)/2}[(A_3\{c_{12}+\left(\xi H_1-\eta H_2\right)c_{22}\}+A_4\left\{\left(\eta H_1+\xi H_2\right)c_{22}\right\})\cos \left(\left(s-1\right){\psi }^{\prime }\right)$$

$$+(-A_3\left\{\left(\eta H_1+\xi H_2\right)c_{22}\right\}+A_4\{c_{12}+\left(\xi H_1-\eta H_2\right)c_{22}\})\sin \left(\left(s-1\right){\psi }^{\prime }\right)]\rangle$$

$$+O\left(r^{s+1}\right),$$

(24b)

$${\tau }_{xy}\left(r,\theta ,z\right)=r^{s-1}{D}_b\left(z\right){\left(ik\right)}^{s}s{c}_{66}\langle {\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{(s-1)/2}[\{A_1\left(\xi +H_1\right)$$

$$+A_2\left(\eta +H_2\right)\}\cos \left(\left(s-1\right)\psi \right)+\left\{-A_1\left(\eta +H_2\right)+A_2\left(\xi +H_1\right)\right\}\sin \left(\left(s-1\right)\psi \right)]$$

$${+\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\}}^{(s-1)/2}[\{-A_3\left(\xi +H_1\right)-A_4\left(\eta +H_2\right)\}\cos \left(\left(s-1\right){\psi }^{\prime }\right)$$

$$+\left\{A_3\left(\eta +H_2\right)-A_4\left(\xi +H_1\right)\right\}\sin \left(\left(s-1\right){\psi }^{\prime }\right)]+O\left(r^{s+1}\right),$$

(24c)

$${\sigma }_z\left(r,\theta ,z\right)=r^{s-1}{D}_b\left(z\right){\left(ik\right)}^{s}s\langle {\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{(s-1)/2}[(A_1\{c_{13}$$

$$+\left(\xi H_1-\eta H_2\right)c_{23}\}+A_2\left\{\left(\eta H_1+\xi H_2\right)c_{23}\right\})\cos \left(\left(s-1\right)\psi \right)+(-A_1\left\{\left(\eta H_1+\xi H_2\right)c_{23}\right\}$$

$$+A_2\{c_{13}+\left(\xi H_1-\eta H_2\right)c_{23}\})\sin \left(\left(s-1\right)\psi \right)]+\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2$$

$${+{\eta }^2{\sin }^2(\theta )\}}^{(s-1)/2}[(A_3\{c_{13}+\left(\xi H_1-\eta H_2\right)c_{23}\}+A_4\left\{\left(\eta H_1+\xi H_2\right)c_{23}\right\})\cos \left(\left(s-1\right){\psi }^{\prime }\right)$$

$$+(-A_3\left\{\left(\eta H_1+\xi H_2\right)c_{23}\right\}+A_4\{c_{13}+\left(\xi H_1-\eta H_2\right)c_{23}\})\sin \left(\left(s-1\right){\psi }^{\prime }\right)]\rangle$$

$$+O\left(r^{s+1}\right),$$

(24d)

$${\tau }_{xz}\left(r,\theta ,z\right)=O\left(r^s\right),\begin{array}{cc}& {\tau }_{yz}\left(r,\theta ,z\right)=O\left(r^s\right)\end{array}.$$

(24e, f)

in which

$$A_1={\overline{A}}_1+{\overline{A}}_2,\begin{array}{cc}& \end{array}{A}_2=i\left({\overline{A}}_1-{\overline{A}}_2\right),$$

(25a)

$$A_3={\overline{A}}_3+{\overline{A}}_4,\begin{array}{cc}& \end{array}{A}_4=i\left({\overline{A}}_3-{\overline{A}}_4\right),$$

(25b)

and

$$D_b\left(z\right)=D_1\sin \left(kz\right)+D_2\cos \left(kz\right),$$

(26)

with

$$D_1=i{\overline{\overline{D}}}_1,\hspace{1em}{D}_2={\overline{\overline{D}}}_2.$$

(27a, b)

It may be noted that since s or Re s (when s is complex) is positive, all the higher order terms in Eqs. (24) vanish as r → 0. The components of displacement can now be expressed in the cylindrical polar coordinate system as follows:

$$u_r(r,\theta ,z)=r^s{D}_b\left(z\right){\left(ik\right)}^s\langle {\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{s/2}[\{A_1\cos \left(\theta \right)+$$

$$\left(H_1{A}_1+H_2{A}_2\right)\sin (\theta )\}\cos \left(s\psi \right)+\left\{A_2\cos \left(\theta \right)+\left(H_1{A}_2-H_2{A}_1\right)\sin (\theta )\right\}\sin \left(s\psi \right)]$$

$$+{\left\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{s/2}[\left\{A_3\cos \left(\theta \right)-\left(H_1{A}_3+H_2{A}_4\right)\sin (\theta )\right\}\cos \left(s{\psi }^{\prime }\right)$$

$$+\left\{A_4\cos \left(\theta \right)-\left(H_1{A}_4-H_2{A}_3\right)\sin (\theta )\right\}\sin \left(s{\psi }^{\prime }\right)]\rangle +O\left(r^{s+2}\right),$$

(28a)

$$u_{\theta }(r,\theta ,z)=r^s{D}_b\left(z\right){\left(ik\right)}^s\langle {\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{s/2}[\{-A_1\sin \left(\theta \right)+$$

$$\left(H_1{A}_1+H_2{A}_2\right)\cos (\theta )\}\cos \left(s\psi \right)+\left\{-A_2\sin \left(\theta \right)+\left(H_1{A}_2-H_2{A}_1\right)\cos (\theta )\right\}\sin \left(s\psi \right)]$$

$$+{\left\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{s/2}[-\{A_3\sin \left(\theta \right)+(H_1{A}_3+$$

$$H_2{A}_4)\cos (\theta )\}\cos \left(s{\psi }^{\prime }\right)-\left\{A_4\sin \left(\theta \right)+\left(H_1{A}_4-H_2{A}_3\right)\cos (\theta )\right\}\sin \left(s{\psi }^{\prime }\right)]\rangle$$

$$+O\left(r^{s+2}\right),$$

(28b)

$$w(r,\theta ,z)=O\left(r^{s+1}\right).$$

(28c)

Similarly, the components of the asymptotic stress field can be conveniently expressed by using standard second-rank tensor transformation rule:

$$\left\{\begin{array}{c}{\sigma }_r\\ {\sigma }_{\theta }\\ {\tau }_{r\theta }\end{array}\right\}=\left[\begin{array}{ccc}{\cos }^2\theta & {\sin }^2\theta & \sin 2\theta \\ {\sin }^2\theta & {\cos }^2\theta & -\sin 2\theta \\ -\frac{1}{2}\sin 2\theta & \frac{1}{2}\sin 2\theta & \cos 2\theta \end{array}\right]\left\{\begin{array}{c}{\sigma }_x\\ {\sigma }_y\\ {\tau }_{xy}\end{array}\right\},$$

(29a)

$$\left\{\begin{array}{c}{\tau }_{rz}\\ {\tau }_{\theta z}\end{array}\right\}=\left[\begin{array}{cc}\cos \theta & \sin \theta \\ -\sin \theta & \cos \theta \end{array}\right]\left\{\begin{array}{c}{\tau }_{xz}\\ {\tau }_{yz}\end{array}\right\}.$$

(29b)

The stress component, σz, is as given in Eq. (24d).

Substitution of Eqs. (24), in conjunction with Eqs. (29), into the boundary conditions on the crack-side surfaces given by Eq. (3b) yield four homogeneous equations, which reduce to

either cos (s-1)π = 0, or sin (s-1)π = 0.

(30a, b)

Equation (30a) contributes to the lowest nonvanishing eigenvalue, s = 1/2, in the range 0 < s < 1, as required by the criterion of locally finite energy. Equation (30b) yields s = 0, 1, which can take care of rigid body translation and rotation, respectively. Interestingly, s = 1 also accounts for the T-stress.

In order to express the asymptotic stress field in terms of the stress intensity factors K_I(z) and K_II(z), it is convenient to obtain from Eq. (10)

$$c_{12}{\overline{A}}_i+c_{11}{p}_i{\overline{B}}_i=-c_{66}\left({\overline{A}}_i+{\overline{p}}_i{\overline{B}}_i\right), \mathrm{i}=1,\dots ,4 (\mathrm{no} \mathrm{sum} \mathrm{on} \mathrm{i})$$

(31)

in which ${\overline{p}}_i$ is the complex conjugate of $p_i$ given by Eqs. (15) and (16), and ${\overline{B}}_i$ is given by Eqs. (18). This will finally yield, with the help of Eqs. (17), (18) and (1), ${\sigma }_y\left(r,\theta ,z\right)$ as follows:

$${\sigma }_y\left(r,\theta ,z\right)=-r^{s-1}{D}_b\left(z\right){\left(ik\right)}^{s}s{c}_{66}\langle {\left\{{\left(\cos (\theta )+\xi \sin (\theta )\right)}^2+{\eta }^2{\sin }^2(\theta )\right\}}^{(s-1)/2}[(A_1\{1$$

$$\xi H_1+\eta H_2\}-A_2\left\{\eta H_1-\xi H_2\right\})\cos \left(\left(s-1\right)\psi \right)-(A_1\left\{\eta H_1-\xi H_2\right\}+A_2\{1+\xi H_1$$

$$+\eta H_2\})\sin \left(\left(s-1\right)\psi \right)]+\{{\left(\cos (\theta )-\xi \sin (\theta )\right)}^2{+{\eta }^2{\sin }^2(\theta )\}}^{(s-1)/2}[(A_3\{1+\xi H_1$$

$$+\eta H_2\}-A_4\left\{\eta H_1-\xi H_2\right\})\cos \left(\left(s-1\right){\psi }^{\prime }\right)-(A_3\left\{\eta H_1-\xi H_2\right\}+A_4\{1+\xi H_1$$

$$+\eta H_2\})\sin \left(\left(s-1\right){\psi }^{\prime }\right)]\rangle +O\left(r^{s+1}\right).$$

(32)

${\sigma }_x\left(r,\theta ,z\right)$ can be similarly obtained from the other relation in Eq. (10).

4.1.1. Symmetric (Mode I) Loading (Extension/Bending)

This relationship among the four unknown coefficients (eigenvectors) under mode I can also be obtained by using the following boundary conditions that correspond to the far-field symmetric loading:

$$\mathsf{\theta }=0:{\tau }_{r\theta }={\tau }_{\theta z}=0,$$

(33a)

$$\mathsf{\theta }=\mathsf{\pi }:{\sigma }_{\theta }={\tau }_{r\theta }={\tau }_{\theta z}=0.$$

(33b)

When s = 1/2, substitution, of Eqs. (29) or equivalently (24) into Eqs. (33), yields the following:

A1 = A3,   A2 = A4,

(34)

$$\frac{A_2}{A_1}=-\frac{\eta }{\xi }\frac{\left(\sqrt{c_{11}{c}_{22}}+c_{12}\right)}{\left(\sqrt{c_{11}{c}_{22}}-c_{12}\right)}.$$

(35)

Finally, on substitution of Eqs. (34) and (35) into the expressions for displacements and stresses, and by defining the mode I stress intensity factor, K_I(z), as

$$K_I(z)=\sqrt{2\pi }{D}_b(z){(ik)}^{1/2}\left(\sqrt{c_{11}{c}_{22}}+c_{12}\right)A_1,$$

(36)

the components of inplane displacements and stresses, in the vicinity of a semi-infinite crack front, under symmetric far-field loading, can be expressed as follows:

$$u(r,\theta ,z)=\frac{K_I(z)}{\left(c_{11}{c}_{22}-c_{12}^2\right)}\sqrt{\frac{r}{2\pi }}{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\{\left(\sqrt{c_{11}{c}_{22}}-c_{12}\right)\cos \left(\psi /2\right)$$

$$-\left(\sqrt{c_{11}{c}_{22}}+c_{12}\right)\frac{\eta }{\xi }\sin \left(\psi /2\right)\}+{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\{\left(\sqrt{c_{11}{c}_{22}}-c_{12}\right)\cos \left({\psi }^{\prime }/2\right)$$

$$+\left(\sqrt{c_{11}{c}_{22}}+c_{12}\right)\frac{\eta }{\xi }\sin \left({\psi }^{\prime }/2\right)\}],$$

(37a)

$$v(r,\theta ,z)=\frac{K_I(z)}{\left(c_{11}{c}_{22}-c_{12}^2\right)}\sqrt{\frac{r}{2\pi }}{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\{-\frac{\left(c_{11}{c}_{22}-c_{12}^2\right)}{2c_{66}\xi }\cos \left(\psi /2\right)$$

$$+2\sqrt{c_{11}{c}_{22}}\eta \sin \left(\psi /2\right)\}+{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\{\frac{\left(c_{11}{c}_{22}-c_{12}^2\right)}{2c_{66}\xi }\cos \left({\psi }^{\prime }/2\right)$$

$$+2\sqrt{c_{11}{c}_{22}}\eta \sin \left({\psi }^{\prime }/2\right)\}],$$

(37b)

$${\sigma }_x(r,\theta ,z)=\frac{K_I(z)}{2\sqrt{2\pi r}}{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\left\{\cos \left(\psi /2\right)-\frac{\eta }{\xi }\sin \left(\psi /2\right)\right\}$$

$$+{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\left\{\cos \left({\psi }^{\prime }/2\right)+\frac{\eta }{\xi }\sin \left({\psi }^{\prime }/2\right)\right\}],$$

(38a)

$${\sigma }_y(r,\theta ,z)=\frac{K_I(z)}{2\sqrt{2\pi r}}{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\left\{\cos \left(\psi /2\right)+\frac{\eta }{\xi }\sin \left(\psi /2\right)\right\}$$

$$+{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\left\{\cos \left({\psi }^{\prime }/2\right)-\frac{\eta }{\xi }\sin \left({\psi }^{\prime }/2\right)\right\}],$$

(38b)

$${\tau }_{xy}(r,\theta ,z)=-\frac{K_I(z)}{2\sqrt{2\pi r}\xi }{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\cos \left(\psi /2\right)$$

$$-{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\cos \left({\psi }^{\prime }/2\right)].$$

(38c)

4.1.2. Skew-symmetric (Mode II) Loading (Sliding Shear/Twisting)

This relationship among the four unknown coefficients (eigenvectors) under mode II can also be obtained by using the following boundary conditions that correspond to the far-field antisymmetric loading:

$$\mathsf{\theta }=0:{\sigma }_{\theta }={\tau }_{\theta z}=0,$$

(39a)

$$\mathsf{\theta }=\mathsf{\pi }:{\sigma }_{\theta }={\tau }_{r\theta }={\tau }_{\theta z}=0.$$

(39b)

When s = 1/2, substitution of Eqs. (29) or equivalently (24) into Eqs. (39) yields the following:

A1 = -A3,   A2 = -A4.

(40)

$$\frac{A_2}{A_1}=\frac{4\sqrt{c_{11}{c}_{22}}{c}_{66}\xi \eta }{\left(c_{11}{c}_{22}-c_{12}^2\right)}.$$

(41)

Finally, on substitution of Eqs. (40) and (41) into the expressions for displacements and stresses, and by defining the mode II stress intensity factor, K_II(z), as

$$K_{II}(z)=2\sqrt{2\pi }{D}_b(z){(ik)}^{1/2}\xi c_{66}{A}_1,$$

(42)

the components of inplane displacements and stresses, in the vicinity of a semi-infinite crack front, under skew-symmetric far-field loading, can be expressed as follows:

$$u(r,\theta ,z)=\frac{K_{II}(z)}{\left(c_{11}{c}_{22}-c_{12}^2\right)}\sqrt{\frac{r}{2\pi }}{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\{\frac{\left(c_{11}{c}_{22}-c_{12}^2\right)}{2c_{66}\xi }\cos \left(\psi /2\right)$$

$$+2\sqrt{c_{11}{c}_{22}}\eta \sin \left(\psi /2\right)\}+{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\{-\frac{\left(c_{11}{c}_{22}-c_{12}^2\right)}{2c_{66}\xi }\cos \left({\psi }^{\prime }/2\right)$$

$$+2\sqrt{c_{11}{c}_{22}}\eta \sin \left({\psi }^{\prime }/2\right)\}],$$

(43a)

$$v(r,\theta ,z)=-\frac{K_{II}(z)}{\left(c_{11}{c}_{22}-c_{12}^2\right)}\sqrt{\frac{r}{2\pi }}{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\{\left(\sqrt{c_{11}{c}_{22}}-c_{12}\right)\cos \left(\psi /2\right)$$

$$+\left(\sqrt{c_{11}{c}_{22}}+c_{12}\right)\frac{\eta }{\xi }\sin \left(\psi /2\right)\}+{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\{\left(\sqrt{c_{11}{c}_{22}}-c_{12}\right)\cos \left({\psi }^{\prime }/2\right)$$

$$-\left(\sqrt{c_{11}{c}_{22}}+c_{12}\right)\frac{\eta }{\xi }\sin \left({\psi }^{\prime }/2\right)\}],$$

(43b)

$${\sigma }_x(r,\theta ,z)=\frac{K_{II}(z)}{2\sqrt{2\pi r}}{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\left\{\frac{\left({\eta }^2-{\xi }^2\right)}{\xi }\cos \left(\psi /2\right)-2\eta \sin \left(\psi /2\right)\right\}$$

$$-{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\left\{\frac{\left({\eta }^2-{\xi }^2\right)}{\xi }\cos \left({\psi }^{\prime }/2\right)+2\eta \sin \left({\psi }^{\prime }/2\right)\right\}],$$

(44a)

$${\sigma }_y(r,\theta ,z)=\frac{K_{II}(z)}{2\sqrt{2\pi r}\xi }{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\cos \left(\psi /2\right)$$

$$+{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\cos \left({\psi }^{\prime }/2\right)],$$

(44b)

$${\tau }_{xy}(r,\theta ,z)=-\frac{K_{II}(z)}{2\sqrt{2\pi r}}{[\left\{{\left(\cos \theta +\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\left\{\cos \left(\psi /2\right)+\frac{\eta }{\xi }\sin \left(\psi /2\right)\right\}$$

$$+{\left\{{\left(\cos \theta -\xi \sin \theta \right)}^2+\eta {\sin }^2\theta \right\}}^{1/4}\left\{\cos \left({\psi }^{\prime }/2\right)-\frac{\eta }{\xi }\sin \left({\psi }^{\prime }/2\right)\right\}].$$

(44c)

It may be noted that in the expressions for displacements and stresses, given by Eqs. (37) and (38), respectively, under mode I loading, there is a coupling between cos(ψ/2) and sin(ψ/2) (and similar coupling between cos(ψ’/2) and sin(ψ’/2)) terms, except for τxy, given by Eq. (38c). Similarly, in the expressions for displacements and stresses, given by Eqs. (43) and (44), respectively, under mode II loading, there is a coupling between cos(ψ/2) and sin(ψ/2) (and similar coupling between cos(ψ’/2) and sin(ψ’/2)) terms, except for σy, given by Eq. (44b). It may be remarked that for $\chi <\sqrt{c_{11}/c_{22}}$, the through-thickness crack considered here will more likely deviate from its plane (010), i.e., x-z plane and its initial propagation direction, i.e., x-direction, because of this coupling, which arises from the complex roots. It can then be inferred that the (010) is not the easy cleavage plane and [100] is not the easy propagation direction, for orthotropic/orthorhombic laminas/single crystals with $\chi <\sqrt{c_{11}/c_{22}}$.

It may, however, be remarked that the Griffith-Irwin approach cannot by itself be considered a sufficient condition for determination of a cleavage system being easy or difficult for crack propagation in single crystals, and an additional conceptual-cum-analytical tool must be developed. This is because Griffith’s criterion is “not really a fracture criterion but only a necessary condition for fracture” [79]. As has been suggested by Sedov [82], similarity analysis is an effective tool to solve complex problems in mechanics. This can be employed to come up with a sorely needed sufficient condition for the problem at hand and is elaborated in Section 7.2 below.

4.2. Case (b): Imaginary Roots

The four imaginary roots of Eq. (11) are given by

$$p_{1,2}=\pm i\left({\xi }^{\prime }+{\eta }^{\prime }\right),\hspace{1em}{p}_{3,4}=\pm i\left({\xi }^{\prime }-{\eta }^{\prime }\right),$$

(45)

where

$${\xi }^{\prime }=\frac{1}{\sqrt{2}}{\left[{\left(\frac{c_{11}}{c_{22}}\right)}^{1/2}+\chi \right]}^{1/2},\hspace{1em}{\eta }^{\prime }=\frac{1}{\sqrt{2}}{\left[-{\left(\frac{c_{11}}{c_{22}}\right)}^{1/2}+\chi \right]}^{1/2},$$

(46a,b)

valid for

$$\chi >\sqrt{c_{11}/c_{22}},$$

(47)

which corresponds to a candidate plane of minimum surface energy.

The final results that satisfy the equilibrium equations (2) can be expressed in the following form:

$$u(x,y,z)=\left({\overline{\overline{D}}}_{1}i\sin \left(kz\right)+{\overline{\overline{D}}}_2\cos \left(kz\right)\right){\left(ik\right)}^s[{\overline{A}}_1{\left(x+i({\xi }^{\prime }+{\eta }^{\prime })y\right)}^s+{\overline{A}}_2{\left(x-i({\xi }^{\prime }+{\eta }^{\prime })y\right)}^s$$

$$+{\overline{A}}_3{\left(x+i({\xi }^{\prime }-{\eta }^{\prime })y\right)}^s+{\overline{A}}_4{\left(x-i({\xi }^{\prime }-{\eta }^{\prime })y\right)}^s],$$

(48a)

$$v(x,y,z)=\left({\overline{\overline{D}}}_{1}i\sin \left(kz\right)+{\overline{\overline{D}}}_2\cos \left(kz\right)\right){\left(ik\right)}^s[{\overline{B}}_1{\left(x+i({\xi }^{\prime }+{\eta }^{\prime })y\right)}^s+{\overline{B}}_2{\left(x-i({\xi }^{\prime }+{\eta }^{\prime })y\right)}^s$$

$$+{\overline{B}}_3{\left(x+i({\xi }^{\prime }-{\eta }^{\prime })y\right)}^s+{\overline{B}}_4{\left(x-i({\xi }^{\prime }-{\eta }^{\prime })y\right)}^s],$$

(48b)

$$w(x,y,z)=\left({\overline{\overline{D}}}_1\cos \left(kz\right)+{\overline{\overline{D}}}_{2}i\sin \left(kz\right)\right){\left(ik\right)}^{s+1}[{\overline{C}}_1{\left(x+i({\xi }^{\prime }+{\eta }^{\prime })y\right)}^{s+1}+{\overline{C}}_2{\left(x-i({\xi }^{\prime }+{\eta }^{\prime })y\right)}^{s+1}$$

$$+{\overline{C}}_3{\left(x+i({\xi }^{\prime }-{\eta }^{\prime })y\right)}^s+{\overline{C}}_4{\left(x-i({\xi }^{\prime }-{\eta }^{\prime })y\right)}^s],$$

(48c)

wherein ξ′ and η′ are as given in Eqs. (46), and ${\overline{A}}_k, {\overline{B}}_k, {\overline{C}}_k$, k = 1,...,4, are undetermined coefficients. It may be noted that ${\overline{B}}_k$ can be expressed in terms of the corresponding ${\overline{A}}_k$, k = 1,...,4, by using Eq. (10) and Eqs. (45).

$${\overline{B}}_1=-i{H}_1^{\prime }{\overline{A}}_1,\hspace{1em}{\overline{B}}_2=i{H}_1^{\prime }{\overline{A}}_2,$$

(49a,b)

$${\overline{B}}_3=-i{H}_2^{\prime }{\overline{A}}_3,\hspace{1em}{\overline{B}}_4=i{H}_2^{\prime }{\overline{A}}_4,$$

(49c,d)

in which

$$H_1^{\prime }=-\frac{\left\{\sqrt{c_{11}{c}_{22}}-c_{66}{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2\right\}}{\left(c_{12}+c_{66}\right)\left({\xi }^{\prime }+{\eta }^{\prime }\right)},\hspace{1em}{H}_2^{\prime }=\frac{\left\{\sqrt{c_{11}{c}_{22}}-c_{66}{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2\right\}}{\left(c_{12}+c_{66}\right)\left({\xi }^{\prime }-{\eta }^{\prime }\right)}.$$

(50a,b)

The corresponding stress field can easily be obtained from Eqs. (48). It is convenient to express the components of the displacement vector and stress tensor, in terms of the cylindrical polar coordinate system (r, θ, z). Expressing

$${\rho }_1\cos ({\psi }_1(\theta ))=r\cos (\theta ),\hspace{1em}{\rho }_1\sin ({\psi }_1(\theta ))=r({\xi }^{\prime }+{\eta }^{\prime })\sin (\theta ),$$

(51a)

$${\rho }_1^{\prime }\cos ({\psi }_1^{\prime }(\theta ))=r\cos (\theta ),\hspace{1em}{{\rho }^{\prime }}_1\sin ({\psi }_1^{\prime }(\theta ))=r({\xi }^{\prime }-{\eta }^{\prime })\sin (\theta ),$$

(51b)

in which

$${\rho }_1=r{\left\{{\cos }^2(\theta )+{({\xi }^{\prime }+{\eta }^{\prime })}^2{\sin }^2(\theta )\right\}}^{1/2},$$

(52a)

$${\rho }_1^{\prime }=r{\left\{{\cos }^2(\theta )+{({\xi }^{\prime }-{\eta }^{\prime })}^2{\sin }^2(\theta )\right\}}^{1/2},$$

(52b)

and

$$\cos \left({\psi }_1(\theta )\right)=\frac{\cos (\theta )}{{\left\{{\cos }^2(\theta )+{({\xi }^{\prime }+{\eta }^{\prime })}^2{\sin }^2(\theta )\right\}}^{1/2}},$$

(53a)

$$\sin \left({\psi }_1(\theta )\right)=\frac{({\xi }^{\prime }+{\eta }^{\prime })\sin (\theta )}{{\left\{{\cos }^2(\theta )+{({\xi }^{\prime }+{\eta }^{\prime })}^2{\sin }^2(\theta )\right\}}^{1/2}},$$

(53b)

$$\cos \left({{\psi }^{\prime }}_1(\theta )\right)=\frac{\cos (\theta )}{{\left\{{\cos }^2(\theta )+{({\xi }^{\prime }-{\eta }^{\prime })}^2{\sin }^2(\theta )\right\}}^{1/2}},$$

(53c)

$$\sin \left({{\psi }^{\prime }}_1(\theta )\right)=\frac{({\xi }^{\prime }-{\eta }^{\prime })\sin (\theta )}{{\left\{{\cos }^2(\theta )+{({\xi }^{\prime }-{\eta }^{\prime })}^2{\sin }^2(\theta )\right\}}^{1/2}},$$

(53d)

the general asymptotic form for the displacement and stress fields can be written as follows:

$$u(r,\theta ,z)=r^s{D}_b\left(z\right){\left(ik\right)}^s[{\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{s/2}\left\{A_1\cos \left(s{\psi }_1\right)+A_2\sin \left(s{\psi }_1\right)\right\}$$

$$+{\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{s/2}\left\{A_3\cos \left(s{\psi }_1^{\prime }\right)+A_4\sin \left(s{\psi }_1^{\prime }\right)\right\}]+O\left(r^{s+2}\right),$$

(54a)

$$v(r,\theta ,z)=r^s{D}_b\left(z\right){\left(ik\right)}^s[{\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{s/2}\left\{-H_1^{\prime }{A}_2\cos \left(s{\psi }_1\right)+H_1^{\prime }{A}_1\sin \left(s{\psi }_1\right)\right\}$$

$$+{\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{s/2}\left\{-{H^{\prime }}_2{A}_4\cos \left(s{\psi }_1^{\prime }\right)+{H^{\prime }}_2{A}_3\sin \left(s{\psi }_1^{\prime }\right)\right\}]+O\left(r^{s+2}\right),$$

(54b)

$$w(r,\theta ,z)=O\left(r^{s+1}\right),$$

(54c)

and

$${\sigma }_x\left(r,\theta ,z\right)=r^{s-1}{D}_b\left(z\right){\left(ik\right)}^{s}s\langle {\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{(s-1)/2}\{c_{11}+$$

$$H_1^{\prime }\left({\xi }^{\prime }+{\eta }^{\prime }\right)c_{12}\}\left\{A_1\cos \left(\left(s-1\right){\psi }_1\right)+A_2\sin \left(\left(s-1\right){\psi }_1\right)\right\}+\{{\cos }^2(\theta )+$$

$${\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2(\theta ){\}}^{(s-1)/2}\{c_{11}+H_2^{\prime }\left({\xi }^{\prime }-{\eta }^{\prime }\right)c_{12}\}\left\{A_3\cos \left(\left(s-1\right){{\psi }^{\prime }}_1\right)+A_4\sin \left(\left(s-1\right){{\psi }^{\prime }}_1\right)\right\}\rangle$$

$$+O\left(r^{s+1}\right),$$

(55a)

$${\sigma }_y\left(r,\theta ,z\right)=r^{s-1}{D}_b\left(z\right){\left(ik\right)}^{s}s\langle {\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{(s-1)/2}\{c_{12}+$$

$$H_1^{\prime }\left({\xi }^{\prime }+{\eta }^{\prime }\right)c_{22}\}\left\{A_1\cos \left(\left(s-1\right){\psi }_1\right)+A_2\sin \left(\left(s-1\right){\psi }_1\right)\right\}+\{{\cos }^2(\theta )+$$

$${\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2(\theta ){\}}^{(s-1)/2}\{c_{12}+H_2^{\prime }\left({\xi }^{\prime }-{\eta }^{\prime }\right)c_{22}\}\left\{A_3\cos \left(\left(s-1\right){{\psi }^{\prime }}_1\right)+A_4\sin \left(\left(s-1\right){{\psi }^{\prime }}_1\right)\right\}\rangle$$

$$+O\left(r^{s+1}\right),$$

(55b)

$${\tau }_{xy}\left(r,\theta ,z\right)=r^{s-1}{D}_b\left(z\right){\left(ik\right)}^{s}s{c}_{66}\langle {\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{(s-1)/2}\{{H^{\prime }}_1$$

$$-\left({\xi }^{\prime }+{\eta }^{\prime }\right)\}\left\{-A_2\cos \left((s-1){\psi }_1\right)+A_1\sin \left((s-1){\psi }_1\right)\right\}+\{{\cos }^2(\theta )+$$

$${{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\}}^{(s-1)/2}\{H_2^{\prime }-\left({\xi }^{\prime }-{\eta }^{\prime }\right)\}\left\{-A_4\cos \left((s-1){\psi }_1\right)+A_3\sin \left((s-1){\psi }_1\right)\right\}\rangle$$

$$+O\left(r^{s+1}\right),$$

(55c)

$${\sigma }_z\left(r,\theta ,z\right)=r^{s-1}{D}_b\left(z\right){\left(ik\right)}^{s}s\langle {\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{(s-1)/2}\{c_{13}+$$

$$H_1^{\prime }\left({\xi }^{\prime }+{\eta }^{\prime }\right)c_{23}\}\left\{A_1\cos \left(\left(s-1\right){\psi }_1\right)+A_2\sin \left(\left(s-1\right){\psi }_1\right)\right\}+\{{\cos }^2(\theta )+$$

$${\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2(\theta ){\}}^{(s-1)/2}\{c_{13}+H_2^{\prime }\left({\xi }^{\prime }-{\eta }^{\prime }\right)c_{23}\}\left\{A_3\cos \left(\left(s-1\right){{\psi }^{\prime }}_1\right)+A_4\sin \left(\left(s-1\right){{\psi }^{\prime }}_1\right)\right\}\rangle$$

$$+O\left(r^{s+1}\right),$$

(55d)

$${\tau }_{xz}\left(r,\theta ,z\right)=O\left(r^s\right),\begin{array}{cc}& {\tau }_{yz}\left(r,\theta ,z\right)=O\left(r^s\right)\end{array}.$$

(55e, f)

in which Ak, k =1, ..., 4, is defined as before in Eqs. (23), while $D_b\left(z\right)$ is same as given earlier in Eqs. (24) and (25).

It may be noted that since s or Re s (when s is complex) is positive, all the higher order terms in Eqs. (55) vanish as r → 0. The components of displacement can now be expressed in the cylindrical polar coordinate system as follows:

$$u_r(r,\theta ,z)=r^s{D}_b\left(z\right){\left(ik\right)}^s\langle {\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{s/2}[\left\{A_1\cos (\theta )-H_1^{\prime }{A}_2\sin (\theta )\right\}\cos \left(s{\psi }_1\right)$$

$$+\left\{A_2\cos (\theta )+H_1^{\prime }{A}_1\sin (\theta )\right\}\sin \left(s{\psi }_1\right)]+{\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{s/2}[\{A_3\cos (\theta )$$

$$-H_2^{\prime }{A}_4\sin (\theta )\}\cos \left(s{\psi }_1^{\prime }\right)+\left\{A_4\cos (\theta )+H_2^{\prime }{A}_3\sin (\theta )\right\}\sin \left(s{\psi }_1^{\prime }\right)]\rangle$$

$$+O\left(r^{s+2}\right),$$

(56a)

$$u_{\theta }(r,\theta ,z)=r^s{D}_b\left(z\right){\left(ik\right)}^s\langle {\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{s/2}[-\left\{A_1\sin (\theta )+H_1^{\prime }{A}_2\cos (\theta )\right\}\cos \left(s{\psi }_1\right)$$

$$+\left\{-A_2\sin (\theta )+H_1^{\prime }{A}_1\cos (\theta )\right\}\sin \left(s{\psi }_1\right)]+{\left\{{\cos }^2(\theta )+{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2(\theta )\right\}}^{s/2}[-\{A_3\sin (\theta )$$

$$+H_2^{\prime }{A}_4\cos (\theta )\}\cos \left(s{\psi }_1^{\prime }\right)+\left\{-A_4\sin (\theta )+H_2^{\prime }{A}_3\cos (\theta )\right\}\sin \left(s{\psi }_1^{\prime }\right)]\rangle$$

$$+O\left(r^{s+2}\right),$$

(56b)

$$w(r,\theta ,z)=O\left(r^{s+1}\right).$$

(56c)

Similarly, the components of the asymptotic stress field can be conveniently expressed by using standard transformation rule, given earlier in Eqs. (27). The stress component, σz, is as given in Eq. (55d).

Substitution of Eqs. (55) in conjunction with Eqs. (27) into the boundary conditions on the crack-side surfaces given by Eq. (3b) yield four homogeneous equations, which reduce to

either cos (s-1)π = 0, or sin (s-1)π = 0.

(57a, b)

Equation (57a) contributes to the lowest nonvanishing eigenvalue, s = 1/2, in the range 0 < s < 1, as required by the criterion of locally finite energy. Equation (57b) yields s = 0, 1, which can take care of rigid body translation and rotation, respectively. Again, s = 1 also accounts for the T-stress.

4.2.1. Symmetric (Mode I) Loading (Extension/Bending)

For s = 1/2,

A2 = A4 = 0;

(58)

and

$$\frac{A_3}{A_1}=-\frac{\left\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }-{\eta }^{\prime }\right)+c_{12}\left({\xi }^{\prime }+{\eta }^{\prime }\right)\right\}}{\left\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }+{\eta }^{\prime }\right)+c_{12}\left({\xi }^{\prime }-{\eta }^{\prime }\right)\right\}}.$$

(59)

Finally, on substitution of Eqs. (58) and (59) into the expressions for displacements and stresses, and by defining the mode I stress intensity factor, K_I(z), as

$$K_I(z)=\sqrt{2\pi }{D}_b(z){(ik)}^{1/2}\frac{c_{66}{\eta }^{\prime }}{\left(c_{12}+c_{66}\right)}\left\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }-{\eta }^{\prime }\right)+c_{12}\left({\xi }^{\prime }+{\eta }^{\prime }\right)\right\}{A}_1,$$

(60)

the components of inplane displacements and stresses, in the vicinity of a semi-infinite crack front, under symmetric far-field loading, can be expressed as follows:

$$u(r,\theta ,z)=\frac{K_I(z)}{\left(c_{11}{c}_{22}-c_{12}^2\right){\eta }^{\prime }}\sqrt{\frac{r}{2\pi }}{[\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }+{\eta }^{\prime }\right)$$

$$+c_{12}\left({\xi }^{\prime }-{\eta }^{\prime }\right)\}\cos \left({\psi }_1/2\right)-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }-{\eta }^{\prime }\right)$$

$$+c_{12}\left({\xi }^{\prime }+{\eta }^{\prime }\right)\}\cos \left({\psi }_1^{\prime }/2\right)],$$

(61a)

$$v(r,\theta ,z)=\frac{K_I(z)}{\left(c_{11}{c}_{22}-c_{12}^2\right){\eta }^{\prime }}\sqrt{\frac{r}{2\pi }}{[\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\{c_{12}$$

$$+\sqrt{c_{11}{c}_{22}}{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2\}\sin \left({\psi }_1/2\right)-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\{c_{12}$$

$$+\sqrt{c_{11}{c}_{22}}{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2\}\sin \left({\psi }_1^{\prime }/2\right)],$$

(61b)

$${\sigma }_x(r,\theta ,z)=\frac{K_I(z)}{2\sqrt{2\pi r}{\eta }^{\prime }}{[\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\left({\xi }^{\prime }+{\eta }^{\prime }\right)\cos \left({\psi }_1/2\right)$$

$$-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\left({\xi }^{\prime }-{\eta }^{\prime }\right)\cos \left({\psi }_1^{\prime }/2\right)],$$

(62a)

$${\sigma }_y(r,\theta ,z)=-\frac{K_I(z)}{2\sqrt{2\pi r}{\eta }^{\prime }}{[\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\left({\xi }^{\prime }-{\eta }^{\prime }\right)\cos \left({\psi }_1/2\right)$$

$$-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\left({\xi }^{\prime }+{\eta }^{\prime }\right)\cos \left({\psi }_1^{\prime }/2\right)],$$

(62b)

$${\tau }_{xy}(r,\theta ,z)=\frac{K_I(z)}{2\sqrt{2\pi r}{\eta }^{\prime }}{[\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\sin \left({\psi }_1/2\right)$$

$$-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\sin \left({\psi }_1^{\prime }/2\right)$$

(62c)

4.2.2. Skew-symmetric (Mode II) Loading (Sliding Shear/Twisting)

For s = 1/2,

A1 = A3 = 0;

(63)

and

$$\frac{A_4}{A_2}=-\frac{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}\frac{\left\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }-{\eta }^{\prime }\right)+c_{12}\left({\xi }^{\prime }+{\eta }^{\prime }\right)\right\}}{\left\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }+{\eta }^{\prime }\right)+c_{12}\left({\xi }^{\prime }-{\eta }^{\prime }\right)\right\}}.$$

(64)

Finally, on substitution of Eqs. (63) and (64) into the expressions for displacements and stresses, and by defining the mode II stress intensity factor, K_II(z), as

$$K_{II}(z)=\sqrt{2\pi }{D}_b(z){(ik)}^{1/2}\frac{c_{66}{\eta }^{\prime }}{\left(c_{12}+c_{66}\right)\left({\xi }^{\prime }+{\eta }^{\prime }\right)}\left\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }-{\eta }^{\prime }\right)+c_{12}\left({\xi }^{\prime }+{\eta }^{\prime }\right)\right\}{A}_2,$$

(65)

the components of in-plane displacements and stresses, in the vicinity of a semi-infinite crack front, under skew-symmetric far-field loading, can be expressed as follows:

$$u(r,\theta ,z)=-\frac{K_{II}(z)}{\left(c_{11}{c}_{22}-c_{12}^2\right){\eta }^{\prime }}\sqrt{\frac{r}{2\pi }}{[\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\{\sqrt{c_{11}{c}_{22}}{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2$$

$$+c_{12}\}\sin \left({\psi }_1/2\right)-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\left\{\sqrt{c_{11}{c}_{22}}{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2+c_{12}\right\}\sin \left({\psi }_1^{\prime }/2\right)],$$

(66a)

$$v(r,\theta ,z)=\frac{K_{II}(z)}{\left(c_{11}{c}_{22}-c_{12}^2\right){\eta }^{\prime }}\sqrt{\frac{r}{2\pi }}{[\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }-{\eta }^{\prime }\right)$$

$$+c_{12}\left({\xi }^{\prime }+{\eta }^{\prime }\right)\}\cos \left({\psi }_1/2\right)-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\{\{\sqrt{c_{11}{c}_{22}}\left({\xi }^{\prime }+{\eta }^{\prime }\right)$$

$$+c_{12}\left({\xi }^{\prime }-{\eta }^{\prime }\right)\}\cos \left({\psi }_1^{\prime }/2\right)],$$

(66b)

$${\sigma }_x(r,\theta ,z)=-\frac{K_{II}(z)}{2\sqrt{2\pi r}{\eta }^{\prime }}[{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}{({\xi }^{\prime }+{\eta }^{\prime })}^2\sin \left({\psi }_1/2\right)$$

$$-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2\sin \left({\psi }_1^{\prime }/2\right)],$$

(67a)

$${\sigma }_y(r,\theta ,z)=\frac{K_{II}(z)}{2\sqrt{2\pi r}{\eta }^{\prime }}{[\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\sin \left({\psi }_1/2\right)$$

$$-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\sin \left({\psi }_1^{\prime }/2\right)],$$

(67b)

$${\tau }_{xy}(r,\theta ,z)=\frac{K_{II}(z)}{2\sqrt{2\pi r}{\eta }^{\prime }}{[\left\{{\cos }^2\theta +{\left({\xi }^{\prime }+{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\left({\xi }^{\prime }+{\eta }^{\prime }\right)\cos \left({\psi }_1/2\right)$$

$$-{\left\{{\cos }^2\theta +{\left({\xi }^{\prime }-{\eta }^{\prime }\right)}^2{\sin }^2\theta \right\}}^{1/4}\left({\xi }^{\prime }-{\eta }^{\prime }\right)\cos \left({\psi }_1^{\prime }/2\right)].$$

(67c)

It may be noted that in the expressions for displacements and stresses, given by Eqs. (61) and (62), respectively, under mode I loading, there is no coupling between cos(ψ/2) and sin(ψ/2) (and similar absence of coupling between cos(ψ’/2) and sin(ψ’/2)) terms. The same is true in the expressions for displacements and stresses, given by (66) and (67), respectively, under mode II loading. It may be remarked that for $\chi >\sqrt{c_{11}/c_{22}}$, the through-thickness crack considered here will more likely propagate in its plane (010), i.e., x-z plane, and in the x-direction, because of this absence of coupling, which arises due to the imaginary roots. Section 7.2 below presents a comparison of the solution involving imaginary roots with its isotropic counterpart and concludes that the expressions for ρ and ψ are similar to their isotropic counterparts. This is especially true for ψ at x = 0. The significance of this, albeit being qualitatively and indirectly inferred, is that for $\chi >\sqrt{c_{11}/c_{22}}$ (or λ >1), the through (010)[001] crack considered here is, like its isotropic counterpart, more likely to propagate in its original plane (010) and direction [100], because of this absence of coupling, which arises due to the imaginary roots. It can then be inferred that the (010) or x-z plane is the easy cleavage plane and x-direction is the easy propagation direction in orthorhombic single crystals with $\chi >\sqrt{c_{11}/c_{22}}$. As has been remarked earlier, determination of a cleavage system being easy or difficult for crack propagation cannot definitively be guaranteed by the Griffith-Irwin based approach, which can serve only as a necessary condition. A similarity-based approach is needed to introduce a sufficient condition (see Section 7.2 below).

Finally, it may be noted that the above expressions for displacements, given by Eqs. (37), (43), (54), (66) and stresses, given by (38), (44), (55), (67), reduce to their two-dimensional counterparts (see, e.g., Sih et al. [17]).

5. Plate Surface Boundary Conditions and Through-Thickness Distribution of Singular Stress Fields

5.1. Satisfaction of traction-free boundary conditions

The stress field in the vicinity of the front of a semi-infinite crack under inplane extension can be recovered if in Eqs. (36), (60) or (42), (65)

$$D_b\left(z\right)=D_{bs}\left(z\right)=D_2\cos \left(kz\right).$$

(68)

By using the boundary condition on the free plate surface, the general form of D_bs(z) can be obtained as

$$D_{bs}\left(z\right)={\displaystyle \sum _{n=0}^{\infty }}{D}_{2n}\cos \left(\frac{\left(2n+1\right)}{2h}\pi z\right).$$

(69)

Hence, K_I = K_Is and K_II = K_IIs would represent symmetric stress intensity factors; see Section 8 and Figure 5(a). If the odd functions are selected from D_b(z), it can yield the out-of-plane bending case given by

$$D_b\left(z\right)=D_{ba}\left(z\right)=D_1\sin \left(kz\right).$$

(70)

D_ba(z) that satisfies the stress-free condition on the plate surface is given by

$$D_{ba}\left(z\right)={\displaystyle \sum _{n=1}^{\infty }}{D}_{1n}\sin \left(\frac{n\pi }{h}z\right),$$

(71)

provided the loading function vanishes at z* = 0, thus eliminating the possibility of discontinuity of the function at z* = 0. In the presence of discontinuity of the function at z* = 0, $D_{ba}\left(z\right)$ can be written as follows:

$$D_{ba}\left(z*\right)=\left|{\displaystyle \sum _{n=0}^{\pm \infty }}{D}_{2n}\cos \left(\frac{\left(2n+1\right)}{2}\pi z*\right)\right|.$$

(72)

As a consequence, K_I = K_Ia and K_II = K_IIa would represent anti-symmetric stress intensity factors; see Section 8 and Figure 5b). If the boundary conditions on the free plate surfaces are satisfied, all the displacements and singular stresses vanish on the plate surfaces in the vicinity of the crack front.

5.2. Hyperbolic Cosine Distributed Far-Field Loading

Hyperbolic cosine distributed far-field loading, which is proportional to cosh(z*), |z*| < 1, is applied. The applied symmetric loading function and the corresponding "stress intensity factors" (valid for |z*| ≤ 1) are proportional to

$$D_{bs}\left(z^{*}\right)=\cosh (z^{*})=\frac{\exp (z^{*})+\exp (-z^{*})}{2}.$$

(73)

The corresponding Fourier series can be derived as follows:

$$D_{bs}\left(z^{*}\right)={\displaystyle \sum _{m=0}^{\infty }\frac{\left\{e^{-1}{\left(-1\right)}^m+e{\left(-1\right)}^m\right\}}{\left\{1+{\left(m+\frac{1}{2}\right)}^2{\pi }^2\right\}}\left(m+\frac{1}{2}\right)\pi \cos \left(\left(\frac{2m+1}{2}\right)\pi z^{*}\right)}.$$

(74)

The applied antisymmetric loading function (valid for |z*| < 1) and the corresponding "stress intensity factors" (valid for |z*| ≤ 1) are proportional to

$$D_{ba}\left(z^{*}\right)=\left|\cosh (-z)\right|=\frac{1}{2}\left|\exp (z^{*})+\exp (-z^{*})\right|.$$

(75)

Since $D_{ba}\left(z^{*}\right)$ has a discontinuity at z* = 0, the corresponding Fourier series can be obtained as given below:

$$D_{bs}\left(z^{*}\right)=\left|{\displaystyle \sum _{m=0}^{\infty }\frac{\left\{e^{-1}{\left(-1\right)}^m+e{\left(-1\right)}^m\right\}}{\left\{1+{\left(m+\frac{1}{2}\right)}^2{\pi }^2\right\}}\left(m+\frac{1}{2}\right)\pi \cos \left(\left(\frac{2m+1}{2}\right)\pi z^{*}\right)}\right|.$$

(76)

6. Stress Intensity Factors and Energy Release Rates for a Through-Thickness Center-Crack (Modes I and II)

6.1. Through-Thickness Distribution of Stress Intensity Factors (Modes I and II)

The stress intensity factors, $K_I(z)$ and $K_{II}(z)$, cannot be determined unless the far- field loading and a characteristic length (e.g., crack geometry) are specified. Sih et al. [17] have shown the applicability of the complex variable approach in conjunction with the eigenfunction expansion approach in the derivation of the two-dimensional stress intensity factors for anisotropic plates. The stress intensity factor for an infinite orthorhombic/tetragonal monocrystalline plate with a central crack of length, 2a, and subjected to far-field mode I/II loading is available for the two-dimensional case [17], and can easily be extended to the present three-dimensional case as follows:

$$K_I\left(z^{*}\right)+\frac{K_{II}\left(z^{*}\right)}{p_3}=2\sqrt{2\pi }\frac{\left(p_1-p_3\right)}{p_3}Li{m}_{\zeta \to {\zeta }_0}{({\zeta }_1-{\zeta }_0)}^{1/2}{\varphi }_1^{\prime }\left({\zeta }_1\right)D_b\left(z^{*}\right),$$

(77)

where

$${\zeta }_1=a+x+p_{1}y,\hspace{1em}{\zeta }_0=a,\hspace{1em}{\varphi }_1\left({\zeta }_1\right)=\frac{a^2}{4\left(p_1-p_3\right)}\left[\frac{2p_3{\sigma }^{\infty }+{\tau }^{\infty }}{{\zeta }_1+\sqrt{{\zeta }_1^2-{\zeta }_0^2}}\right]+C_1{\zeta }_1,$$

(78)

with $C_1$ being a constant. In cylindrical polar coordinates, Eq. (17) can be expressed as follows:

$$K_I\left(z^{*}\right)+\frac{K_{II}\left(z^{*}\right)}{p_3}=\sqrt{2\pi }{a}^2\left({\sigma }_y^{\infty }+\frac{{\tau }_{xy}^{\infty }}{p_3}\right)Li{m}_{r\to 0}\frac{\sqrt{r}}{\sqrt{{\left(a+r\right)}^2-a^2}\left\{a+r+\sqrt{{\left(a+r\right)}^2-a^2}\right\}}{D}_b\left(z^{*}\right),$$

(79)

which finally gives

$$K_I\left(z^{*}\right)={\sigma }_y^{\infty }\sqrt{\pi a}{D}_b\left(z^{*}\right),\hspace{1em}{K}_{II}\left(z^{*}\right)={\tau }_{xy}^{\infty }\sqrt{\pi a}{D}_b\left(z^{*}\right).$$

(80a,b)

for both complex and imaginary roots. Eqs. (20a,b) reduce to their two-dimensional counterparts [53], by taking $D_b\left(z*\right)=1$. It may further be noted that the normalization factor, $K_i\left(z\right)/K_{i,2D}$, I = I, II, is equal to $D_b\left(z*\right)$ for a given far-field loading.

6.2. Through-Thickness Distribution of Energy Release Rates (Modes I and II)

The through-thickness distributions of the energy release rates due far field loadings, ${\sigma }_y^{\infty }$ and ${\tau }_{xy}^{\infty }$, for a center-crack of length 2a, weakening an infinite plate of finite thickness, 2h, can be derived by introducing the thickness-wise partial crack closure method as follows:

$$G_I\left(z^{*}\right)=Li{m}_{\Delta a\to 0}\frac{1}{\left(\Delta a\right)}\left[{\displaystyle \underset{0}{\overset{\Delta a}{\int }}}{\sigma }_y(x,0,z^{*})v(\Delta a-x,\pi ,z^{*})dxd{z}^{*}\right],$$

(81)

which, on substitution of ${\sigma }_y(x,0,z^{*})={{\sigma }_y|}_{\theta =0}$ and $v(\Delta a-x,\pi ,z^{*})={|}_{\theta =\pi }$, obtained from Eqs. (38b) and (37b), respectively, for complex roots, and Eqs. (62b) and (61b), respectively, for imaginary roots, yields

$$G_I\left(z^{*}\right)=Li{m}_{\Delta a\to 0}\frac{1}{\left(\Delta a\right)}\{\begin{array}{cc}\frac{\sqrt{c_{11}{c}_{22}}}{\left(c_{11}{c}_{22}-c_{12}^2\right)}\eta K_I\left(\Delta a,z^{*}\right)K_I\left(0,z^{*}\right){\displaystyle \underset{0}{\overset{\Delta a}{\int }}\sqrt{\frac{\Delta a-x}{x}}dx},& for complex roots\\ \frac{\sqrt{c_{11}{c}_{22}}}{\left(c_{11}{c}_{22}-c_{12}^2\right)}{\xi }^{\prime }{K}_I\left(\Delta a,z^{*}\right)K_I\left(0,z^{*}\right){\displaystyle \underset{0}{\overset{\Delta a}{\int }}\sqrt{\frac{\Delta a-x}{x}}dx},& for imaginary roots\end{array}$$

$$=\{\begin{array}{cc}\frac{\sqrt{c_{11}{c}_{22}}}{\left(c_{11}{c}_{22}-c_{12}^2\right)}\eta {\left({\sigma }_y^{\infty }\right)}^2\pi a{\left(D_b\left(z^{*}\right)\right)}^2& for complex roots\\ \frac{\sqrt{c_{11}{c}_{22}}}{\left(c_{11}{c}_{22}-c_{12}^2\right)}{\xi }^{\prime }{\left({\sigma }_y^{\infty }\right)}^2\pi a{\left(D_b\left(z^{*}\right)\right)}^2,& for imaginary roots\end{array}$$

$$=\frac{{\left({\sigma }_y^{\infty }\right)}^2\pi a\sqrt{c_{11}{c}_{22}}}{\sqrt{2}\left(c_{11}{c}_{22}-c_{12}^2\right)}\sqrt{\sqrt{\left(c_{11}/c_{22}\right)}+\chi }{\left[D_b\left(z*\right)\right]}^2$$

$$=\frac{K_I^2\sqrt{c_{11}{c}_{22}}}{\sqrt{2}\left(c_{11}{c}_{22}-c_{12}^2\right)}\sqrt{\sqrt{\left(c_{11}/c_{22}\right)}+\chi }{\left[D_b\left(z*\right)\right]}^2.$$

(82)

Similarly,

$$G_{II}\left(z^{*}\right)=Li{m}_{\Delta a\to 0}\frac{1}{\left(\Delta a\right)}\left[{\displaystyle \underset{0}{\overset{\Delta a}{\int }}}{\tau }_{xy}(x,0,z^{*})u(\Delta a-x,\pi ,z^{*})dxd{z}^{*}\right],$$

(83)

which, on substitution of ${{\tau }_{xy}(x,0,z^{*})={\tau }_{xy}|}_{\theta =0}$ and $u(\Delta a-x,\pi ,z^{*})={u|}_{\theta =\pi }$, obtained from Eqs. (44c) and (43a), respectively, for complex roots, and Eqs. (67c) and (66a), respectively, for imaginary roots, yields the following:

$$G_{II}\left(z^{*}\right)=Li{m}_{\Delta a\to 0}\frac{1}{\left(\Delta a\right)}\{\begin{array}{cc}\frac{c_{22}}{\left(c_{11}{c}_{22}-c_{12}^2\right)}\eta K_{II}\left(\Delta a,z^{*}\right)K_{II}\left(0,z^{*}\right){\displaystyle \underset{0}{\overset{\Delta a}{\int }}\sqrt{\frac{\Delta a-x}{x}}dx},& for complex roots\\ \frac{c_{22}}{\left(c_{11}{c}_{22}-c_{12}^2\right)}{\xi }^{\prime }{K}_{II}\left(\Delta a,z^{*}\right)K_{II}\left(0,z^{*}\right){\displaystyle \underset{0}{\overset{\Delta a}{\int }}\sqrt{\frac{\Delta a-x}{x}}dx},& for imaginary roots\end{array}$$

$$=\{\begin{array}{cc}\frac{c_{22}}{\left(c_{11}{c}_{22}-c_{12}^2\right)}\eta {\left({\tau }_{xy}^{\infty }\right)}^2\pi a{\left(D_b\left(z^{*}\right)\right)}^2,& for complex roots\\ \frac{c_{22}}{\left(c_{11}{c}_{22}-c_{12}^2\right)}{\xi }^{\prime }{\left({\tau }_{xy}^{\infty }\right)}^2\pi a{\left(D_b\left(z^{*}\right)\right)}^2,& for imaginary roots\end{array}$$

$$=\frac{{\left({\tau }_{xy}^{\infty }\right)}^2\pi a{c}_{22}}{\sqrt{2}\left(c_{11}{c}_{22}-c_{12}^2\right)}\sqrt{\sqrt{\left(c_{11}/c_{22}\right)}+\chi }{\left[D_b\left(z*\right)\right]}^2$$

$$=\frac{K_{II}^2{c}_{22}}{\sqrt{2}\left(c_{11}{c}_{22}-c_{12}^2\right)}\sqrt{\sqrt{\left(c_{11}/c_{22}\right)}+\chi }{\left[D_b\left(z*\right)\right]}^2.$$

(84)

Eqs. (82, 84) reduce to their two-dimensional (plane strain) counterparts [53], by taking $D_b\left(z*\right)$ = 1. It may further be noted that the normalization factor is equal to ${\left(D_b\left(z*\right)\right)}^2$ for a given far-field loading.

For the special case of a tetragonal single crystal, the above energy release rates reduce to

$$G_I\left(z^{*}\right)=\frac{{\left({\sigma }_y^{\infty }\right)}^2\pi a{c}_{11}}{\sqrt{2}\left(c_{11}^2-c_{12}^2\right)}\sqrt{1+\chi }{\left[D_b\left(z*\right)\right]}^2, G_{II}\left(z^{*}\right)=\frac{{\left({\tau }_{xy}^{\infty }\right)}^2\pi a{c}_{11}}{\sqrt{2}\left(c_{11}^2-c_{12}^2\right)}\sqrt{1+\chi }{\left[D_b\left(z*\right)\right]}^2.$$

(85a,b)

7. Necessary and Sufficient Conditions for Easy or Difficult Cleavage Planes

7.1. Crack Deflection Criterion, based on the relative fracture energy

The important issue of a cleavage plane being deemed easy or difficult can be related to a crack deflection criterion, which is based on the relative fracture energy (or the energy release rate) available for possible “fracture paths” [17]. The deflection or kinking of a crack from the cleavage system 1 to the cleavage system 2 is favored if

$$\raisebox{1ex}{$G_1$}\!\left/ \!\raisebox{-1ex}{$\left(2{\mathsf{\Gamma }}_1\right)$}\right.<1<\raisebox{1ex}{$G_2$}\!\left/ \!\raisebox{-1ex}{$\left(2{\mathsf{\Gamma }}_2\right)$}\right.\Rightarrow \raisebox{1ex}{$G_2$}\!\left/ \!\raisebox{-1ex}{$G_1$}\right.>\raisebox{1ex}{${\mathsf{\Gamma }}_2$}\!\left/ \!\raisebox{-1ex}{${\mathsf{\Gamma }}_1$}\right.$$

(86)

in which G_i and Γ_i, i = 1, 2, are energy release rate and surface energy, respectively, of the i^th cleavage system. Atomistic scale modeling of cracks, however, requires consideration of both the long-range elastic interactions and the short-range chemical reactions. The Griffith theory does not take the latter into account [22]. Secondly and more importantly, fracture criteria derived from equilibrium theories such as the Griffith thermodynamics-based) criterion can only be regarded as necessary conditions but not as sufficient [79,83]. The effect of short-range chemical reactions can obviously be encapsulated by atomic scale simulations, such as the investigation of low-speed propagation instabilities in silicon using quantum-mechanical hybrid, multi-scale modelling due to Kermode et al. [84], which, however, entails extensive computational and other resources. Alternatively, and more importantly, such short-range interactions can also be captured by the elastic properties-based parameters (with a few exceptions), such as the inverse anisotropic ratio, λ, or equivalently, the normalized elastic parameter, χ. This is because the elastic properties are controlled by various aspects of the underlying structural chemistry of single crystals, such as the Bravais lattice type, bonding (covalent, ionic, and metallic), bonding (including hybridized) orbitals, electro-negativity of constituent atoms in a compound, polarity, etc. [22]. General theory behind these characteristics pertaining to structural chemistry of crystals are available in well-known treatises (see e.g., [85,86,87]). More specifically, the elastic properties of superconducting YBa₂Cu₃O_7-δ are strongly influenced by oxygen non-stoichiometry (as well as various structural defects). It is known to crystallize in a defect perovskite structure consisting of layers. When δ = 1, the O(1) sites in the Cu(1) layer are vacant and the structure is tetragonal. For δ < 0.65, Cu-O chains along the b-axis of the crystal are formed^‡. Elongation of the b-axis changes the structure from tetragonal (insulator) to orthorhombic (superconductor), with lattice parameters of a = 3.82 Å, b = 3.89 Å, and c = 11.68 Å. Optimum superconducting properties occur when δ ~ 0.07 and all of the O(1) sites are occupied with few vacancies^‡. The coordination geometry of metal centers in YBCO, such as cubic {YO₈}, {BaO₁₀}, square planar {CuO₄} and square pyramidal {CuO₅} as well as structural features such as puckered Cu plane and Cu ribbons were first reported by Williams et al. [88]. Furthermore, Ledbetter [89] and Lin et al. [90] measured the elastic constants of polycrystalline YBCO using ultrasonic methods and found that while the “elastic moduli corresponding mainly to shear modes increase monotonically with oxygen concentration”, their counterparts due to “dilation modes increase up to the values of 6.7 of the oxygen index, after which they begin to decrease”; see also Lubenets et al. [91].

7.2. Comparison of Solutions Involving Complex and Imaginary Roots with Their Isotropic Counterpart

As has been mentioned earlier in Section 4.1, similarity analysis is an effective tool to solve complex problems in fracture mechanics of single crystals [22,82]. In what follows, similarity or dissimilarity of the present asymptotic solutions involving complex and imaginary roots with their isotropic counterparts, which can lead to a sufficient condition for determination of a cleavage system being easy or difficult for crack

------------------------------------------------------------------------------------------------------------

^‡Wikipedia, 2010. http://en.wikipedia.org/wiki/Yttrium_barium_copper_oxide propagation, is presented.

7.2.1. Isotropic Materials:

The inplane displacements for an isotropic material can be rewritten in the form (for n = 0) [25,55]:

$$U(x,y,z)={\left(ik\right)}^s{a}_s\left(z\right){\rho }^s{e}^{ip\psi },\hspace{1em}V(x,y,z)={\left(ik\right)}^s{b}_s\left(z\right){\rho }^s{e}^{ip\psi },$$

(87a,b)

in which

$$p=\pm \left(s\pm 1\right),$$

(88)

$$\rho =\sqrt{x^2+y^2},$$

(89)

and

$$\psi ={\tan }^{-1}\left(\frac{y}{x}\right).$$

(90)

Therefore, for an isotropic material when x = 0, $\psi$ = π/2 for all positive values of y.

7.2.2. Solution Involving Complex Roots:

Going back to (8a,b) and (15), the inplane displacements can be rewritten in the form (for n = 0):

$$U(x,y,z)={\left(ik\right)}^s{a}_s\left(z\right){(x+py)}^s={\left(ik\right)}^s{a}_s\left(z\right){\rho }^s{e}^{is\psi },$$

(91a)

$$V(x,y,z)={\left(ik\right)}^s{b}_s\left(z\right){(x+py)}^s={\left(ik\right)}^s{b}_s\left(z\right){\rho }^s{e}^{is\psi },$$

(91b)

in which $\rho$ and $\psi$ can be rewritten as follows:

$$\rho =\sqrt{{\left(x\pm \xi y\right)}^2+{\eta }^2{y}^2},$$

(92)

$$\psi ={\tan }^{-1}\left(\frac{\pm \eta y}{x\pm \xi y}\right).$$

(93)

Therefore, for an orthorhombic (tetragonal and cubic being special cases) crystal with complex roots when x = 0,

$$\psi ={\tan }^{-1}\left(\frac{\pm \eta }{\pm \xi }\right)$$

(94)

for all positive values of y, which differs from its isotropic counterpart.

7.2.3. Solution Involving Imaginary Roots

Going back to (8a,b) and (45), the inplane displacements can be rewritten in the form (for n = 0):

$$U(x,y,z)={\left(ik\right)}^s{a}_s\left(z\right){(x+py)}^s={\left(ik\right)}^s{a}_s\left(z\right){\rho }^s{e}^{is\psi },$$

(95a)

$$V(x,y,z)={\left(ik\right)}^s{b}_s\left(z\right){(x+py)}^s={\left(ik\right)}^s{b}_s\left(z\right){\rho }^s{e}^{is\psi },$$

(95b)

in which $\rho$ and $\psi$ can be rewritten as follows:

$$\rho =\sqrt{x^2+{\left(\pm {\xi }^{\prime }\pm {\eta }^{\prime }\right)}^2{y}^2},\psi ={\tan }^{-1}\left(\frac{\left(\pm {\xi }^{\prime }\pm {\eta }^{\prime }\right)y}{x}\right).$$

(96)

Therefore, for an orthorhombic (tetragonal and cubic being special cases) crystal with imaginary roots when x = 0, $\psi$ = π/2 for all positive values of y, which is in accord with its isotropic counterpart.

8. Numerical Results and Discussion

Table 1(a) and (b) display the elastic stiffness constants of orthorhombic (superconductor) and tetragonal (insulator) YBCO single crystals. If otherwise not specified, the elastic stiffness constants are measured at the room temperature (Approx. 300^o K). Table 2 shows $\sqrt{c_{11}/c_{22}}$, the normalized elastic parameter, χ, nature of the four roots of characteristic equation (complex or imaginary), and the character of the cleavage plane (easy or not) for a (010)[001] through-thickness crack with [100] being initial propagation direction, while Table 3 exhibits their counterparts for a (010)[100] through-thickness crack with [001] being the initial propagation direction. Table 4, Table 5, Table 6 and Table 7 present similar results for (100)[001]×[010], (100)[010]×[001], (001)[100]×[010] and (001)[010]×[100] through-thickness crack systems.

Next, the effect of elastic constants, c_ij (especially, c₁₂ and to a lesser extent, c₆₆), on the nature (i.e., easy or difficult) of a cleavage system in YBCO (YBa₂Cu₃O_7−δ) is discussed. Only three complete sets of elastic constants are available in the literature accessible to the present author, out of which those due to Ledbetter and Lei [73] are just estimates (marked ** in Table 1(a)), while their experimental counterparts due to Reichard et al. [67] are based on the assumption of tetragonal symmetry; see Table 1(b). This only leaves the experimental measurements (by resonant ultrasound spectroscopy) due to Lei et al. [1], marked * in Table 1(a). However, their $c_{12}$ value appears to be excessively high. This is because, according to these authors themselves, "no wave speed in the crystal depends only on $c_{12}$, it is no way to estimate it directly." It also is well-known that while $c_{12}$ and $c_{66}$ can be measured independently by static tests [76], these constants are always coupled in vibrations-based measurements [77,78]. Therefore, in Table 1 of the present investigation, both $c_{12}$ and $c_{66}$, measured by ultrasound by Saint-Paul and Henry [71], have been utilized (marked ***) in replacement of their counterparts due to Lei et al. [1] in order to assess the fracture characteristics of YBCO, and compare them with experiments by Cook et al. [6] ], Raynes et al. [9] and Goyal et al. [10] among others.

Table 2 shows that the normalized elastic parameter, χ, for YBCO* (measurements reported by Lei et al. [1]) is smaller than $\sqrt{c_{11}/c_{22}}$, giving rise to complex roots (of the characteristic equation) for a (010)[001]×[100] through-crack, weakening a YBCO monocrystalline plate. Same is true for a ($\overline{1}$00)[001]×[010] crack shown in Table 4. These results predict that (010) and ($\overline{1}$00) are difficult cleavage planes, which are in contradiction with the experimental observations by Cook et al. [6], Raynes et al. [9] and Goyal et al. [10] among others. The reason behind this anomaly lies in the excessive values of $c_{12}$ (and $c_{66}$ to a lesser extent) used in the computation of λ and χ. Since the rest of the crack systems are predicted to constitute easy cleavage planes/directions (see Table 3, Table 5, Table 6 and Table 7), the remaining elastic constants measured by Lei et al. [1] are considered to be reasonably accurate. Table 3, Table 5, Table 6 and Table 7 further show that for YBCO** the cleavage systems, (0$\overline{1}$0)[100]×[001], (001)[100]×[010], (100)[010]×[001] and (001)[010]×[100], are predicted to be difficult, thus somewhat invalidating the values of the corresponding elastic constants estimated by Ledbetter and Lei [73]. Furthermore, Table 3, Table 4, Table 5, Table 6 and Table 7 show that for (tetragonal) YBCO^T, all the six cleavage systems investigated here are found to be difficult, thus completely invalidating the values of the corresponding experimentally determined elastic constants reported by Reichard et al. [67]. As can be seen from Table 2, Table 3, Table 4, Table 5, Table 6 and Table 7, only for YBCO***, all the cleavage systems are predicted to be easy, which is in agreement with the experimentally observed fracture characteristics of YBCO due to Cook et al. [6], Raynes et al. [9] and Goyal et al. [10] among others; see also Granozio and di Uccio [14] for a summary of the available experimental results. They [14] have also presented approximate theoretical results of fully oxidized YBCO’s (δ = 0, 1), and concluded that the three lowest surface energies follow the inequality: γ (001) < γ (100) < γ (010). Furthermore, based on experimental results from transmission electron microscopy [92], X-ray photo-emission microscopy [93], low-energy ion scattering spectroscopy [94], and surface polarity [95] analyses performed on fully oxidized YBCO crystals, these authors [14] have shown that the low energy cut is between the Ba=O and Cu=O planes.

Efficacy of the indentation test has extensively been studied in the brittle fracture literature [96,97,98]. Lawn [96] and Anstis et al. [97] have presented the following relationship between fracture toughness and size of a radial crack produced by a Vickers-type sharp indenter:

$$K_c={\chi }_{r}P{c}_0^{-3/2},$$

where ${\chi }_r={\S }_V^R{\left(E/H\right)}^{1/2}$, finally giving rise to the following:

$$K_c={\S }_V^R{\left(E/H\right)}^{1/2}P{c}_0^{-3/2},$$

(97)

in which P, c₀, E and H represent the indentation load, equilibrium half-crack length, Young’s modulus, and hardness (of an isotropic material), respectively, and ${\S }_V^R$ denotes a material-independent constant for the Vickers-produced radial crack. Raynes et al. [9] following the lead of Anstis et al. [97] have determined the fracture toughness of mono-crystalline YBa₂Cu₃O_7−δ, taking into account of its anisotropy. Table 8 presents the critical stress intensity factor or fracture toughness (K_c) and the critical energy release rate or fracture energy (G_c) of the six easy cleavage systems of monocrystalline superconducting YBCO. It is worthwhile to note here that there is some misconception about computation of fracture energy, G_c, from the corresponding measured value of K_c of an anisotropic (e.g., orthorhombic) single crystal in the literature; see e.g., Granozio and di Uccio [14]. The factor $\sqrt{\left(\sqrt{\left(c_{11}/c_{22}\right)}+\chi \right)/2}$ (see Eqn. (82) above) is not accounted for in these authors’ computations. The energy release rate in an anisotropic (e.g., orthorhombic) single crystal not only varies from one cleavage plane to another, but also varies according to propagation direction.

Figure 5(a, b) show that a (010)[001] crack initially propagating in [100] direction would turn into a (100)[001] crack propagating in [010] direction (also in [010] direction because of symmetry), while a (100)[001] crack initially propagating in [010] direction would continue in its original track. This is because, G_c of a (010)[001]×[100] cleavage system, 1.50177 J/m², is higher than its (100)[001]×[010] counterpart, 1.02649 J/m². In a similar vein, as shown in Figure 6(a, b), a (100)[100] crack initially propagating in [001] direction would turn into a (001)[100] crack propagating in [010] direction (also in [010] direction because of symmetry), while a (001)[100] crack initially propagating in [010] direction would continue uninterrupted in its original track. This is because, G_c of a (010)[100]×[001] cleavage system, 1.91945 J/m², is about 2.6 times its (001)[100]×[010] counterpart, 0.73912 J/m². Likewise, as shown in Figure 7(a, b), a (100)[010] crack initially propagating in [001] direction would turn into a (001) [010] crack propagating in [100] direction (also in [100] direction because of symmetry), while a (001)[010] crack initially propagating in [100] direction would continue unhindered in its original track. This is because, G_c of a (100)[010]×[001] cleavage system, 1.02649 J/m², is more than 1.4 times its (001)[010]×[100] counterpart, 0.71163 J/m². Finally, in mono-crystalline superconducting YBCO under triaxial-tension far-field loading, a (010)[100] crack initially propagating in [001] direction would eventually turn into a c-plane cleavage fracture, as shown in Figure 8.

Figure 9(a, b) show variation of the normalized stress intensity factors, $K^{*}\left(z\right)=K\left(z\right)/K_{Plane Strain}$, through the thickness of an orthorhombic mono-crystalline plate, weakened by a through-thickness crack. Variation of the normalized stress intensity factor, $K^{*}\left(z\right)$, through the thickness of the same plate, weakened by any of the six through-cracks investigated here is identical. Figure 9(a) shows the through-thickness variation of $K_S^{*}\left(z\right)=K_S\left(z\right)/K_{Plane Strain}$ for a far-field symmetrically distributed hyperbolic cosine load for mode I (stretching) or mode II (inplane shear), while its skew-symmetric counterpart $K_A^{*}\left(z\right)=K_A\left(z\right)/K_{PlaneStrain}$ for mode I (bending) or mode II (twisting) is displayed in Figure 8(b). Of special significance is the discontinuity in the stress intensity factor at z* = 0 in the skew-symmetric loading case, shown in Figure 9(b). Figure 10 shows the corresponding variation of energy release rate, $G*$, through the top half of the plate thickness. For through-thickness symmetric far-field loading, the crack is expected to grow through thickness in a stable manner till the stress intensity factor or the energy release rate reaches its critical value at the mid-thickness. With further increase of the magnitude of the far-field loading, unstable crack growth is expected to progressively spread throughout the plate thickness. For skew-symmetric loading, as reported on earlier occasions [55], the bottom half will experience crack closure. Such types of results describing the three-dimensional distribution of stress intensity factors and energy release rates are generally unavailable in the fracture mechanics literature.

9. Summary and Conclusions

A modified eigenfunction expansion technique, based partly on separation of the z-variable and in part, on the Eshelby [60]-Stroh [15] type affine transformation, is employed to derive three-dimensional asymptotic displacement and stress fields in the vicinity of the front of a semi-infinite through-thickness crack weakening an infinite orthorhombic single crystal plate. Crack-face boundary conditions and those that are prescribed on the top and bottom (free) surfaces of the orthorhombic plate are exactly satisfied. Explicit expressions for the singular stresses in the vicinity of the front of the through-thickness crack, subjected to far-field inplane mode I and II loadings, are presented.

The present investigation considers six through-crack systems ― (010)[001] with the [100] length direction, (0$\overline{1}$0)[100] with the [001] length direction, ($\overline{1}$00)[001] with the [010] length direction, (100)[010] with the [001] length direction, (001)[0$\overline{1}$0] with the [100] length direction, (001)[100] with the [010] length direction ― weakening orthorhombic YBCO single crystal plates. More important, the present approach predicts whether a crack would propagate in its original plane/direction or deflect to a different one. The present study is unique in the sense that such a fracture mechanics criterion is employed for accurate determination of the full set of elastic constants of mono-crystalline YBCO.

The following interesting conclusions can be drawn the present investigation:

• Atomistic scale modeling of cracks requires consideration of both the long-range elastic interactions and the short-range chemical reactions. The Griffith thermodynamic-based theory does not take the latter into account, and hence must be regarded as a necessary condition but not as sufficient.
• The effect of short-range chemical reactions can be adequately captured by the elastic properties-based parameters, such as the inverse anisotropic ratio, λ, or equivalently, the normalized elastic parameter, χ. This is because the elastic properties are controlled by various aspects of the underlying structural chemistry of single crystals, such as the Bravais lattice type, bonding (covalent, ionic, and metallic), bonding (including hybridized) orbitals, electro-negativity of constituent atoms in a compound, polarity, etc.
• More specifically, the elastic properties of superconducting YBa₂Cu₃O_7-σ are strongly influenced by oxygen non-stoichiometry (as well as various structural defects).
• Similarity or dissimilarity of the present asymptotic solutions involving complex (λ < 1 or equivalently, $\chi <\sqrt{c_{11}/c_{22}}$) and imaginary roots (λ > 1 or equivalently, $\chi >\sqrt{c_{11}/c_{22}}$) with their isotropic (λ = 1) counterparts can lead to a sufficient condition for determination of a cleavage system being easy or difficult for crack propagation.
• The normalized elastic parameter, χ, for YBCO* is smaller than $\sqrt{c_{11}/c_{22}}$, giving rise to complex roots (of the characteristic equation) for a (010)[001] [100] through- crack, weakening a YBCO monocrystalline plate. Same is true for a ($\overline{1}$00)[001]×[010] crack. These results predict that (010) and ($\overline{1}$00) are difficult cleavage planes, which are in contradiction with the experimental observations.
• Only for YBCO***, all the cleavage systems are predicted to be easy, which is in agreement with the experimentally observed fracture characteristics, thus ensuring that a reasonably accurate complete set of nine experimentally determined elastic constants has been arrived at, by employing the present theoretical approach.
• For tetragonal YBCO^T, all the six cleavage systems investigated here are found to be difficult, thus completely invalidating the values of the corresponding experimentally determined elastic constants reported by Reichard et al. [67].
• Finally, generally unavailable results, pertaining to the through-thickness variations of stress intensity factors and energy release rates for a crack corresponding to symmetric and skew-symmetric hyperbolic cosine loads that also satisfy the boundary conditions on the top and bottom surfaces of an orthorhombic monocrystalline plate under investigation, bridge a longstanding gap in the stress singularity/fracture mechanics literature.