Multiscale Failure Analysis and Micromechanical Fem Modeling of Type IV Composite Pressure Vessel Considering Different Fibre Shapes and Dehomogenisation

Composite pressure vessels are widely used in transportation, automotive, and aerospace industries. They are most used for storing gases, which are various types of fuel for internal combustion and electric engines stored under high pressure (1). Therefore, the main task of these structures is to ensure their adequate strength, stiffness and resistance to various external environmental factors such as different temperatures, impacts, etc. When designing structures using composites, it is important to select the appropriate material structure (2). The use of multiscale modelling makes it possible to build a model of the material on a microscopic scale and, based on this, to determine the constants of the resulting composite. Figure 1 describes the scales of analysis for composite materials (3).

Fig. 1.

Composite analysis scales

On a microscopic scale, materials can be chosen, their volume ratio, fibre diameter and the way the fibers are arranged. Depending on the materials used and their proportions, different values of material constants can be obtained. From the microscopic scale one moves to the macro scale, which represents the layer scale, and then to the laminate scale which represents the sequence of layers (4). The last is the component scale, for example the composite overlap in a Type IV pressure vessel. In the case of hydrogen storage vessels, a high operating pressure of 70 MPa is required. Operation under such conditions is provided by this type of tanks design, in which the strength is provided by a carbon fibre composite, and the internal gas diffusion barrier is provided by a plastic liner. The use of a plastic liner allows the weight of such a tank to be reduced compared to the Type III design, where the liner is made of aluminium [5, 6]. Figure 2 shows an example of a Type IV composite pressure vessel. For the composite overlap in these structures, the method of winding the fibre band is used. As a result of using this method, in the dome section of the tank, each layer has a variable thickness and a variable winding angle of the fibre tape, which is an important aspect in the design process [7, 8].

Fig.2.

Construction of composite pressure vessel type IV

Damage analysis of composite pressure vessels is complex, requiring a multiscale approach to accurately capture various damage mechanisms. This research presents a comprehensive framework for damage analysis and finite element micromechanical modeling of Type IV composite pressure vessels, focusing on stress state dehomogenization and damage initiation. Numerical studies of composite pressure vessels necessitate the analysis of strain in homogeneous materials. Multiscale modeling approaches have been employed by numerous researchers, including Sadik L. Omairey (3), Xin Liu (4), to determine the mechanical properties of composite materials based on the properties of their constituent elements. This methodology has been applied to composite vessels by Song Lin (9), P.F. Liu (10) and Nan Zhang (11). In the studies reported here, this approach was implemented using periodic boundary conditions in Abaqus. In the present study, the Material Designer plug-in available in Ansys Workbench and a second novel approach using the Mechanical APDL programming language were used to create equations satisfying periodic boundary conditions to map these conditions for each cubic cell size under different loading conditions. The micromechanical model was constructed using a fibre model based on a scan of its actual structure, in addition to two simplified fibre models with elliptical and circular cross-sections. The resulting material constants were then compared, and the effect of fibre shape on their values was assessed. It was determined that the developed equations will be used in future work for numerical simulations of the resin polymerization process in fiber composites to determine the residual stresses then arising. A detailed process for building a braided composite is presented, using Wang's method, which was also used by Jiao Lin (12), Qian Zhang (13) and Hui Wang (14). A detailed analysis of the composite's strain was conducted using the criteria of Puck, Hashin, and Tsai Wu. The results were compared and the stress distributions in other elements of the tank were shown, which had not been done by previous researchers except for Qian Zhang (13).

The approach used in this study will allow a more accurate prediction of the failure of composites in structures and show the impact that the fibre geometry used has. Novel aspect of the work is the application of a dehomogenisation process to a composite pressure vessel, considering three different fibre geometries, including their actual structure. Dehomogenisation involves transferring the deformation state from the macroscopic scale to the microscopic scale, enabling stress analysis of the constituent materials. A strain tensor using equations describing the three-dimensional strain state obtained from the most critically stressed points was applied to representative volumetric elements. A comparative analysis of the results was performed to assess the influence of fibre shape and orientation on the reinforcement in the composite structure. The equations were created using the Mechanical APDL language and allowed a numerical analysis of the deformation state in the component materials based on the obtained stress distributions in the RVE. This approach makes it possible to observe the behaviour of the fibres and resin with greater accuracy, where and how they will fail or break, which the macro-scale model does not allow. In addition, a comparison was made between the obtained stress distributions of the component materials and the strength criteria of the composites.

The micromechanical model of a composite material made of carbon fibre and epoxy resin (15) was created in Ansys Workbench in two ways, the first using the Material Designer plug-in and the second requiring much more work to create the model using Mechanical APDL language. The composite material was created using Toray T700 carbon fibre and epoxy resin (9). Orthotropic material properties were assigned to carbon fibre and isotropic properties to the resin, as shown in Table 1.

The representative volume element has a cube shape, where the ratio of the materials used is 60% fibre and 40% resin. Three different fibre geometries were used to build a representative volume element. The first uses a circular fibre diameter of 7μm, which is what the fibre manufacturer states. The second geometry is based on a simplified elliptical shape, which is a simplified scan of the actual fibre. The last geometry reproduces the actual scan of the fibre, which is shown in figure 3. Using material designer, only an RVE model with a circular fibre was created, and using Ansys mechanical, models using all three proposed fibre geometries were created. In total, therefore, results from four modelling approaches were obtained and are summarised in Table 6. Figures 4 and 5 show geometric models of RVE with three fibre shapes and a finite element mesh.

Fig. 3.

Scan of the actual fibre structure: a) multifiber view, b) single fibre zoom

Fig. 4.

a) RVE with circular fibre: b) RVE with elliptical fibre

Fig. 5.

RVE with actual fibre: a) geometry, b) mesh

To create a finite element mesh on the RVE in Ansys Mechanical, it is necessary to ensure that there is an equal number of nodes on each opposite face. This is because if the number of nodes is not equal, it will not be possible to solve the model correctly [3, 16, 17, 18], as each node must have a counterpart. The formulation of the equations that describe the dependence of the displacements, the behaviour of the nodes on the individual walls, edges and vertices, and the representation of the periodic boundary conditions was conducted in the mechanical language APDL and is presented in Tables 2-5. Each equation establishes a linkage between the displacements of two nodes located on opposite walls.

Tab. 2.

Surfaces	Edges	Vertices
ux5678−ux1234=0u_x^{5678} - u_x^{1234} = 0	ux56−ux12=0,ux78−ux34=0,ux58−ux14=0,ux67−ux23=0\matrix{ {u_x^{56} - u_x^{12} = 0,} \hfill & {u_x^{78} - u_x^{34} = 0,} \hfill \cr {u_x^{58} - u_x^{14} = 0,} \hfill & {u_x^{67} - u_x^{23} = 0} \hfill \cr }	ux5−ux1=0,ux6−ux2=0,ux8−ux4=0,ux7−ux3=0\matrix{ {u_x^5 - u_x^1 = 0,} \hfill & {u_x^6 - u_x^2 = 0,} \hfill \cr {u_x^8 - u_x^4 = 0,} \hfill & {u_x^7 - u_x^3 = 0} \hfill \cr }
uy2367−uy1458=0u_y^{2367} - u_y^{1458} = 0	uy23−uy14=0,uy67−uy58=0,uy26−uy15=0,uy37−uy48=0\matrix{ {u_y^{23} - u_y^{14} = 0,} \hfill & {u_y^{67} - u_y^{58} = 0,} \hfill \cr {u_y^{26} - u_y^{15} = 0,} \hfill & {u_y^{37} - u_y^{48} = 0} \hfill \cr }	uy2−uy1=0,uy6−uy5=0,uy3−uy4=0,uy7−uy8=0\matrix{ {u_y^2 - u_y^1 = 0,} \hfill & {u_y^6 - u_y^5 = 0,} \hfill \cr {u_y^3 - u_y^4 = 0,} \hfill & {u_y^7 - u_y^8 = 0} \hfill \cr }
uz3478−uz1256=0u_z^{3478} - u_z^{1256} = 0	uz78−uz56=0,uz34−uz12=0,uz48−uz15=0,uz37−uz26=0\matrix{ {u_z^{78} - u_z^{56} = 0,} & {u_z^{34} - u_z^{12} = 0,} \cr {u_z^{48} - u_z^{15} = 0,} & {u_z^{37} - u_z^{26} = 0} \cr }	uz4−uz1=0,uz8−uz5=0,uz7−uz6=0,uz3−uz2=0]\matrix{ {u_z^4 - u_z^1 = 0,} \hfill & {u_z^8 - u_z^5 = 0,} \hfill \cr {u_z^7 - u_z^6 = 0,} \hfill & {u_z^3 - u_z^2 = 0} \hfill \cr }

Tab. 3.

Surfaces	Edges	Vertices
ux3478−ux1256=0u_x^{3478} - u_x^{1256} = 0	ux34−ux12=0,ux78−ux56=0\matrix{ {u_x^{34} - u_x^{12} = 0,} \hfill & {u_x^{78} - u_x^{56} = 0} \hfill \cr }	ux4−ux1=0u_x^4 - u_x^1 = 0
uy2367−uy1458=0u_y^{2367} - u_y^{1458} = 0	uy37−uy48=0,uy26−uy15=0uy23−uy14=0,uy67−uy58=0\matrix{ {u_y^{37} - u_y^{48} = 0,} \hfill & {u_y^{26} - u_y^{15} = 0} \hfill \cr {u_y^{23} - u_y^{14} = 0,} \hfill & {u_y^{67} - u_y^{58} = 0} \hfill \cr }	uy3−uy4=0,uy2−uy1=0uy6−uy5=0,uy7−uy8=0\matrix{ {u_y^3 - u_y^4 = 0,} \hfill & {u_y^2 - u_y^1 = 0} \hfill \cr {u_y^6 - u_y^5 = 0,} \hfill & {u_y^7 - u_y^8 = 0} \hfill \cr }
uz5678−uz1234=0u_z^{5678} - u_z^{1234} = 0	uz78−uz34=0,uz56−uz12=0\matrix{ {u_z^{78} - u_z^{34} = 0,} \hfill & {u_z^{56} - u_z^{12} = 0} \hfill \cr }	uz5−uz1=0u_z^5 - u_z^1 = 0

Tab. 4.

Surfaces	Edges	Vertices
ux2367−ux1458=0u_x^{2367} - u_x^{1458} = 0	ux23−ux14=0,ux67−ux58=0\matrix{ {u_x^{23} - u_x^{14} = 0,} \hfill & {u_x^{67} - u_x^{58} = 0} \hfill \cr }	ux2−ux1=0u_x^2 - u_x^1 = 0
uy5678−uy1234=0u_y^{5678} - u_y^{1234} = 0	uy67−uy23=0,uy58−uy14=0\matrix{ {u_y^{67} - u_y^{23} = 0,} \hfill & {u_y^{58} - u_y^{14} = 0} \hfill \cr }	uy5−uy1=0u_y^5 - u_y^1 = 0
uz3478−uz1256=0u_z^{3478} - u_z^{1256} = 0	uz48−uz15=0,uz37−uz26=0uz34−uz12=0,uz78−uz56=0\matrix{ {u_z^{48} - u_z^{15} = 0,} \hfill & {u_z^{37} - u_z^{26} = 0} \hfill \cr {u_z^{34} - u_z^{12} = 0,} \hfill & {u_z^{78} - u_z^{56} = 0} \hfill \cr }	uz4−uz1=0,uz3−uz2=0uz7−uz6=0,uz8−uz5=0\matrix{ {u_z^4 - u_z^1 = 0,} \hfill & {u_z^3 - u_z^2 = 0} \hfill \cr {u_z^7 - u_z^6 = 0,} \hfill & {u_z^8 - u_z^5 = 0} \hfill \cr }

Tab. 5.

Surfaces	Edges	Vertices
ux5678−ux1234=0u_x^{5678} - u_x^{1234} = 0	ux58−ux14=0,ux67−ux23=0ux78−ux34=0,ux56−ux12=0\matrix{ {u_x^{58} - u_x^{14} = 0,} \hfill & {u_x^{67} - u_x^{23} = 0} \hfill \cr {u_x^{78} - u_x^{34} = 0,} \hfill & {u_x^{56} - u_x^{12} = 0} \hfill \cr }	ux5−ux1=0,ux6−ux2=0ux8−ux4=0,ux7−ux3=0\matrix{ {u_x^5 - u_x^1 = 0,} \hfill & {u_x^6 - u_x^2 = 0} \hfill \cr {u_x^8 - u_x^4 = 0,} \hfill & {u_x^7 - u_x^3 = 0} \hfill \cr }
uy3478−uy1256=0u_y^{3478} - u_y^{1256} = 0	uy48−uy15=0,uy37−uy26=0\matrix{ {u_y^{48} - u_y^{15} = 0,} \hfill & {u_y^{37} - u_y^{26} = 0} \hfill \cr }	uy4−uy1=0u_y^4 - u_y^1 = 0
uz2367−uz1458=0u_z^{2367} - u_z^{1458} = 0	uz26−uz15=0,uz37−uz48=0\matrix{ {u_z^{26} - u_z^{15} = 0,} \hfill & {u_z^{37} - u_z^{48} = 0} \hfill \cr }	uz2−uz1=0u_z^2 - u_z^1 = 0

Fig. 6.

Numbering of RVE vertices

The determination of all nine material constants was achieved by setting six loading conditions, namely tension in the X, Y, Z directions and shear in the XY, XZ and YZ planes. The equations used were divided into three groups, categorised by faces, edges and vertices. Each vertex of the cubic RVE model was assigned a number, as illustrated in Figure 4, and the displacement equations [16, 17, 18, 19] were defined using these numbers of vertices.

The displacement conditions were set using the displacement function and the displacement has the same value for each load case. For example, for tension in the X-direction, the displacement was set on the front surface on X-direction, for shear in the XZ-plane the displacement was set in the X-direction on the top Z-sur-face and in the Z-direction on the front X-surface. Based on the specified loads using periodic boundary conditions, distributions of resultant displacements (Fig.7) and Huber von Mises reduced stresses (Fig. 8). The distributions are shown for RVE with circular cross-section fibre.

Fig. 7.

Displacement distributions of RVE: a) tensile in X-direction, b) tensile in Y-direction, c) tensile in Z-direction

Fig. 8.

Displacement distributions of RVE: a) shear in XY-plane, b) shear in XZ-plane, c) shear in YZ-plane

Fig. 9.

Equivalent stress distributions of RVE: a) tensile X-direction, b) tensile Y-direction, c) tensile Z-direction

Fig. 10.

Equivalent stress distributions of RVE: a) shear in XY-plane, b) shear in XZ-plane, c) shear in YZ-plane

From the displacement distributions it is possible to see how the constituent materials of the composite work and it is evident in the shear and tensile cases that the matrix deforms more due to its much lower stiffness compared to the fibres. From the Huber von Mises reduced stress distributions, higher values are obtained in the fibres, which is also understandable due to their higher stiffness. In the case of tensile and shear in the XY and XZ planes, the highest stress values are obtained at the fibre-resin interface. By far the highest stress level is reached in the case of tension in the X-axis direction, where it is 1354 MPa. Once all the load cases have been simulated, the necessary information for calculating the material constants can be derived from the results. The mechanical properties of the composite represented by them were calculated using the following equations [3, 21]: 1Ei=FiAi2·ΔLiLi,{E_i} = {{{{{F_i}} \over {{A_i}}}} \over {{{2\cdot\Delta {L_i}} \over {{L_i}}}}}, 2vij=ΔLjLjΔLiLi,{v_{ij}} = {{{{\Delta {L_j}} \over {{L_j}}}} \over {{{\Delta {L_i}} \over {{L_i}}}}}, 3Gij=FijAjΔLiΔLj+ΔLjΔLi,{G_{ij}} = {{{{{F_{ij}}} \over {{A_j}}}} \over {{{\Delta {L_i}} \over {\Delta {L_j}}} + {{\Delta {L_j}} \over {\Delta {L_i}}}}},

Where E is Young modulus, v is Poisson ratio, G is Kirchoff modulus and i,j are directions of the local coordinate system (x, y, z). F_i is the force normal to the surface in i direction, F_i,j is the shear force in i direction applied to a surface normal to j direction. L_i,j is dimension of RVE with respect to the local coordinate system, △L_i,j is the displacement under loading conditions. Top and bottom surfaces are normal to Y axis, front and back surfaces are normal to X axis, right and left surfaces are normal to Z axis. From tension in X direction Young’s modulus E₁ and Poisson’s ratio v₁₂, v₁₃ are calculated. Similarly, Young’s modulus E₂ and v₂₃ Poisson’s ratio are calculated from tension in Y direction and E₃ Young’s modulus is obtained from tension in Z direction. Kirchoff modules are calculated from shear loads in XY, YZ and XZ planes. Examples of labelled dimensions and displacements for tension and shear are shown in Figure 11. Table 6 shows the results of the material constants obtained by the two methods and using three different fibre geometry.

Fig.11.

Model displacements to calculate Young's modulus, Poisson's ratio and shear modulus

A comparison of the results obtained reveals that the key significance for the values of the material constants in fibre composites does not depend solely on the ratio of the component materials used, but also on the distribution of the fibres in the resin and the geometry of the fibre. The final column of the table displays the percentage discrepancy between the values obtained using the periodic conditions in mechanical for a fibre with a circular cross-section and its actual scan. The most minute disparities are observed in Young's modulus and Poisson's ratios along the direction of material reinforcement. The maximum discrepancies are observed for the Kirchoff moduli, attaining up to 11% of the total variation. However, the most significant parameter in the case of fibre composites is the Young's modulus in the direction of reinforcement, as it is from this that the strength primarily depends.

A Type IV composite tank with a cylindrical shape was designed, with an internal plastic liner and aluminium bosses at the ends of the tank [13, 21, 22]. Figure 12 shows a cross section of the model with general dimensions.

Fig. 12.

Solid model of composite pressure vessel

The open boss is sealed with a plug to provide more realistic conditions inside the structure. The structure is reinforced with carbon fibre over the entire surface of the tank. The liner thickness is 8mm, the thickness of the composite overlap in the cylindrical section is 26,25mm and the maximum thickness in the dome section resulting from the filament winding method is 41mm. Table 7 shows the material properties of the aluminium and plastic used for the tank bosses and liner. The material constants were chosen based on existing materials. [23, 24]. Only the linear properties of these materials have been used as the simulation includes the working pressure of the tank. The material constants shown in Table 6 were applied sequentially to the composite material, and the effect of the different fibre geometries on the results was checked. Therefore, a total of four simulations containing different material properties will be performed.

Based on the data provided by the Toray composite manufacturer, the composite strength properties necessary to evaluate the stress state of the structure according to the intended criteria were also selected. The properties are given in Table 8. Selected properties are available for T700 carbon-epoxy composite with 60% fibre volume, which is constructed from the fibres used in the micromechanical model.

3.1.

Construction of composite laminate

The ACP Composite module in Ansys Workbench was used to design the overlap composite (25). The entire laminate consists of 42 hoop layers and 30 helical layers (26) with a carbon fibre tape thickness of 0,25mm. The construction of the laminate also used a variable winding polar radius of the helical layers to achieve a smooth decrease in thickness of the laminate away from the bottom of the tank. Table 9 shows all the spiral layer angles used, corresponding to the variable polar radius. Helical layers occur in both the dome and cylindrical sections of the tank; hoop layers occur only in the cylindrical section of the tank and have a constant winding angle of 90°.

Tab. 9.

Winding angleα [°]	Total number of layers	Radius of polar openings r0 [mm]
10	10	30
14,5	4	42
19	4	54
23	4	66
27,5	4	78
31,5	2	90
36,5	2	102

Wang's method was used to calculate the thickness distribution of the helical layers [12, 13, 14, 27]. The winding angels are calculated from equation: 4α=asin(r0R)\alpha = {\mathop{\rm asin}\nolimits} \left( {{{{r_0}} \over R}} \right)

Where R is the radius of the cylinder section, r₀ is, the radius of the polar opening. The thickness distribution can be determined from equations that reflect the actual thickness after filament winding with a high degree of accuracy. In the part of the dome section where the radius of the tank varies from r₀ to r_2b, the following equation is used: 5t(r)=m1+m2·r+m3·r2+m4·r3t(r) = {m_1} + {m_2}\cdotr + {m_3}\cdot{r^2} + {m_4}\cdot{r^3}

The parameters m₁, m₂, m₃, m₄ of the equations are calculated from: 6[ m1m2m3m4 ]=A−1·B\left[ {\matrix{ {{m_1}} \hfill \cr {{m_2}} \hfill \cr {{m_3}} \hfill \cr {{m_4}} \hfill \cr } } \right] = {A^{ - 1}}\cdotB

Where A and B are matrices as shown below: 7A=[ 1r0r02r031r2br2b2r2b3012·r2B3·r2B2π·(r2b2−r02)2·π3·(r2b3−r03)π2·(r2b4−r04)2·π5·(r2b5−r05) ]A = \left[ {\matrix{ 1 & {{r_0}} & {r_0^2} & {r_0^3} \cr 1 & {{r_{2b}}} & {r_{2b}^2} & {r_{2b}^3} \cr 0 & 1 & {2\cdot{r_{2B}}} & {3\cdotr_{2B}^2} \cr {\pi \cdot\left( {r_{2b}^2 - r_0^2} \right)} & {{{2\cdot\pi } \over 3}\cdot\left( {r_{2b}^3 - r_0^3} \right)} & {{\pi \over 2}\cdot\left( {r_{2b}^4 - r_0^4} \right)} & {{{2\cdot\pi } \over 5}\cdot\left( {r_{2b}^5 - r_0^5} \right)} \cr } } \right] 8B=[ tR·π·R·cosα0m0·bmR·nRπ·(acos(r0r2b)−acos(rbr2b))·t0mR·nRπ·(r0r2b·r2b2−r02−r0r2b·r2b2−rb2)·t0vconst ]B = \left[ {\matrix{ {{t_R}\cdot\pi \cdotR\cdot{{\cos {\alpha _0}} \over {{m_0}\cdotb}}} \cr {{{{m_R}\cdot{n_R}} \over \pi }\cdot\left( {{\mathop{\rm acos}\nolimits} \left( {{{{r_0}} \over {{r_{2b}}}}} \right) - {\mathop{\rm acos}\nolimits} \left( {{{{r_b}} \over {{r_{2b}}}}} \right)} \right)\cdot{t_0}} \cr {{{{m_R}\cdot{n_R}} \over \pi }\cdot\left( {{{{r_0}} \over {{r_{2b}}\cdot\sqrt {r_{2b}^2 - r_0^2} }} - {{{r_0}} \over {{r_{2b}}\cdot\sqrt {r_{2b}^2 - r_b^2} }}} \right)\cdot{t_0}} \cr {{v_{const }}} \cr } } \right]

Components of matrices are calculated from equations: 9vconst =∫r0rb2·π·r·mR·nRπ·acos(r0r)·t0dr+∫rbr2b2·π·r·mR·nRπ·(acos(r0r2b)−acos(rbr2b)){v_{const }} = \int_{{r_0}}^{{r_b}} 2 \cdot\pi \cdotr\cdot{{{m_R}\cdot{n_R}} \over \pi }\cdot{\mathop{\rm acos}\nolimits} \left( {{{{r_0}} \over r}} \right)\cdot{t_0}dr + \int_{{r_b}}^{{r_{2b}}} 2 \cdot\pi \cdotr\cdot{{{m_R}\cdot{n_R}} \over \pi }\cdot\left( {{\mathop{\rm acos}\nolimits} \left( {{{{r_0}} \over {{r_{2b}}}}} \right) - {\mathop{\rm acos}\nolimits} \left( {{{{r_b}} \over {{r_{2b}}}}} \right)} \right)

Radius of polar opening with fibre band width: 10rb=r0+b.{r_b} = {r_0} + b

Radius of polar opening with two fibre band width: 11rb=r0+2⋅b{r_b} = {r_0} + 2 \cdot b

Number of fibre bands in the cylindrical section for one layer: 12mR=2·π·R·cosαb{m_R} = {{2\cdot\pi \cdotR\cdot\cos \alpha } \over b}

Number of fibre bands which create one layer in dome section: 13m0=mR{m_0} = {m_R}

Thickness of helical layers in cylindrical section: 14tR=2·t0·nR{t_R} = 2\cdot{t_0}\cdot{n_R} n_R is the total number of helical layers and t₀ is the thickness of fibre band. In the dome section where the radius of the tank varies from r_2B to R, the thickness is calculated using the following equation: 15t(r)=mR·nRπ·(acos(r0r)−acos(rBr))·t0t(r) = {{{m_R}\cdot{n_R}} \over \pi }\cdot\left( {{\mathop{\rm acos}\nolimits} \left( {{{{r_0}} \over r}} \right) - {\mathop{\rm acos}\nolimits} \left( {{{{r_B}} \over r}} \right)} \right)\cdot{t_0}

Based on the calculations using the equations shown, thickness distributions were obtained and used to build the numerical model. The overall lay-up configuration is shown in the Table 10 and the thickness distribution of the composite at the dome part is shown in Figure 13 for the first 8 helical layers.

Fig.13.

Composite thickness distribution at dome section for first 8 helical layers

Tab. 10.

Number of layers	Angle [°]	Quantity
2	10	x2
3	90
1	14,5
1	19
3	90
1	23
1	27,5
3	90
1	31,5
1	36,5
3	90
2	10	x2
3	90
1	14,5
1	19
3	90
1	23
1	27,5
3	90
2	10	x1

A three-dimensional numerical model of the composite pressure vessel was created using Ansys Workbench, and a cross-section of the model with finite element mesh is shown in Figure 14. Once the material had been constructed and the material constants had been determined, a composite overlap was constructed using the ACP Composite Prepost module. The model used an algorithm to create layers with variable thickness and variable fibre angle according to the equations described by Wang's method. Initially, a surface model was constructed upon which a two-dimensional finite element mesh was superimposed. Thereafter, a layer layout was defined as illustrated in Figure 13, incorporating variable thickness and winding angle parameters. Subsequently, a three-dimensional model was generated based on the thickness data of the individual layers and imported into the mechanical environment. The threedimensional models of the liner and the bosses were permanently connected, and a three-dimensional finite element mesh was constructed.

Fig. 14.

Finite element mesh on 3D tank model

In the numerical model of the composite laminate, each layer has a finite element in thickness, and the generation of composite layers allows to capture the structure of the laminate. The thickness distribution of the composite in the finite element part of the dome compared to the data calculated by Wang’s method is shown in Figure 15.

Fig. 15.

Thickness distribution of composite at dome part in 3D model: a) Wang method, b) numerical model

The shape of the composite overlap was smoothed using the extrusion guide function in ACP Prepost. Figure 8 shows diagrams of aspect ratio, element quality, skewness and Jacobian ratio of all finite elements in 3D model of liner and bosses and 2D surface where composite 3D model was built. The maximum value of aspect ratio is 5,07 and most of the elements have a ratio below 3. Minimum value of element quality is 0,24 and more than 90% of all elements have a ratio above 0,4.

The maximum value of skewness is 0,67 and more than 90% of the finite elements have this ratio below 0,4. The last graph is the Jacobian ratio where the minimum value is 0,2 and more than 90% have this ratio above 0,4, so the summary mesh has very good quality in this model.

Fig. 16.

Quality of finite elements mesh graph: a) aspect ratio, b) element quality, c) skewness d) Jacobian ratio

4.1.

Boundary conditions

A fixed support is used in the end face of a closed aluminium boss. The internal pressure in the cylinder is 70MPa. Surfaces of an aluminium boss that are in contact with plastic liner surfaces are shared using the Share Topology function. Frictional contact was defined between the composite inner and outer surfaces of the plastic liner and the aluminium boss with a coefficient of friction μ = 0,2 [4, 5, 26]. On the outer cylindrical surface of the open boss, frictional contact is used to provide more stability for this model during simulation. Due to the use of non-linear frictional contact and high pressure, large deformations were included in the simulation options.

Fig. 17.

Boundary conditions

The Puck, Hashin and Tsai Wu criteria were utilised in the analysis of the composite strain. The initial two criteria enable the evaluation of strain due to fibre or matrix damage, i.e. treating the material as heterogeneous. Conversely, the Tsai Wu criterion does not permit such precise analysis and treats the composite as a homogeneous material. The equations that describe each criterion, incorporating the mechanical properties and components of the stress tensor, are outlined below. Failure modes in Puck criterion are as follows [9, 14, 28]:

Fiber failure 3D: 16σ1≥0;ff=σ1Xt{\sigma _1} \ge 0;{f_f} = {{{\sigma _1}} \over {{X_t}}} 17σ1≤0;ff=σ1Xc.{\sigma _1} \le 0;{f_f} = {{{\sigma _1}} \over {{X_c}}}

Matrix failure 3D: 18σn≥0;fm=((1R⊥+−p⊥θ+R⊥θA)σn)2+(τntR⊥⊥A)2+(τn1R⊥‖A)2+p⊥θ+R⊥θA·σnσn≤0;fm=(τntR⊥⊥A)2+(τn1R⊥‖A)2+(p⊥‖−R⊥θAσn)+p⊥θ−R⊥θAσn\matrix{ {\matrix{ {{\sigma _n} \ge 0;{f_m} = \sqrt {{{\left( {\left( {{1 \over {R_ \bot ^ + }} - {{p_{ \bot \theta }^ + } \over {R_{ \bot \theta }^A}}} \right){\sigma _n}} \right)}^2} + {{\left( {{{{\tau _{nt}}} \over {R_{ \bot \bot }^A}}} \right)}^2} + {{\left( {{{{\tau _{n1}}} \over {R_{ \bot }^A}}} \right)}^2}} + {{p_{ \bot \theta }^ + } \over {R_{ \bot \theta }^A}}\cdot{\sigma _n}} \hfill \cr {{\sigma _n} \le 0;{f_m} = \sqrt {{{\left( {{{{\tau _{nt}}} \over {R_{ \bot \bot }^A}}} \right)}^2} + {{\left( {{{{\tau _{n1}}} \over {R_{ \bot }^A}}} \right)}^2} + \left( {{{p_{ \bot }^ - } \over {R_{ \bot \theta }^A}}{\sigma _n}} \right)} + {{p_{ \bot \theta }^ - } \over {R_{ \bot \theta }^A}}{\sigma _n}} \hfill \cr } } \hfill \cr }

Where σ_n, τ_nt, and τ_n1 are stresses described by equations: 19σn=σ2cos2θ+σ3sin2θ+2τ23sinθcosθ{\sigma _n} = {\sigma _2}{\cos ^2}\theta + {\sigma _3}{\sin ^2}\theta + 2{\tau _{23}}\sin \theta \cos \theta 21τnt=(σ3−σ2)sinθcosθ+τ23(cos2θ−sin2θ){\tau _{nt}} = \left( {{\sigma _3} - {\sigma _2}} \right)\sin \theta \cos \theta + {\tau _{23}}\left( {{{\cos }^2}\theta - {{\sin }^2}\theta } \right) 22τn1=τ13sinθ+τ12cosθ{\tau _{n1}} = {\tau _{13}}\sin \theta + {\tau _{12}}\cos \theta

The parameter R⊥⊥AR_{ \bot \bot }^A is defined as: 23R⊥⊥A=R⊥−2(1+p⊥⊥−)R_{ \bot \bot }^A = {{R_ \bot ^ - } \over {2\left( {1 + p_{ \bot \bot }^ - } \right)}}

The parameters R⊥+,R⊥‖AR_ \bot ^ + ,R_{ \bot }^A are the same as the tensile strength transverse to the fibre direction and the shear strength. Constant parameters for the Puck failure criterion for carbon fibre: 24p⊥‖+=0,35p⊥‖−=0,3p⊥⊥+=0,25p⊥⊥−=0,2p_{ \bot }^ + = 0,35\quad p_{ \bot }^ - = 0,3\quad p_{ \bot \bot }^ + = 0,25\quad p_{ \bot \bot }^ - = 0,2

Hashin failure criteria in Ansys is defined by following equations [29, 30]:

Fiber tensile failure 3D: 25ff=(σ11Xt)2+1S122(τ122+τ132){f_f} = {\left( {{{{\sigma _{11}}} \over {{X_t}}}} \right)^2} + {1 \over {{S_{12}}^2}}\left( {{\tau _{12}}^2 + {\tau _{13}}^2} \right)

Matrix Tensile failure 3D: 26fm=1Yt2(σ22+σ33)2+1s232(τ232−σ22σ33)+1s122(τ122+τ132),σ22+σ33>0{f_m} = {1 \over {Y_t^2}}{\left( {{\sigma _{22}} + {\sigma _{33}}} \right)^2} + {1 \over {s_{23}^2}}\left( {\tau _{23}^2 - {\sigma _{22}}{\sigma _{33}}} \right) + {1 \over {s_{12}^2}}\left( {\tau _{12}^2 + \tau _{13}^2} \right),\quad {\sigma _{22}} + {\sigma _{33}} > 0

Matrix compression failure 3D: 27fm=1Yc[ (Yc2S23)2−1 ](σ22+σ33)+14S232(σ22+σ33)2+1S232( τ232− σ22σ33 )+1s122(τ122+τ132),σ22+σ33<0{f_m} = {1 \over {{Y_c}}}\left[ {{{\left( {{{{Y_c}} \over {2{S_{23}}}}} \right)}^2} - 1} \right]\left( {{\sigma _{22}} + {\sigma _{33}}} \right) + {1 \over {4{S_{23}}^2}}{\left( {{\sigma _{22}} + {\sigma _{33}}} \right)^2} + {1 \over {{S_{23}}^2}}\left( {{\tau _{23}}^2 - } \right.\left. {{\sigma _{22}}{\sigma _{33}}} \right) + {1 \over {{s_{12}}^2}}\left( {{\tau _{12}}^2 + {\tau _{13}}^2} \right),{\sigma _{22}} + {\sigma _{33}} < 0

The Tsai Wu failure criterion in Ansys is defined by the following equation (30): 28f=σ12XtXC+σ22YtYc+σ32ZtZc+τ122 S122+τ132 S132+τ232 S232−σ1σ2XtXcYtYc−σ2σ3YtYcZtZc−σ1σ3XtXcZtZc+σ1(1Xt−1Xc)+σ2(1Yt−1Yc)+σ3(1Zt−1Zc)f = {{\sigma _1^2} \over {{X_t}{X_C}}} + {{\sigma _2^2} \over {{Y_t}{Y_c}}} + {{\sigma _3^2} \over {{Z_t}{Z_c}}} + {{\tau _{12}^2} \over {{\rm{S}}_{12}^2}} + {{\tau _{13}^2} \over {{\rm{S}}_{13}^2}} + {{\tau _{23}^2} \over {{\rm{S}}_{23}^2}} - {{{\sigma _1}{\sigma _2}} \over {\sqrt {{X_t}{X_c}{Y_t}{Y_c}} }} - {{{\sigma _2}{\sigma _3}} \over {\sqrt {{Y_t}{Y_c}{Z_t}{Z_c}} }} - {{{\sigma _1}{\sigma _3}} \over {\sqrt {{X_t}{X_c}{Z_t}{Z_c}} }} + {\sigma _1}\left( {{1 \over {{X_t}}} - {1 \over {{X_c}}}} \right) + {\sigma _2}\left( {{1 \over {{Y_t}}} - {1 \over {{Y_c}}}} \right) + {\sigma _3}\left( {{1 \over {{Z_t}}} - {1 \over {{Z_c}}}} \right)

After the numerical simulation, the following results are presented: the distribution of displacements and deformations of the structure, the stresses in the liner and the bosses, and the distributions of the failure indices of the composite structure using the previously mentioned criteria. In addition, the maximum values of these failure indices for each layer are presented and a comparison is made between the results obtained using these criteria. All distributions are presented for the model, considering the material constants determined for the round fibre in the Material Designer. This is because the differences in results on a macroscopic scale, when using the other composite constants, are negligible.

6.1.

Displacements and strain distributions

The maximum value of the resultant displacement is 3,522mm and it is located in the area of the hole near the aluminium boss, and values of about 3.5mm also occur in the aluminium boss and liner in this section. The displacement value increases along the axis of the vessel away from the fixed end of the closed aluminium boss. The maximum strain occurs in the liner where the shape changes from cylindrical to dome section. The lowest values are found in the composite and ferrule, as the materials used there have a much higher stiffness than the liner material.

Fig. 18.

Distribution of resultant: a) displacement b) strain

6.2.

Stress distribution in aluminium bosses and plastic liner

In the closed aluminium boss, the highest value of von-Mises stress is 255,51MPa and is located at the rounding where the composite laminate ends at a radius of 30mm, which is understandable since there is a sharp change in the stiffness of the model. In the open aluminium boss, the stress distribution is very similar, with the maximum value for this element being 259,32 MPa at the inner hole where the liner has bonded contact with the boss. In both stubs, the yield strength, which is 276 MPa, was not exceeded. The next figure shows the distribution of HMH reduced stress in the plastic liner. The maximum stress value in the liner is attained at the point at which the liner's shape transitions from elliptical to cylindrical. This value falls below the yield strength of the liner material. From the distribution, the values in the cylindrical section remain constant, reaching approximately 40-50 MPa.

Fig.19.

Distribution of Huber von-Mises stress in a) closed boss b) open boss [MPa]

Fig. 20.

Distribution of Huber von-Mises stress in the polyamide liner [MPa]

6.3.

Puck failure criterion

The evaluation of composite laminate was conducted using the puck failure criterion. Figure 21 provides a visual representation of the distribution of fibre and matrix failure. Table 11 presents the maximum values of the puck failure criterion for each layer winding angle. The maximum failure fibre index value of 0,385 is observed in the cylindrical section of the first hoop layer, with values approaching this maximum found in the initial few hoop layers of the cylindrical section. Conversely, the minimum value of 0,02 is attained in the outer layer of the composite layer at an angle of 41° in the dome section. Larger values exceeding 0,3 are attained in the vicinity of the composite's interaction with the aluminium boss and within the composite's inner region, adjacent to the transition from a cylindrical to an elliptical shape. Given a safety factor of 2,5, it can be deduced that the fibre satisfies the stipulated strength criteria, as the permissible value of this factor is 0,4. Additionally, values approaching 0,3 are achieved at the point of contact with the polar hole. The maximum value of the matrix failure index is 0,836, which is located near the polar hole on the external surface of the composite overlap, while a value closer to the maximum is observed at the inner surface of the composite overlap, where contact is made between the aluminium boss and the plastic liner. High values of failure also occur in the section where the shape changes from cylindrical to elliptical. In the cylindrical part, the outer helical layers exhibit higher values of destruction indices ranging from 0,6 to 0,7, while the inner layers demonstrate values of approximately 0,4 to 0,5. For the hoop layers, the maximum value of the matrix failure index is recorded as 0,544, which is observed in the external hoop layer. In contrast, the internal hoop layers exhibit lower values. It is evident that the failure of the matrix is prevented due to the failure index value not exceeding the threshold value of 1. The results presented in Table 11 demonstrate a decrease in the maximum destruction indices of the helical layers for the fibre failure criterion as the winding angle increases. A comparison of the distributions of Puck's fibre and matrix failure indices reveals that for each winding angle, the values of the failure indices of the matrix are higher compared to fibre indices, indicating that the matrix is the primary source of failure in this construction.

Fig. 21.

Puck failure distribution: a) fibre, b) matrix

Tab. 11.

Angle [°]	Fiber failure	Matrix failure
10	0,329	0,836
14,5	0,266	0,695
19	0,238	0,715
23	0,210	0,636
27,5	0,185	0,620
31,5	0,182	0,582
36,5	0,182	0,603
90	0,385	0,544

6.4.

Hashin failure criterion

As with the Puck criterion, the maximum value of the hoop layer failure index is 0,384, which is located in the first hoop layer. Hoop layers that are outside the overlap take smaller values than external hoop layers. A significant difference in this case is that the maximum value of the fibre failure index occurs in the first inner helical layer near the polar hole at 0,774. This value is more than double that of the Puck criterion, because the Hashin criterion also considers the shear stresses for the fibre. For all helical layers, the values of the fibre failure take on values in the prevalent range from 0,38 to 0,51, where for the Puck criterion these values are in the range from 0,18 to 0,33. High values of fibre failure, reaching approximately 0,6, are evident in the area where there is a change of shape from the cylindrical to the elliptical part, which, in the Puck criterion, is evident in the failure of the matrix. The maximum matrix failure value of 0,795 is observed near the polar hole on the outer helical layer of overlap. It is notable that high values, akin to the maximum observed, are also evident at the interface between the liner and aluminium bosses, a phenomenon analogous to the Puck criterion. A comparative analysis of the maximum matrix failure values ascertained in each layer, as illustrated in Tables 11 and 12, reveals a subtle disparity of approximately 0,05. It is observed that the matrix failure values in the cylindrical segment exhibit an upward trend with each subsequent layer, attaining a maximum of approximately 0,6. Concurrently, values of approximately 0,6 are observed in all helical layers within the region undergoing a shape transition from cylindrical to elliptical. A comparative analysis of the distributions and values of matrix failure for both the Puck and Hashin criteria reveals a high degree of similarity, with Hashin failure values exhibiting a slight increase, reaching approximately 0,05, in comparison to Puck failure. However, a more pronounced distinction emerges in the context of fibre failure, where values under the Hashin criterion exceed twice those of the Puck criterion, a consequence of the incorporation of shear stresses.

Fig. 22.

Hashin failure distribution: a) fibre, b) matrix

Tab. 12.

Angle [°]	Fiber failure	Matrix failure
10	0,774	0,795
14,5	0,442	0,629
19	0,448	0,760
23	0,496	0,665
27,5	0,509	0,603
31,5	0,472	0,596
36,5	0,446	0,622
90	0,384	0,598

6.5.

Tsai-Wu failure criterion

Tsai Wu's failure criterion treats composite material as homogeneous, with no distinction made between fibre and matrix. Under this criterion, it is not possible to ascertain which constituent materials are more strained. The maximum value of the composite failure was obtained in the same location as in the previous cases of the Puck and Hashin criterion, around the polar hole in the outer helical layer. Regarding the distribution of failure in the cylindrical component, the values are found to be similar at approximately 0,5-0,6 across the entire thickness of the composite. It is noteworthy that the Tsai Wu criterion values for each layer, as presented in Table 13, are lower than those of the Hashin and Puck criterion for the matrix.

Fig. 23.

Tsai Wu criterion failure distribution

Tab. 13.

Angle [°]	Failure
10	0,759
14,5	0,547
19	0,678
23	0,498
27,5	0,536
31,5	0,566
36,5	0,581
>90	0,580

6.6.

Dehomogenisation and results for RVE

In the multiscale model, a dehomogenisation technique has been employed to transfer the deformation state from the macro to the micro scale. The finite elements were read out for the numerical model of the tank, in which the maximum values of the failure indices were reached according to the strength criteria used. Strain tensors from all simulations containing different fibre geometries for these elements were subsequently measured and transferred to the RVE, utilising the equations shown in the table 15. The maximum strain values for RVEs with circular fibre are shown in Table 14. The equations shown in Table 14 describe the deformations on all RVE faces. However, it should be noted that each RVE edge lies simultaneously on two planes, and each vertex lies simultaneously on three planes. Therefore, there are also separate equations for the edges and vertices themselves, as in the case of periodic boundary conditions. The equation describing the displacement of an edge results from the common part of the two equations describing the displacement of the surfaces on which the edge is located, and the equation describing the displacement of a vertex results from the common part of the three equations describing the displacement of the surfaces on which the vertex is located. The strain tensors were transferred from the most stressed finite elements according to the given failure criteria to the RVE, maintaining the consistency of the fibre coordinate systems. Numerical simulation of the loaded RVEs was performed for all three fibre geometries used, resulting in stress distributions in the constituent materials for different cell types. First, the stress distributions for the circular fibre cell are presented and compared with the failure index values obtained in the macroscale model.

Tab. 14.

Failure criteria/strain	εxx	εyy	εzz	εxy	εxz	εyz
Puck fibre	7,986 · 10−3	–1,011 · 10−2	2,073 · 10−3	–9,570 · 10−5	2,034 · 10−4	6,786 · 10−4
Puck matrix	4,521 · 10–3	2,370 · 10−3	–3,099 · 10−3	–8,194 · 10−3	–8,872 · 10−3	3,176 · 10−3
Hashin fibre	2,230 · 10–3	–6,224 · 10−3	–5,999 · 10−3	1,383 · 10−2	1,201 · 10−3	8,919 · 10−3
Hashin matrix	4,522 · 10−3	2,370 · 10−3	–3,098 · 10−3	–8,194 · 10−3	–8,875 · 10−3	3,171 · 10−3
Tsai Wu	4,521 · 10−3	2,370 · 10−3	–3,099 · 10−3	–8,194 · 10−3	–8,872 · 10−3	3,176 · 10−3

Tab. 15.

Surface X = 0	Surface Y = 0	Surface Z = 0
u(0,y,z) = εxy 1 ε.xzzv(0,y,z) = εyyy + εyzzw(0,y,z) = εyzy + εzzz	u(x, 0, z) = εxxx + εxzzv(x, 0, z) = εxyx + εyzzw(x, 0, z) = εxzx + εzzz	u(x,y,0) = εxxx + εxyyv(x,y,0) = εxyx + εyyyw(x,y,0) = εxzx + εyzy
Surface X = Lx	Surface Y = Lx	Surface Z = Lx
u(Lx,y,z) = εxxLx + εxyy + εxzzv(Lx,y,z) = εxyεx + εyyy + εyzZw(Lx,y,z) = εxzLx + εyzy + εzzz	u(x, Ly, z) = εxxx + εxyLy + εxzzv(x, Ly, z) = εxyx + εyyLy + εyzzw(x, Ly, z) = εxzx + εyzLy + εzzz	u(x,y,Lz) = εxxx + εxyy + εxzLzv(x,y,Lz) = εxyx + εyyy + εyzLzw(x,y,Lz) = εxzx + εyzy + εzzLz

Fig. 24.

Stress distributions for maximum Puck fibre failure criterion strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 25.

Stress distributions for maximum Puck matrix failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 26.

Stress distributions for maximum Hashin fibre failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 27.

Stress distributions for maximum Hashin matrix failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 28.

Stress distributions for maximum Tsai Wu failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

When comparing the equivalent stress values of all the criteria used, the highest values of over 1800 MPa are obtained using the Puck fibre failure criterion. Almost identical values of equivalent stresses around 1100 MPa are obtained for the Puck matrix, Hashin matrix and Tsai Wu failure criteria. It is worth noting that the maximum value of the failure index for the Puck matrix criterion and the Tsai Wu criterion occurs in the same finite element of the composite overlap of the tank, and the maximum value of the Hashin matrix failure criterion occurs close to this finite element. The minimum value of the equivalent stress is found for the Hashin fibre failure criterion and is 676 MPa. There is a characteristic contrast between the values of this stress in the fibres themselves, which is not the case for the other criteria. The tensile strength of the fibres used in the composite is 4900 MPa [9] and the maximum stress in the X direction in the fibre for the Puck criterion is 1833 MPa. Using Equation 19, it can be calculated that the failure index value considering the fibres alone is 0.374, which is very close to the result obtained for the composite treated as a homogeneous material. The Hashin fibre failure criterion also considers the shear stress in the fibre planes, so that the value of the stress in the fibre direction in the most stressed element is only 553 MPa, which is much lower than the Puck fibre criterion because the stress level here is increased by shear stresses in the fibre plane. For the Puck and Hashin matrix criteria and for Tsai Wu, the stress distributions in the grain direction are almost the same, as are the equivalent stresses. The matrix equivalent stress distributions provide interesting results. Based on the literature, the average tensile strength values for the resin used are around 70-90 MPa. For the Puck fibre failure criterion, the maximum stress value is 90 MPa. A much higher level of stress is achieved for the Hashin fibre criterion, up to 154 MPa, which is well above the tensile strength. Despite these high stress values, these characteristic points were not classified as the points with the highest matrix failure indices, so they will require more detailed analysis. On the other hand, for the most stressed components according to the Hashin and Puck matrix criterion and the Tsai Wu criterion, the maximum value of the reduced stress in the resin is 57 MPa, which is still below the allowable value. The fibres are undoubtedly under greater shear stress here, hence the higher failure values. The following images show the stress distributions for RVE with elliptical fibres and with fibres of the actual geometry.

Fig. 29.

Stress distributions for maximum Puck fibre failure criterion strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 30.

Stress distributions for maximum Puck matrix failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 31.

Stress distributions for maximum Hashin fibre failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 32.

Stress distributions for maximum Hashin matrix failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 33.

Stress distributions for maximum Tsai Wu failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Based on the obtained stress distributions for the representative volume element (RVE), it can be observed that the maximum stress values in the constituent materials are very similar to the maximum values obtained for the cell with circular fibres. The largest difference in normal stress within the fiber is observed for the Puck fiber failure criterion, where the values are 1832,9 MPa and 1760 MPa, respectively, corresponding to a difference of approximately 4%. The greatest discrepancy in maximum stress within the resin occurs in the case of the Hashin fiber failure criterion, with values of 154,4 MPa and 177,2 MPa, respectively. Here, the difference is more substantial, reaching nearly 15%. In both considered fiber geometries, the stresses in the resin significantly exceed the allowable stress limits for this material. In the case of the obtained stress distributions for the representative volume element with the actual fiber geometry, significant differences are observed in nearly every instance, which can be attributed to the irregular shape of the fiber. The largest discrepancies occur in the maximum stress values within the resin. The most substantial variation in stress state is found for the Hashin matrix failure criterion, where the stresses are 57 MPa and 138,7 MPa, respectively representing a remarkably large deviation of approximately 143%. For the fiber, the greatest stress difference also occurs under the Hashin matrix criterion, with values of 1128,2 MPa and 1204,5 MPa, corresponding to a difference of about 7%. A surprising result is that the maximum stress values in the resin for the RVE with the actual fiber geometry significantly exceed the allowable stress limit in every case. Additionally, the stress distributions reveal that these elevated values occur exclusively at the points located at the resin–fiber interface. A summary of all maximum values from the stress distributions is presented in Table 16. Red colour indicates resin stress values that exceed the tensile strength of 80MPa.

Fig. 34.

Stress distributions for maximum Puck fibre failure criterion strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 35.

Stress distributions for maximum Puck matrix failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 36.

Stress distributions for maximum Hashin fibre failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 37.

Stress distributions for maximum Hashin matrix failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Fig. 38.

Stress distributions for maximum Tsai Wu failure criteria strains: a) Equivalent Huber von-Mises stress for RVE, b) normal X-direction stress for fibre, c) Equivalent Huber von-Mises stress for matrix

Tab. 16.

Failure criteria	Stress limit [MPa]		Circular fibre			Eliptical fibre			Actual fibre
	fibre	resin	RVE	fibre	resin	RVE	fibre	resin	RVE	fibre	resin
Puck fibre	4900	80	1868	1832,9	89,9	1795,4	1760	90,5	1923	1820	133,4
Puck matrix			1104,8	1128,1	57	1027,4	1067	55,1	1095,6	1133,6	76,3
Hashin fibre			676,1	553,2	154,4	680,9	563,35	177,2	1110,7	580,2	238,3
Hashin matrix			1104,9	1128,2	57	1107,9	1120,3	53,9	1188,9	1204,5	138,7
Tsai wu			1104,8	1128,1	57	1107,8	1120,2	54	1188,8	1204,4	138,7

The use of composite strength evaluation criteria at the microscopic scale provides the opportunity to evaluate the behaviour of the individual materials that make up the composite, as opposed to macroscopic evaluation criteria. Such a criterion, together with the use of micro-modelling technology, makes it possible to modify the structure of the composite material, to change the ratio of materials used and the way they are arranged in the structure, which later makes it possible to assess the influence of these parameters on the results of the obtained failure indices and stresses.

One of the critical points of the structure is where the composite overlaps with the liner and the boss, where there is a sharp change in stiffness, causing an increase in stress in the composite, as can be seen from the failure distributions of the composite for all the criteria used. The values do not exceed the maximum, but this is one of the most stressed areas. For the most part, the maximum values of the failure indices are reached in the outer cylindrical part in the helical layers, where composite cracking would occur at higher pressures. There are also high values in the part of the dome close to the hole (such as in the outer helical layer), but these are small areas around which the values are much smaller. It is possible that increasing the accuracy of the layer distribution in the numerical model, or arranging them differently, would compensate for this phenomenon.

The development of equations written using the mechanical APDL language that represent the periodic boundary conditions enabled the determination of the material constants of the composite. These conditions will be applicable to future studies in which fibre composites will be analysed on a microscopic scale. The equations may prove very useful for analysing temperature-dependent deformation, simulating epoxy resin polymerisation and investigating residual stresses generated during the manufacturing process.

By loading the RVE with a strain tensor, we have the opportunity to check the stress state of the constituent materials of the composite, which, as can be seen from the results obtained, should be left for further analysis, because in the one case where the material is within the safe range according to the Hashin criterion, the maximum equivalent stress in the resin exceeds the allowable value by a considerable margin. Using this approach, we can read off the normal and shear stresses in each direction for the fibres. By analysing the stresses following dehomogenisation, it is possible to enhance the strength criteria of the composites and to compare the results obtained with experimental studies.

Based on the obtained stress distributions in the constituent materials of the composite after dehomogenization, it is evident how significant the impact of fiber geometry is on the results. Discrepancies are also observed between these stress results and the outcomes of the strength criteria, prompting reflection on whether these criteria can be refined to better represent the stress state in fiber-reinforced composites. The exceeding of the allowable stress value in the resin may indicate the formation of microcracks in this material and delamination at the resin-fiber interface, as the highest values occur precisely in these regions of the RVE.

A comprehensive approach based on multiscale modelling with a simultaneous return to microstructure has been developed, which allows a more accurate analysis of the stress state at critical points in the structure. This approach can be employed at a subsequent stage to account for the non-linear material characteristics of the resin. This will facilitate a more profound comprehension and analysis of rheological phenomena in composite structures. The cross-sectional analysis of circular, elliptical and real fibres demonstrate the necessity of investigating responsible structures using a multiscale approach that considers dehomogenisation. This is of particular importance for the correct determination of crack initiation.