Theoretical study of failure in composite pressure vessels subjected to low-velocity impact and internal pressure

Abstract

A theoretical solution is aimed to be developed in this research for predicting the failure in internally pressurized composite pressure vessels exposed to low-velocity impact. Both in-plane and out-of-plane failure modes are taken into account simultaneously and thus all components of the stress and strain fields are derived. For this purpose, layer-wise theory is employed in a composite cylinder under internal pressure and low-velocity impact. Obtained stress/strain components are fed into appropriate failure criteria for investigating the occurrence of failure. In case of experiencing any in-plane failure mode, the evolution of damage is modeled using progressive damage modeling in the context of continuum damage mechanics. Namely, mechanical properties of failed ply are degraded and stress analysis is performed on the updated status of the model. In the event of delamination occurrence, the solution is terminated. The obtained results are validated with available experimental observations in open literature. It is observed that the sequence of in-plane failure and delamination varies by increasing the impact energy.

Keywords

composite pressure vessel / low-velocity impact / failure / theoretical solution / progressive damage modeling

Introduction

The composite pressure vessels with special features such as high strength/stiffness-to-weight ratio, significant fatigue endurance and the superior corrosion resistance were driven forward to be used in widespread applications like aerospace and aviation industries, military, automotive, civil, medical and sport sectors. The technological evolution of pressure vessels can be traced in five distinct generations [¹]. The first generation of composite pressure vessels is full metallic vessels where high weight and low fatigue resistance are their main shortcomings. Over wrapping composite layers on a metallic liner, the second and third generations of pressure vessels have been developed where both weight reduction and prolonged fatigue lifetime goals are achieved. In the second generation, the cylindrical section of the liner is reinforced with composite layers where fibers are oriented along hoop direction. While in the third generation, both end caps and cylindrical section of the liner are fully overwrapped with composite not only in hoop direction but also in helical configuration. Replacing the metallic liner in the third generation with a high-density polymeric liner, a full composite pressure vessel has been developed as the fourth generation [¹]. Recently born fifth generation is a liner-less composite pressure vessel where nano-particles have been incorporated into resin system preventing gas diffusion through the wall thickness in absence of liner.

The safety of composite pressure vessels is required to be examined from various aspects other than internal pressure in accordance with the normative standards [²,³]. As the most challenging safety issue, composite pressure vessels need to withstand the low-velocity impact of external objects. The ability of composite pressure vessels to maintain their structural integrity against low-velocity impact loading is necessary to be verified through a series of costly experiments [²,³]. Composite vessels can be exposed to impact loads not only during production, transportation and installation but also during their operational services. Composite materials are intrinsically weak against transverse impact loads and in some cases the induced failures are not visible. Therefore, simulating the influence of impact on the structural behavior of composite pressure vessels can considerably reduce experimental costs during the early stages of design process. The majority of investigations have analyzed the mechanical response of composite pressure vessels subjected to low-velocity impact either experimentally or numerically via commercial finite element packages [⁴–¹³]. In contrast, very limited studies have theoretically attempted to simulate impact loads on composite vessels using either classical lamination theory (CLT) or first/high-shear deformation theory (FSDT or HSDT) [¹⁴–¹⁶].

Besides the considerable efforts devoted to investigate the failure of composite vessels solely under internal pressure [¹⁷–³¹], some investigations have studied the influence of a low-velocity impact on the vessel when they are not pressurized resembling the impacts prior to operation [⁶–⁸,¹³–¹⁶]. It is of great importance to study the influence of low-velocity impact on composite pressure vessels when they are internally pressurized representing operational conditions. The main purpose of this study is to theoretically predict the occurrence of failure in composite pressure vessels subjected to low-velocity impact loading and internal pressure. Both in-plane and out-of-plane failure modes are considered in parallel and thus Layer-Wise Theory (LWT) [³²] is utilized to extract all in-plane and out-of-plane stress/strain components. The focus of the current research is twofold: 1) both interlaminar and intralaminar failure modes are concurrently taken into account, 2) an internally pressurized composite vessel is exposed to the low-velocity impact.

The current research is an extension of the previously conducted study by the first author and his coworkers [³³] for the purpose of simulating low velocity impact on a composite pressure vessel in absence of internal pressure loading. Unlike the previous study [³³], progressive damage modeling is performed for the purpose of evaluating failure in an internally pressurized composite vessel subjected to the low velocity impact by an external falling ball.

Modeling procedure

The impact loading is assumed to be uniformly applied to a square surface as depicted in Fig. 1 and expanded using Fourier transformation as below [³⁴]:

(1)q(x,\theta ,t)={\displaystyle \sum }_{m=\mathrm{1}}^{\infty }{\displaystyle \sum }_{n=\mathrm{1}}^{\infty }{Q}_{mn}\left(t\right)\cos \left({\displaystyle \frac{m\mathrm{\pi }x}{a}}\right)\cos \left({\displaystyle \frac{n\mathrm{\pi }\theta }{b}}\right),

(2)Q_{mn}\left(t\right)={\displaystyle \frac{\mathrm{4}}{ab}}{{\displaystyle \int }}_{-u/\mathrm{2}}^{u/\mathrm{2}}{{\displaystyle \int }}_{-v/\mathrm{2}}^{v/\mathrm{2}}F\left(t\right)\cos \left({\displaystyle \frac{m\mathrm{\pi }x}{a}}\right)\cos \left({\displaystyle \frac{n\mathrm{\pi }\theta }{b}}\right)\mathrm{d}x\mathrm{d}\theta .

The impact force function F(t) is the result of contact deformation between the striker and the shell. A simple spring-mass model is employed to approximate the contact force. The contact force between the striker and the laminated shell during impact is assumed to be governed by Hertzian theory. By linearizing the Hertzian contact force-deformation, F(t) is determined as below [³⁵]:

(3)F\left(t\right)=K_{\mathrm{2}}\left[A_{\mathrm{1}}\left(C_{\mathrm{1}}-\mathrm{1}\right)\sin {\omega }_{\mathrm{1}}t+A_{\mathrm{2}}\left(C_{\mathrm{2}}-\mathrm{1}\right)\sin {\omega }_{\mathrm{2}}t\right].

The corresponding coefficients of Eq. (3) are expressed in Appendix A.

In the framework of LWT, the displacement fields of a laminated composite are written as below [³²]:

(4)\begin{array}{l}u(x,y,z,t)={\displaystyle \sum }_{I=\mathrm{1}}^N{U}_I(x,y,t){\Phi }^I\left(z\right),\hfill \\ v(x,y,z,t)={\displaystyle \sum }_{I=\mathrm{1}}^N{V}_I(x,y,t){\Phi }^I\left(z\right),\hfill \\ w(x,y,z,t)={\displaystyle \sum }_{I=\mathrm{1}}^N{W}_I(x,y,t){\psi }^I\left(z\right),\hfill \end{array}

where (U_I, V_I, W_I) denotes the nodal values of (u, v, w) and N is the number of nodes as shown in the Fig. 2. {\varphi }^I and {\psi }^I are the global interpolation functions for discretizing in-plane and vertical displacement through the thickness, respectively. Linear variations of aforementioned interpolation functions are expressed as below [³²]:

\begin{array}{l}{\mathrm{\Phi }}^{\mathrm{1}}\left(z\right)={\psi }_{\mathrm{1}}^{\left(\mathrm{1}\right)}\left(z\right),z_{\mathrm{1}}\le z\le z_{\mathrm{2}},\hfill \\ {\mathrm{\Phi }}^I\left(z\right)=\{\begin{array}{cc}{\psi }_{\mathrm{2}}^{(I-\mathrm{1})}\left(z\right),\hfill & z_{I-\mathrm{1}}\le z\le z_I,\hfill \\ {\psi }_{\mathrm{1}}^{\left(I\right)}\left(z\right),\hfill & z_I\le z\le z_{I+\mathrm{1}},\hfill \end{array}\hfill \\ \\ {\mathrm{\Phi }}^N\left(z\right)={\psi }_{\mathrm{2}}^{\left(N_e\right)}\left(z\right),\hspace{1em}{z}_{N-\mathrm{1}}\le z\le z_N,\hfill \end{array}

(5)I=\mathrm{2},\mathrm{3},\leftrightharpoons N_e,

(6)\begin{array}{l}{\psi }_{\mathrm{1}}^{\left(k\right)}=\mathrm{1}-{\displaystyle \frac{\overline{z}}{h_k},}\hfill \\ {\psi }_{\mathrm{2}}^{\left(k\right)}={\displaystyle \frac{\overline{z}}{h_k}},\hfill \end{array}0\le \overline{z}\le h_k,

where N_e indicates the number of layers through the thickness and h_k is the thickness of the k\mathrm{t}\mathrm{h} layer.\overline{z} is taken as z-z_t^k and z_t^k represents the z-coordinate of the top surface of the k\mathrm{t}\mathrm{h} layer [³²].

In parallel with impact loading, the selected square surface undergoes in-plane strain components as below due to internal pressure:

(7)\left\{\begin{array}{l}{\epsilon }_{xx}^{\mathrm{0}}\hfill \\ {\epsilon }_{\theta \theta }^{\mathrm{0}}\hfill \\ {\epsilon }_{\theta x}^{\mathrm{0}}\hfill \end{array}\right\}=\left[\begin{array}{ccc}{a}_{\mathrm{11}}\hfill & a_{\mathrm{12}}\hfill & a_{\mathrm{16}}\hfill \\ a_{\mathrm{21}}\hfill & a_{\mathrm{22}}\hfill & a_{\mathrm{26}}\hfill \\ a_{\mathrm{61}}\hfill & a_{\mathrm{62}}\hfill & a_{\mathrm{66}}\hfill \end{array}\right]\left\{\begin{array}{l}h{\sigma }_{xx}\hfill \\ h{\sigma }_{\theta \theta }\hfill \\ h{\sigma }_{\theta x}\hfill \end{array}\right\},

where matrix [a] is the in-plane compliance matrix of laminated vessel [³⁶]. h stands for the vessel thickness and stress components are calculated as below:

(8){\sigma }_{xx}={\displaystyle \frac{Pr}{\mathrm{2}h}},{\sigma }_{\theta \theta }={\displaystyle \frac{Pr}{h}},{\sigma }_{\theta x}=\mathrm{0},

where P represents internal pressure.

Nonlinear Von-Karman strains [³⁷] are used to express strain fields. Therefore, the strain components arisen from impact loading and internal pressure are combined and described in cylindrical coordinate system as below:

(9)\begin{array}{l}{\epsilon }_{xx}={\displaystyle \frac{\partial u}{\partial x}}+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\left({\displaystyle \frac{\partial w}{\partial x}}\right)}^{\mathrm{2}}+{\epsilon }_{xx}^{\mathrm{0}},\hfill \\ {\epsilon }_{\theta \theta }={\displaystyle \frac{\partial v}{\partial \theta }}+{\displaystyle \frac{w}{R}}+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\left({\displaystyle \frac{\partial w}{\partial \theta }}\right)}^{\mathrm{2}}+{\epsilon }_{\theta \theta }^{\mathrm{0}},\hfill \\ {\epsilon }_{zz}={\displaystyle \frac{\partial w}{\partial z},}\hfill \\ {\epsilon }_{\theta z}={\displaystyle \frac{\partial v}{\partial z}}+{\displaystyle \frac{\partial w}{\partial \theta }}-{\displaystyle \frac{v}{R}},\hfill \\ {\epsilon }_{xz}={\displaystyle \frac{\partial u}{\partial z}}+{\displaystyle \frac{\partial w}{\partial x}},\hfill \\ {\epsilon }_{x\theta }={\displaystyle \frac{\partial u}{\partial \theta }}+{\displaystyle \frac{\partial v}{\partial x}}+{\displaystyle \frac{\partial w}{\partial x}}{\displaystyle \frac{\partial w}{\partial \theta }}+{\epsilon }_{x\theta }^{\mathrm{0}}\mathrm{.}\hfill \end{array}

The governing equations of motion for the present LWT can be derived using the principle of virtual displacements as below [³²]:

(10){{\displaystyle \int }}_{\mathrm{0}}^T\delta L\mathrm{d}t={{\displaystyle \int }}_{\mathrm{0}}^T[\delta K-(\delta U-\delta V\left)\right]\mathrm{d}t=\mathrm{0},

where \mathrm{\delta }U, \mathrm{\delta }V, and \mathrm{\delta }K are the virtual strain energy, virtual work done by applied forces and virtual kinetic energy, respectively. Substituting the assumed displacement field into Eq. (10), virtual energy terms are obtained. The complete form of these terms are summarized in Appendix B.

Finally, the equations of motion are obtained as below [³²]:

(11)\begin{array}{l}\mathrm{\delta }U:{\displaystyle \frac{\partial N_{xx}^I}{\partial x}}+{\displaystyle \frac{\partial N_{x\theta }^I}{\partial \theta }}+{\displaystyle \frac{\partial N_{x\theta }^I}{\partial x}}-\overline{Q}{}_{\theta }^I={\displaystyle \sum _{I,J=\mathrm{1}}^N{I}^{IJ}}{\displaystyle \frac{{\partial }^{\mathrm{2}}{U}_I}{\partial t^{\mathrm{2}}}},\hfill \\ \mathrm{\delta }V:{\displaystyle \frac{\partial N_{\theta \theta }^I}{\partial \theta }}-\overline{Q}{}_x^I+{\displaystyle \frac{Q_x^I}{R}}={\displaystyle \sum _{I,J=\mathrm{1}}^N{I}^{IJ}}{\displaystyle \frac{{\partial }^{\mathrm{2}}{V}_I}{\partial t^{\mathrm{2}}}},\hfill \\ \mathrm{\delta }W:{\displaystyle \frac{\partial Q_{\theta }^I}{\partial x}}+{\displaystyle \frac{\partial Q_x^I}{\partial \theta }}-\overline{Q}{}_z^I+N^I+q{\delta }_{IM}={\displaystyle \sum _{I,J=\mathrm{1}}^N{I}^{IJ}}{\displaystyle \frac{{\partial }^{\mathrm{2}}{W}_I}{\partial t^{\mathrm{2}}}},\hfill \end{array}

where {\delta }_{IM} is Kronecker delta, I^{IJ}={{\displaystyle \int }}_{-\frac{h}{\mathrm{2}}}^{\frac{h}{\mathrm{2}}}{\rho }_{\mathrm{0}}{\Phi }^I{\Phi }^J\mathrm{d}z and

(12)N^I={\displaystyle \sum _{I=\mathrm{1}}^N[{\displaystyle \frac{\partial }{\partial x}}\left(N_{xx}^{IJ}{\displaystyle \frac{\partial W_J}{\partial x}}+N_{x\theta }^{IJ}{\displaystyle \frac{\partial W_I}{\partial \theta }}\right){\displaystyle +{\displaystyle \frac{\partial }{\partial \theta }}\left(N_{\theta \theta }^I{\displaystyle \frac{\partial W_I}{\partial \theta }}+N_{x\theta }^{IJ}{\displaystyle \frac{\partial W_I}{\partial x}}\right)]}\mathrm{.}}

All N's and Q's are outlined in Appendix B. The constitutive equations of a lamina in cylindrical coordinate system are given as below [³²]:

(13)\left\{\begin{array}{c}{\sigma }_{xx}\hfill \\ {\sigma }_{\theta \theta }\hfill \\ {\sigma }_{zz}\hfill \\ {\sigma }_{\theta z}\hfill \\ {\sigma }_{xz}\hfill \\ {\sigma }_{x\theta }\hfill \end{array}\right\}=\left[\begin{array}{cccccc}\hfill {\overline{Q}}_{\mathrm{11}}\hfill & \hfill {\overline{Q}}_{\mathrm{12}}\hfill & \hfill {\overline{Q}}_{\mathrm{13}}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill {\overline{Q}}_{\mathrm{16}}\hfill \\ \hfill {\overline{Q}}_{\mathrm{21}}\hfill & \hfill {\overline{Q}}_{\mathrm{22}}\hfill & \hfill {\overline{Q}}_{\mathrm{23}}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill {\overline{Q}}_{\mathrm{26}}\hfill \\ \hfill {\overline{Q}}_{\mathrm{31}}\hfill & \hfill {\overline{Q}}_{\mathrm{32}}\hfill & \hfill {\overline{Q}}_{\mathrm{33}}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill {\overline{Q}}_{\mathrm{36}}\hfill \\ \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill {\overline{Q}}_{\mathrm{44}}\hfill & \hfill {\overline{Q}}_{\mathrm{45}}\hfill & \hfill \mathrm{0}\hfill \\ \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill {\overline{Q}}_{\mathrm{45}}\hfill & \hfill {\overline{Q}}_{\mathrm{55}}\hfill & \hfill \mathrm{0}\hfill \\ \hfill {\overline{Q}}_{\mathrm{16}}\hfill & \hfill {\overline{Q}}_{\mathrm{26}}\hfill & \hfill {\overline{Q}}_{\mathrm{36}}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill {\overline{Q}}_{\mathrm{66}}\hfill \end{array}\right]\left\{\begin{array}{l}{ϵ}_{xx}\hfill \\ {ϵ}_{\theta \theta }\hfill \\ {ϵ}_{zz}\hfill \\ {ϵ}_{\theta z}\hfill \\ {ϵ}_{xz}\hfill \\ {ϵ}_{x\theta }\hfill \end{array}\right\},

where {\overline{Q}}_{ij}'s are expressed in Appendix C.

Consequently, constitutive equations of laminated composites in the framework of LWT are derived and summarized in Appendix D.

For the case of simply supported element with the length of “a” and the width of “b,” the boundary conditions are expressed as:

(14)\begin{array}{l}u(x,\mathrm{0},t)=\mathrm{0},u(x,b,t)=\mathrm{0},\hspace{1em}v(\mathrm{0},\theta ,t)=\mathrm{0},\hspace{1em}v(a,\theta ,t)=\mathrm{0},\hfill \\ w(x,\mathrm{0},t)=\mathrm{0},w(x,b,t)=\mathrm{0},w(\mathrm{0},\theta ,t)=\mathrm{0},w(a,\theta ,t)=\mathrm{0}.\hfill \end{array}

The simply supported boundary conditions in Eq. (14) are satisfied by the following expansions for the displacement field [³²]:

(15)\{\begin{array}{l}{U}_I(x,\theta ,t)=u\left(x\right)\cos \left(n\theta \right)\cos \left(\omega t\right),\hfill \\ V_I(x,\theta ,t)=v\left(x\right)\sin \left(n\theta \right)\cos \left(\omega t\right),\hfill \\ W_I(x,\theta ,t)=w\left(x\right)\cos \left(n\theta \right)\cos \left(\omega t\right).\hfill \end{array}

Substituting displacement field into Eq. (11), the final governing equations for determining the displacement components are obtained as below:

(16)\begin{array}{l}{L}_{\mathrm{11}}+L_{\mathrm{12}}+L_{\mathrm{13}}=I^{IJ}{\displaystyle \frac{{\partial }^{\mathrm{2}}{U}_I}{\partial t^{\mathrm{2}}}},\hfill \\ L_{\mathrm{21}}+L_{\mathrm{22}}+L_{\mathrm{23}}=I^{IJ}{\displaystyle \frac{{\partial }^{\mathrm{2}}{V}_I}{\partial t^{\mathrm{2}}},}\hfill \\ L_{\mathrm{31}}+L_{\mathrm{32}}+L_{\mathrm{33}}=I^{IJ}{\displaystyle \frac{{\partial }^{\mathrm{2}}{W}_I}{\partial t^{\mathrm{2}}}.}\hfill \end{array}

The reflected differential operators in above equations, i.e., L_ij, are expressed in Appendix E.

The system of nonlinear equations reflected in Eq. (16) can be numerically solved using differential quadratic method (DQM) [³⁷]. The coefficients of partial derivatives are defined as below using DQM [³⁸]:

E_x^n(\ast )={\displaystyle \frac{{\partial }^{\left(n\right)}(\ast )}{\partial x^{\left(n\right)}}}={\displaystyle \sum _{I=\mathrm{1}}^N{C}^{\left(n\right)}}\left(I\right)(\ast ),

E_{\theta }^n(\ast )={\displaystyle \frac{{\partial }^{\left(n\right)}(\ast )}{\partial {\theta }^{\left(n\right)}}}={\displaystyle \sum _{J=\mathrm{1}}^N{C}^{\left(n\right)}}\left(J\right)(\ast ),

E_x^{\mathrm{1}}{E}_{\theta }^{(n-\mathrm{1})}(\ast )={\displaystyle \frac{{\partial }^{\left(n\right)}(\ast )}{\partial x\partial {\theta }^{(n-\mathrm{1})}}}={\displaystyle \sum _{I=\mathrm{1}}^N{C}^{\left(n\right)}}\left(I\right){\displaystyle \sum _{J=\mathrm{1}}^N{C}^{(n-\mathrm{1})}}\left(J\right)(\ast ),

(17)E_x^{(n-\mathrm{1})}{E}_{\theta }^{\mathrm{1}}(\ast )={\displaystyle \frac{{\partial }^{\left(n\right)}(\ast )}{\partial x^{(n-\mathrm{1})}\partial \theta }}={\displaystyle \sum _{J=\mathrm{1}}^N{C}^{(n-\mathrm{1})}}\left(J\right){\displaystyle \sum _{I=\mathrm{1}}^N{C}^{(n-\mathrm{1})}}\left(I\right)(\ast ),

where C_{ijk} are weight coefficients for DQM and defined in Appendix F.

The DQM is a numerical solution technique for initial and/or boundary problems. The method is considered as a promising alternative to the conventional numerical techniques, since higher accuracy is achieved with less computational efforts. Similar to the other numerical techniques, DQM transforms the governing differential equation into a set of analogous algebraic equations in terms of the unknown function values at the preselected sampling points. In the DQM, a partial derivative of a function with respect to a space variable at a discrete point is approximated as a weighted linear sum of the function values at all discrete points in the region of that variable [³⁸,³⁹]. The interested readers are referred to Refs. [³⁸,³⁹] for further study on DQM.

Following the explained procedure, obtained impact loading is substituted into Eq. (11) and corresponding displacement field is derived utilizing outlined DQM. Then, strain components are computed by substituting displacement field into Eq. (4). Strain components are fed into Eq. (13) and finally stress components are computed.

Model validation

Motivated by experimental study performed by Matemilola and Stronge [⁴⁰], the developed modeling procedure is executed for the same impact test on composite pressure vessel under the internal pressure. The investigated vessel is made of graphite/epoxy with diameter of 300 mm and lay-up configuration of [90°/10°/63°/−10°/−90°/10°/63°/−10°/−52°/10°]. The thickness of investigated vessel is 10 mm and the thickness of each layer is taken as 1 mm. The vessel was first subjected to an internal pressure of 400 bar and then a low-velocity impact is applied to its cylindrical section via a dome-nosed impactor with 20-mm nose diameter [⁴⁰]. Mechanical properties of constitutive layers based on Ref. [40] are as follow: Elastic modulus are E_{xx}=\mathrm{60.28}\hspace{0.17em}\mathrm{G}\mathrm{P}\mathrm{a}, E_{yy},\hspace{0.17em}{E}_{zz}=\mathrm{8.963}\hspace{0.17em}\mathrm{G}\mathrm{P}\mathrm{a}, G_{xy},\hspace{0.17em}{G}_{xz}=\mathrm{1.43} GPa, and G_{yz}=\mathrm{2.3}\mathrm{G}\mathrm{P}\mathrm{a}. Poission ratios are {\nu }_{xy},\hspace{0.17em}{\nu }_{xz}=\mathrm{0.174} and {\nu }_{yz}=\mathrm{0.278}. Longitudinal tensile and compressive strength are X_T=\mathrm{837}\mathrm{M}\mathrm{P}\mathrm{a} and X_C=\mathrm{414}\mathrm{M}\mathrm{P}\mathrm{a}, respectively. Transverse tensile and compressive strength are Y_T=\mathrm{25.8}\mathrm{M}\mathrm{P}\mathrm{a} and Y_C=\mathrm{100.2}\mathrm{M}\mathrm{P}\mathrm{a}. Out-of-plane shear strengths are S_{xz}=\mathrm{44.24}\mathrm{M}\mathrm{P}\mathrm{a} and S_{yz}=\mathrm{44.24} MPa and in-plane shear strength is S_{xy}=\mathrm{44.24}\mathrm{M}\mathrm{P}\mathrm{a}. Fiber fracture energies are {\Gamma }_{\mathrm{11}}^T=133 N/mm and {\Gamma }_{\mathrm{11}}^C=40 N/mm, respectively. Matrix fracture energies are {\Gamma }_{\mathrm{22}}^T=\mathrm{0.6}N/mm and {\Gamma }_{\mathrm{22}}^C=\mathrm{2.1}N/mm, respectively.

Following the developed modeling procedure in preceding section, the maximum deflections of the vessel subjected to impact loading are obtained for both cases of with and without internal pressure. The obtained results in comparison with experimental observations are presented in Table 1.

Owning to the minimizing energy through Hamiltonian approach (Eq. (10)), the estimated results are slightly overestimated in comparison with experimental observations. The small different between experimental observation and theoretical modeling established our confidence toward the proper modeling procedure. Therefore, the investigated vessel is exposed to various impact loadings with the velocities of 5, 7 and 10 m/s for more detailed analysis. Contact forces and experienced deflections through the thickness during the impact time with internal pressure of 400 bar are presented in Fig. 3.

As it is evident from Fig. 3, increasing the velocity of impact, both maximum deflection and contact force increase, too. Figure 4 demonstrates variations of strain components during the impact time for the pressurized vessel.

From Fig. 4, it can be observed that in comparison with shear strain component, the other strain components are negligible.

Failure analysis

For the purpose of identifying failure occurrence, proper failure criteria are required to be chosen. 3-D Hashin failure criteria [⁴¹] are utilized for evaluating the occurrence of in-plane failure:

(18)

where X_{\mathrm{T}} and X_{\mathrm{C}} are longitudinal tensile and compressive strength, respectively. Y_{\mathrm{T}} and Y_{\mathrm{C}} are transverse tensile and compressive strength and S_{yz} and S_{xz} indicate out-of-plane shear strength components and S_{xy}is in-plane shear strength. The first two criteria (i.e., “a” and “b”) imply on catastrophic failure modes.

For the case of evaluating out-of-plane failure, following quadratic criterion is employed [⁴²]:

(19){\left({\displaystyle \frac{{\sigma }_{zz}}{Y_{\mathrm{c}}}}\right)}^{\mathrm{2}}+{\left({\displaystyle \frac{{\tau }_{xz}}{S_{xz}}}\right)}^{\mathrm{2}}+{\left({\displaystyle \frac{{\tau }_{yz}}{S_{yz}}}\right)}^{\mathrm{2}}=\mathrm{1}.

Obtained stress components in Eq. (13) are first converted to on-axis stress components for each layer and then calculated on-axis stress components are fed into 3-D Hashin failure criteria to evaluate the occurrence of failure.

Progressive damage modeling

Since, the vessel is subjected to internal pressure and impact loading at the same time, progressive damage modeling is employed for predicting the ultimate failure. The strategy of analyzing failure in pressurized composite pressure vessels subjected to low-velocity impact is presented in Fig. 5 in the framework of progressive damage modeling.

Following to the presented flowchart in Fig. 5, after applying impact to the pre-stressed vessel due to internal pressure, the occurrence of failure is examined. If no failure is observed, the impact energy increases gradually until a failure is observed. In case of delamination failure, the modeling is terminated. Once in-plane failure is identified, mechanical properties of failed ply is degraded on the basis of continuum damage mechanics. In the context of continuum damage mechanic, the failed ply is replaced with an intact ply with reduced level of mechanical properties. Consequently, mechanical properties of the model are updated and the failure is evaluated again using updated status of stress components. This procedure continues until either out-of-plane failure is identified or no more in-plane failure is experienced after degrading mechanical properties of previously failed plies/ply. The latter is also pertinent to the situation when all layers are failed and thus vessel cannot sustain any more loading.

The degradation of mechanical properties is accomplished using damage variables as below [⁴³]:

(20)\left[C_d\right]=\left[\begin{array}{cccccc}\hfill \mathrm{1}-d_{\mathrm{11}}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill \\ \hfill \mathrm{0}\hfill & \hfill \mathrm{1}-d_{\mathrm{22}}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill \\ \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{1}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill \\ \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{1}-d_s\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill \\ \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{1}-d_s\hfill & \hfill \mathrm{0}\hfill \\ \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{0}\hfill & \hfill \mathrm{1}-d_s\hfill \end{array}\right]\left[C\right],

where [C] and [C_d] represent stiffness matrices of undamaged and damaged plies, respectively. d_ij parameters are damage variables defined as below [⁴⁴,⁴⁵]:

(21)\begin{array}{l}{d}_{\mathrm{11}}=\mathrm{1}-\left(\mathrm{1}-d_{\mathrm{11}}^{\mathrm{T}}\right)\left(\mathrm{1}-d_{\mathrm{11}}^{\mathrm{C}}\right),\hfill \\ d_{\mathrm{22}}=\mathrm{1}-\left(\mathrm{1}-d_{\mathrm{22}}^{\mathrm{T}}\right)\left(\mathrm{1}-d_{\mathrm{22}}^{\mathrm{C}}\right),\hfill \\ d_s=\mathrm{1}-\left(\mathrm{1}-d_{\mathrm{11}}^{\mathrm{T}}\right)\left(\mathrm{1}-d_{\mathrm{11}}^{\mathrm{C}}\right)\left(\mathrm{1}-s_{\mathrm{m}\mathrm{t}}{d}_{\mathrm{22}}^{\mathrm{T}}\right)\left(\mathrm{1}-s_{\mathrm{m}\mathrm{c}}{d}_{\mathrm{22}}^{\mathrm{C}}\right)\mathrm{.}\hfill \end{array}

The superscripts “T” and “C” in the damage variables denote the tension and compression, respectively. s_{\mathrm{m}\mathrm{t}} and s_{\mathrm{m}\mathrm{c}} are coefficients to control the loss of shear stiffness due to matrix tension and compression and are assumed as 0.9 and 0.5, respectively [⁴⁴]. d_{\mathrm{11}}^{\mathrm{T}} and d_{\mathrm{11}}^{\mathrm{C}} in Eq. (21) are given by [⁴⁴]:

(22)d_{\mathrm{11}}^{\mathrm{T}\left(\mathrm{C}\right)}={\displaystyle \frac{{\epsilon }_{f,\mathrm{1}}^{\mathrm{T}\left(\mathrm{C}\right)}}{{\epsilon }_{f,\mathrm{1}}^{\mathrm{T}\left(\mathrm{C}\right)}-{\epsilon }_{\mathrm{0},\mathrm{1}}^{\mathrm{T}\left(\mathrm{C}\right)}}}\left(\mathrm{1}-{\displaystyle \frac{{\epsilon }_{\mathrm{0},\mathrm{1}}^{\mathrm{T}\left(\mathrm{C}\right)}}{{\epsilon }_{\mathrm{11}}}}\right),

(23){\epsilon }_{\mathrm{0},\mathrm{1}}^{\mathrm{T}}={\displaystyle \frac{X_{\mathrm{T}}}{E_x}},{\epsilon }_{\mathrm{0},\mathrm{1}}^{\mathrm{C}}={\displaystyle \frac{X_{\mathrm{C}}}{E_x}},

where {\epsilon }_{\mathrm{11}} is the corresponding component of Green-Lagrange strain tensor, {\epsilon }_{\mathrm{0},\mathrm{1}}^{\mathrm{C}} and {\epsilon }_{\mathrm{0},\mathrm{1}}^{\mathrm{T}} are the corresponding initial failure strains. {\epsilon }_{f,\mathrm{1}}^{\mathrm{T}\left(\mathrm{C}\right)} is the strain when corresponding damage variable reaches one and calculated as below [⁴⁴]:

(24){\epsilon }_{f,\mathrm{1}}^{\mathrm{T}\left(\mathrm{C}\right)}={\displaystyle \frac{\mathrm{2}{\Gamma }_{\mathrm{11}}^{\mathrm{T}\left(\mathrm{C}\right)}}{X_{\mathrm{T}\left(\mathrm{C}\right)}}},

where {\Gamma }_{\mathrm{11}}^{\mathrm{T}\left(\mathrm{C}\right)}is the fiber fracture energy and determined experimentally [⁴⁴].

Damage variables for matrix tension/compression (d_{\mathrm{22}}^{\mathrm{T}} and d_{\mathrm{22}}^{\mathrm{C}}) are expressed as below [⁴³]:

(25)d_{\mathrm{22}}^{\mathrm{T}\left(\mathrm{C}\right)}={\displaystyle \frac{{\epsilon }_{f,\mathrm{2}}^{\mathrm{T}\left(\mathrm{C}\right)}}{{\epsilon }_{f,\mathrm{2}}^{\mathrm{T}\left(\mathrm{C}\right)}-{\epsilon }_{\mathrm{0},\mathrm{2}}^{\mathrm{T}\left(\mathrm{C}\right)}}}\left(\mathrm{1}-{\displaystyle \frac{{\epsilon }_{\mathrm{0},\mathrm{2}}^{\mathrm{T}\left(\mathrm{C}\right)}}{{\epsilon }_{\mathrm{22}}}}\right),

where {\epsilon }_{\mathrm{22}} is the corresponding component of Green-Lagrange strain tensor. {\epsilon }_{\mathrm{0},\mathrm{2}}^{\mathrm{T}\left(\mathrm{C}\right)}, {\epsilon }_{f,\mathrm{2}}^{\mathrm{T}\left(\mathrm{C}\right)} are calculated as below [⁴⁴]:

(26){\epsilon }_{\mathrm{0},\mathrm{2}}^{\mathrm{T}\left(\mathrm{C}\right)}={\displaystyle \frac{Y_{\mathrm{T}\left(\mathrm{C}\right)}}{E_y}}\mathrm{;}{\epsilon }_{f,\mathrm{2}}^{\mathrm{T}\left(\mathrm{C}\right)}={\displaystyle \frac{\mathrm{2}{\Gamma }_{\mathrm{22}}^{\mathrm{T}\left(\mathrm{C}\right)}}{Y_{\mathrm{T}\left(\mathrm{C}\right)}l}},

where {\Gamma }_{\mathrm{22}}^{\mathrm{T}\left(\mathrm{C}\right)} is the matrix fracture energy.

Although some advanced techniques are available in literature for predicting crack initiation and propagation in thin shells exposed to dynamic loading [⁴⁶–⁵⁰], the above mentioned failure criteria are employed here in the context of continuum damage mechanics.

Following the explained procedure of progressive damage modeling, the same vessel introduced in Section 3 is subjected to an impactor with various velocities resembling 110, 150, and 270 J energy while it undergoes an internal pressure of 400 bar. The detailed information on experienced failure patterns in various layers of investigated vessels at applied impact load levels is summarized in Table 2.

One of the main advantages of developed technique is placed behind this fact that the sequence of failure can be monitored and analyzed in details. Applying an impact with 110 J energy to the pressurized vessel, it is observed that some layers experience matrix cracking and in-plane shear failure modes while neither fiber failure nor delamination is experienced. Increasing the impact load to 150 J, first matrix failure and in-plane shear failure happen in some layers and then first three layers experience catastrophic mode of fiber failure. Finally, at this level of energy, delamination is experienced before the occurrence of any in-plane failure mode in other layers. At the impact load level of 270 J, all layers experience matrix cracking and in-plane shear failure and then fiber failure is experienced in all layers. As it was explained before, the solution is truncated at this moment, since all layers have already experienced catastrophic failure and thus the internal pressure cannot be accommodated anymore. It can be seen from the results that failure does not occur in layers with the same sequence of lay-up. The observed failure modes are in a good consistency with experimental observations [⁴⁰].

Conclusions

In this research, a theoretical framework is established on the basis of LWT to extract all six components of stress and strain fields in a pressurized composite vessel subjected to low-velocity impact. For this purpose, the strain-displacement field represented by Von-Karman is utilized to determine the strains and then stress components are characterized using constitutive equations. Having in hand stress fields, 3D failure criteria are employed to predict in-plane failure induced by impact loading. Simultaneously, the occurrence of out-of-plane failure in the form of delamination is also examined. After occurrence of in-plane failure, mechanical properties of the failed ply are degraded with respect to the experienced failure mode. Therefore, failed ply is replaced with an intact ply with weaker mechanical properties in the context of progressive damage modeling. In other words, LWT theory is employed in combination with progressive damage modeling to evaluate the damage evolution. The obtained results are validated with available experimental observations in literature. It is realized that the sequence of experiencing either in-plane failure modes or delamination in a pressurized composite vessel depends on the imposed impact energy level.

While numerical modeling through commercial finite element packages is widely used by industrial centers to evaluate the performance of a composite pressure vessel under impact and internal pressure loadings, developed theoretical platform can be conveniently used in the early stages of designing a composite pressure vessel with considerable less computational efforts.

It is noteworthy to mention that experimental procedure as a qualification test on new design schemes is a mandatory task and developed modeling technique does not waive its necessity. However, it can be employed in the design stage as an engineering tool for establishing the confidence of the new design scheme prior to conducting the prototype test.

References

[1] Velosa J C, Nunes J P, Antunes P J, Silva J F, Marques A M. Development of a new generation of filament wound composite pressure cylinders. Composites Science and Technology, 2009, 69(9): 1348–1353

[2] ISO 11439:2000. International Standard, Gas Cylinders-High Pressure Cylinders for the on-Board Storage of Natural Gas as a Fuel for Automotive Vehicles. International Organization for Standardization,2000

[3] EN 14427. Transportable Refillable Fully Wrapped Composite Cylinders for Liquefied Petroleum Gases (LPG). Des Construction, 2004

[4] Alderson K L, Evans K E. Failure mechanisms during the transverse loading of filament-wound pipes under static and low velocity impact conditions. Composites, 1992, 23(3): 167–173

[5] Gning P B, Tarfaoui M, Collombet F, Riou L, Davies P. Damage development in thick composite tubes under impact loading and influence on implosion pressure: Experimental observations. Composites. Part B, Engineering, 2005, 36(4): 306–318

[6] Deniz M E, Karakuzu R, Sari M, Icten B M. On the residual compressive strength of the glass-epoxy tubes subjected to transverse impact loading. Journal of Composite Materials, 2012, 46(6): 737–745

[7] Krishnamurthy K S, Mahajan P, Mittal R K. A parametric study of the impact response and damage of laminated cylindrical composite shells. Composites Science and Technology, 2001, 61(12): 1655–1669

[8] Zhang C, Ren M, Zhao W, Chen H. Delamination prediction of composite filament wound vessel with metal liner under low velocity impact. Composite Structures, 2006, 75(1–4): 387–392

[9] Demir I, Sayman O, Dogan A, Arikan V, Arman Y. The effects of repeated transverse impact load on the burst pressure of composite pressure vessel. Composites. Part B, Engineering, 2015, 68: 121–125

[10] Han M G, Chang S H. Evaluation of structural integrity of type-III hydrogen pressure vessel under low-velocity car-to-car collision using finite element analysis. Composite Structures, 2016, 148: 198–206

[11] Gemi L. Investigation of the effect of stacking sequence on low velocity impact response and damage formation in hybrid composite pipes under internal pressure. A comparative study. Composites. Part B, Engineering, 2018, 153: 217–232

[12] Wu Q, Chen X, Fan Z, Jiang Y, Zhang X, Nie D. Damage behavior of filament-wound composite cylinder under impact by flat-ended impactor. In: ASME 2018 Pressure Vessels and Piping Conference. American Society of Mechanical Engineers, 2018

[13] Abrate S, Ferrero J F, Navarro P. Cohesive zone models and impact damage predictions for composite structures. Meccanica, 2015, 50(10): 2587–2620

[14] Najafi F, Shojaeefard M H, Saeidi Googarchin H. Low-velocity impact response of functionally graded doubly curved panels with winkler-pasternak elastic foundation: An analytical approach. Composite Structures, 2017, 162: 351–364

[15] Khalili S M R, Ardali A. Low-velocity impact response of doubly curved symmetric cross-ply laminated panel with embedded SMA wires. Composite Structures, 2013, 105: 216–226

[16] Arachchige B, Ghasemnejad H, Augousti A T. Theoretical approach to predict transverse impact response of variable-stiffness curved composite plates. Composites. Part B, Engineering, 2016, 89: 34–43

[17] Hwang T K, Hong C S, Kim C G. Probabilistic deformation and strength prediction for a filament wound pressure vessel. Composites B. Engineering (London), 2003, 34(5): 481–497

[18] Hwang T K, Hong C S, Kim C G. Size effect on the fiber strength of composite pressure vessels. Composite Structures, 2003, 59(4): 489–498

[19] Sun X K, Du S Y, Wang G D. Bursting problem of filament wound composite pressure vessels. International Journal of Pressure Vessels and Piping, 1999, 76(1): 55–59

[20] Onder A, Sayman O, Dogan T, Tarakcioglu N. Burst failure load of composite pressure vessels. Composite Structures, 2009, 89(1): 159–166

[21] Liu P F, Zheng J Y. Progressive failure analysis of carbon fiber/epoxy composite laminates using continuum damage mechanics. Materials Science and Engineering A, 2008, 485(1–2): 711–717

[22] Xu P, Zheng J Y, Liu P F. Finite element analysis of burst pressure of composite hydrogen storage vessels. Materials & Design, 2009, 30(7): 2295–2301

[23] Wang L, Zheng C, Luo H, Wei S, Wei Z. Continuum damage modeling and progressive failure analysis of carbon fiber/epoxy composite pressure vessel. Composite Structures, 2015, 134: 475–482

[24] Leh D, Saffré P, Francescato P, Arrieux R, Villalonga S. A progressive failure analysis of a 700-bar type IV hydrogen composite pressure vessel. International Journal of Hydrogen Energy, 2015, 40(38): 13206–13214

[25] Berro Ramirez J P, Halm D, Grandidier J C, Villalonga S. A fixed directions damage model for composite materials dedicated to hyperbaric type IV hydrogen storage vessel—Part I. International Journal of Hydrogen Energy, 2015, 40(38): 13165–13173

[26] Berro Ramirez J P, Halm D, Grandidier J C, Villalonga S. A fixed directions damage model for composite materials dedicated to hyperbaric type IV hydrogen storage vessel—Part II. International Journal of Hydrogen Energy, 2015, 40(38): 13174–13182

[27] Berro Ramirez J P, Halm D, Grandidier J C, Villalonga S, Nony F. 700 bar type IV high pressure hydrogen storage vessel burst— Simulation and experimental validation. International Journal of Hydrogen Energy, 2015, 40(38): 13183–13192

[28] Gentilleau B, Villalonga S, Nony F, Galiano H. A probabilistic damage behavior law for composite material dedicated to composite pressure vessel. International Journal of Hydrogen Energy, 2015, 40(38): 13160–13164

[29] Gentilleau B, Touchard F, Grandidier J C. Numerical study of influence of temperature and matrix cracking on type IV hydrogen high pressure storage vessel behavior. Composite Structures, 2014, 111: 98–110

[30] Rafiee R, Torabi M A. Stochastic prediction of burst pressure in composite pressure vessels. Composite Structures, 2018, 185: 573–583

[31] Rafiee R, Torabi M A, Maleki S. Investigating structural failure of a filament-wound composite tube subjected to internal pressure: Experimental and theoretical evaluation. Polymer Testing, 2018, 67: 322–330

[32] Reddy J N. Mechanics of Laminated Composite Plates and Shells: Theory and Analysis. New York: CRC press, 2004

[33] Rafiee R, Ghorbanhosseini A, Rezaee Sh. Theoretical and numerical analyses of composite cylinders subjected to the low velocity impact. Composite Structures, 2019, 226: 111230

[34] Ramkumar R L, Thakar Y R. Dynamic response of curved laminated plates subjected to low velocity impact. Journal of Engineering Materials and Technology, 1987, 109(1): 67–71

[35] Gong S W, Shim V P W, Toh S L. Impact response of laminated shells with orthogonal curvatures. Composites Engineering, 1995, 5(3): 257–275

[36] Tsai S W, Hahn H T. Introduction to composite materials. Lancaster: Technomic Pub, 1980

[37] Hashiguchi K, Yamakawa Y. Introduction to finite strain theory for continuum elasto-plasticity. New Delhi: John Wiley & Sons, 2012

[38] Civalek O, Ulker M. HDQ-FD integrated methodology for nonlinear static and dynamic response of doubly curved shallow shells. Structural Engineering and Mechanics, 2005, 19(5): 535–550

[39] Bert C W, Malik M. Differntial quadrature method in computational mechanics: A review. Applied Mechanics Reviews, 1996, 49(1): 1–28

[40] Matemilola S A, Stronge W J. Low-speed impact damage in filament-wound CFRP composite pressure vessels. Journal of Pressure Vessel Technology, 1997, 119(4): 435–443

[41] Cheng X, li Z. Damage progressive model of compression of composite laminates after low velocity impact. Applied Mathematics and Mechanics, 2005, 26(5): 618–626

[42] Brewer J C, Lagace P A. Quadratic stress criterion for initiation of delamination. Journal of Composite Materials, 1988, 22(12): 1141–1155

[43] Donadon M V, Iannucci L, Falzon B G, Hodgkinson J M, de Almeida S F M. A progressive failure model for composite laminates subjected to low velocity impact damage. Computers & Structures, 2008, 86(11–12): 1232–1252

[44] Liu P F, Liao B B, Jia L Y, Peng X Q. Finite element analysis of dynamic progressive failure for carbon fiber composite laminates under low velocity impact. Composite Structures, 2016, 149: 408–422

[45] Pederson J. Finite element analysis of carbon fiber composite ripping using ABAQUS. Thesis for the Master’s Degree. Clemson, SC: Clemson University, 2008

[46] Rabczuk T, Zi G, Bordas St, 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

[47] 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

[48] 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

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

[50] 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

RIGHTS & PERMISSIONS

Higher Education Press