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.

1. Unisolvance

Définition

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

est bijective.

Preuve:

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

2. Elément fini de Lagrange

É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\)).

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

Définition

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. Exemples d’éléments finis de Lagrange

3.1. Espaces de polynômes

On notera \(\Pk{k}\) l’espace vectoriel des polynômes de degré total inférieur ou égal à \(k\).

  • Sur , \(\Pk{k}=\hbox{Vect} \{ 1,X,\ldots, X^k\}\quad\) et dim \(\Pk{k}=k+1\).

  • Sur \(^2\), \(\Pk{k}=\hbox{Vect} \{ X^i Y^j ; 0\le i+j \le k\}\quad\) et dim \({ \Pk{k}=\frac{(k+1)(k+2)}{2}}\).

  • Sur \(^3\), \(\Pk{k}=\hbox{Vect} \{ X^i Y^j Z^l ; 0\le i+j+l \le k\}\quad\) et dim \({ \Pk{k}=\frac{(k+1)(k+2)(k+3)}{6}}\).

On notera \(\Qk{k}\) l’espace vectoriel des polynômes de degré inférieur ou égal à \(k\) par rapport à chaque variable.

  • Sur \(\RR\), \(\Qk{k}=\Pk{k}\).

  • Sur \(\RR^2\), \(\Qk{k}=\hbox{Vect} \{ X^i Y^j ; 0\le i,j \le k\}\quad\) et dim \({ \Qk{k}=(k+1)^2}\).

  • Sur \(\RR^3\), \(\Qk{k}=\hbox{Vect} \{ X^i Y^j Z^l ; 0\le i,j,l \le k\}\quad\) et dim \({ \Qk{k}=(k+1)^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}\)
  • \(K=[a,b\)]

  • \({\Sigma=\{a+i\,\frac{b-a}{m},\quad i=0,\ldots,k\}}\)

  • \(P=\Pk{k}\)

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

elements 2D

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

elements 3D

4. Famille affine d’éléments finis

É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\).

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}\).

Définition

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

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++

feelpp mesh

Définition

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}\).

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é.

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].

Définition

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

É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

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

et en particulier on a

\[ T_K(\hat{g}_i) = g^K_i, \quad \forall i \in \set{1,\ldots,\ngeo}\]
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.

Définition: 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<Simplex<d,1> > ou Mesh<Simplex<d> > 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.

5.2. Quelques calculs avec la transformation géométrique

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

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

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

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

5.2.3. 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 \(\hat{{\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\]
\[\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\]
\[\int_{\partial K} f( x ) dx = \int_{\partial \hat{K}} \hat{f}(\xi) \| B_K(\xi) {\hat{\mathbf{n}}}(\xi) \| J_K( \xi ) d \xi\]
\[\int_{\partial K}\ {\mathbf{F}}( x )\ \cdot\ {\mathbf{n}}(x) dx = \int_{\partial \hat{K}} {\hat{\mathbf{F}}}( \xi )\ \cdot \Big(B_K(\xi) {\hat{\mathbf{n}}}(\xi) \Big) J_K( \xi )\ d \xi\]

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.
Définition: 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.

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:23].

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}\]
Notations

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

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

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

Définition: 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\}\]

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

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

Espace \(H^1\)-conforme sur des simplexes (e.g. triangle)
\[ 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\}\]
Espace \(H^1\)-conforme sur des hypercubes (e.g quadrilatère)
\[ 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\}\]
Table 1. Dimension des espaces \(P_{\mathrm{c}, h}^k\) et \(Q_{\mathrm{c}, h}^k\) en dimension 2 et 3 pour \(k \in\{1,2,3\}\).
espace d k=1 k=2 k=3

\(P_{c,h}^k\)

2

\(N_{so}\)

\(N_{so}+N_{ar}\)

\(N_{so}+2N_{ar}+N_{ma}\)

\(Q_{c,h}^k\)

2

\(N_{so}\)

\(N_{so}+N_{ar}+N_{ma}\)

\(N_{so}+2N_{ar}+4N_{ma}\)

\(P_{c,h}^k\)

3

\(N_{so}\)

\(N_{so}+N_{ar}\)

\(N_{so}+2N_{ar}+N_{fa}\)

\(Q_{c,h}^k\)

3

\(N_{so}\)

\(N_{so}+N_{ar}+N_{fa}+N_{ma}\)

\(N_{so}+2N_{ar}+4N_{fa}+8N_{ma}\)

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

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}

Pour calculer \(\Pi^{0,k}_{c,h}\) et \(\Pi^{1,k}_{c,h}\) on a besoin de résoudre les problèmes

Problème: 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}\]
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}\]
L’opérateur \(\Pi_{\mathrm{c}, b}^{1, k}\) est appelé projecteur elliptique. On peut également considérer la projection orthogonale de \(H_0^1(\Omega)\) dans \(P_{\mathrm{c}, b, 0}^k\) avec le produit scalaire \(\int_{\Omega} \nabla v \cdot \nabla w\). Ce dernier satisfait les mêmes propriétés que celles énoncées ci-dessus pour l’opérateur \(\Pi_{\mathrm{c}, b}^{1, k}\).
Des résultats analogues existent pour des opérateurs de projection orthogonale sur les espaces \(Q_{c, b}^k\).

6.2. Interpolant de Lagrange

Définition

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)\]
Théorème

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

De plus, pour tout \(m \in\{2, \ldots, k+1\}\), on a

\[\begin{equation*} \left(\sum_{K \in \mathcal{T}_b}\left|v-\mathcal{I}_b^{\mathrm{Lag}} v\right|_{m, K}^2\right)^{\frac{1}{2}} \leqslant c h^{k+1-m}|v|_{k+1, \Omega_b} . \end{equation*}\]

Ce théorème est une conséquence directe du théorème d’interpolation local. On notera que l’hypothèse portant sur la régularité de la famille \(\left\{\mathcal{T}_h\right\}_{h>0}\) permet de contrôler uniformément les rapports \(\varpi_K\). Enfin, lorsque la fonction à interpoler n’est pas suffisamment régulière, par exemple si \(v \in H^s\left(\Omega_h\right)\) et \(v \notin H^{s+1}\left(\Omega_h\right)\) pour un entier \(s \in\{1, \ldots, k\}\), l’estimation [eq.62] ci-dessus devient

\[\begin{equation*} \left\|v-\mathcal{I}_b^{\mathrm{Lag}} v\right\|_{0, \Omega_b}+h\left|v-\mathcal{I}_b^{\mathrm{Lag}} v\right|_{1, \Omega_b} \leqslant c h^s|v|_{s, \Omega_b} \end{equation*}\]
Example 1. Vérification du théorème précédent sur un exemple

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

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:62] en échelle log-log.

Nous devons obtenir des droites de pentes \(k\) (resp. \(k+1\)) pour la norme \(L^2\) (resp. \(H^1\).)

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

Définition
  • Si \(\hat{P}_{\mathrm{geo}} = \hat{P}\) l’approximation est dite iso-parametrique.

  • Si \(\hat{P}_{\mathrm{geo}} \subset \hat{P}\) l’approximation est dite sub-parametrique.

  • Si \(\hat{P} \subset\hat{P}_{\mathrm{geo}}\) l’approximation est dite sur-parametrique.

Gmsh cite:[gmsh] un mailleur libre permet de générer des maillage d’ordre élevé jusqu’à l’ordre 5 en 2D et 4 en 3D
Théorème

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

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.

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

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

Fonctions de base globales
\[{\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\]

L’espace d’approximation interne est donc alors :

\[V_h = \hbox{Vect }\left\{\varphi_1,\ldots,\varphi_{N_h}\right\}\]
Exemple de fonction de base globale (élément triangulaire \(P_1\))

fonction globale

On remarque 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

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.

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.


1. la transformation et son inverse sont \(\mathcal{C}^1\) et bijectives