Computation of energy release rate for interfacial crack in orthotropic bimaterials

<!-- metadata commented in wiki content


==Computation of energy release rate for interfacial crack in orthotropic bimaterials==

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
'''Salah Bouziane<sup>a,b</sup>, Hamoudi Bouzerd<sup>a,b</sup>'''</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<span style="text-align: center; font-size: 75%;"><sup>a</sup>Civil Engineering department, University of August 20, 1955, Skikda, Algeria. E-mails: [mailto:bouziane_21@yahoo.fr bouziane_21@yahoo.fr], [mailto:bib_ham@yahoo.fr bib_ham@yahoo.fr]</span></div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<span style="text-align: center; font-size: 75%;"><sup>b</sup>Laboratory of Civil Engineering and Hydraulics, University of May 8, 1945, Guelma, Algeria.</span></div>
-->

'''Abstract'''

The interfacial crack in bimaterials is a very interesting problem for composite materials and which has received particular attention from several researchers. In this study, we will propose a numerical modeling of the interfacial crack between two orthotropic materials using a special mixed finite element. For the calculation of the energy release rate, a technique, based on the association of the present mixed finite element with the virtual crack extension method, was used. The numerical model proposed, in this work, was used to study a problem of interfacial crack in bimaterials. Two cases were treated: isotropic and orthotropic bimaterials. The results obtained, using the present element, were compared with the values of the analytical solution and other numerical models found in the literature.

'''Keywords: '''Interfacial crack, Mixed finite element, Virtual crack extension method, Energy release rate, Orthotropic bimaterials

==1. Introduction==

The interfacial fracture is a complex phenomenon which is still badly understood, so it would already be enough to justify its study. Indeed, the interface located between two dissimilar materials is, on the mechanical level, a weak point: when these materials are requested by stresses, of thermal origin for example, the fracture of the interface is a mode usually observed. Moreover, one knows little about the mechanical conditions which lead to this fracture. A comprehension of the interfacial fracture thus represents a significant stake in the field of composite materials.

The problem of the interfacial crack in isotropic bimaterials has been treated by many researchers. We can cite, for example, the work of Williams [1],  Erdogan [2,3], England [4], Rice and Sih [5], Hutchinson et al. [6], Rice [7] and Suo and Hutchinson [8].

The cracks along the interface between two anisotropic plates were initially treated by Gotoh [9]. The case of plane deformation of interfacial crack between two anisotropic materials was studied by Clements [10], Willis [11], Qu and Bassani [12], Suo [13] and Ni and Nemat-Nasser [14]. Bassani and Qu [15] have explicitly resolved the special case of Griffith's problem and the solution of the general problem has been found by Suo [13] and Qu and Li [16]. The crack path in the anisotropic medium was studied theoretically and numerically by Gao et al. [17], a weak plane model was adopted to characterize the anisotropic fracture toughness and the maximum energy release rate criterion was chosen to predict the crack path. The problem of interfacial cracks in anisotropic bimaterials was also treated by Wang et al. [18], Juan and Dingreville [19].

Based on anisotropic elasticity, Tanaka et al. [20] evaluate the energy release rate by the modified crack closure integral of the finite element method, and convert to the stress intensity factor for the cases of cracks on elastic symmetrical planes. Two approaches have been described by Banks-Sills and Ikeda [21] for considering an interface crack between two anisotropic materials. Both approaches have been used for orthotropic and monoclinic materials. The problem of cracked orthotropic bimaterial was also studied by Bouchemella et al. [22]. Fracture analysis of orthotropic cracked media was investigated by applying the recently developed Extended IsoGeometric Analysis (XIGA) [23]  using the T-spline basis functions [24]. The same method XIGA was used by Habib et al. [25](2017) for the analysis of static fracture behaviour for a crack in orthotropic materials.

Khatir and Wahab [26] used an inverse algorithm based on Proper Orthogonal Decomposition (POD) and Radial Basis Functions (RBF) for single and multiple cracks identification in plate structures. The inverse analyses combine experimental fracture mechanics tests with numerical models based on the XIGA method. The eXtended IsoGeometric Analysis combined with Particle Swarm Optimization (PSO) have been  used for crack identification in two-dimensional linear elastic problems (plate) based on inverse problem [27].

In this paper, a numerical modeling has been proposed to study the interfacial crack between two orthotropic materials. This model uses a two-dimensional mixed finite element developed in a natural plane. It is an element with 7 nodes: 5 displacement nodes and 2 stress nodes. The proposed model was used to calculate the energy release rate in a cracked orthotropic bimaterial using a technique that combines the present element with the virtual crack extension method. In this work, two cases of interfacial cracks were treated: an isotropic bimaterial and an orthotropic bimaterial. The results obtained, using the present mixed finite element, were compared with the values of the analytical solution and other numerical models found in the literature.

==2. Numerical modelling of interfacial crack==

The bimaterial has been discredized using a special mixed finite element RMQ-7 ('''R'''eissner '''M'''odified '''Q'''uadrilateral) as shown in Figure 1(a). The present mixed finite element used in this study is two-dimensional element with seven nodes: five displacement nodes and two stress nodes as shown in Figure 1(b). The node 5 coincides with the crack tip. This element was developed by Bouzerd [28], in the physical (x, y) plane, and was reformulated and validated by Bouziane et al. [29] in a natural (ξ, η) plane.

Displacement for the present mixed finite element can be given by

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>u=\sum _{i=1}^{5}{N}_{i}{u}_{i}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(1)</span>
|}


<span style="text-align: center; font-size: 75%;">where </span> <math display="inline">{N}_{i}</math><span style="text-align: center; font-size: 75%;">'' ''are the shape functions and </span> <math display="inline">{u}_{i}</math><span style="text-align: center; font-size: 75%;"> is the nodal displacement corresponding to node ''i''. For the present element, the shape functions are given as follows</span>

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="text-align: center;"|<math display="inline">{N}_{1}=-\frac{1}{4}\left( 1-\right. </math><math>\left. \xi \right) \left( 1-\eta \right) \xi</math> <span style="text-align: center; font-size: 75%;">'',        ''</span> <math display="inline">{N}_{2}=</math><math>\frac{1}{4}\left( 1+\xi \right) \left( 1-\eta \right) \xi</math>

<math display="inline">{N}_{3}=\frac{1}{4}\left( 1+\xi \right) \left( 1+\eta \right)</math> <span style="text-align: center; font-size: 75%;">'',         ''</span> <math display="inline">{N}_{4}=</math><math>\frac{1}{4}\left( 1-\xi \right) \left( 1+\eta \right)</math> <span style="text-align: center; font-size: 75%;">'',     ''</span> <math display="inline">\, \, {N}_{5}=</math><math>\frac{1}{2}\left( 1-{\xi }^{2}\right) \left( 1-\eta \right)</math> 
|  style="text-align: right;"|<span style="text-align: center; font-size: 75%;">(2)</span>
|}


<span style="text-align: center; font-size: 75%;">The element stress component is approximated by</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math display="inline">\left\{ \sigma \right\} =\left[ M\right] \left\{ \tau \right\}</math> <span style="text-align: center; font-size: 75%;">'' ''</span>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(3)</span>
|}


<span style="text-align: center; font-size: 75%;">where ''[M]'' is the matrix of interpolation functions for stresses and ''{τ}'' is the vector of nodal stresses.</span>

For the RMQ-7 element (see Figure 1(b)), the shape functions ''M<sub>i2</sub>'', used to evaluate ''σ<sub>12</sub> ''and ''σ<sub>22</sub>''[29] for nodes 6 and 7 are obtained by

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math display="inline">{M}_{i2}^{6}=\frac{1}{6}\left( 1-2\xi \right) \left( 1-\right. </math><math>\left. 2\eta \right)</math> <span style="text-align: center; font-size: 75%;">'',   ''</span> <math display="inline">{M}_{i2}^{7}=</math><math>\frac{1}{6}\left( 1+2\xi \right) \left( 1-2\eta \right) \quad ,\quad i=1,2</math>'' ''
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(4)</span>
|}


<span style="text-align: center; font-size: 75%;">The element stiffness  matrix ''[K<sub>e</sub>]'' is given by the following expression </span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math display="inline">\left[ {K}_{e}\right] =\left[ \begin{matrix}\left[ {K}_{\sigma \sigma }\right] &\left[ {K}_{\sigma u}\right] \\{\left[ {K}_{\sigma u}\right] }^{T}&\left[ 0\right] \end{matrix}\right]</math> <span style="text-align: center; font-size: 75%;">'' ''</span>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(5)</span>
|}


<span style="text-align: center; font-size: 75%;">where the sub-matrices, </span> <math display="inline">\left[ {K}_{\sigma \sigma }\right]</math> <span style="text-align: center; font-size: 75%;"> and </span> <math display="inline">\left[ {K}_{\sigma u}\right]</math> <span style="text-align: center; font-size: 75%;">, are given by the following relations</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\left[ {K}_{\sigma \sigma }\right] =-t\int_{{A}_{e}}^{}{\left[ M\right] }^{T}\left[ S\right] \left[ M\right] d{A}^{e}</math>

<math>\, \left[ {K}_{\sigma u}\right] =t\int_{{A}_{e}}^{}{\left[ M\right] }^{T}\left[ B\right] d{A}^{e}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(6)</span>
|}


<span style="text-align: center; font-size: 75%;">where ''[S]'' is the compliance matrix, ''[M]'' is the matrix of interpolation functions for stresses, ''[B]'' is the strain-displacement matrix</span> of shape function derivatives<span style="text-align: center; font-size: 75%;">, ''t'' is the thickness,</span> <span style="text-align: center; font-size: 75%;">''A<sup>e</sup>'' is the element area and ''T'' indicate the matrix transpose.</span>

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="text-align: center;width: 69%;"|[[Image:Draft_بوزيان_803129046-image1.png|378px]] 
|  style="text-align: center;width: 30%;"|[[Image:Draft_بوزيان_803129046-image2.png|162px]] 
|-
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">'''(a)''' Discretization of bimaterial</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">'''(b)'''  RMQ-7 element</span>
|-
|  colspan='2'  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">'''Figure 1.'''  Discretization of bimaterial and RMQ-7 element</span>
|}


==3. Computation of energy release rate==

<span style="text-align: center; font-size: 75%;">The virtual crack extension method, associated with the mixed finite element RMQ-7, is used to calculate the energy release rate ''G'' [28]. In this technique, the first calculation of the deformation energy  </span> <math display="inline">{\Pi }_{1}</math><span style="text-align: center; font-size: 75%;"> is carried out in the initial configuration of the crack. The crack is then moved an infinitesimal distance </span> <math display="inline">\delta a</math><span style="text-align: center; font-size: 75%;"> in the direction of its axis. The deformation energy  </span> <math display="inline">{\Pi }_{2}</math><span style="text-align: center; font-size: 75%;"> is evaluated again in the second configuration, the energy released during this crack length variation is</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\delta \Pi ={\Pi }_{2}-\, {\Pi }_{1}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(7)</span>
|}


<span style="text-align: center; font-size: 75%;">The energy release rate ''G'' will be evaluated thereafter starting from the relation</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>G=\frac{\delta \Pi }{\delta a}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(8)</span>
|}


<span style="text-align: center; font-size: 75%;">Calculation by the virtual crack extension method requires two finite element analysis. The use of the RMQ-7 element makes it possible to introduce one mesh for the calculation of the energy release rate, which represents a considerable profit in computing times and setting data compared to the traditional techniques which use two meshes [28].</span>

<span style="text-align: center; font-size: 75%;">Indeed the intermediate displacement node of  the RMQ-7 element is associated to crack tip, and consequently the length of crack "''a''" can be increased by a quantity </span> <math display="inline">\delta a</math><span style="text-align: center; font-size: 75%;"> while acting inside strict of the crack element by translation of the tip crack node without disturbing the remainder of the mesh.</span>

<span style="text-align: center; font-size: 75%;">With the assumption on materials and displacements (linear elastic behaviour and small displacements), the solutions  </span> <math display="inline">u(a)</math><span style="text-align: center; font-size: 75%;"> and  </span> <math display="inline">u(a+</math><math>\delta a)</math><span style="text-align: center; font-size: 75%;"> obtained in the structure with a crack length "''a''" and in the same structure with a crack length </span> <math display="inline">a+</math><math>\delta a</math><span style="text-align: center; font-size: 75%;"> are as close as the disturbance </span> <math display="inline">\delta a</math><span style="text-align: center; font-size: 75%;"> is small compared to dimensions of the crack element. We can thus write with a good approximation</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>u\left( a\right) =u(a+\delta a)</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(9)</span>
|}


<span style="text-align: center; font-size: 75%;">Several calculations on simple examples enabled us to confirm the relation Equation (9), which is theoretically coherent and physically acceptable, considering the assumptions used.</span>

<span style="text-align: center; font-size: 75%;">If we consider that the external loading does not vary during the increase </span> <math display="inline">\delta a</math><span style="text-align: center; font-size: 75%;">, the energy release rate is calculated as follows</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>G=-\frac{\Pi \left( a+\delta a\right) -\Pi (a)}{\delta a}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(10)</span>
|}


<span style="text-align: center; font-size: 75%;">where'' ''</span> <math display="inline">\Pi \left( a+\right. </math><math>\left. \delta a\right)</math> <span style="text-align: center; font-size: 75%;"> and  </span> <math display="inline">\Pi (a)</math><span style="text-align: center; font-size: 75%;"> represent respectively the deformation energy of the cracked structure in the configuration  </span> <math display="inline">a+</math><math>\delta a</math><span style="text-align: center; font-size: 75%;"> and "''a''".  </span>

In its discretized form, the deformation energy is written

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\Pi =\frac{1}{2}\sum _{i=1}^{ne}{\left\{ u\right\} }_{i}^{T}{\left[ K\right] }_{i}{\left\{ u\right\} }_{i}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(11)</span>
|}


<span style="text-align: center; font-size: 75%;">with:</span>

''ne'' : total number of elements in discretized structure,

<math display="inline">{\left\{ u\right\} }_{i}</math><span style="text-align: center; font-size: 75%;">'':'' vertical vector containing the nodal values of element i,</span>

<math display="inline">{\left[ K\right] }_{i}</math><span style="text-align: center; font-size: 75%;">: elementary matrix of element i, and the exponent "''T''" indicates the transposed vector.</span>

By substitution of Equation (11) in Equation (10), the expression of the energy release rate ''G'' becomes

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>G=-\frac{1}{2\delta a}\left[ \sum _{i=1}^{ne}{\left\{ u(a+\delta a)\right\} }_{i}^{T}{\left[ K(a+\delta a)\right] }_{i}{\left\{ u(a+\delta a)\right\} }_{i}-\right. </math><math>\left. \sum _{i=1}^{ne}{\left\{ u(a)\right\} }_{i}^{T}{\left[ K(a)\right] }_{i}{\left\{ u(a)\right\} }_{i}\right]</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(12)</span>
|}


<span style="text-align: center; font-size: 75%;">Taking account of Equation (9), the expression Equation (12) can be written in the following form</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>G=-\frac{1}{2\delta a}\sum _{i=1}^{ne}{\left\{ u(a+\delta a)\right\} }_{i}^{T}\left[ {\left[ K(a+\delta a)\right] }_{i}-{\left[ K(a)\right] }_{i}\right] {\left\{ u(a+\delta a)\right\} }_{i}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(13)</span>
|}


<span style="text-align: center; font-size: 75%;">and as only the crack element is disturbed, then ''G'' results more simply in the relation</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>G=-\frac{1}{2\delta a}{\left\{ u(a+\delta a)\right\} }_{f}^{T}\left[ {\left[ K(a+\delta a)\right] }_{f}-\right. </math><math>\left. {\left[ K(a)\right] }_{fi}\right] {\left\{ u(a+\delta a)\right\} }_{f}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(14)</span>
|}


<span style="text-align: center; font-size: 75%;">where the index "''f''" indicates that the matrix and vector used are those of the crack element.</span>

<span style="text-align: center; font-size: 75%;">The expression Equation (14) shows that only the crack element is concerned, and consequently it is enough to place in the mesh another RMQ-7 element equivalent to that placed on the crack, in other words an element which has the same geometry and made up of same material. The energy release rate is calculated according to the relation Equation (14) with only one discretization starting from the difference of the elementary matrices of the element containing the crack and representing the state  </span> <math display="inline">a+</math><math>\delta a</math><span style="text-align: center; font-size: 75%;">  and its equivalent element representing the state "''a''".  The expression Equation (14) can be written differently as follows</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{G=-\frac{1}{2}\left\{ u\right\} }_{f}^{T}\frac{{\left[ \delta K\right] }_{f}}{\delta a}{\left\{ u\right\} }_{f}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(15)</span>
|}


<span style="text-align: center; font-size: 75%;">In practice, we carry out the discretization of the cracked structure in the configuration </span> <math display="inline">a+</math><math>\Delta a</math><span style="text-align: center; font-size: 75%;">, and we locate the element containing the crack like its equivalent element representing the configuration "''a''", in order to save  their elementary matrices during the assembly operation and before the application of the boundary conditions.</span>

After the resolution phase, the nodal values of the crack element are extracted, and a special module is used to evaluate the energy release rate according to the following formula

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{G=-\frac{1}{2}\left\{ u\right\} }_{f}^{T}\frac{{\left[ \Delta K\right] }_{f}}{\Delta a}{\left\{ u\right\} }_{f}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(16)</span>
|}


==4. Numerical examples ==

<span id='_Hlk536098584'></span><span style="text-align: center; font-size: 75%;">'''4.1 Presentation of the example'''</span>

<span style="text-align: center; font-size: 75%;">The example treated, in this study, is the interfacial crack centered of a bimaterial plate. This example was studied by Chow et al. [30] with plane strain condition. This rectangular bimaterial is made of material #1 and #2 and subjected to a tension σ<sub>22</sub><sup>0</sup> = 1 MPa.</span> <span style="text-align: center; font-size: 75%;">As shown in Figure 2, the dimensions of the bimaterial are the half crack length a=1mm, the width w=20a and the height h=20a. Two cases are treated in this example. In the first case it is assumed that the materials #1 and #2 are isotropic and in the second case the materials are considered to be orthotropic (carbon composites: AS4/3501-6) with lay-up angle of 0 and 90 degree. The material properties of the used materials are defined in Table 1.</span>

{| style="width: 100%;margin: 1em auto 0.1em auto;border-collapse: collapse;" 
|-
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;"> [[Image:Draft_بوزيان_803129046-image3.png|198px]] </span>
|-
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">'''Figure 2.''' Bimaterial plate</span>


|}


<span style="text-align: center; font-size: 75%;">A stress σ<sup>0</sup><sub>11 </sub>is applied to the side of the material #2. In the case of plane strain, this stress is expressed by</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\sigma }_{11}^{0}=\left[ \frac{{\nu }_{12\#2}+{\nu }_{13\#2}{\nu }_{32\#2}}{1+{\nu }_{13\#2}{\nu }_{31\#2}}-\right. </math><math>\left. \frac{{\nu }_{12\#1}+{\nu }_{13\#1}{\nu }_{32\#1}}{1+{\nu }_{13\#2}{\nu }_{31\#2}}\left( \frac{{E}_{1\#2}}{{E}_{1\#1}}\right) \right] {\sigma }_{22}^{0}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(17)</span>
|}


<span style="text-align: center; font-size: 75%;">where E is the Young's modulus and ν is the Poisson's ratio of the material. </span>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<span style="text-align: center; font-size: 75%;">'''Table 1.''' Material property</span></div>

{| style="width: 86%;border-collapse: collapse;" 
|-
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;vertical-align: top;"|'''Isotropic'''
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;vertical-align: top;"|'''Orthotropic (0 degree)'''
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;vertical-align: top;"|'''Orthotropic (90 degree)'''
|-
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;vertical-align: top;"|G<sub>#1</sub> = 1GPa

ν<sub>#1</sub> = ν<sub>#2</sub> = 0.3
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;vertical-align: top;"|E<sub>3</sub> = 142 GPa

E<sub>1</sub>/E<sub>3</sub> =E<sub>2</sub>/E<sub>3</sub> = 6.91x10<sup>-2</sup>

G<sub>12</sub>/E<sub>3</sub> = 2.68x10<sup>-2</sup>

G<sub>13</sub>/E<sub>3</sub>=G<sub>23</sub>/E<sub>3</sub>=4.23x10<sup>-2</sup>

ν<sub>31</sub> = ν<sub>32</sub> = ν<sub>12</sub> = 0.3
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;vertical-align: top;"|E<sub>1</sub> = 142 GPa

E<sub>2</sub>/E<sub>1</sub> =E<sub>3</sub>/E<sub>1</sub> = 6.91x10<sup>-2</sup>

G<sub>23</sub>/E<sub>1</sub> = 2.68x10<sup>-2</sup>

G<sub>13</sub>/E<sub>1</sub>=G<sub>12</sub>/E<sub>1</sub>=4.23x10<sup>-2</sup>

ν<sub>12</sub> = ν<sub>13</sub> = ν<sub>23</sub> = 0.3
|}


<span style="text-align: center; font-size: 75%;">In the example above, the authors (Chow et al. 1995) calculate and compare the stress intensity factors K<sub>1</sub> and K<sub>2</sub>, the energy release rate is calculated according to K<sub>1</sub> and K<sub>2</sub> by the expression given by Qu and Bassani [31]. The results are resumed in Table 2 for the two materials (isotropic and orthotropic).</span>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<span style="text-align: center; font-size: 75%;">'''Table 2.''' Energy release rate in the numerical example</span></div>

{| style="width: 100%;border-collapse: collapse;" 
|-
|  rowspan='2' colspan='2'  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">'''Material'''</span>
|  rowspan='2' style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">'''Exact solution'''</span>
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">'''Hybrid element'''</span>
|  colspan='2'  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">'''Mutual integral'''</span>
|  colspan='2'  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">'''Extrapolation technique'''</span>
|-
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">205 nodes</span>
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">679 nodes</span>
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">237 nodes</span>
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">679 nodes</span>
|  style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">237 nodes</span>
|-
|  rowspan='3' style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|<span style="text-align: center; font-size: 75%;">Isotropic</span>
|  style="border-top: 1pt solid black;"|<span style="text-align: center; font-size: 75%;">G<sub>#2</sub></span><span style="text-align: center; font-size: 75%;">/G<sub>#1</sub></span><span style="text-align: center; font-size: 75%;">=1</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">10,988E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">11,290E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">11,302 E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">11,253 E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">13,132 E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">12,554E-04</span>
|-
| <span style="text-align: center; font-size: 75%;">G<sub>#2</sub></span><span style="text-align: center; font-size: 75%;">/G<sub>#1</sub></span><span style="text-align: center; font-size: 75%;">=5</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">06,453E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">06,606E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">06,614 E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">06,592 E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">07,649 E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">07,326E-04</span>
|-
|  style="border-bottom: 1pt solid black;"|<span style="text-align: center; font-size: 75%;">G<sub>#2</sub></span><span style="text-align: center; font-size: 75%;">/G<sub>#1</sub></span><span style="text-align: center; font-size: 75%;">=50</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">05,353E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">05,460E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">05,461 E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">05,444 E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">06,287 E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">06,026E-04</span>
|-
|  rowspan='3' style="border-top: 1pt solid black;border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">Orthotropic</span>


|  style="border-top: 1pt solid black;"|<span style="text-align: center; font-size: 75%;">[0/0]</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">03,170E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">03,257E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">03,262 E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">03,247 E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">03,793 E-04</span>
|  style="border-top: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">03,540E-04</span>
|-
| <span style="text-align: center; font-size: 75%;">[90/90]</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">02,200E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">02,221E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">02,216 E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">02,221 E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">02,549 E-04</span>
|  style="text-align: center;"|<span style="text-align: center; font-size: 75%;">02,480E-04</span>
|-
|  style="border-bottom: 1pt solid black;"|<span style="text-align: center; font-size: 75%;">[0/90]</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">02,640E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">02,685E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">02,679 E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">02,675 E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">03,094 E-04</span>
|  style="border-bottom: 1pt solid black;text-align: center;"|<span style="text-align: center; font-size: 75%;">03,021E-04</span>
|}


<span style="text-align: center; font-size: 75%;">'''4.2 Results and discussions'''</span>

The mixed finite element RMQ-7 is now used to calculate the energy release rate of the cracked bimaterial plate. For this purpose three meshes (207, 237 and 677 nodes) are used in order to be able to compare the results of RMQ-7 element with the other elements results by using approximately the same number of nodes. The results obtained are resumed in the Table 3.

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<span style="text-align: center; font-size: 75%;">'''Table 3.''' Energy release rate obtained using RMQ-7 element</span></div>

{| style="width: 83%;margin: 1em auto 0.1em auto;border-collapse: collapse;" 
|-
|  rowspan='2' colspan='2'  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|'''Material'''
|  colspan='3'  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|'''RMQ-7 mixed finite element'''
|-
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|207 nodes
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|237 nodes
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|677 nodes
|-
|  rowspan='3' style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|Isotropic
|  style="border-top: 1pt solid black;vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=1
|  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|11,272E-04
|  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|11,205E-04
|  style="border-top: 1pt solid black;text-align: center;vertical-align: top;"|11,126E-04
|-
|  style="vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=5
|  style="text-align: center;vertical-align: top;"|06,393E-04
|  style="text-align: center;vertical-align: top;"|06,486E-04
|  style="text-align: center;vertical-align: top;"|06,438E-04
|-
|  style="vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=50
|  style="text-align: center;vertical-align: top;"|05,274E-04
|  style="text-align: center;vertical-align: top;"|05,278E-04
|  style="text-align: center;vertical-align: top;"|05,297E-04
|-
|  rowspan='3' style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|Orthotropic
|  style="vertical-align: top;"|[0/0]
|  style="text-align: center;vertical-align: top;"|03,225E-04
|  style="text-align: center;vertical-align: top;"|03,237E-04
|  style="text-align: center;vertical-align: top;"|03,167 E-04
|-
|  style="vertical-align: top;"|[90/90]
|  style="text-align: center;vertical-align: top;"|02,260E-04
|  style="text-align: center;vertical-align: top;"|02,293E-04
|  style="text-align: center;vertical-align: top;"|02,168 E-04
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|[0/90]
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|02,691E-04
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|02,764E-04
|  style="border-bottom: 1pt solid black;text-align: center;vertical-align: top;"|02,617 E-04
|}


<span style="text-align: center; font-size: 75%;">According to the number of nodes, the numerical results of the energy release rate for different methods are listed in Table 4, 5 and 6 for both the isotropic bimaterial and anisotropic bimaterial. The difference with exact solution for the different methods are calculated and consigned in Tables 4, 5 and 6. This difference is expressed by the Error (%) calculated as follows</span>

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>Error\, (\%)=\frac{G-{G}_{exact}}{{G}_{exact}}\times 100\, \%</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span style="text-align: center; font-size: 75%;">(18)</span>
|}


<span style="text-align: center; font-size: 75%;">Compared to the exact solution, the numerical results show the accuracy and efficiency of the RMQ-7 element. The difference between the values of exact solution and those of the mixed finite element vary between -0,10% and 4,70%.</span>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<span style="text-align: center; font-size: 75%;">'''Table 4.'''  Energy release rate for crack along bimaterial interface, Mesh 1: 207 nodes</span></div>

{| style="width: 100%;border-collapse: collapse;" 
|-
|  rowspan='2' colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''Material'''
|  rowspan='2' style="border-top: 1pt solid black;vertical-align: top;"|'''Exact solution'''
|  colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''RMQ-7 element'''
|  colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''Hybrid element'''
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|207 nodes
|  style="border-bottom: 1pt solid black;vertical-align: top;"|Error %
|  style="border-bottom: 1pt solid black;vertical-align: top;"|205 nodes
|  style="border-bottom: 1pt solid black;vertical-align: top;"|Error %
|-
|  rowspan='3' style="border-top: 1pt solid black;vertical-align: top;"|Isotropic
|  style="border-top: 1pt solid black;vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=1
|  style="border-top: 1pt solid black;vertical-align: top;"|10,988E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|11,272E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|2,58
|  style="border-top: 1pt solid black;vertical-align: top;"|11,290E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|2,75
|-
|  style="vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=5
|  style="vertical-align: top;"|06,453E-04
|  style="vertical-align: top;"|06,393E-04
|  style="vertical-align: top;"|-0,93
|  style="vertical-align: top;"|06,606E-04
|  style="vertical-align: top;"|2,37
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=50
|  style="border-bottom: 1pt solid black;vertical-align: top;"|05,353E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|05,274E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|-1,48
|  style="border-bottom: 1pt solid black;vertical-align: top;"|05,460E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|2,00
|-
|  rowspan='3' style="border-top: 1pt solid black;border-bottom: 1pt solid black;vertical-align: top;"|Orthotropic


|  style="border-top: 1pt solid black;vertical-align: top;"|[0/0]
|  style="border-top: 1pt solid black;vertical-align: top;"|03,170E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|03,225E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|1,74
|  style="border-top: 1pt solid black;vertical-align: top;"|03,257E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|2,74
|-
|  style="vertical-align: top;"|[90/90]
|  style="vertical-align: top;"|02,200E-04
|  style="vertical-align: top;"|02,260E-04
|  style="vertical-align: top;"|2,73
|  style="vertical-align: top;"|02,221E-04
|  style="vertical-align: top;"|0,95
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|[0/90]
|  style="border-bottom: 1pt solid black;vertical-align: top;"|02,640E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|02,691E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|1,93
|  style="border-bottom: 1pt solid black;vertical-align: top;"|02,685E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|1,70
|}


<span style="text-align: center; font-size: 75%;">For isotropic bimaterials, the RMQ-7 element, for the same number of nodes, shows a clear superiority compared to the eight noded isoparametric displacement finite element (extrapolation technique), and more accurate results than those of the mutual integral method. For example, with the RMQ-7 element, the Error passed from -0,93% to 2,58% with 207 nodes whereas the Error varied from 2,00% to 2,75 using the hybrid element with 205 nodes. For orthotropic bimaterials, the element RMQ-7 shows its performance compared to the classical displacement element. It still gives results clearly closer to the exact solution. Compared to the mutual integral method the RMQ-7 element gives very satisfactory results. Using RMQ-7 element with 677, the difference varied between   -0,10% and -1,45% whereas it is between 0,73% and 2,90% using mutual method with 679 nodes.</span>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<span style="text-align: center; font-size: 75%;">'''Table 5.''' Energy release rate for crack along bimaterial interface, Mesh 2: 237 nodes</span></div>

{| style="width: 100%;border-collapse: collapse;" 
|-
|  rowspan='2' colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''Material'''
|  rowspan='2' style="border-top: 1pt solid black;vertical-align: top;"|'''Exact solution'''
|  colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''RMQ-7 element'''
|  colspan='3'  style="border-top: 1pt solid black;vertical-align: top;"|'''Mutual integral'''
|  colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''Extrapolation technique'''
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|237 nodes
|  style="border-bottom: 1pt solid black;vertical-align: top;"|Error %
|  style="border-bottom: 1pt solid black;vertical-align: top;"|237 nodes
|  colspan='2'  style="vertical-align: top;"|Error %
|  style="border-bottom: 1pt solid black;vertical-align: top;"|237 nodes
|  style="border-bottom: 1pt solid black;vertical-align: top;"|Error %
|-
|  rowspan='3' style="border-top: 1pt solid black;vertical-align: top;"|Isotropic
|  style="border-top: 1pt solid black;vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=1
|  style="border-top: 1pt solid black;vertical-align: top;"|10,988E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|11,205E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|1,98
|  colspan='2'  style="vertical-align: top;"|11,253E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|2,41
|  style="border-top: 1pt solid black;vertical-align: top;"|12,554E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|14,25
|-
|  style="vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=5
|  style="vertical-align: top;"|06,453E-04
|  style="vertical-align: top;"|06,486E-04
|  style="vertical-align: top;"|0,51
|  colspan='2'  style="vertical-align: top;"|06,592E-04
|  style="vertical-align: top;"|2,15
|  style="vertical-align: top;"|07,326E-04
|  style="vertical-align: top;"|13,53
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=50
|  style="border-bottom: 1pt solid black;vertical-align: top;"|05,353E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|05,278E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|-1,40
|  colspan='2'  style="border-bottom: 1pt solid black;vertical-align: top;"|05,444E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|1,70
|  style="border-bottom: 1pt solid black;vertical-align: top;"|06,026E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|12,57
|-
|  rowspan='3' style="border-top: 1pt solid black;border-bottom: 1pt solid black;vertical-align: top;"|Orthotropic


|  style="border-top: 1pt solid black;vertical-align: top;"|[0/0]
|  style="border-top: 1pt solid black;vertical-align: top;"|03,170E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|03,237E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|2,11
|  colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|03,247E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|2,43
|  style="border-top: 1pt solid black;vertical-align: top;"|03,540E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|11,67
|-
|  style="vertical-align: top;"|[90/90]
|  style="vertical-align: top;"|02,200E-04
|  style="vertical-align: top;"|02,293E-04
|  style="vertical-align: top;"|4,23
|  colspan='2'  style="vertical-align: top;"|02,221E-04
|  style="vertical-align: top;"|0,95
|  style="vertical-align: top;"|02,480E-04
|  style="vertical-align: top;"|12,73
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|[0/90]
|  style="border-bottom: 1pt solid black;vertical-align: top;"|02,640E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|02,764E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|4,70
|  colspan='2'  style="border-bottom: 1pt solid black;vertical-align: top;"|02,675E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|1,33
|  style="border-bottom: 1pt solid black;vertical-align: top;"|03,021E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|14,43
|}


<span style="text-align: center; font-size: 75%;">The results obtained, using the present mixed finite element, show the efficiency and accuracy of the proposed numerical model, which can give an acceptable solution with a few degrees of freedom from a unique mesh. It should be noted that during numerical calculation, the choice of the variation of the crack length </span> <math display="inline">\Delta a</math><span style="text-align: center; font-size: 75%;"> is very significant. Indeed, it is necessary that this variation is sufficiently small so that the approximation Equation (9) has a justification, and not too small to avoid problems involved in the precision machine.</span>

The results show also, that the current techniques of the finite elements analysis make it possible to find an effective numerical solution and a high precision to the problems of fracture mechanic.

==Table 6. Energy release rate for crack along bimaterial interface, Mesh 3: 677 nodes==

{| style="width: 100%;border-collapse: collapse;" 
|-
|  rowspan='2' colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''Material'''
|  rowspan='2' style="border-top: 1pt solid black;vertical-align: top;"|'''Exact solution'''
|  colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''RMQ-7 element'''
|  colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''Mutual integral'''
|  colspan='2'  style="border-top: 1pt solid black;vertical-align: top;"|'''Extrapolation technique'''
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|677 nodes
|  style="border-bottom: 1pt solid black;vertical-align: top;"|Error %
|  style="border-bottom: 1pt solid black;vertical-align: top;"|679 nodes
|  style="border-bottom: 1pt solid black;vertical-align: top;"|Error %
|  style="border-bottom: 1pt solid black;vertical-align: top;"|679 nodes
|  style="border-bottom: 1pt solid black;vertical-align: top;"|Error %
|-
|  rowspan='3' style="border-top: 1pt solid black;vertical-align: top;"|Isotropic
|  style="border-top: 1pt solid black;vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=1
|  style="border-top: 1pt solid black;vertical-align: top;"|10,988E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|11,126E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|1,26
|  style="border-top: 1pt solid black;vertical-align: top;"|11,302E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|2,86
|  style="border-top: 1pt solid black;vertical-align: top;"|13,132E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|19,51
|-
|  style="vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=5
|  style="vertical-align: top;"|06,453E-04
|  style="vertical-align: top;"|06,438E-04
|  style="vertical-align: top;"|-0,23
|  style="vertical-align: top;"|06,614E-04
|  style="vertical-align: top;"|2,49
|  style="vertical-align: top;"|07,649E-04
|  style="vertical-align: top;"|18,53
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|G<sub>#2</sub>/G<sub>#1</sub>=50
|  style="border-bottom: 1pt solid black;vertical-align: top;"|05,353E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|05,297E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|-1,05
|  style="border-bottom: 1pt solid black;vertical-align: top;"|05,461E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|2,02
|  style="border-bottom: 1pt solid black;vertical-align: top;"|06,287E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|17,45
|-
|  rowspan='3' style="border-top: 1pt solid black;border-bottom: 1pt solid black;vertical-align: top;"|Orthotropic


|  style="border-top: 1pt solid black;vertical-align: top;"|[0/0]
|  style="border-top: 1pt solid black;vertical-align: top;"|03,170E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|03,167 E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|-0,10
|  style="border-top: 1pt solid black;vertical-align: top;"|03,262E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|2,90
|  style="border-top: 1pt solid black;vertical-align: top;"|03,793E-04
|  style="border-top: 1pt solid black;vertical-align: top;"|19,65
|-
|  style="vertical-align: top;"|[90/90]
|  style="vertical-align: top;"|02,200E-04
|  style="vertical-align: top;"|02,168 E-04
|  style="vertical-align: top;"|-1,45
|  style="vertical-align: top;"|02,216E-04
|  style="vertical-align: top;"|0,73
|  style="vertical-align: top;"|02,549E-04
|  style="vertical-align: top;"|15,86
|-
|  style="border-bottom: 1pt solid black;vertical-align: top;"|[0/90]
|  style="border-bottom: 1pt solid black;vertical-align: top;"|02,640E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|02,617 E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|-0,87
|  style="border-bottom: 1pt solid black;vertical-align: top;"|02,679E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|1,48
|  style="border-bottom: 1pt solid black;vertical-align: top;"|03,094E-04
|  style="border-bottom: 1pt solid black;vertical-align: top;"|17,20
|}


==5. Conclusion==

In this paper, a numerical modeling has been proposed to study the interfacial crack between two orthotropic materials. This model uses a special mixed finite element developed in a natural plane. It is a two-dimensional element with seven nodes: five displacement nodes and two stress nodes.

The proposed numerical model was used to calculate the energy release rate in a cracked orthotropic bimaterial using a technique that combines the present element with the virtual crack extension method. Two cases were treated: isotropic and orthotropic bimaterials.

The accuracy and the efficiency of the proposed model has been evaluated by comparing the numerical solution with an available analytical solution or numerical ones obtained from others methods. Comparisons with existing analytical or numerical solutions for interfacial cracks in orthotropic bimaterials proved that the proposed model provide a good accuracy and efficiency.

==References==

[1] Williams M. L. The stresses around a fault or crack in dissimilar media. Bulletin of seismology society of America, 49: 199-204, 1959.

[2] Erdogan F. Stress distribution in nonhomogeneous elastic plane with cracks. J. appl. Mech., 30: 232-236, 1963.

[3] Erdogan F. Stress distribution in bonded dissimilar materials with cracks.  J. appl. Mech., 32: 403-409, 1965.

[4] England A. H. A crack between dissimilar media. J. appl. Mech., 32: 400-402, 1965.

[5] Rice J. R., Sih G. C. Plane problems of cracks in dissimilar media.  J. appl. Mech., 32: 418-423, 1965.

[6] Hutchinson J. W., Mear M. , Rice J. R. Crack Paralleling an Interface Between Dissimilar Materials. ASME Journal of Applied Mechanics, 54: 828-832, 1987.

[7] Rice J. R. Elastic Fracture Mechanics Concepts for Interfacial Cracks. ASME Journal of Applied Mechanics, 55: 98-103,1988.

[8] Suo Z., Hutchinson J.W. Interface crack between two elastic layers. Int J Fract., 43: 1–18, 1990.

[9] Gotoh  M.  Some problems of bonded anisotropic plates with cracks along the bond. Int. J.  Fract. Mech., 3: 253-265, 1967.

[10] Clements D.L. A crack between dissimilar anisotropic media.  Int.  J.  Engng.  Sci., 9: 257–265, 1971.

[11] Willis J. R. Fracture mechanics of interfacial cracks. J. Mech. Phys. Solids, 19: 353-368, 1971.

[12] Qu J., Bassani J. L. Cracks on bimaterial and bicrystal interfaces.  J. Mech. Phys. Solids, 37(4): 417-433, 1989.

[13] Suo  Z. Singularities, interfaces and cracks in dissimilar anisotropic media. Proc. R.  Soc. Lond., A 427: 331-358, 1990.

[14] Ni L.,  Nemat-Nasser S.  Interface crack in anisotropic dissimilar materials: An analytic solution. J. Mech. Phys. Solids, 39(1): 113-144, 1991.

[15] Bassani J. L., Qu J. Finite crack on bimaterial and bicrystal interfaces. J. Mech. Phys. Solids, 37(4): 435-453, 1989.

[16] Qu J., Li Q. Interfacial dislocation and its applications to interface cracks in anisotropic bimaterials. J. Elasticity, 26: 169-195, 1991.

[17] Gao Y., Liu Z., Zeng Q., Wang T., Zhuang  Z., Hwang K-C. Theoretical and numerical prediction of crack path in the material with anisotropic fracture toughness. Engineering Fracture Mechanics, 180: 330-347, 2017.

<span id='bau010'></span><span id='bau015'>[18] [https://www.sciencedirect.com/science/article/pii/S0020768315000347 Wang] X., [https://www.sciencedirect.com/science/article/pii/S0020768315000347 Zhou], K., [https://www.sciencedirect.com/science/article/pii/S0020768315000347 Wu] M.S. Interface cracks with surface elasticity in anisotropic bimaterials. [https://www.sciencedirect.com/science/journal/00207683 Int. J. of Solids and Structures], 59: 110-120, 2015.</span>

<span id='bau0005'></span><span id='bau0010'>[19] Juan [https://www.sciencedirect.com/science/article/pii/S0022509616305610 P.-A.],  [https://www.sciencedirect.com/science/article/pii/S0022509616305610 Dingreville] R.  Mechanics of finite cracks in dissimilar anisotropic elastic media considering interfacial elasticity.  [https://www.sciencedirect.com/science/journal/00225096 J. of the Mechanics and Physics of Solids], [https://www.sciencedirect.com/science/journal/00225096/99/supp/C 99]: 1-18, 2017.</span>

[20] Tanaka K., Oharada K., Yamada D., Shimizu K. Fatigue crack propagation in short-carbon-fiber reinforced plastics evaluated based on anisotropic fracture mechanics. Int. Journal of Fatigue, 92: 415-425, 2016.

[21] Banks-Sills L. and Ikeda T. Stress intensity factors for interface cracks between orthotropic and monoclinic material. Int. J. Fract., 167(1): 47-56, 2011.

[22] Bouchemella S., Bouzerd, H.,  Khaldi, N. Modélisation des interfaces fissurées des bimatériaux orthotropes. XIXème Congrès Français de Mécanique, Marseille, France, 2009.

[23] Ghorashi S. Sh., Valizadeh  N., Mohammadi S. Extended isogeometric analysis (XIGA) for simulation of stationary and propagating cracks. Int. J. Numer. Methods Eng., 89: 1069 –1101, 2012.

[24] Ghorashi S. Sh., Valizadeh N., Mohammadi S. , Rabczuk T. T-spline based XIGA for fracture analysis of orthotropic media. Computers and Structures, 147: 138–146, 2015.

[25] Habib S. H., Belaidi I., Khatir S., Abdel Wahab M.  Numerical Simulation of cracked orthotropic materials using extended isogeometric analysis.  Journal of Physics: Conf. Series 842 012061, 2017.

[26] Khatir S., Wahab M. A. Fast simulations for solving fracture mechanics inverse problems using POD-RBF XIGA and Jaya algorithm. Engineering Fracture Mechanics, 205: 285-300, 2019.

[27] Khatir S., Wahab M. A., Benaissa B., Köppen  M. Crack identification using eXtended IsoGeometric Analysis and particle swarm optimization. In Fracture, Fatigue and Wear, pp. 210-222, Springer, Singapore, 2018.

[28] Bouzerd H.  Elément fini mixte pour interface cohérente ou fissurée. Thèse de doctorat, Université de Claude Bernard (Lyon I), France, 1992.

[29] Bouziane S., Bouzerd H., Guenfoud M. Mixed finite element for modelling interfaces. European Journal of Computational Mechanics, 18(2): 155-175, 2009.

[30] Chow W. T., Beom H. G., Alturi S. N. Calculation of stress intensity factors  for an interfacial crack between dissimilar anisotropic media, using a hybrid element method and mutual integral:  Computational Mechanics, 15(6): 546-557, 1995.

[31] Qu J., Bassani  J. L. Interfacial fracture mechanics for anisotropic bimaterial.  Journal of Applied Mechanics, 60: 422-431, 1993.