Finite elements for 2D problems of pressurized membranes

This paper presents theoretical and numerical developments of finite elements for axisymmetric and cylindrical bending problems of pressurized membranes. The external loading is mainly a normal pressure to the membrane and the developments are made under the assumptions of follower forces, large displacements and finite strains. The numerical computing is carried out in a different way that those used by the conventional finite element approach which consists in solving the non-linear system of equilibrium equations in which appears the stiffness matrix. The total potential energy is here directly minimized, and the numerical solution is obtained by using optimization algorithms. When the derivatives of the total energy with respect to the nodal displacements are calculated accurately, this approach presents a very good numerical stability in spite of the nil bending rigidity of the membrane. Our numerical models show a very good accuracy by comparisons to analytical solutions and experimental results.


Introduction
The finite element method is widely used to solve a lot of engineering problems. Its formulation leads to the resolution of a non-linear system of equations. In the case of membrane structures, which present a singular stiffness matrix due to the loss of bending stiffness, the numerical solution obtained by using shell elements cannot be straightforward because the convergence of the iterating process is difficult to reach. The numerical method, which is built in this paper, is based on a direct minimization of the total potential energy to compute the state of equilibrium of pressurized membranes.
The literature on the finite element analysis of the pressurized membrane structures is sparse. We can quote the papers of Main et al. [8], Kawabata and Ishii [7]. Bonet [2] presented a finite element analysis of closed membrane structures that contain a constant mass fluid such as air. The membrane formulation presented by the author avoids the need for local co-ordinate axes by using the isoparametric finite element plane as a material reference configuration. Full linearization of the internal pressure forces and the resulting additional term in the tangent operator are derived. General finite element membrane analyses have been presented by Argyris et al. [1], Gruttman and Taylor [5], De Souza Neto et al. [10], Oden and Sato [9].
In the first section of this paper, we will come back on the main analytical results concerned with inflated axisymmetrical pressurized membranes. The second section is devoted to the construction of the axisymmetric finite element. Numerical results are compared with analytical ones and show the good accuracy of our finite elements. In the third section, a cylindrical bending finite element is developed, and the numerical solution is compared to experimental results on inflatable flat panels inflated at high pressure. Once again a good accuracy between experimental and numerical results is obtained.

Analytical solutions for inflated circular membranes
The main bibliographical results are about circular inflated membranes clamped at their rim (Fig. 1). The membrane is supposed to have large deflections and is not pre-tensioned at its initial state (Henchy's problem [6]). This bibliographical recall is mainly based on Fichter's paper [4]. Campbell [3] proposed a solution to Hencky's problem by taking into account the case of a pre-tensioned circular membrane.

Governing equations
The equilibrium equation depends on the assumption made on the effects of the external pressure. If the pressure remains vertical during the inflation, Hencky [6] shows that the radial equilibrium leads to: When the pressure remains orthogonal to the membrane during the inflation, Fichter [4] shows that this equation is rewritten as: The lateral equilibrium leads to the following equation: The behavior is assumed linear and elastic The strain-displacement relation takes into account the non-linear terms: Only the large deflection are here taken into account. The other quadratic terms of the strains are supposed to be small.

Hencky's theory
Hencky [6] searched the solution for the a-dimensional deflection W ðqÞ and stress resultant N ðqÞ in the form of powers series where q represents the a-dimensional radius.
Nomenclature a radius of the circular membrane a 2n , b 2n Hencky's power series coefficients D ij terms of the elasticity matrix e r ! , e z ! unit vectors E, m elasticity modulus and Poisson's ratio E Ã ¼ Eh membrane modulus (product of the elasticity modulus by the membrane thickness) f Xi , f Yi components of the nodal external load F nodal external load h thickness of the membrane E full Green-Lagrange tensor H displacement gradient tensor K linear stiffness tensor i.e. Hessian matrix of the total potential energy L, M third and fourth order displacement derivative of the total potential energy l 0 , l initial length and elongation for the cylindrical bending finite element n 2m , w 2m Fichter's power series coefficients N dimensionless meridian stress resultant (N rr = Eh in Hencky's problem, N rr =pa in Fichter's approach) N rr , N hh radial and lateral stresses resultants, N rr ¼ r rr h, N hh ¼ r hh h p inflating pressure r radius u radial displacement w, W ¼ w=a deflection, a-dimensional deflection U F external force work for the cylindrical bending finite element U p pressure work W strain energy P potential energy q a-dimensional radial co-ordinate, q ¼ r=a r rr , r hh radial and lateral stresses By using Eqs. (1), (4) and (5), the following equations are obtained: Fichter [4] shows that by taking into account the boundary condition wðaÞ ¼ 0 in the preceding equations, it is possible to calculate the Hencky's coefficients a 2n , b 2n by expanding the coefficients of the power series. The values until n ¼ 10 are given in Table 1.

Fichter's theory
Fichter [4] has also searched the solution in the form of powers series: W ðqÞ ¼ By using the same way as Hencky's approach, one can find that [4]: The coefficients n 2m , w 2m are once again obtained from the boundary conditions and can be found in [4]. Note that the undetermined coefficient n 0 depends now both on m and on pa=E Ã . Fig. 2 shows the lateral deflection obtained for each of these analytical methods for three levels of the pressure. The characteristics of the considered membrane are: membrane modulus E Ã ¼ 311 488 Pa m, Poisson ratio m ¼ 0:34, radius a ¼ 0:1425 m. A difference between the deflection shapes appears for medium and high pressure levels (250 and 400 kPa). The differences between these two formulations come from the radial component of the pressure, which was neglected in Hencky's problem [6]. When the pressure is supposed to follow the shape, the deflection is more similar to a spherical one.

An axisymmetrical membrane finite element
In this section, we present the development of an axisymmetrical finite element for the membranes. The developments are made under the assumptions of large deflections, finite strains and follower forces due to the applied pressure. We use the total Lagrangian formulation to represent the kinematics of the deflections. They are described by the mid-surface in the plane containing the meridian curve. We will present the theoretical developments and its numerical implementation.

Theoretical foundations
At its initial position the membrane is circular and flat so that, the geometry is defined by only the radius Table 1 Coefficients of the power series 3 value. The displacement vector between the current position and the initial shape is defined by: The radial u and lateral w components are assumed to depend only on the radius r. Fig. 3 shows the initial and the current positions of the membrane and the displacement vector. The displacement gradient tensor is written in cylindrical co-ordinates and takes into account the simplifications due to the axisymmetrical configuration: The Green-Lagrange tensor in the case of finite strains is defined by: Only three components of the strain tensor have nonnull values One easily identifies the contributions of the three hypothesis made in continuum mechanics: small strains and deflections ðE rr % u ;r ; E hh % u=rÞ, finite strains ðu 2 ;r =2; ðu=rÞ 2 =2Þ, and large deflections ðw ;r =2; w 2 ;r =2Þ. The hypothesis of finite strains and so the use of the strain tensor of Green-Lagrange implies that the stress tensor will be the second tensor of Piola-Kirchoff. The material is supposed to be isotropic and its constitutive law is: The strain energy is given by: We can summarize this result by defining the following energies components: The total strain energy is then: We use linear shape functions for the lateral and radial deflections: wðrÞ ¼ r 2 À r r 2 À r 1 r À r 1 Their derivatives with regards to the r parameter lead to: These two quantities are constant within an elementary volume, so that the elementary strain energies W d ij can be written as (for a finite element between the radii r 1 and r 2 ): 4 The two terms u ;r and w ;r are independent of r, and will be considered as constant during the integration. Therefore we will keep this notation in the following results: W d rh ¼ phD rh ðu ;r þ u 2 ;r =2 þ w 2 ;r =2Þ Â ðu ;r " þ u 2 ;r =2Þðr 2 2 À r 2 1 Þ þ 2ðu 1 r 2 À u 2 r 1 Þð1 þ u ;r Þ þ u 1 r 2 À u 2 r 1 r 2 À r 1 2 ln Note: a singularity arises when one of the radii of the finite element is equal to zero. This singularity disappears when we assume that the strains and stresses cannot be infinite when r ¼ 0. We will therefore simply write that W d hh ¼ W d rr and W d rh ¼ 0 when the radius vanishes.

Pressure work
We have to calculate the work of the applied pressure p: Which can be written as follows: Here n ! is the normal to the surface dS And for an element in its deformed configuration: The formulation is written on the deformed state of membrane and takes into account the real effects of the pressure. The pressure work is exactly written for the complete structure.

Theorem of the potential energy
The potential energy is defined as the difference between the strain energy and the work of the external loads (here the pressure work). We must evaluate the potential energy for each finite element modelling the membrane and then add all these elementary energies to obtain the complete potential energy of the structure: Here i is the finite element index. The displacement solution is obtained by minimizing the total potential energy with respect to each node component. For each element defined between the two nodes 1i and 2i, four non-linear equations have to be solved: We can show that this system of non-linear equations involves until third degrees of displacement components. It can be written in a more usual way as used in finite element method: where K is a second order tensor (Hessian matrix of P, or stiffness matrix), L, M third and fourth order tensors, U T ¼ fu 1i ; w 1i ; u 2i ; w 2i g is the vector of the nodal displacements of one element and F is the vector of the external load applied to the membrane. We give here the components of the stiffness matrix K. L and M tensors are given in Appendix A. The stiffness tensor K can be decomposed as two stiffness tensors, one depending on the mechanical characteristics of the membrane K 0 and the other depending on the pressure K p : With

Mesh sensitivity
Many computations with different meshes are compared in order to evaluate the solution convergence. Figs. 4 and 5 show the mesh sensitivity of the central deflection and the potential energy for a loading pressure equal to 250 kPa. Fig. 4 shows the evolution of the central deflection as a function of the number of axisymmetric finite elements, compared to analytical solutions. The deflection is quite stable beyond 25 elements and the numerical solution seems to be an average of the two analytical solutions.

Comparison between Fichter's analytical model and the finite element solution
Several numerical tests were performed to validate our numerical model. We have not found any analytical reference with a complete model including finite strains, large deflections and follower pressure. We will therefore compare our numerical results with Fichter's theoretical solution [4]. The main difference between these two solutions comes from the hypothesis of finite strains. Fig. 6 shows the different deflections obtained by the two models for three levels of the pressure: 100, 250 and 400 kPa. The mechanical properties of the membrane are E Ã ¼ 311 488 Pa m for the membrane modulus, m ¼ 0:34 for the Poisson's ratio and a ¼ 0:1425 m for the radius. At medium pressure (100 kPa), the deflections given by the two solutions are almost identical because the strains remain small. When the pressure is higher (400 kPa), a difference between the deflection shapes is quite apparent. This pressure value gives high membrane stresses, and the finite strains imply a less spherical shape, but a higher central deflection. The finite element results differ from the analytical model, which is less accurate for   higher pressure levels. The main differences appear on the curvature of the membrane. Fig. 7 shows the radial displacements of the nodes of the discretized structure for pressures varying from 100 to 400 kPa. Even if these displacements are not very high (max of 4.53 mm for 400 kPa) for this size of the membrane, it is important to know their level.

Relative error distribution
We define the relative error between the numerical and Fichter's analytical solution by: where w FE is the deflection obtained with the finite element and w Fichter the deflection given by Fichter's solution. The differences increase with the pressure, as shown in Fig. 8. The maximum difference does not stand in the middle of the membrane but near the rim. Fichter's deflection shows a more spherical shape than the numerical solution according to the error distribution. The maximum relative error reaches 7.3% for a pressure of 400 kPa. These results show that finite strains have an influence when the level of pressure increase and that they cannot be neglected.

Formulation
The formulation of another finite element based on the same energetic approach is now displayed. This element is devoted to cylindrical bending inflatable membranes. It also takes into account finite strains, large deflections and follower forces. The element is defined between the nodes 1 and 2 in the ðx; yÞ plane, and has a unit width along the z-axis. The co-ordinates are named ðx 1 ; y 1 Þ and ðx 2 ; y 2 Þ (Fig. 9).

Strain energy
The strain energy is defined by: Let us denote by l 0 the initial length and by l the length after deformations: The square strain is: The strain energy is be obtained by integration along the element:

Pressure work
The pressure work is equal to the product of the applied pressure by the volume generated by the displacements of the element: The work of the applied forces on the structure nodes is defined by the dot product of the applied forces by the displacement vector:

Potential energy
The total potential energy of one element is the difference between the strain energy and the work of external forces applied to the element: The total potential energy of the structure is calculated by adding the elementary potential energy of each element.

Comparison with experimental results
The validation of this cylindrical bending finite element is done by the comparison with experimental tests made on double-layered pressurized panels.

Double-layered panel
The tested panels are prototypes constructed by TISSAVEL Inc, and made of two parallel-coated woven fabrics. The upper and lower layers are linked with flexible yarns [11]. The yarn density is strong enough to ensure the flatness of the pressurized structure. Fabrics are made with high strength polyester material. The behavior of the panel depends on the inflation pressure, which generates the pre-stress of the fabrics and of the yarns.

Experimental tests
Panels have been tested like beams submitted to bending loads and with various boundary conditions. Fig. 10 shows the experimental device we used to test the panels. Two kinds of boundary conditions have been used: isostatic conditions (simply supported at the two ends), and hyperstatic conditions (clamped at one end and simply supported at the other end). Analytical solutions are available in [12] for the isostatic case.

Numerical results
The cylindrical bending finite element is used to model the upper and lower layer of the panel, and high stiffness elements are used to take into account the vertical yarns linking the two layers. Computations have been done using a membrane stiffness modulus equal to E Ã ¼ 650 000 Pa m evaluated by an inflation test [12]. The length between the supports is 1.6 m for both devices. First, we pressurize the panels to make them strong and then we apply a concentrated load F at the middle of the panel. Fig. 11 shows comparisons for the isostatic case between the numerical model and the experiments for three levels of pressure and for three values of the load F . Fig. 12 shows the comparisons for the hyperstatic case. These two comparisons show the accuracy of the cylindrical bending finite element. Numerical and experimental comparisons for a higher hyperstatic case are available in [13].

Conclusions
We have described two finite elements for 2D problems of inflatable membranes: axisymmetric and cylindrical bending one. The elements are built with the   hypothesis of large deflections, finite strains and with follower pressure loading. A new formulation is used to construct these two elements: the numerical solution is obtained by solving directly the optimization problem formulated by the theorem of the minimum of the total potential energy. The stiffness matrix appears as the Hessian matrix of the energy and depends explicitly on the inflation pressure. These two elements are compared with analytical and experimental results, and in both cases, numerical results are close to reference results, which prove the accuracy of the proposed elements. The energy formulation and the solution techniques for more general cases are straightforward. The numerical technique is only limited by the size of the optimization problem. A further modelling with a 3D triangular finite element based on the same energetic formulation is in progress, and we would be able to predict the deflections of 3D pressurized membrane structures.