UNIT

E DE PROGRAMME : S9UP1
Modelisation de la turbulence
Modelisation de la turbulence
Auteur : Yann MARCHESSE
Departement : Mecanique et

Energetique

Edition : Annee universitaire 2009-2010

ECOLE CATHOLIQUE DARTS ET M

ETIERS
40 Montee Saint-Barthelemy - 69321 Lyon Cedex 05
Tel. : 04 72 77 06 00 - Fax : 04 72 77 06 11
[Link]
1
Mecanique des uides
Modelisation de la turbulence - Application `a la CFD
Yann MARCHESSE
Departement de Mecanique et

Energetique

Ecole Catholique dArts et Metiers - Lyon

Date de compilation du document : 9 mars 2010

ECOLE CATHOLIQUE DARTS ET M

ETIERS
40 Montee Saint-Barthelemy - 69321 Lyon Cedex 05
Tel. : 04 72 77 06 00 - Fax : 04 72 77 06 11
[Link]
2
Sommaire
Avant-propos 7
Nomenclature 9
I Modelisation de la turbulence en vue dune application `a
la CFD 11
1 TURBULENCE 13
1.1 Renseignements sur lagitation turbulente . . . . . . . . . . . . . . . . . . . 14
1.2 Structures coherentes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
1.3 Passage du laminaire vers le turbulent : la transition . . . . . . . . . . . . 17
1.3.1 Cas dune conduite . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
1.3.2 Cas dune couche limite sur une plaque plane . . . . . . . . . . . . 18
1.4 Inuence de la paroi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2

EQUATIONS DE LA M

ECANIQUE DES FLUIDES 23

2.1 Derivee particulaire . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
2.2

Equation de continuite . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.3 Bilan de quantite de mouvement . . . . . . . . . . . . . . . . . . . . . . . . 27
2.4

Equation de lenergie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
2.5 Bilan sur les equations obtenues par lapplication des principes fondamentaux 32
2.6 Conditions limites . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
3 MOD

ELISATION DE LA TURBULENCE 35
3.1 Trois approches numeriques . . . . . . . . . . . . . . . . . . . . . . . . . . 35
3.1.1 Simulation directe numerique (DNS) . . . . . . . . . . . . . . . . . 36
3.1.2 Simulation des grandes echelles (LES) . . . . . . . . . . . . . . . . 36
3.1.3 Modelisation statistique de la turbulence (RANS) . . . . . . . . . . 38
3.1.4 Remarques sur les trois approches . . . . . . . . . . . . . . . . . . . 41
3.2 Mod`ele de turbulence du premier ordre - Concept de la viscosite turbulente 41
3.2.1 Mod`ele algebrique ou mod`ele `a zero equation . . . . . . . . . . . . . 43
3.2.2 Mod`eles de fermeture `a une equation de transport . . . . . . . . . . 44
3.2.3 Mod`eles de fermeture `a deux equations de transport . . . . . . . . . 46
3.2.4

Equation nale . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
3.3 Mod`ele de turbulence du deuxi`eme ordre . . . . . . . . . . . . . . . . . . . 50
3.4 Inuence de la paroi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
3.5 Choisir le mod`ele de turbulence . . . . . . . . . . . . . . . . . . . . . . . . 52
4 SOMMAIRE
3.6 Modelisation instationnaire . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
3.7 Et lindustrie dans tout ca ? . . . . . . . . . . . . . . . . . . . . . . . . . . 53
3.8 Resolution numerique `a partir dun code de calculs . . . . . . . . . . . . . 55
4 Application de la modelisation de la turbulence `a la CFD 57
4.1 Structure dun code CFD . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58
4.1.1 Le pre-processeur . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58
4.1.2 Le solver . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58
4.1.3 Le post-processeur . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
4.2 Le maillage . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
4.3 Les conditions limites . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
4.4 Les conditions initiales . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60
4.5 Les schemas numeriques de convection . . . . . . . . . . . . . . . . . . . . 61
4.6 Convergence dun calcul . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
4.6.1 Les residus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
4.6.2 La sous-relaxation . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
4.6.3 Calcul permanent `a partir dune approche pseudo-transitoire . . . . 63
II Exemples de predictions numeriques decoulements et dechanges
thermiques 65
Presentation des etudes 67
5 Cas test : 2D Model Hill Flow 69
5.1 Approche numerique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
5.1.1 Domaine de calcul et maillage . . . . . . . . . . . . . . . . . . . . . 70
5.1.2 Conditions limites . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71
5.2 Resultats numeriques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71
6

Ecoulement au passage dun obstacle prismatique 73
6.1 Approche numerique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
6.1.1 Domaine de calcul et maillages . . . . . . . . . . . . . . . . . . . . 74
6.1.2 Conditions limites . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
6.1.3 Mod`eles de turbulence . . . . . . . . . . . . . . . . . . . . . . . . . 75
6.1.4 Discretisation temporelle . . . . . . . . . . . . . . . . . . . . . . . . 75
6.2 Resultats numeriques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
6.2.1 Validation du maillage sur le pave . . . . . . . . . . . . . . . . . . . 76
6.2.2 Eorts aerodynamiques . . . . . . . . . . . . . . . . . . . . . . . . . 76
6.2.3 Sillage du pave . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
6.3 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77
7

Etude des param`etres turbulents dans un canal 79
7.1 Approche numerique du calcul annexe . . . . . . . . . . . . . . . . . . . . 79
7.2 Resultats numeriques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81
7.3 Approximation analytique des grandeurs turbulentes . . . . . . . . . . . . 82
7.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84
SOMMAIRE 5
8 Prediction numerique dun echange thermique convectif `a partir de
dierents mod`eles de turbulence 85
8.1 Approche numerique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86
8.1.1 Domaine de calcul et maillages . . . . . . . . . . . . . . . . . . . . 86
8.1.2 Conditions limites . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86
8.1.3 Mod`eles de turbulence . . . . . . . . . . . . . . . . . . . . . . . . . 87
8.2 Resultats numeriques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87
9 Inuence de la nesse du maillage sur la prediction numerique des
echanges thermiques convectifs dune plaque plane 91
9.1 Approche numerique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91
9.1.1 Domaine et maillage . . . . . . . . . . . . . . . . . . . . . . . . . . 91
9.1.2 Conditions limites . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
9.2 Mod`ele propose . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
9.3 Resultats numeriques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
Lectures interessantes 97
6 SOMMAIRE
Avant-propos
Les ecoulements rencontres dans lindustrie ou tout simplement au quotidien sont pour
la plupart domines par des mouvements chaotiques : la turbulence. Au-del`a de la dif-
fusion moleculaire, les echanges de masse, de quantite de mouvement et de chaleur sont
ainsi essentiellement regis par ce mode decoulement. Limpact de la turbulence peut etre
positif ou negatif, lingenieur doit donc etre capable de predire ces eets lors de la concep-
tion de syst`emes. Malheureusement, le mouvement turbulent est tr`es complexe et presente
la plupart du temps des caracteristiques tridimensionnelles et instables. Ce mouvement
consiste en la surperposition de tourbillons dont le spectre de taille est tr`es large, celui-ci
balaye des grosses structures dont la taille depend de la geometrie et correspondant `a des
uctuations de basses frequences jusquaux petites structures associees `a des uctuations
de hautes frequences dont lechelle est celle de la dissipation energetique.
La simulation numerique peut etre dun grand secours `a lingenieur, le calcul analytique
demeurant ici tr`es peu ecace. Il serait toute fois imprudent de considerer cette approche
comme etant lunique voie `a emprunter, en eet la modelisation des ecoulements vient
plutot en complement des essais sur site ou sur maquette.
Lors dune etude, la CFD intervient souvent en amont lors de la conception et evite
alors aux constructeurs de fabriquer de nombreux et surtout co uteux prototypes. Elle peut
alors etre caracterisee de prototypage virtuel. En eet, un des principaux avantages du
calcul numerique est la possibilite de faire varier les param`etres geometriques, dynamiques
ou thermophysiques du probl`eme traite en evitant la repetition dexperiences longues et
lourdes `a gerer. Ensuite, plus en aval de letude, elle peut etre utile `a lanalyse davaries
decelees sur des equipements ou pour ameliorer leur performance.
Le traite de ce cours sinscrit dans cette problematique et est compose de deux parties.
La turbulence est en premier lieu examinee (chapitre 1), suivie des ecritures des
equations mathematiques traduisant trois principes fondamentaux (chapitre 2). On sinte-
resse alors aux consequences des termes uctuants, une fois ceux-ci introduits dans les
equations du mouvement des uides. Cette etape presente essentiellement le tenseur de
Reynolds que lon cherchera ensuite `a estimer `a partir de mod`eles de fermeture (chapitre
3). Ensuite, larticulation dun code de calcul est bri`evement presentee (chapitre 4).
La deuxi`eme partie du document vise dune part `a presenter quelques ecoulements
simules dans des congurations standards, mais surtout contribue `a sensibliser lutilisateur
de codes de calculs aux conditions limites en entree de domaine, au choix des mod`eles de
turbulence, `a la nesse du maillage, et `a la pertincence de lemploi dune loi de paroi.
Enn, cet expose a pour but unique de familiariser le lecteur `a la modelisation de
la turbulence. On ny trouvera donc que le strict minimum des concepts existants. Le
lecteur, plus curieux, est invite quant `a lui `a consulter quelques ouvrages cites en annexe.
8 Avant-propos
Par contre, on trouvera ici quelques informations pratiques pour utiliser tous les logiciels
de calcul.
Nomenclature
La meme notation peut dans certains cas representer des quantites physiques distinctes,
mais replacee dans le contexte du chapitre traite, toute ambiguite devrait etre levee.
Caract`eres usuels
c
v
Chaleur specique `a volume constant
C
k
Constante de Kolmogorov (Eq. 3.6)
D Volume de controle
e

Energie interne par unite de masse
E(k) Spectre de lenergie cinetique turbulente
f
i
Eort volumique
f
v
, f

Coecients damortissement (Eqs. 3.9 et 3.15)

k

Energie cinetique turbulente moyenne k = u

i
u

j
/2
k Conductivite thermique

k Nombre donde

k = 2/
l
m
Longueur de melange
m masse

n Vecteur unitaire sortant dun volume de controle

p Pression
P
k
Production denergie cinetique turbulente P
k
= u

i
u

j
u
i
/x
j
q Densite volumique du taux de chaleur

q Vecteur densite surfacique de chaleur re cue par conduction

r Constante des gaz parfaits
Re Nombre de Reynolds
Re
c
Nombre de Reynolds critique
S Surface dun volume de controle
S
ij
Tenseur des vitesses de deformation
t temps
T Temperature
T Duree dobservation dans loperateur moyenne
T
w
Temperature de paroi
u Composante selon x du vecteur vitesse
u
CL
Vitesse au centre dune conduite
10 Nomenclature
Caract`eres usuels (suite)
u
f
Vitesse de frottement u
f
=
_

p
/
v Composante selon y du vecteur vitesse

V Vecteur vitesse decoulement

V Volume dun volume de controle
w Composante selon z du vecteur vitesse
y
+
Distance adimensionnee `a la paroi
Caract`eres grecs

Epaisseur de couche limite

ij
Symbole de Kronecker
ij
= 1 si i = j et 0 sinon
Taux de dissipation moyen de lenergie cinetique turbulente
Constante de von Karman (Eq. (1.3))
Longueur donde
Coecient de pertes de charge p =
L
D

V
2
2
Viscosite dynamique

e
Viscosite eective

T
Viscosite turbulente
Viscosite cinematique
Masse volumique

k
,

Nombres de Prandtl turbulent (Eqs. 3.7 et 3.14)

ij
Tenseur des contraintes

p
Contrainte parietale

ij
Tenseur des contraintes visqueuses
Taux de dissipation specique moyen
Premi`ere partie
Modelisation de la turbulence en vue
dune application `a la CFD
Chapitre 1
TURBULENCE
La plupart des ecoulements naturels observes quotidiennement sont turbulents (fu-
mee de cigarette, cr`eme de lait versee dans le cafe, nuages). Ils sont dautre part tr`es
diversies : uides biologiques (sang), mouvements de geouides (vents, courants marins),
mouvements de uides stellaires (circulation gazeuse autour des plan`etes). Malgre cela, ces
ecoulements ont des proprietes en commun que nous enoncerons plus loin. Les equations
qui gouvernent le mouvement instantane des uides, quils soient turbulents ou non, ont
ete ecrites par Claude Navier
1
en 1823. Elles sont appelees equations de Navier-Stokes
en raison des perfectionnements apportees ulterieurement par George Stokes
2
. Il sagit ni
plus ni moins des equations de Newton quil faut appliquer `a une particule uide. Ces
equations avaient ete prealablement ecrites par Euler mais Navier eut lingeniosite dy
rajouter un terme de friction entre les diverses couches de uide.
Lequation de Navier-Stokes pour la vitesse instantanee est connue depuis longtemps,
cependant sa resolution reste trop compliquee. Une solution a alors ete dessayer de pro-
poser une solution pour lecoulement moyen. Malheureusement le passage de lequation de
Navier-Stokes `a lequation pour la moyenne fait apparatre un terme qui nest pas connu
de mani`ere exacte. Ce terme represente leet des uctuations sur la vitesse moyenne et
doit etre approche ou modelise.
La turbulence est devenue une science experimentale vers la n du XIX si`ecle quand
1
Claude Navier (1785-1836), physicien fran cais qui, en 1921, proposa les equations du mouvement vi-
bratoire dun solide, generalisees par Cauchy, et resolue par Poisson. Stokes contribua aussi aux equations
de mouvement de uide, connues sous le nom dequations de Navier-Stokes.
2
George Stokes (1819-1903), mathematicien et physicien irlandais qui contribua enormement aux
avancees scientiques en son temps. Il ameliora `a la fois la formulation des equations de mouvement
des uides mais aussi leur comprehension.
14 Turbulence
langlais Osborne Reynolds
3
a pu observer la transition du regime laminaire vers le regime
turbulent. Il mit ainsi en evidence quelques lois assez simples et introduisit un nombre
adimensionnel portant son nom qui caracterise cette transition. Malgre tout, avant les
annees 1950, la turbulence etait un sujet obscur. La seule issue pour lingenieur etait
dexperimenter sur des mod`eles physiques an dameliorer son savoir-faire. Heureuse-
ment apr`es les annees 1960, la situation allait se debloquer avec quelques progr`es ac-
complis en mati`ere de modelisation, en meme temps que la capacite des traitements
numeriques augmentait fortement. Cependant, la predetermination de proprietes statis-
tiques locales etait encore impossible. De plus quelques probl`emes majeurs demeuraient :
limpredicibilite, luniversalite des mod`eles etablis, et la convergence des fermetures. Ac-
tuellement, ces probl`emes ont ete en partie resolus et les calculs numeriques permettent
une bonne estimation de lecoulement moyen en presence dune turbulence developpee.
Des methodes plus recentes permettent aussi destimer le champ uctuant par resolution
des equations de Navier-Stokes instantanees. On parle de fa con abusive dexperience nu-
merique.
1.1 Renseignements sur lagitation turbulente
Un ecoulement represente un glissement de particules uides les unes sur les autres.
Lagitation moleculaire entrane `a son echelle des echanges entre ces couches et une dif-
fusion de lecoulement (non perceptible `a lechelle du milieu continu, i.e. `a lechelle de la
particule uide).
A lechelle du milieu continu, le regime turbulent presente, en plus, des mouvements
desordonnes et tridimensionnels. Notons tout de meme que bien que ce que nous
observons soit tr`es complexe, tr`es desordonne, cest tr`es loin detre le desordre total. En
eet, la turbulence est composee de structures coherentes
4
, en particulier de tourbillons
(que lon observe en aval dun pont par exemple). Il sagit l`a dun melange subtile dordre
et de desordre qui en fait une des principales dicultes vis `a vis de la modelisation. Toute
cette agitation porte le nom dagitation turbulente.
Notons dautre part quil sagit dun mouvement secondaire de lecoulement. La tur-
bulence nest donc pas liee `a la nature du uide mais `a son mode de mouvement. On traite
la plupart du temps de la turbulence pleinement developpee, caracterisee par :
Taille supra-moleculaire : la taille admissible par une structure turbulente (on
dira coherente) est guidee dune part par letendue spatiale disponible (dimensions
dun canal, epaisseur de la couche limite) et dautre part par la viscosite. Il existe
en eet une echelle appelee echelle de Kolmogorov
5
sinterpretant comme etant la
3
Osborne Reynolds, physicien anglais (1842-1912). Ces premiers travaux concernent le magnetisme et
lelectricite. Apr`es 1873, il se concentre principalement sur la dynamique des uides. Il etudie entre autre
les changements de regime dun ecoulement dans une conduite. Son experience restera cel`ebre et portera
son nom. En 1886, il publie The theory of lubrication et invente la tribologie (etude des frottements).
4
Cette coherence nest neanmoins pas observee sur le trajet entier de la structure. Il existe une distance
sur laquelle la structure reste coherente avec elle-meme, apr`es quoi elle sapparie avec une autre pour
former une autre structure de plus grande taille.
5
Andrei Nikolaevich Kolmogorov (1903-1987), physicien et mathematicien russe. Il est `a lorigine de
la theorie de la turbulence aux petites echelles. Il demontra entre autres la loi de puissance en k
5
3
du
Renseignements sur lagitation turbulente 15
plus petite dimension des structures ou tourbillons que lon peut rencontrer dans
un ecoulement turbulent. En-dessous de cette echelle, les eets visqueux font leurs
eets, et lenergie mecanique est totalement transformee en chaleur. Comprenons
bien quil existe alors des structures de taille plus elevee capable de dissiper. Pour
se donner une idee, dans latmosph`ere terrestre par exemple, lechelle de dissipation
est de lordre du millim`etre, alors que les plus grosses structures ont des echelles de
plusieurs milliers de kilom`etres.
Comportement aleatoire : les mesures mettent en evidence laspect chaotique
de toute fonction du champ de lecoulement (temperature, pression, vitesses, masse
volumique). On perd alors toute notion de predictibilite presente dans lecoulement
laminaire. Une approche statistique est alors necessaire (ex : la meteo).
Denombrement inni : la mesure de la turbulence donne des resultats dierents
selon le point considere. Neanmoins, ceux-ci font apparatre des oscillations recou-
vrant un tr`es large spectre en frequence
6
(caracterisees par le nombre donde 2/).
Structures tridimensionnelles : il ny a pas de direction priviligiee meme en
presence dun ecoulement considere comme bidimensionnel.
Intermittence
7
: ce phenom`ene est observe dans des cas precis : lors du passage
laminaire/turbulent, dans les petites echelles (dissipation), et proche des fronti`eres
libres (couche de melange dun jet). Ce phenom`ene est lie `a la penetration massive
de uide non turbulent dans la couche de melange.
Cinematique rotationnelle : les uctuations de vitesse sont porteuses de uc-
tuations de rotationnel. La turbulence ne produit pas de rotationnel mais aura
des eets sur cette quantite tels quun renforcement de sa production proche des
parois (gradient de vitesse important), ou la diusion de vorticite hors des zones de
production sous leet de lagitation moleculaire.
Dynamique non lineaire : interaction entre tous les mouvements presents dans
lecoulement selon une dynamique qui est regie par les equations de Navier-Stokes.

Energie dissipative : la dissipation est donnee par = 2S
ij
S
ij
, avec S
ij
le ten-
seur des vitesses de deformation.

Etant donnees les faibles valeurs de pour leau
ou lair, seules les echelles presentant des vitesses de deformation importantes vont
dissiper : il sagit des petites echelles. Pour se donner une image, les grosses struc-
tures (faible nombre donde) sont porteuses denergie (ex : les bourasques de vent).
Celles-ci se cassent en des structures de plus en plus nes sans production ni
dissipation denergie : cest la cascade denergie (Fig. 1.1). Ce processus continue
jusqu`a ce que la structure soit assez ne (nombre donde important) pour que les
eets de viscosite puissent agir et que la dissipation ait lieu par transfert de lenergie
spectre denergie cinetique representant le mouvement turbulent.
6
La large gamme de frequence est `a mettre en parall`ele avec celle des tailles de structure.
7
Zone o` u lon observe de fa con alternee un melange dagitation et de calme au sein dun ecoulement
turbulent.
16 Turbulence
cinetique turbulente en energie interne. Ce processus, appele Cascade denergie, a
ete enonce par Richardson en 1922 [30] et a ete complete par Kolmogorov en 1941
[14].

k =
2

E(

k)
Dissipation
Fig. 1.1 Spectres
denergie et de dissipa-
tion.
`
A partir de toutes les caracteristiques vues ci-dessus, on propose une denition de la
turbulence :
La turbulence est un mode naturel decoulement dun uide visqueux o` u des
mecanismes internes dechange denergie assurent la creation et le maintien
de toute hierarchie de mouvements chaotiques repartis contin ument sur une
large gamme dechelles macroscopiques.
1.2 Structures coherentes
Nous lavons vu plus haut, la turbulence apparat comme un mode decoulement chao-
tique pour lequel un tourbillon cree est immediatement etire dans toutes les directions,
et nalement detruit par la turbulence. Il existe neanmoins des tourbillons coherents ou
structures coherentes au sein de cette turbulence, cest `a dire des structures ayant un
temps de vie important devant leur temps de retournement
8
selon la denition enoncee
par Marcel Lesieur [20] (Fig. 1.2).
Fig. 1.2 Structures
coherentes dans une couche
de melange air/azote (dapr`es
Brown [3]).
8
Il sagit l` a du temps necessaire au tourbillon pour quil puisse executer une rotation autour de lui-
meme.
Passage du laminaire vers le turbulent : la transition 17
1.3 Passage du laminaire vers le turbulent : la tran-
sition
La transition dun ecoulement laminaire vers un ecoulement turbulent a ete en premier
lieu etudiee de fa con experimentale par O. Reynolds [29]. Cette transition est un probl`eme
tr`es complexe. En consequence, il nexiste pas `a lheure actuelle de mod`ele theorique per-
mettant de predire ecacement le comportement du uide lors de cette phase. Soulignons
tout de meme la premi`ere approche theorique proposee par Landau en 1944 [15] dans
laquelle celui-ci consid`ere la transition comme lapparition de mouvements periodiques
non-stationnaires dont la phase initiale est arbitraire. Landau baptise ce premier degre
dindetermination degre de liberte (introduction historique de ce terme). Lorsque le
rapport des eorts dinertie et de viscosite augmente, le mouvement devient instable
conduisant `a lapparition dun nouveau mouvement periodique (deux degres de liberte).
On atteint alors tr`es vite un nombre tr`es important de degre de liberte, correspondant au
cycle limite de Landau.
Le point de transition est dicilement determine par lobservation. On pref`erera eec-
tuer des mesures et etudier levolution dun param`etre qui depend du nombre de Reynolds
(coecient de perte de charge dans une conduite, epaississement de la couche limite le
long dune plaque plane, etc.).
1.3.1 Cas dune conduite
La transition vers un regime turbulent, dans le cas dun ecoulement en conduite, peut
etre determinee par le changement devolution des pertes de charge. Ces derni`eres sont
caracterisees par un coecient appele coecient de pertes de charge, tel que
p =
L
D

V
2
2
o` u L est la distance entre les deux mesures de pression, D le diam`etre de la conduite, et V
la vitesse moyenne de lecoulement. Le coecient de pertes de charge depend du nombre
de Reynolds et du rapport de la rugosite relative /D caracterisant letat de surface de la
conduite.
Lorsque le nombre de Reynolds est inferieur `a 2000, lecoulement est laminaire, et le
coecient de pertes de charge est le meme quelle que soit la rugosite. Il vaut dans ce cas :
=
64
Re
: Hagen/Poiseuille
Lecoulement quitte le regime laminaire lorsque le nombre de Reynolds depasse 2000, la
region caracterisee par des nombres de Reynolds compris entre 2000 et 4000 etant appelee
regime critique. Le coecient de pertes de charge est tr`es mal deni dans cette region car
la turbulence nest pas compl`ete. En eet, lecoulement alterne entre des regimes laminaire
et turbulent.
Pour des nombres de Reynolds superieurs `a 4000, lecoulement devient turbulent et le
coecient de pertes de charge depend du nombre de Reynolds et de la rugosite relative.
Pour des conduites lisses, et d`es lors que 4000 < Re < 100000, le coecient est estime
par la relation :
=
0, 316
Re
1/4
: Blasius
18 Turbulence
Levolution du coecient de pertes de charge nest pas lineaire decroissante. En eet,
celui-ci voit ses valeurs decroitrent tant que le nombre de Reynolds reste inferieur `a une
premi`ere valeur critique (notee Re
c1
). Le nombre de Reynolds augmente brutalement
au passage de la transition. Pour des nombres de Reynolds superieurs `a une deuxi`eme
valeur critique (notee Re
c2
), le comportement turbulent se traduit par une nouvelle loi de
decroissance dierente de celle observee en laminaire. On pourra alors trouver un crit`ere
de transition, Re
c
tel que
Re
c
=
1
2
[Re
c1
+ Re
c2
] (1.1)
Levolution des pertes de charge avec le nombre de Reynolds est etudiee sur le tube present
au departement de mecanique de lECAM (TP

Etude de lecoulement dans un tube).

Le diam`etre du tube vaut 3 mm, et la distance entre les deux sections de mesure est de
523 mm. Les pertes de charge sont estimees `a partir de tubes piezometriques (eau/air
pour les faibles debits, et eau/mercure pour des debits plus importants). Levolution du
coecient de perte de charge est conforme `a la precedente discussion (Fig. 1.3). Le nombre
de Reynolds critique, estime par la relation (1.1), vaut Re
c
= 2800. On note dautre part
que les valeurs du coecient de perte de charges obtenues `a partir des deux mod`eles
precedents (Hagen/Poiseuille et Blasius) sont tr`es proches des valeurs mesurees.
Fig. 1.3

Evolution du coecient de pertes de charge avec le nombre de Reynolds. Mesures eectuees au
laboratoire de mecanique : , manom`etre Eau/Air ; , manom`etre Eau/Mercure.
1.3.2 Cas dune couche limite sur une plaque plane
Considerons une plaque plongee dans un ecoulement de vitesse U. La vitesse sannule
`a la surface de la plaque, et le passage de cette vitesse nulle `a la vitesse de lecoulement
externe seectue sur une distance tr`es faible appelee couche limite (Region notee 1
sur la gure 1.4). Ainsi la vitesse longitudinale evolue tr`es rapidement avec la distance
transversale dans cette region. Pour des positions proches du bord dattaque (x = 0), la
couche limite est laminaire, et lepaisseur de couche limite, (x) evolue selon x
0,5
. Pour
Inuence de la paroi 19
une distance plus en aval, des instabilites se developpent, caracterisant la transition vers
la turbulence. Au-dela de cette abscisse critique, lepaisseur de couche limite evolue selon
x
0,8
. Le nombre de Reynolds base sur la distance au bord dattaque, Re
x
= Ux/, est
habituellement utilise ici pour caracteriser cette transition. La gure 1.5 rassemble des
mesures de lepaississement de la couche limite en fonction du Reynolds local. On observe
la transition du regime decoulement `a partir de Re = 3, 2 10
5
par le changement
devolution de .
Fig. 1.4 Representation schematique de la transition dune couche limite laminaire vers un regime turbulent
pour une plaque plane sans gradient de pression transversal.
Fig. 1.5

Evolution de lepaississement de la couche limite en
fonction du nombre de Reynolds local. Notez le changement
devolution de lepaississement de couche limite `a partir de
Re = 3, 2 10
5
(dapr`es Schlichting [32]).
1.4 Inuence de la paroi
La turbulence telle que nous lavons denie plus haut est limitee `a une turbulence
pleinement developpee. Lorsquune paroi est proche, les importants gradients modient
fortement cette turbulence. La gure 1.6 presente la repartition des dierents phenom`enes
lies `a lenergie cinetique turbulente prenant place dans la couche limite (i.e. production
et dissipation de lenergie cinetique turbulente, etc.). Ces estimations ont ete eectuees `a
partir de la methode numerique DNS par Moser et al. [26] dans le cas dun ecoulement
dans un canal. Les resultats experimentaux conrment ces tendances, meme sil est tr`es
20 Turbulence
dicile de les obternir par cette voie. On note dapr`es ces resultats quil y a equilibre
entre la production et la dissipation denergie cinetique turbulente pour y
+
> 30 jusqu`a
une distance telle que y/h 0, 5 (cette derni`ere valeur est absente des gures). Ceci tend
`a montrer le caract`ere isole de cette region, car netant pas perturbee par la paroi proche
et lecoulement principal (les autres termes de transfert vus plus haut etant negligeables).
Notons dautre part que dans la region caracterisee par y
+
< 30, la dissipation nest plus
equilibree par le terme de production (les termes uctuants deviennent nuls par la presence
de ladherence) mais par un terme de diusion visqueux. Ces deux regions constituent la
region interne de la couche limite et est denie de fa con plus generale par la relation
y < 0, 1, avec lepaisseur locale de la couche limite.
Fig. 1.6 Bilan de lenergie cinetique
turbulente dans la region interne dune
couche limite (dapr`es Moser [26]). y
+
est
la distance `a la paroi adimensionnee.
Dans le but de decrire la couche limite, on sinteresse generalement `a levolution trans-
versale de la vitesse longitudinale dans la couche limite. On peut baser notre reexion sur
une approche globale pour laquelle on cherche `a decrire notre probl`eme `a partir de pro-
duits sans dimension. Dapr`es ce que nous venons de voir, la vitesse locale u depend
des contraintes locales de frottement exercees par la paroi (
p
), des proprietes du uide
( et ), et de la distance `a la paroi (y). Il est donc possible `a partir du theor`eme de
Vashy-Buckingham de construire deux nombres sans dimension :
u
+
=
u
_

p
/
et y
+
=
y
_

p
/

avec
_

p
/ homog`ene `a une vitesse que lon nomme vitesse de frottement, notee u
f
, et
ne possedant pas de signication physique (on pref`ere neanmoins manipuler une vitesse
au lieu dune tension surfacique). Pour des vitesses moyennes, on aura en eet plutot
interet `a sinteresser au rapport U/u
f
, et non pas `a U seule. De meme, il est preferable
de regarder le rapport u

i
u

j
/u
2
f
au lieu de u

i
u

j
seul. La vitesse de frottement est de lordre
de 5% de la vitesse de debit.
Levolution transversale de la vitesse longitudinale dans la region interne de la couche
limite temoigne de plusieurs zones composant cette region interne de la couche limite :
(1) la sous-couche visqueuse (y
+
< 5) dans laquelle le tenseur de Reynolds est negligeable
devant les contraintes visqueuses (Fig. 1.7) et dans laquelle la vitesse evolue lineairement
avec la distance `a la paroi :
u
+
= y
+
(1.2)
Inuence de la paroi 21
et (2) la zone logarithmique (30 < y
+
< 1000) dominee par les contraintes turbulentes
(Fig. 1.8) et dans laquelle une evolution logarithmique de la vitesse est observee :
u
+
=
1

lny
+
+B (1.3)
et B etant des constantes du mod`ele obtenues experimentalement sur une plaque plane
et lisse (respectivement 0,41 et 5,2). La region situee entre ces deux zones precedentes
(5 < y
+
< 30) est appelee zone tampon. La region interne de la couche limite ne depend
donc que des param`etres y
+
et u
f
, et reste donc isolee de lecoulement principal.
Fig. 1.7 Prol de la vitesse moyenne dans la
sous-couche visqueuse. Estimation numerique `a
partir de la methode DNS (ligne hachuree, Re = 5
600 ; ligne, Re = 13 750 ; ligne hachuree + point,
u
+
= y
+
) (dapr`es Kim [12]).
Fig. 1.8 Prol de la vitesse moyenne dans la
region interne de la couche limite. , Re = 2
970 ; , Re = 4 914 ; , Re = 22 776 ; , Re =
39 582 ; ligne, loi logarithmique, Eq. (1.3) (dapr`es
Weil [36]).
Au-del`a de la zone logarithmique, la relation (1.3) nest plus valable car les eets de
lecoulement principal sont de moins en moins negligeables. La relation de vitesse qui ne
dependait que de y
+
doit tendre de fa con asymptotique vers une relation cette fois fonction
de y/. On entre alors dans une region caracterisee par une loi de vitesse decitaire et
caracterisee par une relation du type :
U
o
u
u
f
= F(
y

)
22 Turbulence
o` u U
o
est la vitesse de lecoulement principal, et F une fonction `a determiner.
Ce chapitre a presente la turbulence et ses caracteristiques an de proposer une
denition. Les mod`eles de turbulence ont pour but de formuler au mieux le compor-
tement des ecoulements turbulents. Nous verrons par la suite que ces mod`eles sont plus
ou moins complexes, et bases la plupart du temps sur lempirisme.
Le prochain chapitre presente bri`evement les equations mathematiques traduisant trois
principes fondamentaux de la mecanique des uides : la conservation de la masse, le
principe fondamental de la dynamique, et la conservation de lenergie.
Chapitre 2

EQUATIONS DE LA M

ECANIQUE
DES FLUIDES
La prediction numerique des ecoulements de uide est basee sur la resolution dequations
de bilan traduisant trois principes fondamentaux : (1) la conservation de la masse, (2)
le principe fondamental de la dynamique, et (3) la conservation de lenergie. Les trois
equations mathematiques resultantes sont lequation de continuite, le bilan de quantite de
mouvement, et lequation de lenergie. Le but de ce chapitre est decrire ces equations.
Ces equations de bilan peuvent etre ecrites sous plusieurs formes (conservatives, non
conservatives, etc.) nentrainant pas de reels changements dans leur comprehension. Nean-
moins, celles-ci deviennent importantes en CFD, car une forme particuli`ere peut etre
ecace dans certains types decoulements, alors quelle peut conduire `a des instabilites
de la solution dans dautres cas.
Lobtention des equations bilan repose sur la denition dun volume de controle de
taille importante et nie (Fig. 2.1) ou alors sur un volume elementaire. Ceux-ci peuvent
etre xes ou en mouvement. Dans le premier cas, lapproche Eulerienne, le uide traverse
les volumes, alors que dans la deuxi`eme approche, les volumes suivent le mouvement du
uide et sont donc toujours constitues des memes particules uides (on parlera dapproche
Lagrangienne). Tout se passe donc comme si un observateur suivait le volume et appliquait
un bilan sur ce volume. Le principal avantage de ces approches est de focaliser notre interet
sur une partie de lecoulement, et non pas sur lecoulement dans son ensemble.
Une fois ces volumes denis, lapplication des principes fondamentaux sur un volume
de controle conduit `a une ecriture integrale, qui peut etre manipulee pour atteindre des
relations locales (aux derivees partielles). Ces derni`eres sont obtenues directement par
24

Equations de la m ecanique des fluides
V
Surface S
V
Surface S

V
(a) (b)
Fig. 2.1 Representation dun volume de controle. Dans le cas (a), le volume de controle est xe et est traverse
par le uide en mouvement ; dans le cas (b), le volume de controle est en mouvement, et est alors
constitue des memes particules uides.
lapplication des trois bilans sur le volume elementaire.
La suite du document se propose de presenter les equations mathematiques caracterisant
lapplication des bilans de masse, de quantite de mouvement, et denergie sur ces volumes.
Avant cette presentation, un rappel sur la derivee particulaire est eectue.
2.1 Derivee particulaire
Les particules etant en mouvement, il est parfois utile de connatre la variation dun
param`etre attache `a un volume (de controle ou elementaire) durant son mouvement. Cette
variation est `a la fois liee au point sur lequel elle est consideree, et aussi quen un point
elle peut evoluer. Dans le premier cas, la variation depend du mouvement (dependance
du param`etre par rapport `a

x , t etant xe), alors que le deuxi`eme est linstationnarite
marquee par la dependance du param`etre par rapport au temps,

x etant xe.
La variation enoncee plus haut est traduit mathematiquement par la derivee particu-
laire
1
, notee D/Dt, et est ecrite :
D
Dt
=

t
+

V .

grad (2.1)
Le terme /t est la derivee locale et rappelle le fait que la grandeur peut varier au
cours du temps pour chacun des points de lespace (par exemple, la temperature evolue
dans la journee pour tout lieu en France). Le terme

V .

grad est une derivee convective

representant la variation dun point `a un autre de lespace pour lesquels les proprietes
de lecoulement sont dierentes (par exemple, `a un instant donne, les temperatures sont
dierentes pour tout lieu en France). Ainsi la variation de la temperature au cours dun
trajet eectuee par un volume de controle peut etre traduite par la relation :
DT
Dt
=
T
t
+u
T
x
+v
T
y
+w
T
z
o` u u, v, et w sont les trois composantes de la vitesse

V exprimees dans un rep`ere cartesien.

Nous notons de suite que si on consid`ere un volume (de controle ou elementaire) xe,
le deuxi`eme terme de la derivee particulaire (2.1) est nul. Ainsi seule la derivee locale sera
1
Une demonstration de cet operateur est proposee dans de nombreux ouvrages, dont la reference [11].

Equation de continuite 25
consideree. Si maintenant les volumes sont en mouvement, la derivee particulaire (enti`ere)
sera utilisee.
2.2

Equation de continuite
Nous allons appliquer dans cette partie le fait que la masse attachee `a un volume
detude soit conservee. Lequation mathematique resultante est lequation de continuite.
Conservation de la masse sur un volume de controle xe
On consid`ere un volume de controle D xe (Fig. 2.1.a) de volume V, borne par une surface
S de normale exterieure unitaire

n en chacun des points de sa surface. Dans ce cas, le
bilan de debit de masse `a travers la surface S
_
S

V .

n dS
est egale `a levolution temporelle de la masse du volume au cours du temps. Cette masse
est donnee par lexpression
_
D
dV
La conservation de la masse pour un volume de controle xe est nalement exprimee par
lequation de continuite sous forme integrale :

t
_
D
dV +
_
S

V .

n dS = 0 (2.2)
Conservation de la masse sur un volume de controle mobile
Considerons cette fois ce meme volume en mouvement avec le uide et pourvu cette
fois dans lecriture dun indice t pour se rappeler quil peut changer de forme au cours
du temps (Fig. 2.1.b). Il est donc constitue des memes particules uides. La masse est
identique `a celle vue precedemment. Le volume etant constitue de volumes elementaires
dotes chacun dune masse elementaire dV , il ny a aucune raison pour que la masse totale
ne change au cours du mouvement meme si le volume peut voir sa forme evoluer dans son
mouvement. Finalement, la conservation de la masse secrit alors dapr`es la denition de
la derivee particulaire :
D
Dt
_
D
t
dV = 0 (2.3)
Il sagit l`a de lequation de continuite integrale faisant apparatre la derivee particulaire
introduite par le fait de considerer le volume de controle en mouvement, et sera denommee
forme non-conservative. Une equation sera dte sous une forme conservative lorsquelle
pourra secrire sous la forme :

t
+ div

F = 0
Cela ne semble pas evident `a premi`ere vue, mais lequation (2.2) est de forme conservative,
lintegrale surfacique pouvant etre transformee en integrale volumique avec lintroduction
26

Equations de la m ecanique des fluides
dun divergent du produit

V . Une discussion sur linteret de la forme de lequation

est menee page 32.
Conservation de la masse sur un volume elementaire xe
On consid`ere ici un volume elementaire dV de dimension dx dy dz (Fig. 2.2) et
xe dans lespace. Une quantite de masse passe `a travers les dierentes surfaces de cet
element. Ainsi la masse passant par la surface de gauche vaut (u)dydz, et celle passant
par la surface de droite secrit
_
(u) +
(u)
x
dx
_
dydz
Le bilan de masse sur ces deux surfaces secrit alors :
_
(u) +
(u)
x
dx
_
dydz (u)dydz =
(u)
x
dxdydz
On peut ensuite eectuer un bilan sur les six surfaces correspondant alors `a levolution
temporelle de la masse dans le volume elementaire. Il vient alors
_
(u)
x
+
(v)
y
+
(w)
z
_
dxdydz =

t
dxdydz
On peut nalement ecrire la forme dierentielle de lequation de continuite :

t
+
_
(u)
x
+
(v)
y
+
(w)
z
_
= 0 (2.4)
que lon peut aussi ecrire

t
+ div
_

V
_
= 0 (2.5)
ou en utilisant lecriture indicee :

t
+

x
j
(u
j
) = 0 (2.6)
Ces equations ont ete obtenues en considerant un volume elementaire xe dans lespace
et traverse par le uide. Cette forme dierentielle est ecrite sous une forme conservative
identiquement `a lexpression (2.2) et `a lopposee de lexpression (2.3).
Conservation de la masse sur un volume elementaire mobile
On consid`ere dans ce cas le volume elementaire de la partie precedente anime dun mou-
vement identique `a lecoulement. Sa masse, notee m, correspond au produit dV .
Lorsquon suit ce volume elementaire dans son mouvement, la variation de la masse nulle :
D(m)
Dt
=
D( dV )
Dt
= 0
que lon peut aussi ecrire
D
Dt
+
_
1
dV
D(dV )
Dt
_
= 0
Bilan de quantite de mouvement 27
dx
dy
dz
x
y
z
Fig. 2.2 Volume elementaire xe.
On note que le deuxi`eme terme correspond `a la variation au cours du temps du volume
de lelement dV considere (par unite de volume). Cette quantite est directement liee au
champ de vitesse, et correspond `a la divergence du champ de vitesse. On obtient nalement
la relation
D
t
+ div

V = 0 (2.7)
ou en utilisant lecriture indicee
D
t
+
u
j
x
j
= 0 (2.8)
Ces relations sont les formes locales (non conservatives) de lequation de continuite at-
teinte en considerant un volume elementaire en mouvement.
Discussion sur les quatre equations obtenues
Lapplication de la conservation de la masse sur des volumes de controle ou elementaire,
xe ou en mouvement, conduit `a lecriture de quatre formes dierentes de lequation de
continuite : deux formes integrales (Eqs. 2.2 et 2.3) et deux formes locales (Eqs. 2.4 et 2.7).
Neanmoins, il sagit l`a de quatre ecritures de la meme relation. Chacune de ces equations
peut etre manipulee an datteindre une autre forme. Le membre de gauche de la relation
locale (2.7) peut etre modiee en utilisant la derivee particulaire :
D
t
+ div

V =

t
+

V .

grad + div

V =

t
+ div
_

V
_
correspondant ainsi au membre de gauche de la relation (2.5).
2.3 Bilan de quantite de mouvement
Nous appliquons dans cette partie le Principe Fondamental de la Dynamique enonce
par Newton et applique sur des volumes (de controle ou elementaire), animes dun mou-
vement aligne sur lecoulement ou xes dans lespace.
`
A lidentique de lequation de conti-
nuite, quatre formes dequation du bilan de quantite de mouvement seront atteintes.
Le PFD stipule le fait que le produit de la masse dun element par son acceleration
est equilibree par lensemble des eorts exterieurs agissant sur cet element. Il existe deux
28

Equations de la m ecanique des fluides
types deorts : (1) les eorts volumiques agissant `a distance comme leort de gravite
et (2) les eorts surfaciques agissant directement sur la surface exterieure `a lelement
considere.
Nous nous placerons dans le cas du volume de controle en mouvement (Fig. 2.3) pour
aboutir au bilan de quantite de mouvement.
V
Surface S

V

n

n
Fig. 2.3 Volume de controle en mouvement sur lequel
sappliquent des eorts volumiques et surfaciques.
Les eorts volumiques par unite de masse sont notes

f , tandis que les eorts surfa-
ciques sont representes par le tenseur des contraintes, note , et evalues en chaque point
de la surface par le produit

n . Le tenseur des contraintes est compose dune contribution

liee `a la pression et dune seconde liee `a la viscosite du uide. Il secrit lorsque lon fait
lhypoth`ese dun uide newtonien :

ij
= p
ij
+
ij
= p
ij
+
_
u
i
x
j
+
u
j
x
i
_

2
3

u
k
x
k

ij
(2.9)
Le Principe Fondamental de la Dynamique applique au volume de controle en mou-
vement secrit alors
D
Dt
_
D
t

VdV =
_
D
t

fdV +
_
S
t

n dS (2.10)
Il sagit l`a du bilan de quantite de mouvement sous forme integrale (et non conservative).
Autres ecritures
Lequation (2.10) a ete obtenue en considerant un volume de controle en mouvement. Il
est possible dobtenir dautres formes de cette equation en considerant cette fois le volume
comme etant xe dans lespace, ou en eectuant le bilan de quantite de mouvement sur
un volume elementaire.
Ainsi lapplication du Principe Fondamental de la Dynamique sur un volume de
controle xe entrane la relation :

t
_
D

VdV +
_
S

V
_

V .

n
_
dS =
_
D

fdV +
_
S

n dS (2.11)
On note labsence de la derivee particulaire, le volume etant xe dans lespace. Lapplica-
tion du PFD sur des volumes elementaires (xe ou en mouvement avec le uide) conduit
naturellement `a des relations locales.
Si on consid`ere un volume elementaire xe, le bilan de quantite de mouvement devient :

_
u
i
t
+u
j
u
i
x
j
_
= f
i
+

ij
x
j
(2.12)

Equation de lenergie 29
Cette ecriture est la forme locale (conservative) du bilan de quantite de mouvement.
Lapplication du PFD sur un volume elementaire mobile conduit `a la forme locale (non
conservative) suivante

Du
i
Dt
= f
i
+

ij
x
j
(2.13)
ou, en faisant intervenir la pression et le tenseur des contraintes visqueuses :

Du
i
Dt
= f
i

p
x
i
+

ij
x
j
(2.14)
Encore une fois, il est important de comprendre que les expressions (2.10), (2.11), (2.12),
et (2.13) ou (2.14) sont des representations dierentes dune unique equation caracterisant
le bilan de quantite de mouvement.
2.4

Equation de lenergie
Nous appliquons dans cette partie le premier principe de la Thermodynamique stipu-
lant que la variation de lenergie associe `a un volume (de controle ou elementaire) est egale
`a lenergie re cue par celui-ci. Cette derni`ere est composee dune partie correspondant au
ux de chaleur `a travers la surface du volume, et dune seconde composante liee au travail
des eorts surfaciques et volumiques.
Considerons un volume elementaire en mouvement (Fig. 2.2) auquel on associe une
energie decomposee en une energie interne par unite de masse, notee e, liee au mouvement
aleatoire des molecules composant le volume, et une energie cinetique par unite de masse,
V
2
/2, consecutive au mouvement du volume lui-meme. Ainsi, levolution au cours du
temps de lenergie du volume en suivant son mouvement vaut

D
Dt
_
e +
V
2
2
_
dxdydz (2.15)
Essayons devaluer `a present le travail des eorts exterieurs, estime par le produit
scalaire de leort par la vitesse. Ainsi le travail des eorts volumiques agissant sur le
volume elementaire anime dune vitesse

V secrit

f .

V dV (2.16)
avec dV = dxdydz le volume de lelement considere. Les eorts surfaciques, quant `a
eux, impliquent les contraintes de pression et les contraintes visqueuses sappliquant sur
les six surfaces. Nous traiterons uniquement des surfaces ayant une normale selon x, et
generaliserons ensuite les formulations obtenues. Ainsi le bilan du travail des eorts de
pression sur les surfaces gauche et droite secrit
_
up
_
up +
(up)
x
dx
__
dydz =
(up)
x
dxdydz
De fa con identique, le bilan du travail des eorts visqueux sur les surfaces inferieure et
superieure secrit
__
u
yx
+
(u
yx
)
y
dy
_
u
yx
_
dxdz =
(u
yx
)
y
dxdydz
30

Equations de la m ecanique des fluides
`
A cette relation, il faut ajouter la contribution des composantes toujours selon x calculees
sur les autres paires de surfaces elementaires. Ainsi, le bilan sur les surfaces ayant une
normale selon x du travail des eorts de pression et de viscosite peut se lire :
_

(up)
x
+
(u
xx
)
x
+
(v
yx
)
y
+
(w
zx
)
z
_
dxdydz
Finalement, en tenant compte des contributions de toutes les composantes cartesiennes,
le travail des eorts volumiques et surfaciques secrit
_

(up)
x

(vp)
y

(wp)
z
+
(u
xx
)
x
+
(u
yx
)
y
+
(u
zx
)
z
+
(v
xy
)
x
+
(v
yy
)
y
+
(v
zy
)
z
+
(w
xz
)
x
+
(w
yz
)
y
+
(w
zz
)
z
_
dxdydz
Interessons-nous maintenant au ux de chaleur, celui-ci provenant `a la fois dune absorp-
tion ou dune emission de rayonnement, mais aussi dun ux de chaleur `a travers sa surface
consecutivement `a un gradient de temperature. Soit q la densite volumique du taux de
chaleur. Le taux de chaleur transfere par rayonnement secrit alors
qdxdydz
Notons de la meme fa con

q le vecteur densite surfacique de chaleur re cue par le volume
par conduction. Comme nous lavons fait precedemment, on sinteresse au bilan de chaleur
transferee par conduction par les surfaces ayant une normale selon x, soit
_
q
x

_
q
x
+
(q
x
)
x
dx
__
dydz =
(q
x
)
x
dxdydz
Le bilan sur lensemble des surfaces, en considerant les eets du rayonnement conduit `a
lexpression
_
q
_
q
x
x
+
q
y
y
+
q
z
z
__
dxdydz
Le ux de chaleur par conduction peut etre modeliser par la loi de Fourier selon laquelle
le ux de chaleur est proportionnel au gradient de temperature :

q = k

gradT
Le bilan precedent secrit donc
_
q +

x
_
k
T
x
_
+

y
_
k
T
y
_
+

z
_
k
T
z
__
dxdydz (2.17)
La forme nale de la conservation de lenergie appliquee sur un volume elementaire
en mouvement secrit `a partir des relations mises en evidence precedemment (Eqs. 2.15,

Equation de lenergie 31
2.16, et 2.17). Il vient

D
Dt
_
e +
V
2
2
_
= q +

x
_
k
T
x
_
+

y
_
k
T
y
_
+

z
_
k
T
z
_

(up)
x

(vp)
y

(wp)
z
+
(u
xx
)
x
+
(u
yx
)
y
+
(u
zx
)
z
+
(v
xy
)
x
+
(v
yy
)
y
+
(v
zy
)
z
+
(w
xz
)
x
+
(w
yz
)
y
+
(w
zz
)
z
+

f .

V (2.18)
La relation (2.18) est la forme non conservative de la conservation de lenergie obtenue
en considerant un volume elementaire en mouvement. On peut evidemment atteindre
dautres formes de cette equation en basant notre reexion sur un volume element xe,
ou en considerant un volume de controle (mobile ou xe).
On pourrait dautre part sinteresser uniquement `a la variation de lenergie interne du
volume au cours du temps. Une demarche possible consiste `a ecrire une equation similaire
`a lequation (2.18) mais adaptee `a lenergie cinetique. Cette derni`ere est facilement atteinte
en multipliant la relation (2.14) par u
i
faisant apparatre alors lenergie cinetique. Cette
relation est ensuite soustraite de lequation (2.18). Il vient alors :

De
Dt
= q +

x
_
k
T
x
_
+

y
_
k
T
y
_
+

z
_
k
T
z
_
p
_
u
x
+
v
y
+
w
z
_
+
xx
u
x
+
yx
u
y
+
zx
u
z
+
xy
v
x
+
yy
v
y
+
zy
v
y
+
xz
w
x
+
yz
w
y
+
zz
w
z
(2.19)
`
A ce stade, le nombre de variables inconnues est de sept (u, v, w, p, , e, et T) pour
cinq equations. An de fermer le probl`eme, on fait souvent lhypoth`ese en aerodynamique
dun gaz parfait traduit par une sixi`eme equation, la relation des gaz parfait :
p = rT
Une septi`eme equation est necessaire, et provient generalement de la thermodynamique
(toutes les equations de laerodynamique ayant ete trouvees). Il sagit l`a de la relation
entre lenergie interne et la temperature. Pour un gaz parfait, il vient :
e = c
v
T
o` u c
v
est la chaleur specique `a volume constant.
32

Equations de la m ecanique des fluides
2.5 Bilan sur les equations obtenues par lapplication
des principes fondamentaux
Les paragraphes precedentes avaient pour but decrire des mod`eles des trois principes
fondamentaux (i.e. conservations de la masse, de la quantite de mouvement, et de lenergie)
appliques sur des volumes detude. Nous avons note que lorsque le volume etait anime
dune vitesse identique `a celle de lecoulement, les equations atteintes etaient de forme
non-conservative. Ce qui nest plus le cas si on consid`ere maintenant un volume xe
dans lespace. Comprenons tout de meme que la distinction forme conservative/forme
non conservative ne porte que tr`es peu dinteret physique. La subtilite devient important
lorsquon aborde la CFD, cest `a dire les techniques numeriques permettant de predire
un ecoulement. En eet, la technique des Volumes Finis, sur laquelle la plupart des codes
commerciaux reposent, est tout particuli`erement adaptee ` a la mecanique des uides car
cette technique est basee sur des bilans de ux sur les surfaces des volumes.
Dautre part, les equations precedentes ont ete etablies sans poser dhypoth`eses im-
portantes sur lecoulement. Ainsi ces equations traitent aussi bien des ecoulements com-
pressibles `a tr`es hautes vitesses que des ecoulements incompressibles visqueux. Neanmoins
dans le deuxi`eme cas, il est preferable de simplier les equations de bilan pour les purger
de termes inutiles car negligeables. Nous ne traiterons ici que des ecritures locales des
equations.
Les equations (2.13) et (2.14) sont les equations de Navier-Stokes, telles quils les
ont trouvees independamment lun de lautre. Cependant, lorsquon aborde la CFD, les
equations de Navier-Stokes representent generalement la serie dequations issues des trois
principes fondamentaux.
2.6 Conditions limites
Les equations de bilan ecrites plus haut decrivent `a la fois des ecoulements superso-
niques, comme des ecoulements tr`es visqueux. Et pourtant, il sagit de la meme equation.
Le type decoulement est alors genere par le type de conditions limites qui sont tr`es
dierentes dans les deux cas precedents. Celles-ci sont guidees par la physique du probl`eme.
Lorsque le uide est reel, les particules uides adh`erent `a toute paroi. Cette condition
de non-glissement est traduit par lexpression :
u = v = w = 0 (en paroi)
Si le uide est denue de viscosite, et devient alors parfait, la condition devient moins
sev`ere :

V .

n = 0 (en paroi)
et permet un glissement des particules uides sur la paroi, seule la composante normale
`a la surface de la vitesse devenant nulle. De meme si cette surface est maintenue `a une
temperature constante, T
w
, on ajoute alors la condition
T = T
w
(en paroi)
Il se peut neanmoins que dans certaines congurations, la temperature parietale reste
inconnue et surtout liee au materiau constituant la surface. La condition limite adaptee
Conditions limites 33
est alors donnee par le ux de chaleur entre la paroi et le uide. La condition limite est
alors ecrite sous forme de gradient `a la paroi :
_
T
n
_
w
=
q
w
k
(en paroi)
o` u n est la direction normale `a la surface, et k la conductivite thermique. Ainsi des
conditions de parois adiabatiques seront caracterisees par un gradient de temperature nul
en paroi.
Ce chapitre nous a permis decrire trois equations mathematiques traduisant le res-
pect de trois principes fondamentaux de la mecanique des uides et de la thermique.
La prediction numerique des ecoulements par les codes de calculs est basee sur ces trois
equations dont lecriture est adaptee `a la discretisation de la methode. Le prochain para-
graphe presente tr`es bri`evement les trois methodes principales de la resolution numerique
des ecoulements, et sattardent tout particuli`erement sur les methodes statistiques.
34

Equations de la m ecanique des fluides
Chapitre 3
MOD

ELISATION DE LA
TURBULENCE
Le chapitre precedent a presente les equations de mouvement dun uide et de trans-
ferts thermiques possibles entre un uide et une paroi. Ces equations sont instantanees
et valables pour tous les types decoulement. Ce chapitre a pour but de presenter le prin-
cipe des trois approches numeriques DNS, LES, et RANS. Nous nous attarderons tout
particuli`erement sur la methode statistique car tr`es employee par les codes du commerce.
Pour des raisons de simplicite decriture, seuls les ecoulements incompressibles seront
developpes ici.
3.1 Trois approches numeriques
Malgre des eorts importants de recherche depuis plus dun si`ecle, la modelisation
des ecoulements turbulents demeure un de `a relever. Il existe principalement trois axes
de recherche, (1) les resolutions numeriques deterministes (DNS, pour Direct Numerical
Simulation), (2) les methodes semi-deterministes (LES, pour Large Eddy Simulation), et
enn (3) les methodes statistiques (RANS pour Reynolds Average Navier-stokes
1
) plus
anciennes et donc largement developpees.
Ces trois methodes ont des objectifs et necessitent des co uts de calcul dierents.
1
Le but de cette methode est en eet decrire une equation de Navier-Stokes moyennee `a partir dune
decomposition de Reynolds (3.1.3)
36 Mod elisation de la turbulence
3.1.1 Simulation directe numerique (DNS)
Cette methode permet de resoudre directement les equations de Navier-Stokes sans
aucune modelisation. Elle presente ainsi lavantage de donner acc`es `a toutes les quantites
instantanees considerees dans lecoulement. Une application possible est donnee sur la
gure 3.1 representant la turbulence creee par un jet dair issu dune tuy`ere.
Tous les mouvements doivent etre resolus par cette methode, la taille de maille doit
donc etre inferieure `a lechelle de dissipation. Le nombre de mailles est alors important.
Ceci a pour consequence des temps de calcul extremement longs, et dautant plus longs
que la vitesse de lecoulement est elevee
2
. La capacite et la performance des calculateurs
actuels ne cessent de progresser mais ne permettent pas encore de sonder des ecoulements
complexes et `a hautes vitesses `a partir de cette methode. Neanmoins, celle-ci permet de
mieux comprendre les comportements turbulents dans des congurations simples, mais
aussi dans certains cas de valider les mod`eles de turbulence issus de la modelisation
statistique (paragraphe 3.1.3 page 38).
Une alternative `a cette methode est de simuler uniquement les grandes echelles et de
modeliser
3
les petites : LES.
Fig. 3.1 Simulation directe numerique dun jet,
M=1.92 (dapr`es Freund [7]).
3.1.2 Simulation des grandes echelles (LES)
Il sagit dune methode numerique intermediaire entre la DNS et les methodes statis-
tiques consistant `a appliquer un ltre spatial en tout point du domaine. Le champ ltre
est obtenu par un produit de convolution dans lespace :
(x) =
_
D
(x

)G(x, x

)dx

2
Pour exemple, la simulation dun jet dair issu dune tuy`ere de rayon r = 1 cm dans un domaine
detendu 0 < x/r < 40, 15 < y/r < 15, 15 < z/r < 15 necessite une grille de 12 millions de points et
un temps de calcul de 100 h sur un super calculateur.
3
On remarquera ici la distinction entre simuler un ecoulement, cest `a dire calculer numeriquement
`a partir des diverses equations, et modeliser, cest `a dire utiliser une loi qui rend compte du phenom`ene
physique.
Trois approches numeriques 37
Fig. 3.2 Simulation des grandes echelles : colli-
sion axiale de deux anneaux tourbillonnaires (dapr`es
Mansel [21]).
o` u D est le domaine de calcul, et G le ltre determinant la taille des tourbillons `a simuler.
Le ltre separe donc les grandes echelles (simulees) des petites structures (modelisees).
On suppose ici que le comportement de ces derni`eres ne depend pas de la geometrie
et est donc isotrope, ce qui nest pas le cas des grandes echelles qui, elles, voient leur
comportement guide par leur environnement. La simulation des grandes echelles reste
neanmoins tr`es delicate car elles ont un comportement anisotrope, elles sont sujettes `a des
eets historiques, et sont fortement dependante du type decoulement et de ses conditions
aux limites. La taille de maille est choisie largement superieure `a lechelle de Kolmogorov
et correspond la plupart du temps `a la taille du ltre. Le ltage precedent secrit alors :
(x) =
1
V
_
V
(x

) dx

, x

V
avec V le volume de la maille.
Malgre la modelisation des petites structures, le temps de calcul demeure important et
les calculs sont limites, comme pour la DNS `a des nombres de Mach faibles. Neanmoins,
alors que la puissance de nos ordinateurs personnels etaient insusantes pour ces types
de calcul, depuis les annees 2000 la tendance sest inversee (Fig. 3.3).
Fig. 3.3 Illustration des
puissances disponibles et
necessaires pour un cal-
cul numerique base sur
la methode LES (Dapr`es
Pope [27]).
1980 2000 2020
besoins
disponibilite
p
u
i
s
s
a
n
c
e
d
e
s
o
r
d
i
n
a
t
e
u
r
s
(
L
o
g
)
La gure 3.2 propose un exemple de calcul possible : la collision de deux anneaux
tourbillonnaires en vue dune meilleure comprehension des interactions entre les structures
turbulentes. La gure 3.4 montre une autre forme dinteraction, ici entre deux tourbillons,
resultant sur lappariement des deux structures. Les etudes fondamentales ont aussi pour
38 Mod elisation de la turbulence
Fig. 3.4 Observation de lappariement de deux tourbillons.
objectif de simuler et predire ces comportements.
La DNS et la LES sont tr`es prometteuses car la perte dinformations est minimale
pour la premi`ere (seules les contributions moleculaires sont ltrees), un peu plus impor-
tante pour la deuxi`eme. Cependant, le co ut dun calcul est inversement proportionnel `a
la precision obtenue. Pour cette raison, dans des congurations industrielles, ces deux
methodes sont inutilisables, et on leur pref`erera les methodes statistiques exposees dans
la suite du document.
3.1.3 Modelisation statistique de la turbulence (RANS)
La strategie adoptee ici consiste `a mettre de cote le mouvement instantane du uide,
dans le but dexprimer les equations du champ moyen. Le soucis de simuler toutes les
petites structures de lecoulement est donc elimine. La taille de maille, ainsi que la valeur
des pas de temps pour des etudes instationnaires, deviennent plus importantes. Cette ap-
proche est donc moins co uteuse en temps de calcul que la DNS et la LES, ceci expliquant
sa grande utilisation dans le monde industriel. Neanmoins, notons tout de suite que cette
approche presente un fort degre dempirisme, rendant alors la methode peu able dans
certaines congurations. Le domaine dapplication est tr`es vaste, avec transfert thermique
(Fig. 3.5) ou sans echange thermique (Fig. 3.6).
Le principe de cette methode repose sur la resolution des equations de Navier-Stokes
moyennees. Celles-ci sont obtenues en introduisant une decomposition des variables du
probl`eme `a traiter :
(x
i
, t) = (x
i
) +

(x
i
, t)
o` u (x
i
) est la moyenne du param`etre (x
i
, t) sur une duree dobservation T :
(x
i
) = lim
T
1
T
_
T
0
(x
i
, t)dt
avec = {u
i
, p, T}. La duree dobservation doit etre importante comparativement `a
lechelle de temps de la turbulence, theoriquement elle doit tendre vers linni. Cependant
si elle est assez grande, la valeur moyenne nen depend plus. Cette decomposition appelee
decomposition de Reynolds (1894), a ete introduite par Boussinesq
4
(1872).
4
Valentin Joseph Boussinesq(1842-1929), physicien et mathematicien fran cais. Boussinesq apporta une
Trois approches numeriques 39
Fig. 3.5 Temperature dans un lm dhuile cisaille par deux disques en rotation symbolisant un embrayage en
position ouverte (dapr`es Changenet [4]).
Fig. 3.6

Ecoulement dair au passage dune roue dentee en mouvement - D
p
= 150 mm, m = 5 mm, largeur
= 24 mm (dapr`es Marchesse [23]).
On suppose dans cette decomposition que les valeurs uctuantes sont centrees, cest
`a dire que leurs valeurs moyennes sont nulles. Les proprietes les plus utilisees sont les
suivantes :
() = +

= +

= et = +

contribution enorme `a la physique mathematique. Et son travail en hydraulique est considerable. Il etudia
les tourbillons, les ondes de surface, la resistance `a lavancement dun obstacle, et les eets refroidissants
dun liquide.
40 Mod elisation de la turbulence
Il existe dautres decompositions, telle que la decomposition triple : = +

+

,
composante moyenne, organisee et aleatoire. Cette decomposition est tr`es utilisee dans la
modelisation de la turbulence atmospherique, cas o` u les tourbillons gardent longtemps une
structure uctuante organisee. Ce formalisme etant plus complexe, on regroupe parfois
la composante moyenne et la composante organisee dans un terme unique. On obtient
alors une decomposition equivalente `a la decomposition de Reynolds. Lorsque les eets
de compressibilite deviennent non-negligeables, la decomposition de Fabre est preferee
`a celle de Reynolds car plus adaptee `a ce type decoulement. On denit dans ce cas la
moyenne massique,

= /, permettant de mieux tenir compte des variations de masse
volumique.
Notons au passage que lapplication de la moyenne de Reynolds conduit `a une gigan-
tesque perte dinformation. Il devient alors illusoire dessayer dobtenir par cette approche
un mod`ele capable de reproduire correctement tous les phenom`enes, dans toutes les simu-
lations. Cette reexion renvoie `a la notion duniversalite des mod`eles. Les decompositions
des variables de notre probl`eme, u
i
= u
i
+ u

i
et p = p + p

, sont introduites dans les

equations precedentes. On obtient alors
u
i
x
i
= 0 (3.1)
et

_
u
i
t
+u
j
u
i
x
j
+u

j
u

i
x
j
_
= f
i

p
x
i
+

ij
x
j
que lon peut mettre sous la forme suivante en utilisant le fait que u

i
/x
i
= 0 :

_
u
i
t
+u
j
u
i
x
j
_
= f
i

p
x
i
+

x
j
[
ij
u

i
u

j
] (3.2)
Les deux relations (3.1) et (3.2) representent les equations de Navier-Stokes moyen-
nees. On remarque que lequation de bilan moyennee (3.2) semble identique `a celle ecrite
pour lecoulement instantane, `a lexception cependant dun terme ajoute, u

i
u

j
, ayant
la dimension dune contrainte. Ce terme porte le nom de tenseur de Reynolds. Ainsi le
champ moyen ne satisfait pas aux equations de Navier-Stokes, telles que nous les avions
ecrites dans le precedent chapitre, le tenseur de Reynolds faisant apparatre le lien etroit
entre les champs moyen et uctuant.
On note ainsi que la turbulence apporte une contrainte supplementaire `a lecoulement.
Neanmoins, les contributions des contraintes de cisaillement et de celles liees `a la turbu-
lence ne sont pas identiques. En eet, les eets visqueux sont dominants dans la region
proche des parois, alors que pour des regions qui en sont eloignees, la contrainte turbulente
a une contribution plus importante (Fig. 3.7).
Limitant la description au champ moyen [u
i
(x
j
, t), p(x
j
, t)], le syst`eme regissant ces
fonctions se reduit en situation isovolumique aux equations (3.1) et (3.2). Le tenseur de
Reynolds fait donc apparatre six termes supplementaires (u
2
, v
2
, w
2
, u

, u

, et
v

) sajoutant aux variables habituelles (u, v, w, et p). Il y a donc 10 inconnues pour 4

equations, et il est alors necessaire de trouver une strategie nous permettant de fermer
ce syst`eme.
Mod`ele de turbulence du premier ordre - Concept de la viscosite turbulente 41
Fig. 3.7 Repartition des contraintes de cisaillement
et turbulentes (normalisees par la contrainte totale,
(y)) dans le cas dun ecoulement dans un canal.
Lignes continues, Re = 13 750 ; lignes discontinues,
Re = 5 600 (dapr`es Moser [26]).
3.1.4 Remarques sur les trois approches
Nous venons de decrire de fa con succinte les trois grandes methodes numeriques ac-
tuelles. Leurs domaines de modelisation (quand il existe) et de simulation sont rassembles
sur la gure 3.8. Dapr`es ce que nous venons de voir, plus le domaine de modelisation est
faible, plus le calcul sera `a meme de donner une description precise de lecoulement. Ce-
pendant, celui-ci necessite alors une taille de maille tr`es ne, entranant irremediablement
des temps de calculs co uteux. Pour cela, la modelisation statistique semble seduisante
car plus rapide. Neanmoins, lutilisateur sembarque d`es lors dans le choix de mod`eles de
turbulence
5
.
Fig. 3.8 Representation des domaines de
simulation et de modelisation pour la DNS,
la LES et la modelisation statistique.
3.2 Mod`ele de turbulence du premier ordre - Concept
de la viscosite turbulente
Lapproche de la modelisation statistique de la turbulence necessite la fermeture du
syst`eme dequations regissant le mouvement moyen. Generalement, les methodes de fer-
meture sont classees en fonction du nombre dequations supplementaires `a resoudre. La
5
Il utilisera dans la plupart des cas le mod`ele k qui, il est vrai predit des ecoulements corrects mais
aussi pour des raisons dignorance dautres mod`eles.
42 Mod elisation de la turbulence
strategie adoptee par les mod`eles du premier ordre et presentes dans la suite du docu-
ment repose sur le concept de viscosite turbulente presente par Boussinesq en 1877
[2]. Son idee est basee sur lobservation tendant `a montrer que le transfert de quantite
de mouvement dans un ecoulement turbulent est fortement domine par le melange des
grosses structures. Cette viscosite est notee
T
et relie lineairement le tenseur de Reynolds
`a lecoulement moyen :
u

i
u

j
=
T
_
u
i
x
j
+
u
j
x
i
_
+
2
3
k
ij
. .
(a)
(3.3)
o` u k est lenergie cinetique moyenne du champ turbulent par unite de masse, appelee de
fa con plus concise lenergie cinetique turbulente :
k =
1
2
u

k
u

k
=
1
2
(u
2
+v
2
+w
2
)
et
ij
le symbole de Kronecker (
ij
= 1 si i = j et
ij
= 0 si i = j). Le terme isotrope
(a) est necessaire, an de ne pas avoir

k = u

i
u

i
/2 0 par contraction des indices
6
. La
viscosite est a priori une fonction locale de lecoulement
T

T
(

x , t) et donc fontion
du mouvement turbulent, contrairemement `a la viscosite moleculaire qui est une propriete
du uide.
La relation (3.3) repose sur des hypoth`eses simplicatrices de lecoulement et sa turbu-
lence : (1) linstantaneite de la reponse de la turbulence ` a une variation du champ moyen,
cest `a dire la non prise en compte de lhistoire de la deformation et de la turbulence ;
(2) la localite, la turbulence etant inuencee que par son voisinage immediat, (3) la faible
inhomogeneite ; (4) la linearite de cette loi de comportement entrainant une surestimation
de la production de la turbulence. La trop forte diusivite turbulente est ainsi `a lorigine
de la dissipation anticipee des structures turbulentes. Ce dernier point est probablement
celui qui fait le plus defaut `a cette approche. En eet, la relation (3.3) est bien adaptee
`a la reproduction des composantes de cisaillement (i.e. u

et v

) qui sont produites

respectivement par le cisaillement d u au sillage et `a la couche limite. En revanche, elle
est incapable de distinguer les trois composantes diagonales, predisant tout simplement
u
2
= v
2
= w
2
= 2k/3. On comprendra mieux que les mod`eles du premier ordre ne
puissent marcher dans toutes les situations. Ceci est grandement ameliore par lintroduc-
tion dune relation non-lineaire de la loi de comportement du materiau turbulent par
lintermediaire de termes quadratiques voire cubiques des gradients de vitesses moyennes.
La viscosite turbulente, inconnue de prime abord, doit donc etre denie. Par analogie
avec la theorie cinetique des gaz, la determination de la viscosite necessite deux echelles,
u l, o` u u est la vitesse quadratique des molecules constituant le gaz et l le libre
parcours moyen de ces particules (i.e. la distance moyenne eectuee par ces particules
entre deux chocs). On pourra alors ecrire dans notre cas :

T
= u
T
l (3.4)
o` u u
T
et l representent respectivement une vitesse et une longueur caracteristique de
la turbulence locale. Ainsi, le propos des mod`eles de turbulence est destimer ces deux
6
En eet, sans ce terme, u

i
u

i
= 2k = 2
T
u
i
x
i
, ce dernier terme etant nul dapr`es la relation 3.1.
Mod`ele de turbulence du premier ordre - Concept de la viscosite turbulente 43
echelles de la turbulence an dapprocher au mieux la valeur de
T
. Le mod`ele de turbu-
lence sera dautant plus complexe quil cherchera `a se rapprocher au plus pr`es de la realite,
`a savoir tenir compte des eets de la convection, de la production et de la dissipation de
la turbulence le long de lecoulement.
3.2.1 Mod`ele algebrique ou mod`ele `a zero equation
Prandtl
7
[28] propose en 1925 une formulation de la viscosite turbulente utilisant le
concept de longueur de melange, l
m
, telle que

T
= l
2
m

u
y

(3.5)
Cette longueur, inspiree du libre parcours moyen dans la theorie cinetique des gaz repose
sur le principe que le transport de la quantite de mouvement seectue sur une distance
egale au libre parcours moyen. Lexpression de cette longueur depend de la conguration
etudiee. Le tableau 3.1 rassemble quelques valeurs de longueurs de melange dans des
congurations frequemment etudiees.
Tableau 3.1 Longueurs de melange dans le cas decoulements bidimensionnels (dapr`es Rodi [31]).

Ecoulement Longueur de melange, l

m
L
Couche de melange 0, 07L epaisseur de la
couche de melange
Jet 0, 09L demi-epaisseur du jet
Sillage 0, 16L demi-epaisseur du sillage
Jet axisymetrique 0, 075L demi-epaisseur du jet
Couche limite :
Sous-couche visqueuse [1 exp(y
+
/26)] epaisseur
et region logarithmique de la
Partie superieure (y/L 0, 22) 0, 09L couche-limite
Ce mod`ele a un degre de generalite tr`es faible puisque dependant du type decoulement
envisage. Cependant lorsquil est applique dans des congurations decoulement adaptees,
les predictions numeriques sont proches de celles mesurees (Fig. 3.9). Ce mod`ele est dt
incomplet car pour decrire la longueur de melange en fonction de la position celle-ci
depend de lecoulement. Dautre part, cette vision de la turbulence est trop simpliste, car
la longueur doit etre dierente que lon se positionne dans lecoulement ou proche des
parois
8
. Il semble en pratique que le concept de viscosite turbulente sous cette forme est
trop imprecis pour permettre des predictions satisfaisantes.
7
Ludwig Prandtl (1875-1953), aerodynamicien allemand. Il est le p`ere de laerodynamique en mettant
en evidence la presence dune couche limite pr`es des parois. Ceci permit dexpliquer les phenom`enes de
portance et de tranee necessaires `a lamelioration de laviation. Il fut le directeur de th`ese de Paul Richard
Heinrich Blasius, qui utilisa le concept de couche limite dans lestimation de la tranee dune plaque plane.
Il montra aussi que les pertes de charge dans une conduite est fonction du nombre de Reynolds.
8
Des modications ont ete apportees au mod`ele de longueur melange pour tenir compte des eets de
parois. Van Driest, Clauser et Klebano proposent par exemple un mod`ele de longueur de melange reduite
dans la region laminaire de la couche limite par lutilisation dun facteur jouant le r ole damortisseur.
44 Mod elisation de la turbulence
De plus, si on se place au centre dune couche de melange caracterisee par un gradient
de vitesse nul, la viscosite turbulente y devient nulle dapr`es la relation (3.5), ce qui est
refute par les mesures. Il devient alors necessaire davoir recours `a dautres theories pour
mettre en evidence les diverses interactions entre les instabilites de tailles dierentes.
Pour toutes ces raisons le mod`ele de longueur de melange nest pas applicable dans
une modelisation generale de la turbulence. De ce fait, il sera insere dans des mod`eles de
turbulence plus complexes pour traiter par exemple le comportement du uide dans les
regions proches des parois.
(a) (b)
Fig. 3.9 Schema dun sillage (a) et predictions numeriques du champ de vitesse dans le sillage dun obstacle
cylindrique circulaire utilisant le concept de longueur de melange (dapr`es Versteeg [35] et Schlichting
[32]).
3.2.2 Mod`eles de fermeture `a une equation de transport
Les mod`eles de fermeture `a une equation de transport reposent sur lhypoth`ese de
Boussinesq, et permettent doter le probl`eme dune viscosite turbulente nulle au milieu
dune couche de melange. On presente dans la suite deux exemples de ces mod`eles `a 1
equation.
Mod`ele de fermeture base sur lequation de transport de lenergie cinetique
turbulente
Prandtl et Kolmogorov [13] proposent en 1940 une relation dans laquelle la viscosite
turbulente est proportionnelle `a la racine carree de lenergie cinetique tubulente :

T
= C
k
l
_
k (3.6)
Le choix du param`etre k vient naturellement de sa presence dans la loi de comportement
(3.3). Cette approche est ainsi basee sur la connaissance en tout point de lecoulement de
lenergie cinetique turbulente :

_
k
t
+u
j
k
x
j
_
. .
[1]
=

x
j
_

k
x
j
_
. .
[2]

x
j
_

2
u

j
u

i
u

i
+p

j
_
. .
[3]
u

i
u

j
u
i
x
j
. .
[4]

i
x
j
u

i
x
j
. .
[5]
La variation totale de lenergie cinetique turbulente [1] est equilibree par :
Mod`ele de turbulence du premier ordre - Concept de la viscosite turbulente 45
[3] la diusion de lenergie cinetique turbulente par eet des uctuations de vitesse
et de pression.
[4] linteraction avec lecoulement moyen. Cette energie, cedee par lecoulement
moyen, est appelee terme de production.
[5] le taux de dissipation moyen (multiplie ici par la masse volumique), note par la
suite , et toujours positif.
Le terme de diusion [3] doit etre modelise, cest generalement le cas `a laide dun terme
de diusion par gradient :

2
u

j
u

i
u

i
+p

j
_

k
k
x
j
(3.7)
o` u
T
est la viscosite turbulente vue plus haut, et
k
est le nombre de Prandtl turbulent
(proche de lunite). On suppose l`a que les deux termes ont des eets diusifs, et peuvent
donc etre modelises de fa con identique. Ceci serait parfaitement justie si lagitation
turbulente avait un comportement de pure agitation (au sens moleculaire), et `a une echelle
tr`es petite devant celle de lecoulement moyen. Il nen est rien. En eet, la presence de
structures coherentes dans la turbulence va `a lencontre de toute separation dechelles.
De la meme mani`ere, le terme de production [4] doit etre modelise, ceci est realise `a
partir de la denition de la viscosite turbulente (3.3) :
P
k
= u

i
u

j
u
i
x
j

T
_
u
i
x
j
+
u
j
x
i
_
u
i
x
j
Lequation de transport de lenergie cinetique turbulente est ecrite `a partir dune analogie
avec lequation de transport de lenergie cinetique (approche phenomenologique), et non
pas `a laide dune approche theorique. On cherche donc `a equilibrer le transport du taux
de dissipation par un terme de convection, un terme source, un terme puits, de la diusion
visqueuse et turbulente. On peut ainsi ecrire :
(k)
t
+u
j
(k)
x
j
=

x
j
__
+

t

k
_
k
x
j
_
+P
k
(3.8)
Cette approche, plus quosee, est tr`es largement utilisee aujourdhui. Elle est loin
detre parfaite, mais personne na propose autre chose pour le moment. Ce mod`ele dispose
nalement de cinq equations pour les cinq fonctions scalaires u
i
, p et k. Il a lavantage
detre peu complexe `a mettre en place, et prend en compte partiellement lhistoire de la
turbulence. Ce mod`ele, comme le precedent, reste incomplet car necessitant de prescrire
lechelle de longueur en tout point, ce qui reste dicile pour des geometries complexes.
Lequation de transport de lenergie cinetique turbulente reste valable pour des ecou-
lements tr`es fortement turbulents. Cest `a dire dans les regions eloignees des eets des
parois (i.e. y
+
superieur `a 30).
Mod`ele de fermeture de Spalart et Allmaras base sur lequation de transport
de la viscosite v
Un des mod`eles de turbulence `a 1 equation de transport, autre que celle de lenergie
cinetique turbulente, le plus cel`ebre est le mod`ele propose par Spalart et Allmaras [34]. Ce
46 Mod elisation de la turbulence
mod`ele a ete adapte `a des congurations de couches de melange 2D et de couches limites
sur des plaques planes. La viscosite turbulente est donnee par la formulation :

T
= f
v1
v (3.9)
o` u v est obtenue par resolution de son equation de transport :

t
( ) +

x
j
( u
j
) = P

+
1

_

x
j
_
( + )

x
j
_
+C
b2

_

x
j
_
2
_
D

(3.10)
o` u P

represente la production de viscosite turbulente et D

sa destruction dans les regions

proches des parois liee `a lamortissement visqueux.

et C
b2
sont des constantes, et est
la viscosite dynamique.
Le mod`ele de turbulence propose par Spalart et Allmaras est un mod`ele dt bas-
Reynolds, cest `a dire construit de telle sorte quil puisse etre utilise dans les regions
de paroi. Ceci est possible `a partir du comportement du coecientf
v1
dans lexpression
(3.9). Aucune loi de paroi nest alors necessaire. Certains codes de calculs commerciaux
proposent neanmoins une modication de ce mod`ele adaptee `a lutilisation dune loi de pa-
roi diminuant alors le nombre de mailles dans la region de paroi. Ce mod`ele est performant
pour des applications aeronautiques car les ecoulements attaches et decolles sont tr`es bien
representes. Cependant, ce mod`ele est bien trop simple dans sa construction pour donner
de bons resultats dans une large gamme decoulements. En eet, son manque duniversa-
lite lui fait defaut. Ce point est ameliore en considerant des mod`eles `a deux equations de
transport.
3.2.3 Mod`eles de fermeture `a deux equations de transport
An de lever ce probl`eme dempirisme concernant lechelle de longueur, il est vite ap-
paru naturel de calculer cette grandeur en resolvant une equation de transport supplemen-
taire. On acc`ede ainsi aux methodes du premier ordre `a deux equations. Nous lavons
vu precedemment, il est tr`es naturel dintroduire lenergie cinetique turbulente dans la
prediction de la viscosite turbulente. Le choix du second param`etre, `a partir duquel une
equation de transport sera ecrite, est plus delicat. Plusieurs suggestions ont ete proposees
basees sur les variables
kl 1/
o` u et sont respectivement le taux de dissipation de lenergie cinetique turbulente et
le taux de dissipation specique ( = /k) representant le taux de dissipation par unite
denergie cinetique turbulente. Dans tous ces cas, la viscosite turbulente est formulee `a
partir de lenergie cinetique turbulente et du second param`etre.
Mod`ele de fermeture k
Le mod`ele de fermeture k est lun des plus utilise dans le monde industriel. La viscosite
turbulente est donnee par lexpression :

T
= C

k
2

(3.11)
Mod`ele de turbulence du premier ordre - Concept de la viscosite turbulente 47
et est donc ecrite par combinaison de lenergie cinetique turbulente et du taux de dissi-
pation
=
u

i
x
j
u

i
x
j
(3.12)
estimes tous deux par leur equation de transport. Respectivement :
(k)
t
+u
j
(k)
x
j
=

x
j
__
+

t

k
_
k
x
j
_
+P
k
(3.13)
()
t
+
(u
i
)
x
i
=

x
j
__
+

t

_

x
j
_
+

k
(C
1
P
k
C
2
) (3.14)
Dans les equations (3.13) et (3.14), P
k
caracterise la generation denergie cinetique turbu-
lente issue de gradient de vitesse moyenne, et des eets de la pesanteur. C
1
, C
2
sont des
constantes.
k
et

sont appelees nombre de Prandtl, respectivement pour k et

9
. Les
constantes impliquees dans ce mod`ele, que lon souhaite les plus universelles possibles,
ont ete determinees dans les congurations de reference suivantes :
C

= 0,09 Loi logarithmique en paroi

C
1
= 1,44 deformation ou cisaillement uniforme
C
2
= 1,92 decroissance turbulente isotrope

k
= 1 comparaison avec experience jet-sillage

= 1,3 comparaison avec experience jet-sillage

Ce mod`ele a ete presente au debut des annees 1970 par W.P. Jones et B.E. Launder [10],
puis les constantes ont ete determinees par B.E. Launder et D. Spalding [19]. Un des
avantages de cette methode est la prise en compte de la variabilite spatiale de lagitation
turbulente, et sa simplicite de mise en uvre. Il sav`ere etre un des mod`eles les plus
repandus dans les applications pratiques `a lusage de lingenieur (ecoulements dans des
conduites par exemple). Utilise en dehors decoulements cisailles simples pour lesquels il
fut initialement con cu, il conduit `a des resultats qui, sans etre toujours quantitativement
corrects, restent le plus souvent qualitativement representatifs.
Rappelons que le mod`ele k presente plus haut, que lon appelle aussi mod`ele k
standard, nest utilisable que dans les ecoulements fortement turbulents, ce qui nest plus
avere dans des regions proches des parois. Lutilisation de lois de paroi est alors une
solution pour tenir compte de la presence des eets visqueux (Cf. paragraphe 3.4). Une
alternative `a lutilisation de ces lois de paroi a ete mise en place par la modication
des equations de transport de k et en integrant un facteur damortissement pour tenir
compte des eets de paroi [18]. Dans ce cas, la viscosite turbulente secrit :

T
= C

k
2

(3.15)
o` u
f

= 1 e
0,0115y
+
Ainsi lorsque lon consid`ere des regions de plus en plus proches de la paroi, ce coecient
tend vers zero, et les eets de la turbulence diminuent.
9
Ces deux nombres caracterisent la correlation pression-vitesse, et restent tr`es diciles `a modeliser.
48 Mod elisation de la turbulence
Le mod`ele k est dautre part base sur une dependance locale et surtout lineaire des
tensions turbulentes vis `a vis du champ moyen (Eq. 3.3). De ce fait, le mod`ele a tendance
`a surestimer les contraintes turbulentes (en leur donnant un aspect diusif quelles nont
generalement pas) et la viscosite turbulente dans les zones de fort cisaillement (couche
limite, couche de melange). Dautre part, cette linearite fait que ce mod`ele a quasiment
aucune chance de bien reproduire des ecoulements complexes, en particulier tridimension-
nels. De plus, lequation du taux de dissipation est une forme approchee, obtenue `a partir
darguments dont certains sont purement intuitifs. Pour cela, il parat moins adapte aux
ecoulements plus complexes (recirculation, anisotropie forte) que le mod`ele RNG, propose
dans la rubrique suivante.
Mod`ele de fermeture k RNG
Le mod`ele de turbulence RNG k [5] est obtenue `a partir de lequation instantanee de
Navier-Stokes en utilisant une technique mathematique appelee groupe de renormalisa-
tion (RNG)
10
. Il sagit donc dun mod`ele standard ameliore. Le mod`ele fait apparatre
des constantes dierentes de celles du mod`ele standard k :
C

=0,0845 ;
k
=

=0,7179 ; C
1
=1,42 ; C
2
=1,68.
et des termes supplementaires dans les equations de transport de k et . Lelimination de
la longueur presente dans ce mod`ele est eectuee pour decrire le plus d`element possible
la variation de la viscosite eective en fonction du nombre de Reynolds local. Ceci permet
ainsi au mod`ele une meilleur manipulation des ecoulements `a faible Reynolds et proches
des parois.
Le principe de ce mod`ele est de quantier les eets des structures turbulentes de petite
echelle. Cette approche sappuie sur lhypoth`ese duniversalite des petites echelles
11
, plus
connue sous le nom dhypoth`ese de Kolmogorov (1941) [14]. Les grosses structures sont
conservees mais `a la dierence de la LES qui resoud de fa con deterministe les equations
de Reynolds pour ces echelles, le mod`ele RNG k introduit deux equations de transport
pour k et similaires au mod`ele standard.
Ce mod`ele a pour interet majeur daboutir sur des resultats plus realistes dans des
regions de fort gradient en ne produisant pas de diusivite excessive qui faisait defaut
au mod`ele k standard. Ceci est tr`es important dans les cas o` u lecoulement subit des
accelerations ou des decelerations importantes (conduite en presence dun convergent ou
dun divergent) ou lorsquil existe des decollements ou des zones de recirculation.
Mod`ele de fermeture k
Ce mod`ele a ete propose par Wilcox en 1988 [37] et repose sur lutilisation dune echelle
caracteristique des gros tourbillons, = /k, et sur la meme linearite entre la viscosite
10
Cette approche mathematique est utilisee dans differentes branches de la physique pour modeliser
des phenom`enes dynamiques complexes.
11
Ici on fait lhypoth`ese que quelle que soit la conguration de lecoulement, les petites structures se
comportent de fa con identique. Il nen est pas de meme pour les plus grandes structures qui elles voient
leur dimension limitee par la geometrie qui les entoure. Leur comportement sen trouve d`es lors tr`es
inuence, on perd ainsi la notion duniversalite. Remarquons dautre part que cette hypoth`ese a dej` a ete
utilisee dans le cadre de la LES (paragraphe 3.1.2 page 36) permettant la creation dun ltre spatiale.
Mod`ele de turbulence du premier ordre - Concept de la viscosite turbulente 49
turbulente et le cisaillement du champ moyen. La viscosite turbulente est donnee dans ce
cas par la relation suivante :

T
=
k

o` u est un coecient permettant de diminuer linuence de la viscosite turbulente en

eectuant une correction bas-Reynolds. Les equations de transport pour k et sont
respectivement :

t
(k) +

x
i
(ku
i
) =

x
j
_

k
k
x
j
_
+G
k
Y
k
(3.16)
et

t
() +

x
i
(u
i
) =

x
j
_

x
j
_
+G

(3.17)
Dans ces deux equations de transport, G represente la generation denergie cinetique
turbulente ou de la dissipation eective, selon lindice utilise, due au gradient de vitesse
moyenne.
k
et

representent la diusivite eective de k et , respectivement. Y

k
et Y

representent la dissipation de k et par la turbulence.

Linteret principal de ce mod`ele est quil est integrable jusqu`a la paroi. Il donne
de meilleurs resultats que le mod`ele k pour des ecoulements caracterises par des gra-
dients de pression adverses en estimant mieux la position des decollements. Dautre
part, lequation de transport est independante de celle de k ce qui rend ce mod`ele
interessant numeriquement. Par contre, ce mod`ele est tr`es sensible au niveau de tur-
bulence de lecoulement externe `a la couche limite, ce niveau etant xe generalement de
fa con arbitraire par lutilisateur. Tel quel, ce mod`ele ne peut pas etre utilise pour des
ecoulements externes. Lecacite de ce mod`ele pour resoudre lecoulement dans la region
proche de la paroi a souleve lidee dun melange de deux mod`eles de turbulence, le mod`ele
k pour resoudre la region aectee par les eets visqueux, et le mod`ele k pour resoudre
lecoulement dans des regions eloignees des parois. On parle alors de mod`eles hybrides, tel
que le mod`ele SST k.
Mod`ele de fermeture SST k
Le mod`ele Shear Stress Transport k, propose par Menter [25] en 1994 est un mod`ele
hybride, cest `a dire utilisant des mod`eles de turbulence dierents pour la region proche
paroi et dans le reste de lecoulement. La transition dun mod`ele vers lautre se fait progres-
sivement `a partir dun melange. Lapproche k est utilisee pour predire lecoulement
dans la region aectee par la viscosite, et poss`ede lavantage de ne pas necessiter de
fonction damortissement, f, alors que lapproche k est utilisee pour lautre partie de
lecoulement. Une fonction permet dassurer la transition dun mod`ele `a lautre.
Les equations de transport du mod`ele SST k sont similaires aux equations (3.16) et
(3.17) utilisees par le mod`ele mod`ele k standard.
3.2.4

Equation nale
Rappelons que la resolution de toutes les equations enoncees dans les paragraphes
prec`edents ont pour but destimer la viscosite turbulente (Eq. 3.3) necessaire dans lob-
tention du champ moyen. Cette estimation introduite dans les equations de Reynolds
permet decrire :
50 Mod elisation de la turbulence

_
u
i
t
+u
j
u
i
x
j
_
= f
i

p
x
i
+

x
j
_

e
_
u
i
x
j
+
u
j
x
i
__
(3.18)
avec
e
= +
T
la viscosite eective. Cette derni`ere est elevee dans les regions com-
portant de forts gradients de vitesse introduisant une production denergie turbulente.
Lecoulement moyen est donc regi par cette equation, en plus de lequation de continuite
moyennee (Eq. 3.1).
3.3 Mod`ele de turbulence du deuxi`eme ordre
Nous lavons vu plus haut la relation utilisant le concept de la viscosite turbulente est
basee sur des hypoth`eses fragiles (page 42). De ce fait, les mod`eles presentes dans les para-
graphes precedents paraissent inadaptes `a des ecoulements complexes. Laissant ce concept
de viscosite turbulente, le mod`ele RSM (pour Reynolds Stress Model) [17, 8, 16] propose
de resoudre lequation de transport du tenseur de Reynolds accompagnee dune autre
equation de transport (pour le taux de dissipation , ou la dissipation specique ) per-
mettant, elle, destimer la longueur caracteristique de turbulence. Ainsi, la modelisation
au second ordre resout les equations du tenseur de Reynolds, alors que les mod`eles du
premier ordre les modelisent. Les mod`eles du second ordre contiennent donc plus de phy-
sique.
Il est possible dobtenir lequation de transport du tenseur de Reynolds en prenant la
moyenne de la quantitee ci-dessous :
u

i
N(u
j
) +u

j
N(u
i
)
o` u N(u
j
) correspond `a loperateur equation instantanee de Navier-Stokes, i.e.,
N(v
j
) =
_
u
i
t
+u
j
u
i
x
j
_
f
i
+
p
x
i

ij
x
j
Lequation de transport du tenseur de Reynolds peut secrire schematiquement :
D
Dt
(u

i
u

j
) +
T
ijk
x
k
= P
ij
+
ij

ij
o` u T
ijk
est un terme en correlation triple, P
ij
est terme de production denergie cinetique
turbulente,
ij
un terme de transport par la pression, et
ij
le tenseur de dissipation (et
non pas le taux de dissipation donne par la relation 3.12) :

ij
= 2
u

i
x
k
u

j
x
k
Le terme T
ijk
secrit
T
ijk
= u

i
u

j
u

k
+p
_

kj
u

i
+
ik
u

j
_

x
k
(u

i
u

j
)
et fait apparatre un tenseur dordre 3 inconnu et necessitant alors une modelisation.
Darly et Harlow [6] proposent par exemple la modelisation suivante :
u

i
u

j
u

k
= C
s
u

k
u

l
u

i
u

j
x
l
Inuence de la paroi 51
Une modelisation est aussi necessaire pour les termes
ij
et
ij
. Le terme de production,
quant `a lui, ne necessite aucune modelisation, ce qui est linteret du mod`ele RSM, car il
ne fait intervenir que des moments dordre inferieur ou egal `a deux :
P
ij
=
_
u

i
u

k
u
j
x
k
+ u

j
u

k
u
i
x
k
_
Le mod`ele RSM est un mod`ele dt Haut Reynolds et necessite alors lutilisation dune loi
de paroi pour traiter la region aectee par les eets visqueux.
3.4 Inuence de la paroi
Les mod`eles de turbulence bases sur lhypoth`ese dune forte turbulence ne sont pas
valables dans la region proche des parois. Il existe alors plusieurs solutions pour tenir
compte de la presence des eets visqueux : (1) ne pas resoudre cette region et utiliser `a
la place une loi de paroi, (2) introduire des fonctions damortissement for cant le compor-
tement du mod`ele, appele d`es lors mod`ele bas-Reynolds, ou (3) resoudre des equations
dierentes pour lecoulement principal et lecoulement proche des parois. Les mod`eles k
standard ou RNG sont bases sur lapproche (1), le mod`ele de Spallart et Allmaras sur
lapproche (2), et le mod`ele SST k sur lapproche (3).
Ainsi dans la premi`ere approche, la region interne nest pas resolue et une loi de paroi
est alors utilisee pour relier la region aectee par les eets visqueux `a la partie de la
couche limite pleinement turbulente (Fig. 3.10.a). Cette loi plus ou moins evoluee permet
destimer les vitesses moyennes ainsi que les quantites turbulentes dans la partie non-
resolue mais repose neanmoins sur des hypoth`eses tr`es contraignantes. Il existe dierents
niveaux de lois de paroi plus ou moins evoluees, la plus simple etant la loi de paroi
standard (relation 1.3) aper cue au paragraphe 1.4 valable dans la region caracterisee par
y
+
= 30. Linconvenient de cette loi provient du fait quelle repose, entre autre, sur un
nombre de Reynolds de lecoulement tr`es important, labsence dun gradient de pression
adverse. Pour satisfaire ces points, dautres lois plus elaborees sont alors disponibles et
construites la plupart du temps sur une decomposition en deux couches de la couche
limite (loi bi-couche). En fonction du point de calcul place dans la couche limite, la loi
de paroi aectera des valeurs adequates. Ainsi lutilisation de ce type de loi permet de
retrouver des informations dans la sous-couche visqueuse. La valeur du crit`ere de maillage
doit donc satisfaire y
+
1. On note donc que lapproche des lois de paroi peut etre moins
consommatrice en temps de calculs et en espace de stockage en fonction du degre de
complexite de la loi, la taille du premier element etant plus importante.
Les approches (2) et (3) necessitent quant `a elles dans tous les cas des tailles de mailles
beaucoup plus petites, assez petites pour chercher des informations dans la sous-couche
visqueuse (Fig. 3.10.b). Le nombre delements de calcul devient alors plus eleve que celui
dun calcul base sur des mod`eles haut-Reynolds utilisant des lois de paroi standards.
Un cas test classique pour tester les capacites des lois de paroi est la couche limite
au-dessus dune plaque plane (Fig. 3.11). On consid`ere dans cet exemple un domaine
dune longueur egale `a 1 m avec les nuds de la ligne inferieure places sur la paroi. La
vitesse de lecoulement exterieur `a la couche limite est egale `a 50 m/s. Le maillage est
construit `a partir de cette valeur ainsi que celle du y
+
qui doit etre en adequation avec
lapproche numerique utilisee dans la region de paroi. On observe quun calcul eectue `a
52 Mod elisation de la turbulence
Loi de paroi
(a) (b)
Fig. 3.10 Traitement numerique de la region proche dune paroi. a, approche des lois de paroi ; b, approche
des mod`eles bas-Reynolds.
laide dun mod`ele bas-Reynolds necessite deux fois plus de volumes quun calcul base sur
des mod`eles haut-Reynolds utilisant une loi de paroi standard.
50 m/s
L = 1 m
L/2
Paroi
S
o
r
t
i
e
E
n
t
r
e
e
Symetrie
y
+
Taille de N
bre
de volume
maille (m)
1 1,25 10
5
14 751
30 3,76 10
4
7 821
Fig. 3.11 Domaine de calculs dune couche limite au-dessus dune plaque plane. Le tableau rassemble les
valeurs de la taille des premiers elements proches de la paroi estimes `a partir de la valeur du crit`ere
y
+
. Les deux maillages utilisent 100 elements le long de la paroi, et une evolution des tailles
delements dans la direction transversale egale `a 5%.
3.5 Choisir le mod`ele de turbulence
Aucun mod`ele de turbulence nest universellement reconnu superieur aux autres pour
toutes les classes de probl`eme. Le choix du mod`ele de turbulence va dependre de consi-
derations telles que la complexite de lecoulement aborde, le niveau de precision requis,
Modelisation instationnaire 53
le temps de calcul et la place memoire disponibles. Pour eectuer le meilleur choix, il
est important de connatre les limites de chacun des mod`eles de turbulence des codes
CFD. Quand cela est possible, une comparaison de resultats numeriques obtenus `a partir
de dierents mod`eles avec des resultats experimentaux reste une aide dans le choix du
mod`ele le plus approprie.
3.6 Modelisation instationnaire
Dans certaines etudes industrielles, il est tr`es important de connatre levolution tem-
porelle du champ instantane. Seules les methodes deterministe et semi-deterministe sont
adaptees `a lobtention de cette information, alors que la methode statistique par essence
ne lest pas. En eet, le champ resolu correspond `a la moyenne temporelle par ergodicite.
Le champ resolu ne contient donc aucune information temporelle, et aucune information
instationnaire nest disponible.
Malgre cela, certaines methodes ont ete mises en place depuis les annees 1990 et sont
basees sur la modication du mod`ele RANS sur des bases plus ou moins empiriques an de
resoudre les instationnarites `a grande echelle. Parmi lensemble des methodes disponibles,
la methode URANS (pour Unsteady Reynolds Average Navier-Stokes) sest repandue avec
succ`es dans le monde industriel. Cette methode consiste `a resoudre les equations RANS
vues plus haut avec une resolution du terme instationnaire (i.e. /t) sans modication
des constantes impliquees dans les mod`eles. Cependant ce terme doit tendre vers zero
dans le cas dun ecoulement statistiquement stationnaire. Or des simulations URANS ont
montre ses capacites `a capturer des structures coherentes `a grande echelle dans diverses
congurations decoulement, avec une bonne estimation du nombre de Strouhal
12
.
Dans le cas o` u on sinteresse au champ aerodynamique moyen (dans le temps) dun
ecoulement caracterise par une generation suivie dune convection dinstabilites, il est
preferable de resoudre le probl`eme de fa con instationnaire au lieu dune resolution per-
manente. Les resultats obtenus sont plus proches des valeurs experimentales. Ceci est le
cas par exemple dun ecoulement au passage dun cube pose sur le sol. Une resolution en
permanent entraine une sous-estimation du coecient de tranee, une surestimation de la
region de recirculation dans le sillage de lobstacle (Fig. 3.12).
3.7 Et lindustrie dans tout ca ?
Dapr`es ce que nous venons de voir, il existe trois principales methodes dapproxima-
tion des eets de la turbulence : la DNS, la LES, et lapproche statistique. Notons que
les deux premi`eres nous permettent de toucher au plus pr`es la turbulence, alors quun
travail statistique au contraire nous eloigne de cette physique pour nous entraner vers
des formulations de fermetures que lon peut dicilement valider.
Cependant, noublions pas qu`a l`ere des ordinateurs de faible puissance
13
, les methodes
statistiques etaient lunique approche pour le petit industriel ne disposant pas, `a la fois
12
Ce nombre sans dimension caracterise les eets instationnaires.
13
Ce nest pas si vieux que ca...
54 Mod elisation de la turbulence
Fig. 3.12 Lignes de courant dans le cas dun ecoulement au passage dun cube (Re
h
= 40 000 - a, donnees
experimentales dapr`es Martinuzzi et al. [24] ; b, LES dapr`es Shah [33] ; c, RANS permanent
dapr`es Iaccarino et al. [9] ; d, RANS transitoire dapr`es Iaccarino et al. [9]).
Resolution numerique ` a partir dun code de calculs 55
de supercalculateurs reserves aux universitaires ou aux grands groupes, mais aussi des
competences adequates.
Aujourdhui, le developpement tr`es important de la puissance des calculateurs met en
avant la methode LES, qui est largement proposee dans la plupart des codes commerciaux
(Fluent, Star CD, CFX...). Il est donc fort `a parier que la simulation
14
sera le nouvel
outil numerique industriel. Cependant, son utilisation recqui`ere des connaissances sur ces
methodes et une experience certaine.
Pour ces raisons, de nos jours le calcul statistique demeure la methode la plus utilisee,
permettant dautre part davoir des resultats corrects pour des temps de calcul satisfai-
sants.
3.8 Resolution numerique `a partir dun code de cal-
culs
La modelisation de la turbulence nest evidemment pas une n en soi, elle est la base
de methodes numeriques proposant des solutions aux equations dierentielles vues plus
haut. Celles-ci sont obtenues `a partir de codes de calculs, commerciaux ou developpes par
des laboratoires de recherche universitaires et industriels.
Le succ`es de ces codes est etroitement lie au developpement des capacites dexecution
et de stockage des ordinateurs, le rapport performance/co ut ayant incroyablement aug-
mente d`es les annees 1950.
`
A cette date, les machines eectuaient quelques centaines
doperations `a la seconde alors que de nos jours les plus rapides op`erent 10
12
operations
pendant cette meme seconde. De meme, ces premi`eres machines, tr`es onereuses et oc-
cupant une pi`ece enti`ere, demandaient une maintenance permanente. Aujourdhui un
ordinateur de bureau poss`ede les memes capacites sinon plus.
On peut alors comprendre aisement linteret grandissant des industriels dutiliser ces
codes de calculs. Cette demande est tellement importante quelle occupe un tiers des
chercheurs en mecanique des uides, ce nombre grandissant encore. Ce champ de recherche
est appele CFD (pour Computational Fluid Dynamics).
Il est important de comprendre que la CFD ne donne quune approximation de la so-
lution recherchee. La dierence entre le calcul et la realite provient derreurs accumulees
par le code, et provenant dorigines multiples : (1) les equations dierentielles generales
sont souvent simpliees en ne gardant que les termes preponderants ; (2) une approxima-
tion est eectuee lors de letape de discretisation ; (3) la resolution est atteinte par des
methodes iteratives et `a moins de ne faire calculer tr`es longtemps le code, les solutions
exactes des equations discretisees ne sont jamais obtenues.
Les composants importants dun code de calculs en mecanique des uides numerique
sont les suivants :
1. Le mod`ele mathematique : le depart de toute resolution numerique est le mod`ele
mathematique, cest `a dire lensemble des equations dierentielles et des conditions
aux limites. Ces premi`eres changent selon le type decoulement traite (incompressible
ou compressible, turbulent, bi- ou tridimensionnel, etc.). Comme il est mentionne
plus haut, ce mod`ele peut inclure des simplications des lois de conservation. En
eet, il est illusoire de mettre en place un mod`ele universel et performant. Les
14
Comparativement `a la modelisation.
56 Mod elisation de la turbulence
equations de bilan ecrites plus haut (Eqs. 2.6 et 2.14) peuvent etre ecrites sous une
forme generalisee (si on ne tient pas compte du terme lie au temps) :
(u
j
)
x
j
=

x
j
_

x
j
_
+q

o` u le terme du membre de gauche est un terme convectif, alors que le premier

terme du membre de droite est un terme diusif. Dans le cas de lequation de
conservation de la masse = 1, et pour le bilan de quantite de mouvement = u
i
.
On supposera les quantites , u
j
, , et q

connues. Ceci nest evidemment pas le

cas car la vitesse `a cette etape nest pas encore calculees, et les proprietes du uide
peuvent changer avec la temperature, et si comme nous lavons vu des mod`eles de
turbulence sont employes, du champ de vitesse aussi. La resolution par methode
iterative consid`ere comme la seule inconnue, les autres etant xees et calcules
dans literation precedente.
2. La methode de discretisation : Il sagit de methodes dapproximation des equations
dierentielles regissant le mouvement du uide par un syst`eme dequations algebriques.
Il existe plusieurs approches, mais les plus utilisees sont les dierences nies, les vo-
lumes nis, et les elements nis. Toutes ces methodes aboutissent `a des solutions
identiques si le maillage est assez n. Neanmoins, certaines methodes sont plus
adaptees `a certains ecoulements.
3. Le maillage : les localisations des points sur lesquels sont resolues les equations
algebriques sont denies par le maillage. Selon lagencement des elements, il peut
etre structure, structure par bloc, non-structure, chacune de ces congurations ayant
ses avantages et ses inconvenients.
4. Lapproximation sur les elements : il faut `a partir des elements mis en place,
denir les approximations `a utiliser. Les solutions sont obtenues sur des nuds,
lapproximation permet datteindre des valeurs sur lelement.
5. La methode de resolution : la discretisation aboutit `a un syst`eme dequations
algebriques (mis sous forme de matrice generalement). Il existe des methodes directes
dinversions et dautres par iterations.
6. Les crit`eres de convergence : les methodes de resolution iteratives demandent
un crit`ere darret. Celui-ci est xe par lutilisateur.
Chapitre 4
Application de la modelisation de la
turbulence `a la CFD
La mecanique des uides numerique, plus connue sous lacronyme CFD pour Computa-
tional Fluid Dynamics, est lanalyse dun syst`eme impliquant un ecoulement de uide, des
echanges thermiques, et des phenom`enes associes, `a partir dune simulation aidee par ordi-
nateur. Cette technique est tr`es puissante et permet detudier des varietes decoulements
du monde industriel ou non tr`es vastes. Ces applications peuvent traiter les eorts exerces
par un ecoulement sur un obstacle (i.e. portance et tranee), les ecoulements dans les turbo-
machines, la combustion dans les moteurs, le refroidissement de composants electroniques,
les transports de polluants par lair, la meteorologie, lingenierie medicale, etc. Loutil CFD
est devenu presque imperatif lors de la conception dun produit. Lindustrie aeronautique
a integre depuis les annees 1960 cet outil dans le dessin, la R&D, et la fabrication des
avions et des moteurs. Certe le co ut de lutilisation de la CFD lie aux licenses, au materiel
informatique, et aux personnes ayant les competences peut etre important. La CFD a ce-
pendant plusieurs avantages tels que la possibilite destimer des champs aerodynamiques
dans des regions inaccessibles experimentalement, des temps de realisation interessants
comparativement au temps necessaire pour realiser des mod`eles reduits pour souerie.
Ainsi les etudes parametriques numeriques deviennent peu couteuses et permettent dat-
teindre une quantite tr`es importante dinformations concernant les ecoulements.
Notons que meme si les codes de calcul deviennent de plus en plus evolues, il demeure
important de valider reguli`erement les resultats numeriques par comparaison avec des
donnees experimentales. De meme, le succ`es dun calcul CFD repose generalement sur
lexperience de lutilisateur et de ses connaissances scientiques.
58 CFD
4.1 Structure dun code CFD
Les codes CFD sont articules autour dalgorithmes de resolution des probl`emes decou-
lements. Ils accueillent generalement des interfaces pratiques et sophistiquees pour le choix
des param`etres de resolution ou de post-traitement. Ainsi les codes contiennent trois
elements importants : (1) un pre-processeur, (2) un solver, et (3) un post-processeur.
4.1.1 Le pre-processeur
Le role du pre-processeur est daccueillir les donnees dentree du probl`eme `a resoudre et
`a les transformer de fa con utile au solver. Ces donnees comprennent le domaine de calcul
maille, les phenom`enes physiques `a modeliser (turbulence, combustion, multiphasique,
etc.), les proprietes physiques du uide, et les conditions limites adaptees au probl`eme.
La solution du probl`eme (i.e. vitesse, pression, temperature) est denie aux nuds de
chaque cellule. La precision de cette solution, mais aussi le temps de calcul, dependent
grandement du nombre de nuds. Un maillage optimise est generalement un maillage
non-uniforme. Ainsi la taille des mailles devra etre petite dans des regions aectees par
des evolutions importantes de la solution (i.e. le sillage dun obstacle, la region proche
dune paroi, etc.) alors que cette taille pourra etre beaucoup plus importante dans les
regions caracterisees par des faibles variations.
Lutilisateur des codes CFD passe plus de 50% de son temps `a dessiner le domaine de
calcul et `a le mailler. En consequence, les codes proposent maintenant des outils propres
de conception ou facilitent limportation de dessins issus doutils de CAO. Dautre part,
les outils de maillage permettent de tenir compte tr`es tot dans le projet de particularites
de maillage liees aux phenom`enes `a traiter par la suite.
4.1.2 Le solver
Les algorithmes du solver ont pour but (1) dintegrer les equations de bilan sur les
volumes de controle, (2) de discretiser les formulations integrales en un syst`eme algebrique
dequations, et (3) destimer la solution de ces equations `a partir de methodes iteratives.
Il existe principalement trois methodes de discretisation : (i) les dierences nis, (ii) les
volumes nis, et (iii) les elements nis. Les codes de CFD sont generalement bases sur
la methode des volumes nis car eectuer des bilans dequations de bilan elles-meme sur
chaque volume rend tr`es attractive cette methode et explique le succ`es de cette approche
aupr`es des fournisseurs de code CFD.
La conservation dune quantite `a travers un volume de controle peut etre exprimee
par lequilibre des processus tendant augmenter ou diminuer cette quantite :
_

_
Taux de
variation de la
quantite
au cours du temps
_

_
=
_

_
Augmentation
de due `a
la convection `a
travers le volume
_

_
+
_

_
Augmentation
de due `a
la diusion `a
travers le volume
_

_
+
_

_
Creation de
`a
linterieur
du domaine
_

_
Les codes de CFD disposent de techniques adaptees pour le traitement des phenom`enes de
transport lies `a la convection (transport consecutif `a lecoulement), `a la diusion (trans-
port du aux variations point-`a-point dune quantite ), et la variation dans le temps de
cette quantite.
Le maillage 59
4.1.3 Le post-processeur
Les post-processeurs ont les capacites dacher le domaine de calcul et son maillage, les
vecteurs-vitesse, des contours, des suivis de particules, des animations, et deectuer des
operations sur lensemble des solutions atteintes par le solver. Lutilisateur a la possibilite
dautre part dexporter ses resultats dans des formats standards pour une utilisation
exterieure.
4.2 Le maillage
Le maillage est une representation discr`ete du domaine reel accueillant le uide en mou-
vement. Il existe dierents types de maillage, structures ou non. Les maillages structures
proposent des mailles alignees les unes aux autres sans que les lignes ne se croisent (Fig.
4.1.a). Ce type de maillage est particuli`erement adapte aux formes simples de geometrie. Il
est neanmoins possible de modier le maillage pour des surfaces courbes (Fig. 4.1.b). Pour
(a) (b)
Fig. 4.1 a, maillage structure de type H; b, maillage structure de type O adapte `a la courbure de la surface
de lobstacle (technique O-Grid).
les geometries plus complexes, il devient impossible de realiser un maillage structure. Les
maillages non-structures permettent alors de mailler convenablement le domaine. Dans ce
cas les elements sont des triangles ou des quadrilat`eres en 2-D (Fig. 4.2), et des tetra`edres
et des hexa`edres en 3-D. Des mailles de ce type peuvent etre generes `a partir dalgo-
rithmes (Delaunay, etc.). Ce type de maillage peut etre rane simplement localement.
Cette simplicite est malheureusement desequilibree par la complexite de la manipulation
des donnees liees au maillage. Les matrices du syst`eme algebrique des equations deviennent
mal-conditionnees, et les solvers adaptes aux maillages irreguliers sont generalement plus
lents que ceux ecrits pour les maillages structures.
4.3 Les conditions limites
Les conditions limites sont appliquees `a tous les nuds places aux fronti`eres du do-
maine et doivent etre coherentes avec le phenom`ene simule. Lapplication de conditions
limites inadaptees conduit `a des resultats imprecis, des instabilites de la resolution, et
60 CFD
Fig. 4.2 Maillage non structure autour de trois cylindres.
de mauvaises convergences du calcul. Les conditions limites suivantes sont generalement
utilisees : la condition de paroi (i.e. adherence pour les uides visqueux et la condition de
glissement pour le uide denue de viscosite), la condition dentree, la condition de champ
lointain, la symetrie, les conditions de periodicite, et les interfaces entre dierents blocs
constituant le domaine.
Il est preferable dimposer les conditions limites dans des regions assez eloignees de
celle dont est porte linteret. Si par exemple une condition de sortie est placee dans une
region trop proche dun obstacle present dans lecoulement, ce dernier ne pourra atteindre
un prol de turbulence parfaitement developpee en sortie de domaine. Ce mauvais condi-
tionnement peut conduire `a des erreurs numeriques. Il est donc necessaire de placer la
condition de sortie dans une region pour laquelle lecoulement nevolue plus. Par exemple
dans le cas de simulations decoulements au passage dune singularite dans une conduite
il est preferable de placer la condition initial trois voire quatre diam`etres hydrauliques
en amont et une vingtaine de diam`etres hydrauliques en aval de la singularite `a etudier.
Pour abaisser le nombre de nuds la condition limite dentree peut etre imposee dans
une region plus proche `a condition dimposer cette fois un prol de turbulence pleinement
developpee au lieu dune vitesse constante.
4.4 Les conditions initiales
Les conditions initiales, pour toutes les variables `a calculer, sont necessaires pour
entamer tout calcul numerique. Dans le cas dun calcul permanent, les valeurs initiales des
variables correspond aux proprietes dun ecoulement `a partir duquel le calcul numerique
va evoluer. La convergence du resultat est plus rapide si celles-ci sont proches des valeurs
nales. Neanmoins, le resultat converge ne doit pas dependre des valeurs initiales. Dans
le cas dun calcul transitoire, ces valeurs correspondent `a lecoulement au temps initial du
calcul.
La plupart des codes de calculs proposent dierentes methodes daectation de ces
valeurs : (i) des valeurs generees automatiquement, (ii) des valeurs explicites, et (iii) la
Les schemas numeriques de convection 61
lecture dun chier de resultats. Dans le premier cas, le code estime les valeurs initiales `a
partir des conditions limites proposees par lutilisateur. Le dernier cas est souvent utilise
lorsque des diculte apparaissent dans la convergence. Ainsi la solution peut etre calculees
`a partir de dierentes etapes dans laquelle un param`etre voit sa valeur augmenter. Chaque
champ initial dun calcul correspondant alors `a la solution du calcul precedent. De meme
on peut entamer une etude instationnaire `a partir de la solution dun calcul permanent
realise dans des conditions similaires.
4.5 Les schemas numeriques de convection
Dans les probl`emes impliquant un ecoulement de uide, les termes de convection jouent
un role tr`es importants. La qualite des resultats, ainsi que la convergence du calcul,
dependent grandement de lestimation de cette grandeur. Il existe plusieurs methodes
dapproximation de ces termes (Upwind, power law, QUICK, etc.) se dierenciant par
lordre des schemas (premier ou second ordre generalement). Par exemple, le schema
Upwind correspond `a un schema du premier ordre, donc moins precis quun schema du
second ordre, mais tr`es stable. Au contraire, lutilisation de schemas du second ordre peut
entrainer dans certains cas des probl`emes de convergence. On pourra alors dans ce cas
entamer le calcul `a partir dun schema du premier ordre, puis dans une deuxi`eme etape,
realiser un calcul base sur les resultats de ce calcul en utilisant cette fois un schema du
second ordre. Les risques doscillations de la solution sont diminues par cette procedure.
4.6 Convergence dun calcul
Une etape importante lors dun calcul est lobservation de la convergence de ce calcul.
Ceci peut etre realise `a partir du trace des residus pour chaque equation `a resoudre par
le code. Le but de cette section est de presenter bri`evement la notion de residu, et les
ameliorations possibles `a apporter pour atteindre une bonne convergence.
4.6.1 Les residus
Lapplication des equations de bilan sur lensemble des volumes nis obtenus par
discretisation spatiale du domaine permet decrire un syst`eme compose de n equations
avec n inconnues :
[A][x] = [b] ou
n

j=1
A
ij
x
j
= b
i
Les methodes iteratives permettant dobtenir les solutions x
j
en ninversant pas directe-
ment la matrice A
ij
isole la solution x
i
telle que
A
ii
x
i
= b
i

j=1,j=i
A
ij
x
j
62 CFD
La solution du membre de gauche est evaluee `a literation (k) `a laide des solutions x
j
du
membre de droite evaluee quant `a elle `a literation precedente (k 1) :
x
(k)
i
=
n

j=1,j=i
_

A
ij
A
ii
_
x
(k1)
j
+
b
i
A
ii
que lon peut manipuler pour faire apparatre le terme x
(k1)
i
dans le membre de droite.
On obtient :
x
(k)
i
= x
(k1)
i
+
n

j=1
_

A
ij
A
ii
_
x
(k1)
j
+
b
i
A
ii
(4.1)
Par le jeu des approximations numeriques, la solution x
j
atteinte `a literation (k) ne
satisfait pas le syst`eme dequations precedents, do` u :
n

j=1
A
ij
x
(k)
j
= b
i
On denit alors le residu r
(k)
i
pour lequation i apr`es literation (k) comme la dierence
entre le membre de droite et le membre de gauche :
r
(k)
i
= b
i

j=1
A
ij
x
(k)
j
Si la resolution converge de fa con satisfaisante, les solutions intermediaires x
(k)
j
se rap-
prochent dune solution telle que les residus tendent vers une valeur nulle. Le residu moyen
considerant les n equations, cest `a dire en considerant toutes les solutions atteintes sur
les volumes nis du domaine, vaut :
r
(k)
=
1
n
n

i=1

r
(k)
i

Il est souvent utile de normaliser cette quantite par le residu de la premi`ere iteration :
R =
r
(k)
r
(1)
Ceci rend plus facile la comparaison des convergences des autres solutions. Ainsi il existe
autant de residus quil y a dequations `a resoudre sur les volumes nis : les equations
selon x, y, et z, les equations de transport des variables caracterisant la turbulence (k, ,
, etc.).
4.6.2 La sous-relaxation
La convergence nest pas toujours satisfaisante et dans ce cas la solution ne doit pas
etre utilisee pour une analyse de lecoulement reel. Ceci peut etre le cas lorsque les condi-
tions limites ne sont pas adaptees, lorsque le maillage est trop grossier et necessite donc
un ranement dans des regions impliquant des gradients importants, ou lorsque des eets
dynamiques sont importants rendant instable la solution. Une solution consiste alors `a
Convergence dun calcul 63
ralentir la marche de la resolution numerique par lintroduction de coecients de re-
laxation, , dans les equations precedentes. Ainsi la solution `a literation (k) est estimee
`a partir dune partie seulement de celle estimee `a literation (k 1) :
x
(k)
i
= x
(k1)
i
+
_
n

j=1
_

A
ij
A
ii
_
x
(k1)
j
+
b
i
A
ii
_
(4.2)
Des valeurs du coecient 0 < < 1 conduit `a la technique de sous-relaxation, et > 1
`a la technique de sur-relaxation.
La convergence des residus normalises de lequation de mouvement selon x est tracee
sur la gure 4.3 pour le cas dun ecoulement tr`es rapide au passage dun nombre inni de
cylindres circulaires alignes (utilisation de condition de periodicite dans ce cas). En ab-
sence de sous-relaxation, le calcul est arrete par le code de calcul alors quune convergence
satisfaisante est atteinte d`es lors quune sous-relaxation est imposee `a levolution de la
solution. On note ainsi que la sous-relaxation permet datteindre des solutions convergees
dans des cas diciles avec cependant un nombre diterations plus importants.
Fig. 4.3 Inuence de la
technique de sous-relaxation sur
la convergence des residus des
equations selon x.
4.6.3 Calcul permanent `a partir dune approche pseudo-transitoire
On note une similitude des equations obtenues pour lapplication de la technique de
sous-relaxation dune resolution dun ecoulement permanent dune part, et celles obtenues
`a laide dune discretisation temporelle dautre part. On peut alors mettre en evidence
une relation analytique entre le coecient de relaxation du calcul permanent et le
pas de temps de lapproche transitoire, dt. Certains codes de calculs utilisent alors cette
possibilite dobtenir les eets de la sous-relaxation sur un calcul permanent au moyen
dun calcul pseudo-transitoire utilisant un pseudo pas de temps.
64 CFD
Deuxi`eme partie
Exemples de predictions numeriques
decoulements et dechanges
thermiques
Interet des etudes de la deuxi`eme
partie du document
La suite du document propose quelques etudes visant `a estimer numeriquement dans
des cas standards les champs aerodynamiques, ainsi que les champs thermiques lorsque
les transferts de chaleur ne sont plus negligeables. Ces calculs ont ete menes `a partir de
deux codes du commerce : Fluent et Ansys-Flotran. Le premier repose sur la methode
numerique des volumes nis tout particuli`erement adaptee `a la prediction des ecoulements
de uide, alors quAnsys-Flotran est base sur la technique des elements nis tr`es utilisee
en calcul de structure.
Il est important de comprendre que linteret principal de ces etudes nest pas de
sarreter sur la forme des champs aerodynamiques et thermiques estimes mais plutot
de simpregner des methodes utilisees pour aboutir `a ces solutions. Ces etudes permet-
tront dautre part dobserver les inuences de maillages ou mod`eles de turbulence sur le
resultats naux. Ainsi le chapitre 5 se propose destimer le champ aerodynamique dun
ecoulement au passage dune colline, le but etant de sinteresser aux resultats atteints par
dierents mod`eles de turbulence. On sinteresse ensuite (chapitre 6) `a un ecoulement au
passage dun obstacle prismatique, le but etant de mettre en evidence la methode desti-
mation dun pas de temps adapte pour lobservation dinstabilites. Le chapitre 7 soul`eve
le probl`eme des conditions limites `a imposer `a lentree dune conduite an detre au plus
pr`es des grandeurs reelles. Un calcul annexe est mene et permet destimer des prols trans-
versaux des grandeurs moyennes et turbulentes, qui seront ensuite considerees comme des
conditions limites dans le calcul principal. Les prols obtenus `a partir de cette methode
seront confrontes `a ceux donnes par une formulation analytique. Le chapitre 8 vise quant
`a lui `a montrer la disparite des resultats obtenus `a partir de mod`eles haut-Reynolds et
un mod`ele bas-Reynolds. Pour cela, on consid`ere un echange thermique convectif entre
un ecoulement et un obstacle elance. La qualite des predictions numeriques sera jugee
par comparaison avec des resultats experimentaux obtenus dans des congurations iden-
tiques. Enn, le chapitre 9 met en evidence linuence de la nesse dun maillage sur la
prediction des echanges thermiques basee sur lutilisation du mod`ele k . Une discussion
est menee `a propos de ladequation entre la taille de la premi`ere maille et la loi de paroi
68
non-equilibree employee dans cette etude.
Chapitre 5
Cas test : 2D Model Hill Flow
On desire dans cette etude confronter les resultats obtenus par le code de calcul FLO-
TRAN avec ceux obtenus dans une conguration reelle. On choisit pour cela un cas test
appele 2D Model Hill Flow
1
(Fig. 5.1), mene par Almeida et. al [1] en 1992.
Fig. 5.1 Cas test 2D Model Hill Flow (dapr`es [1]).
Lauteur met `a la disposition des numericiens le prol de la colline (fonction po-
lynomiale), le prol de vitesse en entree de veine dessai (Fig. 5.2.a), et les prols de
vitesse `a dierentes abscisses (Fig. 5.2.b). Les mesures ont ete eectuees par la methode
L.D.V. (Laser Doppler Velocimeter).
1
Les donnees experimentales sont disponibles sur le site internet http ://[Link]/
70 Cas test : 2D Model Hill Flow
Fig. 5.2 Donnees experimentales du cas test. a, vitesse longitudinale `a lentree ; b, vitesses longitudinales `a
trois abscisses de la veine (dapr`es Almeida [1]).
5.1 Approche numerique
5.1.1 Domaine de calcul et maillage
Lecoulement est calcule sur un domaine representant une partie du tunnel hydrau-
lique dans lequel a ete eectue lensemble des mesures. Lentree et la sortie sont places
respectivement `a des distances egales `a 3H et 30H du centre de la colline, H representant
sa hauteur. La hauteur du domaine est placees `a une longueur de 6, 071H identiquement
au tunnel.
Le maillage est inhomog`ene et constitue delements inclines selon lecoulement (Fig.
5.3).
Fig. 5.3 Maillage de
letude de lecoulement au
passage dune colline.
Le syst`eme dequations de Navier-Stokes moyennees sera fermee `a partir de plusieurs
mod`eles de turbulence, un mod`ele `a zero equation de transport, le mod`ele k standard,
et le mod`ele k RNG.
Le mod`ele de fermeture `a zero equation est le plus simple des trois mod`eles etudies
ici, et dans ce cas la viscosite turbulente est estimee par la relation

T
= L
s

o` u est le taux de dissipation visqueuse, et L

s
une longueur, respectivement
=
_
u
i
x
k
+
u
k
x
i
_
u
i
x
k
et L
s
=
_
L
x
si L
x
> 0
min{0, 4L
n
; 0, 09L
c
} si L
x
0
avec L
x
une echelle representant lobstacle (ici la hauteur de la colline), L
n
la distance
entre le nud considere et la paroi la plus proche, et L
c
etant la valeur la plus elevee de
L
n
.
Les deux autres mod`eles ont ete presentes dans le chapitre 3.
Resultats numeriques 71
5.1.2 Conditions limites
On impose sur les nuds de lentree le prol experimental de vitesse, ainsi quune
intensite turbulente de 3%. Une condition dadherence est imposee sur les nuds des
surfaces inferieure et superieure, et une condition de turbulence pleinement developpee en
sortie.
5.2 Resultats numeriques
Les resultats du code font apparatre une zone de recirculation comme le montre
lexperience, mais plus ou moins grande selon le mod`ele de turbulence employe. Les gures
5.4.a et 5.4.b representent ces recirculations respectivement pour un mod`ele k et rng.
Fig. 5.4

Ecoulement autour de la colline, dapr`es flotran. a, mod`ele k ; b, mod`ele rng.
La confrontation des vitesses experimentales et celles simulees par dierents mod`eles
de turbulence, ceci pour dierentes abscisses, est illustree sur la gure 5.5.
Fig. 5.5 Comparaison des donnees experimentales [1] et celles de flotran.
On remarque que le mod`ele le plus simple, `a zero equation, dune part ne predit
quune tr`es faible zone de recirculation. A labscisse x=70 mm, celle-ci nexiste dej`a plus.
72 Cas test : 2D Model Hill Flow
Ce mod`ele surestime dautre part la vitesse au centre de la veine dessais pour des abscisses
eloignees de la colline. Ce mod`ele est `a proscrire pour cette etude.
Les autres mod`eles, quant `a eux, restent proches des mesures pour cette abscisse
ainsi que celles situees plus en amont. On ne retrouve cependant pas le fort gradient
de vitesse pour des distances proches de la paroi (y/h
max
< 1), alors que celui-ci est
pratiquement restitue sur la paroi superieure (y/h
max
> 5). Ceci nest pas consecutif au
maillage, le crit`ere y
+
est situe entre 8 et 600, lestimation la plus faible etant dans la
zone de recirculation.
Les resultats, plutot satisfaisants, sont equivalents. On ne peut pas dans ce cas estimer
lequel des deux mod`eles `a deux equations de transport est le mieux adapte `a cette etude.
Chapitre 6

Ecoulement au passage dun obstacle

prismatique
On se propose dans cette etude de predire numeriquement le champ aerodynamique
dun ecoulement au passage dun obstacle prismatique, ici un pave de largeur 10 cm et
de hauteur 6 cm (Fig. 6.1). Les calculs seront realises en instationnaires an de mettre en
evidence correctement les instabilites developpees par le decollement des couches limites
laterales. Une des parties importantes ici est lestimation de la valeur du pas de temps de
lapproche instationnaire. Lecoulement est caracterise par un nombre de Reynolds Re =
1,210
5
base sur la largeur du pave.
Pave
x
y
z
O
V

Fig. 6.1 Conguration du calcul visant `a predire le champ aerodynamique au passage dun pave.
74

Ecoulement au passage dun obstacle prismatique
6.1 Approche numerique
On presente dans cette section les maillages sur lesquels les calculs numeriques seront
eectues, ainsi que les mod`eles de turbulence.
6.1.1 Domaine de calcul et maillages
Le domaine de calcul est divise en deux parties, une region proche du pave accueillant
un maillage structure, et une deuxi`eme region caracterisee par un maillage non structure
(Fig. 6.2). La premi`ere region setend sur une longueur 10H par rapport `a la face aval du
cube, et 3H par rapport `a la face amont du cube. Lextension transversale vaut 3H par rap-
port aux faces laterales, et 3H pour lextension verticale. Les extensions supplementaires
de la deuxi`eme region valent respectivement 10H, 5H, 5H, et 5H par rapport aux types
de surfaces enoncees plus haut.
Fig. 6.2 Domaines
detude
Le domaine est maille de fa con structuree dans la premi`ere region `a partir du logiciel
ICEM-CFD. Le maillage est caracterise par le crti`ere y
+
< 1 pour les volumes ayant
en commun une face avec le cube assurant une taille de maille y
m
= 1, 2 10
5
m. La
technique O-Grid permet de respecter ce crit`ere tout en limitant le nombre de nuds total
(Fig. 6.3). La deuxi`eme region est maillee de fa con non-structuree `a partir de CFX-Mesh.
La continuite des tailles de maille est `a peu pr`es assuree au niveau des interfaces des deux
domaines.
6.1.2 Conditions limites
On impose aux nuds dentree une vitesse longitudinale de 18,5 m/s et une intensite
turbulente valant 1%. Les nuds de la sortie accueillent une condition de turbulence
pleinement developpee (i.e. pression statique nulle). La localisation de cette condition est
placee susamment loin pour satisfaire ce crit`ere. Une condition dadherence est imposee
aux nuds du pave et du sol. Une condition de symetrie est imposee sur les surfaces
laterales en sortie de domaine. Des conditions dinterface sont utilisees pour les nuds
places sur les faces communes aux deux domaines.
Approche numerique 75
Fig. 6.3 Maillage du
sol des deux regions
(haut) et de la region
proche du pave utili-
sant la technique O-
Grid (bas)
6.1.3 Mod`eles de turbulence
Les equations de Navier-Stokes moyennees sont fermees `a laide du mod`ele de turbu-
lence SST k . Ce mod`ele bas-Reynolds necessite un maillage optimise dans les regions
proches des parois. Ceci est le cas a priori dans la region du cube. Une loi de paroi est
necessaire neanmoins pour le sol des deux domaines.
Le terme de convection des equations de mouvement est modelise `a partir du schema
du second ordre High Resolution lorsque la taille de maille est ne, et dordre 1 pour des
tailles importantes.
6.1.4 Discretisation temporelle
La valeur du pas de temps du calcul instationnaire est basee sur la periode des lachers
de tourbillon. Dans le cas dobstacle cylindrique circulaire de diam`etre D, cas le plus
etudie, le Strouhal est de lordre de St = fD/V

= 0,2. Cette relation assure dans notre

cas une frequence valant f = 46,25 Hz avec D H = 10 cm, et donc une periode T =
2,210
2
s. On pose dans notre cas une valeur de pas temps egal au cinquanti`eme de la
periode, soit t = T/50 = 4,410
4
s. On sattend neanmoins `a obtenir un nombre de
Strouhal plus faible consecutivement `a la presence de bord saillant.
La convergence du calcul est veriee par levolution des residus des param`etres calcules
sur lensemble des volumes nis des deux domaines, et `a laide de levolution de la tranee
et de la portance (horizontale).
76

Ecoulement au passage dun obstacle prismatique
6.2 Resultats numeriques
6.2.1 Validation du maillage sur le pave
Les valeurs du crit`ere y
+
sont inferieures `a 3 sur la surface du pave. Lutilisation
du mod`ele de turbulence SST k est donc optimale, le code nutilisant pas de loi
de paroi dans ce cas. Les valeurs de ce crit`ere au niveau du sol sont plus etendues. Les
valeurs maximales observees valent 111 dans les regions dacceleration de lecoulement
(par exemple dans la region proche-amont du pave).
6.2.2 Eorts aerodynamiques
Les coecients de portance et tranee valent en moyenne respectivement 0 et 1,1, va-
leurs en accord avec les observations habituelles. On observe dautre part des uctuations
des eorts aerodynamiques consecutivement aux lachers dinstabilites par lobstacle (Fig.
6.4). Le nombre de Strouhal, evalue `a partir des variations temporelles de la portance
et base sur la largeur du pave, vaut St = 0,12. Cette valeur est coherente avec celles
communement observees pour des obstacles non-arrondis (i.e. 0,1 < St < 0,2). Le spectre
associe `a la tranee est caracterise par des frequences plus elevees comparativement `a
celui associe `a la portance. Ceci doit etre du aux eets des instabilites generees par le
decollement de la couche limite inuen cant grandement la tranee.
Fig. 6.4

Evolution
temporelle des eorts
aerodynamiques.
La valeur non-nulle de leort de tranee est lie `a la forte deceleration de lecoulement
dans la region proche-amont du pave entrainant une augmentation de la pression (Fig.
6.5).
6.2.3 Sillage du pave
Le sillage du pave est fortement agite ayant pour consequence des uctuations de
pression parietale importantes dans cette region. Les points places au sol en-dessous dune
instabilites sont caracterises par une depression importante, celle-ci se depla cant ensuite
Conclusion 77
(a) (b)
Fig. 6.5 Champ de pression dans le plan symetrique (a) et dans un plan place `a mi-hauteur du pave (b).
leg`erement dans des regions plus en aval. Dautre part des surpressions quasiment xes
sont observees `a une distance de lordre de 2H en aval de lobstacle.
(a) (b)
Fig. 6.6 Pression parietale dans le sillage du pave (a, t
o
; b, t
o
+ 40t).
Les structures sont observables `a partir du param`etre second invariant, Q :
Q =
1
2
_
S
ij
S
ij

ij

ij
_
(6.1)
avec S
ij
le tenseur de cisaillement et
ij
la partie antisymetrique du tenseur de gradient
de vitesse. On observe dapr`es ce param`etre un tourbillon etire dans la region aval au
pave, et des instabilites convectees ensuite par lecoulement principal (Fig. 6.7).
6.3 Conclusion
Le champ aerodynamique dun ecoulement au passage dun obstacle prismatique a
ete obtenu dans cette etude. La valeur du pas de temps base sur lestimation a priori
78

Ecoulement au passage dun obstacle prismatique
Fig. 6.7 Valeur instantanee de
linvariant Q.
du Strouhal parait adaptee `a lobservation des instabilites generees par le decollement
des couches limites des surfaces laterales du pave. La distribution de la pression sur les
surfaces du pave gen`ere une tranee non-nulle et des eorts de portance uctuant autour
dune valeur nulle. Le sillage est fortement instable et est caracterise par une distribution
de pressions parietales inhomog`ene.
Chapitre 7

Etude des param`etres turbulents dans

un canal
La qualite des predictions decoulements dans un canal repose en partie sur les condi-
tions limites dentree. La methode la plus realiste consiste `a imposer aux nuds de
lentree des prols experimentaux pour les variables moyennes de lecoulement (vitesses et
temperature par exemple), ainsi que pour les param`etres associes `a la turbulence (energie
cinetique turbulente, taux de dissipation, dissipation specique, etc.). Lorsque ces donnees
experimentales sont absentes, il est neanmoins possible dimposer des prols estimes soit
`a partir dun calcul dt annexe, soit `a partir dune formulation analytique.
Le but de cette etude est dobtenir les prols dentree dun ecoulement dans un canal
plan, et de les comparer `a ceux proposes par une formulation analytique. On consid`ere un
canal plan parcouru par un ecoulement dair de vitesse moyenne egale `a 3 m/s (Figure 7.1).
7.1 Approche numerique du calcul annexe
Le calcul annexe est base sur une partie du canal et ne necessite pas la modelisation
de sa longueur totale. On impose aux nuds des surfaces une condition dadherence, et
sur les nuds dentree et de sortie une condition de periodicite. Lecoulement sera
impose par un debit massique par unite de profondeur (le calcul est eectue en 2D) egal
`a 0,735 kg/m.s.
Une fois le calcul eectue par le code, on sinteressera aux prols de la vitesse longi-
tudinale moyenne, denergie cinetique turbulente, et du taux de dissipation. Ces prols
ainsi estimes pourront etre imposes comme condition limite `a lentree dun domaine de
meme dimension, et caracterise par un debit massique identique.
80

Etude des param` etres turbulents dans un canal

V
L
H
Fig. 7.1 Conguration de lecoulement dans le canal.
Le domaine nest pas necessairement long car les prols vont etre identiques quel que
soit lemplacement longitudinal. Dans notre cas, la hauteur vaut 20 cm et la longueur 5 m.
Le maillage est structure avec des elements alignes dans le sens de lecoulement (Figure
7.2). La taille des mailles proches des parois est choisie a priori pour satisfaire une valeur
du y
+
proche de 30. Ceci est realise `a partir de la relation
y
+
= y
u
f

(7.1)
avec une vitesse de frottement, u
f
, supposee etre egale `a 5% de la vitesse exterieure V . On
obtient ainsi la taille de la premi`ere maille dans la region de paroi. Le nombre delement
dans la direction longitudinale est quant `a lui egal `a 100.
Fig. 7.2 Maillage de
letude.
Les equations de Navier-Stokes (incompressible et stationnaire) sont resolues sur le
maillage cite precedemment `a partir du code commercial FLUENT. Toutes les equations
sont resolues au centres des volumes de controle et les valeurs aux fronti`eres de ces volumes
sont estimees par interpolation avec un schema au premier ordre. Le couplage pression-
vitesse est estime au moyen de lalgorithme SIMPLE. La fermeture des equations est
obtenue au moyen dun mod`ele de turbulence k- standard associe `a une loi de paroi
du type standard. Cette derni`ere permet de calculer les caracteristiques de lecoulement
se situant entre la paroi (o` u les eets visqueux sont preponderants) et la premi`ere maille
normale `a la paroi. Enn, le calcul est converge lorsque la somme des residus est inferieure
`a 10
3
pour lensemble des param`etres. La convergence est dautre part veriee par la
stabilite du debit masse en sortie du domaine, et celle de lenergie cinetique turbulente
estimee sur des marqueurs places dans la couche limite superieure.
Resultats numeriques 81
7.2 Resultats numeriques
An de visualiser correctement levolution longitudinale des prols de vitesse longitu-
dinale, denergie cinetique turbulente, et du taux de dissipation, on consid`ere des lignes
verticales placees `a dierentes distances de lentree 0H, 5H, 10H, 15H, 20H, et 25H
appele dans cet ordre Ligne1 `a Ligne6.
Les prols de la vitesse longitudinale, de lenergie cinetique turbulente, et du taux
de dissipation sont identiques, quelle soit la ligne dobservation (Fig. 7.3 et 7.4). Ils cor-
respondent alors `a un prol parfaitement developpe. Levolution transversale de lenergie
cinetique turbulente montre une decroissance quasi lineaire lorsque le point dobservation
seloigne de la paroi. Le taux de dissipation est caracterise, quant `a lui, par une chute
brutale de sa valeur pour atteindre une valeur constante dans une grande majorite de la
conduite.
(a)
(b)
Fig. 7.3

Evolution de la vitesse longitudinale : a, vecteurs ; b, evolution sur les lignes de reference.
82

Etude des param` etres turbulents dans un canal
(a)
(b)
Fig. 7.4

Evolution longitudinale des param`etres turbulents : a, energie cinetique turbulente ; b, taux de dis-
sipation.
7.3 Approximation analytique des grandeurs turbu-
lentes
Nous proposons dans cette partie une formulation de levolution transverse de lenergie
cinetique et du taux de dissipation dans le cas dun ecoulement dans une conduite. Avant
cela, nous nous interesserons `a levolution du prol de la vitesse longitudinale qui peut
etre approchee par la loi en puissance 1/7 :
u(y) = u
CL
_
y

_
1/7
o` u u
CL
est la vitesse au centre du canal, et la demi-hauteur du canal. Les vitesses
longitudinales numeriques et analytiques sont proches (Figure 7.5).
Connaissant la forme des evolutions de ces param`etres dans la region proche de la
paroi (loi de paroi), on suppose que lenergie cinetique turbulente `a lentree varie de la
valeur parietale
k
nw
=
u
2
f
C
1/2

Approximation analytique des grandeurs turbulentes 83

Fig. 7.5 Vitesse
longitudinale (,
resultats numeriques ;
, resultats analy-
tiques).
jusqu`a une valeur au centre du canal
k
CL
=
3
2
(0, 01u
CL
)
2
avec u
f
la vitesse de frottement, C

une constante egale `a 0,09. La vitesse de frottement

est exprimee `a partir de la contrainte parietale et de la masse volumique, u
f
=
_

p
/.
La contrainte parietale peut etre estimee `a partir du mod`ele du coecient de frottement
de Blasius dans les conduites bidimensionnelles :

p
= f
u
2
CL
2
, f = 0, 045
_
u
CL

_
2
Levolution transversale du taux de dissipation peut ensuite etre calculee `a partir de la
valeur locale de lenergie cinetique turbulente et de la longueur de melange l :

i
=
C
3/4

k
3/2
i
l
On pourra estimer la longueur de melange `a partir du mod`ele propose par Prantdl :
l = min(y, 0.085)
avec = 0,41. Dans notre cas, les dierentes valeurs sont rassemblees dans le tableau 7.1.
Les valeurs numeriques et analytiques sont comparees sur la gure 7.6. Lenergie cinetique
turbulente analytiques sont sous-estimees dans le centre du canal, et sont correctement
estimees proche de la paroi, que les taux de dissipation sont comparables pour toute
position transversale.
Tableau 7.1 Estimations des param`etres lies `a la turbulence dans le canal.
V
x
(m/s) u
CL
u
f
k
nw
k
CL
15 16,4 0,589 1,16 0,04
84

Etude des param` etres turbulents dans un canal
(a)
(b)
Fig. 7.6 Approximation des param`etres turbulentes. a, energie cinetique turbulente ; b, taux de dissipation
(, resultats numeriques ; , resultats analytiques).
7.4 Conclusion
Lecoulement dans un canal plan nous a permis de valider des formulations de levolution
transversale de lenergie cinetique turbulente et du taux de dissipation. Celles-ci pourront
ainsi etre utilisees pour simuler de fa con la plus realiste possible un ecoulement amont `a
un obstacle situe dans un canal. La fonction permettant de retrouver les valeurs locales
de lenergie cinetique turbulente pourra etre de la forme :
k(y
m
) = k
nw
+
k
CL
k
nw

y
m
Chapitre 8
Prediction numerique dun echange
thermique convectif `a partir de
dierents mod`eles de turbulence
On etudie ici les echanges thermiques convectifs entre un ecoulement de vitesse V

ayant une temperature T

et un obstacle elance dont la temperature surfacique est notee

T

(Fig. 8.1). Le but de cette etude est de comparer les resultats obtenus `a partir des
mod`eles de turbulence k standard, k Realizable, et SST k . Les resultats
seront atteints `a laide de deux maillages dierents, et lensemble des resultats numeriques
pourront etre compares `a des resultats experimentaux obtenus dans des congurations
identiques.
V

, T

T
p
3 mm
L = {10, 20, 40}
100 mm
Fig. 8.1 Conguration de letude de lechange thermique convectif (dapr`es Marchesse [22]).
Les resultats numeriques et experimentaux seront compares sur dierentes congura-
tions dans lesquelles les largeurs de lobstacle et vitesse de lecoulement pourront varier
(Tab. 8.1), alors que lepaisseur reste xe et egale `a 3 mm, et la temperature dentree
egale `a 383 K quelle que soit la vitesse.
86 Pr ediction num erique dun echange thermique convectif ` a partir de
diff erents mod` eles de turbulence
Tableau 8.1 Valeurs des param`etres variables.
V

(m/s) 3 7 9
L (mm) 10 20 40
8.1 Approche numerique
On presente dans cette section les maillages sur lesquels les calculs numeriques seront
eectues, ainsi que les mod`eles de turbulence.
8.1.1 Domaine de calcul et maillages
Le calcul est eectue sur la moitie du domaine reel (Fig. 8.2), les ecarts des resultats
obtenus sur le domaine entier et la moitie du domaine etant tr`es faibles.
Fig. 8.2 Domaine
du calcul et conditions
aux limites.
Les equations de mouvement et de lenergie sont resolues ensemble sur deux maillages,
les mailles M1 et M2 caracterises respectivement par y
+
= 1 et y
+
= 0,1. La taille de
la premi`ere maille, notees y
M1
et y
M2
, depend ainsi du maillage considere, et est dautre
part dautant plus faible que la vitesse de lecoulement est importante. Les valeurs de la
taille de cette premi`ere maille sont rassemblees dans le tableau 8.2.
Tableau 8.2 Valeurs des param`etres variables.
V

(m/s) y
M1
(m) y
M2
(m)
3 2, 1 10
4
2, 1 10
5
7 9, 0 10
5
9, 0 10
6
9 7, 0 10
5
7, 0 10
6
Tous les maillages sont inhomog`enes et structures avec des elements alignes selon
lecoulement.
8.1.2 Conditions limites
Les conditions limites sont rassemblees sur la gure 8.2. On impose aux nuds dentree
des prols de vitesse longitudinale moyenne et des param`etres associes `a la turbulence qui
dependent alors du mod`ele de turbulence utilise. Ces prols sont obtenus `a partir dun
calcul annexe (Voir le chapitre 7 pour plus de precision). Les temperatures dentree et de
lobstacle valent respectivement 293 K et 383 K, assurant ainsi un echange thermique. Une
condition de symetrie est imposee sur les lignes inferieures et une condition de turbulence
pleinement developpee en sortie de domaine. Les nuds de la ligne superieure accueillent
une condition dadherence et de ux thermique nul (paroi adiabatique).
Resultats numeriques 87
8.1.3 Mod`eles de turbulence
Des mod`eles de turbulence sont necessaires pour fermer le syst`eme dequation RANS.
Trois mod`eles sont etudies ici, les mod`eles k standard, k Realizable, et SST k, o` u
k est lenergie cinetique turbulente, le taux de dissipation, la dissipation specique. Les
deux premiers mod`eles sont valables uniquement dans le cur de lecoulement. Une loi de
paroi amelioree sera donc utilisee pour estimer les caracteristiques de lecoulement dans les
regions proches des parois. Dans cette loi, la region de paroi est scindee en deux regions,
une sous-couche visqueuse et une sous-couche pleinement turbulente, dont la demarcation
seectue `a partir dun nombre de Reynolds de turbulence base sur la distance normale `a
la paroi :
Re
y
=
y

Lutilisation de cette loi necessite un y

+
egal ou inferieur `a lunite.
Lorsque le mod`ele SST k est utilise, aucune loi de paroi est necessaire car ce
mod`ele est adapte aux deux types de regions. Cependant le maillage doit etre construit
an de resoudre lecoulement dans la sous-couche visqueuse (i.e. y
+
< 1).
Ainsi dans tous les cas, la valeur du param`etre y
+
doit etre egale ou inferieure `a lunite.
Les deux maillages M1 et M2 satisfont a priori ce crit`ere.
8.2 Resultats numeriques
Les resultats numeriques et experimentaux sont traces sur les gures 8.3 `a 8.5. On
observe que les valeurs des coecients dechange numeriques sont tr`es proches de celles
mesurees lorsque la largeur de lobstacle est egale `a 40 mm quel que soit le mod`ele. Pour
des valeurs plus faibles de ce param`etre, on note que les mod`eles k surestiment lechange
thermique. Ceci est du `a une surestimation de lenergie cinetique turbulente dans la region
proche du point darret (Fig. 8.6), conduisant alors `a un echange thermique excessif (lecart
avec les mesures vaut 40%). Les resultats numeriques bases sur le mod`ele SST k sont
tr`es proches des valeurs experimentales, lecart nexc`edant pas 15%. Notons dautre part
que les predictions numeriques de ce mod`ele, obtenues cette fois sur le maillage tr`es n
M2, sont ameliorees, les ecarts chutant `a 3%.
Le but de cette etude etait de comparer les echanges thermiques convectifs obte-
nus `a partir de dierents mod`eles de turbulence et calcules sur deux maillages, `a des
resultats experimentaux atteints dans des congurations identiques. Nous avons observes
des dierences notables entre ces resultats. Et seul le mod`ele SST k est `a meme
de predire des echanges convectifs realistes, et ceci tout specialement lorsque le maillage
est tr`es n (y
+
= 0,1). Le mod`ele k surestime de fa con trop importante les echanges
consecutivement `a une mauvaise prediction de lenergie cinetique turbulente dans la region
proche du point darret.
88 Pr ediction num erique dun echange thermique convectif ` a partir de
diff erents mod` eles de turbulence
Fig. 8.3 Estima-
tion du coecient
dechange convectif
(V

= 3 m/s).
Dapr`es Marchesse
[22].
Fig. 8.4 Estima-
tion du coecient
dechange convectif
(V

= 7 m/s).
Dapr`es Marchesse
[22].
Resultats numeriques 89
Fig. 8.5 Estima-
tion du coecient
dechange convectif
(V

= 9 m/s).
Dapr`es Marchesse
[22].
Fig. 8.6 Distribu-
tion du coecient
dechange convectif
le long de lobstacle
de largeur egale `a 40
mm (V

= 3 m/s).
Dapr`es Marchesse
[22].
90 Pr ediction num erique dun echange thermique convectif ` a partir de
diff erents mod` eles de turbulence
Chapitre 9
Inuence de la nesse du maillage sur
la prediction numerique des echanges
thermiques convectifs dune plaque
plane
Le but de cette etude est de mettre en evidence linuence du maillage sur la qualite de
prediction des echanges thermiques entre un ecoulement et une plaque plane non-inclinee
lorsque mod`ele k standard est utilise. Ce mod`ele netant pas valable dans la region
proche des parois, lutilisation dune loi de paroi permet la prediction des caracteristiques
de lecoulement dans cette region. Cependant, la taille de la premi`ere maille doit etre
en adequation avec cette loi. Nous verrons dans cette etude que lutilisation dune loi de
paroi sur un maillage inadapte conduit `a des estimations dechanges thermiques errones.
On consid`ere ici un ecoulement de vitesse
V

egale `a 1 m/s et de temperature T

egale `a 353 K au passage dune plaque plane

de longueur 1 m sans inclinaison ayant une
temperature egale `a 413 K.
(V

,T

)
L = 1 m
T
p
= 413 K
9.1 Approche numerique
9.1.1 Domaine et maillage
Le domaine est construit `a partir dun carre avec la ligne inferieure representant la
plaque plane. De ce fait, la couche limite se developpe d`es lentree du domaine. Les
92 Influence de la finesse du maillage sur la pr ediction num erique des
echanges thermiques convectifs dune plaque plane
equations de mouvement et thermiques sont resolues `a partir du code commercial FLO-
TRAN sur trois maillages (M1, M2, M3) structures et inhomog`enes (Fig. 9.1). Ces maillages
di`erent par la taille des elements proches de la paroi, les elements du maillage M1 sont
beaucoup plus importants que ceux du maillage M3.
Fig. 9.1 Maillages de letude numerique des echanges thermiques entre un ecoulement et une plaque plane.
De la gauche vers la droite : M1, M2, M3.
Le syst`eme dequations RANS est ferme `a partir du mod`ele de turbulence k stan-
dard, et la modelisation de lecoulement dans la couche limite est eectuee `a laide de la
loi de paroi equilibree necessitant une valeur y
+
proche de 30.
9.1.2 Conditions limites
On impose aux nuds de lentree une vitesse longitudinale egale `a 1 m/s et une
temperature egale `a 353 K. Les nuds de la plaque imposent imposent une condition
dadherence et une temperature egale `a 413 K. On impose une condition de symetrie sur
la ligne superieure et une condition de turbulence pleinement developpe en sortie.
9.2 Mod`ele propose
Les echanges thermiques peuvent etre caracterises par un nombre de Nusselt local
(Nu
x
) en fonction du nombre de Reynolds local (Re
x
) :
Nu
x
= 0, 0226 Re
0,8
x
Pr
0,6
_
T

T
p
_
0,4
o` u Pr est le nombre de Prandtl. Ce mod`ele est valable dans la gamme de Reynolds [10
5
;
10
7
]. Ce mod`ele est choisi parmi dautres disponibles car les nombres de Reynolds etant
eleves ceci nous assure la presence dune couche limite tr`es developpee. Lutilisation de la
loi de paroi dans de bonnes conditions (i.e. y
+
proche de 30) sera dautant meilleure.
9.3 Resultats numeriques
Les resultats obtenus `a partir des trois maillages sont rassembles sur la gure 9.2. On
note que les predictions sont tr`es dierentes. Le maillage le plus n dans la region proche
Resultats numeriques 93
de la paroi (i.e. M3) entrane une forte surestimation de lechange consecutivement `a des
valeurs de y
+
tr`es faibles et non-adaptees `a lutilisation de la loi de paroi. Lutilisation
dun maillage grossier (i.e. M1) ne permet pas quant `a lui destimer correctement lechange
thermique, y
+
etant proches de 300 le long de la plaque. Cette valeur nest pas totalement
inadaptee pour lutilisation de la loi de paroi equilibree, mais nest pas optimisee. Ceci
est obtenu `a partir du maillage M2 pour lequel les valeurs du param`etres de maillage
est leg`erement inferieures `a 30 entranant d`es lors une estimation numerique de lechange
tr`es satisfaisante. On notera dautre part que les evolutions du Nusselt pour les faibles
Reynolds sont dierentes pour M1 et M2, les courbes etant soient concaves ou convexes.
(a) (b)
Fig. 9.2 Estimation numerique de lechange thermique entre un ecoulement et une plaque plane.

Evolution
du Nusselt local (a) et du crit`ere de maillage (b).
94 Influence de la finesse du maillage sur la pr ediction num erique des
echanges thermiques convectifs dune plaque plane
Bibliographie
[1] G. P. Almeida, D. F. G. Durao, and M. V. Heitor. Wake ows behind two
dimensional model hills. Exp. Thermal and uid Science, 7 :87, 1992.
[2] J. Boussinesq. Theorie de lecoulement tourbillonnant, 1877. Mem. Pres. Acad.
Sci. XXIII.
[3] G. L. Brown and A. Roshko. On density eects and large structure in turbulent
mixing layers. J. of Fluid Mechanics, 64(4) :775816, 1974.
[4] C. Changenet, Y. Marchesse, and P. Guarat. Numerical investigation for pre-
dicting drag torque and heat transfer in wet disengaged clutches, 2005. International
congress of uid dynamics applications in ground transportation (Lyon).
[5] D. Choudhury. Introduction to the renormalization group method and turbulence
modeling, 1993. Fluent Inc. Technical Memorandum TM-107.
[6] B. J. Daly and F. H. Harlow. Transport Equations in Turbulence. Phys. Fluids,
13 :26342649, 1970.
[7] J. B. Freund and et P. Moin S. K. Lele. Numerical Simulation of a Mach 1.92
Turbulent Jet and Its Sound Field. AIAA Journal, 38 (11) :20232031, 2000.
[8] M. M. Gibson and B. E. Launder. Ground Eects on Pressure Fluctuations in the
Atmospheric Boundary Layer. J. Fluid Mech, 86 :491511, 1978.
[9] G. Iaccarino, A. Ooi, P. A. Durbin, and M. Behnia. Reynolds averaged simula-
tion of unsteady separated ow. Int. J. Heat and Fluid Flow, 24 :147156, 2003.
[10] W. P. Jones and B. E. launder. The prediction of laminarization with a two-
equation model of turbulence. Int. J. Heat Mass Transfer, 15 :301314, 1974.
[11] J. D. Anderson Jr. Computational Fluids Dynamics - The basics with applications.
1995. Edition MacGraw Hill.
[12] J. Kim, P. Moin, and R. Moser. Turbulent statistics in fully developed channel ow
at low Reynolds number. J. of uid mechanics, 177 :133166, 1987.
[13] A. N. Kolmogorov. Curves in Hilbert space invariant with respect to a one-
parameter group of motions. Dokl. Akad. Nauk SSSR, 26 (1) :69, 1940.
[14] A. N. Kolmogorov. Local structure of turbulence in an incompressible uid at very
high Reynolds number. Dokl. Akad. Nauk SSSR, 30 (299), 1940.
[15] L. Landau and E. Lifchitz. Mecanique des uides. 1971. Editions MIR.
[16] B. E. Launder. Second-Moment Closure : Present... and Future ? Inter. J. Heat
Fluid Flow, 10(4) :282300, 1989.
[17] B. E. Launder, G. J. Reece, and W. Rodi. Progress in the Development of a
Reynolds-Stress Turbulence Closure. J. Fluid Mech., 68(3) :537566, 1975.
96 BIBLIOGRAPHIE
[18] B. E. Launder and B. I. Sharma. Application of the energy dissipation model of
turbulence to the calculation of ow near a spinning disc. Letters in Heat and Mass
Transfer, 1 :131138, 1974.
[19] B. E. Launder and D. B. Spalding. The Numerical Computation of Turbulent
Flows. Computer Methods in Applied Mechanics and Engineering, 3 :269289, 1974.
[20] M. Lesieur. La turbulence. DEP Science, 1994.
[21] J.R. Mansfield, O.M. Knio, and C. Meneveau. Dynamic LES of colliding vortex
rings using a 3D vortex method. J Comp Phys, 152 :305345, 1999.
[22] Y. Marchesse and C. Changenet. Forced convective heat transfer over a noncir-
cular slender cylinder. Proceedings of the institution of mechanical engineers, Part
C. Journal of mechanical engineering science., 223(C2) :427437, 2008.
[23] Y. Marchesse, C. Changenet, F. Ville, and P. Velex. Investigations on CFD
simulation for predicting windage power losses in spur gears, 2009. IRF2009 - 3rd
International Conference on Integrity, Reliability, and Failure, Porto (Portugal), 20-
24 Juillet 2009.
[24] R. Martinuzzi and C. Tropea. The ow around surface mounted, prismatic obs-
tacles placed in a fully developped channel ow. Trans. ASME J. Fluid. Eng., 115 :85
93, 1996.
[25] F.R. Menter. Two-equation eddy viscosity models for engineering applications.
A.I.A.A. Journal, 32 :15981605, 1994.
[26] R. D. Moser and et N. Mansour J. Kim. Direct Numerical Simulation of turbulent
channel ow up to Re

= 590. Phys. Fluids, 11 :943946, 1999.

[27] S. B. Pope. Ten questions concerning the large-eddy simulation of turbulent ows.
New Journal of Physics, 6 (35), 2004.
[28] Z. A. Prandtl, 1925. Zs Angew. Math. Mech.
[29] O. Reynolds. An experimental investigation of the circumstances which determine
whether the motion of water shall be direct and sinous, and the law of resistance in
parallel channel. Phil. Trans. Roy. Soc., pages 51105, 1883.
[30] L. F. Richardson. Weather prediction by numerical process, 1922. Cambridge
University Press.
[31] W. Rodi. Experience with two-layer models combining the k with a one-equation
model near the wall. A.I.A.A. Paper, 91-0216, 1991.
[32] H. Schlichting. Boundary-Layer Theory - 7`eme edition. 1960. Edition MacGraw
Hill.
[33] K. B. Shah. Large Eddy Simulation of the ow past a cubic obstacle. Th`ese, 1998.
[34] S. R. Spalart and S. A. Allmaras. A one-equation turbulence model for aero-
dynamic ows, 1992. Technical Report AIAA-92-0439, et aussi dans La Recherche
Aerospatiale, 1(94), pp. 5-21.
[35] H. K. Versteeg and W. Malalasekera. An introduction to computational uid
dynamics - The nite volume method -2nd edition. Pearson Education Limited, 2007.
[36] T. Wei and W. W. Willmarth. Reynolds number eects on the structure of a
turbulent channel ow. J. Fluid Mech, 204 :5795, 1989.
[37] D.C Wilcox. Reassessment of the scale determining equation for advanced turbu-
lence models. A.I.A.A. Journal, 26 (136) :12991310, 1988.
Lectures interessantes
Ouvrages
Computational Fluid Dynamics : The Basics with Applications
John D. Anderson Jr, McGraw-hill International Edition, 1995 ;
Computational Fluid Dynamics : Principles and Applications
J. Blazek, Ed. Elsevier ;
Computational Methods for Fluid Dynamics
Joel H. Ferziger & Milovan Peric, Ed. Springer (3`eme edition) ;
Les ecoulements turbulents, modelisation et simulation
Roland Schiestler,

Edition Hermes ;
Fundamentals of Aerodynamics
John D. Anderson Jr, McGraw-hill International Edition, 2001 ;
Turbulence
Christophe Bailly & Genevi`eve Comte-Bellot, CNRS

Edition;
La turbulence
Marcel Lesieur, Collection Grenoble Sciences, 1994 ;
La turbulence developpee
Marcel Lesieur, La Recherche, n139, pp 1412-1425, decembre 1982 ;
Turbulence en mecanique des uides, Analyse du phenom`ene en vue de
sa modelisation `a lusage de lingenieur
Patrick Chassaing, Cepadu`es-

Editions.
Cours
Modelisation aux moyennes de Reynolds ou Modelisation RANS, Seconde
partie du cours Turbulence et modelisation de la turbulence
Remi Manceau, Master Sciences et technologies - Poitiers, 2004 ;
Sites internet
http ://[Link]/
http ://[Link]/marutgers/