# Formulation théorique et évaluation d'un modèle d'élément fini discrêt pour les plaques composites. Le modèle DMQPml

Sakami Siham, Ayad Rezak, Sabhi Hamid

► To cite this version:

Sakami Siham, Ayad Rezak, Sabhi Hamid. Formulation théorique et évaluation d'un modèle d'élément fini discrêt pour les plaques composites. Le modèle DMQPml. 9e Colloque national en calcul des structures, CSMA, May 2009, Giens, France. hal-01413976

HAL Id: hal-01413976

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

Submitted on 12 Dec 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
# Formulation théorique et évaluation d'un modèle d'élément fini discrêt pour les plaques composites. Le modèle DMQPml

Sakami siham*, Ayad rezak**, Sabhi hamid*

13 mars 2009

(*) Laboratoire d'Analyse des Contraintes Mécaniques-Dynamique des transferts aux interfaces (LACM/DTI). EA3304 URA CEA-Université LRC DSM 0534
UFR Sciences Exactes et Naturelles Moulin de la Housse - BP 1039 51687 REIMS
siham.sakami@univ-reims.fr

(**) Groupe de Recherche en Sciences pour l'Ingénieur (GRESPI)
UFR Sciences Exactes et Naturelles Moulin de la Housse - BP 1039 51687 REIMS
rezak.ayad@univ-reims.fr

**Résumé** – Dans cet article, nous présentons la formulation théorique et l'évaluation d'un nouvel élément fini du premier ordre pour les plaques composites multicouches. Il est basé sur un modèle variationnel en déplacement que nous considérons comme discrêt, dans la mesure où l'on introduit localement et de manière discrête des hypothèses cinématiques et mécaniques. Ce modèle que nous appelons DDM (Discrete Displacement Mindlin) fournit un élément fini géométriquement simple (un quadrilatère à 4 noeuds et 3 ddl par noeud) et efficace par l'existence de courbures de flexion linéaires, issues d'une représentation quadratique des rotations de la normale à la surface moyenne. Nous avons introduit localement des hypothèses modifiées de Mindlin dans l'élément DKMQ proposé par Katili, pour prendre en compte le caractère multicouche des composites. Il prend en compte les effets de CT à travers l'épaisseur et reproduit les résultats de plaques minces, quelque soit l'élancement (absence de verrouillage en CT).

mots clés – Eléments finis spéciaux, plaques multicouches, Modèle DDM

## 1 Présentation de l'élément DMQPml :

Le nouvel élément de plaque, baptisé (Discrete Mindlin Quadrilateral Plate multilayer), est une extension aux cas multicouches du modèle isotrope DKMQ proposé par Katili [2]. L'élément DMQPml (Fig. 1) possède quatre nœuds et 3ddl par nœud ($w_i, \beta_{xi}, \beta_{yi}$) : le déplacement transversal suivant l'axe z et les deux rotations de la normale à la surface moyenne dans les deux plans $x-z$ et $y-z$.

## 2 Champ des déplacements :

Le champ des déplacements virtuels correspondant, supposé petit entre les deux configurations $C^0$ et $C$, est défini par :
<figure><img src="image_2.png" /><figcaption>FIG. 1: Géomètrie de l'élément DMQP</figcaption></figure>

$$ \vec{u}_q^*(x, y, z) = \left\{ \begin{array}{c} \vec{u}^*(x, y, z) \\ \vec{v}^*(x, y, z) \\ \vec{w}^*(x, y, z) \end{array} \right\} = \left\{ \begin{array}{c} 0 \\ 0 \\ \vec{w}^*(x, y) \end{array} \right\} + z \left\{ \begin{array}{c} \vec{\beta}_x^*(x, y) \\ \vec{\beta}_y^*(x, y) \\ 0 \end{array} \right\} \quad (1) $$

$w, \beta_x, \beta_y$ sont respectivement le déplacement transversal et les deux rotations de la normale dans les plans $x-z$ et $y-z$. Nous proposons d'approcher les rotations $\beta_x$ et $\beta_y$ à l'aide d'une interpolation quadratique incomplète qui fera apparaître des accroissements de rotations $\Delta\beta_{sk}$ sur les quatre bords élémentaires. $\Delta\beta_{sk}$ sont les variables associées à la représentation quadratique de $\beta_x$ et $\beta_y$.

$$ \left\{ \begin{array}{c} \beta_x \\ \beta_y \end{array} \right\} = \sum_{i=1}^{4} N_i \left\{ \begin{array}{c} \beta_{xi} \\ \beta_{yi} \end{array} \right\} + \sum_{k=5}^{8} P_k \left\{ \begin{array}{c} C_k \\ S_k \end{array} \right\} \Delta\beta_{sk} \quad (2) $$

## 3 Champs des déformations

On se limitera au cas d'une plaque de Reissner/Mindlin en flexion et CT, sans la prise en compte du comportement membranaire. Seules les courbures de flexion {$\chi$} et les déformations de CT seront considérées dans l'étude. Pour plus de détail sur l'aspect théorique du modèle, l'auteur pourra consulter [5].

### 1. Courbures de flexion :

$$ \{\chi\} = [B_f] \left\{ \begin{array}{c} U_n \\ \Delta\beta_n \end{array} \right\} = [B_{b\beta}] \{U_n\} + [B_{b\Delta\beta}] \{\Delta\beta_n\} \quad (3) $$

### 2. Déformations naturelles de CT. Utilisation de la méthode des déformations de substitution (méthode ANS)

Nous définissons les déformations cartésiennes de CT {$\gamma_0$} en fonction des déformations isoparametriques {$\gamma_\xi$} que nous interpolons linéairement, deux à deux, sur les bords élémentaires (Fig. 2).
<figure><img src="image_3.png" /><figcaption>FIG. 2: Déformations de CT de bords</figcaption></figure>

$$ \{\gamma_0\} = \left\{ \begin{array}{c} \gamma_{xz} \\ \gamma_{yz} \end{array} \right\} = \left[ \begin{array}{cc} j_{11} & j_{12} \\ j_{21} & j_{22} \end{array} \right] \left\{ \begin{array}{c} \gamma_{\xi z} = \frac{1}{2}(1-\eta)\gamma_{\xi} \zeta_5 + \frac{1}{2}(1+\eta)\gamma_{\xi} \zeta_7 \\ \gamma_{\eta z} = \frac{1}{2}(1-\xi)\gamma_{\eta} \zeta_8 + \frac{1}{2}(1+\xi)\gamma_{\eta} \zeta_6 \end{array} \right\} \quad (4) $$

$\gamma_{\xi\zeta 5}, \gamma_{\eta\zeta 6}, \gamma_{\xi\zeta 7}, \gamma_{\eta\zeta 8}$ sont les déformations de CT isoparamétriques définies sur les quatre côtés (5, 6, 7, 8) respectivement. Elles sont reliées aux déformations tangentielles de bord $\gamma_{s5}, \gamma_{s6}, \gamma_{s7}, \gamma_{s8}$ par :

$$ \gamma_{\xi\zeta 5} = +\frac{L_5}{2}\gamma_{s5}; \quad \gamma_{\eta\zeta 6} = +\frac{L_6}{2}\gamma_{s6}; \quad \gamma_{\xi\zeta 7} = -\frac{L_7}{2}\gamma_{s7}; \quad \gamma_{\eta\zeta 8} = -\frac{L_8}{2}\gamma_{s8} \quad (5) $$

Les déformations naturelles ou isoparamétriques de bord $\{\gamma_{sk}\}$ sont projetées sur les degrés de liberté élémentaires $\{U_n\}$ et $\{\Delta\beta_{sk}\}$ sous une forme discrète, définie par une intégrale de contour le long de chaque côté $k$ :

$$ \int_{0}^{L_k} (\gamma_{sk} - \tilde{\gamma}_{sk}) ds = 0 \quad (\text{Hyp. cinématique}) \quad (6) $$

$$ \gamma_{sk} = (w_s + \beta_s)_k \quad (7) $$

$$ \gamma_{sk} L_k = w_j - w_i + \frac{L_k}{2}(\beta_{si} + \beta_{sj}) + \frac{2}{3}L_k \Delta\beta_{sk} \quad (8) $$

L'application de l'hypothèse cinématique (8) sur le côté $k$ conduit à un élément à 8 nœuds avec des degrés de liberté additionnels aux nœuds milieux des côtés. Celui-ci possède encore les rotations $\{\Delta\beta_{sk}\}$ comme ddl dont l'utilisation est peu pratique car elle nécessite des maillages éléments finis adaptés. Ce problème a été signalé dans les travaux de Zienkiewicz et ses collaborateurs [7, 3] ainsi que ceux de Ayad et al.[1]. Nous avons opté pour une méthode heuristique qui consiste à utiliser une deuxième hypothèse modifiée de Mindlin, dans le but d'extraire ces mêmes rotations $\{\Delta\beta_{sk}\}$ et de les éliminer en les substituant dans l'expression des déformations de bords $\{\gamma_{sk}\}$. Nous avons introduit une hypothèse mécanique modifiée de Mindlin qui permet de prendre en compte le caractère multicouche des plaques composites. Elle s'écrit comme suit :

$$ \int_{0}^{L_k} \left( \gamma_{sk} - \langle H_{ck11}^{inv} \quad H_{ck12}^{inv} \rangle \left\{ \begin{array}{c} T_s \\ T_n \end{array} \right\} \right) ds = 0 \quad (\text{hyp. mécanique}) \quad (9) $$

$H_{ck11}^{inv}$ et $H_{ck12}^{inv}$ sont deux termes de comportement en CT de la matrice inverse $[H_{ck}]^{-1}$, avec
$$[H_{ck}] = [P_k] [H_c] [P_k]^{-1} ; \quad ([H_c] : \text{la matrice de comportement de CT})$$

$$[P_k] = [P_k]^{-1} = \begin{bmatrix} \cos\theta_k & \sin\theta_k \\ \sin\theta_k & -\cos\theta_k \end{bmatrix}$$

Nous utilisons les relations d'équilibre sur la surface moyenne pour détérminer les efforts trans-chants {T} en fonction des moments fléchissants $\langle M \rangle = \langle M_s M_n M_{sn} \rangle$. Elles s'écrivent de la façon suivante :

$$T_s = M_{s,s} + M_{sn,n}; T_n = M_{sn,s} + M_{n,n} \quad (10)$$

Nous définissons les déformations {$γ_{sk}$} en fonction des dérivées secondes de $β_s$ et $β_n$ par l'intermédiaire des équations d'équilibre (10), et de deux lois de comportement, en flexion et en cisaillement. Nous obtenons :

$$\gamma_{sk} = (H_{ck11}^{inv} H_{f11} + H_{ck12}^{inv} H_{f13}) \beta_{s,ss} ; \quad \beta_{s,ss} = -\frac{8}{L_k^2} \Delta\beta_{sk} \quad (11)$$

La combinaison des deux équations (8 et 11) donne :

$$\Delta\beta_{sk} = \frac{3}{2 + 3\phi_k^{ml}} \left( \frac{w_i - w_j}{L_k} - \frac{\beta_{si} + \beta_{sj}}{2} \right) \quad (12)$$

$\phi_k^{ml} = 8 \left( H_{f11} \frac{H_{ck11}^{inv}}{L_k^2} + H_{f13} \frac{H_{ck12}^{inv}}{L_k^2} \right)$ est le facteur d'influence du CT. L'équation (12) appliquée aux différents côtés élémentaires k conduit aux résultats suivants :

$$\{\Delta\beta_n\} = [A_n^{ml}] \{U_n\} ; \quad [A_n^{ml}] = [A_\phi^{ml}]^{-1} [A_G] \quad (13)$$

La matrice $[A_\phi^{ml}]$ est de dimension (4x4) et la matrice $[A_G]$ est de dimension (4x20).

# 4 Validation numérique

## Plaque composite à 9 couches simplement supportée sous chargement doublement sinusoidal
$q = q_0 \sin(\pi x/L) \sin(\pi y/L) \text{ MPa}$

Ce problème a été proposé et traité par Pagano et Hatfield [4]. Le matériau de base est un composite unidirectionnel fortement orthotrope. Nous avons étudié une stratification à 9 couches 0°/90°/0°/90°/0°/90°/0°/90°/0°. Les épaisseurs totales de couches à 0° et 90° sont égales, et les couches de même orientation ont toutes la même épaisseur. Pour des raisons de symétrie, seul le quart de la plaque est maillé avec un maillage (6 x 6). Les données physiques du problème sont telles que :

$$(E_1/E_2 = 25 ; v = 0.25 ; G_{12} = G_{13} = 0.5 E_2, G_{23} = 0.25 E_2)$$

Les valeurs des déplacements et des contraintes sont données sous la forme adimensionnelle suivante ($S = L/h$) :

$$(\bar{\sigma}_{xx} \quad \bar{\tau}_{xz}) = \frac{1}{q_0 S^2} (\bar{\sigma}_{xx} S \tau_{xz}) ; \quad \bar{w} = \frac{\pi^4 Q}{12 S^4 h q_0} w ; \quad Q = 4G_{12} + \frac{[E_1 + E_2 (1 + 2v_{23})]}{(1 + v_{12} v_{21})} \quad (14)$$
<figure><img src="image_4.png" /><figcaption>FIG. 3: Plaque carrée simplement supportée sous chargement sinusoïdal 9-couches. Déplacement au centre de la plaque en fonction de L/h</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIG. 4: Plaque carrée simplement supportée sous chargement sinusoïdal 9-couches. Distribution des contraintes à travers l'épaisseur</figcaption></figure>

La figure 3 montre la convergence du déplacement $\bar{w}_C$, vers la solution “*Classical Plate Theory*” [6] dans le cas des plaques minces. Sur la figure 4, nous présentons la distribution à travers l'épaisseur des contraintes planes $\bar{\sigma}_C$ pour L/h=10 et de CT $\bar{\tau}_{xz}$ pour trois élancements différents de L/h. La figure 4-b montre une bonne corrélation des contraintes $\bar{\tau}_{xz}$ pour L/h=10 avec les résultats de Pagano.

## 5 Conclusion :

En conclusion, le modèle DMQPml de plaque multicouche donne globalement des résultats satisfaisants sur la flèche et sur les contraintes. Les travaux actuels consistent, le développement du modèle DMPQml avec une théorie d'ordre supérieure.
Références

[1] R. Ayad, A. Rigolot, and N. Talbi. An improved three-node hybride-mixed element for mindlin-reissner plates. *Int. Jou. Num. Meth. Eng.*, 51:919–942, 2001.

[2] I. Katili. A new discrete kirchhoff-mindlin element based on mindlin-reissner plate theory and assumed shear strain fields - part ii : An extended dkq element for thick-plate bending analysis. *Int. Jou. Num. Meth. Eng.*, 36:1885–1908, 1993.

[3] E. Onate, O.C. Zienkiewicz, B. Suarez, and R.L. Taylor. A general methodology for deriving shear constrained reissner-mindlin plate element. *Int. Jou. Num. Meth. Eng.*, 33:345–367, 1992.

[4] N.J. Pagano and S.J Hatfield. Elastic behaviour of multilayered bidirectional composites. *Ame. Ins. Aer. Ast.*, 10, n °7:931–933, 1972.

[5] S. Sakami. *Modélisation numérique des matériaux composites multicouches à laide d' une approche discrète au sens de Mindlin. Le modèle DDM (Displacement Discrete Mindlin).* Thèse de doctorat, URCA, 3 Octobre 2008.

[6] Y. Stavsky. Bending and stretching of laminated aeolotropic plates. *Jou. Eng. Mech. Division ACSE, 87*:31–56, Dec. 1961.

[7] O.C. Zienkiewicz, R.L. Taylor, P. Papadopoulos, and E. Onate. Plate bending elements with discrete constraints : new triangular elements. *Com. Str.*, 35:505–522, 1990.