Investigating the influence of delamination on the stiffness of composite pipes under compressive transverse loading using cohesive zone method

Abstract

The effect of delamination on the stiffness reduction of composite pipes is studied in this research. The stiffness test of filament wound composite pipes is simulated using cohesive zone method. The modeling is accomplished to study the effect of the geometrical parameters including delamination size and its position with respect to loading direction on stiffness of the composite pipes. At first, finite element results for stiffness test of a perfect pipe without delamination are validated with the experimental results according to ASTM D2412. It is seen that the finite element results agree well with experimental results. Then the finite element model is developed for composite pips with delaminated areas with different primary shapes. Thus, the effect of the size of delaminated region on longitudinal and tangential directions and also its orientation with respect to loading direction on delamination propagation and stiffness reduction of the pipes is assessed.

Keywords

delamination / composite pipes / stiffness test / cohesive zone method

Introduction

Nowadays, the use of composite structures, including reservoirs, pressure vessels, and pipes, has increased dramatically due to their unique characteristics, and thereafter, many studies have been done to identify the properties and predict their behavior in different situations. Various numerical techniques such as finite element and meshfree methods are used for simulating the crack propagation in plates and shells [¹–¹²].

In general, the failure of fiber-reinforced composites can be characterized by two main modes: breakage of resin or fibers (cracking) and separation of layers (delamination) [¹³]. Delamination can be originated from different sources like manufacturing method, applying out-of-plane loading, impact loading, etc., which causes stress concentration on the discontinuities and decreases load carrying capacity of the structure. In addition, since delamination reduces the load carrying capacity of composite layers, it will encourage buckling driven delamination [¹⁴]. One of the most important parameters in composite pipes is stiffness class. The stiffness is defined as the resistance of a pipe to the diametric deformation resulting from applying a compressive transverse loading. The stiffness test of composite pipes is carried out using the ASTM D 2412 standard [¹⁵]. In this test, a Glass-fiber Reinforced Polymers (GRP) pipe with a specific dimension is placed between two pseudo-rigid plates. The upper plate is moved downward at a controlled rate while the bottom plate is completely fixed. During the test, the load and its proportional deformation are recorded and thus the load-deformation diagram can be traced. The choice of the pipe stiffness is directly related to the depth of the burial and the traffic load on the pipeline. Rafiee and Habibagahi [¹⁶] have evaluated experimentally and numerically the stiffness of a GRP pipe with a sand/resin core layer. They also presented a simple analytical method on the basis of solid mechanic approach to estimate the stiffness of a GRP mortar pipe under compressive transverse loading [¹⁶]. Finite element analysis are often successful in predicting failure in composite structures. Such as stress-based methods, fracture mechanics and cohesive elements methods. In general, stress-based methods are suitable to predict the start of delamination, and the fracture mechanics are suitable for predicting the growth of delamination as soon as the initial cracks start. Cohesive zone method is a combination of these two methods [¹⁴]. In recent years, cohesive elements have been used extensively. Wang et al. [¹⁴] used cohesive elements method to analyze delamination in the composite plates. They also investigate the effects of some parameters such as the ratio of length to width of the plane, shape of the delamination area, the delamination size and depth of it on critical load of buckling. Also, Blackman et al. [¹⁷] studied the evolution of mode-one delamination in fiber-reinforced composites using cohesive elements. Their sample is modeled as two-dimensional part. They have also provided an analytical solution for the growth of mode-one delamination using linear traction-separation law [¹⁷]. Alfano [¹⁸] examined the effect of the shape of the cohesive zone on the prediction of crack initiation and growth in the DCB specimen, taking into account different traction-separation models such as parabolic, bilinear, and trapezoidal. Ouyang and Li [¹⁹] used a bilinear traction-separation model to predict the initiation and propagation of mode one delamination in a DCB specimen. They use the classic beam theory and do not consider the effect of shear deformation in the cohesive zone. Hélénon et al. [²⁰] specifically investigated the structural response of the composite layers of a T-shaped sample, and showed that the achievement of actual values of interlaminar strength required a very fine mesh pattern. Rafiee and Habibagahi [²¹] have investigated experimentally and numerically the damage progression in a composite pipe subjected to compressive transverse loading. Both interlaminar and intralaminar failure modes are taken into account, simultaneously. In-plane failure criteria are chosen for identifying the onset of in-plane failure mode while cohesive approach is employed for identifying the initiation of delamination as the out-of-plane failure mode.

Cohesive zone method

Cohesive behavior can be modeled either using surface-based cohesive behavior or cohesive elements. Both methods are very similar and obey the same rules, but they have different applications. The surface-based cohesive behavior is recommended when the thickness of adhesion between layer is negligible compared with the thickness of layers. In this case, the separation take place at the interface between layers. If there is a failure in adhesion which is pertinent to thicker regions, cohesive elements are required to be used.

Figure 1 shows the cohesive zone in mode I, II, and III loading and their corresponding softening behaviors. A high initial stiffness K (penalty stiffness) is used to hold the upper and lower plates of non-adhesive elements together in a linear elastic range (point 1 in Fig. 1). For all loading modes after the normal or shear stresses in the cohesive zone reached the corresponding interlaminar strength (point 2 in Fig. 1), the stiffness gradually decreased to 0. The area under each of the stress-displacement curves (modes I, II, and III) is the corresponding critical fracture energy (G_IC, G_IIC, and G_IIIC, respectively) [²²].

The aforementioned structural behavior shown in Fig. 1 can be defined as:

(1)

(2)d_i={\displaystyle \frac{{\delta }_i^f({\delta }_i^{\max }-{\delta }_i^{\mathrm{0}})}{{\delta }_i^{\max }({\delta }_i^f-{\delta }_i^{\mathrm{0}})}},i=\mathrm{1},\mathrm{2},\mathrm{3};d_i\in [\mathrm{0},\mathrm{1}].

To avoid penetration of layers following equation is considered:

(3){\tau }_{\mathrm{3}}=K{\delta }_{\mathrm{3}},\hspace{1em}{\delta }_{\mathrm{3}}\le \mathrm{0.}

In the above equations, the initial displacements can be calculated as:

(4)\begin{array}{l}{\delta }_{\mathrm{1}}^{\mathrm{0}}=T/K_{\mathrm{t}\mathrm{t}},\\ {\delta }_{\mathrm{2}}^{\mathrm{0}}=S/K_{\mathrm{s}\mathrm{s}},\\ {\delta }_{\mathrm{3}}^{\mathrm{0}}=N/K_{\mathrm{n}\mathrm{n}}\hspace{0.17em},\end{array}

where parameters N, S, and T are the values of normal strength and shear strengths in first and second directions, respectively. The final displacements are also calculated as below:

(5)\{\begin{array}{l}{\displaystyle {\int }_{\mathrm{0}}^{{\delta }_{\mathrm{3}}^f}{\tau }_{\mathrm{3}}d{\delta }_{\mathrm{3}}=G_{\mathrm{I}\mathrm{C}},}\\ {\displaystyle {\int }_{\mathrm{0}}^{{\delta }_{\mathrm{2}}^f}{\tau }_{\mathrm{2}}d\delta \mathrm{2}=G_{\mathrm{I}\mathrm{I}\mathrm{C}},}\\ {\displaystyle {\int }_{\mathrm{0}}^{{\delta }_{\mathrm{1}}^f}{\tau }_{\mathrm{1}}d{\delta }_{\mathrm{1}}=G_{\mathrm{I}\mathrm{I}\mathrm{I}\mathrm{C}}.}\end{array}

As a result, the final displacements can be calculated as:

(6)\begin{array}{l}{\delta }_{\mathrm{1}}^f=\mathrm{2}{G}_{\mathrm{I}\mathrm{I}\mathrm{I}\mathrm{C}}/T,\\ {\delta }_{\mathrm{2}}^f=\mathrm{2}{G}_{\mathrm{I}\mathrm{I}\mathrm{C}}/S,\hspace{0.17em}\hspace{0.17em}\\ {\delta }_{\mathrm{3}}^f=\mathrm{2}{G}_{\mathrm{I}\mathrm{C}}/N.\hspace{0.17em}\end{array}

It should be mention that, D is the scalar damage variable in contact points and initially it has a value of 0. As soon as evolution criterion is defined, the magnitude of D evolves from 0 to 1 [²³]. Mostly, in composite structures delamination propagation occurs in mix-mode loading. By comparing the stress components and their corresponding critical values or strengths, it is easy to determine the beginning of damage at the interface of composite layers under the mode I, II, or III loading. The quadratic stress criterion has been used to calculate the onset of damage [²³]:

(7){\left\{{\displaystyle \frac{\langle t_{\mathrm{3}}\rangle }{N}}\right\}}^{\mathrm{2}}+{\left\{{\displaystyle \frac{\langle t_{\mathrm{2}}\rangle }{S}}\right\}}^{\mathrm{2}}+{\left\{{\displaystyle \frac{\langle t_{\mathrm{1}}\rangle }{T}}\right\}}^{\mathrm{2}}=\mathrm{1},

where t denotes the maximum value of stress in the first, second and third directions. To describe the evolution of damage under a combination of normal and shear separations across the interface, it is useful to introduce an effective separation defined as [²¹]:

(8){\delta }_m=\sqrt{{\delta }_{\mathrm{1}}^{\mathrm{2}}+{\delta }_{\mathrm{2}}^{\mathrm{2}}+{\langle {\delta }_{\mathrm{3}}\rangle }^{\mathrm{2}}}=\sqrt{{\delta }_{\mathrm{s}\mathrm{h}\mathrm{e}\mathrm{a}\mathrm{r}}^{\mathrm{2}}+{\langle {\delta }_{\mathrm{3}}\rangle }^{\mathrm{2}}},

{\delta }_{\mathrm{s}\mathrm{h}\mathrm{e}\mathrm{a}\mathrm{r}}^{} represents the total shear displacement vector. In this study, the power law criterion is used to predict the growth of delamination under mix-mode condition as below [²³]:

(9){\left\{{\displaystyle \frac{G_{\mathrm{I}}}{G_{\mathrm{I}\mathrm{C}}}}\right\}}^{\alpha }+{\left\{{\displaystyle \frac{G_{\mathrm{I}\mathrm{I}}}{G_{\mathrm{I}\mathrm{I}\mathrm{C}}}}\right\}}^{\alpha }+{\left\{{\displaystyle \frac{G_{\mathrm{I}\mathrm{I}\mathrm{I}}}{G_{\mathrm{I}\mathrm{I}\mathrm{I}\mathrm{C}}}}\right\}}^{\alpha }=\mathrm{1.}

Camanho and Devila [²²] studied the effect of a in power law on the prediction of mix mode delamination propagation prediction for epoxy resin and found that assuming a = 1 is conservative enough for accurate prediction of delamination propagation. So, due to complicated procedure for determining this parameter and taking into account that it is intended to compare the results for a same pipe with different delaminated areas, selecting this value for a can be considered as an acceptable compromise in modeling.

Geometrical modeling

Figure 2 illustrates the schematic presentation of a GRP pipe that has been exposed to a compressive transverse loading. Mechanical properties of composite layers and cohesive interface are presented in Tables 1 and 2, respectively. The ply configuration of the investigated pipe is [(±48°)₆]. The overall thickness of the cylindrical shells is H and the inner shell thickness is h.

Experimental study

Figure 3 shows the stiffness testing apparatus. According to ASTM D 2412, a test specimen with the minimum length of 300 mm was prepared and loaded between two rigid parallel flat plates at a controlled rate of approach to one another. Specimens was compressed till the diametric deflection reaches the 5% of the specimen diameter. The test stopped at this stage and applied compressive load was recorded. The recorded compressive load for investigated GFRP pipe at 5% deflection was 9190 N.

Finite element analysis

In this research, ABAQUS commercial finite element software is utilized to simulate the composite pipe stiffness test and to investigate delamination. First, an analysis was carried out without considering the initial delamination to predict the stiffness of the pipe and compare it to the experimental results. Then, in h/H=\mathrm{1}/\mathrm{3}, 15 different delamination modes were simulated. In 11 modes, the delamination was carried out in the form of tangential (circumferential) and longitudinal ellipses at the top of the pipe (zero angle relative to the loading direction) and the next four modes are in 30, 45, 60, and 90 degrees relative to loading direction (Fig. 4). In Fig. 4(b), a and b are oval diameters.

Both internal and external cylindrical parts have been simulated using 4-node quadrilateral conventional shell (S4R) Elements. Due to the very small thickness of entrapped resin between adjacent layers, the adhesion property between the two layers is defined by the surface-based cohesive behavior method and the coefficients of stiffness are considered as follows [²⁴]:

(10)K_{\mathrm{n}\mathrm{n}}=K_{\mathrm{s}\mathrm{s}}=K_{\mathrm{t}\mathrm{t}}={\displaystyle \frac{E_{\mathrm{22}}}{t}},

where t is the thickness of the resin in the above equation. A suitable interaction is defined between the upper plate and the pipe through appropriate coefficient of tangential friction resembling contact status between them.

The size of the elements are chosen small enough to avoid any dependency of the results to the mesh size [²¹].

After defining the criteria for the onset of damage, delamination evolution is defined based on the energy that is dissipated as a result of the damage process. Therefore, gradual degradation process starts on the basis of linear softening law and it continues until the complete separation. It must be noted that cohesive elements is used just in the interface of two pre-delaminated layers.

Results and discussion

Reaction force calculation without considering delamination

To ensure the accuracy of the modeling, two types of 3-D FE models are constructed. In the first model, a cylindrical shell with an internal diameter of 198 mm and a length of 300 mm is constructed using S4R conventional shell element, which consists of 12 layers of Glass/Epoxy composite plies. The magnitude of reaction force at 5% diametric deflection is inserted in Table 3.

In the second model, an internal 4-layer cylindrical part and an external 8-layer cylindrical part are constructed (Fig. 2). In this model, surface-based cohesive behavior is defined at the interface of internal and external parts. It is noteworthy to mention that there is no initial delamination in this model. Obviously, this cohesive behavior is considered between the fourth and fifth layers. The amount of reaction force at 5% diametric deflection is calculated and inserted in Table 3.

It can be seen that the results of finite element simulations are in good agreements with the result of the experiment.

Reaction force calculation with considering delamination

In this part, the effect of delamination size and its location on the stiffness of pipe has been investigated.

Delamination size effect

Eleven models with delamination defects in h/H=\mathrm{1}/\mathrm{3} are simulated with zero angle relative to the loading direction. In all cases, the delamination region has an ellipsoidal shape with diameters a, b and the same area equal to 0.02. By changing the ratio of b/a=\mathrm{1},\mathrm{2},\mathrm{...},\mathrm{6} and b/a=\mathrm{1}/\mathrm{2},\mathrm{1}/\mathrm{3},\mathrm{...},\mathrm{1}/\mathrm{6}, circumferential and longitudinal ellipsoidal delamination shapes are obtained. Nonlinear analysis is taken into account for all models and the magnitudes of reaction force are shown in Fig. 5.

According to Figs. 5(a) and 5(b), the smaller the ratio of b/a, the lower the reaction force and stiffness of the composite pipe subsequently. In addition, Fig. 5(a) shows that the stiffness of composite pipe reduces when the major axis of the elliptical delaminated region is extended perpendicular to the loading direction in comparison with circular delaminated region. In contrast, the pipe stiffness increases when the major axis is placed along circumferential direction compared with circular delaminated region.

Figure 6 shows the damage evolution in longitudinal and circumferential directions. The red zone indicates that the scalar damage parameter d has reached the value of 1 and the cohesive elements are completely damaged.

As it is evident from this Fig. 6, in models where the delamination defect is in tangential direction, the growth of delamination is mostly in longitudinal directions; while in the models accommodating delamination defect in the longitudinal direction, the propagation of the delamination is more in tangential direction.

Effect of delamination location

In this section four models with constant b/a ratio of 1 are constructed where delamination defect is placed along 30, 45, 60, and 90 degrees with respect to the loading direction. The obtained results are compared at h/H=\mathrm{1}/\mathrm{3}, as shown in Fig. 4. The magnitude of reaction force at 5% diametric deflection for different angles can be seen in Fig. 7.

Figure 7 reveals that the existence of delamination damage at the interface of composite layers, decreases the stiffness of composite pipe and the smaller angle of the delamination position with respect to loading direction, leads to more stiffness reduction of the composite pipe. As it is shown in Fig. 8, the growth of damage at 30°is concentrated in the upper part of the delaminated region (in the loading direction), and with increasing angle, separation grows more uniformly around the delamination region. Delamination propagation for different angles are illustrated in Fig. 8.

Conclusions

It is very obvious that the existence of delamination phenomenon at the interface of adjacent layers reduces the stiffness of composite pipes. Naturally, the larger the delamination area, the greater its effect. It is understood that when the major axis of an elliptical delamination defect is oriented along circumferential direction, increasing the length of its major axis will be led to less stiffness reduction of pipe. On the contrary, when the major axis of the elliptical delaminated region is placed along longitudinal direction, the influence of any increase in its length will pronounce the pipe stiffness reduction. Moreover, the propagation of elliptical delaminated region along the minor axis is more than that of along major axis as the elliptical delamination defect tends to be converted to a circular shape.

Examining the different angles of delamination defect with respect to the loading direction, it can also be concluded that, the lower angles of the delamination center with respect to loading direction, leads to more stiffness reduction of the composite pipe. When the center of delaminated region is perpendicular to the loading direction (i.e., 90°), almost no element is experiencing complete failure. It is interpreted that in case of delamination occurrence in a pipe, from practical point of view, the pipe can be rotated along its longitudinal axis for extending the operational lifetime of the pipe.

References

[1] Rabczuk T, Zi G, Bordas S, Nguyen-Xuan H. A simple and robust three-dimensional cracking-particle method without enrichment. Computer Methods in Applied Mechanics and Engineering, 2010, 199(37–40): 2437–2455

[2] Rabczuk T, Belytschko T. Cracking particles: A simplified meshfree method for arbitrary evolving cracks. International Journal for Numerical Methods in Engineering, 2004, 61(13): 2316–2343

[3] Rabczuk T, Areias P M A, Belytschko T. A meshfree thin shell method for non-linear dynamic fracture. International Journal for Numerical Methods in Engineering, 2007, 72(5): 524–548

[4] Rabczuk T, Belytschko T. A three-dimensional large deformation meshfree method for arbitrary evolving cracks. Computer Methods in Applied Mechanics and Engineering, 2007, 196(29–30): 2777–2799

[5] Rabczuk T, Gracie R, Song J, Belytschko T. Immersed particle method for fluid-structure interaction. International Journal for Numerical Methods in Engineering, 2010, 81(1): 48–71

[6] Rabczuk T, Bordas S, Zi G. On three-dimensional modelling of crack growth using partition of unity methods. Computers & Structures, 2010, 88(23–24): 1391–1411

[7] Areias P, Rabczuk T, Msekh M A. Phase-field analysis of finite-strain plates and shells including element subdivision. Computer Methods in Applied Mechanics and Engineering, 2016, 312(C): 322–350

[8] Nguyen-Thanh N, Valizadeh N, Nguyen M N, Nguyen-Xuan H, Rabczuk T. An extended isogeometric thin shell analysis based on Kirchhoff-Love theory. Computer Methods in Applied Mechanics and Engineering, 2015, 284: 265–291

[9] Areias P, Rabczuk T, Camanho P P. Finite strain fracture of 2D problems with injected anisotropic softening elements. Theoretical and Applied Fracture Mechanics, 2014, 72: 50–63

[10] Areiasae P, Rabczukb T, Dias-da-Costacd D. Element-wise fracture algorithm based on rotation of edges. Engineering Fracture Mechanics, 2013, 110: 113–137

[11] Areias P, Rabczuk T. Finite strain fracture of plates and shells with configurational forces and edge rotations. International Journal for Numerical Methods in Engineering, 2013, 94(12): 1099–1122

[12] Amiri F, Millán D, Shen Y, Rabczuk T, Arroyo M. Phase-field modeling of fracture in linear thin shells. Theoretical and Applied Fracture Mechanics, 2014, 69: 102–109

[13] Soteropoulos D, Fetfatsidis K A, Sherwood J A. Using Abaqus to model delamination in fiber reinforced composite materials. In: SIMULIA Community Conference, 2002

[14] Wang R G, Zhang L, Zhang J, Liu W B, He X D. Numerical analysis of delamination buckling and growth in slender laminated composite using cohesive element method. Computational Materials Science, 2010, 50(1): 20–31

[15] Standard AAN. Standard test method for determination of external loading characteristics of plastic/composite pipe by parallel/plate loading. ASTM-D2412, 2004

[16] Rafiee R, Habibagahi M R. On the stiffness prediction of GFRP pipes subjected to transverse loading. KSCE Journal of Civil Engineering, 2018, 22(11): 4564–4572

[17] Blackman B R K, Hadavinia H, Kinlochand A J, Williams J G. The use of a cohesive zone model to study the fracture of fibre composites and adhesively bonded joints. International Journal of Fracture, 2003, 119: 25–46

[18] Alfano G. On the Influence of the shape of the interface law on the application of cohesive-zone models. Composites Science and Technology, 2006, 66: 723–730

[19] Ouyang Z, Li G. Local damage evolution of double cantilever beam specimens during crack initiation process: A natural boundary condition based method. Journal of Applied Mechanic, 2009, 76(5): 051003

[20] Hélénon F, Wisnom M R, Hallett S R, Trask R S. Numerical investigation into failure of laminated composite T-piece specimens under tensile loading. Composites. Part A, Applied Science and Manufacturing, 2012, 43(7): 1017–1027

[21] Rafiee R, Habibagahi M R. Evaluating mechanical performance of GFRP pipes subjected to transverse loading. Thin-walled Structures, 2018, 131: 347–359

[22] Camanho P P, Davila C G. Mixed-Mode Decohesion Finite Elements for The Simulation of Delamination in Composite Materials. NASA/TM-211737. 2002

[23] Alfano G, Crisfield M A. Finite element interface models for the delamination analysis of laminated composites: Mechanical and computational issues. International Journal for Numerical Methods in Engineering, 2001, 50(7): 1701–1736

[24] Diniz Melo J D, Levy Neto F, de Araujo Barros G, de Almeida Mesquita F N. Mechanical behavior of GRP pressure pipes with addition of quartz sand filler. Journal of Composite Materials, 2011, 45(6): 717–726

RIGHTS & PERMISSIONS

Higher Education Press and Springer-Verlag GmbH Germany, part of Springer Nature