r4.10.02 Estimateur d’erreur en résidu#

Résumé

L’estimateur d’erreur en résidu permet d’estimer l’erreur de discrétisation due à la méthode des éléments finis sur les éléments d’un maillage 2D ou 3D. C’est un estimateur d’erreur explicite faisant intervenir les résidus des équations d’équilibre et les sauts des contraintes normales aux interfaces, contrairement à l’estimateur de Zhu‑Zienkiewicz, qui utilise une technique de lissage des contraintes à posteriori [R4.10.01] et [bib5].

Cet estimateur est implanté dans Code_Aster en élasto-plasticité 2D et 3D.

Formulation de l’estimateur en résidu#

Soit \(\Omega\) un ouvert de \({R}^{n}\) , \(n=2\) ou \(3\) , de frontière \(\Gamma\) , et \(T\) une triangulation régulière de \(\Omega\) .

En élasticité linéaire, le problème continu s’écrit:

trouver \((u,\sigma )\) tels que:

\(\lbrace \begin{array}{cc}div\sigma +f=0& \text{dans}\Omega \\ u={u}_{D}& \text{sur}{\Gamma}_{D}\\ \sigma .n={g}_{N}& \text{sur}{\Gamma}_{N}\end{array}\)

\({\Gamma}_{D}\) est la frontière de Dirichlet du maillage

\({u}_{D}\) est le déplacement imposé sur cette frontière

\({\Gamma}_{N}\) est la frontière de Neumann

\(n\) la normale unitaire à \({\Gamma}_{N}\)

\({g}_{N}\) est le chargement appliqué sur cette frontière; il peut être continu ou discrétisé.

\(f\) est une force de type volumique (pesanteur, rotation); elle peut être continue ou discrétisée.

\({\sigma}_{h}\) est la contrainte obtenue par la résolution du problème discret:

\(\lbrace \begin{array}{cc}div{\sigma}_{h}+f=0& \text{dans}\Omega \\ {u}_{h}={u}_{D}& \text{sur}{\Gamma}_{D}\\ {\sigma}_{h}.n={g}_{N}& \text{sur}{\Gamma}_{N}\end{array}\)

avec la relation \({\sigma}_{h}=\text{DB}{u}_{h}\) où:

\(D\) est la matrice de Hooke

\(B\) est l’opérateur linéarisé des déformations

Si \(K\) désigne un élément courant du maillage, l’estimateur d’erreur (noté \(\eta (\Omega )\) ) est défini comme étant la moyenne quadratique des indicateurs d’erreur locaux, notés \(\eta (K)\) :

\(\eta (\Omega )={\left[\sum_{K\in T}\eta {(K)}^{2}\right]}^{1/2}\)

L’indicateur par résidu local

L’indicateur est composé de trois termes; le premier représente le résidu de l’équation d’équilibre sur chaque maille, le deuxième terme le saut des contraintes normales sur les interfaces, le troisième terme l’écart entre les contraintes normales et le chargement imposé sur \({\Gamma}_{N}\) si l’élément intersecte \({\Gamma}_{N}\) .

../../../../_images/10000CC8000069BB000038BE5FC4A4679132172E.svg

Figure 2-a: Eléments internes à un maillage

  • le premier terme de l’estimateur est la norme \({L}^{2}\) du résidu de l’équation d’équilibre sur la maille \(K\) , multiplié par \({h}_{K}\) qui est, soit le diamètre du cercle circonscrit pour un élément fini triangulaire, soit la diagonale maximale pour un quadrangle,

  • le deuxième terme est l’intégrale, sur \(S(K)\) défini [Figure 2-a], des sauts de contraintes normales intégrés sur chaque bord \(F\) de l’élément qui possède un voisin, et multipliés par la racine de \({h}_{F}\) , qui est la longueur de l’arête \(F\) ,

../../../../_images/1000164C000069D500004D9E20F5585200AA3625.svg

Figure 2-b: Eléments situés sur la frontière d’un maillage

  • le troisième terme est l’intégrale, sur l’intersection de chaque bord \(F\) des bords \(\partial K\) de l’élément courant \(K\) avec la frontière de Neumann \({\Gamma}_{N}\) , des sauts entre les contraintes normales de l’élément et la force de Neumann \({g}_{N}\) , multipliée par la racine de \({h}_{F}\) , longueur de l’arête \(F\) .

On a donc la formule suivante pour l’estimateur en résidu:

\(\eta (K)={h}_{K}{\parallel f+div{\sigma}_{h}\parallel }_{{L}^{2}(K)}+\frac{1}{2}\sum_{F\in S(K)}{h}_{F}^{1/2}{\parallel \left[{\sigma}_{h}.n\right]\parallel }_{{L}^{2}(F)}+\sum_{F\subset \partial K\subset {\Gamma}_{N}}{h}_{F}^{1/2}{\parallel {g}_{N}-{\sigma}_{h}.n\parallel }_{{L}^{2}(F)}\) éq 2-1

Pour le choix des différents termes de [éq 2-1], on renvoie à [bib1].

Propriétés de l’estimateur en résidu#

On note

\({\eta}_{\text{EX}}(K)\) l’erreur exacte \({\parallel u-{u}_{h}\parallel }_{{H}^{1}(K)}\) sur l’élément \(K\) (inconnue a priori)

et

\({\eta}_{\text{EX}}(\Omega )\) l’erreur exacte globale \({\parallel u-{u}_{h}\parallel }_{{H}^{1}(O)}\)

On a alors les propriétés suivantes ([bib1]):

  • quelque soit l’élément \(K\) , l’erreur élémentaire \(\eta (K)\) est majorée par l’erreur locale exacte (multipliée par une constante indépendante de la triangulation),

soit \(\forall K\eta (K)\le {C}_{1}\times {\eta}_{\text{EX}}(K)\)

  • l’erreur globale exacte est majorée par l’erreur estimée globale \(\eta (\Omega )\) (multipliée par une constante indépendante de \(T\) )

soit \({\eta}_{\text{EX}}(\Omega )\le {C}_{2}\times \eta (\Omega )\)

Les constantes \({C}_{1}\text{et}{C}_{2}\) dépendent a priori du type d’élément fini et des conditions aux limites du problème. Kelly et Gago [bib3] ont proposé en 2D une constante \({C}_{2}\) ne dépendant que du degré \(p\) du polynôme d’interpolation utilisé:

\(\begin{array}{cc}{C}_{2}={(\frac{1}{24{p}^{2}})}^{1/2}\text{soit}& {C}_{2}=\frac{1}{\mathrm{2p}\sqrt{6}}\text{pour les TRIA3 et QUAD4}(\text{degré 1})\\ & {C}_{2}=\frac{1}{\mathrm{4p}\sqrt{6}}\text{pour les TRIA6 et QUAD8}(\text{degré 2})\end{array}\)

Pour le 3D, on ne dispose pas d’évaluation de la constante. On peut néanmoins dire que l’erreur estimée globale surestime l’erreur exacte globale dans tous les cas. Ce résultat n’est pas forcément vrai au niveau local.

Implantation dans Aster#

L’estimateur en résidu est implanté en 2D et 3D sur tous les types d’éléments. Il est calculé par la commande CALC_ERREUR en activant l’option “ERME_ELEM”.

Cette option calcule sur chaque élément:

  • l’erreur absolue \(\eta (K)\) (voir [éq 2-1]),

  • la norme du tenseur des contraintes \({\parallel {\sigma}_{h}\parallel }_{{L}^{2}(K)}\) qui sert à normer l’erreur absolue,

  • l’erreur relative \({\eta}_{\text{rel}}(K)=100\times \frac{\eta (K)}{\sqrt{\eta {(K)}^{2}+{\parallel {\sigma}_{h}\parallel }_{{L}^{2}(K)}^{2}}}\) .

Nota:

Cette définition de l’erreur relative implique que dans les zones où les contraintes sont très faibles, l’erreur relative peut être importante et non significative. C’est alors l’erreur absolue qu’il faut considérer.

Elle calcule également au niveau global:

  • l’erreur absolue \(\eta (\Omega )={\left[\sum_{K\in T}\eta {(K)}^{2}\right]}^{1/2}\) ,

  • la norme globale du tenseur des contraintes \({\parallel {\sigma}_{h}\parallel }_{{L}^{2}(\Omega )}={\left[\sum_{K\in T}{\parallel {\sigma}_{h}\parallel }_{{L}^{2}(K)}^{2}\right]}^{1/2}\) ,

  • l’erreur relative \({\eta}_{\text{rel}}(\Omega )=100\times \frac{\eta (\Omega )}{\sqrt{\eta {(\Omega )}^{2}+{\parallel {\sigma}_{h}\parallel }_{{L}^{2}(\Omega )}^{2}}}\) .

D’après l’expression [éq 2-1], on voit que pour calculer l’indicateur d’erreur sur la maille \(K\) , on doit connaître:

  • les chargements éventuels \(f\) sur \(K\) et \({g}_{N}\) sur \(\partial K\cap {\Gamma}_{N}\) (ou leur discrétisation \({f}_{h}\) et \({g}_{\text{Nh}}\) ),

  • les quantités \({h}_{K}\) , \({h}_{F}\) et \(n\) liées à la géométrie de l’élément,

  • le champ de contraintes \({\sigma}_{h}\) ,

  • la liste des voisins de \(K\) pour récupérer les contraintes sur ces éléments, nécessaires au calcul du 2ème terme de [éq 2-1].

1 et 2 peuvent être calculés ou récupérés facilement.

3 doit avoir été calculé au préalable par une des options “SIGM_ELNO” ou “SIEF_ELNO”.

Dans le cas contraire, un message d’erreur fatale est émis.

4 nécessite le calcul d’une connectivité particulière maille-mailles, en plus de la connectivité standard maille-nœuds. Ce nouvel objet est stocké dans la structure de données de type maillage.

Pour le détail de l’implantation dans Aster , voir [bib6].

Pour la validation de l’estimateur, voir [bib7].

Bibliographie#

  1. BERNARDI, B. METIVET, R. VERFÜRTH: Groupe de travail maillage adaptatif: analyse numérique d’indicateurs d’erreur. Note HI-72/93/062

  2. BABUSKA, RHEINBOLDT: A posteriori error estimates for the finite element method. Int. J. Num. Meth. Eng., Vol 12, 1597-1659, 1978

  3. J.P. de S.R. GAGO: A posteriori error analysis and adaptivity for the finite element method. Ph. D. Thesis, University of Wales, Swansea, U.K., 1982

  4. D.W. KELLY, J.P. de S.R. GAGO, O. C. ZIENKIEWICZ, I. BABUSKA: A posteriori error analysis and adaptative processes in finite element method: part I - error analysis. Int. J. Num. Meth. Eng., Vol 19, 1593-1619, 1983

    1. DESROCHES: Estimateurs d’erreur en élasticité linéaire. Note HI-75/93/118

    1. CORBEL: Implantation d’un estimateur d’erreur en résidu dans le Code de mécanique Aster. Rapport de fin de stage - Juin 94

    1. NAVAB: Validation d’un estimateur d’erreur en résidu en élasticité bi et tri-dimensionnelle. Rapport de stage - Mars 95

Description des versions du document#

Version Aster

Auteur(s) Organisme(s)

Description des modifications

04/07/09

X.DESROCHES (EDF-R&D/AMA)