MEC8270

ÉLÉMENTS FINIS EN THERMOFLUIDE

Navier-Stokes, Convections Forcée & Naturelle

A. Garon, S. Leclaire

Polytechnique Montréal
Plan du cours

Navier-Stokes

3.0 Les approximations Vitesse-Pression – COMSOL

4.0 Calcul du cisaillement pariétal
4.1 Affichage de ⃗τp
4.2 Affichage de ⃗τp – en 2D Axisymétrique
5.0 Formulation adimensionnelle
5.1 Calcul des forces, moments et puissance
 Chapitre 1 — Navier-Stokes 3.0 Les approximations Vitesse-Pression – COMSOL

Définition (Approximation vitesse-pression Pm Pn )

Pm approximation continue par un polynôme de Lagrange de
degré m des vitesses.
Pn approximation continue par un polynôme de Lagrange de degré
n de la pression.

Définition (Approximation vitesse-pression Pm Pm−1 )

Ces approximations se nomment approximations de Taylor-Hood.

Remarque
Les versions antérieures de COMSOL proposaient toutes les familles
d’approximation, dont les approximations de Crouseix-Raviart à
pression discontinue. Ces approximations sont supérieures en stabilité
et précision par rapport aux approximations de Taylor-Hood. En outre,
elles garantissent la conservation de la masse élément par élément,
contrairement à l’approximation de Taylor-Hood.

A. Garon MEC8270 2 / 27
 Chapitre 1 — Navier-Stokes 3.0 Les approximations Vitesse-Pression – COMSOL

Les approximations de COMSOL

Le logiciel COMSOL propose différentes approximations pour le couple
vitesse-pression.
En particulier, l’approximation P1 P1 : η η

3 3

∥⃗
u − u⃗exact ∥1,Ω = O(h)
∥p − pexact ∥0,Ω = O(h2 )
1 2 ξ 1 2 ξ
Instable en pression – doit être
Vitesse P1 Pression P1
stabilisée

A. Garon MEC8270 3 / 27
 Chapitre 1 — Navier-Stokes 3.0 Les approximations Vitesse-Pression – COMSOL

Remarque
Dans l’étude de l’approximation P1 P1 , nous avons calculé la solution
de Poiseuille. En particulier :
Celle-ci est quadratique en vitesse et linéaire en pression.
L’accélération convective est identiquement nulle : u⃗ · ∇⃗
u = 0.
La solution de Poiseuille est donc également solution de l’équation
de Stokes – écoulement rampant.

Donc les instabilités que nous observons ne proviennent pas du terme

d’accélération convective. Il s’agit d’une instabilité issue d’un choix
incompatible des interpolants des vitesses et de la pression. Ces
interpolants doivent satisfaire une condition théorique dont les divers
acronymes sont :
Condition inf-sup ou
Condition de Ladyzenkaja-Babuska-Brezzi (LBB) ou
Condition de compatibilité vitesse-pression.
A. Garon MEC8270 4 / 27
 Chapitre 1 — Navier-Stokes 3.0 Les approximations Vitesse-Pression – COMSOL

Remarque
L’approximation P1 P1 , ne satisfait pas la condition de compatibilité
(LBB) et doit toujours être stabilisée.

Cependant, la méthode de stabilisation n’est pas la méthode supg que

nous avons étudiée pour contrôler les instabilités convectives. Il s’agit de la
méthode pspg. Nous avons donc :

Pressure Stabilized Contrôler les instabilités des

PSPG
Petrov Galerkin approximations vitesse-pression
Streamline Upwing
SUPG Contrôler les instabilités convectives
Petrov Galerkin

Remarque (COMSOL)
En décochant l’option Streamline diffusion, la stabilisation pspg
devient inopérante.

A. Garon MEC8270 5 / 27
 Chapitre 1 — Navier-Stokes 3.0 Les approximations Vitesse-Pression – COMSOL

Les approximations de COMSOL – Taylor-Hood

Les approximations P2 P1 et P3 P2 appartiennent à la famille
Taylor-Hood. Ces approximations sont stables en pression et ne
requièrent pas de stabilisation pspg. Par contre, à haut nombre de
Reynolds, les instabilités convectives doivent être contrôlées à l’aide de
la méthode supg.
η η

3 3

L’approximation P2 P1 :
u − u⃗exact ∥1,Ω = O(h2 )
∥⃗ 6 5

∥p − pexact ∥0,Ω = O(h2 ) 1 4 2 ξ 1 2 ξ

η Vitesse P2 η Pression P1

3 3

L’approximation P3 P2 : 8 7

u − u⃗exact ∥1,Ω = O(h3 )

∥⃗ 9 10 6
6 5

∥p − pexact ∥0,Ω = O(h3 ) 1 4 5 2 ξ 1 4 2 ξ

Vitesse P3 Pression P2

A. Garon MEC8270 6 / 27
Plan du cours

Navier-Stokes

3.0 Les approximations Vitesse-Pression – COMSOL

4.0 Calcul du cisaillement pariétal
4.1 Affichage de ⃗τp
4.2 Affichage de ⃗τp – en 2D Axisymétrique
5.0 Formulation adimensionnelle
5.1 Calcul des forces, moments et puissance
 Chapitre 1 — Navier-Stokes 4.0 Calcul du cisaillement pariétal

⃗)
Définition (Vecteur-Cisaillement K
⃗ s’obtient de la contraction du tenseur des
Le vecteur-cisaillement K
contraintes avec la normale extérieure à la frontière du fluide.

Définition (Vecteur-Cisaillement-Pariétal ⃗τp )

⃗ en lui retirant sa
Le vecteur-cisaillement ⃗τp s’obtient du vecteur K
composante normale. Le vecteur ⃗τp est donc parallèle à la paroi.

⃗ col = K · n⃗col
K
⃗τp
 
⃗ col − n⃗lig · K
⃗τp = K ⃗ col n⃗col
y
La notation pour les indices sont :
⃗
K n⃗
col ⇒ vecteur colonne et
lig ⇒ vecteur ligne. x

A. Garon MEC8270 7 / 27
 Chapitre 1 — Navier-Stokes 4.0 Calcul du cisaillement pariétal

COMSOL
Construction 2D de τp :

⃗
K = spf.K stressx e⃗x + spf.K stressy e⃗y
⃗τp
Kn = spf.K stressx [Link] + y
spf.K stressy [Link]
⃗τp = (spf.K stressx − Kn [Link]) e⃗x + ⃗
K n⃗
(spf.K stressy − Kn [Link]) e⃗y x
Remarque
Dans cet exemple, nous obtenons
∂u
⃗τp = − µ e⃗x = spf.K stressx e⃗x
∂y
l’expression du cisaillement, si la paroi est parallèle à l’abscisse du
système de coordonnées.
A. Garon MEC8270 8 / 27
 Chapitre 1 — Navier-Stokes 4.0 Calcul du cisaillement pariétal

Définition (Cisaillement pariétal – τp )

Le cisaillement pariétal est la norme euclidienne de ⃗τp .

τpx = spf.K stressx − Kn [Link] ⃗τp

τpy = spf.K stressy − Kn [Link] y
q
τp = τpx 2 + τpy 2
⃗
K n⃗
τp = γ̇ µ x
Définition (taux de cisaillement pariétal – γ̇ (ou shear rate ))
Dans COMSOL nous obtenons le taux de cisaillement comme suit :
[Link]

A. Garon MEC8270 9 / 27
 Chapitre 1 — Navier-Stokes 4.0 Calcul du cisaillement pariétal

Le calcul du vecteur ⃗τp permet de vi-

sualiser le développement de la couche
limite au voisinage d’un point d’arrêt y
ou de décollement. Cette approche τ⃗p
est indépendante du système de co- n⃗
ordonnées et se généralise en 3D. x ⃗
K

Remarque
τp = 0 au point de stagnation de l’écoulement sur une paroi solide.
Cette valeur permet de déterminer la frontière d’une région de
recirculation.

A. Garon MEC8270 10 / 27
 Chapitre 1 — Navier-Stokes 4.0 Calcul du cisaillement pariétal

Considérons l’écoulement dans un conduit avec un changement de section

brusque. Dans cet exemple, le fluide se dirige de la gauche vers la droite.
Au niveau de l’expansion, le fluide se détache brusquement de la paroi
solide et se rattache en aval de ce point sur la paroi inférieure.
Γp

h Γe
y Γs

x
Γp

Norme de la vitesse

A. Garon MEC8270 11 / 27
 Chapitre 1 — Navier-Stokes 4.0 Calcul du cisaillement pariétal

Puisque l’écoulement est stationnaire, nous faisons usage des lignes de

courant pour identifier les zones de recirculation :

On distingue nettement la zone de recirculation après l’expansion.

Nous affichons ⃗τp sur la paroi inférieure.
Le cisaillement décroı̂t et est minimum au point de rattachement.

A. Garon MEC8270 12 / 27
 Chapitre 1 — Navier-Stokes 4.0 Calcul du cisaillement pariétal

Pour identifier avec précision le point de stagnation, on trace sur la

frontière le taux de cisaillement ou (en 2D) la composante-z de la vorticité.

En 2D, le taux de cisaillement ou la composante-z de la vorticité sont

équivalentes.
La composante-z de la vorticité permet de localiser plus facilement le
point de stagnation puisque qu’elle change de signe.

A. Garon MEC8270 13 / 27
Plan du cours

Navier-Stokes

3.0 Les approximations Vitesse-Pression – COMSOL

4.0 Calcul du cisaillement pariétal
4.1 Affichage de ⃗τp
4.2 Affichage de ⃗τp – en 2D Axisymétrique
5.0 Formulation adimensionnelle
5.1 Calcul des forces, moments et puissance
 Chapitre 1 — Navier-Stokes 4.1 Affichage de ⃗
τp

Pour afficher le vecteur ⃗τp sur toutes les parois du domaine de calcul (on
ne peut pas être sélectif !), on doit :
Sélectionner un 2D plot Group
Sélectionner dans ce groupe un Arrow Line
Construire les composantes du vecteur.

À l’exception des coins du domaine

de calcul, les cisaillements aux parois
solides sont en accord avec la solution
de l’écoulement de Poiseuille.
A. Garon MEC8270 14 / 27
Plan du cours

Navier-Stokes

3.0 Les approximations Vitesse-Pression – COMSOL

4.0 Calcul du cisaillement pariétal
4.1 Affichage de ⃗τp
4.2 Affichage de ⃗τp – en 2D Axisymétrique
5.0 Formulation adimensionnelle
5.1 Calcul des forces, moments et puissance
 Chapitre 1 — Navier-Stokes 4.2 Affichage de ⃗
τp – en 2D Axisymétrique

COMSOL
Construction 2D-Axi de τp :

⃗
K = spf.K stressr e⃗r + spf.K stressz e⃗z
⃗τp
Kn = spf.K stressr [Link] + z
spf.K stressz [Link]
⃗τp = (spf.K stressr − Kn [Link]) e⃗r + ⃗
K n⃗
(spf.K stressz − Kn [Link]) e⃗z r
Remarque
Dans cet exemple, nous obtenons
∂w
⃗τp = − µ e⃗r = spf.K stressr e⃗r
∂z
l’expression du cisaillement, si la paroi est parallèle à l’abscisse du
système de coordonnées.
A. Garon MEC8270 15 / 27
Plan du cours

Navier-Stokes

3.0 Les approximations Vitesse-Pression – COMSOL

4.0 Calcul du cisaillement pariétal
4.1 Affichage de ⃗τp
4.2 Affichage de ⃗τp – en 2D Axisymétrique
5.0 Formulation adimensionnelle
5.1 Calcul des forces, moments et puissance
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

L’analyse dimensionnelle permet d’identifier un nombre minimum de

paramètres qui caractérisent la physique du problème étudié. Nous allons
écrire les équations de Navier-Stokes adimensionnelles pour :
Faciliter les études paramétriques
Pour améliorer le conditionnement du système d’équations.

A. Garon MEC8270 16 / 27
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

Considérons l’exemple suivant d’un écoulement en conduite. Nous

adoptons la convention suivante :
Définition (Variable dimensionnelle)
Nous utilisons un astérisque ∗ pour identifier une variable
dimensionnelle.
Par exemple :
p∗ Pression dimensionnelle [Pa]
p Pression adimensionnelle [-]

Γ∗p

h∗ Γ∗e
y∗ Γ∗s

x∗
Γ∗p

A. Garon MEC8270 17 / 27
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

Pour écrire les équations sans dimension, on identifie des grandeurs

caractéristiques pour les variables dépendantes et indépendantes.
Définition (grandeur caractéristique)
Pour ne pas surcharger la notation, les grandeurs caractéristiques ne
sont pas suivies d’un astérisque – même si elles sont dimensionnelles.

L0 Longueur caractéristique
U0 Vitesse caractéristique
P0 Pression caractéristique

Remarque
Nous supposons que les propriétés sont constantes.

A. Garon MEC8270 18 / 27
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

Cas 1 – Le débit est connu à l’entrée

Remarque (L0 )

L0 = h∗
La longueur caractéristique est égale au diamètre de l’entrée.
Γ∗p

h∗ Γ∗e
x ∗ = L0 x y∗ Γ∗s
∗
y = L0 y
∗ x∗
L0 ∇ = ∇ Γ∗p

A. Garon MEC8270 19 / 27
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

Cas 1 – Le débit est connu à l’entrée

Remarque (U0 )
Z
1
U0 = u⃗∗ · n⃗ dΓ∗
L0 Γ∗e

La vitesse caractéristique est égale à la vitesse débitante.

Γ∗p

h∗ Γ∗e
y∗ Γ∗s
u ∗ = U0 u
v ∗ = U0 v
x∗
Γ∗p

A. Garon MEC8270 20 / 27
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

Cas 1 – Le débit est connu à l’entrée

Remarque (P0 )

P0 = ρ U02
La pression caractéristique est proportionnelle au carré de la
vitesse débitante.
Γ∗p

h∗ Γ∗e
y∗ Γ∗s
p ∗ = P0 p
T∗ = P0 T
x∗
Γ∗p

A. Garon MEC8270 21 / 27
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

Nous utilisons ces relations linéaires, entre les variables dimensionnelles

et adimensionnelles, pour transformer les équations de Navier-Stokes.

Ainsi nous transformons

u ∗ · ∇∗ )⃗
ρ(⃗ u ∗ = ∇∗ · T∗
∇∗ · u⃗∗ = 0

sous la forme
U2
 
P0
ρ 0 (⃗
u · ∇)⃗u =∇· T
L0 L0
U0
∇ · u⃗ = 0
L0

A. Garon MEC8270 22 / 27
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

Les facteurs d’échelle

Nous utilisons
U02 U0
ρ et
L0 L0
pour mettre à l’échelle respectivement les équations de Navier-Stokes
et de la conservation de la masse.

u · ∇)⃗
(⃗ u =∇·T
∇ · u⃗ = 0
avec
1  
T = −p I + ∇⃗ u )T
u + (∇⃗
Re
ρU0 L0
Re =
µ
A. Garon MEC8270 23 / 27
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

Si on effectue les simulations adimensionnelles, on doit :

1 Mettre à l’échelle la géométrie.
2 Mettre à l’échelle les conditions limites.
3 Mettre à l’échelle les forces volumiques.
4 Calculer le nombre de Reynolds.
5 Remplacer dans les équations :
1 ρ 7−→ 1
1
2 µ 7−→ Re

Γp

h Γe
y Γs

x
Γp

A. Garon MEC8270 24 / 27
 Chapitre 1 — Navier-Stokes 5.0 Formulation adimensionnelle

Les conditions limites de Dirichlet sur les vitesses doivent être,

naturellement, mises à l’échelle selon la relation

u = u ∗ /U0
v = v ∗ /Uo

et les conditions limites sur les forces de traction au bord du domaine (ce
sont des forces par unité de surface) par la relation

Fn = Fn∗ /P0
Ft = Ft∗ /P0

avec Fn la force de traction normale à la surface et Ft la force de traction

tangente à la surface.
De plus si des forces volumiques Fv∗ sont également présentes, on montre
que celles-ci sont mises à l’échelle par la relation suivante :
U02
Fv∗ = ρ Fv
L0
A. Garon MEC8270 25 / 27
Plan du cours

Navier-Stokes

3.0 Les approximations Vitesse-Pression – COMSOL

4.0 Calcul du cisaillement pariétal
4.1 Affichage de ⃗τp
4.2 Affichage de ⃗τp – en 2D Axisymétrique
5.0 Formulation adimensionnelle
5.1 Calcul des forces, moments et puissance
 Chapitre 1 — Navier-Stokes 5.1 Calcul des forces, moments et puissance

2D – dΓ∗ = dl ∗ [m]

Relation F ∗ et F – avec F⃗ ∗ en [N/m].

Z Z
F⃗ ∗ = ⃗ ∗ dΓ∗ = P0 L0
T ⃗ dΓ = P0 L0 F⃗
T
Γ∗ Γ

⃗ ∗ en [N-m/m].
Relation M ∗ et M – avec M
Z Z
⃗∗ =
M ⃗∗ × T
O ⃗ ∗ dΓ∗ = P0 L2 ⃗ ×T
O ⃗ dΓ = P0 L2 M
⃗
0 0
Γ∗ Γ

Relation Ẇ ∗ et Ẇ – avec Ẇ ∗ en [W/m].

Z Z
∗ ∗ ⃗ ∗ dΓ∗ = U0 P0 L0 ⃗ dΓ = U0 P0 L0 Ẇ
Ẇ = u⃗ · T u⃗ · T
Γ∗ Γ

A. Garon MEC8270 26 / 27
 Chapitre 1 — Navier-Stokes 5.1 Calcul des forces, moments et puissance

2D-Axi – dΓ∗ = 2 π r ∗ dl ∗ [m2 ]

Relation F ∗ et F – avec F⃗ ∗ en [N].

Z Z
F⃗ ∗ = ⃗ ∗ dΓ∗ = P0 L2
T 0
⃗ dΓ = P0 L2 F⃗
T 0
Γ∗ Γ

⃗ ∗ en [N-m].
Relation M ∗ et M – avec M
Z Z
⃗∗ =
M ⃗∗ × T
O ⃗ ∗ dΓ∗ = P0 L3 ⃗ ×T
O ⃗ dΓ = P0 L3 M
⃗
0 0
Γ∗ Γ

Relation Ẇ ∗ et Ẇ – avec Ẇ ∗ en [W].

Z Z
∗ ∗ ⃗ ∗ dΓ∗ = U0 P0 L2 ⃗ dΓ = U0 P0 L2 Ẇ
Ẇ = u⃗ · T 0 u⃗ · T 0
Γ∗ Γ

A. Garon MEC8270 27 / 27