UNIVERSITE DE KINSHASA

Faculté des Sciences et Technologie

Mention Maths, Stat et Informatique

' $

Méthodes des Différences

Finies et des Éléments Finis
Cours de Mathématiques Appliquées
& %

l
Dr Joseph Désiré Bukweli Kyemba
Professeur

ANNEE ACADEMIQUE 2023 - 2024

 i

Objectifs :
Cette UE permet à l’étudiant d’acquérir les notions de mathématiques appliquées.
De manière spécifique, à la fin de cet EC, l’étudiant sera capable de :
– Résoudre les problèmes aux valeurs de bords,
– Résoudre les problèmes de diffusion et du potentiel électrique
à l’aide des méthodes des Différences Finies et/ou des Éléments Finies.

Approches pédagogiques :
Dans cette UE l’enseignement et l’apprentissage utiliseront les approches ci-après :
dans
– Exposés magistraux interactifs ;
– Travaux individuels et en groupes ;
– Travaux dirigés ;
– Utilisation des TIC pour les travaux pratiques.

Modalités d’évaluation :
De préférence, après chaque module (chapitre), l’apprenant sera soumis à un travail à
domicile en groupe suivi d’un exposé et à une interrogation.
Un examen sur toute la matière sera proposé à la fin de l’UE.

Bibliographie
1. M. A. Aziz Alaouti et C. Bertelle, Méthodes numériques appliquées, Paris, 2002.
2. R. L. Berden and J. D. Faires, Numerical Analysis, 9th ed., Books/Cole, Cengage
Learning, Boston, 2011.
3. S. D. Conte and C. de Boor, Elementary numerical Analysisis : An algorithmic
approach, 3rd ed., McGraw-Hill Book Company, New York, 1980.
4. G. Legendre, Méthodes numériques : Introduction à l’analyse et calcul scien-
tifique, Université Paris Dauphine, inédit, 2013.
5. S. Nakamura, Numerical Analysis and Graphic Visualization with MATLAB, 2nd
ed., Prentince-Hall, New Jersey, 2002.
6. W. F. Trench, Elementary Differential equations with boundary value problems,
Brooks/Cole Thomson Learning, 2001.
7. C. J. Zarowski, An introduction to numerical analysis for electrical and computer
engineers, John & Sons Inc., New Jersey, 2004.

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 Table des matières

Table des matières . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1

1 Méthode des différences finies . . . . . . . . . . . . . . . . . . . . . . . . . . 2

1.1 Quelques exemples de différences finies . . . . . . . . . . . . . . . . . . . 2
1.1.1 Différences finies décentrées . . . . . . . . . . . . . . . . . . . . . . . 2
1.1.2 Différences finie centrées . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.2 Méthode des différences finies pour P.V.B linéaires . . . . . . . . . . . . . 3
1.3 Méthode des différences finies pour des EDPs . . . . . . . . . . . . . . . 5
1.3.1 Problème elliptique . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.3.2 Problème parabolique . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
Exercices sur le chapitre . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

2 Méthodes des éléments finis . . . . . . . . . . . . . . . . . . . . . . . . . . . 10

2.1 Méthode des résidus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.2 Méthode des résidus pondérés . . . . . . . . . . . . . . . . . . . . . . . . 13
2.3 Généralisation de la méthode des résidus pondérés . . . . . . . . . . . . 15
2.4 Formulation variationnelle résiduelle pondérée . . . . . . . . . . . . . . . 17
2.5 Eléments finis de Lagrange . . . . . . . . . . . . . . . . . . . . . . . . . . 24
Exercices sur le chapitre . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 Chapitre 1

Méthode des différences finies

La plupart des fonctions ne peuvent pas être représentées par des formules analytiques.
Il est donc impossible de calculer leurs dérivées, même avec un code de calcul formel.
La méthode des différences finies propose un moyen de calculer une approximation
numérique des valeurs des dérivées d’une fonction.

1.1 Quelques exemples de différences finies

On considère une fonction régulière u : [a, b] −→ R. On se donne un entier I > 2, on
pose h = (b − a)/(I) et on définie les points :

t0 = a, t1 = t0 + h . . . , ti = t0 + i h . . . , tI = b.

On dit que l’on discrétisé l’intervalle à l’aide d’une grille définies par se sommets

ti = a + ih, i = 0, 1, · · · , I

avec un pas constant h. On parle d’une grille uniforme.

Les dérivées u0 et u00 sont approximées par leurs différences finies.
Pour un point quelconque ti (1 ≤ i ≤ I − 1 ) à l’intérieur de l’intervalle ]a, b[, on écrit les
développements de Taylor à l’ordre n avec reste autour de ti :

h2 00 hn hn+1 (n+1) +
u(ti+1 ) = u(ti ) + hu0 (ti ) + u (ti ) + · · · + u(n) (ti ) + u (ξi ), (1.1.1)
2 n! (n + 1)!
où ti < ξ+i < ti+1 ; et

h2 00 (−h)n (n) (−h)n+1 (n+1) −

u(ti−1 ) = u(ti ) − hu0 (ti ) + u (ti ) + · · · + u (ti ) + u (ξi ), (1.1.2)
2 n! (n + 1)!
où ti−1 < ξ−i < ti+1 .

1.1.1 Différences finies décentrées

De (1.1.1) on obtient
1 h hn−1 (n) hn
[u(ti+1 ) − u(ti )] = u0 (ti ) + u00 (ti ) + · · · + u (ti ) + u(n+1) (ξ+i ).
h 2 n! (n + 1)!
D’où, on tire une approximation de u0 (ti ) appelée différence décentrée à droite
1
u0 (ti ) = [u(ti+1 ) − u(ti )] (1.1.3)
h

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 3

De même, de (1.1.2) on obtient

1 h hn−1 (n) hn
[u(ti ) − u(ti−1 )] = u0 (ti ) − u00 (ti ) − · · · − u (ti ) + u(n+1) (ξ+i ).
h 2 n! (n + 1)!
D’où, on tire une approximation de u0 (ti ) appelée différence décentrée à gauche
1
u0 (ti ) = [u(ti ) − u(ti−1 )] (1.1.4)
h
Ces deux formules décentrées d’approximation de u0 sont de l’ordre 1 si u00 (ti ) , 0 ou
de l’ordre 2 si u00 (ti ) = 0.

1.1.2 Différences finie centrées

En soustrayant membres à membres les égalités (1.1.1) et (1.1.2) et divisant le résultat
obtenu par 2h, on
1 h3
[u(ti+1 − u(ti−1 ] = u0 (ti ) + u(3) (ti ) + · · ·
2h 6
0
On a ainsi une formule d’approximation de u dite de différence finie centrée :
1
u0 (ti ) = [u(ti+1 ) − u(ti−1 )] (1.1.5)
2h
qui est de l’ordre 2 si u(3) (ti ) , 0 ou de l’ordre au moins 4 sinon.
En additionnant membres à membres les deux dernières égalités et en tirant u00 (ti ),
obtient
1 h2 (4)
u00 (ti ) = [u(ti+1 ) − 2u(ti ) + u(ti−1 )] − u (ti ) − · · ·
h2 12
D’où, une formule d’approximation de u00 (ti ) appelée de différence finie centrée :
1
u00 (ti ) = [u(ti+1 ) − 2u(ti ) + u(ti−1 )] (1.1.6)
h2
qui est de l’ordre 2 si u(4) (ti ) , 0 et de l’ordre au moins 6 sinon.
NB : Dire qu’un formule d’approximation est d’ordre k signifie que l’erreur est une
fonction o(hk ) dont la limite est égale à 0 pour h tendant vers 0. Plus l’ordre est élevé
plus l’approximation est bonne.

1.2 Méthode des différences finies pour P.V.B linéaires

Considérons maintenant le P.V.B. linéaire suivant :
−u00 (t) + p(t)u0 (t) + q(t)u(t) = g(t), a ≤ t ≤ b,




(1.2.1)



 u(a) = α,

u(b) = β.

On suppose que g est de classe C2 , ce qui entraı̂ne u est de classe C4 . On discrétise le

domaine [a, b] à l’aide d’une grille définie par les sommets

t0 = a, t1 = t0 + h . . . , ti = t0 + i h . . . , tI = b.

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 4

On écrit alors l’EDO en chaque sommet internes i = 1, . . . , I − 1

−u00 (ti ) + p(t)u0 (ti ) + q(t)u(ti ) = g(ti )

Puisque u est solution du P.V.B (1.2.1), les différences finies centrés (1.1.5) et (1.1.6) nous
permettent d’obtenir :

1 1
− 2
[u(ti+1 ) − 2u(ti ) + u(ti−1 )] + p(ti ) [u(ti+1 ) − u(ti−1 )] + q(ti ) u(ti ) = g(ti )
h 2h
ou encore
! !
h   h
− 1 + p(ti ) u(ti−1 ) + 2 + h q(ti ) u(ti ) − 1 − p(ti ) u(ti+1 ) = h2 g(ti ).
2
(1.2.2)
2 2

On cherche donc des approximations ui de u(ti ), telle que

! !
h   h
− 1 + p(ti ) ui−1 + 2 + h2 q(ti ) ui − 1 − p(ti ) ui+1 = h2 g(ti ), (1.2.3)
2 2

avec u0 = u(t0 ) = α et uI = u(tI ) = β. Ceci revient à résoudre le système linéaire

Au = b (1.2.4)

où

λ1 µ1 ... ...
 
 0 0  
 u1
 
b1

.. ..
    
ν2 λ2 µ2 . .
 2
 u2  h g(t2 )
    
   
. ... ... ... ..  .. ..
   
0 ..
 
.  . .
    
A =  u =  . et b = 
    
.. . .  ; ..
... ... ...
  
 . . 0

  .. 
  .


 
.. ..
   
 u  h2 g(t )
. νI−2 λI−2 µI−2

.
 



  I−2   I−2 
uI−1 bI−1
   
. . . . . . 0 νI−1 λI−1

0


avec λi = 2 + h2 p(ti ) ; νi = −1 − h2 p(ti ) ; µi= −1 + h2 p(ti ) ; 

b1 = 1 + h2 p(t1 ) α + h2 g(t1 ) et bI−1 = 1 − h2 p(tI−1 )β + h2 g(tI−1 )
On peut alors reconstruire une approximation uh (t) de u(t) sur [a, b] de la manière
suivante :

u0 = α t ∈ [a, a + h/2[








uh (t) =  t ∈ [ti − h/2, ti + h/2[, i = 0, . . . , I − 1

ui (1.2.5)







 u = β t ∈ [b − h/2, b]

I

L’erreur de discrétisation est dans ce cas

Z b !1/2
2
εh = uh (t) − u(t) dt .
a

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 5

Remarque 1.2.1
Si une des conditions de bords du P.V.B sont données par u0 (a) = α ou u0 (b) = β on utilisera
au besoin l’une des formules de différences finies décentrées (1.1.1) ou (1.1.2) et on ajoutera aux
systèmes précédents éventuellement l’une ou les deux équations

1
u0 (a) = [u1 − u0 ] = α
h
1
u0 (b) = [uI − uI−1 ] = β
h
où u0 et uI seront des variables éventuelles. Et dans ce cas la méthode sera d’ordre 1 au moins.
On peut aussi utiliser les différences centrées

1
u0 (a) = [u1 − u−1 ] = α
2h
1
u0 (b) =
[uI+1 − uI−1 ] = β
2h
où u−1 et uI+1 sont appelés sommets fantômes, pour augmenter la performance.

1.3 Méthode des différences finies pour des EDPs

Un exemple d’ EDP elliptique est l’équation de Poisson

1.3.1 Problème elliptique

Considérons, par exemple, le problème elliptique

1

 ∇2 u(x, y) ≡ uxx + u yy = − ρ(x, y), (x, y) ∈ D,


 ε (1.3.1)
 u(x, y) = g(x, y) pour (x, y) ∈ ∂D,


où D = (x, y) | a ≤ x ≤ b, c ≤ y ≤ d avec ρ et g continues dans leurs domaines respec-



tifs.
Nous utilisons la méthode de différences finies pour trouver une approximation de la
solution. Pour cela on se donne deux entiers I et J et on définit les pas de discrétisations
h = (b − a)/I et k = (d − c)/J. Ce qui revient à diviser l’intervalle [a, b] en I parties de
taille h et [c, d] en J parties de taille k ; ou encore à placer une grille sur la surface D en
traçant des lignes verticales et horizontales aux points des cordonnées (xi , yi ), avec

xi = a + i h, i = 0, 1, . . . , I et y j = c + j k, j = 0, 1, . . . , J.

Les lignes x = xi et y = y j sont appelées les lignes de la grille et leurs intersections les
points ou sommets de grillage.
Pour chaque sommet intérieur de la grille c’est-à -dire (xi , y j ), i = 0, 1, . . . , I − 1 et
j = 0, 1, . . . , J − 1 on a les formules de différences finies centrées

∂2 u u(xi+1 , y j ) − 2u(xi , y j ) + u(xi−1 , y j ) h2 ∂4 V

(x i , y j ) = − (ξi , y j ),
∂x2 h2 12 ∂x4

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 6

où ξi ∈ [xi−1 , xi+1 ] et

∂2 u u(xi , y j+1 ) − 2u(xi , y j ) + u(xi , y j−1 ) k2 ∂4 u

(x i , y j ) = − (xi , η j ),
∂y2 k2 12 ∂y4

où η j ∈ [y j−1 , y j+1 ]. En introduisant ces approximations dans l’équation de Poisson, on

obtient
u(xi+1 , y j ) − 2u(xi , y j ) + u(xi−1 , y j ) u(xi , y j+1 ) − 2u(xi , y j ) + u(xi , y j−1 ) ρ(x j , y j )
+ =−
h 2 k 2 ε
h 2 ∂4 u k 2 ∂4 u
+ (ξi , y j ) + (xi , η j ), (1.3.2)
12 ∂x4 12 ∂y4

pour i = 1, . . . , I − 1 et j = 1, . . . , J − 1 ; avec les conditions de bord données par

u(a, y j ) = g(a, y j ) et u(b, y j ) = g(b, y j ) pour j = 0, 1, . . . , J;

u(xi , c) = g(xi , c) et u(xi , d) = g(xi , d) pour i = 1, . . . , I − 1.

Pour chaque couple (i, j) on cherche une approximation ui j de u(xi , x j ) tel que

ui+1 j − 2ui j + ui−1 j ui j+1 − 2ui j + ui j−1 ρ(x j , y j )

+ = −
h2 k2 ε
ou encore
!
h2   h2  
2 1 + 2 ui j − ui+1 j + ui−1 j − 2 ui j+1 − 2ui j + ui j−1 = h2 ρ(xi , x j ), (1.3.3)
k k

pour i = 1, . . . , I − 1 et j = 1, . . . , J − 1, et

u0 j = g(a, y j ) et uI j = g(b, y j ) pour j = 0, 1, . . . , J;

ui 0 = g(xi , c) et ui J = g(xi , d) pour i = 1, . . . , I − 1.

Ceci nous donne un système linéaire de (I − 1)(J − 1) équations à (I − 1)(J − 1) inconnues
ui j avec i = 1, . . . , I − 1 et j = 1, . . . , J − 1. Généralement on pose

ui, j = ui+(j−1)(I−1) , 1 ≤ i ≤ I − 1, 1≤ j≤ J−1

pour avoir le système Au = b, où u = (u1 u2 · · · u(I−1)(J−1) )>

Notons en passant que l’erreur de troncature locale est de l’odre de O(h2 + k2 ).

1.3.2 Problème parabolique

On considère l’EDP parabolique (l’équation de la chaleur ou de diffusion)

α2 uxx = ut (t, x) − f (x), t > 0, a < x < b (1.3.4)

soumise aux conditions typiques de bord données par

u(t, a) = 0, u(t, b) = 0 ∀ t > 0; (1.3.5)

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 7

et la condition initiale

u(0, x) = g(x) pour a < x < b. (1.3.6)

On discrétise le domaine :

Ω = {(t, x) : t ≥ 0, a ≤ x ≤ b}

à l’aide d’une grille de pas quelconque h sur l’axe temporel et k = (b − a)/J

txi = i h, i = 0, 1, . . . et x j = a + j k, j = 0, 1, . . . , J.

En chaque sommet interne de la grille l’EDP s’écrit

∂u ∂u
α2 (ti , x j ) = (ti , x j ) − f (x j ), i = 1, 2, > 0 · · · , j = 1, · · · , J − 1
∂x2 ∂t
En utilisant les différences finies centrées i = 1, 2, · · · , j = 1, · · · , J − 1, on a :

α2 h i 1 h i
u(ti , x j+1 ) − 2u(ti , x j ) + u(ti , x j−1 ) = u(ti+1 , x j ) − u(ti−1 , x j ) − f (x j ),
h2 2h
ou encore pour i = 1, 2, · · · , j = 1, · · · , M − 1
h i hh i
α2 ui, j+1 − 2ui j + ui, j−1 = ui+1, j − ui−1,j − h2 f (x j ),
2
ou encore pour i = 1, 2, · · · , j = 1, · · · , M − 1

h h
ui−1, j + α2 ui,j−1 − 2α2 ui j + α2 ui, j+1 − ui+1, j = h2 f (x j ). (1.3.7)
2 2
La discrétisation des conditions de bords donne

uti ,0 = 0, uti ,J = 0, i = 1, 2, · · ·

et celle de la condition initiale donne

u0, j = g(x j ), j = 0, 1, · · · , J − 1.

Remarque 1.3.1 : Discrétisation partielle

On choisit de discrétiser la fonction inconnue dans la variable d’espace x tout en conservant le
caractère continu du temps t. Pour cela on introduit les fonctions ui (t) qui approchent u(xi , t)
pour tout i.

Les ui (t) ne seront jamais calculés en pratique mais sont des intermédiaires pour dériver
les schémas numériques. En chaque sommet intérieur de la grille

x j = a + j k, j = 0, 1, . . . , J

l’EDP (1.3.4) s’écrit alors

α2 h i
u j+1 (t) − 2u j (t) + u j−1 (t) = u0j (t) − f (x j ), t > 0, j = 1, . . . , J − 1
h2

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 8

ou encore
α2 h i
u0j (t) − u j−1 (t) − 2u j (t) + u j+1 (t) = f (x j ), t > 0, j = 1, . . . , J − 1 (1.3.8)
h2
qui peut se mettre sous forme matricielle suivante :

u0 (t) + Au(t) = b
 >  >
où u(t) = u1 (t) u2 (t) · · · u J−1 (t) et b = f (x1 ) f (x2 ) · · · f (x j−1 ) en tena compte des
conditions aux limites
u0 (t) = 0, u J (t) = 0, t > 0.
Le vecteur condition initiale (1.3.6) et u(0) avec comme composantes :

ui (0) = g(xi ), j = 1, . . . , J − 1 (1.3.9)

Exercices sur le chapitre

1. On donne le PVB u00 = 4(u − t), 0 ≤ t ≤ 1, u(0) = 0, u(1) = 2.
(a) Montrer que ce PVB admet une solution unique que l’on déterminera.
(b) Donner une approximation de cette solution en utilisant le pas de discrétisation
h = 0.25.
N.B. : Utiliser si possible un script MATLAB pour trouver l’approximation
2. Mêmes questions pour les PVB suivants :
(a) u00 = u0 + 2u + 2t + 3, 0 ≤ t ≤ 1, u(0) = 2, u(π/2) = 1 avec h = 0.1.
(b) u00 = u0 +2u+cos t, 0 ≤ t ≤ π/2, u(0) = −0.3, u(π/2) = −0.1 avec h = π/4.
(c) u00 + e−tu + sin u0 = 0, u(1) = u(2) = 0 avec h = 0.5.
(d) u00 = 2u0 − u + t et − t 0 ≤ t ≤ 2 u(0) = 0, u(2) = −4 avec h = 0.5
3. Résoudre les EDPs suivantes par la méthode des différences finies
(a) uxx + u yy = 4, 0 < x < 1, 0<y<2
u(x, 0) = x2 , u(x, 2) = (x − 2)2 , 0≤x≤1
u(0, y) = y2 , u(1, y) = (y − 1)2 , 0≤y≤2
avec h = k = 0.5 ; puis comparer le résultat obtenu avec la solution exacte
u(x, y) = (x − y)2 .
(b) uxx + u yy = 0, 1 < x < 2, 0<y<1
u(x, 0) = 2 ln 2, u(x, 1) = ln(x2 + 1)2 , 1≤x≤2
u(1, y) = ln(y2 + 1), u(2, y) = ln(y1 + 2), 0≤y≤1
avec h = k = 0.25 ; puis comparer le résultat obtenu avec la solution exacte
u(x, y) = ln(x2 + y2 ).
(c) u y − uxx = 0, 0 < x < 2, 0<y<1
u(x, 0) = sin(πx/2), 0≤x≤2

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 9

u(0, y) = u(2, y) = 0, 0≤y

avec N = 4, M = 2 ; puis comparer le résultat obtenu avec la solution exacte
u(x, y) = e−(π /4)y sin(πx/2) au point y = 1.0.
2

(d) u y − uxx = 0, 0 < x < 1, 0<y<1

u(x, 0) = sin(πx), 0≤x≤1
u(0, y) = u(1, y) = 0, 0≤y
u y (x, 0) = 0, 0≤x≤1
avec N = 4, M = 4 ; puis comparer le résultat obtenu avec la solution exacte
u(x, y) = cos(πy) sin(πx) au point y = 0.1.
(e) u y − 1
u
16π2 xx
= 0, 0 < x < 0.5, 0 < y < 0.5
u(x, 0) = 0, 0 ≤ x ≤ 0.5
u(0, y) = u(0.5, y) = 0, 0≤y
u y (x, 0) = sin(4πx), 0 ≤ x ≤ 0.5
avec N = 4, M = 4 ; puis comparer le résultat obtenu avec la solution exacte
u(x, y) = cos(πy) sin(4πx) au point y = 0.5.
(f) uxx + u yy = (x2 + y2 )exy , 0 < x < 2, 0<y<1
u(x, 0) = 1, u(x, 1) = ex , 0≤x≤2
u(0, y) = 1, u(2, y) = e2y , 0≤y≤1
avec h = 0.2, k = 0.1 ; puis comparer le résultat obtenu avec la solution exacte
u(x, y) = exy .
y
(g) uxx + u yy = x
y
+ x, 1 < x < 2, 1<y<2
u(x, 1) = x ln x, u(x, 2) = x ln(4x2 ), 1≤x≤2
u(1, y) = y ln y, u(2, y) = 2 ln(2y), 1≤y≤2
avec h = k = 0.1 ; puis comparer le résultat obtenu avec la solution exacte
u(x, y) = xy ln(xy).

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 Chapitre 2

Méthodes des éléments finis

Nous commençons par la méthode des résidus afin d’introduire d’une manière douce
la méthode des éléments finis. Cette manière de faire permet d’éviter de parler des
espaces fonctionnels (espaces de Sobolev) pour introduire les éléments finis.

2.1 Méthode des résidus

On considère un problème général décrit par une EDO définie sur un domaine Ω avec
les conditions définies sur le bord Γ = ∂Ω. Le schémas de la méthode des résidus est
d’approximer la solution de l’EDO en suivant les étapes suivantes :
1. Supposer une solution d’essai définies avec quelques paramètres. Par exemple,
en dimension 1 on peut supposer que la solution d’essai est de la forme :

û(x) = c0 + c1 x + c2 x2 + c3 x3 + · · · ou
û(x) = a0 + a1 cos x + b1 sin x + a2 cos 2x + b2 sin 2x + · · ·

2. Cette solution doit satisfaire l’EDO sur Ω ainsi que les conditions sur le bord Γ.
En substituant cette solution d’essai dans l’EDO et dans les conditions de bord,
on obtient l’erreur appelée respectivement résidu du domaine et résidu de bord.
3. Déterminer les paramètres inconnues c0 , c1 , c2 , c3 , · · ·, ou a0 , a1 , b1 , a2 , b2 , , · · · de
telle sorte que les résidus soient les plus petits possibles.
L’idéal serait d’annuler les résidus.

Exemple 2.1.1 Soit à résoudre le problème

d2 u
α (x) + q0 = 0, 0 < x < L (2.1.1)
dx2
du
u(0) = 0, (L) = 0. (2.1.2)
dx
où α, q0 et L sont des constantes réelles données avec α , 0 et L > 0.

Ce problème représente la déformation d’une barre de longueur L fixée à une extrémité,

libre à l’autre et soumise à une charge axiale uniforme q0 .
Etape 1 : Choix d’un solution d’essai
Pour résoudre ce problème, on suppose que la solution d’essai est donnée par :

û(x) = c0 + c1 x + c2 x2 . (2.1.3)

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 11

où c0 , c1 et c2 sont des paramètres à déterminer.

Puisque û(0) ' u(0) = 0, on a que c0 = 0. D’autre part, comme û0 (x) = c1 + 2c2 x ' u0 (x)
doit s’annuler en x = L, on a aussi c1 + 2c2 L = 0 =⇒ c1 = 2Lc2 . Ainsi,
 
û(x) = −2Lc2 x + c2 x2 = c2 x2 − 2Lx . (2.1.4)
Etape 2 : Détermination du résidu du domaine
Rd = αû”(x) + q0 = 2αc2 + q0 . (2.1.5)
Etape 3 : Minimisation du résidu du domaine
q0
Rd = 0 ⇐⇒ 2αc2 + q0 = 0 =⇒ c2 = − . (2.1.6)
2α
Donc, la solution approximative est donnée par :
q0  
û(x) = 2Lx − x2 . (2.1.7)
2α
Exemple 2.1.2 La vitesse d’écoulement laminaire d’un fluide visqueux sur une surface plane
inclinée est régie par :
µv00 (x) + ρ g cos θ = 0, x ∈]0, L[, (2.1.8)
où v, la vitesse d’écoulement, µ le coefficient de viscosité, ρ la densité du fluide, g l’accélération
de la pesanteur, θ l’angle entre le plan incliné et la verticale et L la largeur de la planche.
Les conditions de bord sont données par :
v0 (0) = 0 (Pas de contrainte de cisaillement)
v(L) = 0 (Pas de glissement) (2.1.9)
Etape 1 : Choix d’une solution d’essai
Pour trouver une approximation de la distribution de la vitesse v(x) par la méthode des
résidus, on suppose que la solution d’essai est de la forme :
v(x) ' v̂(x) = c0 + c1 x + c2 x2 . (2.1.10)
Il en découle que
v̂0 (x) = c1 + 2c2 x. (2.1.11)
Des conditions de bord on trouve : c1 = 0, c0 = −c2 L2 . Ainsi,
 
v̂(x) = c2 x2 − L2 . (2.1.12)
Etape 2 : Détermination du résidu de domaine
Rd = µ(2c2 ) + ρ g cos θ. (2.1.13)
Etape 3 : Minimisation du résidu du domaine
ρ g cos θ
Rd = 0 ⇐⇒ c2 = − . (2.1.14)
2µ
Ainsi, la vitesse d’écoulement est donnée par
ρ g cos θ  2 
v̂(x) = L − x2 . (2.1.15)
2µ

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 12

Remarque 2.1.3
Les solutions d’essai des exemples précédents sont exactes car on a réussi à avoir Rd = 0 sur
tout le domaine ]0, L[ de l’EDO.

Exemple 2.1.4 Considérons maintenant une poutre de longueur L ayant deux points d’appui
et soumise à une charge transversale uniformément distribuée q0 . Les déformations de cette
poutre est régie par l’EDO

d4 v
β − q0 = 0, x ∈]0, L[. (2.1.16)
dx4
Les conditions de bord sont données par

d2 v
v(0) = 0, (0) = 0
dx2
2
dv
v(L) = 0, (L) = 0. (2.1.17)
dx2
Etape 1 : Choix d’un solution d’essai
On pourrait prendre v̂ comme dans l’exemple précédent, mais pour des raisons de
simplicité prenons, au vu des conditions de bord,

v(x) ' v̂(x) = c1 sin(πx/L). (2.1.18)

On peut remarquer que v̂ ne contient qu’un paramètre c1 et satisfait les conditions de

bord.
Détermination du résidu
d4 v̂
Rd = β 4
− q0 = βc1 (π/L)4 sin(πx/L) − q0 . (2.1.19)
dx
Minimisation du résidu
Contrairement aux 3 exemples précédents, le résidu de domaine Rd est variable sur
]0, L[. Comme il s’agit de déterminer un seul paramètre on peut poser Rd = 0 seule-
ment en seul point pris arbitrairement sur l’intervalle ]0, L[. Cette technique s’appelle
”technique de collocation”.
Pour cela, supposons que Rd = 0 en x = L/4, c’est-à-dire
√
2q0 L4
βc1 (π/L) sin(π/4) − q0 = 0
4
⇐⇒ c1 = . (2.1.20)
βπ4
Par suite,
√
2q0 L4
v̂(x) = sin(πx/L) (2.1.21)
βπ4
On peut répéter le processus en changeant le point où le résidu s’annule et comparer
les résultats. Il est évident que toutes ces solutions ne seront pas exactes.
Comme le problème le problème en question est symétrique par rapport au point
x = L/2, on modifie la solution d’essai comme suit :

v̂(x) = c1 sin(πx/L) + c2 sin(3πx/L). (2.1.22)

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 13

Le résidu du domaine est donné par :

d4 v̂
Rd = β − q0 = βc1 (π/L)4 sin(πx/L) + βc2 (3π/L)4 sin(3πx/L) − q0 . (2.1.23)
dx4
Ayant deux paramètres à déterminer, on peut supposer Rd = 0 en deux points suivants :
x = L/4 et x = L/3. On a
1 1
β(π/L)4 √ c1 + β(3π/L)4 √ c2 − q0 = 0
√2 2
3
β(π/L)4 c1 − q0 = 0.
2
La solution de ce système est

2q0 L4
c1 = √ (2.1.24)
βπ4 3
√ √
(2 − 3) 2q0 l4
c2 = √ . (2.1.25)
β34 π4 3

La solution d’approximation v̂ est alors donnée par

√ √
2q0 L4 (2 − 3) 2q0 l4
v̂(x) = √ sin(πx/L) + √ sin(3πx/L). (2.1.26)
βπ4 3 β34 π4 3

Si on compare les déviations entre la solution exacte et les solutions approchées ci-
dessus déterminées, on constatera que la dernière parait plus bonne que la première.

2.2 Méthode des résidus pondérés

Plutôt que d’essayer d’annuler le résidu du domaine à quelques points sélectionnés au
hasard, on va essayer de minimiser le résidu sur tout l’intervalle ]0, L[. Cette technique
s’appelle ”technique de poids résiduel” qui se formule de la manière suivante
Z L
Wi (x)Rd (x)dx = 0, (2.2.1)
0

où Wi (x) sont appelés fonctions poids et choisies afin d’atteindre cette objectif de min-
imiser le résidu sur tout le domaine ]0, L[ Leur nombre dépend du nombre des
paramètres à déterminer. Mais le choix de ces fonctions reste complètement arbitraire.
La seule exigence est qu’elles soient non nulles et intégrables.
Galerkin (1915) proposa de prendre pour fonctions poids Wi (x) les solution d’essai
elles-mêmes. Ainsi pour l’exemple précédent on prendra les fonctions poids :

W1 (x) = sin(πx/L) et W2 (x) = sin(3πx/L), · · ·

Ainsi, si on choisit v̂(x) = c1 sin(πx/L) alors pour minimiser

Rd = βc1 (π/L)4 sin(πx/L) − q0 ,

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 14

la technique de Galerkin propose

Z L h i
sin(πx/L) β(π/L)4 sin(πx/L) − q0 dx = 0 (2.2.2)
0

c’est-à-dire la technique de Galerkin propose

Z L h i
sin(πx/L) βc1 (π/L)4 sin(πx/L) − q0 dx = 0
0
ou
Z L Z L
c1 β(π/L) 4
sin (πx/L)dx =
2
sin(πx/L)q0 dx
0 0
ou encore
L 2q0 L 4 q0 L4
c1 β(π/L)4 = =⇒ c1 = (2.2.3)
2 π π5 β
Par suite, la solution approchée est
4 q0 L4
v̂(x) = sin(πx/L). (2.2.4)
π5 β
Il est évident que les déviations entre la solution exacte v et cette dernière sont faibles.
Pour résumer, la technique de Galerkin consiste à considérer les fonctions poids comme
les fonction d’essai elles-mêmes. Pour être plus précis, on supposera que la solution
d’essai est donnée par
X
f (x) = Φ(x) + ci Ni (x), (2.2.5)

où Φ(x) et Ni (x) sont totalement connues et présélectionnées de telle sorte que les
conditions de bord sont satisfaites par f (x) et donc, les fonctions poids sont données
par
∂f
Wi (x) = (x) = Ni (x). (2.2.6)
∂ci

Exemple 2.2.1 Considérons le problème de propagation de la chaleur sur un cylindre de rayon

R et de longueur L soumis à une source de chaleur q0 . Cette propagation de chaleur est régie par
l’EDO (en coordonnées cylindriques) :
d2 T 1 dT q0
+ + = 0, (2.2.7)
dr2 r dr k
qui exprime la distribution de la température T en fonction de la localisation radiale r. Les
conditions de bord sont :

T(R) = Tm
dT
q0 πR2 L = (−k)(2πR L) (R), (2.2.8)
dr
où Tm représente la température du mur et la deuxième condition de bord exprime que la chaleur
est produite est égale à la chaleur perdue.

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 15

Etape 1 : Choix de la fonction d’essai

Nous supposons qu’elle est donnée par

T(r) ' T̂(r) = c0 + c1 (r − R) + c2 (r − R)2 . (2.2.9)

La première condition de bord T(R) = Tm implique que c0 = Tm Tandis que la seconde

q0 R
implique que q0 πR2 L = −k(2πR L)c1 . Autrement dit, c1 = − . De sorte que
2k
T̂(r) = Tm − (q0 R/2k)(r − R) + c2 (r − R)2 = Φ(r) + c2 N2 (r),
 
(2.2.10)

avec Φ(r) = Tm − (q0 R/2k)(rR ) et N2 (r) = (r − R)2 . Donc, W(r) = N2 (r)

Etape 2 : Détermination du résidu de domaine
" #
1 q0 R q0
Rd = 2c2 + − + 2c2 (r − R) + . (2.2.11)
r 2k k

Etape 2 : Minimisation du résidu

Nous allons minimiser le résidu sur tout le volume du cylindre en considérant que la
température reste constante sur toute circonférence de rayon r le long du cylindre. En
appliquant la technique de poids résiduel de Galerkin, on a
Z L Z 2π Z R
W(r)[Rd ]r dr dθ dz = 0
0 Z0 R " 0 ! #
1 q0 R q0
⇐⇒ 2πL (r − R) 2c2 + −
2
+ 2c2 (r − R) + r dr = 0
r 2k k
Z 0R h Z R
i q0
⇐⇒ 2c2 r(r − R) + (r − R) dr =
2 3
(r − R)2 (R − 2r) dr
0 2k 0
R4 q0 R4 q0
⇐⇒ − c2 = ⇐⇒ c2 = − .
3 2k 6 4k
Ainsi,
q0 R q0 q0  2 
T̂(r) = Tm − (r − R) − (r − R)2 = Tm + R − r2 . (2.2.12)
2k 4k 4k

2.3 Généralisation de la méthode des résidus pondérés

Pour un champ des scalaires inconnu u, nous supposons une solution d’approximation
de la forme
n
X
u ' û = Φ + ci Ni (2.3.1)
i=1

où les ci sont des paramètres à déterminer par la méthode des résidus pondérés, Φ et Ni
sont des fonction choisies de sorte que û satisfait les conditions de bord. En notant Rd
le résidu du domaine Ω, alors les équations résiduelles pondérées sont données par :
Z
Wi Rd dΩ = 0, i = 1, · · · , n, (2.3.2)
Ω

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 16

où Wi = Ni . On notera que ces équations sont satisfaites par la solution exacte du
problème car son résidu est trivialement nul partout.
Ainsi, l’idée générale de la méthode des résidus pondérés est que si on est capable
de satisfaire l’équation résiduelle pondérée pour un plus grand nombre des fonctions
indépendantes, alors on a la chance que la solution d’essai considérée soit plus proche
de la solution exacte. Plus généralement, si on considère que la fonction d’essai possède
la forme d’une série (entière ou de Fourier), on peut espérer que le résultat obtenu sera
meilleur si on a plus des termes dans la série. Donc, on espère des bonnes propriétés
de convergence.

Exemple 2.3.1 Si on considère une plate forme rectangulaire horizontale de longueur a et de

largeur b simplement supportée sur sa bordure et soumise à un charge verticale q0 uniformément
répartie sur sa surface, alors sa déformation est régie par l’EDP suivante :

∂4 w ∂4 w ∂4 w
!
β + 2 2 2 + 4 − q0 = 0, 0 < x < a, 0 < y < b, (2.3.3)
∂x4 ∂x ∂y ∂y

avec comme conditions de bord :

∂2 w ∂2 w
w(0, y) = 0, w(a, y) = 0, (0, y) = 0, (a, y) = 0, 0 < y < b,
∂x2 ∂x2
(2.3.4)
∂w2
∂w
2
w(x, 0) = 0, w(x, b) = 0, (0, y) = 0, (a, y) = 0, 0 < x < a.
∂x2 ∂y2

Etape 1 : Choix de la solution d’essai

Le conditions de bord poussent à choisir la fonction d’essai

w(x, y) ' ŵ(x, y) = c1 sin(πx/a) sin(πy/b) (2.3.5)

qui, comme on peut le voir, satisfait les conditions de bord du problème.

Etape 2 : Détermination du résidu de domaine
h i2
Rd = β (π/a)2 + (π/b)2 c1 sin(πx/a) sin(πy/b) − q0 (2.3.6)

car
∂4 ŵ
= (π/a)4 c1 sin(πx/a) sin(πy/b)
∂x4
∂4 ŵ
= (π/a)2 (π/b)2 c1 sin(πx/a) sin(πy/b)
∂x2 ∂y2
∂4 ŵ
= (π/b)4 c1 sin(πx/a) sin(πy/b).
∂y 4

Etape 3 : Minimisation du résidu

La technique de Galerkin implique que W(x, y) = sin(πx/a) sin(πy/b) et par suite,
Z bZ a   2 
sin(πx/a) sin(πy/b) β (π/a) + (π/b) c1 sin(πx/a) sin(πy/b) − q0 dxdy = 0,
2 2

0 0

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 17

ou encore
" 2  2 #2 Z b Z a
2 πy
Z bZ a  πy 
π π 2 πx πx
     
βc1 + sin sin dxdy = q0 sin sin dxdy.
a b 0 0 a b 0 0 a b

En intégrant, on obtient l’équation algébrique linéaire suivante :

" 2  2 #2
π π ab 4ab
βc1 + − q0 2 = 0 (2.3.7)
a b 4 π

dont la solution est :

!2
16q0 16q0 a2 b2
c1 = = (2.3.8)
βπ2 ((π/a)2 + (π/b)2 )2 βπ6 a2 + b2

Enfin, la solution d’approximation s’écrit :

!2
16q0 a2 b2
ŵ(x, y) = sin(πx/a) sin(πy/b) (2.3.9)
βπ6 a2 + b2

Remarque 2.3.2 Bien que l’augmentation des termes dans la solution d’essai améliore cette
dernière, elle possède néanmoins un inconvénient sur le plan algorithmique (les calculs devi-
ennent de plus en plus lourds). En lieu et place d’augmenter les termes nous allons, dans un
premier temps, utiliser moins des termes dans la solution d’essai en affaiblissant l’équation
résiduelle pondérée. C’est l’objet de la section suivante

2.4 Formulation variationnelle résiduelle pondérée

En considérant la formulation générale de la méthode des résidus pondérés
Z
W Rd dΩ = 0, (2.4.1)
Ω

il peut être possible, en utilisant la formule d’intégration par parties, de réduire le

degré de dérivabilité de la fonction solution. Ce qui aura pour conséquence de prendre
comme solution d’essai des polynômes de degré faible, par exemple.

Exemple 2.4.1 Considérons encore une fois le problème

d2 u
α (x) + ax = 0, 0 < x < L,
dx2
(2.4.2)
du
u(0) = 0, α (L) = 0.
dx
Prenons û comme solution d’essai. Par conséquent, le résidu du domaine est donné par

d2 û
Rd = α 2 + ax. (2.4.3)
dx

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 18

Avec une fonction poids W on obtient l’équation résiduelle pondérée

Z L !
d2 û
W(x) α 2 + ax dx = 0,
0 dx
(2.4.4)
du
u(0) = 0, α (L) = 0.
dx
Notons que l’intégrale peut se décomposer en
L L
d2 û
Z Z
W(x)α 2 dx + W(x)ax dx = 0, (2.4.5)
0 dx 0

ou encore
Z L ! Z L
dû
W(x)d α + W(x)ax dx = 0. (2.4.6)
0 dx 0

En intégrant par partie le premier terme on obtient

" #L Z L Z L
dû dû dW
W(x)α − α dx + W(x)ax dx = 0. (2.4.7)
dx 0 0 dx dx 0

dû
En posant P = α , appelée aussi force axiale, l’équation précédente s’écrit
dx
Z L Z L
dû dW
W(L)P(L) − W(0)P(0) − α dx + W(x)ax dx = 0. (2.4.8)
0 dx dx 0

C’est la formulation variationnelle ou faible de la méthode des résidus pondérés. On

remarquera dans cette formulation que la fonction û doit simplement être un fois
dérivable au lieu d’être deux fois dérivable.
Puisque û satisfait les conditions de bord, P(L) = 0, û(0) = 0 et on peut exiger que W(x)
satisfasse la condition W(0) = 0. Dès lors, la formulation variationnelle se réduit à
Z L Z L
dû dW
α dx = W(x)ax dx,
0 dx dx 0
(2.4.9)
û(0) = 0, W(0) = 0.

Sa solution est dite solution faible ou affaiblie tandis que les autres solutions sont dites
solution fortes.
Pour illustrer une solution faible notons toujours par û la solution d’essai du problème
variationnelle, mais un polynôme de degré 2 :

û(x) = +c1 x + c2 x2 (2.4.10)

qui satisfait û(0) = 0. Il en découle que

dû
= c1 + 2c2 x. (2.4.11)
dx

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 19

Dès lors, W1 (x) = x et W2 (x) = x2 satisfont W1 (0) = 0 = W2 (0) et par suite

dW1 dW2
= 1, = 2x. (2.4.12)
dx dx
Ainsi, les équations variationnelles (faibles) sont respectivement données par
Z L Z L
α(c1 + 2c2 x)(1)dx = (x)(ax)dx (2.4.13)
0 0
Z L Z L
α(c1 + 2c2 x)(2x)dx = (x2 )(ax)dx. (2.4.14)
0 0

Après intégration on obtient

  L3
α c1 L + c2 L2 = a
3
  L4
α c1 L + c 43 L
2 3
= a
4
ce qui correspond au système
aL2
c1 + Lc2 = (2.4.15)
3α
4L aL2
c1 + = (2.4.16)
3 4α
dont la solution est
7aL2 aL
c1 = et c2 = − . (2.4.17)
12α 4α
Ainsi, la solution faible est donnée par
aL  
û(x) = 7Lx − 3x2 . (2.4.18)
12α
On vérifiera que que sa courbe proche de la courbe de la solution forte
a  2 
u(x) = 3L x − x3 . (2.4.19)
6α
Remarque 2.4.2 Utilisation Fonctions continues et définies par morceaux
Jusqu’à présent les les fonctions d’essai choisies avaient une seule composante sur tout
le domaine. L’idée d’utiliser les fonctions continues mais définies par morceaux vient
du fait qu’il est plus facile d’approcher une courbe par un segment de droite sur un
sous-intervalle de longueur petite.
Nous discrétisons l’intervalle [0, L] en N sous-intervalles Ki = [xi−1 , xi ], 1 ≤ i ≤ N où
x0 = 0 et xN = L comme indiqué dans la figure ci-dessus. Dans chaque intervalle Ki la
fonction u(x) est approchée par le segment de droite d’équation
u(xi ) − u(xi−1 ) x − xi−1 x − xi−1
   
û(x) = u(xi−1 ) + (x − xi−1 ) = u(xi−1 ) 1 − + u(xi )
xi − xi−1 xi − xi−1 xi − xi−1
En posant :
x − xi−1 x − xi−1
ui−1 = u(xi−1 ), ui = u(xi ), λi1 (x) = 1 − et λi2 (x) = ,
xi − xi−1 xi − xi−1

û(x) = ui−1 λi1 (x) + ui λi2 (x), xi−1 ≤ x ≤ xi . (2.4.20)

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 20

1 λ11 λ12 λ21 λ22 λ31 λ32

x
x0 = 0 x1 x2 x3 x4 L = xN
Figure : Maillage du domaine et bases locales

Pour faciliter les calculs on fait une discrétisation régulière avec un pas h = NL . Avec
cette discrétisation du domaine [0, L], les fonctions λi1 et λi2 sont appelées fonctions de
la base locale de l’élément Ki = [xi−1 , xi ]. En posant s = x − xi−1 , ces deux fonctions ainsi
que leurs dérivées sont localement sur définies par :

s dλi1 (s) 1
λi1 (s) = 1− , = − , 0 ≤ s ≤ h, (2.4.21)
h ds h
s dλi2 (s) 1
λi2 (s) = , = − , 0 ≤ s ≤ h. (2.4.22)
h ds h
Ces fonctions de bases locales permettent de définir les fonctions d’interpolation
définies sur tout le domaine [0, L] appelées aussi fonctions chapeaux par :
 x − xi−1
λ (x), xi−1 ≤ x ≤ xi
 i 
 , xi−1 ≤ x ≤ xi
 2i+1

 
 h
φi (x) = 

λ1 (x), xi ≤ x ≤ xi+1 càd φi (x) = 
 x − x i
, xi ≤ x ≤ xi+1 (2.4.23)

 1−
h

 
 0, ailleurs 

 0, ailleurs

pour 1 ≤ i ≤ N − 1 et dont les dérivées sont données par :

1

, xi−1 ≤ x ≤ xi



 h

φ0i (x) =  1

(2.4.24)

 − , xi ≤ x ≤ xi+1
h



 0, ailleurs.


Il est à noter que

 x − x0
λ11 (x), x0 ≤ x ≤ x1 , x0 ≤ x ≤ x1
(
 1−

φ0 (x) = càd φ0 (x) = 

h (2.4.25)
0, ailleurs  0,

ailleurs

 x−x
N−1
λN , xN−1 ≤ x ≤ xN
(
(x), xN−1 ≤ x ≤ xN 
φN (x) = càd φN (x) = 

2
h (2.4.26)

0, ailleurs 
 0, ailleurs

1 1
 
 − , x0 ≤ x ≤ x1  , xN−1 ≤ x ≤ xN
 
φ0 (x) =  φN (x) = 
0
 0

 h  h (2.4.27)
 0, ailleurs  0, ailleurs

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 21

v Ki−1 Ki Ki+1

1 φi−1 φi φi+1

x
x0 = 0 xi−1 xi xi+1 L = xN
Figure : Fonctions de base globale

Chaque fonction chapeaux φi donne la contribution ui dans la valeur de û sur [0, L].
Aussi φi (xi ) = 1 et φi (x j ) = 0 ; φi croit sur Ki−1 et décroit sur Ki , 1 ≤ i ≤ N tandis que φ0
est décroit et φN+1 croit.
Finalement, l’interpolation de la fonction u sur le domaine [0, L] est donnée par :
N
X
û(x) = ui φi (x). (2.4.28)
i=0

Exemple 2.4.3 Soit à résoudre le problème :

d2 u
α (x) + q0 = 0, 0 < x < L,
dx2
(2.4.29)
du
u(0) = u0 , α (L) = PL .
dx
En effet, la formulation résiduelle forte pondérée est :
L !
d2 û
Z
W(x) α 2 + q dx = 0,
0 dx
(2.4.30)
du
u(0) = u0 , α (L) = PL .
dx
En intégrant par partie le premier terme on obtient
" #L Z L Z L
dû dû dW
W(x)α − α dx + W(x)q0 dx = 0. (2.4.31)
dx 0 0 dx dx 0

dû
En posant α = P(x), cette équation peut encore se mettre sous la forme :
dx
Z L Z L
dû dW h iL
α dx = W(x)q0 dx + W(x)P(x) . (2.4.32)
0 dx dx 0
0

Tenant compte de la discrétisation du domaine en N éléments finis, l’équation devient :

N Z h N Z h N h
X dû dW X X ih
α ds = W(s)q ds + W(x)P(s) , (2.4.33)
i=1 0 ds ds i=1 0 i=1
0

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 22

où on a utilisé la règle de Chasles et les cordonnées locales de chaque maille Ki . Sachant
que dans cette maille û est interpolée par
" #
h i ui−1
û(s) = ui−1 λ1 (s) + ui λ2 (s) = λ1 (s) λ2 (s)
i i i i
(2.4.34)
ui

Par conséquent,
" #
dû(s) 1 1 1h i ui−1
= − ui−1 + ui = −1 1 . (2.4.35)
ds h h h ui

On voudrait déterminer la contribution de l’élément Ki en utilisant les fonctions poids

s dW1 1
W1 (s) = λi1 (s) = 1 − =−
h ds h
s dW1 1
W2 (s) = λi2 (s) = = .
h ds h
Le premier membre de l’élément Ki donne :
h h
−ui−1 + ui 1 αh
Z Z " #
dû dW1 ui−1
  i
α ds = α − ds = 1 −1 (2.4.36)
0 ds ds 0 h h h ui

h h
−ui−1 + ui 1 αh
Z Z " #
dû dW2 ui−1
  i
α ds = α ds = −1 1 (2.4.37)
0 ds ds 0 h h h ui

Le second membre donne :

Z h Z h
s
 q0 h 
W1 (s)q0 ds = 1 − q0 ds = (2.4.38)
0 0 h 2
Z h Z h 
s q0 h
W2 (s)q0 ds = q0 ds = (2.4.39)
0 0 h 2

" #h h
dû s
 
W1 (s)α = 1−P(s) = −P0 (2.4.40)
ds 0
h 0
" #h h
dû s
 
W2 (s)α = P(s) = Ph . (2.4.41)
ds 0
h 0

En combinant les résultats des calculs obtenus dans le premier et le deuxième membres
de l’élément Ki , on obtient le système caractéristique de l’élément :
#  q0 h  "
α
" #" #
1 −1 ui−1 −P 0
=  q0 h  + .
 2 
(2.4.42)
h −1 1 ui 2
P h

Dans le premier membre, la quantité αh est connu le nom de coefficient de rigidité

axiale de la route, la matrice s’appelle élément matrice de rigidité tandis que le second
membre s’appelle élément du vecteur force.
Il reste maintenant à assembler ces éléments pour construire le système global.

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 23

Pour illustrer notre propos, on suppose que notre domaine est subdivisé en quatre
sous-domaine. Dans ce cas les équations caractéristiques des quatre mailles sont re-
spectivement données par :
#  q0 h  "
α
" #" #
1 −1 u0 −P 0
=  q0 h  + ,
 2 
(2.4.43)
h −1 1 u1 2
P 1

q0 h
  "
α
" #" # #
1 −1 u1 −P 1
=   + ,
 2

q0 h (2.4.44)
h −1 1 u2 2
 P2

q0 h
  "
α
" #" # #
1 −1 u2 −P 2
=   + ,
 2

q0 h (2.4.45)
h −1 1 u3 2
 P3

q0 h
  "
α
" #" # #
1 −1 u3 −P 3
=   + .
 2

q0 h (2.4.46)
h −1 1 u4 2
 P4
Notons que ces systèmes de huit équations à huit inconnues : les P0 , P1 , P2 , P3 ainsi
que les u1 , u2 , u3 , u4 car u0 et P4 = PL sont données. On peut éliminer les inconnues
P1 , P2 , P3 en additionnant membres à membres la dernière et la première équation
des sous-systèmes consécutifs pour obtenir le système réduit à cinq équations à cinq
inconnues :
  qh  
0   u0   20   −P0 
  
 1 −1 0 0
      
 −1 2 −1 0 0   u1   q0 h   0 
α  
0   u2  =  q0 h  +  0  ,
   
 0 −1 2 −1
 
(2.4.47)
h      
 0 0 −1 2 −1   u3   q0 h   0 
   
   q0 h  
0 0 0 −1 1 u4 PL
 
2

Les quatre dernières équations donnent le système

  αu0 /h 
   
0   u1   q0 h
  
 2 −1 0
α  −1
   
0   u2   q0 h
  
2 −1   0 
=  +  , (2.4.48)
h  0 −1 2 −1   u3   q0 h 0 
   
  q0 h  
0 0 −1 1 u4 PL
   
2

La résolution de ce dernier permet d’obtenir les coefficients u1 , u2 , u3 , u4 et construire

ainsi la fonction d’interpolation de u :
û(x) = u0 φ0 (x) + u1 φ1 (x) + u2 φ2 (x) + u3 φ3 (x) + u4 φ4 (x), (2.4.49)
où
1 − xh , 0 ≤ x ≤ h , 3h ≤ x ≤ 4h
( ( x−3h
φ0 (x) = φ4 (x) = h
0 ailleurs 0 ailleurs
x
, ,
  x−h
h
0≤x≤h h ≤ x ≤ 2h
 h x−2h

 

φ1 (x) =  1 − h , h ≤ x ≤ 2h φ2 (x) =  1 − h , 2h ≤ x ≤ 3h

 x−h


 

 0 ailleurs  0 ailleurs
x−2h
,


 h
2h ≤ x ≤ 3h
φ3 (x) =  1 − h , 3h ≤ x ≤ 4h

 x−3h


 0 ailleurs

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 24

Remarque 2.4.4 Equations nodales

On obtiendrai le même système si on utilisait les fonctions de base globale comme

fonctions poids. En substituant
N
X dφi
dû
(x) = ui (x). (2.4.50)
dx i=0
dx

dans la formulation variationnelle faible (2.4.32), en utilisant W j (x) = φ j (x), on obtient

L N L
dφi dφ j
Z X Z h iL
α ui dx = φ j (x)q0 dx + φ j (x)P(x) .
0 i=0
dx dx 0
0

N L L
dφi dφ j
X Z Z h iL
ui α dx = φ j (x)q0 dx + φ j (x)P(x) . (2.4.51)
i=0 0 dx dx 0
0

C’est l’équation du nœud j = 0, 1, · · · , N. Comme Φ j (x) = 0 en dehors de l’intervalle

[x j−1 , x j ], les intégrales vont se réduire à cet intervalle :
j+1 x j+1 x j+1
dφi dφ j
X Z Z h iL
ui α dx = φ j (x)q0 dx + φ j (x)P(x) , (2.4.52)
i=j−1 x j−1 dx dx x j−1
0

2.5 Eléments finis de Lagrange

Dans cette section nous présentons la plus simple et classique d’éléments finis. Avant
d’aller plus loin nous généralisons quelques idées introduites dans la section précédente.

Définition 2.5.1 Unisolvance

Soient Σ = {a1 , · · · , aN } un ensemble de N points distincts de RN et P un espace vectoriel de
dimension fini de fonctionnelles réelles définies sur RN . On dit que Σ est P-unisolvant si et
seulement si pour tous réels α1 , · · · , αN , il existe une unique fonctionnelle p de P tel que

p(ai ) = αi , i = 1, · · · , N.

Ceci revient encore à dire que l’application L : P −→ RN qui à une fonction p associe
L(p) = p(a1 ), · · · , p(aN ) est bijective.

En pratique on montre que Σ est P-unisolvant en vérifiant dim P = N = CardΣ ensuite

en prouvant l’injectivité ou la surjectivité de l’application L.
L’injectivité de L se démontre en établissant que la seule fonction de P s’annulant sur
tous les points de Σ est la fonction nulle.
La surjectivité de L se démontre en exhibant une famille p1 , · · · , pN d’éléments de P
tels que pi (a j ) = δi j , c’est-à-dire un antécédent pour L d’une base de RN . En effet,
N
X
étant donné des réels α1 , · · · , αN , la fonctionnelle p = αi pi est l’unique satisfaisant
i=1
p(ai ) = αi , i = 1, · · · , N.

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 25

Définitions 2.5.2
Un élément fini de Lagrange est un triplet (K, Σ, P) satisfaisant les 3 conditions suivantes :
1° K est un partie compacte, connexe et d’intérieur non vide de Rn (n = 1, 2, 3) ;
2° Σ = {a1 , · · · , aN } un ensemble de N points distincts K et
3° P un espace vectoriel de dimension fini de fonctionnelles réelles définies sur K et tel que
Σ est P-unisolvant (dim P = N = CardΣ)
On appelle fonctions de bases locales de l’élément (K, Σ, P), les N fonctions p1 , · · · , pN
telles que
pi (a j ) = δi j , 1 ≤ i, j ≤ N.
On appelle opérateur de P-interpolation sur Σ l’opérateur πK qui, à toute fonctionnelle v
définie sur K, associe la fonctionnelle πK v de P définie par
N
X
πK v = v(ai )pi .
i=1

πK v est donc l’unique élément de P qui prend les mêmes valeur que v sur les points de Σ.

Avant de donner quelques exemples d’éléments finis de Lagrange rappelons que si on

note Pk l’espace de polyn
n ômes deo degré inférieur ou égal à k, alors sur :
– R est Pk = Vect 1, x, · · · , xk et dim Pk = k + 1 ;
n o (k + 1)(k + 2)
– R2 est Pk = Vect xi y j , 0 ≤ i + j ≤ k et dim Pk = ;
2
n o (k + 1)(k + 2)(k + 3)
– R3 est Pk = Vect xi y j zl , 0 ≤ i + j + l ≤ k et dim Pk = .
6
Si on note Qk l’espace de polynômes de degré inférieur ou égal à k par rapport à une
seule des variables, alors sur :
– R est Qk = Pk ; n o
– R2 est Qk = Vect xi y j , 0 ≤ i, j ≤ k et dim Qk = (k + 1)2 ;
n o
– R3 est Qk = Vect xi y j zl , 0 ≤ i, j, l ≤ k et dim Qk = (k + 1)3 .

Exemples 2.5.3 En dimension 1-D

– Elément P1 : K = [a, b], Σ = {a, b} et P = P1 .

– Elément P2 : K = [a, b], Σ = {a, a+b
2
, b} et P = P2 .
– Elément PN : K = [a, b], Σ = {xi = a + i b−a
N
, i = 0, · · · , N} et P = PN .
On note que les polynômes de Lagrange :
N
Y x − xj
pi (x) = , i = 0, · · · , N
j=0
xi − x j
j,i

constituent aussi une base  de PN , tel que pi (x j ) = δi j . Elles sont

 donc les fonctions de
base locale de l’élément : [a, b], {xi = a + i N , 1 ≤ i ≤ N}, PN . . En particulier,
b−a

– p1 (x) = x−b
a−b
, p2 (x) = x−a
b−a
sont les fonctions de base locale de l’élément ([a, b], {a, b}, P1 )
(x−b)(2x−a−b) −4(x−a)(x−b) (x−a)(2x−a−b)
– p1 (x) = (b−a)2
, p2 (x) = (b−a)2 , p3 (x) = (b−a)2
sont les fonctions de base
 n o 
locale de l’élément [a, b], a, 2 , b , P2 .
a+b

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 26

Exemples 2.5.4 En dimension 2-D

y
a3 (x3 , y3 )
P(λ1 , λ2 , λ3 )

a2 (x2 , y2 )

a1 (x1 , y1 )

x
Fig. Cordonnées barycentriques :
P est le barycentre du triangle

Soit K un triangle de sommets a1 , a2 , a3 .

On appelle cordonnées barycentriques de K les fonctions affines λ1 , λ2 , λ3 de K dans R
définies par

λ j (ai ) = δ ji , 1 ≤ i, j ≤ 3. (2.5.1)

On peut remarquer que la somme λ1 + λ2 + λ3 est une fonction affine qui vaut 1 sur
chacun des 3 sommets. C’est donc la fonction constante égale à 1.
Si on pose (xi , yi ) les coordonnées du sommet ai et λ j (x, y) = α j x + β j y + γ j , alors le
système linéaire (2.5.1) s’écrit :

α j xi + β j yi + γ j = 0



α jx j + β j y j + γ j = 1

(2.5.2)



 α x +β y +γ =0

j k j k j

où {i, j, k} est une permutation de {1, 2, 3}. La solution du système ci-dessus est :
yk − yi xi − xk xk yi − xi yk
αj = , βj = , γj = ,
∆ ∆ ∆
où ∆ = (xi − x j )(y j − yk ) − (x j − xk )(yi − y j ) = 2 l’aire du triangle des sommets ai , a j , ak .
Ainsi,
yi − yk xk − xi xi yk − xk yi
λ j (x, y) = x+ y+ (2.5.3)
∆ ∆ ∆
– Elément P1 : K le triangle de sommets a1 , a2 , a3
Σ = {a1 , a2 , a3 } et P = P1 = Vect {λ1 , λ2 , λ3 }.
Les fonctions de bases locales les coordonnées barycentriques : p j = λ j , (1 ≤ j ≤ 3)
définies-ci-dessus car p j (ai ) = λ j (ai ) = δi j .
– Elément P2 : K le triangle de sommets a1 , a2 , a3
a +a
Σ = {a1 , a2 , a3 , na12 , a13 , a23 } , où ai j = i 2 j o
P = P2 = Vect λ1 λ2 , λ1 λ3 , λ2 , λ3 , λ21 , λ22 , λ23 .

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 27

Les fonctions de bases locales sont pi = λi (2λi − 1), (1 ≤ i ≤ 3) et pi j = 4λi λ j ,

(1 ≤ i < j ≤ 3). En effet, 
pi (a j ) = λi (a j ) 2λi (a j ) − 1 = δi j (2δi j − 1) = δi j .
h i h i
Comme λi (ak j ) = 12 λi (ak ) + λi (a j ) = 12 δik + δi j ,
    h  i
pi (ak j ) = λi (ak j ) 2λi (ak j ) − 1 = 12 δik + δi j δik + δi j − 1 = 0 pour 1 ≤ k < j ≤ 3,
1 si k = i < j = l
(
pi j (akl ) = 4λi (akl )λ j (akl ) = [δik + δil ] [δ jk + δ jl ] =
0 sinon,
pi j (ak ) = 4λi (ak )λ j (ak ) = 4δik δ jk = 0 pour 1 ≤ i < j ≤ 3.

Exemples 2.5.5 En dimension 3-D

– Elément tétraédrique P1 : K = tétraèdre le triangle de sommets a1 , a2 , a3 , a4
Σ = {a1 , a2 , a3 , a4 } et P = P1 = Vect {λ1 , λ2 , λ3 , λ4 }.
où les λ j sont les fonctions coordonnées barycentriques du tétraèdre K définies par
λ j (ai ) = δi j pour 1 ≤ i, j ≤ 4.
Les fonctions de base locale sont les p j = λ j .
– Elément tétraédrique P2 : K = tétraèdre le triangle de sommets a1 , a2 , a3 , a4
a +a
Σ = {a1 , a2 , a3 , a4 , a12 , a13 , a14 , a23 , a24 , a34 }, avec ai j = i 2 j pour 1 ≤ i < j ≤ 4
P = P2  
Les fonctions de base locale sont : p j = λ j 2λ j − 1 pour 1 ≤ i, j ≤ 4. et pi j = 4λi λ j pour
1 ≤ i < j ≤ 4.

Définition 2.5.6  
Deux élément finis K̂, Σ̂, P̂ et (K, Σ, P) sont dits affine-équivalents si et seulement si CardΣ =
N = CardΣ̂ et il existe une fonction affine inversible F (F(x̂) = Bx̂ + b) telle que
– K = F(K̂),
– ai = F( n âi ) pour i = 1,o · · · , N,
– P = p̂ ◦ F : p̂ ∈ P̂ . Autrement dit les fonctions de base locale de K et celles de K̂ sont
liées par : pi = p̂i ◦ F pour i = 1, · · · , N.

Notons en passant que si K ⊂ Rn , alors F est inversible si B est une matrice inversible
régulière de dimension n et b un vecteur de Rn .

Définition 2.5.7 Famille affine d’éléments finis

On appelle famille affine
 d’éléments
 finis une famille d’éléments finis tous affine-équivalents à
un même élément fini K̂, Σ̂, 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.
Voici quelques éléments de référence :
– En 1-D : l’intervalle [0, 1] ;
– En 2 − D : le triangle unité de sommets (0, 0), (1, 0) et (0, 1) ;
– En 3 − D : le tétraèdre unité de sommets (0, 0, 0), (1, 0, 0), (0, 1, 0) et (0, 0, 1).

Remarque 2.5.8 Lien entre le problème global et les éléments locaux

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 28

Supposons qu’on voudrait résoudre une EDO ou une EDP sur un domaine Ω et soit H
l’espace (de Hilbert) dans lequel on cherche une solution de la formulation variation-
nelle du problème. On réalise un maillage de Ω par une famille affine de Ne éléments
finis (Ki , Σi , Pi )1≤i≤Ne . Il est à noter qu’un nœud sera en général commun à plusieurs
éléments adjacents. Le nombre total de nœuds Nh est donc inférieur à Ne × CardΣi et on
a alors dim Hh = Nh . Notons a1 , · · · , aNh les nœuds du maillage. Le problème se ramène
à la détermination des valeurs de uh aux points ai : ce sont les degrés de liberté (ddl) du
problème approché.
On construit alors une base de Hh en associant à chaque ddl ai une fonction de base
globale φi satisfaisant :

φi |K j ∈ P j , j = 1, · · · , Ne et φi (a j ) = δi j , 1 ≤ i, j ≤ Nh (2.5.4)

L’espace d’approximation Hh est ainsi :

n o
Hh = Vect φ1 , · · · , φNh . (2.5.5)

Il est facile alors de remarquer qu’une telle fonction φi est nulle partout, sauf sur les
éléments dont ai est nœud. En effet, si ai n’appartient pas à un élément K, φi est nulle
sur tous les nœuds de K, et donc sur K tout entier par unisolvance.
De plus, sur un élément K dont ai est un nœud, φi vaut 1 sur ai et 0 sur les autres nœuds
de K. Donc φi |K est une fonction de base locale de K. On voit donc que la fonction
de base globale φi est construite comme réunion des fonctions de base locales sur les
éléments du maillage dont ai est un nœud.
C’est à ce niveau que se situe le lien entre les définitions locales introduites ci-haut 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ée le domaine avec une famille d’éléments finis affine-équivalents).
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
en dimension 1, ou encore à un sommet ou une arête en dimension 2 (ou à un sommet,
une arête ou une face en dimension 3).

Exemple 2.5.9 Soit à résoudre le PVB

−∆u(x, y) = f (x, y), (x, y) ∈ Ω,

(
(2.5.6)
u(x, y) = 0, (x, y) ∈ ∂Ω

où f est fonction carrée intégrable sur Ω =]0, 1[×]0, 1[.

Si on note Hh l’espace d’approximation de la solution , alors ce problème a pour

formulation résiduelle :

Trouver û ∈ Hh tel que

" "
∂2 û ∂2 û
!
− W + dxdy = W. f dxdy ∀W ∈ Hh .
Ω ∂x2 ∂y2 Ω

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 29

En intégrant par parties et en utilisant la condition de bord on obtient la formulation

variationnelle (ou formulation faible) suivante :

Trouver û ∈ H01 (Ω) tel que

" "
∂û ∂W ∂û ∂v
!
+ dxdy = W. f dxdy ∀W ∈ Hh .
Ω ∂x ∂x ∂y ∂y Ω

Nh
X
Sachant que û(x, y) = Ui φi (x, y) et utilisant la technique de Galerkin (prendre W =
i=1
φ j 1 ≤ j ≤ Nh ), le problème se réduit à :

Trouver U = U1 , · · · , UNh ∈ RNh tel que

Nh " "
∂φi ∂φ j ∂φi ∂φ j
X !
Ui + dxdy = f.φ j dxdy ∀j = 1, · · · , Nh .
i=1 Ω ∂x ∂x ∂y ∂y Ω

√
Nous discrétisons Ω en triangles de même diamètre h = 2/N pour obtenir un maillage
triangulaire (ou triangulation) conforme comme indiqué sur la figure ci-dessous :

yN = 1

ij
.. T2
ij ij
. T3 T1
..
. ij
T4
ij
T6
ij
T5
y1

y0 = 0 Triangles entourant un
x0 = 0 x1 .... .... xN = 1 x nœud intérieur (xi , y j )
Fig. 2.1 : Triangulation du domaine
 
Ce qui nous permet d’avoir une famille affine d’éléments finis Ki , Σ j , Pi où
1≤i≤N2 −2
– chaque Ki est un triangle rectangle dont les dimensions des deux côtés de l’angle
droit sont égaux à 1/N ;
– chaque Σi est l’ensemble de trois sommets de Ki ;
– Pi est l’espace des fonctions engendré par les coordonnées barycentriques du
triangle Ki .
Notons que le ddl est le nombre des nœuds intérieurs (N − 1)2 et chaque nœud intérieur
est entouré par six triangles. En notant φi j la fonction de base globale associée au nœud

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 30

intérieur des coordonnées (xi , x j ), on peut montrer que

 y−y ij


 1 − x−x
h
i
− h j , (x, y) ∈ T1
y−y ij
1 − h j,

(x, y) ∈ T2




ij
1 + x−x ,

i
(x, y) ∈ T3


h


y−y
φi j (x, y) =  ij

1 + x−xi
+ h j , (x, y) ∈ T4

 h
y−y ij
1 + h j,

(x, y) ∈ T5




ij
1 − x−x ,

i
(x, y) ∈ T6


h



 0,

ailleurs
ij ij ij ij ij ij
où T1 , T2 , T3 , T4 , T5 , T6 sont les triangles entourant le nœud (xi , y j ). Aussi,
ij ij
 


 − 1h , (x, y) ∈ T1 

 − h1 , (x, y) ∈ T1
ij ij
− h , (x, y) ∈ T2
 
 1
0, (x, y) ∈ T2

 

 

ij ij
,

 1 
(x, y) ∈ T3 0, (x, y) ∈ T3

∂φi j ∂φ
 
 h1

 

i j
(x, y) =  ij
(x, y) =  ij
 
, , (x, y) ∈ T4
 1
(x, y) ∈ T4 et
∂x  h
ij ∂y  h
ij
, (x, y) ∈ T5
  1
0, (x, y) ∈ T5

 

h

 

ij ij
− h1 ,
 
(x, y) ∈ T6 0, (x, y) ∈ T6

 


 

 
 0,

ailleurs  0, ailleurs


Ainsi, la fonction d’approximation du problème est interpolée comme suit :

N−1 X
X N−1
Uh (x, y) = Ui, j φi j (x, y). (2.5.7)
i=1 j=1

Le problème d’approximation dévient alors :

Trouver Ui, j ∈ R, 1 ≤, i, j ≤ N − 1 tel que

" "
N−1 X
N−1
∂φi j ∂φkl ∂φi j ∂φkl
X !
Ui, j + dxdy = f.φkl dxdy 1 ≤ k, l, ≤ N − 1.
i=1 j=1 Ω ∂x ∂x ∂y ∂y Ω

Nous allons déterminer l’équation de chaque nœud (xk , yl ), 1 ≤ k, l ≤ N − 1. Chacune

de ces équations est appelée équation nodale.
Comme la fonction de base globale φk,l s’annule en dehors des triangles entourant le
nœud (xk , yl ) cette double somme se réduit aux indices i = k, j = l et aux indices i
et j des nœuds des triangles entourant le nœud (xk , yl ) comme indiqué sur la figure
ci-dessous :

(xk−1 , yl+1 ) (xk , yl+1 )

T2kl
T3kl T1kl
(xk−1 , yl ) (xk+1 , yl )
T4kl T6kl
T5kl
(xk , yl−1 ) (xk+1 , yl−1 )
Figure : Nœuds des triangles entourant le nœud (xk , yl )

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 31

Autrement dit,
"
N−1 X
N−1
∂φi j ∂φkl ∂φi j ∂φkl
X !
Ui,j + dxdy = Uk,l · Ik,l + Uk−1,l · Ik−1,l
i=1 j=1 Ω ∂x ∂x ∂y ∂y
+Uk+1,l · Ik+1,l + Uk,l+1 · Ik,l+1 + Uk,l−1 · Ik,l−1 + Uk−1,l+1 · Ik−1,l+1 + Uk+1,l−1 · Ik+1,l−1 ,

avec
" 6 "
∂φkl ∂φkl ∂φkl ∂φkl ∂φkl ∂φkl ∂φkl ∂φkl
! X !
Ik,l = + dxdy = + dxdy
Ω ∂x ∂x ∂y ∂y kl ∂x ∂x ∂y ∂y
" " "α=1 T α
" "
2 1 1 2 1
= 2 dxdy + 2 dxdy + 2 dxdy + 2 dxdy + 2 dxdy
h T1kl h T2kl h T3kl h T4kl h T5kl
"
1 8 8 h2
+ 2 dxdy = 2 (Aire triangle) = 2 · = 4;
h Tkl h h 2
6

"
∂φk−1,l ∂φkl ∂φk−1,l ∂φkl
!
Ik−1,l = + dxdy
∂x ∂x ∂y ∂y
" Ω
∂φk−1,l ∂φkl ∂φk−1,l ∂φkl
!
= + dxdy
T3kl ∪T4kl ∂x ∂x ∂y ∂y
" "
−1 1 2
= dxdy − dxdy = − · (Aire triangle) = −1;
h2 Tkl h2 Tkl h2
3 4

"
∂φk+1,l ∂φk+1,l ∂φkl ∂φkl
!
Ik+1,l = + dxdy
∂x ∂x ∂y ∂y
"Ω
∂φk+1,l ∂φk+1,l ∂φkl ∂φkl
!
= + dxdy
T1kl ∪T6kl ∂x ∂x ∂y ∂y
" "
1 1 2
= − 2 dxdy − 2 dxdy = − 2 · (Aire triangle) = −1;
h Tkl h Tkl h
1 6

"
∂φk,l−1 ∂φkl ∂φk,l−1 ∂φkl
!
Ik,l−1 = + dxdy
∂x ∂x ∂y ∂y
" Ω
∂φk,l−1 ∂φkl ∂φk,l−1 ∂φkl
!
= + dxdy
T4kl ∪T5kl ∂x ∂x ∂y ∂y
" "
−1 1 2
= 2
dxdy − 2 dxdy = − 2 · (Aire triangle) = −1;
h Tkl h Tkl h
4 5

"
∂φk,l+1 ∂φkl ∂φk,l+1 ∂φkl
!
Ik,l+1 = + dxdy
∂x ∂x ∂y ∂y
"Ω
∂φk,l+1 ∂φkl ∂φk,l+1 ∂φkl
!
= + dxdy
T1kl ∪T2kl ∂x ∂x ∂y ∂y
" "
−1 1 2
= 2
dxdy − 2 dxdy = − 2 · (Aire triangle) = −1;
h Tkl h Tkl h
1 2

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024

 32

"
∂φk−1,l+1 ∂φkl ∂φk−1,l+1 ∂φkl
!
Ik−1,l+1 = + dxdy
∂x ∂x ∂y ∂y
"Ω
∂φk−1,l+1 ∂φkl ∂φk−1,l+1 ∂φkl
!
= + dxdy = 0
Tkl ∪Tkl ∂x ∂x ∂y ∂y
2 3

"
∂φk+1,l−1 ∂φkl ∂φk+1,l−1 ∂φkl
!
Ik+1,l−1 = + dxdy
∂x ∂x ∂y ∂y
"Ω
∂φk+1,l−1 ∂φkl ∂φk+1,l−1 ∂φkl
!
= + dxdy = 0
Tkl ∪Tkl ∂x ∂x ∂y ∂y
2 3

Ainsi,
"
N−1 X
N−1
∂φi j ∂φkl ∂φi j ∂φkl
X !
Ui, j + dxdy = 4Ukl − Uk−1,l − Uk+1,l − Uk,l−1 − Uk,l+1 .
i=1 j=1 Ω ∂x ∂x ∂y ∂y

Le système peut encore se mettre sous la forme :

6 "
Uk−1,l − 2Uk,l + Uk+1,l Uk,l−1 − 2Uk,l + Uk,l+1 1 X
− − = 2 f φkl dxdy, 1 ≤ k, l ≤ N − 1,
h2 h2 h α=1 Tαkl

avec Uk,l = 0 sur ∂Ω, c’est-à-dire U0,l = UN,l = 0, 0 ≤ l ≤ N et Uk,0 = Uk,N = 0, 0 ≤ k ≤ N.

Pour la mise en équation, on pose
6 "
X
Ukl = Uk+(l−1)(N−1) et Fk+(l−1)(N−1) = f · φkl dxdy, 1 ≤ k, l ≤ N − 1,
α=1 Tαkl

En prenant N = 4, on obtient le système suivant :

    
 4 −1 0 −1 0 0 0 0 0   U1   F1 
 −1 4 −1 0 −1 0 0 0 0   U2   F2
    

 
 0 −1 4 0 0 −1 0 0 0   U3   F3
   
 
 −1 0 0 4 −1 0 −1 0 0   U4   F4
    

  U5  = 
     
 0 −1
 0 −1 4 −1 0 −1 0     F5 

 0
 0 −1 0 −1 4 0 0 −1   U6   F6 

    
 0
 0 0 −1 0 0 4 −1 0   U7   F7 

  
 0
 0 0 0 −1 0 −1 4 −1   U8   F8 

0 0 0 0 0 −1 0 −1 4 U9 F9
    

Exercices sur le chapitre

Résoudre par la méthode des éléments finis les exercices du chapitre précédent.
Indication : Utiliser les éléments de Lagrange d’ordre un.

Mathématique Appliquée: Notes de cours Pr J.D. Bukweli ©2024