r4.10.06 Estimateurs d’erreur en quantités d’intérêt

Résumé:

Dans ce document, on présente le cadre de l’estimation d’erreur en quantité d’intérêt. Ce cadre est utilisé avec des estimateurs en résidu explicites ainsi qu’avec des estimateurs en lissage. Un certain nombre de quantités d’intérêt sont détaillées, notamment leur écriture en terme de chargement pour le problème dual.

Table des Matières

Estimateurs d’erreurs a posteriori

Définition des estimateurs a posteriori

En estimation a posteriori , l’objectif n’est pas de trouver une estimation de la fonction erreur \(e\) mais de déterminer une estimation d’une mesure de l’erreur. Ainsi les estimateurs qui ont été développés ces trente dernières années fournissent une estimation fiable et précise de l’erreur dans des normes globales, telles que les normes \({H}^{1}\) , normes \({L}^{2}\) ou encore normes en énergie, en post-traitement de la solution éléments finis. Ce choix de normes globales est imposé par la forme bilinéaire. Quand la forme est définie positive et symétrique, comme c’est le cas ici, cela induit un produit scalaire pour lequel la norme associée est la norme en énergie; il est ainsi naturel d’estimer l’erreur dans cette norme.

Le principe de base de ce type d’estimation est d’utiliser la solution approchée pour estimer l’erreur de discrétisation. Par opposition à l’estimation a priori , les estimations a posteriori ne peuvent se faire qu’une fois la solution approchée calculée.

Une fonction \(\eta (h,{u}^{h},d)\) est une estimation d’erreur a posteriori si:

(1581)\[\parallel u-{u}^{h}\parallel \le \eta (h,{u}^{h},d)\]

où \(\parallel \cdot \parallel\) est une norme sur les champs de déplacements, \(h\) la taille des éléments, \(d\) un ensemble de données du problème et \({u}^{h}\) la solution approchée. De plus, si \(\eta (h,d,{u}^{h})\) peut être localisée sous la forme:

(1582)\[\eta (h,{u}^{h},d)={(\sum_{E}{\eta}_{E}{({u}^{h},d)}^{2})}^{\frac{1}{2}}\]

alors les quantités \({\eta}_{E}({u}^{h},d)\) , contributions élémentaires de l’estimation de l’erreur globale \(\eta (h,{u}^{h},d)\) , sont appelées indicateurs locaux d’erreur. Ils fournissent une base pour l’adaptation de maillages.

De nombreux estimateurs existent et peuvent être classés en trois catégories:

• Les estimateurs d’erreurs basés sur les défauts d’équilibre issus des travaux de Babuška et Rheinboldt [bib3];

• Les indicateurs d’erreurs construits à partir de contraintes lissées issus des travaux de Zienkiewicz et Zhu [bib4];

• Les mesures d’erreurs fondées sur le concept d’erreur en relation de comportement issues des travaux de Ladevèze [bib5].

Ces différents estimateurs, aujourd’hui bien maîtrisés, seront présentés de manière plus ou moins détaillée. Une revue plus détaillée pourra être trouvée dans [bib6] à [bib10].

Qualité des estimateurs

Lorsqu’un estimateur est défini, il faut s’interroger sur les critères et les moyens permettant de juger des performances de cet estimateur. De manière générale, il faut chercher à construire un estimateur dont le comportement asymptotique (quand la taille \(h\) des éléments tend vers zéro) suit celui de l’erreur. Ce comportement se traduit par l’existence de deux constantes \({C}_{1}\) et \({C}_{2}\) , dépendant des données du problème et de la discrétisation mais pas de la taille des éléments, vérifiant la relation suivante (Ladevèze & Pelle, 2004):

(1583)\[{C}_{1}\eta \le \parallel e\parallel \le {C}_{2}\eta\]

où \(\eta\) représente l’estimation de la mesure \(\parallel e\parallel\) de l’erreur \(e\) sur le domaine considéré. Ainsi pour juger des performances des estimateurs, il existe différents critères intrinsèques.

Indice d’efficacité

L’indice d’efficacité \(\gamma\) est défini comme le rapport entre l’erreur calculée par un estimateur \({e}_{\text{estimée}}\) et l’erreur vraie \({e}_{\text{vraie}}\) :

(1584)\[\gamma =\frac{{e}_{\text{estimée}}}{{e}_{\text{vraie}}}\]

A moins de disposer d’une solution analytique, l’erreur vraie est calculée comme la différence entre une solution obtenue sur un maillage très fin (solution overkill ) et la solution obtenue sur un maillage donné. Un indice d’efficacité proche de l’unité caractérise un bon estimateur. Si cette propriété est atteinte quand la taille des éléments tend vers zéro, l’estimateur est dit asymptotiquement exact. Néanmoins, l’estimation doit être suffisamment précise (\(\gamma\) proche de 1) pour des maillages assez grossiers afin d’être utilisable pour les discrétisations courantes en ingénierie. Enfin, il est souhaitable de surestimer l’erreur(\(\gamma >1\) ) pour pouvoir être utilisable comme critère d’arrêt dans un processus adaptatif.

Des études systématiques ont été effectuées par Strouboulis et Haque [bib11] et par Babuška [bib12] sur différents cas-tests ayant une solution analytique ou pour lesquels une solution de référence peut être déterminée sur un maillage très fin. Globalement, il en ressort que la qualité d’un estimateur dépend de la topologie du maillage, de la régularité de la solution et de la régularité (aplatissement) des éléments.

Le seul indice global \(\gamma\) ne permet pas de rendre compte de la manière dont l’erreur locale se comporte: il se peut, que dans quelques cas, l’estimateur soit satisfaisant pour ce critère sans qu’il ne réussisse à localiser les zones où l’erreur en énergie est importante [bib13].

Indice de robustesse

Afin de s’affranchir du comportement local de l’estimateur, un autre critère a été défini par Babuška et Rodriguez [bib14]. Pour une triangulation donnée \(T\) , l’indice d’efficacité local \({\gamma}_{\omega}\) relatif au sous-domaine \(\omega\) est défini par:

(1585)\[{\gamma}_{\omega}=\frac{{e}_{\text{estimée}}^{\omega}}{{e}_{\text{vraie}}^{\omega}}\]

Babuška et Strouboulis ont montré que l’on pouvait accéder numériquement à la plage de variation asymptotique de \({\gamma}_{\omega}\) [bib15]:

(1586)\[0<{C}_{\inf}^{\omega}\le {\gamma}_{\omega}\le {C}_{\sup}^{\omega}<\infty\]

L’écart à 1 de \({\gamma}_{\omega}\) est mesuré par l’indice \({R}_{\omega}\) défini par:

(1587)\[{R}_{\omega}=\max\left\lbrace ∣1-{C}_{\sup}^{\omega}∣+∣1-{C}_{\inf}^{\omega}∣;∣1-\frac{1}{{C}_{\sup}^{\omega}}∣+∣1-\frac{1}{{C}_{\inf}^{\omega}}∣\right\rbrace\]

L’indice de robustesse \(R\) est le plus grand \({R}_{\omega}\) obtenu en variant la position de la cellule \(\omega\) de son voisinage sur les maillages appartenant à une classe de triangulation. Cette approche permet de mettre en évidence l’erreur de discrétisation sur un sous-domaine mais également la distorsion des éléments du maillage et la régularité de la solution. Bien que possédant des bases mathématiques, cette technique semble compliquée à mettre en œuvre pour des maillages complexes utilisés en ingénierie.

Erreur locale

Erreur de pollution

Définition de l’erreur de pollution

Pour définir l’erreur de pollution, nous allons fournir une définition pragmatique en essayant de comprendre quel peut être son comportement et donc son influence au cours d’un processus adaptatif. Pour cela, nous allons étudier la convergence de l’erreur en norme de l’énergie dans différentes situations:

• Convergence de l’erreur globale pour un raffinement uniforme global;

• Convergence de l’erreur dans la zone d’intérêt pour un raffinement global;

• Convergence de l’erreur globale pour un raffinement local;

• Convergence de l’erreur dans la zone d’intérêt pour un raffinement local.

Fig. 164 Géométrie et maillage (192 éléments)

Le problème à résoudre est le problème dont la géométrie et le maillage sont donnés par la figure Fig. 164. L’effort imposé est unitaire, le module d’Young est égal à \(E=1.0\) Pa et le coefficient de Poisson vaut \(\nu =0.3\) . Enfin la modélisation employée est une modélisation en contraintes planes. Une série d’approximations éléments finis sera déterminée sur une succession de maillages linéaires obtenus par raffinement uniforme de la structure ou par raffinement uniforme de la zone d’intérêt \(\omega\) . L’erreur en norme de l’énergie est alors calculée (donc sans estimation) à partir d’une solution obtenue sur un maillage très fin (solution overkill ) qui sera notre solution de référence. Enfin l’erreur absolue est normée par la norme en énergie de la solution sur la structure ou sur la zone d’intérêt selon le cas. Les résultats sont présentés dans la figure Fig. 165.

Pour un raffinement global, l’erreur globale converge «un peu moins vite» que \(p=1\) ; en divisant la taille des mailles par 10, l’erreur n’est pas divisée par 10 (même en divisant la taille par 16 l’erreur n’est pas divisée par 10) (marqueurs carrés). Toutefois cette vitesse de convergence n’est atteinte qu’asymptotiquement. Il faut noter que la courbe représentée semble se rapprocher de la convergence théorique à partir du dernier segment (pour des tailles mailles entre 1/8 et 1/16). Ce n’est donc que la première partie d’une courbe de convergence qui est tracée, quand la convergence en \(p=1\) n’est pas encore atteinte (régime pré-asymptotique).

Lorsque seule la zone d’intérêt est raffinée, la solution n’est pas améliorée dans le reste de la structure. L’erreur globale reste pratiquement inchangée (marqueurs ronds). Ceci indique que l’erreur dans la zone d’intérêt contribue très peu à l’erreur globale.

Fig. 165 Mise en évidence de l’erreur de pollution: convergence de l’erreur

Par contre en étudiant les courbes relatives à l’erreur dans la zone d’intérêt, il est visible que pour améliorer de manière significative la solution dans cette zone il ne suffit pas de raffiner uniquement cette zone. En effet les deux courbes (marqueurs losanges et marqueurs triangulaires) montrent une décroissance de l’erreur locale mais celle obtenue par le raffinement local montre une décroissance plus lente: les deux courbes s’écartent. Il semble donc que l’erreur dans le reste de la structure a aussi une influence sur la qualité de la solution dans la zone d’intérêt. Ce qui signifie que la somme des contributions locales des estimateurs globaux sur le sous domaine \(\omega\) ne représente pas exactement l’erreur locale. Il semble donc que l’erreur locale vraie est composée de deux parties: une composante locale et une composante complémentaire issue du reste de la structure. La conséquence de ceci est qu’en ne raffinant que la zone d’intérêt, l’erreur dans le reste de la structure ne diminue pas (i.e. la composante complémentaire reste constante). On peut alors imaginer que la courbe représentant l’erreur locale avec un raffinement local (marqueurs losanges) converge vers une valeur non nulle si le raffinement est poursuivi: cette valeur non nulle représentera la composante complémentaire (la composante locale tendant vers zéro). La valeur de la composante complémentaire est alors appelée «erreur de pollution».

Une définition plus formelle peut être donnée. L’erreur de pollution a été introduite par Babuška et al. [bib16], [bib17], [bib18]. Ils présument que l’erreur de discrétisation dans un sous domaine \(\omega\) de la structure est composée:

• d’une erreur locale sur \(\omega\) , faisant abstraction du reste de la structure;

• d’une erreur de pollution ou transportée, issue de la discrétisation dans le reste de la structure.

La relation fondamentale entre l’erreur globale \(e\) et la solution \({\Phi}_{E}\) du problème local de l’erreur (solution issue de l’ Equilibrated Residual Method par exemple) est rappelée:

(1588)\[\begin{array}{cc}a(e,v)=\sum_{E\in \Omega }{a}_{E}({\Phi}_{E},v)& \forall v\in V\end{array}\]

Un erreur locale \({e}_{\omega}^{\text{loc}}\) sur un sous-domaine \(\omega\) peut être définie par:

(1589)\[\begin{array}{cc}a({e}_{\omega}^{\text{loc}},v)=\sum_{E\subset \omega }{a}_{E}({\varphi}_{E},v)& \forall v\in V\end{array}\]

Comme la fonction indicatrice d’erreur locale \({\varphi}_{E}\) ne dépend que des résidus locaux de l’approximation éléments finis évaluée sur le patch d’éléments \(E\text{*}\) , l’erreur locale \({e}^{\text{loc}}\) ne dépend que des résidus locaux sur le patch \(\tilde{\omega}\) , constitué des éléments de \(\omega\) et des éléments voisins. La composante complémentaire de l’erreur vraie est l’erreur de pollution \({e}_{\omega}^{\text{pol}}\) :

(1590)\[\begin{array}{cc}a({e}_{\omega}^{\text{pol}},v)=\sum_{E\not\subset \omega }{a}_{E}({\varphi}_{E},v)& \forall v\in V\end{array}\]

Cette partie de l’erreur vraie provient des résidus situés en dehors du sous domaine et est transportée jusqu’au sous domaine, affectant ainsi la précision. L’erreur \({e}_{\omega}\) dans un sous-domaine \(\omega\) est donc la somme des deux composantes définies précédemment:

(1591)\[{e}_{\omega}={e}_{\omega}^{\text{loc}}+{e}_{\omega}^{\text{pol}}\]

Au vu de la figure Fig. 165 et de la définition précédente, il devient évident que la maîtrise de l’erreur de pollution est primordiale. En effet, en dehors de la zone d’intérêt, lorsque le maillage est grossier, il peut-être une source importante d’erreur de pollution. En conséquence, un raffinement local sur un sous domaine doit être contrebalancé par un raffinement adéquat en dehors du sous domaine afin de contrôler l’erreur de pollution. Nous verrons dans la suite que les estimateurs d’erreur en quantité d’intérêt permettent d’assurer cette équilibre entre raffinement local et raffinement en dehors du sous domaine.

Estimation de l’erreur de pollution

Pour estimer l’erreur de pollution, différentes méthodes ont été proposées. La méthode proposée par Babuška (Babuška et al. , 1995) est liée aux fonctions de Green qui décrivent l’interaction entre différents points du domaine. L’erreur de pollution est estimée avec l’erreur d’approximation de la fonction de Green associée à certains points. Une autre méthode a été utilisée par Huerta et Díez (Huerta & Díez, 2000) et consiste à trouver une approximation de \({e}_{\omega}^{\text{pol}}\) par résolution du problème définissant cette erreur. D’autres travaux peuvent être également cités (Ainsworth, 1999), (Oden & Feng, 1996). L’intérêt industriel de l’évaluation de cette erreur étant limité, nous ne développerons pas davantage ces méthodes d’estimation. Les nombreux travaux sur le calcul de bornes de l’erreur rendent inutiles l’estimation de l’erreur de pollution car elle est comprise dans les bornes. De plus, pour l’estimation d’une erreur locale, comme nous l’avons déjà indiqué, les estimateurs en quantité d’intérêt permettent de s’affranchir de l’erreur de pollution grâce au problème dual.

Erreur en quantité d’intérêt

Quantités d’intérêt

Les estimations d’erreurs en norme de l’énergie sont trop abstraites pour fournir aux utilisateurs de codes de calculs suffisamment d’informations sur les différents aspects de la solution cherchée. Il est alors plus utile d’estimer une erreur en termes de quantités d’intérêt; quantités ayant un sens physique (une moyenne de déplacements dans une sous-région du domaine physique ou une moyenne de contraintes sur une interface ou bien encore, dans le cas de vibrations, l’une des fréquences propres). Mathématiquement, elles sont caractérisées par des fonctionnelles linéaires ou non dans l’espace de fonctions auxquelles appartiennent les solutions.

Plusieurs exemples sont donnés: la moyenne d’une composante du déplacement sur un domaine \(\omega\) :

(1592)\[Q(v)=\frac{1}{\mid \omega \mid }{\int}_{\omega}{v}_{x}d\Omega\]

la moyenne d’une composante du tenseur des contraintes sur un domaine \(\omega\) :

(1593)\[Q(v)=\frac{1}{\mid \omega \mid }{\int}_{\omega}{\sigma}_{xx}d\Omega\]

On cherche alors à calculer la quantité \({\varepsilon}^{Q}\) qui peut être exprimée sous la forme suivante:

(1594)\[{\varepsilon}^{Q}=Q(u)-Q({u}^{h})\]

Dans le cas où \(Q\) est linéaire, l’égalité suivante est vérifiée:

(1595)\[{\varepsilon}^{Q}=Q(u-{u}^{h})=Q(e)\]

Dans la suite, seul le cas des fonctionnelles linéaires sera considéré car cela correspond au cadre de l’approche simplifiée présentée dans la suite. De manière générale, le traitement utilisé pour les quantités non linéaires (comme la contrainte de Von Mises ou l’intégrale \(J\) en mécanique de la rupture) est la linéarisation [bib19] et [bib20].

Problème dual et relation fondamentale

L’erreur en quantité d’intérêt peut être estimée en utilisant des propriétés de dualité des opérateurs intervenant dans la formulation faible du problème de référence. D’un point de vue mécanique, cela correspond à l’existence d’une fonction solution d’un problème dual, dite fonction d’influence, représentant l’influence des efforts du problème de référence sur l’énergie totale. La fonction d’influence est ensuite utilisée comme pondération des résidus d’équilibre du problème de référence.

Une approche simplifiée est présentée, dans la mesure où le cadre est limité à l’élasticité linéaire et aux quantités d’intérêt linéaires. Ceci se traduit par le fait que \(Q(\cdot )\) est une forme linéaire et que \(a(\cdot ,\cdot )\) est une forme bilinéaire. Cette approche permet de mieux comprendre la nécessité d’utiliser l’erreur sur la fonction d’influence et pas seulement la fonction d’influence.

La relation entre l’erreur en quantité d’intérêt et le résidu peut être exprimée par l’introduction d’une fonction inconnue \(\omega\) , appelée fonction d’influence:

(1596)\[Q(e)=\omega ({R}_{h}^{u})\]

De plus, on montre dans [bib21] que la relation précédente devient:

(1597)\[Q(e)={R}_{h}^{u}(\omega )\]

En combinant la relation précédente avec l’équation des résidus:

(1598)\[a(e,v)={R}_{h}^{u}(v)\]

une relation entre l’erreur sur le problème de référence, la fonction d’influence et l’erreur en quantité d’intérêt est obtenue:

(1599)\[a(e,\omega )=Q(e)\]

Cette égalité est nécessairement vérifiée lorsque \(\omega \in V\) est solution du problème suivant:

(1600)\[\begin{split}\lbrace \begin{array}{cc}\text{Trouver}\omega \in V\text{tel}\text{que}& \\ a(v,\omega )=Q(v)& \forall v\in V\end{array}\end{split}\]

Ce problème est appelé problème dual ou adjoint. On remarque que si \(\omega\) pouvait être calculée de manière exacte, alors \(Q(u)\) pourrait être déterminée directement à partir des données du problème car:

(1601)\[Q(u)=a(u,\omega )=l(\omega )\]

Mais malheureusement, il est aussi difficile de déterminer la fonction d’influence \(\omega\) que de déterminer \(u\) dans le problème primal. C’est pourquoi une approximation éléments finis \({\omega}^{h}\in {V}^{h}\) de \(\omega\) qui satisfait

(1602)\[\begin{array}{cc}a({v}^{h},{\omega}^{h})=Q({v}^{h})& \forall {v}^{h}\in {V}^{h}\end{array}\]

est introduite. On note également qu’en substituant \(e\) à \({v}^{h}\) dans le problème dual, on obtient:

(1603)\[Q(e)=a(e,\omega )\]

Remplacer \(\omega\) par \({\omega}^{h}\) dans la relation précédente n’est pas suffisant pour obtenir une estimation de \(Q(e)\) car \(a(e,{\omega}^{h})=0\) (propriété d’orthogonalité). C’est pourquoi il faut évaluer l’erreur sur la fonction d’influence notée \(\varepsilon =\omega -{\omega}^{h}\) . En introduisant cette erreur dans la relation précédente, une représentation de l’erreur est obtenue:

(1604)\[Q(e)=a(e,\varepsilon )\]

Dans le cas des problèmes non linéaires et/ou des quantités d’intérêt non linéaires, ce qui précède ne s’applique pas. Dans ce cas, le problème primal s’écrit:

(1605)\[\begin{split}\lbrace \begin{array}{cc}\text{Trouver}u\in V\text{tel que}& \\ A(u;v)=L(v)& \forall v\in V\end{array}\end{split}\]

où \(A(\cdot ;\cdot )\) est un forme semi-linéaire telle que les arguments qui suivent le point-virgule sont linéaires et \(L(\cdot )\) une forme linéaire continue sur \(V\) .

La démarche consiste à considérer le problème primal comme un problème de minimisation sous contrainte [bib19], [bib20]:

(1606)\[\begin{split}\lbrace \begin{array}{c}\mathrm{Trouver}u\in V\mathrm{tel}\mathrm{que}\\ Q(u)=\underset{z\in S}{\inf}Q(z)\\ \text{où}\\ S=\left\lbrace z\in V;A(z;v)=L(v),\forall v\in \in V\right\rbrace \end{array}\end{split}\]

Pour résoudre le problème de minimisation, le couple \((u,\omega )\in V\times V\) est cherché tel qu’il satisfait le problème primal et le problème dual:

(1607)\[\begin{split}\lbrace \begin{array}{cc}A(u;\tilde{v})=L(\tilde{v})& \forall \tilde{v}\in V\\ {A}^{'}(u;v,\omega )={Q}^{'}(u;v)& \forall v\in V\end{array}\end{split}\]

Finalement, la relation entre l’erreur en quantité d’intérêt \({\varepsilon}^{Q}\) et le résidu est donnée par:

(1608)\[{\epsilon}^{Q}=Q(u)-Q({u}^{h})=R({u}^{h};\omega )+\mathrm{\Delta A}+\mathrm{\Delta Q}\]

où \(\Delta A\) et \(\Delta Q\) font intervenir les dérivées secondes et troisièmes de \(A\) et \(Q\) .

Estimation de l’erreur en quantité d’intérêt

A partir de la relation fondamentale , différentes techniques ont été développées afin d’estimer l’erreur en quantité d’intérêt.

La technique la plus simple et la plus directe est d’utiliser les estimateurs globaux en norme de l’énergie. En appliquant le théorème de Cauchy-Schwartz, une relation entre l’erreur en quantité d’intérêt et les normes en énergie des erreurs sur le problème primal et le problème dual peut être trouvée [bib22]:

(1609)\[\mid Q(e)\mid \le {\parallel e\parallel }_{e}{\parallel \varepsilon \parallel }_{e}\]

Cette expression nous permet ainsi d’estimer facilement l’erreur en quantité d’intérêt en utilisant n’importe quel estimateur d’erreur en norme en énergie, ceci sans développement supplémentaire dans un code de calcul disposant d’estimateurs globaux. Bien que cette estimation est la plus simple, elle est aussi la plus grossière: elle ne prend pas en compte les interactions locales entre l’erreur primale et l’erreur duale (schématiquement cela veut dire qu’une erreur sur le déplacement ne se traduit pas par une erreur en quantité d’intérêt). Ceci peut mener à de lourdes surestimations de l’erreur en quantité d’intérêt.

Une approche un peu différente, appelée Dual Weighted Residual (DWR) [bib19], consiste à considérer l’expression exacte de l’erreur en quantité d’intérêt:

(1610)\[Q(e)=a(e,\omega -{\omega}^{h})={R}_{h}^{u}(\omega -{\omega}^{h})\]

Le problème dual est résolu en utilisant une méthode d’approximation d’ordre élevée. Une autre méthode, moins coûteuse, consiste à utiliser des fonctions d’interpolation d’ordre élevées définies sur des patchs élémentaires du domaine. Cette approche conduit à une borne supérieure garantie de l’erreur.

Une approche « bornes exactes » ( exact bounds approach ) [bib23] utilise les propriétés de l’estimation basée sur l’analyse duale. Ainsi grâce à une formulation en déplacements de la méthode des éléments finis, une borne inférieure de l’énergie de déformation exacte est obtenue; grâce à une formulation en contraintes, une borne supérieure de l’énergie complémentaire exacte est obtenue. D’un point de vue pratique, cette méthode est très coûteuse puisqu’elle implique de résoudre deux fois le problème primal et le problème dual.

D’autres méthodes existent pour estimer l’erreur en quantité d’intérêt mais d’un point de vue industriel, seule l’approche utilisant des estimateurs globaux en norme de l’énergie semble être convaincante. Elle ne nécessite pas de développement supplémentaire pour accéder à une indication de l’erreur et elle est aussi beaucoup moins coûteuse.

Bornes de l’erreur

Une nouvelle expression pour l’erreur

Bien que très simple à obtenir, la relation produit une surestimation conséquente de l’erreur. Différentes techniques ont été proposées pour remédier à ce problème. Prudhomme et Oden emploient la relation du parallélogramme pour les problèmes avec des formes bilinéaires symétriques donnant accès à une borne supérieure et une borne inférieure de l’erreur [bib24], [bib21]].

Si \(e\) et \(\varepsilon\) sont les erreurs des approximations éléments finis, respectivement des solutions des problèmes primal et dual alors:

(1611)\[Q(e)=a(e,\varepsilon )=\frac{1}{4}{\parallel e+\varepsilon \parallel }_{e}^{2}-\frac{1}{4}{\parallel e-\varepsilon \parallel }_{e}^{2}\]

On introduit le facteur scalaire \(s\in R\)

(1612)\[Q(e)=a(e,\varepsilon )=\frac{1}{4}{\parallel se+{s}^{-1}\varepsilon \parallel }_{e}^{2}-\frac{1}{4}{\parallel se-{s}^{-1}\varepsilon \parallel }_{e}^{2}\]

Le scalaire \(s\) est choisi de manière à ce que les quantités \({\parallel se\parallel }_{e}\) et \({\parallel {s}^{-1}\varepsilon \parallel }_{e}\) aient la même amplitude, ce qui implique de choisir \(s=\sqrt{{\parallel \varepsilon \parallel }_{e}/{\parallel e\parallel }_{e}}\) .

En effet \(e\) et \(\varepsilon\) peuvent être d’ordre de grandeur très différent. Le scalaire \(s\) vise à «normer» les deux termes, afin d’éviter de faire une différence de deux termes très voisins car cela peut induire une erreur numérique importante.

Formellement, n’importe quel estimateur d’erreurs existant peut être utilisé pour évaluer les normes en énergie de la relation du parallélogramme; plus cette estimation sera «juste» meilleure sera la qualité de l’estimateur d’erreur en quantité d’intérêt.

Construction des bornes

Les estimateurs globaux \({\eta}_{\inf}^{+}\) , \({\eta}_{\sup}^{+}\) , \({\eta}_{\inf}^{-}\) et \({\eta}_{\sup}^{-}\) sont définis tels que:

(1613)\[{\eta}_{\inf}^{+}\le {\parallel se+{s}^{-1}\varepsilon \parallel }_{e}\le {\eta}_{\sup}^{+}\]

(1614)\[{\eta}_{\inf}^{-}\le {\parallel se-{s}^{-1}\varepsilon \parallel }_{e}\le {\eta}_{\sup}^{-}\]

En réécrivant les deux inégalités précédentes ainsi:

(1615)\[\frac{1}{4}{({\eta}_{\inf}^{+})}^{2}\le \frac{1}{4}{\parallel se+{s}^{-1}\varepsilon \parallel }_{e}^{2}\le \frac{1}{4}{({\eta}_{\sup}^{+})}^{2}\]

(1616)\[\frac{1}{4}{({\eta}_{\inf}^{-})}^{2}\le \frac{1}{4}{\parallel se-{s}^{-1}\varepsilon \parallel }_{e}^{2}\le \frac{1}{4}{({\eta}_{\sup}^{-})}^{2}\]

et en les additionnant:

(1617)\[\frac{1}{4}{({\eta}_{\inf}^{+})}^{2}-\frac{1}{4}{({\eta}_{\sup}^{-})}^{2}\le \frac{1}{4}{\parallel se+{s}^{-1}\varepsilon \parallel }_{e}^{2}-\frac{1}{4}{\parallel se-{s}^{-1}\varepsilon \parallel }_{e}^{2}\le \frac{1}{4}{({\eta}_{\sup}^{+})}^{2}-\frac{1}{4}{({\eta}_{\inf}^{-})}^{2}\]

un encadrement de l’erreur en quantité d’intérêt est obtenu.

Ensuite les bornes \({\eta}_{\inf}^{Q}\) , \({\eta}_{\sup}^{Q}\) sont introduites telles que:

(1618)\[{\eta}_{\inf}^{Q}=\frac{1}{4}{({\eta}_{\inf}^{+})}^{2}-\frac{1}{4}{({\eta}_{\sup}^{-})}^{2}\]

(1619)\[{\eta}_{\sup}^{Q}=\frac{1}{4}{({\eta}_{\sup}^{+})}^{2}-\frac{1}{4}{({\eta}_{\inf}^{-})}^{2}\]

L’erreur en quantité d’intérêt \(Q(e)\) est bornée par \({\eta}_{\inf}^{Q}\) et \({\eta}_{\sup}^{Q}\) :

(1620)\[{\eta}_{\inf}^{Q}\le Q(e)\le {\eta}_{\sup}^{Q}\]

Il est également possible d’utiliser les estimations \({\eta}_{\text{eei}}^{Q}\) et \({\eta}_{\mathrm{ees}}^{Q}\) définies telles que:

(1621)\[{\eta}_{\text{eei}}^{Q}=\frac{1}{4}{({\eta}_{\inf}^{+})}^{2}-\frac{1}{4}{({\eta}_{\inf}^{-})}^{2}\]

(1622)\[{\eta}_{\text{ees}}^{Q}=\frac{1}{4}{({\eta}_{\sup}^{+})}^{2}-\frac{1}{4}{({\eta}_{\sup}^{-})}^{2}\]

comme estimation de l’erreur en quantité d’intérêt \(Q(e)\) .

Enfin, une dernière estimation peut être obtenue, en moyennant les bornes \({\eta}_{\inf}^{Q}\) et \({\eta}_{\sup}^{Q}\) ou les estimations \({\eta}_{\mathrm{eei}}^{Q}\) et \({\eta}_{\mathrm{ees}}^{Q}\) , notée \({\eta}_{\mathrm{moy}}^{Q}\) :

(1623)\[\begin{split}\begin{array}{ccc}{\eta}_{\mathrm{moy}}^{Q}& =& \frac{1}{2}({\eta}_{\mathrm{eei}}^{Q}+{\eta}_{\mathrm{ees}}^{Q})\\ & =& \frac{1}{2}({\eta}_{\inf}^{Q}+{\eta}_{\sup}^{Q})\\ & =& \frac{1}{8}({({\eta}_{\inf}^{+})}^{2}+{({\eta}_{\sup}^{+})}^{2})-\frac{1}{8}({({\eta}_{\inf}^{-})}^{2}+{({\eta}_{\sup}^{-})}^{2})\end{array}\end{split}\]

La figure suivante résume les différentes combinaisons possibles de l’estimation d’erreur en norme de l’énergie pour obtenir une estimation ou des bornes de l’erreur en quantité d’intérêt. Ceci montre qui si l’on dispose de bornes pour un estimateur en norme de l’énergie alors on dispose de bornes pour l’erreur en quantité d’intérêt.

Fig. 166 Estimations et bornes de l’erreur en quantité d’intérêt (Prudhomme et al. , 2003)

Estimations des bornes de l’erreur

Dans le chapitre précédent, il a été montré que des bornes de l’erreur globale étaient accessibles en utilisant des estimateurs en résidus implicites. Les grandes lignes sont données ici pour obtenir des bornes de l’erreur en quantité d’intérêt en utilisant les estimateurs implicites. Il faut noter que ceci n’est qu’un exemple et que d’autres méthodes d’estimation peuvent fournir des bornes de l’erreur en quantité d’intérêt.

Dans le but d’obtenir les bornes supérieures précédemment définies, les fonctions \({\psi}_{E}^{u}\) et \({\psi}_{E}^{\omega}\) , correspondant respectivement aux erreurs \(e\) et \(\varepsilon\) , sont calculées en introduisant les problèmes locaux relatifs aux problèmes primal et dual. Les estimations suivantes sont alors obtenues:

(1624)\[{\eta}_{\sup}^{+}={\parallel s{\psi}_{E}^{u}+{s}^{-1}{\psi}_{E}^{\omega}\parallel }_{e}\]

Pour l’estimation des bornes inférieures, il faut procéder de manière analogue que pour l’estimation des bornes supérieures en construisant des fonctions continues aux interfaces à partir de \({\psi}_{E}^{u}\) et \({\psi}_{E}^{\omega}\) :

(1625)\[{\eta}_{\inf}^{+}=\frac{∣(s{R}_{h}^{u}+{s}^{-1}{R}_{h}^{\omega})(s{\chi }^{u}+{s}^{-1}{\chi }^{\omega})∣}{{\parallel s{\chi }_{E}^{u}+{s}^{-1}{\chi }_{E}^{\omega}\parallel }_{e}}\]

Estimation de l’erreur en quantité d’intérêt

Estimateur d’erreur en quantité d’intérêt basé sur les résidus explicites

L’erreur en quantité d’intérêt peut être estimée à partir de n’importe quel estimateur de l’erreur globale. Ainsi, un estimateur d’erreur en quantité d’intérêt est construit à partir de la relation du parallélogramme et de l’estimateur d’erreur basé sur les résidus explicites \(\eta\) dont l’expression sur chaque élément \({\eta}_{E}\) est [R4.10.02]:

(1626)\[{\eta}_{E}^{2}={h}_{E}^{2}{\parallel {r}_{E}\parallel }_{{L}^{2}(E)}^{2}+\frac{1}{2}\sum_{\Gamma \not\subset {\Gamma}_{F}}{l}_{\Gamma}{\parallel {t}_{\Gamma}\parallel }_{{L}^{2}(\Gamma )}^{2}+\sum_{\Gamma \subset {\Gamma}_{F}}{l}_{\Gamma}{\parallel {t}_{\Gamma}\parallel }_{{L}^{2}(\Gamma )}^{2}\]

L’utilisation de l’expression donnée par l’inégalité de Cauchy-Schwartz:

(1627)\[\mid Q(e)\mid \le {\parallel e\parallel }_{e}{\parallel \varepsilon \parallel }_{e}\]

avec l’estimateur de type résidu explicite donne facilement une estimation de l’erreur en quantité d’intérêt en estimant \({\parallel e\parallel }_{e}\) et \({\parallel \varepsilon \parallel }_{e}\) mais cette estimation est grossière.

La relation du parallélogramme donnée par:

(1628)\[Q(e)=a(e,\varepsilon )=\frac{1}{4}{\parallel e+\varepsilon \parallel }_{e}^{2}-\frac{1}{4}{\parallel e-\varepsilon \parallel }_{e}^{2}=\frac{1}{4}{\parallel se+{s}^{-1}\varepsilon \parallel }_{e}^{2}-\frac{1}{4}{\parallel se-{s}^{-1}\varepsilon \parallel }_{e}^{2}\]

avec \(s=\sqrt{{\parallel \varepsilon \parallel }_{e}/{\parallel e\parallel }_{e}}\) (cf. §3.2.4), fournit une estimation plus précise car elle permet un réel couplage entre l’erreur primale et l’erreur duale. L’utilisation de l’estimation de l’erreur de chacun des deux problèmes n’est pas possible, la norme d’une somme n’est pas égale à la somme des normes:

(1629)\[{\parallel e+\varepsilon \parallel }_{e}^{2}=a(e+\varepsilon ,e+\varepsilon )={\parallel e\parallel }_{e}^{2}+{\parallel \varepsilon \parallel }_{e}^{2}+2a(e,\varepsilon )\]

(1630)\[{\parallel e-\varepsilon \parallel }_{e}^{2}=a(e-\varepsilon ,e-\varepsilon )={\parallel e\parallel }_{e}^{2}+{\parallel \varepsilon \parallel }_{e}^{2}-2a(e,\varepsilon )\]

C’est pourquoi il faut repartir de la définition de l’estimateur d’erreur en résidu afin d’évaluer directement \({\parallel se+{s}^{-1}\varepsilon \parallel }_{e}^{2}\) et \({\parallel se-{s}^{-1}\varepsilon \parallel }_{e}^{2}\) . L’erreur en solution \(e=u-{u}^{h}\) étant solution d’un problème d’élasticité dont les données en effort sont les résidus d’équilibre, il s’agit donc de formuler les problèmes dont \((se+{s}^{-1}\epsilon )\) et \((se-{s}^{-1}\varepsilon )\) sont solutions.

Le problème primal et le problème dual sont donnés par les relations suivantes:

(1631)\[\begin{split}\lbrace \begin{array}{cc}a(u,v)=l(v)& \forall v\in V\\ a(v,\varepsilon )=Q(v)& \forall v\in V\end{array}\end{split}\]

Les résidus associés à ces problèmes sont respectivement donnés par:

(1632)\[\begin{split}\lbrace \begin{array}{c}{R}_{h}^{u}(v)=l(v)-a({u}^{h},v)\\ {R}_{h}^{\omega}(v)=Q(v)-a(v,{\omega}^{h})\end{array}\end{split}\]

En remplaçant \(u\) par \({u}^{h}+e\) et \(\omega\) par \({\omega}^{h}+\varepsilon\) et en utilisant la symétrie de \(a(\cdot ,\cdot )\) , on obtient les équations gouvernant les erreurs \(e\) et \(\varepsilon\) :

(1633)\[\begin{split}\lbrace \begin{array}{c}a(e,v)={R}_{h}^{u}(v)\\ a(\varepsilon ,v)={R}_{h}^{\omega}(v)\end{array}\end{split}\]

Par combinaison linéaire, on a les deux problèmes suivants:

(1634)\[\begin{split}\lbrace \begin{array}{c}a(se+{s}^{-1}\varepsilon ,v)=(s{R}_{h}^{u}+{s}^{-1}{R}_{h}^{\omega})(v)\\ a(se-{s}^{-1}\varepsilon ,v)=(s{R}_{h}^{u}-{s}^{-1}{R}_{h}^{\omega})(v)\end{array}\end{split}\]

Ainsi l’erreur numérique \((se+{s}^{-1}\epsilon )\) (et \((se-{s}^{-1}\varepsilon )\) ) est générée par le résidu total \((s{R}_{h}^{u}+{s}^{-1}{R}_{h}^{\omega})\) (respectivement \((s{R}_{h}^{u}-{s}^{-1}{R}_{h}^{\omega})\) ). De manière analogue à l’estimateur en résidu, il existe une constante \({C}^{+}\) (et \({C}^{-}\) ), indépendante de la taille des éléments, telle que:

(1635)\[{\parallel se+{s}^{-1}\varepsilon \parallel }_{e}^{2}\le {C}^{+}\left[\sum_{E}{h}_{E}^{2}{\parallel s{r}_{E}^{u}+{s}^{-1}{r}_{E}^{\omega}\parallel }_{{L}^{2}(E)}^{2}+\sum_{\Gamma \not\subset {\Gamma}_{U}}{l}_{\Gamma}{\parallel s{t}_{\Gamma}^{u}+{s}^{-1}{t}_{\Gamma}^{\omega}\parallel }_{{L}^{2}(\Gamma )}^{2}\right]\]

(1636)\[{\parallel se-{s}^{-1}\varepsilon \parallel }_{e}^{2}\le {C}^{-}\left[\sum_{E}{h}_{E}^{2}{\parallel s{r}_{E}^{u}-{s}^{-1}{r}_{E}^{\omega}\parallel }_{{L}^{2}(E)}^{2}+\sum_{\Gamma \not\subset {\Gamma}_{U}}{l}_{\Gamma}{\parallel s{t}_{\Gamma}^{u}-{s}^{-1}{t}_{\Gamma}^{\omega}\parallel }_{{L}^{2}(\Gamma )}^{2}\right]\]

Enfin l’estimateur d’erreur en quantité d’intérêt \({}^{\mathrm{REQI}}\eta _{E}\) est défini tel que:

(1637)\[{}^{\mathrm{REQI}}\eta _{E}=\frac{1}{4}{({\eta}_{E}^{+})}^{2}-\frac{1}{4}{({\eta}_{E}^{-})}^{2}\]

avec \({\eta}_{E}^{+}\) et \({\eta}_{E}^{-}\) définit par:

(1638)\[({\eta}_{E}^{+}{)}^{2}={c}^{+}\left[{h}_{E}^{2}{\parallel s{r}_{E}^{u}+{s}^{-1}{r}_{E}^{\omega}\parallel }_{{L}^{2}(E)}^{2}+\frac{1}{2}\sum_{\Gamma \not\subset {\Gamma}_{F}}{l}_{\Gamma}{\parallel s{t}_{\Gamma}^{u}+{s}^{-1}{t}_{\Gamma}^{\omega}\parallel }_{{L}^{2}(\Gamma )}^{2}+\sum_{\Gamma \subset {\Gamma}_{F}}{l}_{\Gamma}{\parallel s{t}_{\Gamma}^{u}+{s}^{-1}{t}_{\Gamma}^{\omega}\parallel }_{{L}^{2}(\Gamma )}^{2}\right]\]

(1639)\[({\eta}_{E}^{-}{)}^{2}={c}^{-}\left[{h}_{E}^{2}{\parallel s{r}_{E}^{u}-{s}^{-1}{r}_{E}^{\omega}\parallel }_{{L}^{2}(E)}^{2}+\frac{1}{2}\sum_{\Gamma \not\subset {\Gamma}_{F}}{l}_{\Gamma}{\parallel s{t}_{\Gamma}^{u}-{s}^{-1}{t}_{\Gamma}^{\omega}\parallel }_{{L}^{2}(\Gamma )}^{2}+\sum_{\Gamma \subset {\Gamma}_{F}}{l}_{\Gamma}{\parallel s{t}_{\Gamma}^{u}-{s}^{-1}{t}_{\Gamma}^{\omega}\parallel }_{{L}^{2}(\Gamma )}^{2}\right]\]

Compte tenu des définitions des bornes de l’erreur en quantité d’intérêt données dans la partie 1, l’estimateur que nous venons de définir ne borne pas l’erreur mais en donne une valeur approchée:

(1640)\[Q(e)\approx \sum_{E}{}^{\mathrm{REQI}}\eta _{E}\]

L’estimateur ainsi défini s’apparente donc davantage à un indicateur d’erreur permettant de guider une stratégie de maillage. Cet indicateur est programmé dans Code_Aster , il est disponible pour tous les types d’éléments 2D et 3D (sauf les pyramides) ainsi que pour les éléments linéaires et quadratiques.

Estimateur en quantité d’intérêt basé sur le lissages de champ

Un estimateur d’erreur est construit directement à partir des techniques de lissage [R4.10.01] et de la relation fondamentale entre l’erreur en quantité d’intérêt, l’erreur sur la solution du problème primal et celle sur la solution du problème dual:

(1641)\[Q(e)=a(e,\varepsilon )={\int}_{\Omega}\mathrm{Tr}\left[{}^{T}(\sigma (u)-\sigma ({u}^{h})){K}^{-1}(\sigma (\omega )-\sigma ({\omega}^{h}))\right]d\Omega\]

La démarche consiste à construire à partir des solutions éléments finis \({u}^{h}\) et \({\omega}^{h}\) des champs de contraintes continues \(\tilde{\sigma}({u}^{h})\) et \(\tilde{\sigma}({\omega}^{h})\) le plus proche possible des solutions exactes. Cette opération de lissage est réalisée par les différentes techniques décrites dans la première partie. Une borne supérieure est ainsi obtenue [bib9]:

(1642)\[a(e,\varepsilon )\le C{\int}_{\Omega}\mathrm{Tr}\left[(\tilde{\varepsilon}(u)-\varepsilon ({u}^{h}))K(\tilde{\varepsilon}(\omega )-\varepsilon ({\omega}^{h}))\right]d\Omega\]

Les contributions locales de d’erreur en quantité d’intérêt peuvent être définies ainsi:

(1643)\[{}^{\mathrm{ZZQI}}\eta _{E}={\int}_{E}\mathrm{Tr}\left[(\tilde{\varepsilon}(u)-\varepsilon ({u}^{h}))K(\tilde{\varepsilon}(\omega )-\varepsilon ({\omega}^{h}))\right]dE\]

L’estimateur est la somme sur les éléments des contributions et s’écrit:

(1644)\[{}^{\mathrm{ZZQI}}\eta =\sum_{E}{}^{\mathrm{ZZQI}}\eta _{E}\]

Cet estimateur est très simple à mettre en œuvre et ne nécessite pratiquement pas de développement pour un code pourvu de techniques de lissage. Cet indicateur est programmé dans Code_Aster , il est disponible pour tous les types d’éléments 2D.

Estimateur en quantité d’intérêt basé sur les bases hiérarchiques

Le point de départ pour construire cet estimateur est encore la relation fondamentale de l’erreur en quantité d’intérêt:

(1645)\[Q(e)=a(e,\varepsilon )\]

La résolution du problème primal et du problème dual avec des maillages plus «riches» (pratiquement ce sont des maillages raffinés uniformément) permet d’obtenir la solution améliorée \(u\text{*}\) du problème primal et la solution améliorée du problème dual \(\omega \text{*}\) .

En considérant que les approximations \(u\text{*}\) et \(\omega \text{*}\) sont suffisamment précises par rapport aux approximations originales \({u}^{h}\) et \({\omega}^{h}\) alors l’erreur en quantité d’intérêt est donnée par:

(1646)\[Q(e)=a(e,\varepsilon )=a(u-{u}^{h},\omega -{\omega}^{h})\approx a(u\text{*}-{u}^{h},\omega \text{*}-{\omega}^{h})=Q(e\text{*})\]

Un estimateur local d’erreur en quantité d’intérêt peut être défini ainsi:

(1647)\[{}^{\mathrm{BHQI}}\eta _{E}={\int}_{E}\mathrm{Tr}\left[(\varepsilon (u\text{*})-\varepsilon ({u}^{h}))K(\varepsilon (\omega \text{*})-\varepsilon ({\omega}^{h}))\right]dE\]

L’estimateur global est la somme sur les éléments et s’écrit ainsi:

(1648)\[{}^{\mathrm{BHQI}}\eta =\sum_{E}{}^{\mathrm{BHQI}}\eta _{E}\]

Cet estimateur peut être d’un intérêt limité pour une utilisation industrielle compte tenu du coût de calcul engendré. En effet le problème primal et le problème dual devront être chacun résolus deux fois: une fois sur le maillage initial et une fois sur un maillage plus fin donc avec un surcoût de calcul. Toutefois cet estimateur peut rendre l’erreur en quantité d’intérêt disponible pour des classes de problèmes pour lesquels aucun estimateur n’est disponible ou implémenté (problèmes aux valeurs propres ou de dynamique dans Code_Aster par exemple).

Formulation de la quantité d’intérêt

Comme nous l’avons déjà vu dans la première partie, au lieu d’évaluer une mesure de l’erreur telle que la norme en énergie par exemple, il est plus utile d’estimer l’erreur sur une quantité d’intérêt ayant un sens physique. Cette quantité d’intérêt ( quantity of interest ou output ) peut être représentée par une fonctionnelle linéaire \(Q(\cdot )\) définie sur l’espace des fonctions tests. Le but de l’estimation en quantité d’intérêt est donc de s’assurer de la qualité de \(Q({u}^{h})\) en estimant la quantité:

(1649)\[{\varepsilon}^{Q}=Q(u-{u}^{h})=Q(e)\]

Pour obtenir l’estimation, il faut utiliser la fonctionnelle \(Q(\cdot )\) comme le chargement du problème dual défini par:

(1650)\[\begin{split}\lbrace \begin{array}{cc}\text{Trouver}\omega \in V\text{tel que}& \\ a(v,\omega )=Q(v)& \forall v\in V\end{array}\end{split}\]

De nombreuses quantités physiques sont utilisées en conception et doivent être contrôlées. Ces grandeurs peuvent être, par exemples, la moyenne d’une composante du déplacement ou la moyenne d’une composante des contraintes sur un sous-domaine (zone d’intérêt):

(1651)\[Q(v)=\frac{1}{\mid \omega \mid }{\int}_{\omega}{v}_{x}d\Omega\]

(1652)\[Q(v)=\frac{1}{\mid \omega \mid }{\int}_{\omega}{\sigma}_{xx}d\Omega\]

Ces quantités ne sont pas utilisables directement comme chargement du problème dual. Il faut les exprimer de manière adéquate afin de retrouver le chargement à imposer sur le problème local. La forme générale de ces quantités est la suivante [bib25]:

(1653)\[Q(v)={\int}_{\Omega}{f}^{Q}\cdot vd\Omega +{\int}_{{\Gamma}_{F}}{F}^{Q}\cdot v\mathrm{d\Gamma }-a({u}^{Q},v)\]

où \({f}^{Q}\) , \({F}^{Q}\) et \({u}^{Q}\) sont données. Dans la suite, nous verrons que des quantités d’intérêt simples se déduisent directement de cette forme générale pour en constituer des cas particuliers dans lesquels un ou plusieurs termes sont nuls.

Grandeurs associées aux déplacements

Composante du déplacement dans un domaine

La première quantité liée au déplacement est la valeur moyenne d’une composante du déplacement dans un sous‑domaine \(\omega\) :

(1654)\[\begin{array}{cc}{Q}^{\text{depl}}(v)=\frac{1}{\mid \omega \mid }{\int}_{\omega}{v}_{i}d\Omega & i=1,2\text{ou}3\end{array}\]

\({Q}^{\text{depl}}(v)\) peut s’écrire:

(1655)\[{\int}_{\Omega}{f}^{Q}\cdot vd\Omega\]

qui correspond au cas particulier où \({F}^{Q}=0\) et \({u}^{Q}=0\) et où \({f}^{Q}\) est un vecteur constant sur \(\omega\) et nul ailleurs.

Le problème dual à résoudre s’écrit de la manière suivante:

(1656)\[\begin{array}{cc}{\int}_{\Omega}{\sigma}_{ij}(v){\varepsilon}_{ij}(\omega )d\Omega ={\int}_{\Omega}{f}_{i}^{Q}{v}_{i}d\Omega & \forall v\in V\end{array}\]

En utilisant la propriété de symétrie du tenseur des contraintes et en dérivant, on obtient:

(1657)\[{\int}_{\Omega}{({\sigma}_{ij}(\omega ){v}_{i})}_{,j}d\Omega -{\int}_{\Omega}{\sigma}_{ij,j}(\omega ){v}_{i}d\Omega ={\int}_{\Omega}{f}_{i}^{Q}{v}_{i}d\Omega\]

Enfin en utilisant la formule de Stokes, on a:

(1658)\[{\int}_{\Omega}({\sigma}_{ij,j}(\omega )+{f}^{Q}){v}_{i}d\Omega -{\int}_{\partial \Omega }{\sigma}_{ij}(\omega ){n}_{j}{v}_{i}d\Omega =0\]

Ceci nous permet d’écrire le problème dual local:

(1659)\[\begin{split}\lbrace \begin{array}{cccc}{\sigma}_{ij,j}(\omega )+{f}_{i}^{Q}& =& 0& \mathrm{dans}\Omega \\ {u}_{i}(x)& =& 0& \mathrm{sur}{\Gamma}_{U}\\ {\sigma}_{ij}(\omega ){n}_{j}(x)& =& 0& \mathrm{sur}{\Gamma}_{F}\end{array}\end{split}\]

On note ainsi qu’imposer \({Q}^{\mathrm{depl}}(v)\) comme chargement du problème dual revient à imposer un chargement volumique constant sur \(\omega\) . Il faut noter que si on recherche l’erreur sur le vecteur déplacement complet, il faudra résoudre autant de problème dual que de composantes (2 en dimension 2 et 3 en dimension 3).

Composante du déplacement sur un bord

Une autre quantité utile liée au déplacement est la valeur moyenne d’une composante du déplacement sur un bord \(\gamma \subset {\Gamma}_{F}\) :

(1660)\[\begin{array}{cc}{Q}^{\mathrm{depl}}(v)=\frac{1}{\mid \gamma \mid }{\int}_{\gamma}{v}_{i}d\Gamma & i=1,2\mathrm{ou}3\end{array}\]

\({Q}^{\mathrm{depl}}(v)\) peut s’écrire:

(1661)\[{\int}_{{\Gamma}_{F}}{F}^{Q}\cdot vd\Gamma\]

qui correspond au cas particulier où \({f}^{Q}=0\) et \({u}^{Q}=0\) et où \({F}^{Q}\) est un vecteur constant sur \(\gamma\) et nul ailleurs et avec pour chaque \({v}_{i}\) , \({F}_{i}^{Q}=\frac{1}{\mid \gamma \mid }\) ; les autres composantes de \({F}^{Q}\) étant nulles.

Le problème dual local s’écrit de la manière suivante:

(1662)\[\begin{split}\lbrace \begin{array}{cccc}{\sigma}_{ij,j}(\omega )& =& 0& \mathrm{dans}\Omega \\ {u}_{i}(x)& =& 0& \mathrm{sur}{\Gamma}_{U}\\ {\sigma}_{ij}(\omega ){n}_{j}(x)& =& {F}_{i}^{Q}& \mathrm{sur}{\Gamma}_{F}\end{array}\end{split}\]

Imposer \({Q}^{\mathrm{depl}}(v)\) comme chargement du problème dual revient à imposer un chargement linéique en dimension 2 et surfacique en dimension 3 restreint à \({\Gamma}_{F}\) .

Déplacement normal à un bord

Enfin la dernière quantité liée au déplacement est la moyenne du déplacement normal à un bord \(\gamma \subset {\Gamma}_{F}\) :

(1663)\[{Q}^{\mathrm{depl}}(v)=\frac{1}{\mid \gamma \mid }{\int}_{\gamma}{v}_{n}d\Gamma\]

\({Q}^{\mathrm{depl}}(v)\) peut s’écrire:

(1664)\[{\int}_{{\Gamma}_{F}}\mathrm{pn}\cdot vd\Gamma\]

où \(n\) est le vecteur normal au bord \(\gamma\) et \(p=\frac{1}{\mid \gamma \mid }\) sur \(\gamma\) et nul ailleurs.

Le problème dual local s’écrit de la manière suivante:

(1665)\[\begin{split}\lbrace \begin{array}{cccc}{\sigma}_{ij,j}(\omega )& =& 0& \mathrm{dans}\Omega \\ {u}_{i}(x)& =& 0& \mathrm{sur}{\Gamma}_{U}\\ {\sigma}_{ij}(\omega ){n}_{j}(x)& =& {\mathrm{pn}}_{i}& \mathrm{sur}{\Gamma}_{F}\end{array}\end{split}\]

Imposer \({Q}^{\mathrm{depl}}(v)\) comme chargement du problème dual revient à imposer un chargement de pression linéique en dimension 2 et surfacique en dimension 3.

Grandeurs associées aux contraintes

Composante des contraintes dans un domaine

Nous nous intéressons maintenant à l’estimation de l’erreur sur la valeur moyenne d’une composante des contraintes dans un sous‑domaine \(\omega\) :

(1666)\[\begin{array}{cc}{Q}^{sigma}(v)=\frac{1}{\mid \omega \mid }{\int}_{\omega}{\sigma}_{ij}d\Omega & i,j=1,2\mathrm{ou}3\end{array}\]

Pour déterminer \({Q}^{sigma}(v)\) , il suffit de déterminer:

(1667)\[-a({u}^{Q},v)={\int}_{\Omega}\sigma (v):\Sigma ({u}^{Q})d\Omega\]

qui correspond au cas particulier où \({f}^{Q}=0\) et \({F}^{Q}=0\) et où \(\Sigma ({u}^{Q})\) est un opérateur symétrique constant sur \(\omega\) et nul ailleurs.

Pour cela, il suffit d’écrire que la composante des contraintes est égale à la trace du produit entre des opérateurs dont l’un est non nul dans \(\omega\) :

(1668)\[\begin{array}{cc}{Q}^{sigma}(v)=\frac{1}{\mid \omega \mid }{\int}_{\omega}{\sigma}_{ij}d\Omega ={\int}_{\omega}\sigma :\Sigma d\Omega +\underset{=0}{\underset{\underbrace{}}{{\int}_{\Omega \setminus \omega }\sigma :\Sigma d\Omega }}={\int}_{\Omega}\sigma :\Sigma d\Omega & i,j=1,2\mathrm{ou}3\end{array}\]

avec pour chaque \({\sigma}_{ij}\) , \({\Sigma}_{ij}=\frac{1}{\mid \omega \mid }\) pour \(i=j\) ; \({\Sigma}_{ij}=\frac{1}{2\mid \omega \mid }\) pour \(i\ne j\) .

Le problème dual à résoudre s’écrit de la manière suivante:

(1669)\[{\int}_{\Omega}{\sigma}_{i}(v){\varepsilon}_{ij}(\omega )d\Omega ={\int}_{\Omega}{\sigma}_{ij}(v){\Sigma}_{ij}d\Omega\]

En dérivant, on obtient:

(1670)\[{\int}_{\Omega}{(({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl}){v}_{i})}_{,j}d\Omega -{\int}_{\Omega}{({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl})}_{,j}{v}_{i}d\Omega =0\]

Enfin en utilisant la formule de Stokes, on a:

(1671)\[{\int}_{\Omega}{(({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl}){v}_{i})}_{,j}d\Omega -{\int}_{\Omega}({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl}){n}_{j}{v}_{i}d\Omega =0\]

Ceci nous permet d’écrire le problème dual local:

(1672)\[\begin{split}\lbrace \begin{array}{cccc}{({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl})}_{,j}& =& 0& \mathrm{dans}\Omega \\ {u}_{i}(x)& =& 0& \mathrm{sur}{\Gamma}_{U}\\ ({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl}){n}_{j}(x)& =& 0& \mathrm{sur}{\Gamma}_{F}\end{array}\end{split}\]

On note ainsi qu’imposer \({Q}^{sigma}(v)\) comme chargement du problème dual revient à imposer une déformation initiale \(\Sigma\) . Il faut noter que si on recherche l’erreur sur le tenseur des contraintes complet, il faudra résoudre autant de problème dual que de composantes (3 en dimension 2 et 6 en dimension 3).

Contrainte normale sur un bord

Il peut être intéressant de considérer l’erreur sur la contrainte normale sur un bord \(\gamma\) :

(1673)\[\begin{array}{c}{Q}^{sigma}(v)={\int}_{\gamma}n\cdot \sigma (v)nd\Omega \end{array}\]

Cette quantité ne se présente pas sous la forme générale. L’introduction d’une fonction auxiliaire est nécessaire [bib25]. Cette fonction \(\chi\) est définie telle que \(\chi {\mid}_{\gamma}=1\) et égale à zéro sur le reste de \(\partial \Omega\) .

Si \({n}^{\gamma}=n{\mid}_{\gamma}\) , la fonction \(\chi\) vérifie l’équation suivante:

(1674)\[\begin{split}\begin{array}{cc}a(u,\chi {n}^{\gamma})& ={\int}_{\Omega}{\sigma}_{ij}(u){\varepsilon}_{ij}(\chi {n}^{\gamma})d\Omega \\ & ={\int}_{\partial \Omega }{\sigma}_{ij}(u){n}_{j}{(\chi {n}^{\gamma})}_{i}d\Omega -\underset{=0}{\underset{\underbrace{}}{{\int}_{\Omega}{({\sigma}_{ij}(u))}_{,j}{(\chi {n}^{\gamma})}_{i}d\Omega }}\\ & ={\int}_{\gamma}{\sigma}_{ij}(u){n}_{j}{(\chi {n}^{\gamma})}_{i}d\Omega \end{array}\end{split}\]

Ce qui montre que la quantité d’intérêt peut s’écrire de manière équivalente:

(1675)\[{Q}^{sigma}(u)={\int}_{\gamma}n\cdot \sigma (u)nd\Omega =a(u,\chi {n}^{\gamma})=:\tilde{Q}(u)\]

Mais comme l’égalité précédente n’est vraie que pour un \(u\) particulier, il faut noter que:

(1676)\[{Q}^{sigma}(v)={\int}_{\gamma}n\cdot \sigma (v)nd\Omega \ne a(v,\chi {n}^{\gamma})=:\tilde{Q}(v)\]

Pour plus de facilité, la quantité \(\tilde{Q}(v)\) sera utilisée à la place de la quantité \({Q}^{sigma}(v)\) . En effet elle correspond à utiliser \({u}^{Q}=-\chi {n}^{\gamma}\) .

La quantité d’intérêt à utiliser s’écrit donc de la manière suivante:

(1677)\[\tilde{Q}(v)=-a({u}^{Q},v)={\int}_{\Omega}\sigma (v):\Sigma (-\chi {n}^{\gamma})d\Omega\]

qui correspond au cas particulier où \({f}^{Q}=0\) et \({F}^{Q}=0\) .

De manière analogue au cas de la moyenne des contraintes, le problème dual local s’écrit:

(1678)\[\begin{split}\lbrace \begin{array}{cccc}{({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl})}_{,j}& =& 0& \mathrm{dans}\Omega \\ {u}_{i}(x)& =& 0& \mathrm{sur}{\Gamma}_{U}\\ ({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl}){n}_{j}(x)& =& 0& \mathrm{sur}{\Gamma}_{F}\end{array}\end{split}\]

Imposer \(\tilde{Q}(v)\) comme chargement du problème dual revient à imposer un champ de déformation initiale calculée à partir du champ de déplacement \({u}^{Q}=-\chi {n}^{\gamma}\) . Cette quantité d’intérêt n’a pas été implémentée dans Code_Aster .

Contrainte de Von Mises dans un domaine

La contrainte de Von Mises est une quantité mécanique très utile en terme de dimensionnement pour l’ingénieur. La difficulté pour estimer l’erreur sur cette quantité réside dans le fait que son expression n’est pas linéaire par rapport au déplacement et donc inutilisable avec la technique présentée. La difficulté est levée grâce à une linéarisation de l’expression et sous certaines hypothèses [bib21].

Le tenseur des contraintes \(\sigma\) peut se décomposer en la somme d’un tenseur déviatorique \({\sigma}^{d}\) (dont la trace est nulle) et d’un tenseur sphérique \({\sigma}^{s}\) (dont les termes extra diagonaux sont nuls):

(1679)\[\sigma ={\sigma}^{d}+{\sigma}^{s}\]

où \({\sigma}^{d}\) et \({\sigma}^{s}\) sont définis par:

(1680)\[{\sigma}^{d}=\sigma -\frac{1}{3}\mathrm{tr}(\sigma )I {\sigma}^{s}=\frac{1}{3}\mathrm{tr}(\sigma )I\]

La contrainte de Von Mises, notée \({\sigma}_{\mathrm{vm}}\) , est définie par:

(1681)\[{\sigma}_{\mathrm{vm}}=\sqrt{\frac{3}{2}{\sigma}^{d}:{\sigma}^{d}}\]

En élasticité linéaire, la contrainte dépend linéairement du déplacement, ce qui nous permet d’écrire:

(1682)\[{\sigma}^{d}(u)={\sigma}^{d}({u}^{h})+{\sigma}^{d}(e)\]

En utilisant ceci dans la définition de la contrainte de Von Mises, on obtient:

(1683)\[{\sigma}_{\mathrm{vm}}(u)=\sqrt{\frac{3}{2}}(({\sigma}^{d}({u}^{h})+{\sigma}^{d}(e)):({\sigma}^{d}({u}^{h})+{\sigma}^{d}(e)))\]

L’expression précédente est réécrite afin de faire apparaître une fonction dont le développement limité est connu:

(1684)\[\begin{split}\begin{array}{ccc}{\sigma}_{\mathrm{vm}}(u)& =& \sqrt{\frac{3}{2}({\sigma}^{d}({u}^{h}):{\sigma}^{d}({u}^{h})+2{\sigma}^{d}({u}^{h}):{\sigma}^{d}(e)+{\sigma}^{d}(e):{\sigma}^{d}(e))}\\ & =& {\sigma}_{\mathrm{vm}}({u}^{h})\sqrt{1+\frac{3}{2}\frac{2{\sigma}^{d}({u}^{h}):{\sigma}^{d}(e)+{\sigma}^{d}(e):{\sigma}^{d}(e)}{{\sigma}_{\mathrm{vm}}^{2}({u}^{h})}}\end{array}\end{split}\]

Le développement limité à l’ordre 1 quand \(e\) tend vers zéro permet d’obtenir l’approximation suivante:

(1685)\[\begin{split}\begin{array}{ccc}{\sigma}_{\mathrm{vm}}(u)& \approx & {\sigma}_{\mathrm{vm}}({u}^{h})(1+\frac{3}{2}\frac{{\sigma}^{d}({u}^{h}):{\sigma}^{d}(e)}{{\sigma}_{\mathrm{vm}}^{2}({u}^{h})})\\ & \approx & {\sigma}_{\mathrm{vm}}({u}^{h})+\frac{3}{2}\frac{{\sigma}^{d}({u}^{h}):{\sigma}^{d}(e)}{{\sigma}_{\mathrm{vm}}({u}^{h})}\end{array}\end{split}\]

Finalement, une approximation de l’erreur sur la contrainte de Von Mises est déterminée:

(1686)\[{\sigma}_{\mathrm{vm}}(u)-{\sigma}_{\mathrm{vm}}({u}^{h})\approx \frac{3}{2}\frac{{\sigma}^{d}({u}^{h}):{\sigma}^{d}(e)}{{\sigma}_{\mathrm{vm}}({u}^{h})}\]

Ainsi, pour estimer l’erreur sur la moyenne de la contrainte de Von Mises dans un domaine:

(1687)\[{\varepsilon}^{Q}=\frac{1}{\mid \omega \mid }{\int}_{\omega}({\sigma}_{\mathrm{vm}}(u)-{\sigma}_{\mathrm{vm}}({u}^{h}))d\Omega\]

la fonctionnelle suivante peut-être utilisée:

(1688)\[{Q}^{\mathrm{vm}}(v)=\frac{3}{2\mid \omega \mid }{\int}_{\omega}\frac{{\sigma}^{d}({u}^{h}):{\sigma}^{d}(v)}{{\sigma}_{\mathrm{vm}}({u}^{h})}d\Omega\]

Comme pour la composante des contraintes sur un domaine, pour déterminer \({Q}^{\mathrm{vm}}(v)\) , il suffit de déterminer:

(1689)\[-a({u}^{Q},v)={\int}_{\Omega}\Sigma ({u}^{Q}):\sigma (v)d\Omega\]

Ce qui correspond au cas particulier où \({f}^{Q}=0\) et \({F}^{Q}=0\) et où \(\Sigma ({u}^{Q})\) est un opérateur symétrique constant par élément sur \(\omega\) et nul ailleurs. Le fait que dans l’équation apparaisse \({\sigma}^{d}(v)\) et non \(\sigma (v)\) nécessite de réécrire la fonctionnelle \({Q}^{\mathrm{vm}}(v)\) pour faire apparaître \(\sigma (v)\) :

(1690)\[{Q}^{\mathrm{vm}}(v)=\frac{3}{2\mid \omega \mid }{\int}_{\omega}\frac{{\sigma}^{d}({u}^{h}):(\sigma (v)-{\sigma}^{s}(v))}{{\sigma}_{\text{vm}}({u}^{h})}d\Omega\]

(1691)\[{Q}^{\mathrm{vm}}(v)=\frac{3}{2\mid \omega \mid }{\int}_{\omega}\frac{{\sigma}^{d}({u}^{h}):\sigma (v)-{\sigma}^{d}({u}^{h}):{\sigma}^{s}(v)}{{\sigma}_{\mathrm{vm}}({u}^{h})}d\Omega\]

Enfin, la double contraction entre un tenseur déviatorique et un tenseur sphérique étant nulle, la fonctionnelle s’écrit:

(1692)\[{Q}^{\mathrm{vm}}(v)=\frac{1}{\mid \omega \mid }{\int}_{\omega}(\frac{3}{2{\sigma}_{\mathrm{vm}}({u}^{h})}{\sigma}^{d}({u}^{h})):\sigma (v)d\Omega\]

avec \(\Sigma ({u}^{Q})=\frac{3}{2{\sigma}_{\mathrm{vm}}({u}^{Q})}{\sigma}^{d}({u}^{Q})\) et \({u}^{Q}={u}^{h}\) .

De manière analogue au cas de la moyenne des contraintes, le problème dual local s’écrit:

(1693)\[\begin{split}\lbrace \begin{array}{cccc}{({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl})}_{,j}& =& 0& \mathrm{dans}\Omega \\ {u}_{i}(x)& =& 0& \mathrm{sur}{\Gamma}_{U}\\ ({\sigma}_{ij}(\omega )-{a}_{ijkl}{\Sigma}_{kl}){n}_{j}(x)& =& 0& \mathrm{sur}{\Gamma}_{F}\end{array}\end{split}\]

Imposer \({Q}^{\mathrm{vm}}(v)\) comme chargement du problème dual revient à imposer un champ de déformation initiale.

Grandeurs liées à la mécanique de la rupture

Les facteurs d’intensité des contraintes

Une méthode de calcul de \({K}_{I}\) , \({K}_{\mathrm{II}}\) et \({K}_{\mathrm{III}}\) en 2D (plan et axisymétrique) et 3D est rappelée; elle est basée sur l’extrapolation des sauts de déplacements sur les lèvres de la fissure. Cette méthode n’est applicable qu’au cas des fissures planes, dans des matériaux homogènes et isotropes ; elle est utilisée dans la commande POST_K1_K2_K3.

Les facteurs d’intensité de contraintes sont identifiés à partir du saut des déplacements \(\left[U\right]\) par une méthode des moindres carrés sur un segment de longueur \({r}_{\max}\) . Ceci s’exprime par le problème de minimisation suivant:

(1694)\[K\text{minimise}F(k)=\frac{1}{2}{\int}_{0}^{{r}_{\max}}{(C\left[U(r)\right]-k\sqrt{r})}^{2}\mathrm{dr}\]

Pour résoudre ce problème, une valeur de \(k\) est cherchée telle que la dérivée de \(F\) s’annule:

(1695)\[F'(K)={\int}_{0}^{{r}_{\max}}\sqrt{r}(K\sqrt{r}-C\left[U(r)\right])\mathrm{dr}=0\]

La séparation de l’intégrale donne:

(1696)\[{\int}_{0}^{{r}_{\max}}rK\mathrm{dr}={\int}_{0}^{{r}_{\max}}\sqrt{r}C\left[U(r)\right]\mathrm{dr}\]

(1697)\[K{\left[\frac{{r}^{2}}{2}\right]}_{0}^{{r}_{\max}}={\int}_{0}^{{r}_{\max}}\sqrt{r}C\left[U(r)\right]\mathrm{dr}\]

La formule explicite pour calculer \(K\) se déduit directement:

(1698)\[K=\frac{\mathrm{2C}}{{r}_{\max}^{2}}{\int}_{0}^{{r}_{\max}}\sqrt{r}\left[U(r)\right]\mathrm{dr}\]

En réalisant l’intégration par une méthode des trapèzes, la relation précédente devient:

(1699)\[K=\frac{C}{{r}_{\max}^{2}}\sum_{i=1}^{N}({\left[U\right]}_{i-1}\sqrt{{r}_{i-1}}+{\left[U\right]}_{i}\sqrt{{r}_{i}})({r}_{i}-{r}_{i-1})\]

où \(N\) est le nombre de nœuds sur le segment \(\left[0,{r}_{\max}\right]\) . On note que \(K\) est une forme linéaire du champ de déplacement pour un \({r}_{\max}\) fixé, typiquement cela correspond à 4 ou 5 éléments.

Pour exprimer les facteurs d’intensité des contraintes en terme de chargement pour le problème dual, nous allons utiliser l’équation .

En faisant apparaître les composantes du déplacement sur la lèvre supérieure \({U}_{i}^{\sup}(r)\) et les composantes du déplacement sur la lèvre inférieure \({U}_{i}^{\inf}(r)\) , on obtient:

(1700)\[{K}_{i}=\frac{2C}{{r}_{\max}^{2}}{\int}_{\Gamma}({U}_{i}^{\sup}(r)-{U}_{i}^{\inf}(r))\sqrt{r}\mathrm{dr}\]

Le contour \(\Gamma\) est constitué des deux lèvres de la fissure; \(\Gamma\) est la lèvre supérieure et \({\Gamma}^{\inf}\) est la lèvre inférieure. Ceci permet de séparer l’intégrale précédente en deux parties:

(1701)\[{K}_{i}=\frac{\mathrm{2C}}{{r}_{\max}^{2}}\left[{\int}_{{\Gamma}^{\sup}}{U}_{i}(r)\sqrt{r}\mathrm{dr}-{\int}_{{\Gamma}^{\inf}}{U}_{i}(r)\sqrt{r}\text{dr}\right]\]

Afin de faire apparaître la fonction test \(v\) , les composantes des déplacements sont écrites comme le produit scalaire de deux vecteurs:

(1702)\[{K}_{i}=\frac{\mathrm{2C}}{{r}_{\max}^{2}}\left[{\int}_{{\Gamma}^{\sup}}v\cdot {g}_{i}(r)\sqrt{r}\mathrm{dr}-{\int}_{{\Gamma}^{\inf}}v\cdot {g}_{i}(r)\sqrt{r}\mathrm{dr}\right]\]

Avec \(i=I,\mathrm{II},\mathrm{III}\) ; \({g}_{I}\) le vecteur normal au plan de fissure, \({g}_{\mathrm{II}}\) le vecteur normal au fond de fissure et \({g}_{\mathrm{III}}\) le vecteur tangent au fond de fissure. Ceci signifie donc que l’expression de la quantité d’intérêt est valide en dimension 3 et permet donc l’estimation d’erreur des facteurs d’intensité des contraintes pour les problèmes en dimension 3.

Enfin en regroupant tous les termes, l’expression suivante est obtenue:

(1703)\[Q(v)={K}_{i}(v)={\int}_{{\Gamma}^{\sup}}{f}^{\sup}\cdot v\mathrm{dr}+{\int}_{{\Gamma}^{\inf}}{f}^{\inf}\cdot v\mathrm{dr}\]

avec:

(1704)\[{f}^{\sup}=\frac{2C}{{r}_{\max}^{2}}{g}_{i}(r)\sqrt{r}\]

Cette expression correspond à un chargement linéique en dimension2 (surfacique en dimension3) à imposer sur une partition des lèvres de la fissure partant de la pointe et de longueur \({r}_{\max}\) . Le choix du rayon \({r}_{\max}\) doit être fait de manière à ce qu’il soit dans la zone singulière et contiennent les éventuels éléments de Barsoum.

Mise en œuvre et utilisation dans Code_Aster

L’estimateur d’erreur en quantité d’intérêt est implémenté dans Code_Aster en mécanique pour tous les éléments, pour les modélisations C_PLAN, D_PLAN, AXIS et 3D.

Pour estimer l’erreur en quantité d’intérêt, nous avons vu qu’il fallait résoudre un problème adjoint dont le chargement en effort est la quantité d’intérêt. Le tableau suivant résume et regroupe les chargements à imposer dans Code_Aster , par le biais de la commande AFFE_CHAR_MECA, pour chacune des quantités d’intérêt. Le tableau est donné pour du 2D mais le 3D se déduit facilement. Pour les quantités composantes des contraintes, la valeur à indiquer est 1.0 pour les termes diagonaux et 0.5 pour les autres.

	\(Q(v)\)	Mot-clé de AFFE_CHAR_MECA	Composante
Moyenne d’une composante du déplacement sur un contour	\(\frac{1}{\mid \gamma \mid }{\int}_{\gamma}{v}_{x}\mathrm{d\Gamma }\)	FORCE_CONTOUR	FX=1.
	\(\frac{1}{\mid \gamma \mid }{\int}_{\gamma}{v}_{y}d\Gamma\)	FORCE_CONTOUR	FY=1.
Moyenne d’une composante du déplacement sur un sous‑domaine	\(\frac{1}{\mid \omega \mid }{\int}_{\omega}{v}_{x}d\Omega\)	FORCE_INTERNE	FX=1.
	\(\frac{1}{\mid \omega \mid }{\int}_{\omega}{v}_{y}d\Omega\)	FORCE_INTERNE	FY=1.
Moyenne du déplacement normal sur un contour	\(\frac{1}{\mid \gamma \mid }{\int}_{\gamma}{v}_{n}d\Gamma\)	PRES_REP	PRES=1.
Moyenne d’une composante des contraintes sur un sous‑domaine	\(\frac{1}{\mid \omega \mid }{\int}_{\omega}{\sigma}_{xx}d\Omega\)	PRE_EPSI	EPXX=1.
	\(\frac{1}{\mid \omega \mid }{\int}_{\omega}{\sigma}_{yy}d\Omega\)	PRE_EPSI	EPYY=1.
	\(\frac{1}{\mid \omega \mid }{\int}_{\omega}{\sigma}_{xy}d\Omega\)	PRE_EPSI	EPXY=0.5

Lorsque les quantités d’intérêt ne sont pas linéaires ou lorsque leur expression ne permet pas de les exprimer facilement en terme de chargement un traitement particulier est nécessaire. Pour les quantités non-linéaires, la théorie ne s’applique pas; c’est pourquoi il faut les linéariser. La contrainte de Von Mises en est l’exemple typique. La linéarisation nous amène à la formule suivante:

(1705)\[{Q}^{\mathrm{vm}}(v)=\frac{1}{\mid \omega \mid }{\int}_{\omega}(\frac{3}{2{\sigma}_{\mathrm{vm}}({u}^{h})}{\sigma}^{d}({u}^{h})):\sigma (v)d\Omega\]

Ainsi on approche la quantité d’intérêt Von Mises par une autre quantité linéaire. Le chargement à imposer dans Code_Aster est un champ de contrainte composé:

(1706)\[\frac{\mathrm{3K}}{2}\frac{{\sigma}^{d}({u}^{h})}{{\sigma}_{\mathrm{vm}}({u}^{h})}\]

Lorsque la définition de la quantité n’est pas adaptée à l’expression d’un chargement, il faut, là aussi, utiliser un traitement qui sera propre à chaque quantité d’intérêt. Pour nos applications, les facteurs d’intensité des contraintes (FIC) sont des grandeurs intéressantes. En utilisant une méthode basée sur l’extrapolation des sauts des déplacements sur les lèvres de la fissure, on peut exprimer les FIC de la manière suivante, avec \(C=1\) en contrainte plane et \(C=1-{\nu}^{2}\) en déformation plane:

(1707)\[K=\frac{\mathrm{2C}}{{r}_{\max}^{2}}{\int}_{0}^{{r}_{\max}}\sqrt{r}\left[U(r)\right]\mathrm{dr}\]

Ce qui permet d’obtenir l’expression suivante:

(1708)\[Q(v)={K}_{i}(v)={\int}_{{\Gamma}^{\sup}}{f}^{\sup}\cdot v\mathrm{dr}+{\int}_{{\Gamma}^{\inf}}{f}^{\inf}\cdot v\mathrm{dr}\]

avec:

(1709)\[{f}^{\sup}=\frac{\mathrm{2C}}{{r}_{\max}^{2}}g(r)\sqrt{r}\]

avec \(g(r)\) vecteur qui oriente la fissure.

Dans Code_Aster , ce vecteur est connu sous le nom de VECT_K1. On voit ainsi qu’il faut imposer un effort sur chacune des lèvres de la fissure (:math:`` et \({f}^{\inf}\) ). Le tableau suivant reprend ces deux résultats.

	\(Q(v)\)	Mot-clé de AFFE_CHAR_MECA_F	Composante
Moyenne de la contrainte de Von Mises	\(\frac{1}{\mid \omega \mid }{\int}_{\omega}{\sigma}_{\mathrm{vm}}d\Omega\)	PRE_SIGM	Champ composé
Facteur d’intensité de contrainte	\({K}_{I}\) \({K}_{\mathrm{II}}\) \({K}_{\mathrm{III}}\)	FORCE_CONTOUR	F+ sur lèvre sup. F- sur lèvre inf.

Le calcul de l’erreur en quantité d’intérêt passe par la résolution d’un calcul élastique (correspondant au problème dual) dont le chargement est d’une part, le même que celui du problème principal (ou primal) sur les conditions de Dirichlet (donc en déplacement imposé), et d’autre part, celui représentant la quantité d’intérêt choisie . Il faut ensuite fournir ce résultat dans CALC_ERREUR via le mot-clé RESU_DUAL.

Enfin, les différentes options possibles dans CALC_ERREUR sont:

• “QIZ1_ELEM” (respectivement “QIZ2_ELEM”): estimateur d’erreur en quantité d’intérêt basé sur la méthode de Zhu-Zienkiewicz. Il faut calculer préalablement l’option “ERZ1_ELEM” ou “ERZ2_ELEM”.

• “QIRE_ELEM” (par élément) ou “QIRE_ELNO” (par élément aux nœuds): estimateur d’erreur en quantité d’intérêt basé sur les résidus en mécanique. Il faut calculer préalablement l’option “ERME_ELEM”.

Des exemples d’application sur des structures industrielles pourront être trouvé dans la thèse de J. DELMAS [bib26].

Liste de tests mettant en oeuvre les estimateurs d’erreur en quantité d’intérêt:

Test de l’estimateur d’erreur en quantite d’interet	sslv113c
test de l’estimateur d’erreur en quantite d’interet	sslv113d
Indicateur d’erreur mailles volumiques et fonctions	sslv153a
Détection des singularités dans une plaque fissurée (TRIA3)	ssnp130a
Estimateurs d’erreur en quantités d’intérêt pour le calcul des fic	zzzz257a
Estimateurs d’erreur en quantités d’intérêt pour le calcul des fic	zzzz257b
Estimateurs d’erreur en quantités d’intérêt pour le calcul des fic	zzzz257c
Estimateurs d’erreur en quantités d’intérêt pour le calcul des fic	zzzz257d

Bibliographie

[bib1]

16. 7. Ciarlet: The finite element methods for elliptic problems. North-Holland, 1978.

[bib2]

7. Strang & G. J. Fix: An analysis of the finite element method. Prentice-Hall, 1973.

[bib3]

1. Babuška and W. Rheinboldt: « A posteriori error estimates for the finite element method », International Journal For Numerical Methods In Engineering, (12), 1978, pp.1597-1615.

[bib4]

15. 3. Zienkiewicz & J. Z. Zhu: A simple error estimator and adaptive procedure for practical engineerng analysis. International Journal for Numerical Methods in Engineering, 24, pp.337-357, 1987.

[bib5]

16. Ladevèze: « Comparaison de modèles de milieux continus », Université Pierre et Marie Curie - Paris VI, 1975.

[bib6]

18. Verfürth: A review of a posteriori error estimation techniques for elasticity problems. Computer Methods in Applied Mechanics and Engineering, 176, pp.419-440, 1999.

[bib7] M. Ainsworth & J. T. Oden: A posteriori error estimation in finite element analysis. Computer Methods in Applied Mechanics and Engineering, 142, pp.1-88, 1997.

[bib8] P. Ladevèze & J. P. Pelle: Mastering Calculations in Linear and Nonlinear Mechanics, Springer NY, 2004.

[bib9]

5. Ramm, E. Rank, R. Rannacher, S. K., E. Stein, G. Wendland, W. Wittum, P. Wriggers & W. Wunderlich: Error-controlled adaptive finite elements in solid mechanics. Stein (eds.), Wiley, 2003.

[bib10]

20. Grätsch & K. Bathe: A posteriori error estimation techniques in practical finite element analysis. Computers and Structures, 83, pp.235-265, 2005.

[bib11]

20. Strouboulis & K. A. Haque: Recent experiences with error estimation and adaptivity, Part II: Error estimation for h-adaptive approximations on grids of triangles and quadrilaterals. Computer Methods in Applied Mechanics and Engineering, 100, pp.359-430, 1992.

[bib12]

1. Babuška, T. Strouboulis, A. Mathur & C. Upadhyay: Pollution-error in the h-version of the finite-element method and the local quality of a-posteriori error estimators. Finite Elements in Analysis and Design, 17, pp.273-321, 1994.

[bib13]

14. Wiberg & F. Abdulwahab: Patch recovery based on superconvergent derivatives and equilibrium. International Journal for Numerical Methods in Engineering, 36, pp.2703-2724, 1993.

[bib14]

1. Babuška & R. Rodriguez: Problem of the selection of an a posteriori error indicator based on smoothening techniques. International Journal for Numerical Methods in Engineering, 36, pp.539-567, 1993.

[bib15]

1. Babuška, T. Strouboulis, C. S. Upadhyay, S. K. Gangaraj & K. Copps: Validation of a posterior error estimators by numerical approach. International Journal for Numerical Methods In Engineering, 37, pp.1073-1123, 1994.

[bib16]

1. Babuška, T. Strouboulis, A. Mathur & C. Upadhyay: Pollution-error in the h-version of the finite-element method and the local quality of a-posteriori error estimators. Finite Elements in Analysis and Design, 17, pp.273-321, 1994.

[bib17]

1. Babuška, T. Strouboulis, C. S. Upadhyay & S. K. Gangaraj: A posteriori estimation and adaptive control of the pollution error in the h-version of the finite element method. International Journal for Numerical Methods in Engineering, 38, pp.4207-4235, 1995.

[bib18]

1. Babuška, T. Strouboulis, S. Gangaraj & C. Upadhyay: Pollution error in the h-version of the finite element method and the local quality of the recovered derivatives. Computer Methods in Applied Mechanics and Engineering, 140, pp.1-37, 1997.

[bib19] (1,2,3)

18. Becker & R. Rannacher: An optimal control approach to shape a posteriori error estimation in finite element methods. Acta Numerica, 10, pp.1-120, 2001.

[bib20] (1,2)

10. 20. Oden & S. Prudhomme: Estimation of modeling error in computational mechanics. Journal of Computational Physics, 182, pp.496-515, 2002.

[bib21] (1,2,3)

10. 20. Oden & S. Prudhomme: Goal-oriented error estimation and adaptivity for the finite element method. Computers and Mathematics with Applications, 41, pp.735-56, 2001.

[bib22]

13. Ainsworth & J. T. Oden: A posteriori error estimation in finite element analysis. Wiley-Interscience, 2000.

[bib23]

1. Sauer-Budge, J. Bonet, A. Huerta & J. Peraire: Computing bounds for linear functionals of exact weak solutions to poisson’s equation. SIAM Journal on Numerical Analysis, 42, pp.1610-1630, 2004.

[bib24]

19. Prudhomme & J. T. Oden: On goal-oriented error estimation for elliptics problems: application to the control of pointwise errors. Computer Methods in Applied Mechanics and Engineering, 176, pp.313-331, 1999.

[bib25] (1,2)

14. Parés, P. Díez & A. Huerta: Subdomain-based flux-free a posteriori error estimators. Computer Methods in Applied Mechanics and Engineering, 195, pp.297-323, 2006.

[bib26]

J. Delmas: Stratégies de contrôle d’erreur en calcul de structures industrielles. Mise en oeuvre d’estimation d’erreur en quantité d’intérêt et d’adaptation de maillage, thèse de doctorat de l’université de Picardie Jules Verne, 2008. Disponible au téléchargement à cette adresse: http://tel.archives-ouvertes.fr/tel-00311947

Description des versions du document

Indice doc	Version Aster	Auteur(s) ou contributeur(s), organisme	Description des modifications
A	8.3	J.DELMAS EDF/R&D/AMA	Texte initial