Introduction

Composite overwrapped pressure vessels (COPV) are used in the rocket and spacecraft industry due to their high strength and lightweight. Consisting of a thin, non-structural liner wrapped with a structural fiber composite COPV are produced to hold the inner pressure of tens and hundreds of atmospheres. COPV have been one of the most actual and perspective directions of research, supported especially by NASA [1,2].

Designing of a highly reliable and efficient COPV requires a technology for analisys of its deformation behaviour and strength assessment. This technology should allow one to obtain target COPV parameters through changing vessel's geometry, structural and mechanical material parameters while keeping its useful load.

Combination of mathematical modeling and numerical optimization makes it possible to reduce the cost and the duration of identifying the best parameters for a COPV. However, this approach is characterized by a number of hurdles. Overcoming these hurdles determine success of an optimum designing of such structures.

So far, there have been two main approaches in optimization of composite structures: analytical and numerical ones.

In the first approach the problems are solved basing on their simplified statement, for example using the momentless (membrane) shell theory and the netting model of composite material (CM) [3][4][5][6]. The obtained results may be far from reality, however they are of value for testing of numerical optimization methods.

Application of the numerical approach in designing, on the other hand, produces a number of challenges that must be overcome, e.g. lack of reliable methods for global optimization; nonconvexity and nonlinearity of constraint functions; ill-conditioned boundary value problems; different scaling of optimization criteria represent just some of the obstacles that prevent reliable optimization of COPV.

Numerical analysis is usually a computation-intensive process and takes considerable time. One way to solve this problem is approximation of the objective function using different approaches, such as response surface method [7] and neural network [8]. Some kinds of numerical analyses use a small number of design variables, functions and/or corresponding set of their discrete values (analytical geometry parametrization [9], finite set of feasible winding angles [10]). It leads to reduction in the number of objective function calculations.

Another way is reasonable simplification of the elasticity problem statement, for example by using the membrane theory or other shell theories [9,11,12], that leaves the question of results validity. And this is the approach we have applied in our study. For validation we have used the Timoshenko [13] and Andreev-Nemirovskii [14] shell theories, acconting transverse shears with different degrees of accuracy.

Of course, it should be taken into account that the computed solutions are not optimum in the strict mathematical sense. However, these solutions could provide the considerable economy of the weight while keeping the required strength, and, therefore, they are important from the engineering point of view.

The Problem Statement and the Mathematical Models

Let's consider a multilayer composite pressure vessel at a state of equilibrium under equidistributed inner pressure. We need to determine the parameters of structure and CM meeting the following requirements:

𝑉 ≥ 𝑉 0 , 𝑃 ≥ 𝑃 0 , 𝑀 ≤ 𝑀 0 ,(1)

where 𝑉 is the volume of the vessel, 𝑃 is inner pressure and 𝑀 is the vessel's mass and they are constrained by some preset values 𝑉 0 , 𝑃 0 , 𝑀 0 .

We define the optimization problems the following way: to find extremum of one functional from (1) under other constraints.

The structures optimization problem statement includes selection of objective functional, formulation of constitutive equations and constraints on performance and design variables.

The mathematical models describing the vessel's state are based on the following assumptions:

1) the vessel is a multilayer thin-walled structure; 2) the vessel's layers can have different mechanical characteristics; 3) the reinforced layer's material is quasi-homogeneous; 4) the vessel's main loading is high inner pressure, whose alteration happens rather slowly during operation.

These assumptions allow us to reduce dimension of the corresponding mathematical problem and to build the mathematical vessel's models based on the different theories of multilayer non-isotropic shells.

Let's consider the vessel as a shell rigidly compressed on the edges. Taking into account a symmetry plane in the middle of the vessel, it is enough to calculate and design only it's half. The type of loading and boundary conditions allow considering the axisymmetric problem statement.

The half of shell is set by rotation of the generatrix 𝑟 = 𝑟(𝜃) around axis 0𝑧 (fig. 1) where 𝑟 is the current point of the shell radius, 𝜃 is the angle between the normal to the shell surface and axis 0𝑧 changing within [𝜃 0 ; 𝜃 1 = 90 ∘ ]. The full shell is set by reflection the shell's half about plane 0𝑥𝑦. The Kirchhoff-Love shell theory [15] (KLST) and the theories with shear terms (Timoshenko [13] (TiST) and Andreev-Nemirovskii [14] (ANST)) are used to solve the direct calculation problems of multilayer composite vessels, to analyse their behavior and to verify optimization problem solutions. The used coordinate system is (𝜃, 𝜙, 𝜁), where 𝜙 denotes polar angle, 𝜁 -normal to the surface. The load is equidistributed inner pressure 𝑞 = (0, 0, 𝑞 3 ). On the top fixed edge (𝜃 = 𝜃 0 ) all displacement components, angle of normal rotation and additional shear term (ANST) are equal to zero; on the edge (𝜃 = 90 ∘ ) we use the symmetry conditions: transverse force, the first displacement component, angle of normal rotation and additional shear term (ANST) are equal to zero. One could find the full systems of equations in the paper [16].

Relations between stresses and strains are defined by the structural models [17]. The main idea of these models is that CM characteristics are calculated through matrix and fibers mechanical characteristics, fibers volume content and winding angles. The stress-strain state of matrix and fibers are evaluated through stresses and strains of the composite shell. A failure criterion is applied for every component of CM. Here we use the Mises criterion to determine the first stage of failure.

The objective function whose minimum is required is the minimum mass:

𝑀 = 2𝜋 ∫︁ 𝜃1 𝜃0 𝑟𝑅 1 ℎ𝑑𝜃 [𝜌 𝑚 (1 − 𝜔 𝑟 ) + 𝜌 𝑟 𝜔 𝑟 ] → min(2)

where 𝜌 𝑚 , 𝜌 𝑟 are the densities of matrix and reinforcing fibers, 𝜔 𝑟 is the volume content of reinforcement. We chose the following design functions: the curvature radius 𝑅 1 (𝜃) to define the generatrix; the thickness of the shell ℎ(𝜃); the winding angle 𝜓(𝜃) (fig. 1).

The solution has to satisfy the constraints on the shell's inner volume:

𝜋 ∫︁ 𝜃1 𝜃0 𝑟 2 𝑅 1 sin 𝜃𝑑𝜃 = 𝑉 0 ,(3)

and the strength requirement:

max{𝑏𝑠 𝑟 , 𝑏𝑠 𝑚 } ≤ 1,(4)

where 𝑏𝑠 𝑟 , 𝑏𝑠 𝑚 are the normalized von Mises stresses in the matrix and fibers [17].

Note that the factor of safety is widely used while solving engineering problems. It can be considered by correction of the right part of the inequality (4). We used the following constraints on the design functions:

0 ≤ 𝜓 ≤ 90, ℎ * 0 ≤ ℎ ≤ ℎ * 1 , 𝑅 * 0 ≤ 𝑅 1 ≤ 𝑅 * 1 .(5)

The method of the continuous geodesic winding have been widely used in the manufacturing of composite shells of revolutions. In this case the winding angles are defined by the Clairaut's formula:

𝑟 sin 𝜓(𝑟) = 𝐶,(6)

Mathematical and Information Technologies, MIT-2016 -Mathematical modeling where the constant 𝐶 is defined, as a rule, from the condition at the shell's equator. The thickness equation is

ℎ(𝑟) = ℎ 𝑅 𝑅 cos 𝜓 𝑅 𝑟 cos 𝜓(𝑟) ,(7)

which has the singularity at the edge where the winding angle has to be equal to 90 ∘ . The formula ( 7) is applied into practice at 𝑟 ≥ 𝑟 0 + 𝑟 𝜔 , where 𝑟 𝜔 is equal to the width of the reinforcement tape. As a result the equation defining the vessel's thickness takes the form:

ℎ(𝑟) = ⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ ℎ 𝑅 𝑅 cos 𝜓 𝑅 𝑟 𝜔 cos 𝜓(𝑟 0 + 𝑟 𝜔 ) , 𝑟 ≤ 𝑟 0 + 𝑟 𝜔 ; ℎ 𝑅 𝑅 cos 𝜓 𝑅 𝑟 cos 𝜓(𝑟) , 𝑟 ≥ 𝑟 0 + 𝑟 𝜔 .(8)

We did not consider the problem of fibers slippage. The main goal of the study was to demonstrate the potentials of using CM in one COPV design approach.

Direct Problems. Analysis of the Shell Theories

Estimation of composite vessel stress-strain state using the offered models leads to the solution of boundary value problems for stiff systems of differential equations. These problems are ill-conditioned, and their solutions have big gradients near the edges. Numerical analysis was performed by the spline collocation and discrete orthogonalization methods, implemented in the COLSYS [18] and GMDO [19] software. These computing tools have proved to be effective in numerical solving of wide range of composite shell mechanics problems [20].

We investigated the vessel's deformations by computing of its stress-strain state based on the different shell theories. The vessel's shape was a part of a toroid: 𝑅 1 = 2.46 m, 𝜃 0 = 0.108 ∘ , 𝜃 1 = 90 ∘ (the computed half), 𝑟(𝜃 0 ) =0.04 m. The CM parameters were: 𝐸 𝑚 = 3 • 10 9 Pa, 𝜈 𝑚 = 0.34, 𝐸 𝑟 = 300 • 10 9 Pa, 𝜈 𝑟 = 0.3, 𝜔 𝑟 = 0.55, 𝑉 0 =350 liters where 𝐸 𝑚 , 𝐸 𝑟 are the Young's modulus of the matrix and fibers, 𝜈 𝑚 , 𝜈 𝑟 -their Poisson's ratio.

Fig. 2 shows the stress-strain state characteristics of the vessel with the thickness ℎ = 0.6 cm, reinforced in the circumferential direction (𝜓 = 90 ∘ ) under the load of 170 atm. On the left are the displacements of the reference surface along the generatrix 𝑢 1 (𝑟) (dashed curves) and the normal displacement of this surface 𝑤(𝑟) (solid curves). On the right are the distribution of normalized von Mises stresses (nVMS) along the thickness in the matrix 𝑏𝑠 𝑚 (𝑟). The solid curves correspond to a slice at the shell edge, the dashed curves -to a slice at 𝜃 = 0.1. The curves without symbols correspond to KLST simulations, the curves marked with △ -to those using TiST, and -to ANST.

It's easy to see that the basic kinematic characteristics coincide both qualitatively and quantitatively. Small differences are observed only for the stresses w, mm q, rad u, mm bs m z, mm Fig. 2. The stress-strain state characteristics of the composite vessel computed using different shell theories and deformations near the compressed edge. The maximum results and qualitative difference were obtained for ANST. This is due to accounting for the transverse shears by non-linear distribution in a thickness of a shell. Earlier it was shown [20] that ANST's based results were closest to the ones of 3D elastic theory in most cases.

The winding angle's influence on the COPV performance was investigated using parametric analysis.

Dependence of the maximum nVMS in the matrix 𝑏𝑠 𝑚 (dashed curves) and the fibers 𝑏𝑠 𝑟 (dash-dotted curves), and the maximum size of the displacement vector ||𝑣|| (solid curves) are shown in fig. 3. KLST's results are drawn without marks, TiST -with symbols △, ANST -with . It was revealed that the winding angles corresponding to minimum stresses values were almost insensitive to the thickness variation. The change of ℎ from 0.6 to 1.6 cm corresponded to the angle's change about 0.2 ∘ .

Additionally we investigated stress-strain state of the vessel (the thickness ℎ = 0.6 cm, the winding angles 𝜓 = ±43.2), when nVMS in the matrix and the fibers were near their minimum (fig. 4). The adopted notation is the same as in fig. 2. And again the difference is visible only in a very small region near the edge but now this difference is small enough to be neglected. Moreover the displacement values of the reference surface, the efforts and the moments completely coincide for all the theories.

All the theories (KLST, TiST, ANST) provided similar estimated characteristics of stress-strain state. This vessel was characterized not only by essential decrease of the maximal nVMS in the matrix and fibers, but also by their almost uniform distribution along the generatrix. At the same time the values of bending moments significantly reduced bringing vessel's stress-strain state close to momentless.

The performed analysis showed that the optimizing problem can be solved using rather simple shell theories (KLST, TiST). These theories are characterized by lower computational complexity of corresponding boundary value problem if compared to ANST. It takes from 10 to 20 times less resources.

One can see that the winding angle as a design parameter gives an opportunity to increase the vessel's strength significantly. The difference between the "best" and "worst" designs can reach 20 -35 times comparing their nVMS in the matrix and fibers. The "worst" designs have the winding angle close to 90 ∘ . In this case are considerable transverse shears near the compressed edge, and the loading is redistributed to a rather weak matrix while the fibers remain unloaded.

Inverse problems. Optimization of the Vessel

Inverse problems involve not only numerical methods for fast and reliable solving of direct boundary value problems, but also require numerical optimization methods for identifying design parameters.

Here we considered conditional optimization problem, including direct constraints on design functions and trajectory constraints on the solution imposed at the end of the interval. The sequential unconstrained optimization is one of the most widespread approaches to solution of such problems. The main idea of the method is terminal functional convolution and multiple solutions of one-criteria problem using different optimization methods [21]. In our study the modified Lagrange function was used for the convolution.

Hence we sought for solution of a nonconvex problem of finite-dimensional optimization [22] by discretization of design functions. The methods implemented in the OPTCON-A software [23] were used to get the corresponding solution.

In our study several vessels with different type of design parameters were investigated. The parameters were either functions or constants. Additionally a design of continuous winding on the geodesic path was considered, where the only design function was its curvature radius.

Uniform mesh for design functions discretization included 7 points, except for the geodesic winding design with nonuniform mesh of 17 points. The distance between points was also the solution of corresponding optimization problem. Approximation of the design functions was carried out using the 3rd degree natural splines.

The masses of these vessels are shown in tab. 1. Numbers after "F" denote the design parameters-functions, after "C" -the design parameters-constants with notations: 1 -𝜓, 2 -ℎ, 3 -𝑅 1 . "Geod" denotes the design with continuous geodesic winding.

We used the mass of C123 design as a basis for further comparisons. It was about 22 kg (tab. 1).

The considered design with the continuous geodesic winding has been one of COPV widely used in practice [24]. Its mass was 19 kg while the design F3C12 (the constant winding angles and thickness, the variable geometry) was about 16 kg. The main reason was the circumferential winding near the opening, which did not allow using all the fiber resources. Therefore it was necessary to increase the vessel thickness near the holes.

Comparison of the designs with two constants and one function (F1C23, F2C13, F3C12) shows that the possibility to change the value of parameter along the radius (parameter's variability) is the most critical for the shell geometry (𝑅 1 ) and its thickness.

Fig. 5 presents the design parameters and the vessel's generatrix for several designs.

Important additional design characteristic is its "adaptability in manufacturing". For example, the 5 -10 times difference of thickness along the meridian would become a serious obstacle for vessels manufacturing. Thus, designs of nearly minimum mass possessing good properties and satisfying to the given technological constraints could be of great value, than optimum without them.

Let's consider the following characteristics of the obtained solutions (tab. 1): the difference between the maximal and minimum values of angles (𝜓 max −𝜓 min ), the maximum thickness ratio (to ℎ max /ℎ min ) and the maximum curvature radius ratio (𝑅 1max /𝑅 1min ), which characterizes deviation of the generatrix from a circle arch.

According to tab. 1 the design with the geodesic continuous winding has the thickness ratio about 10 and large gradient near the edge. We found the design F123, where the parameter variance is small -no more than 10% for the thickness and 4 ∘ for the winding angle with low gradients. The designs with the constant thickness (F13C2, F3C12) showed that it was possible to receive the vessels with weight close to minimum, having on hand only such design functions as 𝑅 1 and 𝜓.

Inverse Problem. Verification of the Solutions

We verified the solutions of optimization problem by substituting the obtained design parameters into the direct problem. Fig. 6 shows key characteristics of the stress-strain state for the designs F123, F23C1, C123 and Geod. It is noteworthy that the stress-strain state of design F123 is almost momentless, and the fibers are equally stressed. The influence of transverse shear is minimum.

Let's substitute this solution in the direct problem. All the three considered bs r bs m q q, rad T , 11 MN/m M , 11 N Fig. 7. The stress-strain state characteristics of the vessel with the optimized design functions based on the three shell theories shell theories have yielded close results (fig. 7). The difference is noticeable only for ANST in narrow zones (less than 1% of all the area of calculation) at the edges, where non-linear accounting for transverse shear gives difference of about 5%. At the same time the estimated efforts and bending moments are very close for all the theories, and the bending moments are very small. Thus, it is possible to use the simplest shell theory to solve such optimization problem and the estimation of stress-strain state will be close to those obtained using more complex theories.

Conclusions

A technology of COPV optimization has been developed. It makes possible to obtain high pressure vessel designs that not only meet such requirements as minimum mass, preset volume and strength, but also possess a number of additional valuable engineering characteristics including stress-strain state close to momentless and almost equally stressed fibers.

Non-constant design parameters, such as thickness, winding angles and curvature radius of composite shell give a possibility for additional reduction of COPV mass while keeping its strength. The obtained design with the variable design parameters are up to 27% lighter if compared to the best design with the constant parameters.

The optimization problem solutions have been verified by solving the direct problems with obtained design parameters using the classical shell theory and theories with shear terms.

Our study has demonstrated acceptability and convenience of using simple mathematical models based on the Kirchhoff-Love and Timoshenko shell theories for numerical solving of the optimization problems.