Stochastic MMSIF of multiple edge cracks FGMs plates subjected to combined loading using XFEM

Abstract The second order statistics of multiple edge crack functionally graded materials (FGMs) under tensile, shear and combined loading assuming uncertain system parameters is presented in this paper. The uncertain parameters used under the present study are the material properties, and crack parameters such as crack length and crack angle. In this present analysis extended finite element method (XFEM) is used. The stochastic analysis is carried out using second order perturbation technique (SOPT) for the evaluation of mean and coefficient of variance (COV) of mixed mode stress intensity factor (MMSIF).


Introduction
FGMs are the perfect blend of two (or more) material phases of different materials and functional performance that vary in definite direction within the geometry. FGMs were at first intended to resist in thermal environment as in aerospace structural applications and fusion reactors however now daily it is utilized as imperative composite materials all through the world for various modern applications with ease of debonding issue. Safety and performance are the prime components for any structure and the presence of cracks in structure diminish its strength, toughness and lastly performance. Subsequently, it is important to look at the fracture parameters, for example, MMSIF against crack development [1].
Structures are subjected to various types of loading during many applications. Presence of these loadings, the service life of such structures decreases extensively. The impact is more articulated with the presence of discontinuities like cracks. The impact of these loading on the MM-SIF is one of the vital researches for ideal execution. Researchers are doing impressive work on deterministic investigation for examination of different discontinuities of FGMs and composite panels using conventional finite element method (FEM). In conventional FEM it is required to update the mesh during evolving of crack is very costly and time consuming. To overcome this problem, the XFEM are being used by the researches where, instead of remeshing enrichment functions are used.
In this direction, Noor [2] studied nonlinear analysis of different applications by global-local methodology to predict nonlinear and post buckling behaviour in structures, where major characteristics of these problems are identified. Swenson and Ingraffea [3] presented finite element model to study mix mode dynamic crack propagation where he removed the constraint that crack propagation is along mesh line. In this present work analysis of curve crack under biaxial loading is also presented. Beissel et al. [4] presented an algorithm of finite element analysis for the fracture mechanics on dynamic elastic-plastic behaviour which is suitable for crack propagation in any direction. It is based on tracking the path of crack tip and to indicate propagation, T* fracture parameter is employed. Song et al. [5] studied and compared three finite element methods for the dynamic crack advancement in brittle material. In this work XFEM, inter element crack method and element deletion method are utilised where, XFEM and inter element method show good result for crack speed and crack path and element deletion method can predict crack branching. Kim et al. [6] presented an analysis of interacting cracks by a generalized finite element method (GFEM) and global-local functions are used for enrichment, which is computable to analyse fracture mechanics problems in three-dimensional geometries. It is demonstrated that GFEM with enrichment functions is more robust as compared to global-local FEM. Srivastava and Lal [7] estimated stress at crack tip by utilising predictive method, where FEM is used to predict the stress intensity factor (SIF) of a multi edge cracked plate. Wang et al. [8] presented numerical study on two phase particles by utilising the XFEM in a reinforced composite material. The proposed method is used without tip-enriched functions which can avoid expensive meshing strategies and considerable flexibility is obtained.
Due to presence of uncertainties or variations in material properties, loading and crack parameters, the prediction of system behaviour using deterministic analysis is not enough. For reliable design, proposed response may not vary from actual response stochastic/ probabilistic analysis is widely used. The stochastic analysis is a numerical tool widely used for quantification of uncertainties at various levels and their effect on structural response. To overcome this difficulty, Lin and Yang [9] presented the use of Marcov random processes and utilized the basics of fracture mechanics to study of fatigue crack growth advancement. Besterfield et al. [10] applied probabilistic FEM for different uncertainties which are present in material properties, geometry, applied loads and crack geometries to estimate fatigue crack growth and probability of fatigue failure. Liu et al. [11] presented first and second order reliability method for solving curvilinear fatigue growth problem, where Monte Carlo simulation (MCS) is used. Rahman [12] evaluated the probabilistic characteristics of reliability analysis in a body with various crack geometry by applying the new dimensional decomposition method. In this work MCS technique is used for the reliability analysis of the cracked structure for random load and material properties. Xu and Rahman [13] presented dimension-reduction method to calculate statistical moments in a body for different material properties, geometry and loads, where Taylor expansion is used to estimate precise statistical moments. Xu and Rahman [14] presented computational methods to predict the failure probability in a system having random geometry, loads and material properties. In this study MCS is used which is very effective to predict probability of failure. Chakraborty and Rahman [15] proposed multiscale models for fracture analysis in FGM composite. The proposed models are computationally inexpensive models but cannot find exact probability of fracture initiation. Rahman and Rao [16] presented reliability analysis in a linear elastic body with different material properties by element-free Galerkin method, which shows good agreement with MCS. Rao and Rahman [17] proposed Galerkin-based meshless method to find SIF. In this present work it is confirmed that the result from proposed method and finite difference method show good agreement with each other. Reddy and Rao [18] utilized fractal FEM to present stochastic fracture mechanics anal-ysis and this method is compared with the MCS. Evangelatos and Spanos [19] utilized the discretization of the partial integrodifferential equation by the collocation approach. The result of utilisation of the discretization, is the formation of peridynamic stiffness obtained by finite element method. Nobile and Gentilini [20,21] studied structural behaviour of three-dimensional cracked structure by using virtual distortion method approach through Neumann expansion. This study is done for random crack location and crack depth.
Many researchers have done work based on stochastic XFEM. Nouy et al. [22,23] presented XFEM based stochastic analysis for random multi-phase materials. This method is combination of XFEM and spectral stochastic methods. Lang et al. [24] utilized the XFEM and spectral stochastic FEM to study the geometry which contains inclusion for the thermal analysis, where XFEM with Polynomial chaos on level set method are employed. Nespurek [25] presented fracture reliability and fatigue analysis of twodimensional structural. The proposed method is used by extended FEM with FORM and MCS. Gope et al. [26] employed photo elastic and FEM for studding two inline cracks and two eccentric edge cracks where, they determined different modes of SIF by FEM and the results are compared with experimental value. Chen teal [27] evaluated SIF and T-stress of cracked finite plate by the spline fictitious method, this method is accurate and efficient where, meshing is not required near crack face. Khanna et al. [28] utilized the strip yield model with plate theory to find out the plasticity effect. The results obtained by the first order plate theory utilized in advanced fatigue and fracture analysis of cracks at free surface. Kotousov et al. [29] developed an approach to find displacement in transverse direction at the tip of elastic plate. In the proposed approach first order plate theory can be useful in various problems, which can be reduced to a modified Helmholtz. Monfared et al. [30] studied mixed mode fracture behaviour for the distribution of dislocations on functionally graded orthotropic strip, where Fourier transform is utilised to develop a system of singular integral equations and these equations are solved numerically by Lobatto-Chebyshev integration formula. Muthu et al. [31] proposed a level set method for modelling multiple cracks by using coarse mesh free nodal discretization. The proposed method can handle multiple cracks and their propagation by maximum tangential principal stress criterion. Kim and Paulino [32] presented MMSIF ofcracked FGMs plate where material properties vary with respect to geometry. MMSIF is calculated by using various techniques with FEM and the present results are compared with experimental and theoretical results. Shi et al. [33] eval-uated the first mode and mixed mode SIF for the twodimensional cracked FGMs plate by utilising weight function method, where it is confirmed that present method is valid for the derivation of weight function in FGMs. In conventional FEM due to limitation in updating of mesh, many researchers have started working on XFEM, where fracture problem is solved by utilising enrichment functions without using re-meshing during evaluation of crack. Kumar [34] presented fracture analysis of components with different discontinuities by analytical methods.
In this present work MMSIF of various edges cracked FGMs plate by stochastic based XFEM under the action of tensile, shear and combined loadings in MATLAB environment is presented. In present paper, an effective stochastic method is also utilized for the prediction of fracture response regarding mean and COV of MMSIF with assumed input uncertain parameters. Many researchers have worked on cracked FGMs plate with uniaxial loading and fracture analysis of FGMs plate with shear and combined loadings is unaccounted from literature. Apart from this stochastic based fracture analysis of FGMs plate under different loadings is utmost important for reliability analysis of crack structure plate and it is also presented in this present work which makes it unique. Figure 1 shows a body of area Ω, outer boundary Γ with crack boundary (Γc). The body has uniform body loads b, and the surface load at the boundary Γt. Boundary conditions (BC) are applied at Γu, where Γ = Γu ∪ Γ t ∪ Γc.The parameters u is the recommended displacement and tis the tractions. It is assumed crack surface is traction free.

Governing general equation
The equilibrium equations and BC are written as With Where, Γc is the traction free crack The BC from equation (2-4) can be written as

The XFEM formulation
The linear equation in the discrete system can be composed as [35] [K] Here, {u h }, [K] are the vector of degrees of freedom for nodes and stiffness matrix and {F} is the external force vector. The kinematic condition for plane pressure condition can be composed as The C ij can be composed as For orthotropic materials, E and υ can be additionally characterized as The expression for crack tip enrichment functions f i (r,θ) which incorporates all displacement states can be given as sin and θ j = arctan The parameters α 1 , α 2 , β 1 and β 2 are defined as 11 , β 2 = c 12 + c 66 2c 66

Formulation of MMSIF
The classical form of the interaction integral (M-intergal) is the combination of auxiliary and actual fields for evalution of MMSIF can be composed as (Kim and Paulino [32], Mohammadi [35]) The M-interaction integral is represented as (16) and t 12 = − c 22 2 lm where s 1 and s 2 are the roots of above equation. The MMSIF can be found by considering the two states; as first state (︀ K aux I = 1 and K aux II = 0 )︀ and second state as. From Eq. (15), the actual MSIF associated with state 1 and 2 are communicated as (Shrivastava and Lal [36]), For first state: For second state: Here, K I and K II are first and second mode SIF.

Crack growth analysis with MMSIF propagation path
There are numerous criteria being broadly utilized for determination and continuity of crack path amid advancing crack for linear elastic crack issue. Out of these, global tracking algorithm given by (Oliver et al. [37]) is one the best reliable crack growth criteria. Out of all these tracking criteria, the global crack growth tracking criteria presented and finding the path of discontinuity simultaneously and also now every time it is not required to track individual crack propagation step. The mentioned criteria can be carefully and effectively utilized in various crack problems with XFEM, as it can gives great result for the prediction of crack paths. This can be done by finding extra global system equations which containing single DOF per node. An iso-linecan be composed as The scalar function ψ can be found at the crack tip, crack propagation criteria utilized for the demonstration of the crack growth. So, for present study the energy based crack growth criteria is used and from this crack propagation condition can be defined as, Where Ac is new crack part segment.

Stochastic analysis of mixed mode SIF for multiple random variables
The stochastic problem can be characterized by mean and random counterpart where random part is little when contrasted with mean part. The classified problem is substituted in expansion of Taylor series and neglecting the higher order terms since second order approximation is adequate to measure MMSIF. (Nouy et al. [38,39], Nespurek [40]). The zeroth-order equation gives mean response whereas first and second order equations give standard deviation (SD) response regarding input random variables. The proportion of mean and standard deviation is characterized as COV of reaction. A response variable (K) is the function of k where r = {︃ 1 for rst mode On expending the above equation Where x i and µ b i are the random variable inputs and expected mean. From Eq. (22), the first-order mean of K, denoted by E(K ′ ), can be expressed as From Eq. (22), the variance of response variable K for the first order is written as Where, number of variables are denoted by n and cov is resolved from the autocorrelation function of the underlying stochastic field of b, which can be represented as where, Here [C] is the covariance matrix. The parameters σ b1 , σ b2 , . . . , σ bn are termed as Where cov (b1, b 2 , . . . , b i ) and Kr (︀ µ b1 , µ b2 , . . . , µ b i )︀ are the coefficient of variance of random input variables and their mean values. The parameters ρ b i ,b j and cov are termed as the COC and coefficient of variance from different input random variables. For uncorrelated random variables, the value of cov Similarly, for correlated random variables, ρ b i ,b j is assumed Perturbation method of second order is represented as.
The variance of response variable K for second order is written as var The proportion of expected mean and relating standard deviation with random parameters is characterized as COV of MMSIF.

Results and discussion
Here, MATLAB code is utilized find MMSIF of mean and COV of edge cracked FGMs plate under different loadings by XFEM with SOPT.
The essential system random variables (b i ) used here for edge cracks plates are defined as, b 1 Where E 1 , E 2 are the longitudinal Young's modulus and a, α are the crack length and crack angle respectively.
The normalized first and second mode SIF (K I and K II ) used here for tensile (uniaxial), shear and combined loadings (tensile and shear) are represented as, The domain for computing, interaction integral is taken to exist within radius of size r d = 3 √ ae. The total number of elements considered over here are 30 × 60 considered in this analysis. The material properties and geometry are used in the present investigation as (unless otherwise stated) For double edge crack: E 1 = 1.0, E 2 = E(W) and  Figure 2 shows geometry of FGMs plate with (a) double edge cracked plate under tensile loading. (b) crack enrichments. Two cracks are eccentrically placed from the x-axis at ey. Figure 3 shows the validation of normalized first mode SIF with respect to a/W for double edge cracked isotropic plate under uniaxial tensile loading. Here, analytical study for normalized first mode SIF is compared with result of present work using XFEM. The results from the present study are near to analytical results. The normalized mean  SIF using analytically can be calculated by following expression.
Where, α = a W for, 0 < α < 0.7 f (α) = 1.12 − 0.23α + 10.55α 2 − 21.72 α 3 + 30.39α 4 Table 1 shows the validation study of normalized mean first mode SIF results using present method are validated with the results available in literature for various values of modulus ratio and crack length to width ratio. The present approach using XFEM with partition of unity approach is in very good agreement with Kim and Paulino [32] using conventional FEM and Shi et al. [33] using weighted function approach. Table 2 shows the comparative observations of normalized mean and COV of MMSIF by utilising first order perturbation technique (FOPT), SOPT and MCS for E 2 /E 1 = 0.45, a =0.45, α = 0, double crack at ey (+1 and −1). To validate the present result obtained by stochastic XFEM algorithm, results are obtained by FOPT and SOPT method for normalized mean and COV MMSIF and resembled it with MCS technique. This is because results of stochastic analysis of     Table 2, it is concluded that present FOPT and SOPT results are in very good agreement with MCS using direct use of computer code in MATLAB. The mean K I and K II is highest for combined loading while COV of first mode is highest for tensile loading. Similarly, the mean and COV of second mode SIF is highest for combined loading. With increase in the COV of input random variables, the mean and corresponding COV of first and second mode increases. With the increase of COV in crack length (b 4 ), the mean and corresponding COV of MMSIF increases in MCS, FOPT and SOPT. It is observed that higher the randomness in input parameter, higher will be the sensitiveness in MMSIF.
In MCS method, it is necessary to take large sample size. So, 3000 number of samples points are utilized for satisfactory convergence of result, which leads to high computational cost as compared to perturbation technique. Table 4 shows the effect of eccentricity (ey) where ey is a distance of crack along y-axis from centre of the plate and crack length to width ratio (a/W) with random system properties {bi, (i=4 and 5) = 0.1} on the mean and corresponding COV of mixed mode SIF of multiple edge cracked plate subjected to uniaxial tensile, shear, combine (Tensile and shear) stress for a = 0.45, α = 0, and E 2 /E 1 = 0.65 using SOPT. Here, two cracks on the left-hand side of the plate have been considered at different values of ey. It is observed that for the constant eccentricity, the mean of first mode SIF is maximum for combined loading whereas corresponding COV of first mode SIF is maximum for tensile loading. Also, for a given crack length to width ratio, with the increase of eccentricity, the mean of first mode SIF decreases for tensile loading and increases for shear and combined loading, while corresponding COV decreases. Table 5 shows the effect of crack length on the normalized mean and COV (bi= {i=4 and 5} = 0.1) of K I and K II of edge cracked FGMs plate experiencing to uniaxial tensile shear and combined stresses for α = 45 and E 2 /E 1 = 0.45.   Here, two cracks on the left-hand side of the plate have been considered at ey = (+ 1 and −1 from centre). It is observed that for constant value of crack length the mean of first mode (K I ) SIF is maximum for combined loading whereas corresponding COV of first mode (K I ) SIF is maximum for tensile loading. Here, it is also observed that when crack length to width ratio increases from 0.4 to 0.6 the mean of first mode SIF increases for all loading condition whereas corresponding COV of first mode SIF shows random nature. Table 6 shows the effect of crack angle and modulus ratios on the normalized mean and COV (bi= {i=4 and 5} = 0.1) of K I and K II of edge cracked FGMs plate under uniaxial tensile shear and combined stresses for a = 0.45 and E 2 /E 1 = 0.45. Here, two cracks on the left-hand side of the plate has been considered at ey= (+ 1 and −1) It is observed that for constant value of crack angle the mean of first mode SIF is maximum for combined loading whereas corresponding COV of first mode (K I ) SIF is maximum with tensile loading. It is also observed that when angle of the crack increases from 0 ∘ to 45 ∘ the mean of first mode SIF decreases for all loading condition whereas corresponding COV of first mode SIF shows random nature. Table 7 shows the effect of E 2 /E 1 on the normalized mean and COV (bi= {i=4 and 5} = 0.1) of K I and K II of edge cracked FGMs plate experiencing tensile (uniaxial), shear and combined stresses for a = 0.5, α = 0. Here, two cracks on the left-hand side of the plate has been considered at ey = (+ 1 and -1). It is observed that for constant value of E 2 /E 1 the mean of first mode (K I ) SIF is maximum for combined loading whereas corresponding COV of first mode SIF is maximum for tensile loading. Here, it is also observed that when E 2 /E 1 increases from 0.1 to 5 the mean of first mode SIF decreases for all loading condition whereas corresponding COV of first mode SIF shows random nature. a) When plate is under tensile stress, crack propagating path is almost parallel to crack and normalized first mode SIF increases for step (1 to 2 and 3 to 4) and increases for step (2 to 3) because of variation in stress concentration, while K II shows very less change. Also, normalized first mode SIF of first crack is more as compared to second crack. b) When plate is under shear stress, crack propagating path is going down and normalized first mode SIF increases for step (1 to 2 and 3 to 4) and decreases for step (2 to 3) while K II shows very less change. Also, the normalized first mode SIF of first crack is less as compared to second crack. c) When plate is under combined stress, crack propagating path is going down and normalized first mode SIF increases for step (1 to 4) and decreases for step (2 to 3) while K II shows very less change. Here it is also observed that normalized first mode SIF for first crack is less as compared to first crack. ing path decreases till number step 2 and again shows very minor increment till number step 3 for normalized first mode SIF. While second mode SIF shows very minor variation. Also, the normalized first mode SIF for third crack is more as compared to other cracks. b) When plate is under combined stress, normalized K I increases. While K II follows the almost constant nature. Also, the normalized first mode SIF for third crack is more as compared to other cracks.

Conclusions
The stochastic XFEM using SOPT and MCS are used to calculate the mean and COV of normalized (first and second mode) MMSIF of multiple edge crack FGMs plate under uniaxial tensile, shear and combined loadings. The following observations are listed from the limited study, 1. Random change in crack length and crack angle are most dominants as compared to other random variables, also, higher modulus ratios would be more preferable for higher reliability and safety of crack structures. 2. For the constant eccentricity ratios and E 2 /E 1 , the mean of K I is maximum for combined loading whereas corresponding COV of K I SIF is maximum for tensile loading. With the increase of eccentricity ratio, the mean of first mode SIF decreases for tensile loading and increases for combined as well as shear loading, while corresponding COV decreases. When E 2 /E 1 increases from 0.1 to 5 the mean of first mode SIF decreases for all loading condition whereas corresponding COV of first mode SIF shows random nature. 3. When crack angle increases from 0 ∘ to 45 ∘ the mean of K I decreases for all loading condition whereas corresponding COV of K I shows random nature. 4. For constant value of crack length the mean related to K I SIF is maximum for combined loading whereas corresponding COV of K I SIF is maximum for tensile loading. When crack length to width ratio increases from 0.4 to 0.6 the mean of first mode SIF increases for all loading condition whereas corresponding COV of first mode SIF shows random nature. 5. When the number of steps of crack increases then normalized first mode SIF shows random nature for double crack and increases for triple crack. While second mode SIF follows the almost constant nature.