This document is under active development and discussion! If you find errors or omissions in this document, please don’t hesitate to submit an issue or open a pull request with a fix. We also encourage you to ask questions and discuss any aspects of the project on the Feel++ Gitter forum. New contributors are always welcome!

# Functional analysis

## 1. Outils d’analyse fonctionnelle

$\newcommand{\Feel}{\textsc{Feel++}} \newcommand{\eg}{\emph{e.g.}} \newcommand{\ie}{\emph{i.e.}} \newcommand{\ds}{\displaystyle} \newcommand{\cf}{{\it cf \hspace*{1 mm}}} \newcommand{\bdm}{\begin{displaymath}} \newcommand{\edm}{\end{displaymath}} \newcommand{\bi}{\begin{itemize}} \newcommand{\ei}{\end{itemize}} \newcommand{\be}{$$} \newcommand{\ee}{$$} \newcommand{\diver}{\hbox{div}} \newcommand{\hun}{\hbox{H}^1(\Omega)} \newcommand{\czero}{\hbox{C}^0(\bar{\Omega})} \newcommand{\dpa}[2]{\frac{\partial #1}{\partial #2}} \newcommand{\NN}{\hbox{I\hspace*{-2 mm} N}} \newcommand{\CC}{\mathbf{C}} \newcommand{\RR}{{\mathbb{R}}} \newcommand{\Rn}{\RR^n} \newcommand{\ZZ}{{\bf Z}} \newcommand{\dd}[2]{\frac{\partial #1}{\partial #2}} \newcommand{\saut}{\vspace*{5 mm}\\} \newcommand{\Th}{\mathcal{T}_h} \newcommand{\Thi}{\mathcal{T}_h^i} \newcommand{\Thb}{\mathcal{T}_h^b} \newcommand{\Rh}{\mathcal{R}_h} % simplicial conformal refinement \newcommand{\Sh}{\mathcal{S}_h} % stencil \newcommand{\Nh}{\mathcal{N}_h} % neighbours of neighbours \newcommand{\Fh}{\mathcal{F}_h} \newcommand{\Fhi}{\mathcal{F}_h^i} \newcommand{\Fhb}{\mathcal{F}_h^b} \newcommand{\jump}[1]{[\![ #1 ]\!]} \newcommand{\Nma}{{{N_{\mathrm{ma}}}}} \newcommand{\Nso}{{{N_{\mathrm{so}}}}} \newcommand{\Nno}{{{N_{\mathrm{no}}}}} \newcommand{\Ne}{{{N_{\mathrm{e}}}}} \newcommand{\Ilag}[1]{{{\mathcal{I}^{\mathrm{lag}}_{#1}}}} \newcommand{\diam}{{{\operatorname{diam}}}} \renewcommand{\max}{{{\operatorname{max}}}} \newcommand{\essinf}{{{\operatorname{ess}\, \operatorname{inf}}}} \newcommand{\calTh}{{{\mathcal{T}_h}}} \newcommand{\Ck}[1]{{{\mathcal{C}^#1}}} \newcommand{\Ck}[1]{{{\mathcal{C}^#1}}} \newcommand{\Pk}[1]{{{{\mathbb{P}_#1}}}} \newcommand{\Qk}[1]{{{{\mathbb{Q}_#1}}}} \newcommand{\Pch}[1]{{{{P^#1_{c,h}}}}} \newcommand{\Pcho}[1]{{{{P^#1_{c,h,0}}}}} \newcommand{\Qch}[1]{{{{Q^#1_{c,h}}}}} \newcommand{\Ich}[2]{{{{\mathcal{I}^#1_{c,h}}#2}}} \newcommand{\set}[1]{\left\{#1\right\}} \newcommand{\Next}{{{\mathrm{n}}}} \newcommand{\nf}{{{n_f}}} \newcommand{\ngeo}{{{n_{\mathrm{geo}}}}} \newcommand{\geo}{{{\mathrm{geo}}}} \renewcommand{\dim}[1]{{{\operatorname{dim}(#1)}}} \newcommand{\card}[1]{{{\operatorname{card}(#1)}}} \newcommand{\poly}[1]{{{\mathbb{#1}}}} \newcommand{\R}[1][d]{{\mathbb{R}^{#1}}} \newcommand{\tr}{{{\operatorname{tr}}}} \newcommand{\Id}{\mathcal{I}} \newcommand{\disp}[1][u]{{\mathbf{#1}}} \newcommand{\stresst}[1][(\mathbf{u})]{{\mathbf{\sigma#1}}} \newcommand{\deformt}[1][(\mathbf{u})]{{\mathbf{\varepsilon#1}}} \newcommand{\domain}{\Omega} \newcommand{\prect}[1]{\left$$#1\right$} %%% Local Variables: %%% mode: latex %%% TeX-master: "mef-intro" %%% End:$$

### 1.1. Quelques rappels

#### 1.1.1. Normes et produits scalaires

Soit $E$ un espace vectoriel.

Definition 1.
$\|.\|$ : $E \rightarrow \RR$ est une norme sur $E$ ssi elle vérifie : (N1):: $\left( \| x \| = 0 \right) \Longrightarrow (x=0)$ (N2):: $\forall\, \lambda\in\RR,\; \forall x\in E, \quad \| \lambda x \| = |\lambda| \; \| x \|$ (N3):: $\forall\, x,y \in E, \quad \| x+ y \| \le \|x \| + \|y\|$ (inégalité triangulaire)

Pour $E=\RR^n$ et $x=(x_1,\ldots,x_n) \in\RR^n$, on définit les normes

$\| x \|_1 = \sum_{i=1}^n |x_i| \qquad \| x \|_2 = \left( \sum_{i=1}^n x_i^2 \right)^{1/2} \qquad \| x \|_\infty = \sup_{i} |x_i|$

Definition 2: Produit scalaire
On appelle produit scalaire sur $E$ toute forme bilinéaire symétrique définie positive. $\quad<.,.>$ : $E\times E \rightarrow \RR$ est donc un produit scalaire sur $E$ ssi il vérifie :
S1
$\forall\; x,y \in E, \quad <x,y> = <y,x>$
S2
$\forall\; x_1,x_2,y \in E, \quad <x_1+x_2,y> = <x_1,y>+ <x_2,y>$
S3
$\forall\; x,y \in E, \, \forall\, \lambda\in\RR,\quad <\lambda x,y> = \lambda <x,y>$
S4

$\forall\; x \in E, x\ne 0, \quad <x,x>\; > 0$

A partir d’un produit scalaire, on peut définir une norme induite : $\| x \| = \sqrt{<x,x>}$
On a alors, d’après (N3), l’inégalité de Cauchy-Schwarz : ${ | <x,y> | \le \| x \| \; \| y \| }$

Example 1.
Pour $E=\RR^n$, on définit le produit scalaire $<x,y> = \sum_{i=1}^n x_i \, y_i.$ Sa norme induite est $\| . \|_2$ définie précédemment.
Definition 3: Espace normé
Un espace vectoriel muni d’une norme est appelé espace normé.
Definition 4: Espace préhilbertien
Un espace vectoriel muni d’un produit scalaire est appelé espace préhilbertien. En particulier, c’est donc un espace normé pour la norme induite.

#### 1.1.2. Suites de Cauchy - espaces complets

Definition 5.
Soit $E$ un espace vectoriel et $(x_n)_n$ une suite de $E$. $(x_n)_n$ est une suite de Cauchy si et seulement si
$\forall \varepsilon > 0,\;\; \exists N / \forall p>N, \forall q>N, \quad \|x_p - x_q \| < \varepsilon$
Theorem 1.
Toute suite convergente est de Cauchy.
 La réciproque de Theorem 1 est fausse.
Definition 6: Espace complet
Un espace vectoriel est complet ssi toute suite de Cauchy y est convergente.
Definition 7: Espace de Banach
Un espace normé complet est un espace de Banach.
Definition 8: Espace de Hilbert
Un espace préhilbertien complet est un espace de Hilbert.
Definition 9: Espace euclidien
Un espace de Hilbert de dimension finie est appelé espace euclidien.

#### 1.1.3. Espaces fonctionnels

Definition 10.
Un espace fonctionnel est un espace vectoriel dont les éléments sont des fonctions.
Example 2: ${\cal C}^p([a;b$)]
${\cal C}^p([a;b$)] désigne l’espace des fonctions définies sur l’intervalle $[a,b$], dont toutes les dérivées jusqu’à l’ordre $p$ existent et sont continues sur $[a,b$].

Dans la suite, les fonctions seront définies sur un sous-ensemble de $\RR^n$ (le plus souvent un ouvert noté $\Omega$), à valeurs dans $\RR$ ou $\RR^p$.

Example 3.
La température $T(x,y,z,t)$ en tout point d’un objet $\bar{\Omega}\subset \RR^3$ est une fonction de $\bar{\Omega} \times \RR \longrightarrow \RR$.

Les normes usuelles les plus simples sur les espaces fonctionnels sont les normes $\bf L^p$ définies par :

$\| u \|_{L^p} = \left ( \int_{\Omega } |u|^p \right) ^{1/p} \quad ,\; p\in [1,+ \infty[ , \qquad \hbox{et}\qquad \| u \|_{L^\infty} = {\hbox{Sup}}_{\Omega } |u|$

Comme on va le voir, ces formes $L^p$ ne sont pas nécessairement des normes. Et lorsqu’elles le sont, les espaces fonctionnels munis de ces normes ne sont pas nécessairement des espaces de Banach. Par exemple, les formes $L^\infty$ et $L^1$ sont bien des normes sur l’espace ${\cal C}^0([a;b$)], et cet espace est complet si on le munit de la norme $L^\infty$, mais ne l’est pas si on le munit de la norme $L^1$.

Pour cette raison, on va définir les espaces ${\cal L}^p(\Omega)(p\in [1,+ \infty[$) par

$[{\cal L}^p(\Omega) = \left\{ u : \Omega \rightarrow \RR, \hbox{ mesurable, et telle que }\int_\Omega |u|^p<\infty \right\}$
 on rappelle qu’une fonction $u$ est mesurable ssi $\{ x / |u(x)|0$.

Sur ces espaces ${\cal L}^p(\Omega)$, les formes $L^p$ ne sont pas des normes. En effet, $\| u \|_{L^p} = 0$ implique que $u$ est nulle presque partout dans ${\cal L}^p(\Omega)$, et non pas $u=0$. C’est pourquoi on va définir les espaces $\bf L^p(\Omega)$ :

Definition 11: Égalité presque partout
$L^p(\Omega)$ est la classe d’équivalence des fonctions de ${\cal L}^p(\Omega)$ pour la relation d’équivalence égalité presque partout. Autrement dit, on confondra deux fonctions dès lors qu’elles sont égales presque partout, c’est à dire qu’elles ne diffèrent que sur un ensemble de mesure nulle.
Theorem 2: Les espaces $L^p(\Omega)$ sont complets
La forme $L^p$ est une norme sur $L^p(\Omega)$, et $L^p(\Omega)$ muni de la norme $L^p$ est un espace de Banach (c.a.d. est complet).

Un cas particulier très important est $p=2$. On obtient alors l’espace fonctionnel $L^2(\Omega)$, c’est à dire l’espace des fonctions de carré sommable sur $\Omega$ (à la relation d’équivalence égalité presque partout près). A la norme $L^2$ : $\| u \|_{L^2} = \left( \int_\Omega u^2 \right)^{1/2}$, on peut associer la forme bilinéaire $(u,v)_{L^2} = \int_\Omega u\, v$. Il s’agit d’un produit scalaire, dont dérive la norme $L^2$.

D’où le théorème suivant

Theorem 3.
$L^2(\Omega)$ est un espace de Hilbert.

### 1.2. Notion de dérivée généralisée

Nous venons de définir des espaces fonctionnels complets, ce qui sera un bon cadre pour démontrer l’existence et l’unicité de solutions d’équations aux dérivées partielles, comme on le verra plus loin notamment avec le théorème de Lax-Milgram.

Toutefois, on a vu que les éléments de ces espaces $L^p$ ne sont pas nécessairement des fonctions très régulières.

Dès lors, les dérivées partielles de telles fonctions ne sont pas forcément définies partout.

Pour s’affranchir de ce problème, on va étendre la notion de dérivation.

Le véritable outil à introduire pour cela est la notion de distribution, due à L. Schwartz (1950).

Par manque de temps dans ce cours, on se contentera ici d’en donner une idée très simplifiée, avec la notion de dérivée généralisée.

Cette dernière a des propriétés beaucoup plus limitées que les distributions, mais permet de “sentir" les aspects nécessaires pour mener à la formulation variationnelle.

Dans la suite, $\Omega$ sera un ouvert (pas nécessairement borné) de $\RR^n$.

#### 1.2.1. Fonctions tests

Soit $\varphi : \Omega \rightarrow \RR$.

Definition 12.
On appelle support de $\bf \varphi$ l’adhérence de $\{ x \in \Omega / \varphi(x) \ne 0 \}$.
Example 4.
Pour $\Omega = ]-1,1$$, et $\varphi$ la fonction constante égale à 1, $\hbox{Supp}\, \varphi = [-1,1$]. Definition 13: Espace des fonctions tests On note ${\cal D}(\Omega)$ l’espace des fonctions de $\Omega$ vers $\RR$, de classe ${\cal C}^\infty$, et à support compact inclus dans $\Omega$. ${\cal D}(\Omega)$ est parfois appelé espace des fonctions-tests. Example 5. L’exemple le plus classique dans le cas 1-D est la fonction \varphi(x) = \left\{ \begin{array}{ll} { e^{- \frac{1}{1-x^2}} } & \hbox{si } |x|<1\\ 0 & \hbox{si } |x|\ge 1\\ \end{array} \right.$ $\varphi$ est une fonction de ${\cal D}(]a,b$)$ pour tous $a < -1 < 1 < b$. Cet exemple s’étend aisément au cas multi-dimensionnel ($n>1$). Soit $a\in\Omega$ et $r>0$ tel que la boule fermée de centre $a$ et de rayon $r$ soit incluse dans $\Omega$. On pose alors : \[ \varphi(x) = \left\{ \begin{array}{ll} { e^{- \frac{1}{r^2-|x-a|^2}} } & \hbox{si } |x-a|<r\\ 0 & \hbox{sinon }\\ \end{array} \right.$ $\varphi$ ainsi définie est un élément de ${\cal D}(\Omega)$. Theorem 4: Adhérence de $\overline{{\cal D}(\Omega)$ $\overline{{\cal D}(\Omega) } = L^2(\Omega)$ #### 1.2.2. Dérivée généralisée Soit $u\in {\cal C}^1(\Omega)$ et $\varphi \in {\cal D}(\Omega)$. Par intégration par parties (annexe [sec:green]), on a : $\int_\Omega \partial_i u\; \varphi = - \int_\Omega u \; \partial_i\varphi + \int_{\partial \Omega} u \; \varphi \; {\bf e}_i.{\bf n}$ Ce dernier terme (intégrale sur le bord de $\Omega$) est nul car $\varphi$ est à support compact (donc nul sur $\partial \Omega$). Or $\int_\Omega u \; \partial_i\varphi$ a un sens par exemple dès que $u\in L^2(\Omega)$. Donc le terme $\int_\Omega \partial_i u\; \varphi$ a aussi du sens, sans que $u$ ne soit nécessairement de classe ${\cal C}^1$. Ceci permet de définir $\partial_i u$ même dans ce cas. Definition 14. cas 1-D $\quad$ Soit $I$ un intervalle de $\RR$, pas forcément borné. On dit que $u\in L^2(I)$ admet une dérivée généralisée dans $L^2(I)$ ssi $\exists u_1\in L^2(I)$ telle que $\forall \varphi\in {\cal D}(I), \quad \int_I u_1\;\varphi = - \int_I u \varphi'$ Example 6. Soit $I=]a,b[$ un intervalle borné, et $c$ un point de $I$. On considère une fonction $u$ formée de deux branches de classe ${\cal C}^1$, l’une sur $]a,c[$, l’autre sur $]c,b[$, et se raccordant de façon continue mais non dérivable en $c$. Alors $u$ admet une dérivée généralisée définie par $u_1(x)=u'(x)\quad \forall x\ne c$. En effet : $\forall \varphi\in {\cal D}(]a,b[)\qquad \int_a^b u \varphi' = \int_a^c + \int_c^b = - \int_a^c u' \varphi - \int_c^b u'\varphi + \underbrace{(u(c^-)-u(c^+))}_{=0} \, \varphi(c)$ par intégration par parties. La valeur $u_1(c)$ n’a pas d’importance: on a de toute façon au final la même fonction de $L^2(I)$, puisqu’elle est définie comme classe d’équivalence de la relation d’équivalence égalité presque partout. Definition 15. En itérant, on dit que $u$ admet une dérivée généralisée d’ordre $\bf k$ dans $L^2(I)$, notée $u_k$, ssi ${\forall \varphi\in {\cal D}(I), \quad \int_I u_k\;\varphi = (- 1)^k \; \int_I u \varphi^{(k)} }$ Ces définitions s’étendent naturellement pour la définition de dérivées partielles généralisées, dans le cas $n>1$. Theorem 5: Unicité de la dérivée généralisée Quand elle existe, la dérivée généralisée est unique. Theorem 6. Quand $u$ est de classe ${\cal C}^1(\bar{\Omega})$, la dérivée généralisée est égale à la dérivée classique. ### 1.3. Espaces de Sobolev #### 1.3.1. Les espaces $H^m$ Definition 16. ${ H^1(\Omega) = \left\{ u \in L^2(\Omega)\; / \; \partial_i u \; \in L^2(\Omega), \quad 1 \le i \le n \right\} }$ où $\partial_i u$ est définie au sens de la dérivée généralisée. $H^1(\Omega)$ est appelé espace de Sobolev d’ordre 1. Definition 17. Pour tout entier $m\ge 1$, $H^m(\Omega) = \left\{ u \in L^2(\Omega) \; / \; \partial^\alpha u \; \in L^2(\Omega) \quad \forall \alpha =(\alpha_1,\ldots,\alpha_n) \in \NN^n\hbox{ tel que}\; |\alpha|= \alpha_1+\cdots+\alpha_n \le m \right\}$ </div>+ $H^m(\Omega)$ est appelé *espace de Sobolev d’ordre $\bf m$. Par extension, on voit aussi que $H^0(\Omega)=L^2(\Omega)$. Dans le cas de la dimension 1, on écrit plus simplement pour $I$ ouvert de $\RR$ : + $H^m(I) = \left\{ u \in L^2(I) \; / \; u', \ldots, u^{(m)} \in L^2(I) \right\}$ Theorem 7. $H^1(\Omega)$ est un espace de Hilbert pour le produit scalaire $(u,v)_1 = \int_\Omega u \, v\, + \sum_{i=1}^n \; \int_\Omega \partial_i u \; \partial_i v = (u,v)_0 + \sum_{i=1}^n (\partial_i u, \partial_i v )_0$ en notant stem:[(.,.)_0] le produit scalaire stem:[L^2]. On notera stem:[\|.\|_1] la norme associée à stem:[(.,.)_1]. </div>+ On définit de même un produit scalaire et une norme sur $H^m(\Omega)$ par $$(u,v)_m = \sum_{|\alpha| \le m} ( \partial^\alpha u , \partial^\alpha v )_0 \qquad \hbox{ et }\qquad \| u \|_m = (u,u)_m^{1/2}]$ Theorem 8. $H^m(\Omega)$ muni du produit scalaire $(.,.)_m$ est un espace de Hilbert.[thr:8] Theorem 9. Si $\Omega$ est un ouvert de $\RR^n$ de frontière $\partial\Omega$ “suffisamment régulière" (par exemple ${\cal C}^1$), on a l’inclusion : $H^m(\Omega) \subset {\cal C}^k(\bar{\Omega})$ pour ${ k < m-\frac{n}{2} }$ Example 7. En particulier, on voit que pour un intervalle $I$ de $\RR$, on a $H^1(I) \subset {\cal C}^0(\bar{I})$, c’est à dire que, en 1-D, toute fonction $H^1$ est continue. L’exemple de $u(x) = x\, \sin\frac{1}{x}$ pour $x\in$0,1]] et $u(0)=0$ montre que la réciproque est fausse. L’exemple de $u(x,y) = | \ln (x^2+y^2) |^k$ pour $0<k<1/2$ montre qu’en dimension supérieure à 1 il existe des fonctions $H^1$ discontinues. #### 1.3.2. Trace d’une fonction Pour pouvoir faire les intégrations par parties qui seront utiles par exemple pour la formulation variationnelle, il faut pouvoir définir le prolongement (la trace) d’une fonction sur le bord de l’ouvert $\Omega$. $n=1$ (cas 1-D) on considère un intervalle ouvert $I=]a,b[$ borné. On a vu que $H^1(I) \subset {\cal C}^0(\bar{I})$. Donc, pour $u\in H^1(I)$, $u$ est continue sur $[a,b$], et $u(a)$ et $u(b)$ sont bien définies. $n>1$ on n’a plus $H^1(\Omega) \subset {\cal C}^0(\bar{\Omega})$. Comment alors définir la trace ? La démarche est la suivante : • On définit l’espace ${\cal C}^1(\bar{\Omega}) = \left\{ \varphi : \Omega \rightarrow \RR \;/\; \exists O \hbox{ ouvert contenant } \bar{\Omega},\; \exists \psi \in {\cal C}^1(O),\; \psi_{|\Omega} = \varphi \right\}$ Autrement dit, ${\cal C}^1(\bar{\Omega})$ est l’espace des fonctions ${\cal C}^1$ sur $\Omega$, prolongeables par continuité sur $\partial\Omega$ et dont le gradient est lui-aussi prolongeable par continuité. Il n’y a donc pas de problème pour définir la trace de telles fonctions. • On montre que, si $\Omega$ est un ouvert borné de frontière $\partial\Omega$ “assez régulière", alors ${\cal C}^1(\bar{\Omega})$ est dense dans $H^1(\Omega)$. • L’application linéaire continue, qui à toute fonction $u$ de ${\cal C}^1(\bar{\Omega})$ associe sa trace sur $\partial\Omega$, se prolonge alors en une application linéaire continue de $H^1(\Omega)$ dans $L^2(\partial\Omega)$, notée $\gamma_0$, qu’on appelle application trace. On dit que $\gamma_0(u)$ est la trace de $u$ sur $\partial\Omega$.  Pour une fonction $u$ de $H^1(\Omega)$ qui soit en même temps continue sur $\bar{\Omega}$, on a évidemment $\gamma_0(u) = u_{|\partial\Omega}$. C’est pourquoi on note souvent par abus simplement $u_{|\partial\Omega}$ plutôt que $\gamma_0(u)$. On peut de façon analogue définir $\gamma_1$, application trace qui permet de prolonger la définition usuelle de la dérivée normale sur $\partial\Omega$. Pour $u\in H^2(\Omega)$, on a $\partial_i u \in H^1(\Omega)$, $\forall i=1,\ldots,n$, et on peut donc définir $\gamma_0(\partial_i u)$. La frontière $\partial\Omega$ étant “assez régulière" (par exemple, idéalement, de classe ${\cal C}^1$), on peut définir la normale $n=\left( \begin{array}{l} n_1 \\ \vdots \\ n_n \end{array} \right)$ en tout point de $\partial\Omega$. On pose alors ${\gamma_1(u) = \sum_{i=1}^n \gamma_0(\partial_i u) n_i}$. Cette application continue $\gamma_1$ de $H^2(\Omega)$ dans $L^2(\partial\Omega)$ permet donc bien de prolonger la définition usuelle de la dérivée normale. Dans le cas où $u$ est une fonction de $H^2(\Omega)$ qui soit en même temps dans ${\cal C}^1(\bar{\Omega})$, la dérivée normale au sens usuel de $u$ existe, et $\gamma_1(u)$ lui est évidemment égal. C’est pourquoi on note souvent, par abus, $\partial_n u$ plutôt que $\gamma_1(u)$. #### 1.3.3. Espace $H^1_0(\Omega)$ Definition 18. Soit $\Omega$ ouvert de $\RR^n$. L’espace $H^1_0(\Omega)$ est défini comme l’adhérence de ${\cal D}(\Omega)$ pour la norme $\|.\|_1$ de $H^1(\Omega)$. (on rappelle que ${\cal D}(\Omega)$ est l’espace des fonctions ${\cal C}^\infty$ sur $\Omega$ à support compact, encore appelé espace des fonctions tests) Theorem 10. Par construction $H^1_0(\Omega)$ est un espace complet. C’est un espace de Hilbert pour la norme $\|.\|_1$ Si $n=1$ (cas 1-D)} on considère un intervalle ouvert $I=$a,b[] borné. Alors \[ H^1_0(]a,b[) = \left\{ u \in H^1(]a,b[),\; u(a)=u(b)=0 \right\}$ Si $n>1$ Si $\Omega$ est un ouvert borné de frontière“assez régulière" (par exemple ${\cal C}^1$ par morceaux), alors $H^1_0(\Omega) = \ker \gamma_0$. $H^1_0(\Omega)$ est donc le sous-espace des fonctions de $H^1(\Omega)$ de trace nulle sur la frontière $\partial\Omega$. Definition 19. Pour toute fonction $u$ de $H^1(\Omega)$, on peut définir : ${ |u|_1 = \left( \sum_{i=1}^n \| \partial_i u \|_0^2 \right)^{1/2} = \left( \int_\Omega \sum_{i=1}^n \left( \partial_i u \right)^2 dx \right)^{1/2} }$ Theorem 11: Inégalité de Poincaré Si $\Omega$ est borné dans au moins une direction, alors il existe une constante $C(\Omega)$ telle que $\forall u \in H^1_0(\Omega), \; \|u\|_0 \le C(\Omega)\; |u|_1.$ On en déduit que $|.|_1$ est une norme sur $H^1_0(\Omega)$, équivalente à la norme $\|.\|_1$. Le résultat précédent s’étend au cas où l’on a une condition de Dirichlet nulle seulement sur une partie de $\partial\Omega$, si $\Omega$ est connexe. On suppose que $\Omega$ est un ouvert borné connexe, de frontière ${\cal C}^1$ par morceaux. Soit $V=\left\{ v\in H^1(\Omega),\, v=0 \hbox{ sur }\Gamma_0 \right\}$ où $\Gamma_0$ est une partie de $\partial\Omega$ de mesure non-nulle. Alors il existe une constante $C(\Omega)$ telle que $\forall u \in V, \; \|u\|_{0,V} \le C(\Omega)\; |u|_{1,V}$, où $\|.\|_{0,V}$ et $|.|_{1,V}$ désignent les norme et semi-norme induites sur $V$. On en déduit que $|.|_{1,V}$ est une norme sur $V$, équivalente à la norme $\|.\|_{1,V}$. ### 1.4. Exercices 1. Montrer que les fonctions définies par ([eq:fonction-test1]) et ([eq:fonction-test2]) sont bien ${\cal C}^\infty$ à support compact. 2. Montrer que ${\cal C}^0([a,b$)] est un espace complet pour la norme $L^\infty$. 3. Montrer que ce n’est pas le cas pour la norme $L^1$ (exhiber une suite de Cauchy non convergente dans ${\cal C}^0([a,b$)]). 4. Démontrer que, lorsqu’elle existe, la dérivée généralisée est unique. 5. Démontrer que, pour une fonction de classe ${\cal C}^1$, la dérivée généralisée est égale à la dérivée classique. 6. Soit une fonction de $[a,b$] vers $\RR$, formée de deux branches de classe ${\cal C}^1$ sur $[a,c[$ et $]c,b$], et discontinue en $c$. Montrer qu’elle n’admet pas de dérivée généralisée. (il faudrait alors avoir recours à la notion de distribution pour dériver cette fonction). 7. Montrer que $|.|_1$ est une norme sur $H^1_0(\Omega)$, équivalente à la norme $\|.\|_1$ # Finite Element Method ## 2. Introduction à la méthode des éléments finis ### 2.1. Formulation variationnelle #### 2.1.1. Exemple 1-D Soit à résoudre le problème (  $({\cal P}) \left\{ \begin{array}{l} -u"(x)+c(x)u(x) = f(x),\quad a où $f$ et $c$ sont des fonctions données continues sur $[a,b$]. On supposera de plus que la fonction $c$ est strictement positive sur $[a,b$]. Un tel problème est appelé problème aux limites. Definition 20. Une solution classique (ou solution forte) de $({\cal P})$ est une fonction de ${\cal C}^2([a,b$)] telle que $u(a)=u(b)=0$ et $\forall x \in$a,b[, \; -u"(x)+c(x)u(x) = f(x)]. En faisant le produit scalaire $L^2(]a,b[)$ de l’équation différentielle avec une fonction-test $v \in {\cal D}(]a,b[)$ (c’est à dire en intégrant sur $[a,b$]), on a :  \[ - \int_a^b u"(x)v(x)\; dx + \int_a^b c(x)u(x)v(x)\; dx = \int_a^b f(x)v(x)\; dx$ soit, en intégrant par parties le premier terme :  $\int_a^b u'(x)v'(x)\; dx + \int_a^b c(x)u(x)v(x)\; dx = \int_a^b f(x)v(x)\; dx$ car $v(a)=v(b)=0$ puisque $v\in{\cal D}(]a,b[)$. Chaque terme de cette équation a en fait un sens dès lors que $v\in H^1_0(]a,b[)$. De plus, ${\cal D}(]a,b[)$ étant dense dans $H^1_0(]a,b[)$ (§[sec:H10]), cette équation est vérifiée pour tout $v\in H^1_0(]a,b[)$. On peut donc définir le nouveau problème :  $({\cal Q})\quad \left\{ \begin{array}{l} \hbox{Trouver }u\in H^1_0(]a,b[) \hbox{ tel que }\\ \ds{ \int_a^b u'(x)v'(x)\; dx + \int_a^b c(x)u(x)v(x)\; dx = \int_a^b f(x)v(x)\; dx \quad\forall v\in H^1_0(]a,b[)} \end{array} \right.$ (2) Ce problème est la formulation variationnelle (ou formulation faible) du problème $({\cal P})$. Toute solution de $({\cal Q})$ est appelée solution faible. Il est immédiat que toute solution forte de $({\cal P})$ est aussi une solution faible. #### 2.1.2. Exemple 2-D Soit $\Omega$ ouvert borné de $\RR^n$. On veut résoudre le problème  $({\cal P})\qquad \left\{ \begin{array}{ll} -\Delta u+u = f & \hbox{ dans }\Omega\\ \partial_n u=0& \hbox{ sur }\partial\Omega\\ \end{array}\right.$ (3) [Une solution classique de ce problème est une fonction de ${\cal C}^2(\bar{\Omega})$ vérifiant ([eq:modele-2D]) en tout point de $\Omega$.  Au passage, on voit que ceci impose que $f$ soit ${\cal C}^0(\bar{\Omega})$. Toute solution classique vérifie donc :  $\forall v\in {\cal C}^1(\bar{\Omega}) \quad \int_\Omega -\Delta u\; v + \int_\Omega uv = \int_\Omega fv$ soit par intégration par parties :  $\ds{ \int_\Omega \nabla u \cdot \nabla v + \int_\Omega uv = \int_\Omega fv }$ puisque $\partial_n u = 0$ sur $\partial\Omega$. Or, $\overline{{\cal C}^1(\bar{\Omega})} = H^1(\Omega)$. On peut donc définir le nouveau problème :  $({\cal Q})\quad \left\{ \begin{array}{l} \hbox{Trouver }u\in H^1(\Omega) \hbox{ tel que }\\ \qquad \ds{ \int_\Omega \nabla u \cdot \nabla v + \int_\Omega u v = \int_\Omega f v \quad\forall v\in H^1(\Omega)} \end{array} \right.$ (4) C’est la formulation variationnelle de $({\cal P})$. On voit aussi que ce problème est défini dès lors que $f\in L^2(\Omega)$. ### 2.2. Formulation générale Les exemples précédents montre que, d’une façon générale, la formulation variationnelle sera obtenue en faisant le produit scalaire $L^2(\Omega)$ de l’équation avec une fonction $v$ appartenant à un espace à préciser (c’est à dire en multipliant par $v$ et en intégrant sur $\Omega$), et en intégrant par parties les termes d’ordre les plus élevés en tenant compte des conditions aux limites du problème. On arrive alors à une formulation du type :  $\hbox{Trouver }u\in V \hbox{ tel que } a(u, v) =l(v) \quad\forall v\in V$ où $a(.,.)$ est une forme sur $V\times V$ (bilinéaire si l’EDP de départ est linéaire) et $l(.)$ est une forme sur $V$ (linéaire si les conditions aux limites de l’EDP de départ le sont). ### 2.3. Théorème de Lax-Milgram #### 2.3.1. Définitions et théorèmes On va introduire ici un outil important pour assurer l’existence et l’unicité de solutions à la formulation variationnelle de problèmes aux limites de type elliptique. Soit $V$ un espace de Hilbert. Definition 21: Continuité d’une forme linéaire Une forme linéaire $l(u)$ sur $V$ est continue ssi il existe une constante $K$ telle que $|l(u)| \le K\, \|u\| \quad \forall u \in V$ Definition 22: Continuité d’une forme bilinéaire Une forme bilinéaire $a(u,v)$ sur $V\times V$ est continue ssi il existe une constante $M$ telle que $|a(u,v)| \le M\, \|u\| \|v\| \quad \forall (u,v) \in V\times V.$ Definition 23: Coercivité d’une forme bilinéaire Une forme bilinéaire $a(u,v)$ sur $V\times V$ est coercive (ou V-elliptique) ssi il existe une constante $\alpha >0$ telle que $a(u,u) \ge \alpha \, \|u\|^2, \; \forall u \in V.$ Theorem 12: Lax-Milgram Soit $V$ un espace de Hilbert. Soit $a$ une forme bilinéaire continue coercive sur $V$. Soit $l$ une forme linéaire continue sur $V$. Alors il existe un unique $u\in V$ tel que $a(u,v)=l(v),\; \forall v\in V.$ De plus, l’application linéaire $l \rightarrow u$ est continue. Proof 1: Lax-Milgram La démonstration générale de ce théorème peut être trouvée par exemple dans cite:[raviart-thomas-1983]. Theorem 13. On prend les mêmes hypothèses que pour le théorème de Lax-Milgram, et on suppose de plus que $a$ est symétrique, c’est à dire que $a(u,v)=a(v,u)\quad\forall u,v$. On définit alors la fonctionnelle $J(v)=\frac{1}{2}\, a(v,v)-l(v)$, et on considère le problème de minimisation : $\hbox{Trouver } u\in V \hbox{ tel que } J(u) = \min_{v\in V} J(v)$ Alors ce problème admet une solution unique, qui est également la solution du problème variationnel précédent. Proof 2. La démonstration de ce théorème vient du fait que $J$ est une fonctionnelle quadratique, et que l’on a $\nabla J[u$(v) = a(u,v) - l(v)].  C’est de cette propriété que vient l’utilisation du terme “variationnel", puisqu’elle montre le lien avec le “calcul des variations". #### 2.3.2. Retour à l’exemple 1-D En reprenant l’exemple 1-D précédent, on peut poser :  $a(u,v) = \int_a^b u'(x)v'(x)\; dx + \int_a^b c(x)u(x)v(x)\; dx$ et  $l(v) = \int_a^b f(x)v(x)\; dx$ $a$ ainsi définie est une forme bilinéaire symétrique continue coercive sur $H^1_0(a,b) \times H^1_0(a,b)$, et $l$ est une forme linéaire continue sur $H^1_0(a,b)$. Donc le problème (2) admet une solution unique d’après le théorème de Lax-Milgram. Cherchons maintenant à interpréter cette solution $u$ de ([eq:FV]). Prenons $v=\varphi \in {\cal D}(]a,b[)$. Alors  $\int_a^b u'(x)\varphi'(x)\; dx + \int_a^b c(x)u(x)\varphi(x)\; dx = \int_a^b f(x)\varphi(x)\; dx$ soit, en intégrant par parties :  $- \int_a^b u"(x)\varphi(x)\; dx + \int_a^b c(x)u(x)\varphi(x)\; dx = \int_a^b f(x)\varphi(x)\; dx$ c’est à dire  $(-u"+cu,\varphi)_0 = (f,\varphi)_0\; \forall \varphi \in {\cal D}(]a,b[).$ ${\cal D}(]a,b[)$ étant dense dans $L^2(]a,b[)$, on a donc : $-u"+cu=f$ dans $L^2(]a,b[)$. $u$ étant dans $L^2(]a,b[)$, et $f$ et $c$ étant dans ${\cal C}^0([a,b])$, donc également dans $L^2(]a,b[)$, on en déduit que $u"=cu-f$ est aussi dans $L^2(]a,b[)$. Puisque $u$ est dans $H^1_0(]a,b[)$ et que $u"$ est dans $L^2(]a,b[)$, on en déduit que $u$ est dans $H^2(]a,b[)$. Donc $u$ est dans ${\cal C}^1([a,b])$ (§[sec:sobolev]). De ce fait, $cu-f$, c’est à dire $u"$, est dans ${\cal C}^0([a,b])$. Donc $u'$ est dans${\cal C}^1([a,b])$, donc $u$ est dans ${\cal C}^2([a,b])$. La solution faible $u$ est donc aussi solution forte du problème de départ. En résumé : • On est parti d’un problème $({\cal P})$ et on a introduit sa formulation variationnelle $({\cal Q})$. • On a montré l’existence et l’unicité d’une solution faible (en utilisant le théorème de Lax-Milgram). Toute solution forte étant aussi solution faible, ceci prouve qu’il y a au plus une solution forte pour $({\cal P})$. • On a prouvé que cette solution faible est bien une solution forte. Le problème de départ $({\cal P})$ a donc une solution unique. L’intérêt de cette démarche est . la formulation variationnelle se prête bien à l’étude de l’existence et de l’unicité de solutions, . on travaille dans des espaces de Hilbert, ce qui va permettre de faire de l’approximation interne. #### 2.3.3. Équations elliptiques d’ordre 2 Soit $\Omega$ un ouvert borné de $\RR^n$, de frontière $\partial\Omega$ assez régulière. Soient des fonctions $\alpha_{ij}$ ($1\le i,j \le n$) dans ${\cal C}^1(\bar{\Omega})$ et $\beta$ dans ${\cal C}^0(\bar{\Omega})$. On considère le problème :  $({\cal P})\qquad \left\{ \begin{array}{rl} {\ds -\sum_{i,j=1}^n \partial_i (\alpha_{ij} \, \partial_j\, u) + \beta\, u = f } & \hbox{ dans }\Omega \\ u= 0 & \hbox{ sur }\Gamma_0 \\ {\ds \sum_{i,j=1}^n \alpha_{ij} \, \partial_j u\; n_i = g } & \hbox{ sur }\Gamma_1 % \end{array}\right.$ (5) où $\Gamma_0$ et $\Gamma_1$ forment une partition de $\partial\Omega$ ($\Gamma_0 \cap\Gamma_1 = \emptyset$ et $\Gamma_0 \cup\Gamma_1 = \partial\Omega$). Une solution classique de $({\cal P})$, sous l’hypothèse que $f\in{\cal C}^0(\bar{\Omega})$ et $g\in{\cal C}^0(\Gamma_1)$, sera une fonction de ${\cal C}^2(\bar{\Omega})$ vérifiant l’équation en chaque point de $\Omega$.La formulation variationnelle de $({\cal P})$ est obtenue par intégration par parties. Elle s’écrit :  $({\cal Q})\quad \left\{ \begin{array}{l} \hbox{Trouver }u\in V \hbox{ tel que }\\ \qquad \ds{ \int_\Omega \left( \sum_{i,j=1}^n \alpha_{ij} \, \partial_j u\; \, \partial_i v + \beta\, u v \right) = \int_\Omega f v + \int_{\Gamma_1} gv \qquad\forall v\in V} \end{array} \right.$ (6) avec $\ds{ V = \left\{ v \in H^1(\Omega) \; , \; v=0 \hbox{ sur }\Gamma_0 \right\} }$. Cette formulation est en fait définie dès lors que $\beta$ et les $\alpha_{ij}$ sont dans $L^\infty(\Omega)$, $f$ dans $L^2(\Omega)$ et $g$ dans $L^2(\Gamma_1)$. Posons  $\ds{a(u,v) = \int_\Omega \left( \sum_{i,j=1}^n \alpha_{ij} \, \partial_j u\; \partial_i v + \beta\, u v \right) }, \quad \ds{l(v) = \int_\Omega f v + \int_{\Gamma_1} gv. }$ Il est immédiat que $a$ est une forme bilinéaire continue et $l$ une forme linéaire continue sur $V$. Si l’EDP de départ 5 vérifie les deux hypothèses d’ellipticité : • il existe $\alpha >0$ tel que $\forall \xi=(\xi_1, \ldots , \xi_n)\in\RR^n$, ${\sum_{i,j=1}^n \alpha_{ij}(x) \, \xi_i \, \xi_j \ge \alpha \, \| \xi \|^2 }$ presque pour tout $x\in\Omega$ • il existe $\beta_0$ tel que $\beta(x) \ge \beta_0$ presque pour tout $x\in\Omega$ alors $a$ est coercive : • sur $H^1_0(\Omega)$ dès que $\ds{\alpha_0 > \frac{-\alpha}{C(\Omega)^2}}$ (et donc en particulier si $\beta\ge 0$) où $C(\Omega)$ est la constante de l’inégalité de Poincaré, voir le théorème [thr:11]. • sur $H^1(\Omega)$ si $\beta > 0$ Par application du théorème de Lax-Milgram, on a donc existence et unicité d’une solution à la formulation variationnelle $({\cal Q})$ : • si $\Gamma_0 = \partial\Omega$ (c’est à dire $\Gamma_1=\emptyset$) et si $\ds{\beta > \frac{-\alpha}{C(\Omega)^2}}$ • si $\Gamma_1\ne \emptyset$ et si $\beta > 0$ ### 2.4. Approximation interne #### 2.4.1. Principe général Soit $\Omega$ un domaine ouvert de $\RR^n$ ($n=1,2$ ou 3 en pratique), de frontière $\partial\Omega$, et sur lequel on cherche à résoudre une équation aux dérivées partielles, munie de conditions aux limites. En écrivant la formulation variationnelle, on obtient un problème de la forme  $({\cal Q})\qquad \hbox{Trouver } u\in V \hbox{ tel que } a(u,v)=l(v), \quad\forall v\in V$ où $V$ est un espace de Hilbert. Sous réserve que l’équation de départ ait de bonnes propriétés, c’est à dire par exemple qu’on soit dans les hypothèses du théorème de Lax-Milgram, $({\cal Q})$ admet une solution unique $u$. Pour obtenir une approximation numérique de $u$, on va maintenant remplacer l’espace $V$ qui est en général de dimension infinie par un sous-espace $V_h$ de dimension finie, et on va chercher à résoudre le problème approché  $\label{eq:6} ({\cal Q}_h)\qquad \hbox{Trouver } u_h\in V_h \hbox{ tel que } a(u_h,v_h)=l(v_h), \quad\forall v_h\in V_h$ $V_h$ étant de dimension finie, c’est un fermé de $V$. $V$ étant un espace de Hilbert, $V_h$ l’est donc aussi. D’où l’existence et l’unicité de $u_h$, à nouveau par exemple d’après le théorème de Lax-Milgram. L’espace $V_h$ sera en pratique construit à partir d’un maillage du domaine $\Omega$, l’indice $h$ désignant la ``taille typique'' des mailles. Lorsque l’on construit des maillages de plus en plus fins, la suite de sous-espaces $(V_h)_h$ formera une approximation interne de $V$, c’est à dire que, pour tout élément $\varphi$ de $V$, il existe une suite de $\varphi_h\in V_h$ telle que $\|\varphi-\varphi_h\|\longrightarrow 0$ quand $h\longrightarrow 0$.  Cette méthode d’approximation interne est également appelée méthode de Galerkin. #### 2.4.2. Interprétation de $u_h$ On a $a(u,v)=l(v), \forall v\in V$, donc en particulier $a(u,v_h)=l(v_h), \forall v_h\in V_h$, car $V_h\subset V$. Par ailleurs, $a(u_h,v_h)=l(v_h), \forall v_h\in V_h$. Par différence, on en déduit que $a(u-u_h,v_h)=0,\quad \forall v_h\in V_h \label{eq:ortho}$  Dans le cas où $a(.,.)$ est symétrique, il s’agit d’un produit scalaire sur $V$. $u_h$ peut alors être interprétée comme la projection orthogonale de $u$ sur $V_h$ au sens de $a(.,.)$. #### 2.4.3. Estimation d’erreur On a :  $\begin{array}{ll} a(u-u_h,u-u_h) & = a(u-u_h,u-v_h+v_h-u_h) \quad\forall v_h\in V_h\\ & =a(u-u_h,u-v_h) + a(u-u_h,v_h-u_h) \end{array}$ Or $v_h-u_h \in V_h$. Donc $a(u-u_h,v_h-u_h)=0$ d’après ([eq:ortho]). On a donc :  $a(u-u_h,u-u_h) = a(u-u_h,u-v_h) \quad\forall v_h\in V_h$ (7) $a$ étant coercive, il existe $\alpha > 0$ tel que $a(u-u_h,u-u_h) \ge \alpha \|u-u_h\|^2$, où $\|.\|$ est une norme sur $V$. Par ailleurs, $a$ étant continue, il existe $M > 0$ tel que $a(u-u_h,u-v_h)\le M \|u-u_h\| \, \|u-v_h\|$. En réinjectant ces deux inégalités de part et d’autre de 7 et en simplifiant par $\|u-u_h\|$ on obtient  $\|u-u_h\| \le \frac{M}{\alpha}\; \|u-v_h\| \quad \forall v_h\in V_h$ (8) c’est à dire  $\|u-u_h\| \le \frac{M}{\alpha}\; d(u,V_h)$ (9) où $d$ est la distance induite par $\|.\|$. Cette majoration est appelée lemme de Céa. Elle ramène l’étude de l’erreur d’approximation $u-u_h$ à l’étude de l’erreur d’interpolation $d(u,V_h)$. ### 2.5. Principe général de la méthode des éléments finis La démarche générale de la méthode des éléments finis est la suivante. On a une EDP à résoudre sur un domaine $\Omega$. On écrit la formulation variationnelle de cette EDP, et on se ramène donc à un problème du type  $({\cal Q})\qquad \hbox{Trouver } u\in V \hbox{ tel que } a(u,v)=l(v), \quad\forall v\in V$ On va chercher une approximation de $u$ par approximation interne. Pour cela, on définit un maillage du domaine $\Omega$, grâce auquel on va définir un espace d’approximation $V_h$, s.e.v. de $V$ de dimension finie $N_h$ (par exemple $V_h$ sera l’ensemble des fonctions continues sur $\Omega$ et affines sur chaque maille). Le problème approché est alors  $({\cal Q}_h)\qquad \hbox{Trouver } u_h\in V_h \hbox{ tel que } a(u_h,v_h)=l(v_h), \quad\forall v_h\in V_h$ Soit $(\varphi_1,\ldots,\varphi_{N_h})$ une base de $V_h$. En décomposant $u_h$ sur cette base sous la forme  $u_h = \sum_{i=1}^{N_h} \mu_i \; \varphi_i$ le problème $({\cal Q}_h)$ devient  $\hbox{Trouver } \mu_1,\ldots,\mu_{N_h} \hbox{ tels que } \sum_{i=1}^{N_h} \mu_i \; a(\varphi_i,v_h)=l(v_h), \quad\forall v_h\in V_h$ ou encore par linéarité de $a$ et $l$ :  $\hbox{Trouver } \mu_1,\ldots,\mu_{N_h} \hbox{ tels que } \sum_{i=1}^{N_h} \mu_i \; a(\varphi_i,\varphi_j)=l(\varphi_j), \quad\forall j=1,\ldots,N_h$ c’est à dire résoudre le système linéaire  $\left( \begin{array}{ccc} a(\varphi_1,\varphi_1) & \cdots & a(\varphi_{N_h},\varphi_1)\\ \vdots & & \vdots\\ a(\varphi_1,\varphi_{N_h}) & \cdots & a(\varphi_{N_h},\varphi_{N_h})\\ \end{array}\right) \left( \begin{array} \mu_1\\ \vdots\\ \mu_{N_h}\\ \end{array}\right) = \left( \begin{array} l(\varphi_1)\\ \vdots\\ l(\varphi_{N_h})\\ \end{array}\right)$ soit  $A\mu = b$ (10)  La matrice $A$ est a priori pleine. Toutefois, pour limiter le volume de calculs, on va définir des fonctions de base $\varphi_i$ dont le support sera petit, c’est à dire que chaque fonction $\varphi_i$ sera nulle partout sauf sur quelques mailles. Ainsi les termes $a(\varphi_i,\varphi_j)$ seront le plus souvent nuls, car correspondant à des fonctions $\varphi_i$ et $\varphi_j$ de supports disjoints. La matrice $A$ sera donc une matrice creuse, et on ordonnera les $\varphi_i$ de telle sorte que $A$ soit à structure bande, avec une largeur de bande la plus faible possible. A ce niveau, les difficultés majeures en pratique sont de trouver les $\varphi_i$ et de les manipuler pour les calculs d’intégrales nécessaires à la construction de $A$. Sans rentrer pour le moment dans les détails, on peut toutefois indiquer que la plupart de ces difficultés seront levées grâce à trois idées principales : Le principe d’unisolvance On s’attachera à trouver des degrés de liberté (ou ddl) tels que la donnée de ces ddl détermine de façon univoque toute fonction de $V_h$. Il pourra s’agir par exemple des valeurs de la fonction en quelques points. Déterminer une fonction reviendra alors à déterminer ses valeurs sur ces ddl. Définition des $\varphi_i$ On définira les fonctions de base par $\varphi_i=1$ sur le $i^{\hbox{\tiny ème}}$ ddl, et $\varphi_i=0$ sur les autres ddl. La manipulation des $\varphi_i$ sera alors très simplifiée, et les $\varphi_i$ auront par ailleurs un support réduit à quelques mailles. La notion de famille affine d’éléments Le maillage sera tel que toutes les mailles soient identiques à une transformation affine près. De ce fait, tous les calculs d’intégrales pourront se ramener à des calculs sur une seule maille de référence, par un simple changement de variable. ### 2.6. Retour à l’exemple 1-D On reprend le problème 1-D 1. On a écrit sa formulation variationnelle [sec:modele-1D> et montré <<sec:modele-1D2] qu’elle admet une solution unique. On s’intéresse à présent à la construction de l’espace d’approximation $V_h$. #### 2.6.1. Construction du maillage La première étape consiste à construire un maillage de $\Omega = ]a,b[$ en définissant une subdivision (pas nécessairement régulière) $a=x_0 < x_1 < \ldots < x_N < x_{N+1}=b$. Le maillage est donc une collection indexée de ($=N$) intervalles \{I_i=[x_{i,1},x_{i,2}]\}_{i=1,…​\Nma} et on a  $[a,b]=\cup_{i=1}^\Nma [x_{1,i},x_{2,i}] \quad \mbox{et} \quad ]x_{1,i},x_{2,i}[ \cap ]x_{1,j},x_{2,j}[ = \emptyset \quad \mbox{ pour } i\neq j$ (11) Definition 24. Les intervalles $I_i$ sont appelées de mailles ou éléments ou cellules du maillage, on a noté $\Nma$ le nombre de maillage Definition 25. Les points $x_i$ sont appelés les sommets du maillage, on note $\Nso=N+1$ le nombre de sommets. On note $h_i = x_{i+1}-x_i$ et $h = \max_{1\leq i \leq \Nma} h_i$. Le maillage est dit uniforme si $h_i=h$ pour tout $i=\{1,...,\Nma\}$. Enfin on note $\calTh=\{I_i\}_{i=\{1,...,\Nma\}}$, $h$ représentant la finesse globale du maillage.  En 1D on a $\Nso = \Nma+1$, en dimension supérieure des relations existent entre le nombre de sommets et de mailles en fonction du type de maille, ce sont les relations d’Euler. #### 2.6.2. Construction de l’espace d’approximation L’étape suivante est de choisir les fonctions de forme ou fonctions de base sur chaque maillage. On choisit les fonctions de $V_h$ telle que leur restriction sur chaque maillage soit un espace polynomial. Definition 26: Espaces $\Pk$ Soit un entier $k \leq 1$. En dimension 1, on appelle l’espace vectoriel des polynômes à coefficients réels de degré inférieur ou égal à $k$. On pose alors  $W_h = \{w_h \in L^2(\Omega); \forall i \in \{ 1,…​,\Nma\}, {w_h}_{|I_i} \in \Pk\}$ (12) $W_h$ est un espace de dimension finie égale à $(k+1)*\Nma$ mais il n’est pas inclus dans $H^1_0(\Omega)$ et ne peut donc pas être utilisé pour l’approximation du problème ([eq:FV]). En effet les fonctions de $w_h \in W_h$ peuvent être discontinues aux interfaces entre les maillages et un résultat d’analyse fonctionnelle montre que dans ces conditions $w_h \ni H^1(\Omega)$. Par ailleurs les fonctions de $W_h$ ne sont pas nécessairement nulles au bord de $\Omega$. On pose donc  $\label{eq:3} V_h = W_h \cap H^1_0(\Omega).$ en d’autres termes, en dimension, on a  $\label{} V_h = \left\{ v_h \in {\cal C}^0 (a,b) \; ; \; {v_h}_{|I_i} \in \Pk \hbox{ et } v_h(a)=v_h(b)=0 \right\}$ (13) Le problème approché sur $V_h$ est :  $({\cal Q}_h)\qquad \hbox{Trouver } u_h\in V_h \hbox{ tel que } a(u_h,v_h)=l(v_h), \quad\forall v_h\in V_h]$ (14) On s’intéresse à présent à des exemples concrets d’espaces d’approximations dans les deux sections suivantes Element fini de Lagrange et Element fini de Lagrange. #### 2.6.3. Element fini de Lagrange On introduit les espaces vectoriels suivants:  $\Pch{1} = \{ v_h \in C^0(\Omega);\; \forall i \in \{ 1,…​,\Nma\} {v_h}_{|I_i} \in \Pk[1] \}$ (15) et  $\Pcho{1} = \{ v_h \in \Pch{1};\; v_h(a)=v_h(b)=0 \}$ (16) Les éléments de ces espaces sont des fonctions continues et affines par morceaux. Ils sont dérivables par morceaux sur chaque maille et ils sont continus aux interfaces entre les mailles. On a le résultat d’analyse fonctionnelle suivant: Theorem 14. $\Pch{1} \subset H^1(\Omega)$ et $\Pcho{1} \subset H^1_0(\Omega)$. On introduit la famille de fonctions $\{\varphi_1,...,\varphi_\Nso\}$ que l’on définit sur chaque maille de la manière suivante, pour tout $i \in \{2,...,\Nso-1\}$,  $\varphi_i(x) = \left\{ \begin{split} \ds{\frac{1}{h_{i-1}} (x-x_{i-1})} & \mbox{ si } x \in I_{i-1}\\ \ds{\frac{1}{h_{i}} (x_{i+1}-x)} & \mbox{ si } x \in I_{i}\\ 0 & \mbox{ sinon}, \end{split} \right.$ (17) et  $\begin{split} \varphi_1(x) &= \left\{ \begin{split} \ds{\frac{1}{h_{1}} (x_2-x)} & \mbox{ si } x \in I_{1}\\ 0 & \mbox{ sinon}, \end{split} \right.\\ \varphi_\Nso(x) &= \left\{ \begin{split} \ds{\frac{1}{h_{\Nso-1}} (x-x_{\Nso-1})} & \mbox{ si } x \in I_{\Nso-1}\\ 0 & \mbox{ sinon}, \end{split} \right. \end{split}$ (18)  Les fonctions $(\varphi_i)_{i=1,...,\Nso}$ sont dans $\Pch{1}$ et $(\varphi_i)_{i=2,...,\Nso-1}$ sont dans $\Pcho{1}$.  Les fonctions $(\varphi_i)_{i=1,...,\Nso}$ satisfont les relations  $\varphi_i(x_j) = \delta_{ij},\quad i,j \in \{1,…​,\Nso\},$ où $\delta_{ij}$ désigne le symbole de Kronecker tel que $\delta_{ij} = 1$ si $i=j$ et $\delta_{ij}=0$ si $i \neq j$. Les fonctions $\varphi_i$ sont appelées fonctions chapeau du fait de leur graphe, voir figure [fig:chapeau]. Proposition 1. . La famille $\{\varphi_1,...,\varphi_\Nso\}$ est une base de $\Pch{1}$. . La famille $\{\varphi_2,...,\varphi_{\Nso-1}\}$ est une base de $\Pcho{1}$. Corollary 1. $\dim \Pch{1} = \Nso = \Nma+1$ et $\dim \Pcho{1} = \Nso-2 = \Nma-1$. On introduit l’opérateur d’interpolation suivant:  $\Ich{1} : \Ck{0}(\bar{\Omega}) \ni v \mapsto \sum_{i=1}^\Nso v(x_i) \varphi_i \in \Pch{1}.$ (19) Pour toute fonction $v \in \Ck{0}(\bar{\Omega})$, $\Ich{1}{v}$ est la seule fonction continue affine par morceaux prenant les mêmes valeurs que $v$ aux sommets $x_i, i=1,...,\Nso$. $\Ich{1}{v}$ est appelée l’interpolé de Lagrange de $v$ de degré $1$.  En dimension 1, les fonctions de $H^1(\Omega)$ sont continues, on peut donc voir comme un opérateur de $H^1(\Omega)$ dans $H^1(\Omega)$. On montre qu’il est continu et que sa norme $\|\Ich{1}\|_{\mathcal{L}(H^1(\Omega),H^1(\Omega))}$ est uniformément bornée en $h$, c’est à dire qu’il existe une constante $c$, indépendante de $h$, telle que pour tout $v \in H^1(\Omega)$  $\|\Ich{1} v \|{1,\Omega} \leq c \|v\|{1,\Omega}]$ (20) #### 2.6.4. Estimation de l’erreur d’interpolation Proposition 2. Pour tout $h$, et tout $v \in H^2(\Omega)$, on a $\|v - \Ich{1} v\|_{0,\Omega} \leq h^2 |v|_{2,\Omega}\quad \mbox{ et }\quad |v - \Ich{1} v|_{1,\Omega} \leq h |v|_{2,\Omega}$ On dit que l’erreur d’interpolation est d’ordre 2 en norme $L^2$ et d’ordre 1 en semi-norme $H^1$ et donc en norme $H^1$. #### 2.6.5. Element fini de Lagrange On introduit les espaces vectoriels suivants: stem:[\label{eq:9}  $\Pch{k} = \{ v_h \in C^0(\Omega);\; \forall i \in \{ 1,…​,\Nma\}, {v_h}_{|I_i} \in \Pk\}$ (21) et  $\Pcho{k} = \{ v_h \in \Pch;\; v_h(a)=v_h(b)=0 \}$ (22) #### 2.6.6. Operateur d’interpolation On introduit l’opérateur d’interpolation suivant: $\label{eq:24} \Ich : \Ck{0}(\bar{\Omega}) \ni v \mapsto \sum_{i=1}^\Nno v(x_i) \varphi_i \in \Pch.$ Pour toute fonction $v \in \Ck{0}(\bar{\Omega})$, $\Ich v$ est la seule fonction continue polynomial de degré $k$ par morceaux prenant les mêmes valeurs que $v$ aux sommets $x_i, i=1,...,\Nso$. $\Ich v$ est appelée l’interpolé de Lagrange de $v$ de degré $k$.  En dimension 1, les fonctions de $H^1(\Omega)$ sont continues, on peut donc voir comme un opérateur de $H^1(\Omega)$ dans $H^1(\Omega)$. On montre qu’il est continu et que sa norme $\|\Ich\|_{\mathcal{L}(H^1(\Omega),H^1(\Omega))}$ est uniformément bornée en $h$, c’est à dire qu’il existe une constante $c$, indépendante de $h$ mais dépendante de $k$, telle que pour tout  $v \in H^1(\Omega) \|\Ich{k}{v} \|{1,\Omega} \leq c \|v\|{1,\Omega}$ (23) Le résultat suivant permet d’estimer la précision de l’opérateur d’interpolation, Proposition 3. Il existe une constante $c$, indépendante de $h$ mais dépendante de $k$, telle que pour tout $h$ et pour tout $v \in H^{k+1}(\Omega)$, on a $\|v - \Ich{k}{v}\|_{0,\Omega} + h |v - \Ich{k}{v}|_{1,\Omega} \leq c\; h^{k+1}\; |v|_{k+1,\Omega}$ et $\sum_{m=2}^{k+1} h^m \left( \sum_{i=0}^N |v - \Ich{k}{v}|^2_{m,I_i}\right)^{1/2} \leq c\; h^{k+1}\; |v|_{k+1,\Omega}]$  L’estimation [eq montre que l’erreur d’interpolation est d’ordre $k+1$ en norme $\|\cdot\|_{0,\Omega}$ et qu’elle est d’ordre $k$ en norme $|\cdot|_{1,\Omega}$. Elle est donc d’ordre $k$ en norme $\|\cdot\|_{1,\Omega}$. #### 2.6.7. Analyse de convergence Nous nous intéressons à présent à l’analyse de la convergence de $u_h$ du problème approché de 14 vers la solution $u$ du problème exact [eq:F] lorsque $V_h=\Pcho{1}$ ou plus généralement $V_h=\Pcho{k},\; k\geq 1$. ##### Estimation en norme $H^1$ Il s’agit dans un premier temps d’estimer l’erreur $u-u_h$ en norme $H^1$. Pour cela on part de l’estimation 8, on a  \begin{align} \|u-u_h\|{1,\Omega} &\leq c\; \inf{v_h \in \Pcho} \|u-v_h\|{1,\Omega}\\ & \leq c\; \|u-\Ich{u}\|{1,\Omega}\\ & \leq c\; h^k |u|_{k+1,\Omega} \end{align} (24) pourvu que la solution exacte soit suffisamment régulière, c’est à dire $u \in H^{k+1}(\Omega)$.  On notera que $\Ich{k}{u} \in \Pcho{k}$ puisque $u \in H^1_0(\Omega)$ et donc que $\Ich{u}(a)=\Ich{u}(b)=0$. On a donc le résultat suivant Proposition 4. Soit un entier $k\geq 1$. On suppose que la solution du problème 2 est dans $H^{k+1}(\Omega)$. On note $u_h$ la solution du problème approché 14 avec l’espace d’approximation $V_h = \Pcho{k}$. Alors, il existe une constante $c$, indépendante de $h$, telle que $\|u-u_h\|_{1,\Omega} \leq c\; h^k |u|_{k+1,\Omega}$  On dit que l’estimation d’erreur [eq est optimale car elle est du même ordre que l’erreur d’interpolation en norme $H^1$, voir la proposition Proposition 3. #### 2.6.8. Estimation en norme $L^2$ Proposition 5. Avec les hypothèses de la proposition [prop:1] et en supposant que $\alpha \in \Ck{1}(\bar{\Omega})$. Alors, il existe une constante $c$, indépendante de $h$, telle que $\|u-u_h\|_{0,\Omega} \leq c\; h^{k+1} |u|_{k+1,\Omega}$  On dit que l’estimation d’erreur [eq est optimale car elle est du même ordre que l’erreur d’interpolation en norme $L^2$, voir la proposition Proposition 3. ### 2.7. Formulation algébrique $V_h=P_{c,h}^1$ En décomposant la solution approchée $u_h$ sur cette base sous la forme ${u_h = \sum_{i=1}^N \mu_i \; \varphi_i}$, on obtient, comme au paragraphe Principe général de la méthode des éléments finis, le système linéaire $A\mu=b$, avec :  $\begin{array}{rcl} A_{ji}=a(\varphi_i,\varphi_j) & = & \ds{\int_a^b \left[ \varphi_i'(x) \varphi_j'(x) + c(x) \varphi_i(x) \varphi_j(x)\right]\; dx }\\ & = & \ds{ \sum_{k=0}^N \int_{x_k}^{x_{k+1}} \left[\varphi_i'(x) \varphi_j'(x) + c(x) \varphi_i(x) \varphi_j(x)\right]\; dx } \end{array}$ Le support de $\varphi_i$ étant réduit à $[x_{i-1},x_{i+1}$], on en déduit que  $\left\{ \begin{array}{lll} a(\varphi_i,\varphi_j) & = & 0 \qquad \hbox{si }|i-j|\ge 2\\ & & \\ a(\varphi_i,\varphi_{i+1}) & = & \ds{ \int_{x_i}^{x_{i+1}} \left[ \varphi_i'(x) \varphi_{i+1}'(x) + c(x) \varphi_i(x) \varphi_{i+1}(x)\right] \; dx}\\ & & \\ a(\varphi_i,\varphi_{i-1}) & = & \ds{ \int_{x_{i-1}}^{x_i} \left[ \varphi_i'(x) \varphi_{i-1}'(x) + c(x) \varphi_i(x) \varphi_{i-1}(x)\right] \; dx}\\ & & \\ a(\varphi_i,\varphi_{i}) & = & \ds{ \int_{x_{i-1}}^{x_{i+1}} \left[ \varphi_i'^2(x) + c(x) \varphi_i^2(x)\right] \; dx}\\ \end{array}\right.$ $A$ est donc tridiagonale. #### 2.7.1. Exercices 1. Dans le Définitions et théorèmes, montrer que, dans le cas où $a$ est symétrique, si $u$ est solution du problème variationnel, alors elle est solution du problème de minimisation. 2. Montrer que $\nabla J[u$(v) = a(u,v) - l(v)]. 3. Montrer que, si $a$ est coercive, la matrice $A$ de ([eq:lin]) est inversible. (C’est donc la démonstration du théorème de Lax-Milgram en dimension finie.) 4. Pour l’exemple 1-D traité dans ce chapitre, démontrer qu’on est bien dans les hypothèses du théorème de Lax-milgram 5. Calculer explicitement la matrice $A$ pour cet exemple. 6. Pour le problème 2-D du §[sec:modele-2D], montrer que la formulation variationelle ([eq:FV2]) admet une solution unique, qui est aussi solution classique si $f \in H^2(\Omega)$. 7. Démontrer les résultats du §[sec:elliptique] ## 3. Eléments finis de Lagrange On va présenter dans ce chapitre le type le plus simple et le plus classique d’éléments finis. ### 3.1. Unisolvance Definition 27. Soit $\Sigma=\{ a_1,\ldots,a_N\}$ un ensemble de $N$ points distincts de $\RR^n$. Soit $P$ un espace vectoriel de dimension finie de fonctions de $\RR^n$ à valeurs dans $\RR$.On dit que $\Sigma$ est $P$-unisolvant ssi pour tous réels $\alpha_1,\ldots,\alpha_N$, il existe un unique élément $p$ de $P$ tel que $p(a_i)=\alpha_i,\; i=1,\ldots,N$. Ceci revient à dire que la fonction :  $\begin{array}{rcl} {\cal L} : P & \longrightarrow & \RR^N\\ p & \longrightarrow & (p(a_1),\ldots,p(a_N)) \end{array}$ (25) est bijective. Proof 3. En pratique, on montrera que $\Sigma$ est $P$-unisolvant en vérifiant que $\opdim P= Card \Sigma$, puis en montrant l’injectivité ou la surjectivité de ${\cal L}$. L’injectivité de ${\cal L}$ se démontre en établissant que la seule fonction de $P$ s’annulant sur tous les points de $\Sigma$ est la fonction nulle. La surjectivité de ${\cal L}$ se démontre en exhibant une famille $p_1,\ldots,p_N$ d’éléments de $P$ tels que $p_i(a_j)=\delta_{ij}$, c’est à dire un antécédent pour ${\cal L}$ de la base canonique de $\RR^N$. En effet, étant donnés des réels $\alpha_1,\ldots,\alpha_N$, la fonction ${p=\sum_{i=1}^N \alpha_i p_i}$ vérifie alors $p(a_j)=\alpha_j,\; j=1,\ldots,N$. ### 3.2. Elément fini de Lagrange Definition 28: Élément fini de Lagange Un élément fini de Lagrange est un triplet $(K,\Sigma,P)$ tel que * $K$ est un élément géométrique de $\RR^n$ ($n=$ 1, 2 ou 3), compact, connexe, et d’intérieur non vide. * $\Sigma=\{a_1,\ldots,a_N\}$ est un ensemble fini de $N$ points distincts de $K$. * $P$ est un espace vectoriel de dimension finie de fonctions réelles définies sur $K$, et tel que $\Sigma$ soit $P$-unisolvant (donc dim $P = N$). Definition 29: Fonctions de base locales Soit $(K,\Sigma,P)$ un élément fini de Lagrange. On appelle fonctions de base locales de l’élément les $N$ fonctions $p_i$ ($i=1,\ldots,N$) de $P$ telles que $p_i(a_j)=\delta_{ij}\qquad 1\le i,j \le N.$ On vérifie aisément que $(p_1,\ldots,p_N)$ ainsi définie forme bien une base de $P$. Definition 30. On appelle opérateur de interpolation (ou encore $P-$interpolation) sur $\Sigma$ l’opérateur $\pi_K$ qui, à toute fonction $v$ définie sur $K$, associe la fonction $\pi_K v$ de $P$ définie par ${\pi_K v = \sum_{i=1}^N v(a_i)\, p_i}$. $\pi_K v$ est donc l’unique élément de $P$ qui prend les mêmes valeurs que $v$ sur les points de $\Sigma$. ### 3.3. Exemples d’éléments finis de Lagrange #### 3.3.1. Espaces de polynômes On notera $\Pk$ l’espace vectoriel des polynômes de degré total inférieur ou égal à $k$. • Sur , $\Pk=\hbox{Vect} \{ 1,X,\ldots, X^k\}\quad$ et dim $\Pk=k+1$. • Sur $^2$, $\Pk=\hbox{Vect} \{ X^i Y^j ; 0\le i+j \le k\}\quad$ et dim ${ \Pk=\frac{(k+1)(k+2)}{2}}$. • Sur $^3$, $\Pk=\hbox{Vect} \{ X^i Y^j Z^l ; 0\le i+j+l \le k\}\quad$ et dim ${ \Pk=\frac{(k+1)(k+2)(k+3)}{6}}$. On notera $\Qk$ l’espace vectoriel des polynômes de degré inférieur ou égal à $k$ par rapport à chaque variable. • Sur $\RR$, $\Qk=\Pk$. • Sur $\RR^2$, $\Qk=\hbox{Vect} \{ X^i Y^j ; 0\le i,j \le k\}\quad$ et dim ${ \Qk=(k+1)^2}$. • Sur $\RR^3$, $\Qk=\hbox{Vect} \{ X^i Y^j Z^l ; 0\le i,j,l \le k\}\quad$ et dim ${ \Qk=(k+1)^3}$. #### 3.3.2. Exemples 1-D Elément $P_1$ • $K=[a,b$] • $\Sigma=\{a,b\}$ • $P=\Pk[1$] Elément $P_2$ • $K=[a,b$] • ${\Sigma=\{a,\frac{a+b}{2},b\}}$ • $P=\Pk[2$] Elément $\Pk$ • $K=[a,b$] • ${\Sigma=\{a+i\,\frac{b-a}{m},\quad i=0,\ldots,k\}}$ • $P=\Pk$ #### 3.3.3. Exemples 2-D triangulaires Elément $\Pk[1$] • $K$=triangle de sommets $a_1, a_2, a_3$ • $\Sigma=\{a_1,a_2,a_3\}$ • $P=\Pk[1$] Les fonctions de base sont définies par $p_i(a_j)=\delta_{ij}$. Ce sont donc les coordonnées barycentriques : $p_i=\lambda_i$ (cf annexe [ch:bary]). Elément $\Pk[2$] • $K$=triangle de sommets $a_1, a_2, a_3$ • $\Sigma=\{a_1,a_2,a_3, a_{12}, a_{13}, a_{23}\}$, où ${a_{ij}=\frac{a_i+a_j}{2}}$. • $P=\Pk[2$] Les fonctions de base sont $p_i=\lambda_i (2\lambda_i -1)$ et $p_{ij}=4\lambda_i\lambda_j$. Un exemple de calcul de ces fonctions de base est donné en annexe [ch:bary]. Eléments finis triangulaire $P_1$, triangulaire $P_2$ et rectangulaire $Q_1$ #### 3.3.4. Exemples 2-D rectangulaires Elément $\Qk[1$] • $K$=rectangle de sommets $a_1, a_2, a_3, a_4$, de côtés parallèles aux axes • $\Sigma=\{a_1,a_2,a_3,a_4\}$ • $P=\Qk[1$] Les fonctions de base sont ${p_i(X,Y)=\frac{(X-x_j)(Y-y_j)}{(x_i-x_j)(y_i-y_j)}}$, où $(x_i,y_i)$ sont les coordonnées de $a_i$, et où $a_j$, de coordonnées $(x_j,y_j)$ est le coin opposé à $a_i$. #### 3.3.5. Exemples 3-D Elément tétraèdrique $\Pk[1$] • $K$=tétraèdre de sommets $a_1, a_2, a_3, a_4$ • $\Sigma=\{a_1,a_2,a_3, a_4\}$ • $P=\Pk[1$] Elément tétraèdrique $\Pk[2$] • $K$=tétraèdre de sommets $a_1, a_2, a_3, a_4$ • ${ \Sigma=\{a_i\}_{1\le i\le 4} \cup \{a_{ij}\}_{1\le i < j \le 4} }$ • $P=\Pk[2$] Les fonctions de base sont $p_i=\lambda_i (2\lambda_i -1)$ et $p_{ij}=4\lambda_i\lambda_j$. Elément parallélépipèdique $Q_1$ • $K$=parallélépipède de sommets $a_1, \ldots , a_8$ de côtés parallèles aux axes • $\Sigma=\{a_i\}_{1\le i\le 8}$ • $P=\Qk[1$] Elément prismatique • $K$=prisme droit de sommets $a_1, \ldots , a_6$ • $\Sigma=\{a_i\}_{1\le i\le 6}$ • $P=\{p(X,Y,Z)=(a+bX+cY)+Z(d+eX+fY), \;\; a,b,c,d,e,f \in \RR\}$ Eléments finis tétraèdriques $P_1$ et $P_2$, parallélépipèdique $Q_1$, et prismatique ### 3.4. Famille affine d’éléments finis Definition 31: Équivalence affine Deux éléments finis $(\hat{K},\hat{\Sigma},\hat{P})$ et $(K,\Sigma,P)$ sont affine-équivalents ssi il existe une fonction affine $F$ inversible ($F: \hat{x} \longrightarrow B\hat{x}+b$) telle que (i) $K=F(\hat{K})$ (ii) $a_i=F(\hat{a}_i) \qquad i=1,\ldots,N$ et (iii) $P=\{ \hat{p}\circ F^{-1} , \quad \hat{p}\in \hat{P} \}$  Si l’on est dans $\RR^n$, $B$ est donc une matrice $n\times n$ inversible, et $b$ est un vecteur de $\RR^n$. Property 1. Soient $(\hat{K},\hat{\Sigma},\hat{P})$ et $(K,\Sigma,P)$ deux éléments finis affine-équivalents, via une transformation $F$. On note $\hat{p}_i \; (i=1,\ldots,N)$ les fonctions de base locales de $\hat{K}$. Alors les fonctions de base locales de $K$ sont les $p_i=\hat{p}_i\circ F^{-1}$. Definition 32. On appelle famille affine d’éléments finis une famille d’éléments finis tous affine-équivalents à un même élément fini $(\hat{K},\hat{\Sigma},\hat{P})$, appelé élément de référence. D’un point de vue pratique, le fait de travailler avec une famille affine d’éléments finis permet de ramener tous les calculs d’intégrales à des calculs sur l’élément de référence. Les éléments de référence sont : • En 1-D : le segment $[0,1]$ • En 2-D triangulaire : le triangle unité, de sommets $(0,0)$, $(0,1)$ et $(1,0)$. • En 2-D rectangulaire : le carré unité $[0,1]\times[0,1]$. • En 3-D tétraèdrique : le tétraèdre unité, de sommets $(0,0,0)$, $(1,0,0)$, $(0,1,0)$ et $(0,0,1)$. • En 3-D parallélépipèdique : le cube unité $[0,1]\times[0,1]\times[0,1]$. • En 3-D prismatique : le prisme unité de sommets $(0,0,0)$, $(0,1,0)$, $(1,0,0)$, et $(0,0,1)$, $(0,1,1)$, $(1,0,1)$. ### 3.5. Maillages Nous étendons ici aux dimension 2 et 3 les notions élémentaires de maillage vues en 1D, voir la figure Maillage Feel++. Maillage Feel++ Definition 33. Un maillage est constituée d’une famille d’éléments(ou mailles ou cellules) $\{K_e\}_{e=1,...,N_e}$ où $N_e$ est le nombre d’éléments, nous noterons $\calTh = \{K_m\}_{m=1,...,N_e}$ avec $h=\max_{1\le e\le N_e} h_{K_e}$ et $h_{K_e} = \diam(K_e) = \max_{x_1,x_2 \in K_e} \|x_1-x_2\|,\, e \in \{1,...,\Ne\}$ On travaille par la suite avec des familles de maillage et on les note $\set{\mathcal{T}_h}_{h > 0}$. Definition 34: Famille de maillage quasi-uniforme On dira qu’une famille de maillage $\set{\mathcal{T}_h}_{h > 0}$ est quasi-uniforme s’il existe une constante $c$ telle que $\forall h,\ \forall K \in \calTh,\ h_K \geq c h$  Cela veut dire que les élements sont tous de la même taille pour $h$ donné. #### 3.5.1. Transformation géométrique Un maillage est généré par 1. un élément de reference noté $\hat{K}$ 2. une famille de transformations géométriques mappant $\hat{K}$ vers les éléments $K_e, e=1,\ldots,\Ne$ dans le maillage Nous supposerons que les transformations sont des $\mathcal{C}^1-$ diffeomorphismes [1]. Definition 35. Pour une cellule $K \in \mathcal{T}_h$, on note $T_K$ la transformation géométrique $T_K: \hat{K} \mapsto K$ Afin de spécifier la transformation géométrique, on considère l’élément fini de Lagrange, noté $(\hat{K},\hat{P}_{\mathrm{geo}}, \hat{\Sigma}_{\mathrm{geo}})$, tel que • $\ngeo = \card{\hat{\Sigma}_{\mathrm{geo}}}$ • $\set{\hat{g}_1,\ldots,\hat{g}_{\ngeo}}$ les noeuds de $\hat{K}$ • $\set{\hat{\psi}_1,\ldots,\hat{\psi}_{\ngeo}}$ les fonctions de forme Definition 36: Élément fini géométrique On dit que $(\hat{K},\hat{P}_{\mathrm{geo}}, \hat{\Sigma}_{\geo})$ est l’élément fini géométrique, $\set{\hat{g}_1,\ldots,\hat{g}_{\ngeo}}$ sont les noeuds géométriques et $\set{\hat{\psi}_1,\ldots,\hat{\psi}_{\ngeo}}$ sont les fonctions de formes géométriques Transformation géométrique associée à un triangle Pour chaque $K \in \mathcal{T}_h$, on a un $\ngeo$-uplet $\set{g^K_1,\ldots,g^K_\ngeo}$. La transformation géométrique est définie comme suit  $T_K: \hat{x} \in \hat{K} \mapsto \sum_{i=1}^\ngeo\ g^K_i \hat{\psi}_i(\hat{x})$ (26) et en particulier on a  $T_K(\hat{g}_i) = g^K_i, \quad \forall i \in \set{1,\ldots,\ngeo}$ (27)  On a $T_K \in [\hat{P}_\geo(\hat{K})$^d] et que $\set{g^K_1,\ldots,g^K_\ngeo}$ sont les noeuds géométriques de $K$. $T_K$ est un $\mathcal{C}^1$-diffeomorphism donc la numérotation des noeuds $\set{g^K_1,\ldots,g^K_\ngeo}$ doit être compatible avec les noeuds de l’élément finit géométrique.  La numérotation locale des entités géométriques dans doit être consistente avec la numérotation locale des générateurs de maillage. voir $\triangleright$ format de fichier Gmsh pour une numérotation locale. Un cas particulier est la transformation géométrique affine. Definition 37: Maillage affine Quand toutes les transformations géométriques $\set{T_K}_{K \in \mathcal{T}_h}$ sont affines, cela veut dire que pour tout $K \in \mathcal{T}_h$, il existe un vecteur $b_K \in \RR^d$ et une matrice $J_K \in \RR^{d\times d}$ tels que $T_K : \hat{x} \in \hat{K} \mapsto b_K + J_K \hat{x} \in K$ On dit que le maillage est affine.  Si l’élément fini géométrique est $(\hat{K},\poly{P}_1,\Sigma_\ngeo)$ alors les éléments $K$ sont soit des triangles soit des tétrahèdres.  `Mesh >` ou `Mesh >` est le type pour les maillages affines formés de simplexes dans $\RR^d$. `1` indique l’ordre de l’élément fini géométrique et est la valeur par défaut. #### 3.5.2. Quelques calculs avec la transformation géométrique ##### Gradient, Inverse et Jacobien On note $\xi$ un ensemble de $n$ points dans $\hat{K}$ et on note $\nabla T_K(\xi)$ le gradient de $T_K$ aux points $\xi$  $\nabla T_K( \xi )\ =\ \sum_{i=0}^{\ngeo}\ g^K_i\ \nabla \psi_i (\xi)$ (28) et $B_K(\xi) = \nabla T_K^{-1}(\xi)$ l’inverse $\xi$ et finalement $J_K(\xi)$ le jacobien de $T_K$ en $\xi$  $J_K(\xi)\ =\ |\det( \nabla T_K(\xi) )|$ (29) ##### Dérivation dans l’élément de référence Afin de dériver un polynome dans l’élément réel $K$, grâce à la transformation géométrique et la règle de différentiation des fonctions composées, nous dérivons seulement dans l’élément de référence $\hat{K}$. Soit $f: K \mapsto \RR$ et $\hat{f}: \hat{K} \mapsto \RR$ telle que $\hat{f} = f \circ T_K$  $\nabla f\ =\ \hat{\nabla} \underbrace{\hat{f}(\xi)}_{f \circ T_K(\xi)} B_K(\xi)$ (30) en 2D, on a, en notant $x=T_K(\xi)$,  $\nabla f(x) = \begin{pmatrix} \frac{\hat{\partial} \hat{f} (\xi)}{\partial \xi_1} & \frac{\hat{\partial} \hat{f} (\xi)}{\partial \xi_2} \end{pmatrix} \begin{pmatrix} B_{K_{11}}(\xi) & B_{K_{12}}(\xi)\\ B_{K_{21}}(\xi) & B_{K_{22}}(\xi)\\ \end{pmatrix}$ (31) ##### Intégration dans l’élément de référence De manière similaire, au lieu de calculer les intégrales sur l’élément réel $K$, nous appliquons un changement de variables et calculons les intégrales sur l’élément de réference $\hat{K}$. Soit $f: K \mapsto \RR$ et $\hat{f}: \hat{K} \mapsto \RR$ telle que $\hat{f} = f \circ T_K$, et ${\mathbf{F}}: K \mapsto \RR^d$ et ${\hat{\mathbf{F}}}: \hat{K} \mapsto \RR^d$ telle que $\hatglossterm::[\mathbf{F$} = {\mathbf{F}} \circ T_K], On a alors les relations suivantes:  $\int_{K} \ f\ dx\ =\ \int_{\hat{K}} f(T_K(\xi) ) J_K( \xi )\ d \xi \ =\ \int_{\hat{K}} \hat{f}(\xi) J_K( \xi )\ d \xi]$ (32)  $\int_{K}\ \nabla f\ dx\ =\ \int_{\hat{K}} \Big(\hat{\nabla} \underbrace{\hat{f}(\xi)}_{f \circ T_K(\xi)} B_K(\xi)\Big) J_K( \xi )\ d \xi$ (33)  $\int_{\partial K} f( x ) dx = \int_{\partial \hat{K}} \hat{f}(\xi) \| B_K(\xi) {\hat{\mathbf}}(\xi) \| J_K( \xi ) d \xi$ (34)  $\int_{\partial K}\ {\mathbf{F}}( x )\ \cdot\ {\mathbf}(x) dx = \int_{\partial \hat{K}} {\hat{\mathbf{F}}}( \xi )\ \cdot \Big(B_K(\xi) {\hat{\mathbf}}(\xi) \Big) J_K( \xi )\ d \xi$ (35) où ${\mathbf{n}}(x)$ est la normale extérieure unitaire à $\partial K$ évaluée en $x \in \partial K$, et ${\hat{\mathbf{n}}}(\xi)$ la normale unitaire extérieure à $\hat{K}$ évaluée en $\xi \in \partial \hat{K}$.  Feel++ effectue automatiquement pour vous les changements de variables dans les intégrales. Definition 38: Maillage conforme Un maillage est dit conforme si l’intersection de deux éléments est soit vide, un sommet, une arête ou une face.  On ne manipule que des maillages conformes dans le cours mais Feel++ peut traiter des maillages non-conformes par exemple dans le contexte de méthode de décomposition de domaines. ### 3.6. Espaces élément fini de Lagrange Soit $\mathcal{T}_h$ un maillage généré par $(\hat{K}, \hat{P}_{\mathrm{geo}}, \hat{\Sigma}_{\mathrm{geo}})$, une cellule $K \in \mathcal{T}_h$ est alors l’image de $\hat{K}$ par la transformation géométrique $T_K$ défini par [eq. L’objectif à présent est de générer la famille d’éléments finis de Lagrange grâce à l’élément fini de référence $(\hat{K},\hat{P}, \hat{\Sigma})$  $\{(K,P_K,\Sigma_K)\}_{K \in \mathcal{T}_h}$ (36) Notation 1. On note $\{\hat{x}_1,...,\hat{x}_{\nf}\}$ les noeuds de l’élément fini. On note $\{\hat{\Psi}_1,...,\hat{\Psi}_{\nf}\}$ les fonctions de forme élément fini. Proposition 6. Soit $K \in \mathcal{T}_h$ et $P_K=\{\hat{p} \circ T_K^{-1}; \hat{p} \in \hat{P}\}$. Pour tout $i \in \{1,...,\nf\}$, on pose $x_{K,i} = T_K(\hat{x}_i)$ On note $\Sigma_K$ l’ensemble des degrés de liberté associé à $\{x_{K,1},...,x_{K,\nf}\}$ Alors $(K,P_K,\Sigma_K)$ est un élément fini de Lagrange. Les fonctions de forme sont définies de la façon suivante $\Psi_{K,i} = \hat{\Psi}_i \circ T_K^{-1}, \quad i=1,...,\nf$ et l’opérateur d’interpolation local comme $\Ilag{K}: v \in \mathcal{C}^0(K) \mapsto \sum_{i=1}^\nf\ v(x_{K,i})\ \Psi_{K,i}\ \in P_K$ Une propriété importante de $\Ilag{K}$ est que $\forall v \in \mathcal{C}^0(K),\quad \Ilag{K}( v \circ T_K ) = \Ilag{K}( v ) \circ T_K.$ Theorem 15. Soit $T_K$ une transformation affine. Soit $\mathbb{P}_k \subset \hat{P}$ et $k+1 > \frac{d}{2}$. Soit $h_K$ le diamètre de $K$ et $\rho_K$ le diamètre de la plus grande boule inscrite dans $K$ et $\omega_K = \frac{h_K}{\rho_K}$. Alors il existe une constante $c$ independente de $K$ telle que $\forall v \in H^{k+1}(K)$ et pour tout $m \in \{0,...,k+1\}$, $|v - \Ilag{K}(v)|_{m,K} \leq c h^{k+1-m} \omega_K^m |v|_{k+1,K}$  $\omega_K$ devrait être aussi proche de $1$ que possible La deuxième hypothèse technique permet d’assurer que $H^{k+1}(K) \subset C^0(K)$. on obtient des résultats similaires si $v$ n’est pas suffisamment régulière Definition 39: Espace $H^1$-conforme Un espace vectoriel $V_h$ de fonctions définies sur un domaine $\Omega_h$ est dit être $H^1$-conforme si $V_h \subset H^1(\Omega_h)$. Afin de construire un tel espace on introduit tout d’abord  $W_h = \{v_h \in L^2(\Omega_h); \forall K \in \mathcal{T}_h, v_h|_K \in P_K\}$ (37) mais ce n’est pas suffisant: les fonctions $W_h$ peuvent avoir des sauts entre les éléments du maillage. Nous avons donc besoin d’assurer la continuité de ces fonctions  $V_h = W_h \cap C^0(\Omega_h) = \{ v_h \in W_h; \forall F \in \mathcal{F}^i_h. \jump{v_h}_F = 0\}$ (38) Concernant l’implémentation, nous avons besoin de d’indentifier les degrés de liberté communs entre les éléments quand nous construisons la tables des degrés de liberté. Voici deux exemples d’espace $H^1$-conforme  $P^k_{c,h} = \{ v_h \in C^0(\Omega_h); \forall K \in \mathcal{T}_h, v_h \circ T_K \in \mathbb{P}_k\}$ (39)  $Q^k_{c,h} = \{ v_h \in C^0(\Omega_h); \forall K \in \mathcal{T}_h, v_h \circ T_K \in \mathbb{Q}_k\}$ (40) #### 3.6.1. Projections orthogonales On note  \begin{aligned} \Pi^{0,k}_{c,h} : Lˆ2(\Omega) \rightarrow P^k_{c,h}\\ \Pi^{1,k}_{c,h} : H^1(\Omega) \rightarrow P^k_{c,h} \end{aligned} (41) les projections orthogonales associées respectivement aux produits scalaires $(u,v)_{0,\Omega} = \int_\Omega u v$ and $(u,v)_{1,\Omega} = \int_\Omega u v + \int_\Omega \nabla u \cdot \nabla v$ On a pour $l=1,...,k$ et si $v \in H^{l+1}(\Omega)$  $\begin{array}{rl} \|v - \Pi^{0,k}_{c,h}(v)\|{0,\Omega} & \leq c h^{l+1} |v|{k+1,\Omega} \\ \|v - \Pi^{1,k}_{c,h}(v)\|{1,\Omega} & \leq c h^{l} |v|{k+1,\Omega} \end{array}$ (42) Pour calculer $\Pi^{0,k}_{c,h}$ et $\Pi^{1,k}_{c,h}$ on a besoin de résoudre les problèmes Problem 1: Projection $l^2$ Soit $v$ une fonction de $L^2$, calculer $\Pi^{0,k}_{c,h}(v) \in P^k_{c,h}$ tel que $\forall v_h \in P^k_{c,h}$ on a $(\Pi^{0,k}_{c,h}(v), v_h)_{0,\Omega} = (v, v_h)_{0,\Omega}$ Problem 2: Projection $H^1$ Soit $v$ une fonction de $H^1$, calculer $\Pi^{1,k}_{c,h}(v) \in P^k_{c,h}$ tel que $\forall v_h \in P^k_{c,h}$ on a $(\Pi^{1,k}_{c,h}(v), v_h)_{1,\Omega} = (v, v_h)_{1,\Omega}$ #### 3.6.2. Interpolant de Lagrange Definition 40. Notons $V_h$ un espace $H^1$ conforme, $\{\Psi_i\}_{1,...,N}$ une base nodale de $V_h$ et $\{x_1,...,x_N\}$ les noeuds associés alors l’interpolant de Lagrange est défini par $\Ilag{h}: v \in C^0(\Omega_h) \mapsto \sum_{i=1}^N v(x_{i}) \Psi_i \in V_h$  La restriction de l’interpolant de Lagrange à une cellule $K$ coincide avec l’interpolant de Lagrange appliqué à la fonction dans la cellule $K$: $\Ilag{h}(v)|_K = \Ilag{h}(v|_K)$ Theorem 16. Supposons $\{\mathcal{T}_h\}_{h>0}$ une famille de maillage quasi-uniforme et conformes, $\mathbb{P}_k \subset \hat{P}$ et $k+1 > \frac{d}{2}$. Alors il existe une constante $c$ telle que pour tout $h$ et $v \in H^{k+1}(\Omega_h)$ $\|v - \Ilag{h}(v)\|_{0,\Omega_h} + h |v - \Ilag{h}(v)|_{1,\Omega_h} \leq c h^{k+1} |v|_{k+1,\Omega_h}$ Nous allons vérifier sur un exemple ce théorème. Nous considérons pour cela $\alpha$ un réel et $\mathbf{x}=(x_1,...,x_d) \in \mathbb{R}^d$ un point de $\Omega = [0,1]^d, d=1,2,3$ et $v$ la fonction définie par  $\begin{array}{rl} v : \Omega &\rightarrow \mathbb{R}\\ \mathbf{x} &\rightarrow ( \mathbf{x} \cdot \mathbf{x} )^{\alpha/2}\ \Pi_{i=1}^d( 1-x_i^2) \end{array}$ (43) Nous construisons l’interpolant de Lagrange de $v$ dans $P^k_{c,h}$ avec $k=1,...,5$ et $d=1,2,3$ et étudions l’erreur d’interpolation $L^2$ et $H^1$ du théorème [eq en échelle log-log. Nous devons obtenir des droites de pentes $k$ (resp. $k+1$) pour la norme $L^2$ (resp. $H^1$.) ### 3.7. Approximation iso-paramétrique Quand le domaine est courbe, si nous désirons obtenir des propriétés de convergence optimale nous avons besoin de discretiser le bord du domaine avec suffisamment de précision Notons $(\hat{K},\hat{P}_{\mathrm{geo}}, \hat{\Sigma}_{\geo})$ l’élément fini géométrique et $(\hat{K},\hat{P}, \hat{\Sigma})$ l’élément fini de référence pour $V_h$. Definition 41. Si $\hat{P}_{\mathrm{geo}} = \hat{P}$ l’approximation est dite iso-parametrique. Definition 42. Si $\hat{P}_{\mathrm{geo}} \subset \hat{P}$ l’approximation est dite sub-parametrique. Definition 43. Si $\hat{P} \subset\hat{P}_{\mathrm{geo}}$ l’approximation est dite sur-parametrique.  Gmsh [1] un mailleur libre permet de générer des maillage d’ordre élevé jusqu’à l’ordre 5 en 2D et 4 en 3D Theorem 17. Supposons que $\{\mathcal{T}_h\}_{h>0}$ une famille de maillage quasi-uniformes et conformes, $\mathbb{P}_k \subset \hat{P}$ et $k+1 > \frac{d}{2}$ et $k_{\mathrm{geo}} = k$. Alors il existe une constante $c$ telle que pour tout $h$ et $v \in H^{k+1}(\Omega_h)$ $\|v - \Ilag{h}(v)\|_{0,\Omega_h} + h |v - \Ilag{h}(v)|_{1,\Omega_h} \leq c h^{k+1} |v|_{k+1,\Omega_h}$  Les résultats sont identiques à ceux du theorème Theorem 16. Nous allons vérifier sur un exemple à l’erreur d’approximation de l’élément géométrique. Considérons les cercles unité généré par une transformation affine, noté $\Omega^1_h$, et d’ordre 2, noté $\Omega^2_h$, et calculons leur aire respective. Construisons une famille de maillage $\{\calTh\}_{h >0}$, par exemple $h\in \{0.4, 0.2, 0.1, 0.05\}$ et calculons l’erreur entre le calcul exact de l’aire $\pi$ et le calcul numérique $\int_{\Omega^1_h} 1$ et $\int_{\Omega^2_h} 1$ respectivement. Le listing suivant présente le code C++ pour effectuer cela Calcul de l’aire d’un cercle ``````// genere le maillage pass:[$\Omega^1_h$] associe une transformation affine auto mesh1 = unitCircle<1>(); auto area1 = integrate( _range=elements(mesh1), _expr=cst(1.)).evaluate().norm(); std::cout << "|pi-surface<1>| = " << math::abs(pi-area1) << "\n"; // genere le maillage pass:[$\Omega^2_h$] associe une transformation d'ordre 2 auto mesh2 = unitCircle<2>(); auto area2 = integrate( _range=elements(mesh2), _expr=cst(1.)).evaluate().norm(); std::cout << "|pi-surface<2>| = " << math::abs(pi-area2) << "\n";`````` La table Convergence de l’approximation $P_1$ et $P_2$ présente les erreurs d’approximation et la figure [fig:circle] présente les courbes de convergence en échelle log-log ainsi que les pentes associées à ces courbes. On s’attend d’après le théorème [eq appliqué à des pentes à l’élément fini géométrique. Table 1. Convergence de l’approximation $P_1$ et $P_2$ h error $P_1$ error $P_2$ 0.4000 1.41593e-01 4.87932e-04 0.2000 2.62996e-02 1.65709e-05 0.1000 5.73876e-03 7.86831e-07 0.0500 1.34419e-03 4.31428e-08 0.0250 3.36079e-04 2.69659e-09 La table Convergence de l’approximation $P_1$ et $P_2$ présente les résultats de convergence. On observe un phénomène de super-convergence pour le cas $\Omega^2_h$, on obtient un ordre de convergence $4$ et nous devrions obtenir $3$. ### 3.8. Feel++ En Feel++, un maillage `Mesh` est décomposé en un ensemble d’éléments décomposés en sous entités (volume,face,arête,point): • faces(decomposés en sous entités), • arêtes(decomposés en sous entités) et • points. et à chaque élément $K$ est associé une transformation géométrique $T_K$. L’ordre polynomial de la transformation géométrique est donné par le second argument template Type d’un object maillage en Feel++ ``````Mesh<Simplex<d, k_geo>> (1) Mesh<Hypercube<d, k_geo>> (2)``````  1 Un maillage de simplexes en dimension $d$ avec transformation géométrique d’ordre $k_{\mathrm{geo}}$. 2 Un maillage d’hypercubes en dimension $d$ avec transformation géométrique d’ordre $k_{\mathrm{geo}}$. Ain de parcourir les éléments et faces du maillage, Feel++ fournit des fonctions renvoyant des itérateurs (début et fin) sur ces ensembles `elements(mesh)` retourne 2 itérateurs sur l’ensemble des éléments du maillage `markedelements(mesh,<int>)` et `markedelements(mesh,<string>)` retourne 2 itérateurs sur les éléments marqués par l’entier `<int>` et la chaîne des caractères `<string>` respectivement, cela correspondra typiquement à des propriétés de matériau `boundaryfaces(mesh)` retourne 2 itérateurs sur les faces au bord du maillage `markedfaces(mesh,<int>)` et `markedelements(mesh,<string>)` retourne 2 itérateurs sur les faces marquées par l’entier `<int>` et la chaîne des caractères `<string>` respectivement, ca correspondra typiquement aux conditions aux limites  link pending to Feel++ user manual mesh iterators section. L’espace d’approximation $V_h$ $H^1$ conforme (espaces de functions continues sur $\Omega$ polynomiales par morceaux de degré $\leq k$) est défini comme suit ``````FunctionSpace<Mesh<Simplex<d, k_geo> >, bases<Lagrange<k> > > V_h; (1) FunctionSpace<Mesh<Hypercube<d, k_geo> >, bases<Lagrange<k> > > V_h; (2)``````  1 $\Pch{k}$ 2 $\Qch{k}$ ### 3.9. Du problème global aux éléments locaux On va maintenant faire le lien entre la résolution d’un problème par méthode d’éléments finis et les notions qui viennent d’être introduites. Soit une EDP à résoudre sur un domaine $\Omega$, et $V$ l’espace de Hilbert dans lequel on cherche une solution de la formulation variationnelle du problème. On réalise un maillage de $\Omega$ par une famille affine de $N_e$ éléments finis $(K_i,\Sigma_i,P_i)_{i=1,\ldots,N_e}$. Par unisolvance, la solution approchée $u_h$ sera entièrement définie sur chaque élément $(K_i,\Sigma_i,P_i)$ par ses valeurs sur les points de $\Sigma_i$, qu’on appellera les noeuds du maillage. Il est à noter qu’un noeud sera en général commun à plusieurs éléments adjacents. Le nombre total de noeuds $N_h$ est donc inférieur à $N_e\times\hbox{Card} \Sigma_i$, et on a dim $V_h = N_h$. Notons $a_1,\ldots,a_{N_h}$ les noeuds du maillage. Le problème approché se ramène donc à la détermination des valeurs de $u_h$ aux points $a_i$: ce sont les degrés de liberté du problème approché. On va construire une base de $V_h$ en associant à chaque ddl $a_i$ un vecteur de la base. On définit ainsi les fonctions de base globales $\varphi_i$ ($i=1,\ldots,N_h$) par  ${\varphi_i}_{|K_j} \in P_j, \quad j=1,\ldots,N_e \mbox{ et } \varphi_i(a_j)=\delta_{ij}, 1\le i,j \le N_h$ (44) L’espace d’approximation interne est donc alors :  $V_h = \hbox{Vect }\left\{\varphi_1,\ldots,\varphi_{N_h}\right\}$ (45) Exemple de fonction de base globale (élément triangulaire $P_1$) Il est facile de remarquer qu’une telle fonction $\varphi_i$ est nulle partout, sauf sur les éléments dont $a_i$ est un noeud. En effet, si $a_i$ n’appartient pas à un élément $K$, $\varphi_i$ est nulle sur tous les noeuds de $K$, et donc sur $K$ tout entier par unisolvance. De plus, sur un élément $K$ dont $a_i$ est un noeud, $\varphi_i$ vaut 1 sur $a_i$ et 0 sur les autres noeuds de $K$. Donc $\varphi_i\, _{|K}$est une fonction de base locale de $K$.  la fonction de base globale $\varphi_i$ est donc construite comme la réunion des fonctions de base locales sur les éléments du maillage dont $a_i$ est un noeud. C’est à ce niveau que se situe le lien entre les définitions locales introduites au Elément fini de Lagrange et le problème global approché à résoudre. Par ailleurs, ceci implique que tous les calculs à effectuer sur les fonctions de base globales peuvent se ramener à des calculs sur les fonctions de base locales, et donc simplement à des calculs sur l’élément de référence (car on a maillé le domaine avec une famille d’éléments finis affine-équivalents).  Maillage non conforme Ce type de définition des fonctions de base n’est possible que si le maillage est conforme, c’est à dire si l’intersection entre deux éléments est soit vide, soit réduite à un sommet ou une arête en dimension 2 (ou à un sommet, une arête ou une face en dimension 3). On interdit ainsi les situations du type de celle de la figure. ### 3.10. Exercices 1. Calculer les fonctions de base locales des éléments finis de Lagrange introduits dans ce chapitre. 2. Donner l’allure des fonctions de base globales correspondantes. Sont-elles continues ? dérivables ? 3. Pour les éléments finis de Lagrange introduits dans ce chapitre, écrire le changement de variable affine entre élément quelconque et élément de référence. # Approximation de problèmes coercifs Dans ce chapitre, on s’intéresse à l’approximation de problèmes coercifs: • le Laplacien et des variantes sur les conditions aux limites et le l’équation elle-même; • l’élasticité linéaire qui permet de modéliser de petites déformations d’un milieu continu déformable. ## 4. Le Laplacien On s’instéresse dans cette section à l’approximation élément fini conforme du problème suivant: Problem 3: Formulation forte du Laplacien On cherche $u$ telle que $\begin{split} -\Delta u &= f \mbox{ dans } \Omega\\ u &= 0 \mbox{ sur } \partial \Omega \end{split}$ ### 4.1. Cadre Mathematique On suppose que $f \in L^2(\Omega)$. La formulation faible du problème Problem 3 est la suivante: Problem 4: Formulation faible pour des conditions de Dirichlet homogènes On cherche $u \in H^1_0(\Omega)$ telle que $\label{eq:65} \int_\Omega \nabla u \cdot \nabla v = \int_\Omega f\, v,\quad \forall v \in H^1_0(\Omega)$ ### 4.2. Problème bien posé Introduisons • $V = H^1_0(\Omega)$ doté de la norme $\|\cdot\|_{1,\Omega}$ telle que $\|v\|_{1,\Omega} = (\|v\|^2_{0,\Omega} + \|\nabla v\|^2_{0,\Omega})^{1/2}$ • la forme bilinéaire $a \in \mathcal{L}(V \times V, \R)$ telle que $a(u,v) = \int_\Omega \nabla u \cdot \nabla v$ • la forme linéaire $\ell \in \mathcal{L}(V, \R)$ telle que $l(v) = \int_\Omega f \nabla v$ Le problème Problem 4 s’écrit sous forme abstraite Problem 5. On cherche $u \in V$ telle que $a(u, v) = \ell(v), \quad \forall v \in V$ L’espace $V$ est un espace de Hilbert et les formes $a$ et $\ell$ sont continues sur $V\times V$ et $V$ respectivement. Il ne reste plus qu’à vérifier si le problème est bien posé (existence d’une solution unique). Pour cela on utilise démontre la coercivité de la forme bilinéaire $a$ sur $V \equiv H^1_0(\Omega)$. Ceci se fait grâce au lemme suivant: Lemma 1. Soit $\Omega$ un ouvert borné de $\R$. Il existe une constante $c_\Omega$ (dépendente de $\Omega$ telle que $\forall v \in H^1_0(\Omega),\quad \|v\|_{0,\Omega} \le c_\Omega \|\nabla v\|_{0,\Omega}$  $c_\Omega$ est homogène à une longeur et peut être interprétée comme une mesure caractéristique de $\Omega$. Grâce à l’inégalité de Poincaré, on a le résultat suivant Proposition 7. La forme bilinéaire $a$ du problème Problem 4 est coercive Proof 4. On note tout d’abord que par l’inégalité de Poincaré et la définition de $\|\cdot\|_{1,\Omega}$ $\|v\|^2_{1,\Omega} \le (1 + c^2_\Omega) \|\nabla v\|^2_{0,\Omega}$ On en déduit que $\forall v \in H^1_0(\Omega),\quad a(v,v) = \|\nabla v\|^2_{0,\Omega} \ge \frac{1}{1+c^2_\Omega} \|v\|^2_{1,\Omega}$ Le Lemme de Lax-Milgram Theorem 12 permet alors de conclure sur l’existence d’une solution unique pour le problème Problem 4. ### 4.3. Approximation conforme On utilise une approximation conforme par éléments finis de Lagrange. On considère $\Omega$ un polygone ou polyhèdre régulier de $\RR^2$ ou $\RR^3$ respectivement et un maillage $\calTh = \{K_e\}_{e=1...\Ne}$ de $\Omega$. On considère un élément fini de référence $(\hat{K},\hat{P},\hat{\Sigma})$ tel que $\Pk \subset \hat{P}$ et $k+1 > \frac{d}{2}$, voir le théorème Theorem 17. On note $L^k_{c,h} = \{ v_h \in C^0(\bar{\Omega}); \forall K \in \mathcal{T}_h,\ v_h \circ T_K \in \hat{P}\}$ où $T_K$ est la tranformation géométrique de $\hat{K}$ dans $K$. • Si on utilise $\hat{P}=\Pk{k}$ on a $L^k_{c,h} = P^k_{c,h}$. • Si on utilise $\hat{P}=\Qk{k}$ on a $L^k_{c,h} = Q^k_{c,h}$. Afin de construire un espace d’approximation conforme $V_h \subset V =H^1_0(\Omega)$ on prend $\label{eq:71} V_h = L^k_{c,h} \cap H^1_0(\Omega)$ c’est à dire que les fonctions de $V_h$ satisfont les conditions aux limites outre le fait d’être dans $L^k_{c,h}$. Le problème discret s’écrit alors Problem 6. Trouver $u_h \in V_h$ telle que $a(u_h,v_h) = \ell(v_h),\, \forall v_h \in V_h$ qui est bien posé (existence et unicité de $u_h$) car $a$ est coercive sur $V$ et que $V_h \subset V$. On a le résultat suivant: Theorem 18. On suppose que $u$ solution de Problem 4 est dans $H^{k+1}(\Omega) \cap H^1_0(\Omega)$, alors il existe une constante $c_1$ telle que pour tout $h$ $\|u-u_v\|_{1,\Omega} \le c_1 h^k |u|_{k+1,\Omega}$ et il existe une constante $c_2$ telle que pour tout $h$ $\|u-u_v\|_{0,\Omega} \le c_2 h^{k+1} |u|_{k+1,\Omega}$ Proof 5. La preuve de ([eq) est obtenue grâce au lemme de Cea (8) et du théorème d’interpolation ([eq). La preuve de ([eq) est obtenue grâce au lemme d’Aubin-Nitsche qui permet d’affirmer qu’il existe une constante $c$ telle que $\|u-u_h\|_{0,\Omega} \leq c h |u-u_h|_{1,\Omega}$ et donc que ([eq) se déduit de ([eq). ### 4.4. Implémentation avec Feel++ Avec Feel++, on ne construit pas explicitement l’espace $V_h$ mais $L^k_{c,h}$. Le traitement des conditions aux limites de Dirichlet du problème ([eq) peut être effectué de diverses façons, nous en verrons une. Commencons par le maillage, dans un premier temps nous définissons le type du maillage contenant soit des éléments de type simplexe (segment,triangle, tetrahèdre) ou de type hypercube (segment, quadrangle, hexahèdre). `````` // un maillage de simplexe dans pass:[$\R$] telle que la transformation // géométrique pass:[$T_K,\ K \in \calTh$] soit affine typedef Mesh<Simplex<d,1> > mesh_type; // un maillage d'hypercube dans pass:[$\R$] telle que la transformation // géométrique pass:[$T_K,\ K \in \calTh$] soit affine en chacune des variables // typedef Mesh<Hypercube<d,1> > mesh_type; // generate the mesh associated to the unit square pass:[$[0,1]^2$] using triangles auto mesh = unitSquare();``````  Le mot clé `auto` permet de faire de l’inférence de type, pour plus de détails consultez la page C++11 de Wikipedia. Ensuite nous pouvons définir l’espace $L^k_{c,h}$, `````` // Vh est une structure de donnée allouée dynamiquement auto Vh = Pch<1>( mesh ); // u est un élément de Vh auto u = Vh->element(); // u est un autre élément de Vh auto u = Vh->element();`````` À présent, nous définissons les formes bilinéaires $a$ et $\ell$ qui sont respectivement des formes bilinéaires et linéaires. `````` auto a = form2( _test=Vh, _trial=Vh ); (1) a = integrate( _range=elements(mesh), _expr=gradt(u)*trans(grad(v)) ); (2) auto l = form1( _test=Vh ); (3) l = integrate( _range=elements(mesh), _expr=f*id(v) ); (4)``````  1 $a \in \mathcal{L}(V_h \times V_h,\ \RR)$ 2 $a = \sum_{e=1...\Ne} \int_{K_e} \nabla \varphi_j \cdot \nabla \varphi_i,\quad i,j=1...,\dim{V_h}$ 3 $\ell \in \mathcal{L}(V_h,\ \RR)$ 4 $\ell = \sum_{e=1...\Ne} \int_{K_e} f \varphi_i,\quad i=1...,\dim{V_h}$ Afin de traiter les conditions aux limites de Dirichlet homogènes, on peut utiliser le mot-clé `on` qui permet de les imposer de manière forte. `` a += on(_range=boundaryfaces(mesh), _element=u, _rhs=l, _expr=constant(0.) );``  Le mot-clé `constant` permet de transformer une type numérique ( entier, flottant) en une expression utilisable par le langage de Feel++. Notez également l’opération `+=` qui permet de rajouter le traitement des conditions de Dirichlet tout en gardant les contributions précédentes. L’opération `=` aurait d’abord remis à $0$ les entrées de la matrice associée à $a$. Enfin nous pouvons résoudre le problème Problem 6 `` a.solve( _rhs=l, _solution=u );`` Le listing complet ### 4.5. Conditions aux limites #### 4.5.1. Conditions aux limites de Dirichlet non homogène On suppose toujours $f \in L^2(\Omega)$ et on se donne une fonction $g \in C^{0,1}(\partial \Omega)$ [$g$ est Lipschitzienne sur $\partial \Omega$]. On s’intéresse au problème suivant: Problem 7. On cherche $u : \Omega \rightarrow \RR$ telle que $\begin{split} -\Delta u &= f \mbox{ dans } \Omega\\ u &= g \mbox{ sur } \partial \Omega \end{split}$  L’hypothèse $g \in C^{0,1}(\partial \Omega)$ permet d’affirmer qu’il existe $u_g \in H^1(\Omega)$ telle que $u_{g_{|\partial \Omega}} = g$. On se ramène au problème avec conditions de Dirichlet homogène en faisant le change d’inconnue $u_0=u-u_g$ et on s’intéresse au problème suivant: Problem 8. On cherche $u_0 \in H^1_0(\Omega)$ telle que $a(u_0,v) = \ell(v) - a(u_g,v),\quad \forall v \in H^1_0(\Omega)$ Ce problème est bien posé d’après Lax-Milgram, voir section précédente. Theorem 19. On suppose que $u$ solution de [prob:8] est dans $H^{k+1}(\Omega) \cap H^1_0(\Omega)$, alors il existe une constante $c_1$ telle que pour tout $h$ $\|u-u_v\|_{1,\Omega} \le c_1 h^k |u|_{k+1,\Omega}$ et il existe une constante $c_2$ telle que pour tout $h$ $\|u-u_v\|_{0,\Omega} \le c_2 h^{k+1} |u|_{k+1,\Omega}$ Avec Feel++, les conditions Dirichlet non-homogènes sont traitées par exemple avec le mot-clé `on`. Conditions de Dirichlet non homogènes `````` auto g = sin(2*pi*Px() ); (1) (2) ... a += on( _range=boundaryfaces(mesh), _expr=g ); (3)``````  1 définition de la fonction, p.ex $g=sin(2 \pi x)$ 2 définition de $a$ 3 ajout des conditions de Dirichlet non-homogènes  Il n’y a pas besoin de rajouter le terme $a(u_g,v)$ au second membre $\ell(v)$, il est pris en compte automatiquement par `on`. Voici le listing complet de l’exemple du laplacien avec conditions de Dirichlet non-homogène ### 4.6. Condition aux limites de Neumann Étant donnés un réel $\mu$ strictement positif, $f \in L^2(\Omega)$ et $g \in L^2(\partial \Omega)$, on s’intéresse au problème suivant: Problem 9. On cherche $u : \Omega \rightarrow \RR$ telle que $\begin{split} -\Delta u + \mu u &= f, \mbox{ dans } \Omega\\ \partial_\Next u &= g, \mbox{ sur } \partial\Omega \end{split}$ où $\partial_\Next u = \nabla u \cdot \Next = \sum_{i=1}^d n_i \partial_i u$ dénote la dérivée normale de $u$ avec $\Next=(n_1,...,n_d) \in \RR^d$ la normale extérieure unitaire en un point du bord de $\Omega$. La formulation faible s’écrit Problem 10. On cherche $u \in H^1(\Omega)$ telle que $a( u, v ) = \ell(v),\ \forall v \in H^1(\Omega)$ avec $a( u, v ) = \int_\Omega \nabla u \cdot \nabla v + \mu u v$ et $\ell( v ) = \int_\Omega f v + \int_{\partial\Omega} g v$ On a $a(v, v) = \int_\Omega \nabla v \cdot \nabla v + \mu v v \ge \min(1,\mu) \int_\Omega \nabla v \cdot \nabla v + v v = \min(1,\mu) \|v\|_{1,\Omega}, \quad \forall v \in H^1(\Omega)$ ce qui nous permet d’affirmer que $a$ est coercive sur $H^1(\Omega)$ et que le problème Problem 10 est bien posé grâce à Lax-Milgram. On peut approcher le problème Problem 10 par des éléments finis de Lagrange. On utilise la formulation développée dans la section Approximation conforme Problem 11. On cherche $u_h \in P^k_{c,h}$ tel que $a_h(u_h,w_h)=\ell(w_h), \quad \forall w_h \in P^k_{c,h}$  Par rapport aux conditions de Dirichet, les conditions de Neumann sont directement (naturellement) traitées par la formulation faible. Elles ne sont pas directement imposées dans l’espace et les fonctions tests peuvent prendre des valeurs non-nulles au bord. Ces conditions sont traitées de manière approchée et non pas exacte. L’analyse du problème Problem 11 est identique par le théorème Theorem 18 aux mêmes estimations que dans le cas Dirichlet homogène. #### 4.6.1. Cas $\lambda=0$ Le cas $\lambda=0$ présente quelques difficultés techniques. On a Problem 12. On cherche $u : \Omega \rightarrow \RR$ telle que $\begin{split} -\Delta u &= f, \mbox{ dans } \Omega\\ \partial_\Next u &= g, \mbox{ sur } \partial\Omega \end{split}$  On observe ici qu’une condition nécessaire d’existence de solution est que $\int_\Omega f + \int_{\partial \Omega} g = 0$ de par le théorème de la divergence.  la deuxième observation est que la solution de Problem 12 est connue à une constante additive près: si $u$ est solution et $c\in \RR$ alors $u+c$ est également solution. Il convient alors de chercher la solution dans l’espace fonctionel suivant  $H^1_*(\Omega) = \{ v \in H^1(\Omega):\quad \int_\Omega v = 0 \}$ La formulation faible du problème Problem 12 s’écrit alors Problem 13. On cherche $u\in H^1_*(\Omega)$, telle que $a(u,v) = \ell(v) \quad \forall v\in H^1_*(\Omega)$ avec $a(u,v)=\int_\Omega \nabla u \cdot \nabla v$. Le caractère bien posé du problème Problem 13 résulte de la coercivité de $a$ sur $H^1_*(\Omega)$ qui elle-même résulte du Lemme suivante Lemma 2: Lemme d’inégalité de Poincaré-Wirtinger Soit $\Omega$ un ouvert borné de $\RR^d$, il existe une constante $C_\Omega$ telle que $\|v\| \leq C_\Omega \|\nabla v\|_{0,\Omega},\quad \forall v \in H^1_*(\Omega)$ Le problème Problem 13 peut être approché par des éléments finis de Lagrange ce qui conduit au problème discret suivant Problem 14. On cherche $u_h \in P^k_{c,h}$ telle que a_h(u_h,v_h)=\ell(v_h) \quad \forall v_h\in P^k_{c,h}  L’espace d’approximation $P^k{c,h}$ n’est pas conforme dans $H^1_*(\Omega)$  le problème Problem 14 est singulier, il a une infinité de solution. L’une d’entre elle peut être approchée par une méthode itérative de type gradient conjugué. Notons $u^*_h$ cette solution alors la solution à moyenne nulle est $u_h=u^*_h-\frac{1}{|\Omega|}\int_\Omega u^*_h$ Il s’agit donc d’effectuer a posteriori du calcul un post-processing pour se ramener à la solution à moyenne nulle. Dans Feel++ on utilisera la fonction `mean` ``````// ... // u est la solution du problème discret a.solve(_rhs=l,_solution=u); // calcul de la valeur moyenne de u double meanu = mean(_range=elements(mesh),_expr=idv(u))(0,0); // calcul de la solution à moyenne nulle u = u - meanu u.add(-meanu);`````` ### 4.7. Conditions aux limites de Robin Étant donnés un réel $\mu$ strictement positif, $f \in L^2(\Omega)$ et $g \in L^2(\partial \Omega)$, on s’intéresse au problème suivant: Problem 15. On cherche $u : \Omega \rightarrow \RR$ telle que $\begin{split} -\Delta u &= f, \mbox{ dans } \Omega\\ \mu u + \partial_\Next u &= g, \mbox{ sur } \partial\Omega. \end{split}$ La formulation faible s’écrit Problem 16. On cherche $u \in H^1(\Omega)$ telle que $\label{eq:84} a( u, v ) = \ell(v),\ \forall v \in H^1(\Omega)$ avec $a( u, v ) = \int_\Omega \nabla u \cdot \nabla v + \int_{\partial \Omega} \mu u v$ et $\ell( v ) = \int_\Omega f v + \int_{\partial\Omega} g v$ On a $\begin{split} a(v, v) & = \int_\Omega \nabla v \cdot \nabla v + \int_{\partial\Omega} \mu v v \\ & \geq \min(1,\mu)\left( \int_\Omega \nabla v \cdot \nabla v + \int_{\partial\Omega} v v\right) \\ &\geq \min(1,\mu) \|v\|_{1,\Omega}, \quad \forall v \in H^1(\Omega) \end{split}$ La forme bilinéaire $a$ est donc coercive et le problème Problem 16 est bien posé grâce à Lax-Milgram. On considère le problème discret suivant Problem 17. On cherche $u_h \in P^k_{c,h}(\Omega)$ telle que $a_h(u_h,v_h) = \ell(v_h)\quad \forall v_h \in P^k_{c,h}$ Le problème est bien posé ($P^k_{c,h} \subset H^1(\Omega)$).  Comme pour le problème avec conditions de Neumann, les fonctions tests peuvent prendre des valeurs non nulles au bord. Les conditions de Robin(ou Fourier) ne sont satisfaites que de manière approchée et non pas exacte. La convergence de $u_h$ est donnée par le théorème Theorem 18. Considérons $\Omega=[0,1$^2] et les données $\mu=0.01$, $f=1$ et $g=0$. ## 5. Advection-diffusion-réaction avec diffusion dominante On s’intéresse au problème suivant: Problem 18. On cherche $u : \Omega \rightarrow \RR$ telle que $\begin{split} -\nabla \cdot ( \mathbf{\alpha} \nabla u ) + \mathbf{\beta} \cdot \nabla u + \mu u &= f \mbox{ dans } \Omega\\ u &= 0 \mbox{ sur } \partial \Omega\\ \end{split}$ La variation sur les conditions aux limites de la section [sec:cond-aux-limit] s’appliquent. • $-\nabla \cdot ( \mathbf{\alpha} \nabla u )$ est un terme de diffusion, • $\mathbf{\beta} \cdot \nabla u$ est un terme de convection, • $\mu u$ est un terme de réaction Ce type d’équation est très fréquente en ingéniérie, biologie ou encore finance. On suppose que $\mathbf{\alpha} \in [L^{\infty}(\Omega)]^{d\times d}$, $\mathbf{\beta} \in [L^{\infty}(\Omega)]^d$ et $\mu \in L^\infty(\Omega)$. On suppose que l’opérateur $\mathcal{L}$ tel que $\mathcal{L} u = -\nabla \cdot ( \alpha \nabla u ) + \beta \cdot \nabla u + \mu u$ est elliptique au sens suivant: Definition 44. L’opérateur $\mathcal{L}$ est elliptique si il existe une constante $\alpha_0$ telle que presque pour tout $x\in\Omega$ $\forall \xi=(\xi_1, \ldots , \xi_n)\in\RR^n,\quad { \sum_{i,j=1}^n \alpha_{ij}(x) \, \xi_i \, \xi_j \ge \alpha_0 \, \| \xi \|^2 }$  Le Laplacien est dans la catégorie des opérateurs elliptiques, il correspond à $\mathbf{\beta} = 0$, $\mu = 0$ et $\mathbf{\alpha}=\mathcal{I}_d$ avec $\mathcal{I}_d$ la matrice identité de $\RR^{d\times d}$. La formulation faible s’écrit sous la forme suivante: Problem 19. On cherche $u \in H^1_0(\Omega)$ telle que $a(u,v) = \ell(v) \quad \forall v \in H^1_0(\Omega)$ avec $a(u,v)=\int_\Omega (\mathbf{\alpha} \nabla u ) \cdot \nabla v + (\mathbf{\beta} \cdot \nabla u) v + \mu u v$ et $\ell(v) = \int_\Omega f v$ On note $\gamma = \essinf_{x \in \Omega} (\mu -\frac{1}{2} \nabla \cdot \mathbf{\beta})$, on peut alors montrer que sous la condition $\alpha_0 > - \min( 0, \gamma c_\Omega )$ où $c_\Omega$ est la constante de l’inégalité de Poincaré. La forme bilinéaire $a$ est coercive et donc que, grâce au théorème de Lax-Milgram, le problème Problem 19 est bien posé. La coercivité est garantie si $\alpha_0$ est suffisamment grand c’est à dire que si la diffusion est dominante. L’approximation élément fini est similaire à celle du Laplacian, de plus les variantes sur les conditions aux limites s’appliquent également: condition de Dirichlet non homogène, de Neumann ou de Robin. ## 6. Élasticité Linéaire On s’intéresse dans cette section à l’approximation par éléments finis de problèmes de mécanique des milieux continus en 3D. Soit $\mathbf{f}: \Omega \rightarrow \mathbb{R}^3$ la charge extérieur s’appliquant au domaine $\Omega$. On note $\disp: \Omega \rightarrow \mathbb{R}^3$ le déplacement de la structure induit par cette charge $\mathbf{f}$. En supposant que les déformations soient petites pour être modélisées dans le cadre de l’elasticité linéaire, on a la relation suivante à l’équilibre: $\nabla \cdot \stresst + \mathbf{f} = \mathbf{0} \mbox{ dans } \Omega$ où $\stresst: \Omega \rightarrow \RR^{d\times d}$ est le tenseur des contraintes défini par la relation $\stresst = \lambda \tr(\deformt) \Id + 2 \mu \deformt$ et $\deformt : \Omega \rightarrow \RR^{d\times d}$ est le tenseur des déformations défini par $\deformt = \frac{1}{2} \left( \nabla \disp + \nabla \disp^T \right),$ $\lambda$ et $\mu$ les coefficients de Lamé et $\Id$ la matrice identité de $\RR^{d\times d}$. On a alors $\stresst = \lambda( \nabla \cdot \disp ) \Id + \mu( \nabla \disp + \nabla \disp^T)$ ### 6.1. Coefficients de Lamé Ils sont des coefficients phénoménologiques contraints par les relations suivants: • $\mu >0$ • $\lambda + \frac{2}{3} \mu \ge 0$ Dans ce qui suit, on supposera que $\lambda \ge 0$ et ces coefficients constants. D’un point de vue pratique, ces coefficients sont obtenus par les module d’Young $E$ et coefficient de Poisson $\nu$ tels que $\lambda = \frac{E \nu}{( 1+\nu )*( 1-2 \nu )} , \quad \mu =\frac{E}{2 ( 1+\nu )}$ ## 7. Approximation de problèmes mixtes ### 7.1. Model Problems We consider now model problems as systems of PDEs where several functions are unknowns and which don’t play the same roles mathematically and physically. Stokes $\left\{\begin{array}[c]{rl} -\Delta u + \nabla p & = f\ \mbox{ in } \Omega\\ \nabla \cdot u & = 0\ \mbox{ in } \Omega \end{array}\right.$ where $u: \Omega \mapsto \RR^d$ is a velocity and $p: \Omega \mapsto \RR$ is a pressure. Darcy $\left\{\begin{array}[c]{rl} \sigma + \nabla u & = f\ \mbox{ in } \Omega\\ \nabla \cdot \sigma & = g\ \mbox{ in } \Omega \end{array}\right.$ where $\sigma: \Omega \mapsto \RR^d$ is a velocity and $u: \Omega \mapsto \RR$ is a hydraulic charge(pressure). ### 7.2. Applications We shall focus on Stokes, but the abstract setting of the next section is the same for Stokes and Darcy. Stokes and incompressible Navier-Stokes for Newtonian fluids The Stokes model is the basis for fluid mechanics models and is a simplication of the Navier-Stokes equations where the viscous effects/terms are much bigger than the convective ones $\left\{\begin{array}[c]{rl} \rho( \frac{\partial u}{\partial t} + u \cdot \nabla u) - \nu \Delta u + \nabla p & = f\ \mbox{ in } \Omega\\ \nabla \cdot u & = 0\ \mbox{ in } \Omega \end{array}\right.$ The first equation results from the conservation of momentum and the second from the conservation of mass. The well-posedness of these problems results from a so-called which is not automatically transfered at the discrete level. In practice In order to ensure that the finite element approximation is well-posed, we will need to choose approximation spaces that satisfy a compatibility condition that ensures that a discrete inf-sup condition is satisfied. ### 7.3. Saddle point problems #### 7.3.1. Abstract Continuous Setting Denote • $X$ and $M$ two Hilbert spaces.[2] • two linear forms $f \in X'=\mathcal{L}(X, \RR)$ and $g \in M'=\mathcal{L}(M, \RR)$ • $a \in \mathcal{L}(X\times X, \RR)$ and $b \in \mathcal{L}(X\times M, \RR)$ two bilinear forms We are interested in the following abstract problem: Problem 20. Look for $(u,p) \in X \times M$ such that $\left\{ \begin{array}[c]{rl} a(u,v) + b(v,p) & = f(v), \quad \forall v \in X\\ b(u,q) & = g(q), \quad \forall q \in M \end{array} \right.$ ##### Definition of a saddle point problem Definiton 1. If the bilinear form $a$ is symmetric and positive on $X\times X$, we say that problem [prob:chmixte:1] is a saddle point problem. The structure of the problem is as follows • the space of solution is the same of the test space • the unknown $p$ does not appear in the second equation • the unknown functions $u$ and $p$ are coupled via the same bilinear form $b$ is the first and second equation. The next question is : #### 7.3.2. Well posed problem ##### Reformulation Let’s rewrite Problem Problem 20. Denote $V=X\times M$ and introduce $c \in \mathcal{L}(V\times V, \RR)$ such that $c((u,p),(v,q)) = a(u,v)+b(v,p)+b(u,q)$ and $h\in \mathcal{L}(V,\RR)$ such that $h(v,q) = f(v)+g(q)$ then problem Problem 20 reads Problem 21. Look for $(u,p) \in V$ such that $\begin{array}[c]{rl} c((u,p), (v,q)) & = h(v,q), \quad \forall (v,q) \in V \end{array}$ [env.theorem#thr:chmixte:1] We suppose that $a$ is coercive over $X$, the problem [prob:chmixte:2] is well-posed if and only if the bilinear form $b$ satisfies the following inf-sup condition: there exists $\beta > 0$ such that $\inf_{q \in M} \sup_{v \in X} \frac{b(v,q)}{||v||_X ||q||_M} \geq \beta$  Lax-Milgram provides only a sufficient condition for well-posedness  The form $c$ in Problem [prob:chmixte:2] does not satisfy Lax-Milgram. Let’s introduce the so-called Lagrangian $l \in \mathcal{L}(X\times M, \RR)$ defined by $l(v,q) = \frac{1}{2} a(v,v) + b(v,q) - f(v) - g(q)$ We say that the point $(u,p)\in X\times M$ is a saddle point of $l$ if $\forall (v,q) \in X\times M, \quad l(u,q) \leq l(u,p) \leq l(v,p)$ Proposition 8. Under the hypothesys of theorem [thr:chmixte:1], the Lagrangian $l$ defined by has a unique saddle point. Moreover this saddle point is the unique solution of problem Problem 20. ### 7.4. Finite element approximation #### 7.4.1. Abstract Discrete Problem We now turn to the approximation of the problem Problem 20 by a standard Galerkin method in a conforming way. Denote the two spaces $X_h \subset X$ and $M_h \subset M$, we consider the following problem: Problem 22. Look for $(u_h,p_h) \in X_h \times M_h$ such that $\left\{ \begin{array}[c]{rl} a(u_h,v_h) + b(v_h,p_h) & = f(v_h), \quad \forall v_h \in X_h\\ b(u_h,q_h) & = g(q_h), \quad \forall q_h \in M_h \end{array} \right.$ [env.theorem#thr:chmixte:2] We suppose that $a$ is coercive over $X$ and that $X_h \subset X$ and $M_h \subset M$. Then the problem [prob:chmixte:3] is well-posed if and only if the following discrete inf-sup condition is satisfied: there exists $\beta_h > 0$ such that $\inf_{q_h \in M_h} \sup_{v_h \in X_h} \frac{b(v_h,q_h)}{||v_h||_{X_h} ||q_h||_{M_h}} \geq \beta_h$ The compatibility condition problem Problem 22, to be well posed, requires that the spaces $X_h$ and $M_h$ satisfy the condition. This is known as the Babuska-Brezzi (BB) or Ladyhenskaya-Babuska-Brezzi (LBB). Regarding error analysis, we have the following lemma Thanks to the Lemma of Céa applied to Saddle-Point Problems, the unique solution $(u,p)$ of problem Problem 22 satisfies $\begin{array}[c]{rl} ||u-u_h||_X & \leq c_{1h} \inf_{v_h \in X_h} ||u-v_h||_X + c_{2} \inf_{q_h \in M_h} ||q-q_h||_M\\ ||p-p_h||_X & \leq c_{3h} \inf_{v_h \in X_h} ||u-v_h||_X + c_{4h} \inf_{q_h \in M_h} ||q-q_h||_M \end{array}$ where • $c_{1h} = (1+\frac{||a||_{X,X}}{\alpha})(1+\frac{||b||_{X,M}}{\beta_h})$ with $\alpha$ the coercivity constant of $a$ over X. • $c_{2} = \frac{||b||_{X,M}}{\alpha}$ • $c_{3h} = c_{1h} \frac{||a||_{X,X}}{\beta_h}$, $c_{4h} = 1+ \frac{||b||_{X,M}}{\beta_h}+\frac{||a||_{X,X}}{\beta_h}$  The constants $c_{1h}, c_{3h}, c_{4h}$ are as large as $\beta_h$ is small. #### 7.4.2. Linear system associated The discretisation process leads to a linear system. We denote * $N_u = \dim {X_h}$ * $N_p = \dim {M_h}$ * $\{\phi_i\}_{i=1,...,N_u}$ a basis of $X_h$ * $\{\psi_k\}_{k=1,...,N_p}$ a basis of $M_h$ * for all $u_h = \sum_{i=1}^{N_u} u_i \phi_i$, we associate $U \in \R{N_u}$, $U=(u_1,\ldots,u_{N_u})^T$, the component vector of $u_h$ is $\{\phi_i\}_{i=1,\ldots,N_u}$ * for all $p_h = \sum_{k=1}^{N_p} u_k \psi_k$, we associate $P \in \R{N_p}$, $P=(p_1,\ldots,p_{N_p})^T$, the component vector of $p_h$ is $\{\psi_k\}_{k=1,\ldots,N_p}$ The matricial form of problem Problem 22 reads $\begin{bmatrix} \mathcal{A} & \mathcal{B}^T\\ \mathcal{B} & 0 \end{bmatrix} \begin{bmatrix} U \\ P \end{bmatrix} = \begin{bmatrix} F\\ G \end{bmatrix}$ where the matrix $\mathcal{A} \in \R{N_u,N_u}$ and $\mathcal{B} \in \R{N_p,N_u}$ have the coefficients $\mathcal{A}_{ij} = a(\phi_j,\phi_i), \quad \mathcal{B}_{ki} = b(\phi_i,\psi_k)$ and the vectors $\mathcal{F} \in \R{N_u}$ and $\mathcal{G} \in \R{N_p}$ have the coefficients * $F_i=f(\phi_i)$ * $G_k=g(\psi_k)$ NOTE: 1. Since $a$ is symmetric and coercive, $\mathcal{A}$ is symmetric positive definite 2. The matrix of the system is symmetric but not positive 3. The inf-sup condition is equivalent to the fact that $\mathcal{B}$ is of maximum rank, i.e. $\ker(\mathcal{B}^T) = \{0 \}$. 4. From theorem [thr:chmixte:2], the matrix of the system is invertible When the Inf-sup condition is not satisfied, the counter examples when the inf-sup condition is not satisfied(e.g. $\mathcal{B}$ is not maximum rank ) occur usually in two cases: Locking $\dim {M_h} > \dim {X_h}$: the space of pressure is too large for the matrix $\mathcal{B}$ to be maximum rank. In that case $\mathcal{B}$ is injective ($\ker(\mathcal{B}) = \{0\})$. We call this locking. Spurious modes there exists a vector $Q^* \neq 0$ in $\ker(\mathcal{B}^T)$. The discrete field$q^*_h$ in $M_h$, $q^*_h=\sum_{k=1}^{N_p} Q^*_k \psi_k$, associated is called a spurious mode. $q^*_H$ is such that $b(v_h,q^*_h)=0.$ We now introduce the Uzawa matrix as follows Definition 45. The matrix $\mathcal{U} = \mathcal{B} \mathcal{A}^{-1} \mathcal{B}^T$ is called the Uzawa matrix. It is symmetric positive definite from the properties of $\mathcal{A}$, $\mathcal{B}$ Applications The Uzawa matrix occurs when eliminating the velocity in system and get a linear system on $P$: $\mathcal{U} P = \mathcal{B} \mathcal{A}^{-1} F - G$ then one application is to solve by solving iteratively and compute the velocity afterwards. ### 7.5. Mixed finite element for Stokes #### 7.5.1. Variational formulation We start with the Well-posedness at the continuous level • We consider the model problem with homogeneous Dirichlet condition on velocity $u = 0$ on $\partial \Omega$ • We suppose the $f \in [L^2(\Omega)$^d] and $g \in L^2(\Omega)$ with $\int_\Omega g = 0$ Introduce $L^2_0(\Omega) = \Big\{ q \in L^2(\Omega): \int_\Omega q = 0 \Big\}$ The condition comes from the divergence theorem applied to the divergence equation and the fact that $u=0$ on the boundary $\int_\Omega g = \int_\Omega \nabla \cdot u = \int_{\partial \Omega} u \cdot n = 0$ This is a necessary condition for the existence of a solution $(u,p)$ for the Stokes equations with these boundary conditions. We turn now to the variational formulation. The Stokes problem reads Problem 23. Look for $(u,p) \in [H^1_0(\Omega)$^d \times L^2_0(\Omega)] such that $\left\{ \begin{array}[c]{rl} \int_\Omega \nabla u : \nabla v -\int_\Omega p \nabla \cdot v & = \int_\Omega f \cdot v, \quad \forall v \in [H^1_0(\Omega)]^d\\ \int_\Omega q \nabla \cdot u & = - \int_\Omega g q, \quad \forall q \in L^2_0(\Omega) \end{array} \right.$ We retrieve the problem Problem 20 with $X=[H^1_0(\Omega)$^d] and $M=L^2_0(\Omega)$ and $\begin{array}[c]{rlrl} a(u,v) &= \int_\Omega \nabla u : \nabla v,& \quad b(v,p) &= -\int_\Omega p \nabla \cdot v,\\ \quad f(v) &= \int_\Omega f \cdot v,& \quad g(q) &= - \int_\Omega g q \end{array}$  Pressure up to a constant The pressure is known up to a constant, that’s why we look for them in $L^2_0(\Omega)$ to ensure uniqueness. #### 7.5.2. Finite element approximation Denote $X_h \subset [H^1_0(\Omega)$^d] and $M_h \subset L^2_0(\Omega)$ Problem 24. Look for $(u_h,p_h) \in X_h \times M_h$ such that $\left\{ \begin{array}[c]{rl} \int_\Omega \nabla u_h : \nabla v_h + \int_\Omega p_h \nabla \cdot v_h & = \int_\Omega f \cdot v_h, \quad \forall v_h \in X_h\\ \int_\Omega q_h \nabla \cdot u_h & = -\int_\Omega g q_h, \quad \forall q_h \in M_h \end{array} \right.$  This problem, thanks to theorem [thr:chmixte:2] is well-posed if and only if $X_h$ and $M_h$ are such that there exists $\beta_h > 0$ $\inf_{q_h \in M_h} \sup_{v_h \in X_h} \frac{\int_\Omega q_h \nabla \cdot v_h}{||v_h||_{X_h} ||q_h||_{M_h}} \geq \beta_h$ #### 7.5.3. Some counter examples: bad finite element for Stokes In this section, we present two classical bad finite element approximations. ##### Finite element $\poly{P}_1/\poly{P}_0$: locking Thanks to the Euler relations, we have $\begin{array}[c]{rl} N_{\mathrm{cells}} - N_{\mathrm{edges}} + N_{vertices} &= 1-I\\ N^\partial_{\mathrm{vertices}} - N^\partial_{\mathrm{edges}} &= 0 \end{array}$ where $I$ is the number of holes in $\Omega$. We have that $\dim {M_h} = N_{\mathrm{cells}}$,$\dim {X_h} = 2 N^i_{\mathrm{vertices}}$ and so $\dim {M_h} - \dim {X_h} = N_{\mathrm{cells}} - 2 N^i_{\mathrm{vertices}} = N^\partial_{\mathrm{edges}} - 2 > 0$ so $M_h$ is too rich for the condition and we have $\ker(\mathcal{B}) = \{0\}$ such that the only discrete $u_h^*$, with components $U^*$, satisfying $\mathcal{B} U^*$ is the null field, $U^*=0$. ##### Finite element $\poly{Q}_1/\poly{P}_0$: spurious mode We can construct in that case a function $q_h^*$ on a uniform grid which is equal alternatively -1, +1 (chessboard) in the cells of the mesh, then $\forall v_h \in [Q^1_{c,h}]^d, \quad \int_\Omega q^*_h \nabla \cdot v_h = 0$ and thus the associated $X_h$, $M_h$ do not satisfy the condition. ##### Finite element $\poly{P}_1/\poly{P}_1$: spurious mode We can construct in that case a function $q_h^*$ on a uniform grid which is equal alternatively -1, 0, +1 at the vertices of the mesh, then $\forall v_h \in [P^1_{c,h}]^d, \quad \int_\Omega q^*_h \nabla \cdot v_h = 0$ and thus the associated $X_h$, $M_h$ do not satisfy the condition. #### 7.5.4. Mini-Element The problem with the $\poly{P}_1/\poly{P}_1$ mixed finite element is that the velocity is not rich enough. A cure is to add a function $v_h^*$ in the velocity approximation space to ensure that $\int_\Omega q^*_h \nabla \cdot v_h^* \neq 0$ where $q_h^*$ is the spurious mode. To do that we add the bubble function to the $\poly{P}_1$ velocity space. Definition 46. Recall the construction of finite elements on a reference convex $\hat{K}$. We say that $\hat{b}: \hat{K} \mapsto \RR$ is a bubble function if: * $\hat{b} \in H^1_0(\hat{K})$ * $0 \leq \hat{b}(\hat{x}) \leq 1, \quad \forall \hat{x} \in \hat{K}$ * $\hat{b}(\hat{C}) = 1, \quad \mbox{where} \hat{C}$ is the barycenter of $\hat{K}$ Example The function $\hat{b} = (d+1)^{d+1} \Pi_{i=0}^d\ \hat{\lambda}_i$ where $(\hat{\lambda}_0, \ldots, \hat{\lambda}_d)$ denote the barycentric coordinates on $\hat{K}$ Denote now $\hat{b}$ a bubble fonction on $\hat{K}$, we set $\hat{P} = [\poly{P}_1(\hat{K}) \oplus \mathrm{span} (\hat{b})]^d,$ and introduce Aligned 1. X_h &=& \Big\{ v_h \in [C0(\bar{\Omega})]d : \forall K \in \mathcal{T}_h, v_h \circ T_K \in \hat{P}; v_{h_|{\partial \Omega}} = 0 \Big\}\\ M_h &=& P^1_{c,h} Lemma 3. The spaces $X_h$ and $M_h \cap L^2_0(\Omega)$ satisfy the compatibility condition uniformly in $h$. Theorem 20. Suppose that $(u,p)$, solution of problem [prob:chmixte:1], is smooth enough, ie. $u \in [H^2(\Omega)$^d \cap [H1_0(\Omega)]d] and $p\in H^1(\Omega) \cap L^2_0(\Omega)$. Then there exists a constant $c$ such that for all $h >0$ $\| u- u_h \|_{1,\Omega} + \|p-p_h\|_{0,\Omega} \leq c h (\|u\|_{2,\Omega} + \|p\|_{1,\Omega})$ and if the Stokes problem is stabilizing then $\|u-u_h\|_{0,\Omega} \leq c h^2 ( \|u\|_{2,\Omega} +\|p\|_{1,\Omega}).$ Definition 47: Stabilizing Stokes problem We say that the Stokes problem is stabilizing if there exists a constant $c_S$ such that for all $f \in [L^2(\Omega)$^d], the unique solution $(u,p)$ of with $g=0$ is such that: $\|u\|_{2,\Omega} + \|p\|_{1,\Omega} \leq c_S \|f\|_{0,\Omega}$ A sufficient condition for stabilizing Stokes problem is that the $\Omega$ is a polygonal convex in 2D or of class $C^1$ in $\RR^d, d=2,3$. #### 7.5.5. Taylor-Hood Element The mini-element solved the compatibility condition problem, but the error estimation in equation is not optimal in the sense that 1. the pressure space is sufficiently rich to enable a $h^2$ convergence in the pressure error, 2. but the velocity space is not rich enough to ensure a $h^2$ convergence in the velocity error. The idea of the Taylor-Hood element is to enrich even more the velocity space to ensure optimal convergence in $h$. Here we will take $[\poly{P}_2$^d] for the velocity and $\poly{P}_1$ for the pressure. Introduce \begin{aligned} \label{eq:chmixte:39} X_h &=& [P^2_{c,h}]^d\\ M_h &=& P^1_{c,h} \end{aligned} Lemma 4. The spaces $X_h$ and $M_h \cap L^2_0(\Omega)$ satisfy the compatibility condition uniformly in $h$. Theorem 21. Suppose that $(u,p)$, solution of problem Problem 20, is smooth enough, ie. $u \in [H^3(\Omega)$^d \cap [H1_0(\Omega)]d] and $p\in H^2(\Omega) \cap L^2_0(\Omega)$. Then there exists a constant $c$ such that for all $h >0$ $\| u- u_h \|_{1,\Omega} + \|p-p_h\|_{0,\Omega} \leq c h^2 (\|u\|_{3,\Omega} + \|p\|_{2,\Omega})$ and if the Stokes problem is stabilizing then $\|u-u_h\|_{0,\Omega} \leq c h^3 ( \|u\|_{3,\Omega} +\|p\|_{2,\Omega}).$ Generalized Taylor-Hood element We consider the mixed finite elements $\poly{P}_k/\poly{P}_{k-1}$ and $\poly{Q}_k/\poly{Q}_{k-1}$ which allows to approximate the velocity and pressure respectively with, on Simplices \begin{aligned} \label{eq:chmixte:42} X_h &=& [P^{k}_{c,h}]^d\\ M_h &=& P^{k-1}_{c,h} \end{aligned}] On Hypercubes latexmath:[\begin{aligned} \label{eq:chmixte:43} X_h &=& [Q^{k}_{c,h}]^d\\ M_h &=& Q^{k-1}_{c,h} \end{aligned} We then have $\|u-u_h\|_{0,\Omega} + h ( \| u- u_h \|_{1,\Omega} + \|p-p_h\|_{0,\Omega} ) \leq c h^{k+1} (\|u\|_{k+1,\Omega} +\|p\|_{k,\Omega})$ There are other stable discretization spaces • Discrete inf-sup condition: dictates the choice of spaces • Inf-sup stables spaces: • $\mathbb Q_k$-$\mathbb Q_{k-2}$, $\mathbb Q_k$-$\mathbb Q^{disc}_{k-2}$ • $\mathbb P_k$-$\mathbb P_{k-1}$, $\mathbb P_k$-$\mathbb P_{k-2}$, $\mathbb P_k$-$\mathbb P^{disc}_{k-2}$ • Discrete inf-sup constant independent of $h$, but dependent on $k$ ##### Numerical validation: Test case We consider the Kovasznay solution of the steady Stokes equations, see Kovasznay . The exact solution is  $\begin{array}{r c l} \mathbf{u}(x,y) & = & \left(1 - e^{\lambda x } \cos (2 \pi y), \frac{\lambda}{2 \pi} e^{\lambda x } \sin (2 \pi y)\right)^T \\ p(x,y) & = & -\frac{e^{2 \lambda x}}{2} \\ \lambda & = & \frac{1}{2 \nu} - \sqrt{\frac{1}{4\nu^2} + 4\pi^2}. \end{array} — The domain is defined as $\domain = (-0.5,1) \times (-0.5,1.5)$ and $\nu = 0.035$. The forcing term for the momentum equation is obtained from the solution and is \[ \mathbf{f} = \left( e^{\lambda x} \left( \left( \lambda^2 - 4\pi^2 \right) \nu \cos (2\pi y) - \lambda e^{\lambda x} \right), e^{\lambda x} \nu \sin (2 \pi y) (-\lambda^2 + 4 \pi^2) \right)^T$ Dirichlet boundary conditions are derived from the exact solution. :leveloffset!: = Solving Linear Systems NOTE: documentation pending [appendix] = Bibliography bibliography::[]
(46)
+ 1. la transformation et son inverse sont $\mathcal{C}^1$ et bijectives 2. An euclidian space which is complete for the norm induced by the scalar product Version 1.0 Last updated 2017-03-29 15:09:20 CEST