# Modélisation du plissage dans les structures membranaires

Adama Diaby, Anh Le Van, Christian Wielgosz

► To cite this version:

Adama Diaby, Anh Le Van, Christian Wielgosz. Modélisation du plissage dans les structures membranaires. Revue Européenne de Mécanique Numérique/European Journal of Computational Mechanics, 2006, 15 (1-2-3), pp.143 - 154. 10.3166/remn.15.143-154. hal-01007071

HAL Id: hal-01007071

https://hal.science/hal-01007071v1

Submitted on 5 Oct 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Modélisation du plissage dans les structures membranaires

Adama Diaby — Anh Le van — Christian Wielgosz

GéM - Faculté des Sciences de Nantes
2, rue de la Houssinière, BP 92208
F-Nantes 44322 cedex 3

{adama.diaby, anh.le-van, christian.wielgosz}@univ-nantes.fr

RÉSUMÉ. Dans cet article on traite de l'apparition des plis dans les membranes. La mise en œuvre numérique est basée sur la formulation lagrangienne totale et les éléments finis utilisés sont des quadrilatères à 8 nœuds ou des triangles à 6 nœuds, sans rigidité à la flexion. Le matériau constitutif est modélisé par une loi hyperélastique compressible. Le plissage est simulé à l'aide des techniques de calcul de bifurcation pure, sans recours aux imperfections. La méthode classique de longueur d'arc quadratique est modifiée au moyen d'un procédé spécifique pour faire face aux racines complexes qui apparaissent dans l'équation du second degré de longueur d'arc. L'application de la formulation proposée à un ensemble d'exemples numériques typiques montre sa capacité à modéliser correctement les plis dans les membranes.

ABSTRACT. The problem of the wrinkling of the membranes is addressed in this paper within the total Lagrangian formulation. The material obeys a compressible hyperelastic constitutive equation and the resulting finite elements are either 8-node quadrilateral or 6-node triangular membrane elements without flexural stiffness. The appearance of wrinkles in the membrane structures is dealt with using the bifurcation analysis without any imperfections. The standard spherical arclength method is modified by means of a specific solution procedure to cope with the occurrence of complex roots when solving the quadratic constraint equation. Applying the proposed formulation to a set of typical numerical examples shows its ability to correctly predict the wrinkling behaviour in membranes structures.

MOTS-CLÉS : éléments finis, membranes, hyperélasticité, bifurcation, plis.

KEYWORDS: finite element, membranes, hyperelasticity, bifurcation, wrinkles.
# 1. Introduction

L'étude des plis des membranes hyperélastiques non linéaires a été abordée par plusieurs auteurs à l'aide de différentes approches. Les plus utilisées sont les techniques du calcul de la bifurcation, communément désignées par tensile field theory dans la littérature. La théorie du tensile field est généralement appliquée aux éléments membranes sans rigidité à la flexion. Dans cette approche l'état de chaque élément (plié, inactif ou tendu) est déterminé en utilisant un critère d'apparition de plis. Il existe en général trois types de critères : le premier est basé sur les contraintes principales (Contri et al., 1988 ; Tabarrok et al., 1992 et Fujikake et al., 1989), le second sur les déformations principales (Miller et al., 1982 ; 1985), et le troisième sur les contraintes et les déformations principales (Kang et al., 1997 ; Stanuszek, 2003 ; Adler, 2000). Durant le calcul numérique, les propriétés matérielles de la structure sont ajustées pour tenir compte de l'état de chaque élément. Bien que cette approche prédise la formation des plis, elle ne donne aucun renseignement sur l'amplitude, la longueur d'onde ou le nombre d'ondes. Pour prendre en compte les détails des plis, une autre simulation du plissage dans les membranes peut être faite à l'aide des éléments finis de coque. Ces éléments ayant une rigidité non nulle à la flexion, ils facilitent l'étude des plis par l'analyse de la bifurcation. On peut citer entre autres modélisations utilisant cette approche, celle de Wong et al. (2002a ; 2002b).

Très peu d'articles concernent l'étude des plis dans les membranes en couplant technique de la bifurcation et éléments finis membranes sans rigidité de flexion. Les plis sont des déformées qui apparaissent à cause d'un flambement local dans les régions de la structure où se développent des contraintes compressives. L'étude présentée dans ce papier entre dans la lignée de Miyamura (2000) : formulation d'un élément fini membrane sans rigidité à la flexion combinée à l'analyse de la bifurcation pour l'étude des plis dans les structures membranaires. Cependant, la théorie standard de longueur d'arc quadratique (Crisfield, 1991) n'est pas adéquate pour résoudre le problème des plis dans les membranes. En effet, l'apparition des plis dans ces structures est un phénomène très critique qui donne naissance à des racines complexes récurrentes dans l'équation du second degré de longueur d'arc. Puisqu'une simple réduction de la longueur d'arc $\Delta l$ ne résout pas le problème, Lam et al. (1992) ont proposé dans ce cas une reformulation de la longueur d'arc quadratique pour faire aboutir les calculs. La procédure de ces auteurs a été greffée à notre algorithme rendant ainsi possible le traitement des plis.

# 2. Mise en œuvre numérique

## 2.1. Formulation de l'élément fini membrane 3D

L'élément membrane est développé en formulation lagrangienne totale, il a une rigidité nulle à la flexion et il satisfait à la condition de contrainte plane. Le matériau obéit à une loi de comportement hyperélastique compressible.
Considérons une membrane dont la position de référence est définie par la surface S₀ et l'épaisseur h₀. Sa position actuelle est définie par la surface S et l'épaisseur h. La surface de référence S₀ est discrétisée en éléments finis isoparamétriques, on notera e l'élément courant et NNE son nombre de nœuds. Soit eξ l'élément de référence de e, les positions de référence X et actuelle x d'une particule courante de la surface de la membrane sont paramétrées par les coordonnées (ξ₁, ξ₂) selon la chaîne (ξ₁, ξ₂) ∈ eξ ↦→ X ∈ e ⊂ S₀ ↦→ x ∈ S.

La discrétisation du principe des travaux virtuels conduit à un système d'équations non linéaires d'inconnue {U} – vecteur des déplacements nodaux de la membrane – système que l'on résout par la méthode de Newton-Raphson. La matrice raideur tangente élémentaire due aux forces internes dans chaque élément e est calculée par :

$$
\forall i, j \in [1,3_{\text{NNE}}], \quad K_{ij}^{e} = \int_{e} N_{a,\alpha} N_{b,\beta} \left[ \delta_{pq} \boldsymbol{\Sigma}^{\alpha\beta} + x_{p,\gamma} x_{q,\delta} \left( \frac{\partial \boldsymbol{\Sigma}}{\partial \mathbf{E}} \right)^{\gamma\alpha\beta\delta} \right] h_o dS_o \quad [1]
$$

Dans la relation [1], l'intégrale est effectuée sur l'élément *e* de la surface de référence S₀. Les entiers a,b∈[1,NNE] et p,q∈[1,3] sont reliés aux indices i,j∈[1,3NNE] par i = 3(a-1)+p, j = 3(b-1)+q, et il y a sommation sur α,β,γ,δ∈{1,2}.
La notation N<sub>a,α</sub> désigne la dérivée de la fonction de forme N<sub>a</sub> par rapport à ξ<sub>α</sub>, Σ le tenseur des contraintes de Piola-Kirchhoff de seconde espèce, E le tenseur des déformations de Green-Lagrange.

En présence d'une pression p qui est une force suiveuse appliquée sur la
membrane, on ajoutera la matrice suivante à la précédente : ∀i,j ∈ [1,3 NNE],

$$
\begin{equation}
\begin{split}
(K_{\text{pression}}^e)_{ij} = \frac{1}{2} \int_{e_{\xi}} p \left[ \begin{aligned}[t]
& (N_a N_{b,l} - N_b N_{a,l}) x_{p \wedge q,2} - (N_a N_{b,2} - N_b N_{a,2}) x_{p \wedge q,1} \end{aligned} \right] d\xi_1 d\xi_2 \\
& + \frac{1}{2} \oint_{\partial e_{\xi}} p N_a N_b (x_{p \wedge q,2} v_1 - x_{p \wedge q,1} v_2) ds
\end{split}
\tag{2}
\end{equation}
$$

Les entiers a, b∈[1, NNE] et p,q∈[1,3] sont reliés aux indices i,j∈[1,3NNE] de la même manière que ci-dessus. La notation x<sub>p∧q</sub> est définie par x<sub>p∧q</sub> = ε<sub>pqm</sub> x<sub>m</sub> (p et q étant fixés, il y a au plus une valeur de m telle que ε<sub>pqm</sub> soit non nul), v est la normale unitaire extérieure à l'élément de référence e<sub>ξ</sub>, en un point de son frontière ∂e<sub>ξ</sub>. L'intégrale de surface sur e<sub>ξ</sub> est symétrique, alors que l'intégrale curviligne sur ∂e<sub>ξ</sub> est antisymétrique. Si la pression est appliquée sur une portion de surface S<sub>p</sub>⊂S, la contribution de l'intégrale curviligne en un point situé sur la frontière d'un élément intérieur à S<sub>p</sub> s'annule avec celle d'un élément adjacent, si bien qu'il ne reste en fin de compte que la contribution aux points situés sur la frontière de S<sub>p</sub>.

Dans notre étude, d'une part la pression p est uniforme sur S_p, et, d'autre part, la surface S_p est soit à bord fixe, soit coupée en deux par un plan de symétrie. Par conséquent, la contribution des intégrales curvilignes dans la matrice tangente est nulle.
## 2.2. Paramètres de contrôle et méthode de pilotage

On partitionne le vecteur {$\tilde{U}$} des déplacements nodaux de la structure en deux parties : l'une notée {$\tilde{U}$} comprenant les degrés de liberté inconnus, l'autre notée {$\tilde{U}$} comprenant des degrés de liberté imposés. Le vecteur force externe {$\tilde{\Phi}$} de la structure est décomposé sur le même modèle : une partie {$\tilde{\Phi}$} correspondant à {$\tilde{U}$} comprend les composantes de forces imposées, et une autre partie {$\tilde{\Phi}$} correspondant à {$\tilde{U}$} comprend les composantes de forces inconnues (forces de réactions). Dans les exemples numériques, on suppose soit un déplacement imposé proportionnel, soit un chargement proportionnel. Ainsi,

$$ \text{soit } \{\tilde{U}\} = \{\tilde{U}\}^0 + \lambda \{\tilde{U}\}^{\text{ref}} \quad \text{soit } \{\tilde{\Phi}\} = \lambda \{\tilde{\Phi}\}^{\text{ref}} \quad [3] $$

où $\lambda$ est le paramètre de contrôle, {$\tilde{U}$}^ref et {$\tilde{\Phi}$}^ref désignent les quantités de référence données. Le vecteur {$\tilde{U}$}^0 relatif aux déplacements imposés nuls ne change pas la valeur de {$\tilde{U}$}. Dans certains cas, on peut avoir besoin de plusieurs paramètres de contrôle, par exemple dans le cas examiné plus loin de la membrane carrée soumis à un glissement imposé, il faut deux paramètres de déplacement, l'un pour étirer au préalable la membrane et l'autre pour imposer le glissement. Les relations [3a] et [3b] se généralisent facilement à de tels cas.

Pour traiter les points critiques (points limites, points de bifurcation), nous utilisons la méthode de longueur d'arc qui consiste à avancer sur une courbe de réponse d'une longueur d'arc $\Delta l$ donnée, (Wempner, 1971, Riks, 1972 ; 1979). L'équation supplémentaire à celle d'équilibre s'écrit, respectivement dans les cas [3a] et [3b] :

$$ \text{soit } \|\Delta \tilde{U}\|^2 + \Delta\lambda^2 \|\overline{\tilde{U}}^{\text{ref}}\|^2 = \Delta l^2 \quad \text{soit } \|\Delta \tilde{U}\|^2 + \Delta\lambda^2 C_{\text{ref}}^2 = \Delta l^2 \quad [4] $$

Dans la relation [4b], le scalaire $C_{\text{ref}}$ est un facteur d'échelle ayant la dimension d'une longueur; lorsqu'il est pris égal à zéro, on parle de méthode de longueur d'arc cylindrique. L'équation [4] combinée à l'équation d'équilibre de la structure conduit à une équation de second degré en $\Delta\lambda$ (Crisfield, 1991). Dans le calcul des plis, on rencontre très fréquemment des difficultés pour passer sur une branche bifurquée souhaitée, dans la mesure où cette équation de second degré a de façon récurrente des racines complexes. Un remède efficace pour surpasser cette difficulté a été proposé par (Lam et al., 1992). Pour fixer les idées, supposons que nous sommes dans le cas du chargement proportionnel, la méthode de Lam et Morley consiste alors à projeter le résidu sur le vecteur force externe : à l'itération où apparaissent les racines complexes, le résidu est scindé en une composante selon la direction de la
force et en une autre selon la perpendiculaire à cette force. C'est cette dernière qui est responsable des racines complexes et qui doit être éliminée. L'algorithme modifié est résumé ci-dessous :

Boucle des itérations

(...)

Résoudre l'équation quadratique de la méthode standard de la longueur d'arc.

Si les racines sont complexes, alors

(i) Réécrire la correction $\delta U$ du déplacement comme $\delta U = x \delta U_H + \eta \delta U_\Phi$, où
$$ \delta U_H = -K_T^{-1} H, $$ (K_T est la matrice raideur tangente, H est la composante orthogonale de la projection du résidu selon la force de référence $\{\tilde{\Phi}\}^{\text{ref}}$),
$$ \delta U_\Phi = -K_T^{-1} \Phi, $$ et x et $\eta$ sont des scalaires.

(ii) Le nouveau déplacement $\delta U$ conduit à une nouvelle équation de second degré en x. On choisit parmi les deux racines de cette nouvelle équation celle qui fait avancer sur la courbe dans le sens souhaité.

Fin si.

(...)

Fin de boucle des itérations.

Pour être sûr que l'équation quadratique en x fournit des racines réelles, le scalaire $\eta$ doit être pris à 5 % de $|\eta_2 - \eta_1|$, où $\eta_1$ et $\eta_2$ sont les racines d'une nouvelle équation quadratique en $\eta$.

## 2.3. Analyse de bifurcation

Lorsqu'on parcourt une courbe d'équilibre donnée, on peut traverser un ou plusieurs points critiques qui sont des points singuliers de la matrice raideur tangente $K_T$. Une détection pratique et à moindre coût des points critiques est basée sur la factorisation de Crout $K_T = LDL^T$, où la matrice L est triangulaire inférieure et où la matrice D est diagonale. Le nombre de valeurs propres négatives de $K_T$ est égal au nombre d'éléments diagonaux négatifs, appelés « pivots », de D.

Lorsqu'un point critique sur la branche actuelle est détecté, on revient en arrière et on examine sa nature pour savoir s'il s'agit d'un point limite ou d'un point de bifurcation. Pour cela, la longueur d'arc $\Delta l$ est recalculée plusieurs fois par une méthode de type de dichotomie comme celle proposée dans (Shi et al., 1996). Une manière efficace pour distinguer un point limite d'un point de bifurcation consiste à calculer le paramètre k de (Bergan et al., 1978), appelé *current stiffness parameter* et défini dans le cas d'un chargement proportionnel par $k = (\tilde{\Phi}^{\text{ref}} \cdot \tilde{U}_{\phi \text{ref}})/(\tilde{U}_{\phi \text{ref}} \cdot \tilde{U}_{\phi \text{ref}})$
où $\tilde{\mathbf{U}}_{\phi\text{ref}} = \tilde{\mathbf{K}}_T^{-1} \tilde{\mathbf{\Phi}}^{\text{ref}}$, $\tilde{\mathbf{K}}_T$ est la matrice carrée extraite de la matrice tangente $\mathbf{K}_T$ selon le mode de partition $\mathbf{U} = (\tilde{\mathbf{U}}, \bar{\mathbf{U}})$ décrit plus haut. Le signe du paramètre k change au passage d'un point limite, alors qu'il reste inchangé au passage d'un point de bifurcation.

Le branchement sur une nouvelle courbe d'équilibre s'effectue par injection des vecteurs propres (Wagner et al., 1988). Au début de l'incrément de la nouvelle branche, on calcule le vecteur propre $\tilde{\mathbf{Z}}$, solution de $\tilde{\mathbf{K}}_T \cdot \tilde{\mathbf{Z}} = \mathbf{0}$ et on effectue la première correction $\delta\lambda = 0$ et $\delta \tilde{\mathbf{U}} = \pm\Delta\ell \frac{\tilde{\mathbf{Z}}}{\|\tilde{\mathbf{Z}}\|}$.

En couplant les techniques de calcul de la bifurcation à la procédure de Lam et Morley pour le traitement des racines complexes, on arrive à simuler l'apparition des plis dans les membranes, comme montrent les exemples numériques proposés ci-dessous.

# 3. Exemples numériques

## 3.1. Hémisphère gonflé puis appuyé

Le premier exemple traité porte sur une membrane hémisphérique, de rayon de référence $R_o = 15.5$ cm et d'épaisseur de référence $h_o = 25$ µm ($R_o/h_o = 6102$), gonflée puis appuyée au pôle par une force verticale. Vu la symétrie du problème, seul un quart de l'hémisphère est modélisé avec 25×25 éléments. La loi de comportement utilisée ici comme dans tous les autres exemples est la loi hyperélastique néo-hookéenne compressible définie par le potentiel volumique
$$ \rho_o \psi = \frac{\mu}{2} (\text{tr} \mathbf{C} - 3) - \mu \ln \mathbf{J} + \frac{\lambda}{2} (\ln \mathbf{J})^2, \text{ où } \mathbf{C} \text{ est le tenseur de dilatation lagrangien} $$
$$ \mathbf{C} = 2\mathbf{E} + \mathbf{I}, \quad \mathbf{J} = \sqrt{\det \mathbf{C}}, \quad \lambda \text{ et } \mu \text{ sont des coefficients matériels. Pour les applications numériques, on prendra } \mathbf{E} = 270 \text{ MPa, } \nu = 0.4, \mathbf{E} \text{ et } \nu \text{ étant reliés à } \lambda \text{ et } \mu \text{ par } \lambda = \frac{\text{Ev}}{(1+\nu)/(1-2\nu)}, \mu = \frac{\text{E}}{2/(1+\nu)}. $$

Le chargement est appliqué en deux étapes : d'abord la membrane est gonflée jusqu'à une pression donnée, ensuite on applique une force verticale au pôle.

Dans un premier temps nous décidons de traiter ce problème sans bifurcation, c'est-à-dire de rester sur la branche fondamentale. Comme le montre la figure 1, le pôle s'enfonce dans l'hémisphère sans exhiber le moindre pli.

Ainsi seule la prise en compte de la bifurcation révèle la présence des plis comme le montrent les figures 2a et 2b. La figure 2a présente la déformée de la membrane pour une pression de $p = 120$ Pa (donnée selon (Szyszkowski et al., 1987) et la figure 2b pour une pression de $p = 5000$ Pa (donnée selon (Stanuszek, 2003)).
On remarque sur les deux figures 2a et 2b que dans le voisinage du point d'application de la force, les plis ont une distribution périodique et sont symétriques. Sur la figure 2a, où $p = 120$ Pa, on voit que les plis sont très fins comme dans les tests expérimentaux menés par (Szyszkowski *et al.*, 1987). Par contre, sur la figure 2b, avec une pression intérieure plus importante, $p = 5\,000$ Pa, les plis sont plus accentués, comme cela a été observé dans (Stanuszek, 2003).

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Déformée de la membrane hémisphérique sans calcul de la bifurcation</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Déformées de l'hémisphère pincé : (a) $p = 120$ Pa ; (b) $p = 5\,000$ Pa</figcaption></figure>

Enfin, il convient de remarquer que mailler un quart de la structure mène à une sélection des modes et qu'on ne peut donc détecter que les modes qui respectent la symétrie exploitée. Cette remarque vaut aussi pour l'exemple de l'airbag plus loin où nous maillons un huitième de la structure.
## 3.2. Membranes carrées ou rectangulaires soumises à un glissement imposé

Les exemples présentés dans cette section concernent les membranes carrées ou rectangulaires, fixées sur un côté et soumises à un glissement imposé sur le côté opposé. Commençons par une membrane carrée qui a, dans la configuration initiale, un côté de 25 cm et une épaisseur de 25 µm. Les propriétés matérielles sont $E = 2\,500\ \text{MPa}$ et $v = 0.34$. Ici, le pilotage de la longueur d'arc en force (relation [3b]) entraîne de sérieux problèmes de convergence et doit être remplacé par le pilotage de la longueur d'arc en déplacement (relation [3a]). La membrane entière est maillée avec $40 \times 40$ éléments. Dans la première étape de calcul, on applique un étirement selon la direction y de 0.5 % du côté, puis dans la deuxième étape on impose le glissement selon x jusqu'à 5 % du côté. Le processus de formation des plis dans la membrane carrée cisaillée est décrit sur la figure 3. Au début, pour un glissement assez faible, plusieurs plis apparaissent simultanément dans la structure et s'étendent sur toute la surface quasiment. Au fur et à mesure que le glissement augmente, le nombre de plis diminue et les plis se réorganisent selon la diagonale du carré.

<figure><img src="image_4.png" /><figcaption>**Figure 3.** Déformées successives de la membrane carrée sous glissement imposé</figcaption></figure>

Les plis obtenus pour un glissement de 1.4 cm sont très proches de ceux obtenus par Tessler *et al.*, 2003) : selon la diagonale, on observe la présence de quatre gros plis, et la naissance de deux petits plis à peine perceptibles. Nous avons représenté sur la figure 4 les branches fondamentale et bifurquée de la membrane carrée sous glissement imposé. Le point choisi pour tracer la flèche w est de coordonnées (0 cm, 12.5 cm, 0 cm), l'origine étant située en bas à gauche de la membrane.
La branche fondamentale correspond à la déformation de la membrane dans son plan, en l'absence des plis. Sur la branche bifurquée on remarque que lorsqu'on soumet la membrane à un glissement imposé, le point de bifurcation n'apparaît pas avant que le déplacement horizontal du bord supérieur ne dépasse une certaine quantité. Ceci est dû à l'étirement initial apporté à la membrane avant le glissement.

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Courbes d'équilibre d'une membrane carrée sous glissement imposé (branche fondamentale et branche bifurquée)</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** Déformées successives de la membrane rectangulaire sous glissement imposé</figcaption></figure>

Passons maintenant au cas de la membrane rectangulaire dont le processus de formation des plis est, comme on va le voir, différent de celui de la membrane carrée. Les dimensions de référence retenues ici sont : longueur L = 40 cm, largeur $l$ = 16 cm et épaisseur = 50 µm. Le module d'Young du matériau est E = 6.46 GPa
et le coefficient de Poisson $v = 0.35$. Comme pour la membrane carrée, on doit appliquer un étirement initial qui vaut ici 3 % de la largeur $l$.

La figure 5 montre qu'avec la membrane rectangulaire, les plis prennent d'abord naissance sur les bords de la structure. Au fur et à mesure que le glissement augmente, les plis convergent vers le centre de la membrane.

### 3.3. Airbag circulaire

Le dernier exemple traité est la modélisation des plis dans un airbag circulaire. Dans la configuration de référence, l'airbag est une membrane d'épaisseur $h_o$ composée de deux faces circulaires accolées, de rayon $R_o$ et se trouvant dans un même plan. Les valeurs numériques utilisées sont celles de Cirak *et al.*, 2001) et reprises par (Flores *et al.*, 2003) : module d'Young $E = 60$ MPa, coefficient de Poisson $v = 0.3$, $h_o = 400$ µm et $R_o = 35$ cm. Bien que le gonflement d'un airbag soit en pratique un problème dynamique, nous le résolvons ici dans le cadre quasi statique.

En tenant compte des symétries, seul un huitième de l'airbag est maillé avec 30×30 éléments. L'airbag est gonflé jusqu'à une pression de 5 000 Pa. Sous l'effet de la pression, l'airbag se déforme en prenant des flèches hors plan, pendant que son bord initialement circulaire se déplace librement dans le plan initial. La forte singularité numérique au démarrage du calcul est surmontée en étirant au préalable la membrane par une force morte appliquée sur le bord circulaire de l'airbag et éliminée au fur et à mesure de la progression du calcul. Ce procédé fait apparaître naturellement les plis sans avoir recours aux techniques de calcul de la bifurcation. Le résultat obtenu est montré sur la figure 6.

<figure><img src="image_7.png" /><figcaption>**Figure 6.** Déformée de l'airbag circulaire</figcaption></figure>

On observe une distribution régulière des plis le long de la grande circonférence de la membrane. Il est évident qu'il faut utiliser un maillage plus fin si on veut représenter correctement le nombre des plis, comme dans les autres exemples d'ailleurs.
## 4. Conclusion

L'objectif de ce travail a été de proposer un outil numérique robuste pour étudier l'apparition des plis dans les membranes. A cette fin, un élément fini membrane sans rigidité à la flexion a été développé en formulation lagrangienne totale et utilisé en analyse de bifurcation pure sans ajout d'aucunes imperfections dans la structure. Les exemples sélectionnés ont montré la capacité de la formulation à prédire la bifurcation et les régions de plis dans la membrane. La singularité de la matrice dans la configuration de référence est surpassée par une précontrainte, qui est soit une précontrainte introduite dans la loi de comportement, soit une force ou un déplacement imposé sur les bords de la structure. Une technique simple proposée par (Lam et al., 1992) a été incorporée dans le programme de calcul pour traiter l'apparition éventuelle de racines complexes lorsqu'on suit la méthode standard de longueur d'arc quadratique.

## 5. Bibliographie

Adler A.L., Finite element approaches for static and dynamic analysis of partially wrinkled membrane structures, PhD Thesis, Departement of Aerospace Engineering, University of Colorado, Boulder, CO, USA, 2000.

Bergan P.G., Horrigmoe G., Krakeland B., Soreide B., « Solution techniques for non-linear finite element problems », *International Journal for Numerical Methods in Engineering*, vol. 12, 1978, p. 1677-1696.

Cirak F., Ortiz M., « Fully $C^1$-conforming subdivision elements for finite deformations thin-shell analysis », *International Journal for Numerical Methods in Engineering*, vol. 51, 2001, p. 813-833.

Contri P., Schrefler B.A., « A geometrically nonlinear finite element analysis of wrinkled membrane surface by no-compression material model », *Communication in Applied Numerical Methods*, vol. 4, 1988, p. 5-15.

Crisfield M.A., *Non-linear Finite Element Analysis of Solids and Structures*, Volume 1, John Wiley & Sons, 1991.

Flores F.G., Onate E., « Rotation-free element for finite strain analysis of elastic-plastic thin shell and membrane structure », *Textile composites and inflatable structures*, Onate E and Kröplin B (Eds.), CIMNE: Barcelona, 2003, p. 148-153.

Fujikake M., Kojima O., Fukushima S., « Analysis of fabric tension structures », *Computers & Structures*, vol. 32, 1989, p. 537-547.

Kang S., Im S., « Finite element analysis of wrinkling membranes », *Journal of Applied Mechanics*, vol. 64, 1997, p. 263-269.

Lam W.F., Morley C.T., « Arc-length method for passing limit points in structural calculation », *Journal of Structural Engineering*, vol. 118, 1992, p. 19-185.
Miller R.K., Hedgepeth J.M., « An algorithm for finite element analysis of partly wrinkled membranes », *AIAA*, vol. 20, 1982, p. 1761B1763.

Miller R.K., Hedgepeth J.M., Weingarten V.I., Das P., Kahyai S., « Finite element analysis of partly wrinkled membranes », *Computers & Structures*, vol. 20, 1985, p. 631B639.

Miyamura T., « Wrinkling on stretched circular membrane under in-plane torsion: birfurcation analyses and experiments », *Engineering Structures*, vol. 23, 2000, p. 1407-1425.

Riks E., « The application of Newton's method to the problem of elastic stability », *Journal of Applied Mechanics*, vol. 39, 1972, p. 1060-1066.

Riks E., « An incremental approach to the solution of snapping and buckling problems », *International Journal of Solids and Structures*, vol. 15, 1979, p. 529-551.

Shi J., Moita G.F., « The post-critical analysis of axisymmetric hyper-elastic membranes by the finite element method », *Computer Methods in Applied Mechanics and Engineering*, vol. 135, 1996, p. 265-281.

Stanuszek M., « FE analysis of large deformations of membranes with wrinkling », *Finite Elements in Analysis and Design*, vol. 39, 2003, p. 599-618.

Szyszkowski W., Glockner G., « Spherical membranes subjected to vertical concentrated loads: an experimental study », *International Journal of Engineering Structures*, vol. 9, 1987, p. 183-192.

Tabarrok B., Qin Z., « Nonlinear analysis of tension structures », *Computers & Structures*, vol. 45, 1992, p. 973B984.

Tessler A., Sleight D.W., Wang JT., « Nonlinear shell modelling of thin membranes with emphasis on structural wrinkling », *44th AIAA/ASME/ASCE/AHS Structures, Structural dynamics, and Materials Conference*, Nortfolk, Virginia, USA, AIAA 2003-1931.

Wagner W., Wriggers P., « A simple method for the calculation of post-critical branches », *Engineering and Computations*, vol. 5, 1988, p. 103-109.

Wempner GA., « Discrete approximations related to nonlinear theories of solids », *International Journal of Solids and Structures*, vol. 7, 1971, p. 1581-1599.

Wong Y.W., Pellegrino S., « Amplitude of wrinkles in thin membranes », *In New approaches to structural mechanics shells and biological structures*, (Dordrecht, The Netherlands, 2002a), Kluwer Academic Publishers, 257-270.

Wong Y.W., Pellegrino S., « Computation of wrinkle amplitudes in thin membranes », *43rd AIAA Structures, Structural dynamics, and Materials Conference*, Denver, Co, USA, AIAA-2002b-1369.