1. Introduction

The concentration of stresses in the vicinity of multiple interacting cracks represents one of the most common fracture sources in structural elements. A notable example is the stress analysis of structural composites containing microcracks. The microcracking process has been used often as a mechanism to enhance the fracture toughness of materials. Microcracks can either enhance or shield a macrocrack depending on their positions and orientations. Also the knowledge about crack interaction is of fundamental importance in predicting the effective moduli of cracked solids. Hence, the subject of interaction effect between multiple cracks in fracture mechanics has attracted the attention of many researchers. Due to complexities in analysis of multiple cracks, researchers have often studied two-dimensional statement of the problem, when cracks have infinite length in one direction (tunnel form) ^[1,2,3]. The research of three-dimensional crack interactions has been in progress for many years, and a variety of methods have been developed to tackle this problem ^[4]. As a powerful analysis tool for linear elasticity problems the boundary integral equation method has been widely employed in the study of three-dimensional crack interactions ^[5]. The main advantages of this method are its more direct approach in formulating the problems and reduction the problem dimension to the two-dimensional work domain of boundary integrals. A closed form of boundary integral equations solutions obtained by analytical methods, is limited only to three-dimensional problems with simple configuration ^[6,7,8]. However, it is difficult to analytically solve the problem for the arbitrary array of cracks, and hence different numerical approaches have been developed. Among these numerical methods, e.g., the collocation method ^[6] or the boundary-elementmethod ^[9,10] are outstanding and have demonstrated their flexibility and versatility in many crack-related problems. But, computational difficulty in the above numerical methods exists in the case of a large number of interacting cracks. Since in these methods of computing the costs are mostly spent on solving the systems of linear algebraic equations, to which the system of boundary integral equation in the discretization process is reduced. Hence, in a problem with a large number of cracks, the number of computing operations increases significantly, since it is generally proportional to the square of the number of cracks. In order to overcome the above computational difficulties encountered in the analysis of crack interactions, a method of effective stress field is developed ^[11,12]. By decomposing the problem into several sub-problems and treating the traction on the surface of crack in each sub-problem as an unknown function, the stress intensity factors for problems containing irregularly distributed cracks can be calculated, because each sub-problem is reduced to solution of a system of only three integral equations.

Of particular interest is the consideration of three-dimensional problems of elasticity theory for limited bodies of complex shape with internal interacting cracks. Hitherto in the literature dominate the solutions of three-dimensional problems for interaction of cracks with simple topological forms of boundary surface of structural elements ^[13,14]. Mostly these are flat or cylindrical surfaces. The boundary element method proves to be especially efficient at the same time to satisfy of the boundary conditions on the body surface and on surfaces of cracks. In this paper we use the method of boundary integral equations to obtain approximate values of the SIF on the cracks contours near the boundary of the body of complex shape under static loading. Application of the proposed effective stress field method in the multiple interaction cracks problems for the limited body makes it possible to take into account the complexity of the boundary geometry.

Here the boundary integral equation method is modified and applied to analyze three-dimensional elastostatic problem for a limited body containing interacting penny-shaped cracks near the boundary surface. The problem is reduced to a system of coupled boundary integral equations consisting of three displacement integral equations on the boundary surface and 3 N traction integral equations on the crack-surface, where N is the number of interacting cracks. This modification allows the direct determination of the stress intensity factors on the crack contour. The regularization procedures for the boundary integral equations with different order singularities rest on the implementation of singularity subtraction technique and the change of variables technique. The discrete analogue of boundary integral equations as a system of linear algebraic equations is obtained by means of collocation and boundary element scheme. Numerical calculations are carried out for the problem with uniform tensile loading along the prismatic rod axis, two penny-shaped cracks are perpendicular to its axis and lying on the plane of middle cross-section of the rod. The effects of the cracks location and distance relative to the boundary surface of the rod on the mode-Ⅰ SIF on cracks contour are investigated.

---

2. Materials and Method

---

2.1. Direct Method Investigation of Interacting Cracks in an Infinite Solid

Let's consider a three-dimensional infinite isotropic solid which contains N flat arbitrarily oriented cracks (Figure 1). Assume that the ${{n}^{th}}$ crack occupies a region ${{S}_{n}}\,\,\left( n=\overline{1,N} \right)$ and is bounded by a smooth contour.

The opposite crack faces are not in contact and are load-free. At the centers of all cracks, we choose local coordinate systems ${{O}_{k}}\,{{x}_{1k}}\,{{x}_{2k}}\,{{x}_{3k}}$ such that ${{k}^{th}}$ crack lies in the plane ${{x}_{3k}}=0$. The matrix material is specified by the shear modulus G and Poisson's ratio v. The strain and stress state in the cracked material is caused by the remote loading $\left\{ \sigma _{ij}^{0} \right\}\,\,\left( i,j=\overline{1,3} \right)$ which is characterized by the displacement components $u_{1}^{0},\,u_{2}^{0},\,u_{3}^{0}$ in a certain principal coordinate system $O\,{{x}_{10}}\,{{x}_{20}}\,{{x}_{30}}$ in the continuous solid under the assumption of crack absence. We characterize the relative location of cracks in the body by the following parameters (Figure 1): ${{d}_{nk}}$ is the distance between the centers of the ${{n}^{th}}$ and ${{k}^{th}}$ defects, ${{e}_{jnk}}$ $j=\overline{1,3}$ are the direction cosines of the vector ${{d}_{nk}}$, and ${{l}_{ijnk}}$, $\left( n,k=\overline{0,N} \right)$ are the cosines of the angles between ${{i}^{th}}$ and ${{j}^{th}}$ axes of the ${{k}^{th}}$ and ${{n}^{th}}$ coordinate systems, respectively.

After representation the solutions of the Lamé equilibrium equations in the form of the combination of elastic potentials of simple and double layers, satisfaction the boundary conditions for all cracks, and taking into account the principle of superposition, we obtain a system of $3N$ boundary integral equations ^[8] in the form

$\sum\limits_{i=1}^{3}{\sum\limits_{p=1}^{3}{\sum\limits_{k=1}^{N}{\iint\limits_{{{S}_{k}}}{\,{{K}_{ijpkm}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{km}} \right){{l}_{3pkm}}\,{{l}_{ijkm}}\,\Delta {{u}_{pk}}\left( \mathbf{\xi } \right){{d}_{\mathbf{\xi }}}{{S}_{k}}}=\frac{1-\nu }{G}{{N}_{jm}}\left( {{\mathbf{x}}_{mm}} \right)\,\,\,\, \\ i,j=\overline{1,3}\,,\,\,\,k,m=\overline{1,N}}}}\,,$

(1)

Here, ${{\mathbf{x}}_{km}}$ and ${{\mathbf{x}}_{mm}}$ are the same points of view on the surface ${{S}_{m}}$ in ${{k}^{th}}$ and ${{m}^{th}}$ local coordinate systems, respectively; the relation between the coordinates of this point in different local coordinate systems will have the form: ${{x}_{jkm}}={{e}_{jkm}}\,{{d}_{km}}+\sum\limits_{q=1}^{2}{{{l}_{qjkm}}\,{{x}_{qmm}}}$; unknown functions $\Delta {{u}_{jn}}\left( {{\mathbf{x}}_{nn}} \right)\frac{u_{jn}^{-}\left( {{\mathbf{x}}_{nn}} \right)-u_{jn}^{+}\left( {{\mathbf{x}}_{nn}} \right)}{4\pi }$, characterizing ${{n}^{th}}$ crack opening in the direction of axes ${{O}_{n}}\,{{x}_{jn}}$; ${{N}_{jm}}\left( {{x}_{mm}} \right)=\sum\limits_{i=1}^{3}{\sum\limits_{p=1}^{3}{\sigma _{ip}^{0}\,{{l}_{3pm0}}\,{{l}_{ijm0}}}}$ are the stress functions in the plane that coincides with the surface of the ${{n}^{th}}$ crack; kernels of integral equations (1) have the form

$\begin{align} & K_{ijpkm}^{{}}\,\left( \mathbf{\xi },\,{{\mathbf{x}}_{km}} \right)=-\frac{G}{4\,\pi \,\left( 1-\nu \right)\,r_{km}^{3}}\sum\limits_{s=1}^{3}{\left\{ 3\,{{r}_{km,s}} \right.{{n}_{kms}}} \\ & \left[ \left( 1-2\,\nu \right)\,{{\delta }_{ij}}{{r}_{km,p}} \right.+\nu \left( {{\delta }_{jp}}{{r}_{km,i}}+ \right.\left. {{\delta }_{ip}}{{r}_{km,j}} \right)-\left. 5\,{{r}_{km,i}}\,{{r}_{km,j}}\,{{r}_{km,p}} \right]+ \\ & 3\,\left[ \nu \,({{n}_{kmi}}\,{{r}_{km,j}}\,{{r}_{km,p}}+\,{{n}_{kmj}}\,{{r}_{km,i}}\,{{r}_{km,p}})+\left. \left( 1-2\,\nu \right)\,{{n}_{kmp}}\,{{r}_{km,i}}\,{{r}_{km,j}} \right] \right. \\ & +\left( 1-2\,\nu \right)\left. \left[ {{\delta }_{jp}}\,{{n}_{kmi}}+{{\delta }_{ip}}\,{{n}_{kmj}} \right]-\left( 1-4\,\nu \right)\,{{\delta }_{ij}}\,v\,{{n}_{kmp}} \right\}. \\ \end{align}$

(2)

In the Eq. (2) the notations ${{r}_{km}}=\left| {{\mathbf{x}}_{km}}-\mathbf{\xi } \right|$, ${{r}_{km,j}}={({{x}_{kmj}}-{{\xi }_{j}})}/{{{r}_{km}}}\;$ are used, ${{\mathbf{n}}_{km}}$, with the components $n_{kmj}^{{}}$, is the normal vector to ${{m}^{th}}$ the crack-surface in ${{k}^{th}}$ local coordinate system, ${{\delta }_{ij}}$ is Kronecker symbol.

Each equation of the system (1) contains a hypersingular integral provided that $k=m$. These integrals should be regularized. The collocation method of numerical solution is reduced to the discretization of crack surfaces and determination of the unknown functions of jumps of displacements at the nodes of the region of integration. Finally, the problem is reduced to solution of a linear algebraic system of equations ^[6]. Dimension of this system is $3\,N\,M$, where $M$ is the number of collocation nodes. For multiple interaction crack problems, it is necessary to have a simplified procedure for determination the stress intensity factors on the contour of a given crack with regard for the influence of neighboring cracks.

---

2.2. Effective Stress Field Method for Investigation of Interacting Cracks in an Infinite Solid

Consider a crack, with occupied region ${{S}_{q}}$. In the vicinity of its contour we determine the degree of stress concentration. In the system of boundary integral equations (1), we select three equations obtained by satisfying the boundary conditions on this crack surface. We isolate hypersingular integrals, whose densities are functions characterizing the jump of displacements of the opposite surfaces of the ${{q}^{th}}$ crack on the left side ofthese equations. Next we transfer all other regular integrals, which have a physical meaning of stresses induced on the ${{q}^{th}}$ crack surfaces by the crack opening displacements of neighboring cracks to the right-hand side of Eq. (1). If opening of neighboring cracks considered are known, we obtain a system of three integral equations for the functions of jumps of displacements of the isolated crack whose surfaces are subjected to the action of fictitious forces $N_{jq}^{*}$:

$\begin{align} & \sum\limits_{i=1}^{3}{\iint\limits_{{{S}_{q}}}{\Delta {{u}_{iq}}\left( \mathbf{\xi } \right)\,{{K}_{ij3qq}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{qq}} \right)\,{{d}_{\mathbf{\xi }}}{{S}_{q}}}}= \\ & \frac{1-\nu }{G}N_{jq}^{*}\left( {{\mathbf{x}}_{qq}} \right)=\frac{1-\nu }{G}{{N}_{jq}}\left( {{\mathbf{x}}_{qq}} \right)\,- \\ & \sum\limits_{i=1}^{3}{\sum\limits_{p=1}^{3}{\sum\limits_{\begin{smallmatrix} k=1 \\ k\ne q \end{smallmatrix}}^{N}{\iint\limits_{{{S}_{k}}}{{{K}_{ijpkq}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{kq}} \right){{l}_{3pkq}}\,{{l}_{ijkq}}\,\Delta {{u}_{pk}}\left( \mathbf{\xi } \right)\,{{d}_{\mathbf{\xi }}}{{S}_{k}}}}}},\,\,\, \\ & i,j=\overline{1,3}\,,\,\,\,k=\overline{1,N}\, \\ \end{align}$

(3)

The kernels of hypersingular integrals can be simplified to the form ^[12]: (given that $n_{qqj}^{{}}={{\delta }_{3j}}$)

$\begin{align} & {{K}_{ij3qq}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{qq}} \right)=\frac{1+\nu \left( {{\delta }_{1i}}+{{\delta }_{2i}} \right)}{{{\left| {{\mathbf{x}}_{qq}}-\mathbf{\xi } \right|}^{3}}}{{\delta }_{ij}}+{{\left( -1 \right)}^{{{\delta }_{ij}}}} \\ & \frac{3\,\nu \,{{\left( {{x}_{qq1}}-{{\xi }_{1}} \right)}^{\left( {{\delta }_{i2}}+{{\delta }_{j2}} \right)}}{{\left( {{x}_{qq2}}-{{\xi }_{2}} \right)}^{\left( {{\delta }_{i1}}+{{\delta }_{j1}} \right)}}}{{{\left| {{\mathbf{x}}_{qq}}-\mathbf{\xi } \right|}^{5}}}\,\left( 1-{{\delta }_{i3}} \right)\left( 1-{{\delta }_{3j}} \right). \\ \end{align}$

For the approximate solution of the problem, we propose to replace the unknown functions of the crack opening on the right side of Eq. (3) by the corresponding functions for isolated cracks whose surfaces are subjected to the action of given loads. Therefore, we have the possibility of independent determination of the effective stress fields created by each neighboring crack separately. And the problem does not become more complicated in the case of increasing the number of interacting cracks because we actually find the numerical solutions only of a system of three integral equations for the functions of crack opening. For the numerical realization of this procedure, it is necessary to have closed solutions for isolated cracks whose surfaces are subjected to the action of a given external load.

---

2.3. Direct Method for Investigation of Interacting Cracks in a Limited Body

Let the homogeneous elastic body, which is limited by a smooth surface S, contain inside a system of N interacting cracks (Figure 2). The surface of the body S can be loaded by the known static surface forces $N_{0}^{{}}$. The crack surfaces are load-free.

The method of potentials defines the stress-strain state of the body under the external load in the following form ^[9]:

$\begin{align} & u_{i}^{{}}\left( {{\mathbf{x}}_{0}} \right)=\sum\limits_{j=1}^{3}{\iint\limits_{S}{\left[ U_{ij}^{{}}\left( \mathbf{\xi },{{\mathbf{x}}_{0}} \right)t_{j}^{{}}\left( \mathbf{\xi } \right)-T_{ij}^{{}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{0}} \right)u_{j}^{{}}\left( \mathbf{\xi } \right) \right]{{d}_{\mathbf{\xi }}}S}}+ \\ & \sum\limits_{k=1}^{N}{\sum\limits_{s=1}^{3}{\sum\limits_{j=1}^{3}{\iint\limits_{{{S}_{k}}}{T_{sjk0}^{{}}\left( \mathbf{\xi },{{\mathbf{x}}_{k0}} \right)l_{sik0}^{{}}\Delta u_{jk}^{{}}\left( \mathbf{\xi } \right){{d}_{\mathbf{\xi }}}{{S}_{k}}}}}} \\ \end{align}$

in principal coordinate system;

$\begin{align} & \sigma _{ijm}^{{}}\left( {{\mathbf{x}}_{0}} \right)=\sum\limits_{s=1}^{3}{\sum\limits_{w=1}^{3}{\sum\limits_{p=1}^{3}{\iint\limits_{S}{\left[ S_{swp0m}^{{}}\left( \mathbf{\xi },{{\mathbf{x}}_{0m}} \right){{l}_{sim0}}\,{{l}_{wjm0}}t_{p}^{{}}\left( \mathbf{\xi } \right)-K_{swp0m}^{{}}\left( \xi ,\,{{\mathbf{x}}_{0m}} \right){{l}_{sim0}}\,{{l}_{wjm0}}u_{p}^{{}}\left( \mathbf{\xi } \right) \right]d{{S}_{\mathbf{y}}}}}}}+ \\ & \sum\limits_{k=1}^{N}{\sum\limits_{s=1}^{3}{\sum\limits_{w=1}^{3}{\sum\limits_{p=1}^{3}{\iint\limits_{{{S}_{k}}}{K_{swpkm}^{{}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{km}} \right){{l}_{sikm}}\,{{l}_{wjkm}}\,\Delta u_{pk}^{{}}\left( \mathbf{\xi } \right){{d}_{\mathbf{\xi }}}{{S}_{k}}}}}}} \\ \end{align}$ in ${{m}^{th}}$ local coordinate system.

where $\,i,j=\overline{1,3}\,,\,\,\,k,m=\overline{1,N}$ ${{t}_{j}}$ and ${{u}_{j}}$ are the components of the force vectors and the displacement vectors on the boundary surface of the body;

$\begin{align} & S_{ijpkm}^{{}}\,\left( \mathbf{\xi },\,{{\mathbf{x}}_{km}} \right)=\frac{1}{8\,\pi \,\left( 1-\nu \right)\,r_{km}^{2}} \\ & \left\{ \left( 1-2\,\nu \right)\left. \left( {{r}_{km,i}}\,{{\delta }_{jp}} \right.+{{r}_{km,j}}\,{{\delta }_{pi}}-{{r}_{km,p}}\,{{\delta }_{ij}} \right)+3\,{{r}_{km,i}}\,{{r}_{km,j}}\,{{r}_{km,p}} \right\}; \\ \end{align}$

$U_{ij}^{{}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{0}} \right)=\frac{1}{16\pi G(1-\nu ){{r}_{{}}}}\left[ (3-4\nu ){{\delta }_{ij}}+{{r}_{,j}}{{r}_{,i}} \right];$

$T_{ij}^{{}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{0}} \right)=\frac{1}{8\pi (1-\nu )r_{{}}^{2}}\sum\limits_{s=1}^{3}{\left\{ \left[ (1-2\nu ){{\delta }_{ij}}+3{{r}_{,i}}{{r}_{,j}} \right]{{r}_{,s}}{{n}_{s}}\left. -(1-2\nu )\left[ {{r}_{,i}}{{n}_{j}}-{{r}_{,j}}{{n}_{i}} \right] \right\} \right.};$

${{n}_{i}}\left( \mathbf{\xi } \right)$, $\left( i=\overline{1,3} \right)$ are the components of the vector normal to the body surface at the point $\mathbf{\xi }\in S$.

Satisfying the boundary conditions on the boundary surface and on the cracks surfaces we obtain a system of $3\left( N+1 \right)$ boundary integral equations for the crack opening functions and boundary forces and displacements on the body surface:

$\begin{align} & \sum\limits_{i=1}^{3}{\sum\limits_{p=1}^{3}{\sum\limits_{k=1}^{N}{\iint\limits_{{{S}_{k}}}{\,{{K}_{ijpkm}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{km}} \right){{l}_{3pkm}}\,{{l}_{ijkm}}\,\Delta {{u}_{pk}}\left( \mathbf{\xi } \right){{d}_{\mathbf{\xi }}}{{S}_{k}}}=\frac{G}{1-\nu }\sum\limits_{i=1}^{3}{\sum\limits_{p=1}^{3}{{{P}_{j}}\,{{l}_{3pm0}}\,{{l}_{ijm0}}}}}}}- \\ & \sum\limits_{i=1}^{3}{\sum\limits_{p=1}^{3}{\sum\limits_{k=1}^{N}{\iint\limits_{S}{\left[ S_{ijp0m}^{{}}\left( \mathbf{\xi },{{\mathbf{x}}_{0m}} \right){{l}_{3p0m}}\,{{l}_{ij0m}}t_{p}^{{}}\left( \mathbf{\xi } \right)-K_{ijp0m}^{{}}\left( \xi ,\,{{\mathbf{x}}_{0m}} \right){{l}_{3p0m}}\,{{l}_{ij0m}}u_{p}^{{}}\left( \mathbf{\xi } \right) \right]{{d}_{\mathbf{\xi }}}S}}}},\,\,\, j=\overline{1,3}\,,\,\,\,m=\overline{1,N}\, \\ \end{align}$

(4)

$\begin{align} & \sum\limits_{j=1}^{3}{\iint\limits_{S}{T_{\,ij}^{{}}\left( \mathbf{x},\,\mathbf{\xi } \right)\,\left[ u_{\,j}^{{}}\left( \mathbf{\xi } \right)-u_{\,j}^{{}}\left( \mathbf{x} \right) \right]\,{{d}_{\mathbf{\xi }}}S}}- \\ & \sum\limits_{j=1}^{3}{\iint\limits_{S}{U_{ij}^{{}}\left( \mathbf{x},\,\mathbf{\xi } \right)\,P_{\,j}^{{}}\left( \mathbf{\xi } \right)\,d\,{{S}_{\mathbf{\xi }}}}}=\sum\limits_{k=1}^{N}{\sum\limits_{s=1}^{3}{\sum\limits_{j=1}^{3}{\iint\limits_{{{S}_{k}}}{T_{sjk0}^{{}}\left( \mathbf{\xi },{{\mathbf{x}}_{k0}} \right)l_{sik0}^{{}}\Delta u_{jk}^{{}}\left( \mathbf{\xi } \right){{d}_{\mathbf{\xi }}}{{S}_{k}}}}}},\,\,\,j=\overline{1,3}\,. \\ \end{align}$

(5)

---

2.4. Effective Stress Field Method for Investigation of Interacting Cracks in a Limited Body

For the approximate determination of the ${{q}^{th}}$ crack opening we propose to replace the actual stress field acting upon the crack by a certain effective field that can readily be found without solving the complicated systems of integral equations (4) and (5). For this we assume that all crack openings do not influence the displacements and forces of the boundary surface nodes. To this end, we replace the unknown functions ${{t}_{j}}$ and ${{u}_{j}}$ $\left( j=\overline{1,3} \right)$ on the right side of Eq. (4) by the corresponding functions $t_{j}^{*}$ and $u_{j}^{*}$ for a solid body without cracks whose boundary is subjected to the action of given loads. All functions of other crack openings we also replaced by the corresponding functions for isolated cracks. Therefore, we determine the effective stress fields created by each neighboring crack and by boundary load separately and obtain a system of three integral equations:

$\begin{align} & \sum\limits_{i=1}^{3}{\iint\limits_{{{S}_{q}}}{\,{{K}_{ij3qq}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{qq}} \right)\,\Delta {{u}_{iq}}\left( \mathbf{\xi } \right){{d}_{\mathbf{\xi }}}{{S}_{q}}}=\frac{1-\nu }{G}\sum\limits_{i=1}^{3}{\sum\limits_{p=1}^{3}{{{P}_{j}}\,{{l}_{3pq0}}\,{{l}_{ijq0}}}}}- \\ & \sum\limits_{i=1}^{3}{\sum\limits_{p=1}^{3}{\sum\limits_{\begin{smallmatrix} k=1 \\ k\ne q \end{smallmatrix}}^{N}{\iint\limits_{{{S}_{k}}}{{{K}_{ijpkq}}\left( \mathbf{\xi },\,{{\mathbf{x}}_{kq}} \right){{l}_{3pkq}}\,{{l}_{ijkq}}\,\Delta {{u}_{pk}}\left( \mathbf{\xi } \right)\,{{d}_{\mathbf{\xi }}}{{S}_{k}}}}}}- \\ & \sum\limits_{i=1}^{3}{\sum\limits_{p=1}^{3}{\iint\limits_{S}{\left[ S_{ijp0q}^{{}}\left( \mathbf{\xi },{{\mathbf{x}}_{0}} \right){{l}_{3p0q}}\,{{l}_{ij0q}}t_{p}^{*}\left( \mathbf{\xi } \right)-K_{ijp0q}^{{}}\left( \xi ,\,{{\mathbf{x}}_{0}} \right){{l}_{3p0q}}\,{{l}_{ij0q}}u_{p}^{*}\left( \mathbf{\xi } \right) \right]{{d}_{\mathbf{\xi }}}S}}},\,\,\,j=\overline{1,3}\,. \\ \end{align}$

(6)

---

2.5. Adaptation of BIEs to Numerical Solution

The system of singular integral equations (4) and (5) is proposed to be solved by the hybridboundary element method and collocation method.

The BIEs (5) contain the integrals over the body surface $S$ with kernels $T_{ij}^{{}}$ and $U_{ij}^{{}}$ in weakly singular form, and the integrals over the crack-surfaces ${{S}_{k}}$ and they are regular because the source point and the integration point do not coincide with each other. The BIEs (4) for direct method and the BIEs (6) for effective stress field method contain hypersingular kernels $K_{ij3qq}^{{}}$ in the integral over the crack-surface ${{S}_{q}}$, and the integrals over the fiber interface $S$ describe the influence of the body surface on the crack and they are regular for the same reason as mentioned above. Because obtained boundary integral formulation of the problem involves integrals with singularities of different orders, subsequent numerical solution of the BIEs requires special techniques, which are presented in the following section.

From the physical content of the crack-opening-displacements (CODs) it follows that these functions should vanish at the contour of the crack domain ${{S}_{k}}$. To incorporate this condition into BIEs implicitly, the CODs in the case of penny-shaped crack with radius ${{a}_{k}}$ as the complex functions are presented as

$\Delta u_{jk}^{{}}\left( \mathbf{\xi } \right)=\sqrt{a_{k}^{2}-\xi _{1}^{2}-\xi _{2}^{2}}\,\alpha _{jk}^{{}}\left( \mathbf{\xi } \right)$

(7)

where $\alpha _{jk}^{{}}$ are new unknown functions, which are defined over the domain ${{S}_{k}}$ and allotted by the differentiability properties due to their sufficient smoothness.

The representations of CODs in the form (7) simplifies the determination of SIFs in the crack vicinity, because they are expressed directly through the solutions of BIEs obtained. In particular, the mode-Ⅰ SIF $K_{\text{I}}^{{}}$ as the function of polar angle φ of the crack front point is given by the formula ^[6]:

$K_{\text{I}}^{{}}\left( \varphi \right)=-2\,\pi \,\sqrt{{{a}_{k}}\,\pi }\frac{G}{1-\nu }{{\alpha }_{3k}}\left( {{a}_{k}}cos\varphi ,\,{{a}_{k}}sin\varphi \right)$

(8)

At before discretization stage, the regularization procedure concerning the BIEs (1) or (3) or (4) or (6) defined on the crack-surface ${{S}_{k}}$ bases on the singularity subtraction technique with the following interpretation of hypersingular integrals due to the kernels ${{K}_{ij3kk}}\,\left( i,j,p=\overline{1,3};\,\,k=\overline{1,N} \right)$:

$\begin{align} & \iint\limits_{{{S}_{k}}}{\frac{\sqrt{a_{k}^{2}-\xi _{1}^{2}-\xi _{2}^{2}}}{{{\left| \mathbf{x}-\mathbf{\xi } \right|}^{3}}}{{\overline{K}}_{ijk}}(\mathbf{\xi },\mathbf{x}){{\alpha }_{j3}}(\mathbf{\xi })d{{S}_{\mathbf{\xi }}}}= \\ & I_{00}^{{}}(\mathbf{x}){{\alpha }_{j3}}(\mathbf{x})+I_{10}^{{}}(\mathbf{x})\frac{\partial {{\alpha }_{j3}}(\mathbf{x})}{\partial x_{1}^{{}}}+I_{01}^{{}}(\mathbf{x})\frac{\partial {{\alpha }_{j3}}(\mathbf{x})}{\partial {{x}_{2}}}+ \\ & \frac{1}{2}{{I}_{20}}(\mathbf{x})\frac{{{\partial }^{2}}{{\alpha }_{j3}}(\mathbf{x})}{\partial x_{1}^{2}}+I_{11}^{{}}(\mathbf{x})\frac{{{\partial }^{2}}{{\alpha }_{j3}}(\mathbf{x})}{\partial x_{1}^{{}}\partial x_{2}^{{}}}+\frac{1}{2}{{I}_{02}}(\mathbf{x})\frac{{{\partial }^{2}}{{\alpha }_{j3}}(\mathbf{x})}{\partial x_{2}^{2}}+ \\ & \iint\limits_{{{S}_{k}}}{\frac{\sqrt{a_{k}^{2}-\xi _{1}^{2}-\xi _{2}^{2}}}{{{\left| \mathbf{x}-\mathbf{\xi } \right|}^{3}}}{{\overline{K}}_{ij3k}}(\mathbf{\xi },\mathbf{x})}\left[ {{\alpha }_{j3}}(\mathbf{\xi })-{{\alpha }_{j3}}(\mathbf{x})- \right. \\ & ({{\xi }_{1}}-x_{1}^{{}})\frac{\partial {{\alpha }_{j3}}(\mathbf{x})}{\partial x_{1}^{{}}}-({{\xi }_{2}}-x_{2}^{{}})\frac{\partial {{\alpha }_{j3}}(\mathbf{x})}{\partial x_{2}^{{}}}-\frac{1}{2}{{({{\xi }_{1}}-x_{1}^{{}})}^{2}} \\ & \frac{{{\partial }^{2}}{{\alpha }_{j3}}(\mathbf{x})}{\partial x_{1}^{2}}-({{\xi }_{1}}-x_{1}^{{}})({{\xi }_{2}}-x_{2}^{{}})\frac{{{\partial }^{2}}{{\alpha }_{j3}}(\mathbf{x})}{\partial x_{1}^{{}}\partial x_{2}^{{}}}- \\ & \left. \frac{1}{2}{{({{\xi }_{2}}-x_{2}^{{}})}^{2}}\frac{{{\partial }^{2}}{{\alpha }_{j3}}(\mathbf{x})}{\partial x_{2}^{2}} \right]d{{S}_{\mathbf{\xi }}} \\ \end{align}$

(9)

where $\mathbf{x}\in {{S}_{k}}$, the regular function ${{\overline{K}}_{ij3kk}}$ describes the relation of the distances between the points $\mathbf{x}$ and $\mathbf{\xi }$ and their preimages, namely ${{\overline{K}}_{ijk}}={{K}_{ij3kk}}{{\left| \mathbf{x}-\mathbf{\xi } \right|}^{3}}$, the last integral in the right part exists in the ordinary sense, the integrals

${{I}_{ij}}\left( \mathbf{x} \right)=\iint\limits_{{{S}_{k}}}{\frac{\sqrt{a_{k}^{2}-\xi _{1}^{2}-\xi _{2}^{2}}}{{{\left| \mathbf{x}-\mathbf{\xi } \right|}^{3}}}{{\overline{K}}_{ijk}}(\mathbf{\xi },\mathbf{x}){{\left( \xi _{1}^{{}}-x_{1}^{{}} \right)}^{i}}{{\left( \xi _{2}^{{}}-x_{2}^{{}} \right)}^{j}}d{{S}_{\mathbf{\xi }}}},\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\mathbf{x}\in {{S}_{k}}$

(10)

are evaluated analytically by the integration by parts ^[6].

Applying the presentations (9) with the analytical results of the integrals (10) yields the regular version of the BIEs (6). The discretization of BIEs (6) is realized by a similar collocation scheme, where the unified boundary element mesh is formed by a uniform division of the circular domain ${{S}_{k}}$ (mapping of the crack-surface ${{S}_{k}}$) in the direction of the polar coordinates $r$ and $\varphi$(Figure 3).

The surface of the body is meshed by the boundary elements. Two-dimensional integrals over these surfaces we replace by the sum of integrals on superparametric boundary elements (Figure 2).

In this analysis, linear approximation within each element is used for the boundary quantities ${{u}_{j}}$ and ${{t}_{j}}$, while constant approximation is adopted for the crack-surface functions $\alpha _{jk}^{{}}$ by using the standard Gaussian quadrature's for the integration. For calculating weakly singular integrals, we propose the method of regularized mappings ^[9]. The main point of this method consists in that the quadrangle boundary element by using the form function is transformed to the square. Then this square divides diagonally into two triangles. Triangles transform by not one-to-one mapping to the squares so that the singular point is transformed to the side of the square. In this case the Jacobian of this mapping is equal to zero on this side and singularity is liquidated. Hypersingular integrals, the density of which is a function of cracks opening, is regularized numerically ^[6]. Discrete analogue of integral equations (4), (5) or (6) is constructed by replacing the unknown functions and kernels of integral equations on their values at the nodes of boundary elements and at the collocation nodes. A crack opening function determine the SIF ^[6].

---

3. Results

In the numerical examples the uniaxial tensile loading of constant amplitude in the normal direction relative to the cracks surfaces is applied at remote places of an infinite solid containing two interacting cracks. To verify the validity of the present effective stress field method, we examine three configurations with two interacted coplanar penny-shaped cracks of the same radii $a$, and varying distance between their centers $d$. The vector between centers of cracks is: a) in the crack plane; b) at 45 degrees to the crack plane; c) perpendicular to the crack plane. Remote loading acts perpendicular to the cracks surfaces. The Poisson's ratios of the solid are postulated as $\nu =0.33$. The mode-Ⅰ SIF was calculated for the considered cases and normalized by the mode-Ⅰ SIF $K_{\text{I}}^{*}=2{{N}_{0}}\sqrt{{a}/{\pi }\;}$ for asingle penny-shaped crack of radius $a$ in an infinite solid under the same loading conditions, so that $\bar{K}_{\text{I}}^{{}}={K_{\text{I}}^{{}}}/{K_{\text{I}}^{*}}\;$. For numerical analysis by collocation Scheme 160 nodes on each surface of the crack were used(Figure 4).

Numerical investigation was carried out also for Configuration Ⅰ of cracks interaction in a limited body. Two penny-shaped cracks of the same radius $a$ are asymmetrically located in asquare middle cross-section of the rod with respect to its center (Figure 5). The distances between the cracks centres were chosen as $d=2.2\,a$ and the distance from the boundary surface to the nearest crack is $0.1\,a$. More details about locations of cracks are presented in Figure 5. Assume that tensile load ${{N}_{0}}$ is applied to the ends of the rod in the normal direction to its cross-section (perpendicular to the cracks surfaces). In numerical simulations of limited body we use the same materialand normalization procedure of a SIF, that for an infinite solid. It is assumed also, that the length of the rod () is large enough to influence deformation of it ends to the cracks opening was absent. In our analysis for the mesh of the rod boundary surface 770 boundary elements were used.

---

4. Discussion

Due to the violation of problem symmetry mode Ⅱ and mode Ⅲ SIF also will not be zero but they have insignificant value. We present mode-Ⅰ SIF on the charts for easy comparison with the case of interaction of cracks in an infinite solid. Variation of the mode-Ⅰ SIF with the distance between the cracks centers is plotted in the Figure 4a for Configuration Ⅰ of location of cracks and in the Figure 4b and Figure 4c for Configuration Ⅱ and Ⅲ, respectively. As it is expected, a neighboring crack in the same plane provides an amplification effect for the crack, what is exhibited by increasing of the mode-Ⅰ SIF in comparison with that for a single crack. So the extremal SIF increasing is fixed at the crack front point nearest to the neighboring crack. Opposite tendency takes place for interacting penny-shaped cracks in parallel planesof the infinite solid. The results of the calculation by approximate method effective stress field for interacting of two penny-shaped cracks with a radius R can considered sufficiently accurate when the distance between the nearest points of interacting cracks is not less than: 0 .05R for Configuration Ⅰ; 1.5R for Configuration Ⅱ; 1.4R for Configuration Ⅲ. In the case where discrepancy between the proposed method and direct method is a significant, problem of inaccuracy reducing can be resolved by multiple using of effective stress field method for each interacting cracks alternately.

Numerical results for the limited body describe the SIF-behavior depending on angular coordinates of the contour of the first crack (Figure 6) and of the second crack (Figure 7). Qualitatively interaction of penny-shaped cracks in an infinite solid and in a limited body differ by more complicated distribution of the SIF along the crack front in the second case due to the appearance of additional extremums and changing their positions. So, the crack is located near the boundary surface and near the other cracks is characterized by two maximal values of the SIF. Maximal ones are observed at the front point nearest to the boundary of the body and to the neighboring crack. Location of the point with the absolute maximum of the SIF depends on the distance between the interacting objects and boundary surface of the body. Boundary surface of complex shape leads to occurrence of two non-equivalent minimums of the function of SIF along the crack front. The absolute minimum is at the front point equidistant from the boundary surface and from the neighboring crack.

---

5. Conclusion

The improved BIEM-based numerical approach has been presented in this paper to analyze the SIFs in the vicinity of two interacted penny-shaped cracks embedded in an elastic infinite solid and in a limited elastic body having rounded edges. Numerical results are presented and discussed for uniform remote loading along normal direction to the cracks surfaces, when defects are on the same plane and on the coplanar planes. The influence of the boundary surface and neighboring defectson the crack are estimated for the dependences of mode-Ⅰ SIF versus the crack-crack distance and angular coordinate of the contour of crack.

Comparison between the results obtained by direct method and effective stress field method forsolution of the problem on interaction of two penny-shaped cracks shows perfect agreement in the case of location of cracks in the same plane. In the case ofinteraction oftwo cracks in coplanar planes, we cannot use the method of effective stress field for small distances between their centers. In the case of occurrence of "not deep" cracks in the body (no more than 3 diameter of a crack from the boundary surface), influence of the body surface on the stress state at the crack front cannot be neglected. Approaching the crack to the surface of the body or to other crack in the same plane increases the SIF over all crack front. In the local sense the boundary surface influence leads to the complication of variation of the SIF along the crack front. It is expressed in increasing the number of the points with minimums and maximums of the SIF, position of which at the crack front depends on the crack-crack location and distance to the surface. Inaccuracy of the effective stress field method in the case closely located cracks is reduced three times already in the second iteration in the scheme which is shown in Figure 8.

Although our attention in the numerical analysis is focused on the mode-Ⅰ SIF, the proposed boundary integral formulation and solution procedure give the possibility for the determination of mixed-mode SIFs also for the arbitrary problem geometry and loading conditions. Proposed method can be easily generalized on the cases with multiple cracks having complicated shapes.

---

Conflict of Interest

The author declares that there is no conflict of interest regarding the publication of this manuscript.

---