Le Laplacien

On s’instéresse dans cette section à l’approximation élément fini conforme du problème suivant:

Problème: 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}\]

1. Cadre Mathematique

On suppose que \(f \in L^2(\Omega)\).

La formulation faible du problème Problème: Formulation forte du Laplacien est la suivante:

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)\]

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, \RR)\) telle que \(a(u,v) = \int_\Omega \nabla u \cdot \nabla v \)

  • la forme linéaire \(\ell \in \mathcal{L}(V, \RR)\) telle que \(l(v) = \int_\Omega f \nabla v \)

Le problème Formulation faible pour des conditions de Dirichlet homogènes s’écrit sous forme abstraite

Problème:

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

Soit \(\Omega\) un ouvert borné de \(\R{d}\). 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

La forme bilinéaire \(a\) du problème Formulation faible pour des conditions de Dirichlet homogènes est coercive

Preuve

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 [thr:12] permet alors de conclure sur l’existence d’une solution unique pour le problème Formulation faible pour des conditions de Dirichlet homogènes.

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 [thr:16].

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

Problème

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:

Théorème

On suppose que \(u\) solution de Formulation faible pour des conditions de Dirichlet homogènes 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}\]
Preuve

La preuve de ([eq:73]) est obtenue grâce au lemme de Cea ([eq:cea]) et du théorème d’interpolation ([eq:62]).

La preuve de ([eq:74]) 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:74]) se déduit de ([eq:74]).

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:64]) 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 $\R{d}$ telle que la transformation
// géométrique $T_K,\ K \in \calTh$  soit affine
typedef Mesh<Simplex<d,1> > mesh_type;
// un maillage d'hypercube dans $\R{d}$ telle que la transformation
// géométrique $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 $[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 Problème

  a.solve( _rhs=l, _solution=u );

Le listing complet

5. Conditions aux limites

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:

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:

Problème

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.

Théorème

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

5.2. 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:

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

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 [prob:13] est bien posé grâce à Lax-Milgram.

On peut approcher le problème [prob:13] par des éléments finis de Lagrange.

On utilise la formulation développée dans la section Approximation conforme

Problème

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 Problème est identique par le théorème Théorème aux mêmes estimations que dans le cas Dirichlet homogène.

5.2.1. Cas \(\mu=0\)

Le cas \(\mu=0\) présente quelques difficultés techniques.

On a

Problème

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 Problème 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 Problème s’écrit alors

Problème

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 Problème résulte de la coercivité de \(a\) sur \(H^1_*(\Omega)\) qui elle-même résulte du Lemme suivante

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 Problème peut être approché par des éléments finis de Lagrange ce qui conduit au problème discret suivant

Problème

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 Problème 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

a.solve(_rhs=l,_solution=u);

double meanu = mean(_range=elements(mesh),_expr=idv(u))(0,0);

u.add(-meanu);

5.3. 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:

Problème

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

Problème

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 Problème est bien posé grâce à Lax-Milgram.

On considère le problème discret suivant

Problème

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 Théorème.

Considérons \(\Omega=[0,1\)^2] et les données \(\mu=0.01\), \(f=1\) et \(g=0\).

6. Et après ?

Les étapes suivantes sont par exemple: