DEPARTEMENT GENIE MECANIQUE

TP: MODELISATION ET CALCUL PAR LA METHODE DES

ELEMENTS FINIS DE STRUCTURES FORMEES DE PLAQUES

M1
H6 Génie
des matériaux

_AZZOUG HOURIA
_BENBEKKA MASSINISSA

Année universitaire :
2023/2024

1
 Liste des figures
Figure 1 Efforts intérieurs dans une plaque. ................................................................................................. 9

Figure 2(Description d'une plaque géométrie et variables) ......................................................................... 11

Figure 3(varaibles cinématiques utilisées pour l'étudedes plaques en flexion) ............................................ 13

Figure 4(Déformation d'une plaque mince) ................................................................................................ 16

Figure 5(cas de l'état plan de déformation) ................................................................................................... 9

Figure 6(cas de l’état plan de contrainte) ................................................................................................... 12

Figure 7.2(Elément triangulaire simple d’une structure) ............................................................................. 35

Figure 8(: Organigramme pour l’étude de la flexion des plaques) .............................................................. 46

2
 sommair

2 Introduction ..................................................................................................................... 5

3 Différents types de plaques ................................................................................................ 6

4 Différents types de comportements des plaques .................................................................... 6

4.1 Etat membranaire et état flexionnel ......................................................................................... 8
4.2 Hypothèses pour l’étude des plaques ....................................................................................... 9
4.2.1 Quelques approches d'analyse des plaques ............................................................................. 11
4.2.2 L'approche multicouche ........................................................................................................ 13
4.2.3 Approche monocouche équivalente........................................................................................ 14

5 Plaques minces et épaisses et critère de classification ........................................................... 14

5.2 Théorie des plaques minces ................................................................................................... 15
5.2.1 Hypothèses pour l’étude des plaques minces ......................................................................... 15
5.2.2 Champs de déplacement 𝑤𝛼.................................................................................................. 16

6 Introduction à la théorie d’élasticité .................................................................................... 7

6.1 Définition des Forces .............................................................................................................. 7
6.2 Etude de la théorie d’élasticité plane ....................................................................................... 7
6.2.1 Etat plan de déformation ......................................................................................................... 7
6.2.2 Matrice des constantes élastiques (matrice contrainte-déformation) ....................................... 10
6.2.3 Etat plan de contraintes ......................................................................................................... 12
6.2.4 Matrice des constantes élastiques .......................................................................................... 13
6.2.5 L‘énergie de déformation de l’état plan de contrainte ............................................................. 14

7 Principe de la méthode des éléments finis ........................................................................... 28

7.1 Caractéristique d'un élément fini ........................................................................................... 29
7.2 Formulation générale de la méthode des éléments finis ......................................................... 30
7.2.1 Etude d’un élément, relation fondamentale ............................................................................ 30
7.2.2 Calcul des déplacements des nœuds....................................................................................... 32
7.3 Structures élastiques plans ............................................................................................... 33
7.3.1 Elément triangulaire où quadrilatère simple pour l’analyse des plaques ................. 33
7.3.2 Déformation.......................................................................................................................... 36
7.3.3 Contrainte ............................................................................................................................. 36
7.3.4 Matrice de rigidité de l'élément ............................................................................................. 36
7.3.5 Forces nodales équilibrant une densité de force extérieure ..................................................... 37
7.3.6 Effet d'une variation linéaire.................................................................................................. 37

3
 7.4 Plaques minces ..................................................................................................................... 37
7.4.1 Généralités ............................................................................................................................ 37
7.4.2 Utilisation d'éléments triangulaires et quadrilatères ................................................................ 39
7.4.3 Étude de l'élément 𝐴𝑖 𝐴j 𝐴𝑘 simplement appuyé aux nœuds................................................. 39
7.4.4 Matrice de rigidité généralisée de l'élément 𝐴𝑖 𝐴j 𝐴𝑘 ............................................................. 42

8 APLCATIONS NUMERIQUES .......................................................................................... 44

9 Conclusion ...................................................................................................................... 47

4
1 Introduction
La méthode des éléments finis (MEF) est une approche puissante et largement utilisée pour analyser le
comportement des structures formées de plaques. Les structures formées de plaques, telles que les coques
minces, les panneaux de revêtement et les coques composites, sont couramment utilisées dans des
applications telles que l'aérospatiale, l'automobile, la construction navale et la construction civile en raison
de leur légèreté et de leur efficacité structurelle.
La MEF consiste à diviser la structure en un grand nombre d'éléments finis interconnectés qui représentent
la géométrie de la structure. Pour les structures formées de plaques, les éléments finis sont généralement des
éléments à quatre côtés (quadrilatères) ou à trois côtés (triangles) qui sont disposés sur la surface de la plaque.
Chaque élément fini a des propriétés mécaniques qui lui sont propres, telles que la rigidité et l'élasticité, qui
sont utilisées pour décrire le comportement de la structure dans son ensemble.
Une fois la structure discrétisée en éléments finis, les équations d'équilibre et les conditions aux limites sont
appliquées à chaque élément pour former un système d'équations linéaires. Ce système est ensuite résolu
numériquement pour obtenir les déformations, les contraintes et autres caractéristiques de la structure sous
les charges spécifiées. La MEF permet de modéliser des structures complexes et de prendre en compte
différents types de charges et de conditions aux limites, ce qui en fait un outil précieux pour la conception et
l'analyse des structures formées de plaques.
En résumé, la méthode des éléments finis est une approche puissante pour l'analyse des structures formées de
plaques, offrant aux ingénieurs et aux concepteurs la possibilité de concevoir des structures plus efficaces et
plus sûres dans une variété d'applications. Son utilisation est essentielle pour comprendre le comportement
des structures sous des charges diverses et pour optimiser leur conception en fonction des exigences
spécifiques de chaque application
L'analyse numérique des plaques est se faire par différences finis, éléments finis ou éléments de frontières.
Les avantages et potentialité de la méthode des éléments sont tels que les autres méthodes sont peu utilisées,
sauf dans certains cas particuliers. La méthode des éléments finis est extrêmement puissante puisqu’elle
permet d’étudier des structures continues ayant des propriétés géométriques et des conditions de charges
appliquées. Donc, elle est bien adaptée à l'analyse des structures tridimensionnelles à parois minces,
puisqu'elle seule aujourd'hui est capable de résoudre ce type de problème et toute généralité. Elle est donc
auxiliaire aux plus précieux pour l'ingénieur.

5
Les plaques sont des structures très utilisées dans l'industrie sous marine, aérospatial, le génie civil et dans
des constructions courantes (ponts, toits de bâtiments,...), dans le domaine de l'énergie, et dans la
conception industrielle (turbines, pièces de mécanique, carrosserie de voiture,...), et même dans le monde
du vivant (artères, bronches,...), etc. Tous ces domaines sont stratégiques et économiquement très
importants. C'est pour cette raison que les plaques ont fait l'objet de très grand nombre de travaux depuis
plus d'un siècle.
Une plaque est un solide bidimensionnel, défini par une surface de référence plane (plan 𝑥 𝑦) et par une
petite épaisseur (notée h (x , y)) par rapport aux autres dimensions (longueur et largeur, figure II.4), elle
peut être constituée d'un matériau homogène ou peut être obtenue par l'empilement de différentes
couches de matériaux orthotropes. Une analyse élastique des plaques consiste à une analyse d'un modèle
mathématique bidimensionnel ce qui nous conduit à une analyse d'un problème en élasticité plane. [22]

2 Différents types de plaques

En fonction de la nature des matériaux qui les constituent et de la géométrie de leur section
transversale, les plaques peuvent être classées en trois catégories;
1) Les plaques isotropes: elles sont constituées d'un matériau isotrope (acier, béton) et leur
section transversale est homogène. Elles sont définies par deux paramètres élastiques
(𝐸 𝑒𝑡 : respectivement le module d'élasticité et le coefficient de Poisson). On les retrouve
dans les constructions civiles courantes (bâtiments, ouvrages d'art,...).
2) Les plaques orthotropes: leurs propriétés élastiques sont différentes dans deux directions
perpendiculaires. L'orthotropie peut être naturelle (bois) ou techniques (dalles rédies). Le
comportement de ces dalles est défini par quatre paramètres élastiques et on les retrouve dans
les constructions navales, aerononavales, de réservoirs de l'industrie chimique, des bâtiments
et d'ouvrages d'art. [23]

3) Les plaques anisotropes: leurs propriétés élastiques sont différentes dans toutes les directions.
Neuf paramètres élastiques sont suffisants pour les devenir. Elles sont souvent constituées de
matériaux composites et sont surtout utilisées dans l'industrie aéronavale.

3 Différents types de
comportements des plaques
L'énergie de déformation d'une plaque peut être décomposée en trois termes principaux : la

6
flexion, la membrane et le cisaillement. Lorsqu'on fait tendre l'épaisseur vers zéro, l'énergie de
cisaillement devienne négligeable et la déformation subie par la plaque appartient à 1'une des trois
catégories asymptotiques suivantes en fonction de la géométrie, des conditions aux limites et des
forces appliquées :

7
 1) La flexion dominante;

2) La membrane dominante;

3) Le cas mixte;

Dans les deux premières catégories, c'est la partie correspondante de l'énergie de déformation qui
est dominante, tandis que dans les cas mixtes aucune partie n'est négligeable par rapport à l'autre.
En considérant la même force, matériau et conditions aux limites, une plaque (plaque) peut exhiber
des comportements asymptotiques complètement différents en fonction de la nature géométrique
de la surface moyenne. [24] [25]

L'espace continu des déplacements inextensionnels « avec énergies de membrane et cisaillement

nul, est nommé 𝑉0 » qui joue un rôle fondamental dans le classement évoqué ci-dessus : dans le
cas flexion-dominante, le problème limite est formulé dans cet espace qui doit être différent de
zéro, tandis que dans le cas membrane-dominante, cet espace est réduit au champ nul.

3.1 Etat membranaire et état

flexionnel
L'état de contrainte par le quel l'élément structural résiste aux actions extérieures, est caractérisé
par les efforts intérieurs définis au niveau de la surface moyenne. La trace de cette surface dans
une section droite s'appelle la ligne moyenne. Les efforts intérieurs sont décrits par unité de
longueur de ligne moyenne dans les sections droites.

L'état membranaire s'associe aux efforts intérieurs de type force agissant dans la surface
moyenne, à savoir les efforts normaux et les efforts tangentiels.

L'état flexionnel regroupe les efforts intérieurs de caractère flexionnel, soit les moments de
flexion, les moments de torsion et les efforts tranchants.

 Selon le mode du comportement, on peut distinguer quatre types d'éléments structuraux :

L'élément de paroi est défini par la géométrie plane de surface moyenne (plan moyen) et par
son épaisseur : sollicitée par des charges agissant dans son plan moyen, il résiste par un état
membranaire (figure. II. 3a); les efforts normaux et tangentiels résultent d'ailleurs de l'état plan
de contrainte.
 L'élément de plaque est défini par la géométrie plane de sa surface moyenne (plan ou feuillet
moyen); il résiste aux charges agissantes normalement à son plan moyen par un état flexionnel
(figure. II.3b).
 L'élément de plaque-membrane est la superposition des deux cas précédents et réunit donc
l'état membranaire de paroi et l'état flexionnel de plaque (figure. II.3c); bien que plan, il se
comporte de manière spatiale, pouvant être soumis à des charges quelconques, tant parallèles
que perpendiculaires à son plan moyen; il constitue la base des structures plissées.
 En fin l'élément structural de plaque est, par nature, la courbe spatiale; il utilise les deux états
d'effort, membranaire et flexionnel, pour s'opposer aux actions arbitraires pouvant le solliciter
(figure. II.3d); exceptionnellement grâce à sa courbure, une plaque peut ne

8
 résister aux charges que par l'état membranaire (structure gonflable, textile, peau, etc…).
[26] 𝑧

𝑉𝑦
𝑀𝑦
𝑀𝑦𝑥
𝑁𝑥 𝑉𝑥

𝑥 𝑥 𝑀𝑥
𝑧
𝑀𝑥𝑦

(a) Etat membranaire de paroi (b) Etat flexionnel de plaque

x
(c) Plaque membrane : superposition des états (d) Plaque cinq efforts
membranaire et flexionnel intérieurs par section
droit
Figure 1 Efforts intérieurs dans une plaque.

3.2 Hypothèses pour l’étude des

plaques
La théorie des plaques, qu'on considère dans ce qui suit, est basée sur :

- L'hypothèse dite des sections droites (ou planes) : les points matériels situés sur une normale à la surface
moyenne non déformée restent sur une droite dans la configuration déformée. les déplacements 𝑢 et 𝑣

9
(suivant 𝑥 et 𝑦) d'un point quelconque , 𝑦 ,𝑧 varient alors linéairement en 𝑧 et le déplacement transversal
𝑤 (suivant𝑧) n'est fonction que de 𝑥 et 𝑦 .cette hypothèse a été proposée par plusieurs auteurs, dont les
plus connus sont Reissner, Hencky, Bollé, Mindlin et permet de prendre en compte l'influence des
déformations de CT. Elle a été initialement proposée pour les plaques isotropes mais ensuite a été
appliquée aux plaques orthotropes et multicouches.

10
 - L'hypothèse d'une déformation transversale 𝜀𝑧 nulle (pas de variation d'épaisseur).
- L'hypothèse des contraintes planes : dans les relations de comportement, la contrainte 𝜎𝑧 est
négligeable par rapport aux autres composantes du tenseur de contraintes.
- L'hypothèse d'anisotropie plane pour chaque couche dans le cas d'une plaque composite (l'axe
𝑧 est l'axe d'orthotropie de toutes les couches).
- La prise en compte d'une loi de comportement de cisaillement transversal CT dépendant de
facteurs de correction déterminés a priori par des considérations d'équivalences énergétiques, en
généralisant la démarche utilisée pour les poutres et en s'inspirant de la théorie mixte de Reissner.

Cette théorie des plaques, dite du premier ordre (relative aux effets de CT) est appelée, dans la
suite de ce chapitre, théorie (ou modèle) de Reissner/Mindlin. Elle fait intervenir (figure II.5) :

- Cinq variables cinématiques indépendantes (les déplacements de membrane 𝑢 et 𝑣 dans le plan

de référence 𝑧 = 0, le déplacement transversal 𝑤 et les rotations 𝛽𝑥 et 𝛽𝑦 de la normale à la surface
moyenne dans les plans 𝑥𝑧 et 𝑦𝑧 respectivement) ;
- Trois efforts résultants de membrane notés 𝑁𝑥, 𝑁𝑦, 𝑁𝑥𝑦; trois moments de flexion (ou couples)
notés 𝑀𝑥, 𝑀𝑦, 𝑀𝑥𝑦 ; deux efforts tranchants notés 𝑇𝑥, et .
z

h
x, y
P
h Ly
y

Lx

z, w
Plan de référence

u, v
Nx, Ny, Nxy
P
y, v W, βx, βy
Mx, My, Mxy
Tx, Ty
x, u
Figure 2(Description d'une plaque géométrie et variables)

3.2.1 Quelques approches d'analyse des plaques

La théorie de la flexion des plaques la plus simple est la théorie de Kirchhoff-Love, appelée aussi
théorie classique de la flexion (Classical Plate Theory). Elle est basée sur le champ de
déplacements suivant:

11
 Dans lequel (𝑢, 𝑣, 𝑤) sont les composantes du déplacement, dans un repère de coordonnées
cartésiennes, d'un point appartenant au feuillet moyen (plan de la plaque de coordonnée z = 0).
Le champ de déplacement, composé des équations (II.37), (II.38) et (II.42), repose sur les
hypothèses classiques selon lesquelles une section droite, normale au feuillet moyen avant
déformation, reste droite et normale après déformation (Figure. II.3.a). Sous ces hypothèses, les
déformations à l'effort tranchant et les contraintes normales au feuillet moyen sont négligées : l'état
de déformation est généré uniquement par la flexion et les extensions du feuillet moyen (tensions
de membrane).
Les bases d'une théorie prenant en compte les déformations au cisaillement ont été posées par
Reissner (1945). Puis la théorie complète a été publiée par Mindlin (1951). Cette théorie est connue
sous le nom de théorie améliorée (Improved Plate Theory -IPT).

Elle est basée sur le champ de déplacement suivant:

𝑢 =
−𝑧Ø(
𝑥, 𝑦) (II.38)
{ = −𝑧Ø𝑥(𝑥, 𝑦)
𝑤 = 𝑤0(𝑥, 𝑦)

Où les termes Øy et Øx désignent les rotations des sections droites autour des axes y et x,
respectivement. Selon cette théorie, une section droite, normale au feuillet moyen avant
déformation, reste droite mais n'est plus normale au feuillet moyen après déformation (Figure.
II.5). De plus, il est possible de prendre en compte, en plus des déformations à l'effort tranchant,
l'inertie de rotation des sections droites.
Un modèle de plaque plus complexe peut être développé en utilisant comme composantes du
champ de déplacement des fonctions de la variable z d'ordre deux (Reddy (1984)). Ce modèle
introduit des nouvelles inconnues dont l'interprétation physique n'est pas évidente. De plus, l'apport
de ce type de modèle sur les résultats n'est pas justifié en regard de la complexité introduite dans
les calculs (Reddy (1984)).

12
 z

M
x z
h
y Feuillet (ligne) moyen

Configuration non déformée

x
y

𝚿𝑦
Configuration déformée – théorie améliorée

Figure 3(varaibles cinématiques utilisées pour l'étudedes plaques en flexion)

3.2.2 L'approche multicouche

Une structure composite multicouche peut être considérée comme un corps hétérogène constitué
d'un nombre fini de couches homogènes anisotropes collées. La caractéristique géométrique d'une
plaque est une épaisseur faible par rapport aux autres dimensions. La modélisation des structures
multicouches modernes avec une forte anisotropie (par exemple : faible rapport du module de
cisaillement transverse de l'âme par rapport au module d'élasticité longitudinal des peaux dans le
cas des structures sandwich) exige des théories raffinées qui prennent en compte une bonne
description des cisaillements transverses. On trouve dans (Carrera) des revues complètes sur les
différents modèles existants de type élasticité tridimensionnelle ou de type plaque.

L'intérêt d'une approche tridimensionnelle réside dans l'obtention de résultats satisfaisants

tridimensionnels, utilisés notamment comme une référence. L'adoption d'une approche
tridimensionnelle ne présente toutefois d'utilité que dans la mesure où les équations différentielles
finalement obtenues peuvent être résolues. L'approche tridimensionnelle (3D) est donc limitée à
certains cas de géométrie, empilement et chargement simple. De même, la prise en compte des
endommagements spécifiques aux stratifies (délaminage, fissure transverse,. )
exige une bonne description des champs au voisinage des interfaces, sites privilégies de ces
endommagements. La encore une approche 3D peut fournir des informations dans certains cas

13
 (bord libre droit, plaque trouée...) mais il est préférable d'envisager des démarches
bidimensionnelles (2D) de type plaque plus opératoire. Dans la suite ce sont ces approches 2D
qui seront mentionnées.
Durant ces dernières années, plusieurs approches bidimensionnelles ont été développées pour la
modélisation des structures multicouches tenant compte les cisaillements transverses ou des
endommagements. Ils peuvent être regroupés en fonction du type d'approche adopté :

 L'approche monocouche équivalente

 L'approche par couche
 L'approche développement asymptotique
3.2.3 Approche monocouche équivalente
Dans l'approche monocouche équivalente, le nombre d'équations ne dépend pas du nombre de
couche, la plaque multicouche est homogénéisée et est donc considérée comme une seule couche.

Depuis le premier travail de Sophie Germain en 1815 sur les plaques minces en passant par les
modèles du premier ordre de Love - Kirchhoff et de Reissner-Mindlin, de nombreux chercheurs
ont d'enveloppé des théories de plaques à partir de cinématiques ou champs de contraintes
plus raffinés.

Nous passons en revue, dans ce qui suit, les principaux modèles.

4 Plaques minces et épaisses et

critère de classification
- Le plan 𝑥𝑦 noté A est généralement le plan moyen de la plaque :
ℎ
−𝑡 ≤ 𝑧 ≤ 𝑡 Avec 𝑡 =
2

C'est également le plan neutre si les propriétés matérielles sont symétriques par rapport au
plan 𝑥𝑦.
- Pour les plaques homogènes isotropes, la validité de la théorie de plaque retenue dépend des
caractéristiques géométriques. On admet généralement les hypothèses de Mindlin si
4 ≤ 𝐿/ℎ ≤ 20 et celles de Kirchhoff si 𝐿/ℎ > 20 où 𝐿 est une dimension caractéristique dans
le plan 𝑥𝑦.

- Pour les plaques composites le rôle des déformations de CT dépend non seulement des
caractéristique géométriques (l'élancement 𝐿/ℎ ), mais également des caractéristiques mécaniques
représentées par le rapport 𝐸/𝐾𝐺 (où 𝐸 est un module caractéristique intervenant dans la flexion,
𝐺 un module de CT et 𝑘 un facteur de correction de CT). comme pour les poutres on peut utiliser
le coefficient Ø pour évaluer l'influence du CT :
Ø = (ℎ/𝑙)2(𝐸/𝐾𝐺)
)

Par exemple, pour une plaque carrée homogène isotrope avec 𝐸/𝐾𝐺 = 2 :

14
 𝐿/ℎ = 4 Ø = 0.125 ; 𝐿/ℎ = 20 Ø = 0.005; 𝐿/ℎ = 100 Ø = 0.002 (II.40)

Pour une plaque composite avec 𝐸/𝐾𝐺 = 50 :

L/h = 4 Ø = 3.125 ; L/h = 20 Ø = 0.125; L/h = 100 Ø = 0.005 [Link].[Link].1 (II.
41)
On observe que Ø est identique pour une plaque isotrope avec 𝐿/ℎ = 4
et pour une plaque
composite (𝐸/𝐾𝐺 = 50) avec 𝐿/ℎ = 20 (l'influence du TC sera sensiblement la même pour les
deux cas).
- La théorie de Kirchhoff peut être interprétée comme un cas particulier de la théorie de Reissner
/ Mindlin. (Ainsi, un ″bon″ modèle élément fini basé sur la théorie de Reissner/Mindlin devra
donner des résultats en accord avec la théorie de Kirchhoff si l'influence du CT est faible).
- La détermination de la loi de comportement en cisaillement (relation entre efforts tranchants et
cisaillements moyens) peut se faire expérimentalement ou par calculs à partir de différentes
hypothèses à caractère statique ou dynamique.
- Pour les plaques composites, la précision des résultats de la théorie du premier ordre en CT (par
rapport aux solutions tridimensionnelles) dépend beaucoup de la loi de comportement modifiée
par les facteurs correction de CT. Pour les plaques homogènes, le facteur de correction de CT est
généralement admis 𝑘 = 5/6 (comme pour les poutres rectangulaires).
- Les théories de plaques dites d'ordre supérieur sont basées sur des distributions plus générales
des champs de déplacements (et éventuellement des champs de contraintes , 𝜎𝑦𝑧, 𝜎𝑧) : les
déplacements 𝑢, 𝑣, 𝑤 sont représentés par les termes 1, 𝑧, 𝑧2, 𝑧3 ou trigonométriques. Ces théories
d'ordre supérieur permettent d'éviter l'introduction des facteurs de correction de CT mais font
généralement intervenir plus de cinq variables cinématiques .signalant qu'une des premières
contributions à ces théories fut celle de Lévy en 1877 (ainsi que le soulignent Bollé et Reissner).
[22]

4.2 Théorie des plaques minces

Les modèles sont basés sur une distribution linéaire des déplacements dans l'épaisseur. L'hypothèse
adoptée est celle de Love-Kirchhoff [27] de contraintes planes, les déformations dues aux
cisaillements transverses sont négligées. La normale reste droite et est perpendiculaire à la surface
moyenne après avoir déformée (Figure. II.6 et II.7)
4.2.1 Hypothèses pour l’étude des plaques minces
Les hypothèses cinématiques adoptées pour les plaques minces par Kirchhoff de contrainte plane,
généralisent à deux dimensions celles adoptées pour les poutres sans déformation à l'effort
tranchant.

1) La plaque est d'épaisseur petite devant les autres dimensions. Elle possède un plan moyen
aussi appelé plan neutre.
2) L'hypothèse des sections droites inclut celle dite de conservation des normales (les points
matériels situés sur une normale à la surface moyenne non déformée restent sur la normale à

15
 la surface moyenne déformée) ce qui conduit au modèle classique généralement associé au
nom de Kirchhoff (mais qui doit également beaucoup aux scientifiques français du 19 e siècle
: Sophie Germain, Lagrange, Poisson, Boussinesq, Levy, Saint-Venant).
3) La déformation en cisaillement transverse est donc négligée.
4) Les termes non linéaires du déplacement sont négligés, en particulier, l'inertie de rotation est
négligée. Seul le déplacement transversal 𝑤 est considéré.
5) La contrainte 𝜎𝑧 dans la direction transversale est nulle. Elle doit en effet s'annuler sur les
faces extérieures et, du fait que la plaque est mince, il est naturel d'admettre qu'elle est nulle
en tout z.
6) Les éléments basés sur la théorie sont les mieux adaptés à la modélisation des plaques minces.
Leur formulation ne nécessite que l'approximation d'une seule composante de déplacement
(𝑥, 𝑦). Les dérivées secondes des déplacements interviennent dans l'expression de l'énergie de
déformation. Ainsi, la condition de compatibilité inter-éléments correspond à la continuité 𝐶1
du déplacement traverse c.-à-d. à la continuité de 𝑤 et de la
pente normale 𝑑w.
𝑑𝑥

Figure 4(Déformation d'une plaque mince)

4.2.2 Champs de déplacement 𝑤𝛼

𝑤0

𝑢0

(Figure. 5) : Cinématique de love - Kirchhoff.

16
 Les composantes des champs de déplacements pour cette théorie des plaques minces s'expriment comme
suit :

𝑢𝛼 (𝑥1, 𝑥2, 𝑥3 = 𝑧) = 𝑢0(𝑥1, 𝑥2) − 𝑧. (𝑥1, 𝑥2) 𝛼 = 1, 2

𝛼
{
𝑢3 (𝑥1, 𝑥2, 𝑥3 = 𝑧) = (𝑥1, 𝑥2)

Avec :

𝑢𝛼∶ 𝑙𝑒 𝑑é𝑝𝑙𝑎𝑐𝑒𝑚𝑒𝑛𝑡 𝑑𝑒 𝑚𝑒𝑚𝑏𝑟𝑎𝑛𝑒 𝑑𝑎𝑛𝑠 𝑙𝑎 𝑑𝑖𝑟é𝑐𝑡𝑖𝑜𝑛

𝑤∶ 𝑙𝑎 𝑓𝑙è𝑐ℎ𝑒 𝑑𝑒 𝑙𝑎 𝑝𝑙𝑎𝑞𝑢𝑒

𝑤𝛼∶ 𝑙𝑎 𝑟𝑜𝑡𝑎𝑡𝑖𝑜𝑛 𝑑𝑢𝑒 à 𝑙𝑎 𝑓𝑙𝑒𝑥𝑖𝑜𝑛 (𝑠𝑎𝑛𝑠 𝑐𝑖𝑠𝑎𝑖𝑙𝑙𝑒𝑚𝑒𝑛𝑡)

D'autre part, les rotations sont données par :

Dans ce cas, le champ de déplacement est défini uniquement par une variable, c'est le déplacement
transverse 𝑤.

17
5 Introduction à la théorie
d’élasticité
Un comportement élastique signifie que l’application des contraintes mène à des déformations
réversibles quand on retire les contraintes, L’élasticité linéaire implique que la relation entre les
contraintes et les déformations obéit à une équation linéaire, cette notion est valable seulement
dans le cas des petites déformations.

Dans ce chapitre une brève exposition des équations de base de la théorie de l’élasticité et des
relations entre contraintes et déformations dans le cas tridimensionnel, sera présentée
ultérieurement. De ce cas général, dérive le cas particulier de l’élasticité plane (État plan de
déformation et état plan de contrainte). Ce chapitre est consacré aussi à une présentation détaillée
de la théorie des plaques.

5.1 Définition des Forces

Il y a deux sortes de forces qui agissent sur un solide, les forces de volume et les forces de surface.
Les forces de volume sont celles qui agissent dans le volume (forces de gravité,
électrostatique,…etc.). Les forces surfaciques sont celles qui agissent sur la surface extérieure ou
celles de séparation des deux parties d’un seul corps et il résulte de ces forces des déplacements et
des déformations. [20]

Equations d’équilibre

Tous les problèmes de calcul de structures élastiques se ramènent à un système d'équations aux
dérivées partielles, qu'il est souvent impossible de résoudre analytiquement.
Les équations d’équilibre pour un corps élastique déformable par rapport aux axes sont résumées
comme suit :

5.2 Etude de la théorie d’élasticité

plane
5.2.1 Etat plan de déformation
On peut dire qu’on est en présence d’un problème de déformation plane lorsqu’on a un corps dont
l’une des dimensions est très importante par rapport aux deux autres, et chargé par une force
perpendiculaire à la direction de la longueur (Figure. II.1). [21]
Un tel problème est défini par :
𝑢 = 𝑢 (𝑥, 𝑦)

7
{𝑣 = 𝑣 (𝑥, 𝑦)
𝑤=0 [Link].[Link].1 (
I
I
.
2
)

8
 x
z

Figure 5(cas de l'état plan de déformation)

Les relations entre les déformations et les déplacements sont :

Dans le cas où la déformation est plane, les contraintes agissant sur un tel corps sont obtenues à
partir des équations :
𝜎𝑥 = 2𝐺. 𝜀𝑥 + 𝜆(𝜀𝑥 + 𝜀𝑦)
𝜎𝑦 = 2𝐺. 𝜀𝑦 + 𝜆(𝜀 𝑥 + 𝜀𝑦)
𝜎𝑧 = 𝜆(𝜀𝑥 + 𝜀𝑦) (II.4)
𝜏𝑥𝑦 = 𝐺. 𝛾𝑥𝑦
Inversement 𝜀𝑥, 𝜀𝑦, 𝛾𝑥𝑦 s’exprime en fonction de 𝜎𝑥, 𝜎𝑦, 𝜏𝑥𝑦 par les relations :

𝑍 𝑒𝑡 𝑍 Étant nécessairement nuls, les équations d’équilibres deviennent :

9
 Et les conditions aux limites deviennent :

Les conditions d’intégrations du tenseur de déformation se réduisent à la seule condition :

La MEF est toujours associée à l’analyse approximative de la forme faible des équations
différentielles, on peut déterminer u qui satisfait :

∫Ω𝜎𝑖j𝜀𝑖j(𝛿𝑢)𝑑Ω = ∫𝛛Ω𝛿𝑢. 𝑡𝑑𝑠 + ∫Ω𝛿𝑢. 𝑓𝑑Ω (II.8)

Après intégration dans le premier terme de la partie gauche, on obtient l’équation suivante :
𝜎𝑖j𝜀𝑖j(𝛿𝑢) = 𝜎𝑥𝜀𝑥(𝛿𝑢) + 𝜏𝑥𝑦𝜀𝑥𝑦(𝛿𝑢) + 𝜏𝑦𝑥𝜀𝑦𝑥(𝛿𝑢) + 𝜎𝑦𝜀𝑦(𝛿𝑢) (II.9)
En utilisant la propriété de symétrie de σij(𝑢) et 𝜀𝑖j(𝛿𝑢) par rapport à 𝑖 et j .
𝜎𝑖j𝜀𝑖j(𝛿𝑢) = 𝜎𝑥𝜀 𝑥(𝛿𝑢) + 𝜏𝑥𝑦𝛾𝑥𝑦(𝛿𝑢) + 𝜎𝑦𝜀 𝑦(𝛿𝑢) = {𝜀(𝛿𝑢)}𝑇{𝜎} (II.10)

{(𝛿𝑢)}, {𝜎} Sont définies par :

𝜀𝑥(𝛿𝑢) 𝜎𝑥

{𝜀(𝛿𝑢)} = { (𝛿𝑢) } ; {𝜎} = { 𝜎𝑦 } (II.11)

𝛾𝑥𝑦(𝛿𝑢) 𝜏𝑥𝑦

5.2.2 Matrice des constantes élastiques (matrice contrainte-

déformation)
Basée sur la relation entre σij et 𝜀𝑖j, on lie {(𝑢)} et {(𝑢)} par la matrice [𝐷]

{𝜎(𝑢)} = [𝐷]{𝜀(𝑢)}

[𝐷] Étant la matrice des constantes élastique.

Les composantes de σij dans {(𝑢)} sont données par,

10
11
On peut les mettre sous forme matricielle :

Les intégrants σij et 𝜀𝑖j(𝛿𝑢) dans la partie gauche de la relation (II.8) peuvent être exprimés par :
𝜎𝑖j𝜀𝑖j(𝛿𝑢) = {𝜀(𝛿𝑢)}𝑇[𝐷]{𝜀(𝑢)}
5.2.3 Etat plan de contraintes
Contrairement à la condition de déformation plane, qui stipule que la dimension longitudinale dans
la direction z est importante par rapport aux dimensions suivant x et y ; la condition de contrainte
plane est caractérisée par le fait que la dimension suivant la direction z est très petite par rapport
aux deux autres dimensions figure (II.2).

x z
z

Figure 6(cas de l’état plan de contrainte)

Ainsi les contraintes , 𝜏𝑥𝑧, 𝜏𝑧𝑥, 𝜏𝑧𝑦, 𝜏𝑦𝑧, deviennent nulles. La déformation 𝜀𝑧𝑧 n’est pas nulle par
effet de poisson.
Dans les deux cas de figures, états plans de contraintes et de déformations, on peut tirer profit
des situations qui les caractérisent pour rendre les calculs plus faciles .De la même manière que
pour l’état plan de déformation, on aboutit ainsi à :

𝜎𝑖j𝜀𝑖j(𝛿𝑢) = 𝜎𝑥𝜀 𝑥(𝛿𝑢) + 𝜏𝑥𝑦𝛾𝑥𝑦(𝛿𝑢) + 𝜎𝑦𝜀𝑦(𝛿𝑢) = {𝜀(𝛿𝑢)}𝑇{𝜎} (II.18)

{(𝛿𝑢)}, {𝜎} Sont définies par

𝜀𝑥(𝛿𝑢) 𝜎𝑥

{𝜀(𝛿𝑢)} = { (𝛿𝑢) } ; {𝜎} = { 𝜎𝑦 } (II.19)

𝛾𝑥𝑦(𝛿𝑢) 𝜏𝑥𝑦

Dans le cas où la déformation est plane, les relations entre les déformations et les contraintes
seront :

12
 Inversement 𝜀1, 𝜀2, 𝛾3 s'expriment en fonction de 𝜎1, 𝜎2, 𝜎3 par les relations :

Comme dans le cas de la déformation plane, on a nécessairement 𝑍 = 0; 𝑍 = 0, les équations

d'équilibre et les conditions aux limites se réduisent aux équations (II.6) et (II.7). Les conditions
d'intégrabilité du tenseur d e déformation se réduisent aux quatre conditions:

5.2.4 Matrice des constantes élastiques

Basée sur la relation entre σij et 𝜀𝑖j , on lie {𝜎(𝑢)} et {𝜀(𝑢)} par la matrice [𝐷]

{𝜎(𝑢)} = [𝐷]{𝜀(𝑢)}

Les composantes de σij dans {(𝑢)} sont données par,

Si on conduit en plus des transformations, on peut donner la relation suivante :

13
 La relation (II.22, 23) devienne:

On peut les mettre sous forme matricielle

De la même façon que précédemment, les intégrants 𝜎𝑖j et 𝜀𝑖j(𝛿𝑢) dans la partie gauche de la relation
(II.8) peuvent être exprimés par :
𝜎𝑖j(𝑢)𝜀𝑖j(𝛿𝑢) = {𝜎(𝛿𝑢)}𝑇[𝐷]{𝜀(𝑢)} (II.31)
5.2.5 L‘énergie de déformation de l’état plan de contrainte
L'expression de l'énergie de déformation dans l'état de contrainte tridimensionnel peut être écrite
sous la forme.
[Link].[Link].1
(II.32)
Où 𝜎 est le vecteur pour les six composants de contrainte et 𝜀 contient les six composants
correspondants de déformation.
D'après la loi de Hooke, l'équation d'énergie de déformation peut être exprimée en fonction de contrainte
comme suit :

Dans le cas de l'état plan de contrainte on a : 𝜎𝑧 = 𝜏𝑥𝑧 = 𝜏𝑦𝑧 = 0 , donc l'équation (II.34)
devient :

14
 L'analyse numérique des plaques est se faire par différences finis, éléments finis ou éléments de
frontières. Les avantages et potentialité de la méthode des éléments sont tels que les autres
méthodes sont peu utilisées, sauf dans certains cas particuliers. La méthode des éléments finis est
extrêmement puissante puisqu’elle permet d’étudier des structures continues ayant des propriétés
géométriques et des conditions de charges appliquées. Donc, elle est bienne adaptée à l'analyse des
structures tridimensionnelles à parois minces, puisqu'elle seule aujourd'hui est capable de résoudre
ce type de problème et toute généralité. Elle est donc auxiliaire aux plus précieux pour l'ingénieur.

6 Principe de la méthode des

éléments finis
La méthode des éléments finis est une méthode approchée de calcul numérique, permettant de
déterminer l'équilibre élastique des structures continues à deux ou trois dimensions : structures
planes, solides élastiques, plaques minces, membranes et coques.

Grâce au principe de d'Alembert, elle permet d'étudier les vibrations de ces structures. Elle
peut également être étendue aux structures visco-élastiques ou plastiques.
La méthode des éléments finis consiste à déterminer les déplacements en un certain nombre de
points Av(v = 1, 2, . . . , N) de la structure; ces déplacements sont définis par des matrices
colonnes ∆v qui sont les inconnues qu'il faut calculer. La méthode repose sur deux hypothèses.

La première hypothèse consiste à diviser, au moyen de lignes ou de surfaces imaginaire

passant par les nœuds 𝐴𝑣, la structure en un nombre fini d'éléments (𝜇 = 1, 2, . . . , 𝑀)
et à supposer que les élément 𝐸𝜇ne sont mécaniquement assemblés entre eux qu'aux nœuds
𝐴 𝑣.

La seconde hypothèse consiste à représenter le déplacement en un point 𝑀 de l’élément 𝐸𝜇

par une fonction d'interpolation (𝑀) prenant les valeurs ∆𝛼, ∆𝛽, ∆𝛾 lorsque le point M vient aux
nœuds ∆𝛼, ∆𝛽, ∆𝛾,……. de l'élément.

Le théorème des travaux virtuels permet de déterminer quelles forces il faut appliquer
aux nœuds ∆𝝰, ∆β, ∆ 𝗒 ,…….pour assurer l’équilibre de Eµ connaissant les déplacements ∆𝝰,
∆β, ∆𝗒,……des nœuds, les forces extérieures appliquées à l'élément et les déformations
éventuelles (dilatation thermique, retrait) imposées à l'élément. On obtient ainsi les
matrices de rigidité généralisées des éléments. Il est donc possible de calculer les
déplacements des nœuds.
La première hypothèse conduit à n'écrire que l'équilibre d'ensemble des éléments; localement, les
équations de l'équilibre ne sont pas vérifiées, puisqu'on néglige les réactions mutuelles qui
s'exercent le long de la frontière commune à deux éléments.

28
 La deuxième hypothèse peut donner lieu à des discontinuités du déplacement (ou de ses dérivées)
interpolé à partir des valeurs aux nœuds lorsqu'on traverse la frontière commune à deux éléments.

Il est souhaitable de choisir des fonctions d'interpolation ne conduisant pas à des discontinuités du
déplacement, mais ce n'est pas toujours possible.
Pour que la méthode des éléments finis converge vers la solution exacte lorsque le nombre
des éléments tend vers l'infini, la plus grande dimension de chaque élément tendant vers zéro,
il est nécessaire que les deux critères suivants soient respectés.
- La fonction d'interpolation doit être choisie de façon à n'entraîner aucune déformation d'un
élément non sollicité dont les déplacements des nœuds résultent d'un déplacement d'ensemble de
l’élément.
- La fonction d'interpolation doit être choisie de façon à obtenir un état de déformation uniforme
d'un élément non sollicité lorsque les déplacements des nœuds sont compatibles avec cet état
de déformation uniforme.

Enfin, pour obtenir une bonne approximation, il faux choisir des éléments dont les différentes
dimensions sont sensiblement égales : triangles voisins d'un triangle équilatéral, rectangles voisins
d'un carré, etc.

La méthode des éléments finis peut être considérée comme un cas particulier de la méthode de
Ritz pour la résolution approchée des problèmes d'élasticité. La méthode de Ritz consiste à
représenter le déplacement par une combinaison linéaire de fonctions vérifiant les conditions
aux limites, et à adopter pour solution approchée la combinaison linéaire qui rend minimale
l’énergie potentielle totale.

Ce dernier aspect de la méthode des éléments finis a permis de l'appliquer à de nombreux

problèmes de champ, lorsqu’une méthode variationnelle est possible, par exemple déterminer dans
un domaine donné les valeurs numériques d’une fonction scalaire 𝑈 (𝑀) vérifiant une équation
aux dérivées partielles linéaire du second ordre, et satisfaisant à des conditions aux limites linéaires
sur la frontière du domaine. De tels problèmes, dont le problème de Dirichlet et le problème de
Neumann ne sont que des cas particuliers, se rencontrent souvent. [34]

6.1 Caractéristique d'un élément fini

Attributs d'un élément fini

Un découpage, artificiel (par l’imagination) d’un milieu continu, en éléments finis permet, d’isoler
un élément fini pour l’étudier et en établir les caractéristiques. L’identification d’un élément fini
comprend les points suivants.

- Géométrie : un élément fini peut être segment de droit ou de courbe, triangle ou quadrilatère
(plan ou courbe), tétraèdre, prismes ou hexaèdre (𝑏𝑟𝑖𝑞𝑢𝑒). Les frontières entre éléments peuvent
être respectivement des points, des segments de droit ou de courbes, des faces planes ou courbes.

29
 - Matériau: le matériau d'élément est défini par une loi de comportement (loi de Hooke isotrope,
etc.)
- Nœuds: les nœuds définissent la géométrie et assurent la connexion des éléments les uns aux
autres. Ils occupent les sommets, les extrémités, les milieux des arêtes et faces, etc. on y choisit
et définis les degrés de liberté du problème.
- Degrés de liberté: la fonction d'approximation choisie (en général le champ des déplacements)
est exprime en fonction des valeurs particulières qu'elles prennent aux nœuds communs des
différents éléments adjacents, permet de reconstituer la solution complète (assemblage) tout
veillant à respecter certains règles, dites critère de convergence.
- Forces nodales: à travers les nœuds transitent des forces associées aux degrés de liberté. elle
sont dues aux charges appliquées à l'élément (poids propre, charge uniforme, température..)
Ces paramètre d'identification permettent de construire les deux caractéristique clés d'un élément
fini que sont sa matrice de rigidité et son vecteur force. [35]

6.2 Formulation générale de la

méthode des éléments finis
6.2.1 Etude d’un élément, relation fondamentale
Considérons les éléments 𝐸𝜇 où aux appuis les nœuds sont nommés par : ∆𝛼, ∆𝛽, ∆𝛾,…Un
repère orthonormé étant choisi, nous représentons le déplacement 𝑈 (𝑀) d'un point 𝑀 de l'élément
par la formule matricielle : [36]

Les conditions :

𝑈(𝐴𝛼) = ∆𝛼, 𝑈(𝐴𝛽) = ∆𝛽, …

Montrent que la fonction d'interpolation Ф('𝑀) doit vérifier, quels que soient les nœuds ∆𝛼
et ∆𝛽, les relations :

𝐼 : Désignant la matrice unité.

La déformation (𝑀) se déduit de 𝑈 (𝑀) par un opérateur différentiel linéaire particulier à

chaque type de structure :

Remplaçons (𝑀) par son expression (III.1); en posant :

(𝑀) = {Φ(𝑀)} 𝑒𝑡 𝐵𝜇(𝑀) = {Φ𝜇(𝑀)}

𝛼 𝛼

30
 Nous obtenons la formule qui donne la déformation (𝑀) en fonction des déplacements des
nœuds de l’élément :

Si 𝑟 désigne la matrice d'élasticité du matériau constitutif de la structure, qui peut éventuellement

être fonction de 𝑀, la contrainte au point 𝑀 est définie par la matrice colonne :

(𝑀) = 𝑟𝐷0(𝑀)

Plus généralement, supposons l’élément soumis, du fait de la dilatation thermique ou du retrait, à une
déformation imposée 𝐷0(𝑀), la contrainte en point 𝑀 est donnée par la relation :

(𝑀) = [(𝑀) − 𝐷0(𝑀)] (III.3)

Proposons nous de déterminer quelles forces, définies par la matrice colonne 𝐹𝜇 formée par la
juxtaposition des matrices colonnes 𝐹𝜇, 𝐹𝜇, 𝐹𝜇 … …, il faut appliquer aux nœuds ∆𝛼, ∆𝛽, ∆𝛾,…….
𝛼 𝛽 𝛾

pour assurer l'équilibre de l'élément sous l'effet :

- des déplacements des nœuds définis par la matrice colonne ∆𝜇;

- d'une densité de force extérieure (𝑀) appliquée à l'élément;
- d'une déformation 𝐷0(𝑀) imposée à 1'élément.
Appliquons le théorème des travaux virtuels :
(𝛿∆𝜇)𝑡𝐹𝜇 + ∫ 𝜇[𝛿𝑈(𝑀)]𝑡𝑓(𝑀)𝑑𝑀 = ∫ 𝜇[𝛿𝐷(𝑀)]𝑡𝐶(𝑀)𝑑𝑀 (III.4)
𝐸 𝐸
Or il résulte des relations (III.1) et (III.2) que :
(𝑀) = Ф(𝑀)𝛿∆𝜇 𝑒𝑡 𝛿𝐷(𝑀) = 𝐵𝜇(𝑀)𝛿∆𝜇 (III.5)

En tenant compte de l’expression (III.7) de (𝑀), nous obtenons ainsi la relation ∶

(𝛿∆𝜇)𝑡𝐹𝜇 + (𝛿∆𝜇)𝑡∫ 𝜇[Ф(𝑀)]

𝐸
𝑡𝑓(𝑀)𝑑𝑀 = (𝛿∆𝜇 )𝑡∫ [𝐵𝜇(𝑀)]𝑡𝑟[𝐷(𝑀) − 𝐷 (𝑀)]𝑑𝑀 (III.6)
𝜇
𝐸 0

Nous déduisons de la relation précédente, exacte quelle que soit la variation virtuelle 𝛿∆𝜇 la formule
fondamentale :

𝐹𝜇 = 𝑅𝜇∆𝜇 + 𝐹'𝜇 + 𝐹''𝜇 (III.7)

Dans laquelle :

31
 𝑅𝜇 : est la matrice de rigidité généralisée de l’élément 𝐸𝜇.

𝐹'𝜇 : est la matrice qui définit les forces qu’il faut appliquer aux nœuds pour équilibrer la densité
de force (𝑀) lorsque les nœuds sont fixes (∆𝜇= 0).

𝐹''𝜇 : est la matrice qui définit les forces qu’il faut appliquer aux nœuds pour que la déformation imposée
𝐷0(𝑀) n’entraîne aucun déplacement des nœuds (∆𝜇= 0).
6.2.2 Calcul des déplacements des nœuds
La formule (III.7) donne la matrice définissant la force qu’il faux appliquer au nœud 𝐴𝑣 :

La sommation étant étendue aux nœuds de l’élément 𝐸𝜇 et :

En supposant que des forces concentrées définies par la matrice colonne 𝐹𝑣 soient appliquées au
nœud 𝐴𝑣, l’équation d’équilibre du nœud 𝐴𝑣 s’écrit :

𝐹𝑣 = ∑𝛼 𝑅𝑣𝛼∆𝛼 + 𝐹' + 𝐹'' (III.11)

𝑣 𝑣
Avec :

𝑅 = ∑ 𝑅𝜇 𝐹' = ∑ 𝐹'𝜇 𝐹'' = ∑ 𝐹''𝜇

𝑣𝛼 𝜇 𝑣𝛼 𝑣 𝜇 𝑣 𝑣 𝜇 𝑣

La première sommation étant étendue à l'ensemble des éléments ayant en commun les nœuds
𝐴𝑣 et 𝐴𝛼 et les deux dernières à l'ensemble des éléments ayant en commun le nœud 𝐴𝑣, Si l'on
pose :
∆= 𝐶𝑜𝑙(∆𝑣), 𝐹 = 𝐶𝑜𝑙(𝐹𝑣), 𝐹' = 𝐶𝑜𝑙(𝐹'), 𝐹'' = 𝐶𝑜𝑙(𝐹'')
𝑣 𝑣
L’ensemble des équations d'équilibre des nœuds équivaut à l'unique équation matricielle :

𝐹 = 𝑅∆ + 𝐹' + 𝐹'' (III.12)

Dont la résolution exige l'inversion de la matrice de raideur 𝑅. Cette inversion est facilitée parce
que la matrice 𝑅 a beaucoup d'éléments nuls; en effet, toute sous-matrice 𝑅𝛼𝛽 relative à deux
nœuds 𝐴𝛼 et 𝐴𝛽 qui n'appartient pas à un ou plusieurs éléments est nulle. En numérotant
convenablement les éléments, la matrice 𝑅 prend la forme :

32
 Et il est facile de montrer que 1'inversion de la matrice R revient à inverser un certain nombre
de matrices dont la somme des ordres est égale à l'ordre de la matrice 𝑅.

6.3 Structures élastiques plans

6.3.1 Elément triangulaire où quadrilatère simple pour
l’analyse des plaques
La structure est divisée en triangles où rectangle. Un élément 𝐸𝜇 comporte trois nœuds pour un
élément triangulaire, et quatre nœuds pour un élément quadrilatère, dans cette paragraphe, on
prend le cas le plus simple c’est le cas du triangle de nœud 𝐴𝑖, 𝐴j, 𝐴𝑘 rangés dans l'ordre où on les
rencontre en parcourant le périmètre du triangle dans le sens trigonométrique ([Link].1). Nous
désignons par (𝑥𝑖, 𝑦𝑖), (𝑥j, 𝑦j) et (𝑥𝑘, 𝑦𝑘) les coordonnées des nœuds dans un repère orthonormé
fixe 𝑂𝑥𝑦.
Si l'on représente les composantes (𝑢, 𝑣) du déplacement au point (𝑥, 𝑦) de l'élément 𝐸𝜇 par des
fonctions linéaires :
𝑢(𝑥, 𝑦) = 𝑎 + 𝑏𝑥 + 𝑐𝑦
{ (III.14)
(𝑥, 𝑦) = 𝑎′ + 𝑏′𝑥 + 𝑐′𝑦

On dispose de six constantes que l'on peut déterminer en fonction des

déplacements (𝑢𝑖, 𝑣𝑖), (𝑥j, 𝑦j), (𝑥𝑘, 𝑦𝑘), des nœuds 𝐴𝑖, 𝐴j, 𝐴𝑘.

Il est intéressant d'utiliser les coordonnées trilinéaires liées au triangle 𝐴𝑖 𝐴j 𝐴𝑘; si 𝑑𝑖, 𝑑j 𝑒𝑡 𝑑𝑘
désignent les distances du point 𝑀 aux côtés du triangle dont les hauteurs sont ℎ𝑖, ℎj 𝑒𝑡 ℎ𝑘 :

Les coordonnées trilinéaires s'expriment en fonction des coordonnées (x, y ) par les formules :

33
Dans lesquelles :

34
 𝑎𝑖 = 𝑥j𝑦𝑘 − 𝑥𝑘𝑦j 𝑏𝑖 = 𝑦j − 𝑦𝑘 𝑐𝑖 = 𝑥𝑘 − 𝑥j
{𝑎j = 𝑥𝑘𝑦𝑖 − 𝑥𝑖𝑦j {𝑏j = 𝑦𝑘 − 𝑦𝑖 {𝑐j = 𝑥𝑖 − 𝑥𝑘
𝑎𝑘 = 𝑥𝑖𝑦j − 𝑥j𝑦j 𝑏𝑘 = 𝑦𝑖 − 𝑦j 𝑐j = 𝑥𝑖 − 𝑥𝑘

2𝐴 = 𝑎𝑖 + 𝑎j + 𝑎𝑘

𝐴 est l′airedu triangle.

Inversement (𝑥, 𝑦) s’exprime en fonction de (𝜉𝑖, 𝜉j, 𝜉𝑘) par :

y 𝐴𝑘 y

0 x 0 x
Figure 7.2(Elément triangulaire simple
Figure 7.1(Elément quadrilatère simple
d’une structure)
d'une structure)

Les coordonnées trilinéaires vérifient la relation :

𝜉𝑖 + 𝜉j + 𝜉𝑘 = 1 (III.18)

Les coordonnées trilinéaires des nœuds 𝐴𝑖, 𝐴j 𝑒𝑡 𝐴𝑘 sont (1, 0, 0), (0, 1, 0), (0, 0, 1). Les
formules d’interpolations s’écrivent donc :

Ou sous forme matricielle :

(𝑀) = Φ𝜇 ∆𝜇= 𝜉𝑖∆𝑖 + 𝜉j∆j + 𝜉𝑘∆𝑘 (III.20)

La fonction d'interpolation choisie entraîne la continuité de 𝑢 et de 𝑣 à la traversée de la frontière

commune 𝐴𝑖, 𝐴j à deux éléments, puisque 𝑢 𝑒𝑡 𝑣 varient linéairement le long de
𝐴𝑖 𝐴j et que les déplacements (𝑢𝑖, ) et (𝑢j, 𝑣j ) sont les mêmes pour les deux éléments [37].

35
 6.3.2 Déformation
En un point 𝑀 de l'élément, la déformation est définie par la matrice colonne :

formule ion :
Soit, compte d’interpolat
tenu de la

6.3.3 Contrainte
Elle est définie par la matrice colonne [37]:

(𝑀) = [𝜎1 𝜎2 𝜏3] (III.23)

Les fonctions d’élasticité plane montrent que la contrainte (𝑀) se déduit de la déformation
(𝑀) par la relation :

(𝑀) = (𝑀) (III.24)

Dans laquelle la matrice d’élasticité 𝑟 a pour expression, dans le cas de la contrainte plane :

Et dans le cas de la déformation plane :

Avec : 𝐸 module d'Young

coefficient de Poisson

6.3.4 Matrice de rigidité de l'élément

Elle est de la forme :

Les sous-matrices sont données par la première formule (III.10). La fonction à intégrer se réduisant

36
 à une constante, on trouve, dans le cas de la contrainte plane, en désignant par 𝑒 l'épaisseur de
l'élément : et, dans le cas de la déformation plane, pour une tranche d'épaisseur égale à l'unité.

6.3.5 Forces nodales équilibrant une densité de force extérieure

Lorsque les composantes (𝑝, 𝑞) de la densité de force extérieure sont constantes, la seconde
formule donne, moyennant des calculs simples :

6.3.6 Effet d'une variation linéaire

Une dilatation uniforme 𝜏 causée par la température ou le retrait impose à l'élément une
déformation :

La dernière formule (III.10) donne, moyennant des calculs simples :

Nous avons ainsi calculé tous les termes qui figurent dans la formule fondamentale (III.7). Nous
pouvons donc calculer les déplacements des nœuds et en déduire la contrainte constante
𝐶𝜇 = 𝑟𝐵𝜇∆𝜇 (III.29)

de l'élément . Il est naturel de supposer que la contrainte ainsi calculée est la contrainte au centre
de gravité de 1'élément, et que la contrainte en un nœud est la moyenne des contraintes calculées
pour les éléments qui aboutissent à ce nœud.

6.4 Plaques minces

6.4.1 Généralités
Considérons une plaque mince dont le plan moyen 𝑂𝑥𝑦 devient après déformation une surface
𝑆 d'équation 𝑍 = 𝑤 (𝑥, 𝑦); le plan tangent à 𝑆 au point 𝑀 (𝑥, 𝑦) a pour équation [38] :

37
 Donc le déplacement d'un élément de plaque de centre M est le produit d'une translation de
composante 𝑤 suivant 𝑂𝑧, et d'une rotation de composantes 𝛼 𝑒𝑡 𝛽 suivant 𝑂𝑥 𝑒𝑡 𝑂𝑦 :

Le déplacement de l'élément de plaque de centre 𝑀 est défini par la matrice colonne :

(𝑀) = [𝑤 𝛼 𝛽] (III.31)

La plaque est sollicitée par une densité de force généralisée définie par la matrice colonne :

𝑡
(𝑀) = [(𝑀) 𝛾(𝑀) 𝛾𝑦(𝑀)] (III.32)

𝑞 (𝑀) est la composante de la densité de force suivant 𝑂𝑧 ; (𝑀) 𝑒𝑡 𝛾(𝑀) sont les composantes
de la densité de couple suivant 𝑂𝑥 𝑒𝑡 𝑂𝑦 .

La méthode des éléments finis consiste à remplacer la plaque par une structure composée
d’éléments polygonaux 𝐸𝜇 assemblés entre eux de façon rigide aux nœuds 𝐴𝑣. Les déplacements
des nœuds et les forces nodales appliquées à l'élément 𝐸𝜇 sont définis par les matrices colonnes
[38]:

La contrainte de la plaque est définie par la matrice colonne :

𝑡
(𝑀) = [𝑀𝑥 𝑀𝑦 𝑀𝑥𝑦]
(III.34)

Dans laquelle 𝑀𝑥, 𝑀𝑦 𝑒𝑡 𝑀𝑥𝑦 désignent les moments fléchissant et le moment de torsion. La
déformation de la plaque est définie par la matrice colonne :

La théorie des plaques montre que :

[Link].[Link].1 (𝑀) = (𝑀) (III.36)

38
La matrice d’élasticité 𝑟 ayant pour expression :

Avec : 𝐸 module d′Young

ℎ épaisseur de la plaque
coefficient de Poisson

On notera que la connaissance de la fonction 𝑤 (𝑥, 𝑦) entraîne celle des matrices (𝑀) et (𝑀).

6.4.2 Utilisation d'éléments triangulaires et quadrilatères

Nous supposons la plaque divisée en éléments triangulaires 𝐴𝑖 𝐴j 𝐴𝑘, ([Link].2) et en élément
quadrilatère Ai Aj Ak Al ([Link].1), l'ordre Ai Aj Ak Al correspondant au parcours du périmètre du
rectangle ou du triangle dans le sens trigonométrique.

Si nous représentons le déplacement transversal w de l’élément Eµ par un polynôme du troisième

degré en x et y, ce polynôme contient dix coefficients. Comme il n’y a que neuf composantes de
déplacement des nœuds pour l’élément triangulaire, le polynôme dépendra linéairement d'un
paramètre que l'on pourrait déterminer en écrivant que l'énergie, de déformation est minimale.

Nous procéderons de façon différente. Nous étudierons d'abord l'élément Ai Aj Ak simplement

appuyé aux nœuds ; dans ce cas wi = wj = wk = 0 et les réactions verticales des appuis se
déterminent à l’aide des seules équations d’équilibre de la statique. Nous passerons ensuite au cas
général en donnant un déplacement d’ensemble à l'élément supposé sur appuis simples.

Pour conserver la symétrie de l'élément, nous utiliserons les coordonnées trilinéaires ((i, (j, (k)
définies au paragraphe III.4.1.1, et nous conserverons les notations définies dans ce paragraphe.

6.4.3 Étude de l'élément 𝐴𝑖 𝐴j 𝐴𝑘 simplement appuyé aux nœuds

Nous accentuons toutes les quantités relatives à l’élément 𝐴𝑖 𝐴j 𝐴𝑘, simplement appuyé. Les
déplacements des nœuds sont définis par les matrices [39]:

Et les forces nodales correspondantes sont définies par les matrices :

Représentons le déplacement transversal w′ par un polynôme du troisième degré en 𝑥 et 𝑦; ce polynôme

39
doit vérifier les conditions :

Ce polynôme dépend donc d'un paramètre λ; en observant que la fonction (i (j (k ainsi que ses
dérivées par rapport à x et à y sont nulles aux nœuds Ai Aj et Ak, on met sans difficultés
w′(x, y)sous la forme :

Les coefficients 𝑏𝑖, 𝑏j, 𝑏𝑘, et 𝑐𝑖, 𝑐j, 𝑐𝑘 sont définis au paragraphe III.3.1.1.

Lorsque l’élément est soumis à un état de déformation uniforme, les dérivées secondes de w′ sont
des constantes, et w′ est un polynôme du second degré qui s'annule aux nœuds, donc de la forme:

Si nous voulons que l'expression (III.37) puisse représenter un état de déformation uniforme, il
faut que cette expression se réduise à l'expression (III.38) lorsqu'on donne aux rotations les
valeurs correspondant à un état de déformation uniforme défini par (III.38). Un calcul simple
montre que cela n'est possible que lorsque 𝜆 = 1/2, valeur que nous adopterons désormais
dans la fonction d'interpolation (III.37). Cette fonction d'interpolation n'entraîne pas la continuité
du déplacement (𝑀)à la traversée de la frontière commune à deux éléments.
La déformation 𝐷 (𝑀) se déduit de l'expression de 𝑤′.

𝛛2w′ 𝛛2w′ 𝛛 2w′ (III.39)

(𝑀) = [ 2 ]=∑ 𝛽′(𝑀)∆′

𝛛𝑥2 𝛛𝑦2 𝑖,j,𝑘 𝑖 𝑖

𝛛𝑥𝛛𝑦

β′µ
i (M) Étant une matrice de dimension (3,2) dont il serait facile de donner 1'expression.
Les sous-matrices de la matrice de rigidité 𝑅′𝜇 de l'élément 𝐸𝜇 simplement appuyé aux nœuds
sont alors données par la première formule (III.10) :

Le calcul de 𝑅𝑚𝑛
′𝜇 peut être achevé en utilisant la valeur des intégrales :

40
41
 6.4.4 Matrice de rigidité généralisée de l'élément 𝐴𝑖 𝐴j 𝐴𝑘
Au cours de la déformation de la plaque, les nœuds 𝐴𝑖, 𝐴j 𝑒𝑡 𝐴𝑘 se déplacent normalement au
plan moyen et viennent en 𝐴′, 𝐴′, 𝑒𝑡 𝐴′ [39] :

L’équation du plan Ai Aj Ak est :

𝑧 = 𝑤𝑖𝜉𝑖 + 𝑤j𝜉j + 𝑤𝑘𝜉𝑘 (III.42)

Il ne résulte que le déplacement normal 𝑤 a pour valeur :

𝑤 = 𝑤′ + 𝑤𝑖𝜉𝑖 + 𝑤j𝜉j + 𝑤𝑘𝜉𝑘 (III.43)

𝑤′ Étant donné par la formule (III.38) dans laquelle 𝜆 = 1/2. En dérivant l'expression précédente
par rapport à 𝑥 𝑒𝑡 à 𝑦, nous trouvons :

relations équivalentes à relation matricielle :

En
plaça
nt
succe
ssive
ment
𝑀 en
𝐴𝑖, 𝐴j
𝑒𝑡 𝐴𝑘,
nous 𝑇 Désignant
obten la matrice de
ons dimension
ainsi (6.9) :
six

Les équations de la statique permettent de déterminer les réactions verticales d'appui de l'élément
simplement appuyé. Nous savons donc calculer les forces nodales 𝐹𝑢 connaissant les forces

42
 nodales 𝐹′𝑢. Nous trouvons la relation:

𝐹𝜇 = 𝑇𝑡𝐹′𝜇 (III.45)

Que l'on peut également démontrer au moyen du théorème des travaux virtuels
(𝛿∆𝜇)𝑡 𝐹𝜇 = (𝛿∆′𝜇)′𝜇 (III.46)

Nous déduisons successivement des relations précédentes :

𝐹𝜇 = 𝑇𝑡𝐹′𝜇 = 𝑇𝑡𝑅′𝜇∆′𝜇= (𝑇𝑡𝑅′𝜇∆′𝜇)∆𝜇 (III.47)

Donc la matrice de rigidité généralisée 𝑅𝜇 se déduit de la matrice de rigidité de l'élément

supposé sur appuis simples par la formule :

𝑅𝜇 = 𝑇𝑡𝑅′𝜇𝑇 (III.48)

43
 CATIONS
A

7 APL NUMERIQUES
Après la présentation de la théorie relative au comportement des plaques, et la formulation de la nouvelle
variante de l’élément SBRPS. Il est évident que la solution des équations différentielles gouvernant, le
comportement des plaques est impossible par une des méthodes analytiques. Pour cette raison, le recours à
une méthode numérique est devenu une obligation. Pour cela, la méthode des éléments finis est sélectionnée
comme une méthode numérique de cette étude. La structure d’un programme de calcul en FORTRAN77 est
formée d'un ensemble de blocs, ces derniers sont constitués par les différentes soubroutines du programme
source. En plus en sa forme, il est totalement identique aux programmes développés dans le domaine du
calcul numérique.
La validation de l’élément formulé, est représentée par la comparaison entre les résultats analytiques, indiqués
dans la littérature et les résultats donnés par les autres éléments à travers des tests similaires.
A cet effet une étude comparative entre les résultats obtenus à l’aide du programme FORTRAN77 a été
réalisée avec ceux des logiciels de calcul tel que (ANSYS), (SAP2000) et (RDM6).
Le programme de FORTRAN77, est structuré par un ensemble de blocs. Sa composition passe par des
opérations très diverses, à savoir : (1) Organisation des données, (2) Intégration analytique, (3) Résolution du
système d'équations. Le traitement numérique s'effectue au niveau des éléments finis, pour cela ce programme
inclut quelques blocs fonctionnels caractérisés par les étapes suivantes :
Lecture, vérification et organisation des données décrivant le maillage (nœuds et éléments), les paramètres
physiques (conductivités, module d'élasticité, etc....), les sollicitations et conditions aux limites ;
Construction des matrices et vecteurs élémentaires, puis assemblage de celle ci pour former la matrice globale
et le vecteur global des sollicitations;
Résolution du système d'équation après prise en compte des conditions aux limites;
Impression des résultats après calcul éventuel des variables additionnelles (déplacements, moments)
Développement du programme pour le calcul des plaques
On considère une plaque homogène soumise à une charge uniformément repartie ou concentrée, qu’on
discrétise en éléments rectangulaires égaux, (𝑛) éléments suivant l’axe 𝑥 et (𝑛) éléments suivant l’axe 𝑦, on
calcule les déplacements et les rotations (𝑤, , 𝜃𝑦) au niveau de chaque nœud, en fin on procède au calcul des
réactions et des contraintes dans le milieu de chaque élément, avec la détermination des valeurs maximales.

44
Les processus de calcul dans le programme se font comme suit :

Discrétisation de la plaque en éléments finis, puis définition du nombre d’éléments selon les axes(𝑥, 𝑦).
Détermination des caractéristiques physiques et géométriques de la plaque.
Définition du chargement de la plaque.
Numérotation des nœuds, puis définition des conditions aux limites.
Assemblage des éléments pour la détermination de la matrice de rigidité de la plaque.
Résolution du système d’équations.
Détermination des déplacements et des rotations des nœuds.
En fin, vient le calcul des contraintes au niveau de chaque élément.

La structure du programme de l’analyse des plaques se présente selon l’ordre suivant :

Un programme principal qui regroupe la grande partie des opérations effectuées, allant de la lecture des
données dans un fichier texte, dont le nom (donné[Link]), et finissant par trouver les résultats dans un autre
fichier texte, dont le nom (ré[Link]).
Un fichier d’introduction des données dans lequel l’utilisateur devra mettre les valeurs numériques
concernant l’application voulue, selon l’ordre suivant :

La dimension de la plaque suivant l’axe x (a).

la dimension de la plaque suivant l’axe y (b).
Le nombre des nœuds par élément.
Le nombre de degrés de liberté par élément.
Le nombre de degrés de liberté par Nœud.
L’épaisseur de la plaque (h), Module de Young (E), coefficient de poisson ( ).
Type de chargement et sa valeur.

Un fichier de résultat, où sont détaillés tous les résultats, à savoir :

Valeurs des déplacements et des rotations au niveau de chaque nœud.

Valeurs des contraintes.
Valeurs des réactions au niveau des appuis.

Pour plus de détails, voir ce programme sous forme d’organigramme en (figure. V.1).

45
 Début

Introduction du nombre d’éléments

suivant l’axe x (n) et suivant l’axe y (n)
de la plaque.

Maillage (Numérotation des nœuds pour chaque

élément de la plaque)

Introduction des caractéristiques

géométriques et physiques de la plaque

Calcul de la matrice de rigidité élémentaire [Ke]

pour chaque élément

Assemblage des matrices élémentaires [Ke] pour

obtenir la matrice globale [K°] globale

Introduction des conditions aux limites

et des forces appliquées aux

Nœuds

Résolution du système d’équations

[OK] {𝛿} = {f}

Détermination des déplacements

de chaque nœud, et des contraintes
dans le milieu de chaque élément

Fin

Figure 8(: Organigramme pour l’étude de la flexion des plaques)

46
8 Conclusion
Comme on a vue ci-dessus, il existe donc de différents types de plaques (isotropes,
anisotropes et orthotropes) qui se comportent avec différents comportements (flexionnel,
membrane ou mixte), dont les études se faits dans différents états (statique, dynamique). Ces
études se basent généralement sur différents théories, principalement on utilise celle de Love -
Kirchhoff dans le cas des plaques minces, et celle de Reissner - Mindlin dans le cas des
plaques épaisses.

47