v3.02.306 SSLP306 - Plaque encastrée soumise à une charge transverse polynomiale#

Résumé:

Le test a pour but de valider l’élément créé automatiquement à partir de Fenicsx basé sur la théorie de Reissner-Mindlin et la procédure MITC (en anglais, Mixed Interpolation of Tensorial Component).

La solution de référence est analytique.

Le nouvel élément fini PLAQ_MITC (discrétisé avec QUAD9 et TRI6) est comparé avec les éléments DKT et DST.

Solution de référence#

Méthode de calcul utilisée pour la solution de référence#

Pour les modélisations A, B, C, la valeur du déplacement vertical et les valeurs des rotations autour des axes \(X\) et \(Y\) sont données par:

\[\begin{split}\begin{align} w(x, y) &= \frac{1}{3} x^{3} (x - 1)^{3} y^{3} (y - 1)^{3} - \frac{2 h^{2}}{5 (1 - \nu)} \Big[ y^{3} (y - 1)^{3} x (x - 1) (5 x^{2} - 5 x + 1) \\ & \quad + x^{3} (x - 1)^{3} y (y - 1) (5 y^{2} - 5 y + 1) \Big] \\ \theta_x(x, y) &= x^{3} (x - 1)^{3} y^{2} (y - 1)^{2} (2 y - 1) \\ \theta_y(x, y) &= -y^{3} (y - 1)^{3} x^{2} (x - 1)^{2} (2 x - 1) \end{align}\end{split}\]

Fig 1. Solution de référence pour les modélisations A, B, C.

Pour la modélisation D, à cause de la rotation de la géométrie, l’expression du déplacement vertical \(w\) s’obtient en faisant le changement de variable :

\[\begin{split}\begin{aligned} x \rightarrow \cos(\theta) x + \sin(\theta) y \\ y \rightarrow -\sin(\theta) x + \cos(\theta) y \\ \textrm{où} \; \theta = \pi / 6 \end{aligned}\end{split}\]

Alors que pour les rotations \(\theta_x,\,\theta_y\) il faut de plus inclure la transformation :

\[\begin{split}\begin{aligned} \theta_x \rightarrow \cos(\theta) \theta_x - \sin(\theta) \theta_y \\ \theta_y \rightarrow \sin(\theta) \theta_x + \cos(\theta) \theta_y \\ \textrm{où} \; \theta = \pi / 6 \end{aligned}\end{split}\]

Fig 2. Solution de référence pour la modélisation D.

Résultats de référence#

Pour les modélisations A, B, C, E :

Déplacement au point \((0.5,\,0.5)\) avec \(h = 0.1\), \({w} = 9.254092261904761\cdot{10}^{-5} \mathrm{m}\)
Déplacement au point \((0.5,\,0.5)\) avec \(h = 0.001\), \({w} = 8.13813244047619\cdot{10}^{-5} \mathrm{m}\)

Pour la modélisation D :

Déplacement au point \((\frac{1}{\sqrt{2}}\cos(\frac{5\pi}{12}),\,\frac{1}{\sqrt{2}}\sin(\frac{5\pi}{12}))\) avec \(h = 0.1\), \({w} = 9.254092261904761\cdot{10}^{-5} \mathrm{m}\)
Déplacement au point \((\frac{1}{\sqrt{2}}\cos(\frac{5\pi}{12}),\,\frac{1}{\sqrt{2}}\sin(\frac{5\pi}{12}))\) avec \(h = 0.001\), \({w} = 8.13813244047619\cdot{10}^{-5} \mathrm{m}\)

Incertitude sur la solution#

Solution analytique.

Références bibliographiques#

Hale JS, Brunetti M, Bordas SP, Maurini C. Simple and extensible plate and shell finite element models through automatic code generation tools. Computers & Structures. 2018 Oct 15;209:163-81.
Chinosi C, Lovadina C. Numerical analysis of some mixed finite element methods for Reissner-Mindlin plates. Computational Mechanics. 1995 Apr;16(1):36-44.

Modélisation A#

Caractéristiques de la modélisation#

C’est une modélisation PLAQ_MITC fait via fenicsx.

Le maillage de \(\Omega=[0,1]^2\) est directement construit par:

mesh = CA.Mesh.buildSquare(refine=FINE_REFINEMENT)

Conditions limites:

sur les bords \({0, 1}\times[0, 1] \cup [0, 1]\times{0, 1}\):

MECA_IMPO=_F(DZ=0, DRX=0, DRY=0, GROUP_MA=("LEFT", "RIGHT", "TOP", "BOTTOM"))

Remarque:: A présent, PLAQ_MITC ne considère que la mouvement de type COQUE et on n’impose donc pas le déplacements DX, DY, DRZ.

Chargement:

PRES_REP=_F(GROUP_MA="SURFACE", PRES=PRES)

où PRES est une FORMULE qui représente l’équation dans le problème de référence.

Caractéristiques du maillage#

Nombre de nœuds : 16641

Nombre de mailles et types : 4352, 256 SEG3, 4096 QUAD9

Valeurs testées#

Localisation	Épaisseur	Type de valeur	Référence	Aster	% différence
Point \((0.5, 0.5)\)	0.1	DZ	\(9.254092261904761\cdot{10}^{-5}\)	\(9.225572268464905\cdot{10}^{-5}\)	< 1
Point \((0.5, 0.5)\)	0.001	DZ	\(8.13813244047619\cdot{10}^{-5}\)	\(8.113126374037738\cdot{10}^{-5}\)	< 1

Modélisation B#

Caractéristiques de la modélisation#

C’est une modélisation DST.

Conditions limites:

sur les bords \({0, 1}\times[0, 1] \cup [0, 1]\times{0, 1}\):

MECA_IMPO=_F(DX=0, DY=0, DZ=0, DRX=0, DRY=0, DRZ=0, GROUP_MA=("LEFT", "RIGHT", "TOP", "BOTTOM"))

Chargement:

FORCE_COQUE=_F(GROUP_MA="SURFACE", PRES=PRES)

où PRES est une FORMULE qui représente l’équation dans le problème de référence.

Caractéristiques du maillage#

Nombre de nœuds : 4225

Nombre de mailles et types : 4352, 256 SEG2, 4096 QUAD4

Valeurs testées#

Localisation	Épaisseur	Type de valeur	Référence	Aster	% différence
Point \((0.5, 0.5)\)	0.1	DZ	\(9.254092261904761\cdot{10}^{-5}\)	\(9.220016849373356\cdot{10}^{-5}\)	< 1
Point \((0.5, 0.5)\)	0.001	DZ	\(8.13813244047619\cdot{10}^{-5}\)	\(8.10641555340742\cdot{10}^{-5}\)	< 1

Modélisation C#

Caractéristiques de la modélisation#

C’est une modélisation DKT.

Conditions limites :

Idem que modélisation B.

Chargement:

Idem que modélisation B.

Caractéristiques du maillage#

Idem que modélisation B.

Valeurs testées#

Localisation	Épaisseur	Type de valeur	Référence	Aster	% différence
Point \((0.5, 0.5)\)	0.1	DZ	\(9.254092261904761\cdot{10}^{-5}\)	\(8.106306974024127\cdot{10}^{-5}\)	< 20
Point \((0.5, 0.5)\)	0.001	DZ	\(8.13813244047619\cdot{10}^{-5}\)	\(8.106306973545047\cdot{10}^{-5}\)	< 1

Remarque :: Le modèle DKT fonctionne moins bien lorsque l’épaisseur de la plaque augmente. Une erreur importante est attendue pour le cas \(h=0,1\).

Modélisation D#

Caractéristiques de la modélisation#

C’est une modélisation PLAQ_MITC fait via fenicsx.

Le géométrie est définie par \(\Omega=\{(x,y) | 0 \le \cos(\pi/6) x + \sin(\pi/6) y \le 1, 0 \le -\sin(\pi/6) x + \cos(\pi/6) y \le 1\}\).

Conditions limites :

sur tous les bords, on impose

MECA_IMPO=_F(DZ=0, DRX=0, DRY=0, GROUP_MA=("LEFT", "RIGHT", "TOP", "BOTTOM"))

Remarque :: A présent, PLAQ_MITC ne considère que le mouvement de type COQUE et on n’impose donc pas le déplacements DX, DY, DRZ.

Chargement :

PRES_REP=_F(GROUP_MA="SURFACE", PRES=PRES)

où PRES est une FORMULE qui représente l’équation dans le problème de référence.

Caractéristiques du maillage#

Nombre de nœuds : 16641.

Nombre de mailles et types : 4352, 256 SEG3, 4096 QUAD9.

Valeurs testées#

Localisation	Épaisseur	Type de valeur	Référence	Aster	% différence
Point \((0.5, 0.5)\)	0.1	DZ	\(9.254092261904761\cdot{10}^{-5}\)	\(9.225572268464905\cdot{10}^{-5}\)	< 1
Point \((0.5, 0.5)\)	0.001	DZ	\(8.13813244047619\cdot{10}^{-5}\)	\(8.113126374037738\cdot{10}^{-5}\)	< 1

Modélisation E#

Caractéristiques de la modélisation#

Idem modélisation A.

Caractéristiques du maillage#

Nombre de nœuds : 16641

Nombre de mailles et types : 8448, 256 SEG3, 8192 TRIA6

Valeurs testées#

Localisation	Épaisseur	Type de valeur	Référence	Aster	% différence
Point \((0.5, 0.5)\)	0.1	DZ	\(9.254092261904761\cdot{10}^{-5}\)	\(9.22439938\cdot{10}^{-5}\)	< 1
Point \((0.5, 0.5)\)	0.001	DZ	\(8.13813244047619\cdot{10}^{-5}\)	\(8.11265671\cdot{10}^{-5}\)	< 1

Synthèse des résultats#

Les modélisations A, B, C, E présentent de bons résultats pour des plaques minces (à faibles épaisseurs par rapport aux autres dimensions). Le modèle C avec DKT fonctionne moins bien lorsque l’épaisseur de la plaque augmente. Ce comportement attendu s’explique à cause des déformations en cisaillement plus importantes qu’un élément fini DKT n’arrive pas à modéliser. A noter que l’erreur accroît lors que l’on discrétise davantage (voir figure ci-dessous).

../../../../_images/mitc_erreur_mince.svg

../../../../_images/mitc_erreur_epais.svg

Fig 1. Évolution de l’erreur de DZ en fonction du niveau de discrétisation pour une plaque mince, \(h=0.001\), (à gauche) et épaisse, \(h=0.1\) (à droite).

Il n’y a pas de différences entre les résultats des modélisations A et D (pas d’influence dû au changement du repère).