Etude de la convection naturelle dans une
cavité par une méthode de réanalyse
Jean Félix Durastanti, Raouf Khelalfa, Youssef Sfaxi
IUT de Sénart- Université de PARIS EST
Département Génie Industriel et Maintenance
Avenue Pierre Point - 77567 Lieusaint cedex
durastanti@[Link] ; [Link]@[Link]; sfaxi@[Link]
Section de rattachement : 62
Secteur : Secondaire
RESUME. On s’intéresse dans cet article à la convection naturelle en régime stationnaire dans
une cavité carrée, pour laquelle l’influence du nombre de Rayleigh est déterminante. On traite le
problème par une méthode numérique originale basée sur la réanalyse des équations de Navier
Stokes. Le principe de cette méthode est de découpler les équations phénoménologiques des
conditions aux limites par l’intermédiaire de variables d’interfaces afin d’obtenir un système de
taille réduite permettant de traiter facilement une modification des conditions aux limites. Le
nombre d’opérations nécessaire se trouve alors considérablement diminué par rapport à une
résolution classique.
MOTS-CLES : Convection naturelle, nombre de Rayleigh, réanalyse, Navier Stokes
1. Introduction
Après avoir validé la méthode de réanalyse aux transferts de chaleur (Peureux 1995),
(Peureux et al. 1995), nous présentons ici son application aux équations de Navier
Stokes, dont la particularité est d’introduire une non-linéarité pour l’obtention du champ
des vitesses inconnu. On rappelle qu’elle consiste à formuler une solution générale à un
problème puis de particulariser cette solution en prenant en compte les conditions aux
limites (Rabearivelo 2005). En fait il s’agit de découpler lors de la modélisation du
système physique (V) de surface (S), la phénoménologie du problème et les conditions
aux limites.
1
2. Réanalyse des équations de Navier Stokes
Dans un milieu continu incompressible de domaine volumique (V) et de frontière
(S), l’évolution dynamique et thermique est régie par les équations classiques de Navier
Stokes :
∇ .V = 0
DV
ρ = − ∇p + ∇.( µ ∇V ) + ρ f
[1]
Dt
DT Dp
ρ cp − ∇.( k V T ) = β T + 2 µ ( D : D)
Dt Dt
1 ∂ρ
avec β=−
ρ ∂ T
p
Le traitement du problème stationnaire [1], par la méthode des résidus
pondérés conduit classiquement à la résolution d’un système de taille importante de la
forme :
[ K ]{ Yn } − { F } = 0 [2]
avec : [K ]= ∑ [K
e
e
] { F } = ∑ { F e}
e
{ }
Y n est le vecteur des valeurs discrètes de Y aux différents nœuds de l’élément
obtenu après discrétisation du problème consistant à diviser le milieu continu ( V ) en
élémentaires ( Ve ) de frontière ( Se ) et à évaluer
un certain nombre d’éléments finis
chaque grandeur physique de Y et ses variations sur la base d’une approximation
nodale linéaire.
e
La matrice [ K ] représente la matrice élémentaire équivalente à la matrice des
e
conductances en thermique. Le vecteur { F } correspond aux sollicitations sur un
élément volumique ( Ve ).
La matrice [ K ] est décomposée dans un premier temps sous la forme :
[ K ] = [ K V ] + [ KS ] [3]
2
où [ K V ] est liée à l’élément de volume et traduit la phénoménologie du problème.
La matrice [ K S ] correspond quant à elle à l’incidence des conditions aux limites sur
la frontière ( S ) du volume ( V ). Le système [2] devient :
[ K V ]{ Yn } = { FG } avec : { FG } = { F } − [ K S ]{ Yn }
On découple les termes non linéaires de la matrice [K] en posant :
[ KV ] = [ K V l ] + [ KV n l ] [4]
Soit :
[ K V l ]{ Yn } = { FG } − { H G } [5]
avec : { H G } = [ K V nl ]{ Yn }
l
Le système obtenu ci-dessus est singulier puisque la matrice [ K V ] est singulière.
Soit r l’ordre de la singularité et donc n − r le rang du système. Pour que le
problème soit bien posé, p conditions aux limites sont imposées avec p ≥ n sous la
forme linéaire :
[ CL ]{ Y n } = { δ } [6]
avec [ CL ] la matrice des conditions aux limites de taille p × n et { δ } le vecteur
dont les p composantes peuvent être du point de vue physique des vitesses, des
températures, des pressions voire des flux imposés au milieu continu ( V ) à travers sa
frontière ( S ).
Le vecteur { FG } des sollicitations globales comporte des termes connus { Fd } et
des termes inconnus { Fr } au nombre maximal de p.
{ FG } = { Fd } + { Fr } [7]
{ Fr } = [ L ] . { φ } [8]
La matrice [ L ] est obtenue en condensant { Fr } en un vecteur { φ } qui permet de
réduire la taille de n à p sans changer la nature physique des composantes par
suppression des composantes nulles de { Fr } .
3
La matrice [ K V ] est régularisée par une méthode (Loredo 1995) faisant intervenir
[ K α ] qui est symétrique et définie positive en posant :
[ Kα ] = [ K V l ] + α { R }. { R }t [9]
où α est un réel non nul strictement positif et { R } la matrice des vecteurs propres
associés à la valeur propre nulle.
Le traitement des p conditions aux limites permet de former le système de réanalyse non
linéaire de dimension réduite p + r << n
[ CL ] [ K α ] − 1 [ L ] [ CL ]{ R } { φ }
{ }
Z1
= [10]
{ R } [ L ]
t
[ 0] { ω } { }
Z 2
avec :
{ Z } = { δ } − [ CL ] [ K ] { F } + [ CL ][ K
1 α
−1
d α ] − 1 { HG }
{ Z } = − { R} ( { F }+ {H } )
2
t
d G
et d’obtenir la solution :
{ Yn } = [ Kα ] − 1 ( { Fd } + [ L ] {φ } )
[11]
+ { R }{ ω } − [ Kα ] − 1 { HG }
avec : {ω } = { R } t . {Yn }
−1
{R} et [ K α ] sont indépendantes des conditions aux limites et calculées une
−1
seule fois. [ K α ] est obtenue par la méthode de décomposition issue de la
modification des méthodes de Gauss ou de Cholesky (Rabearivelo 2005) afin de
conserver l’avantage du stockage profil malgré la structure pleine. Lors d’une
modification d’une ou plusieurs conditions aux limites, seule la matrice CL est [ ]
modifiée.
4
3. Application à la convection naturelle dans une cavité carrée
3.1. Ecriture des équations
Afin de valider la théorie exposée précédemment, on traite par la réanalyse le cas
connu en mécanique des fluides de l’écoulement de convection naturelle dans une
cavité carrée en régime stationnaire engendré par un faible écart de températures, ce qui
justifie l’hypothèse de Boussinesq.
Considérons une cavité carrée de coté L , remplie d’un fluide newtonien où un gradient
de température TH − TC > 0 est imposé entre les deux parois verticales. Les parois
haute et basse sont adiabatiques ( q = 0 ) (Fig.1).
u= v = 0 q=0
u= v = 0 u= v = 0
L y
T = TH T = TC
u= v = 0 q=0
Fig.1: Convection naturelle dans une cavité
En négligeant partout les variations de la masse volumique sauf dans le terme de
gravité, on obtient en variable sans dimension [12] :
5
x y u v p ρ
=x = y =u =v = p =ρ
L L U0 U0 p0 ρ0
[12]
T − TC µ cp k
= T = 1 = 1 = 1
T H − TC µ0 cp 0 k0
les équations de Boussinesq [13] modélisant l’écoulement :
∂u ∂v
+ =0
∂x ∂y
∂u ∂u ∂P ∂2 u ∂2 u
Re v + u = − + +
∂y ∂x ∂x ∂x2 ∂y2
[13]
∂v ∂v ∂P ∂2 v ∂2 v
Re v + u = − + + +
Gr
( T −T )
∂y ∂x ∂y ∂x ∂y
2 2 0
Re
∂T ∂T 1 ∂2 T ∂2 T
v + u − + = 0
∂y ∂x Pe ∂ x 2 ∂y2
avec les conditions aux limites suivantes :
u ( x ;1) = v ( x ;1) = 0 u (1 ; y ) = v (1 ; y ) = 0 T (0 ; y ) = 1
∂ T ( x ; 1) ∂ T ( x ; 0) [14]
T (1 ; y ) = 0 = = 0
∂y ∂y
On a pris comme échelle respective pour la vitesse et la pression :
α U0
U0 = Ra Pr p0 = µ
L L
telles que :
Pe g β ( TH − TC ) L 3
Pr = Ra = Gr Pr = Pr
Re ν2
6
3.2. Résultats
On a résolu le problème de la cavité [13]-[14] à partir de la méthode des éléments
finis en résolvant le système de réanalyse [10] avec des éléments rectangulaires à 8
nœuds pour des Rayleigh compris entre 103 et 105, pour un nombre de Prandtl
Pr = 0.72 . On compare nos résultats à ceux proposés par différents auteurs (Vahl
davis 1983), (Gresho 1980), (Winters 1980) dans le tableau suivant où sont reportées les
valeurs de umax (maximum de la vitesse horizontale sans dimension sur la section
verticale située au milieu de la cavité) à la hauteur adimensionnelle où cette vitesse est
obtenue ainsi que les valeurs de vmax
Nombre
de Rayleigh
Auteur umax ( 1/ 2 ; y ) vmax ( x ; 1/ 2 )
Gresho 3.656 à y = 0.812 3.704 à x = 0.166
Winters 3.64 à y = 0.81 3.69 à x = 0.1
10 3
Vahl Davis 3.649 à y = 0.813 3.697 à x = 0.178
Réanalyse 3.93 à y = 0.8 3.696 à x = 0.1
19.675 à
Gresho 16.193 à y = 0.822
x = 0.1187
4 Winters 16.2 à y = 0.82 19.7 à x = 0.12
10
Vahl Davis 16.718 à y = 0.823 19.167 à x = 0.119
Réanalyse 16.95 à y = 0.8 20.49 à x = 0.1
68.896 à
Gresho 34.62 à y = 0.856
x = 0.0663
5 Winters 34.8 à y = 0.86 68.6 à x = 0.066
10
Vahl Davis 34.73 à y = 0.855 68.59 à x = 0.066
Réanalyse 39.09 à y = 0.8 65.13 à x = 0.1
Les résultats obtenus permettent de valider la méthode. Les différences observées
sont dues en grande partie à l’option qui a été prise de conserver le même nombre de
mailles compte tenu des faibles nombres de Rayleigh envisagés ici.
4. Conclusion
On a montré dans cet article, qu’une méthode de réanalyse basée sur le découplage
des termes phénoménologique et des conditions aux limites ainsi qu’une prise en
7
compte des non linéarités, était adaptable aux équations de Navier Stokes en régime
stationnaire et conduisait à la résolution d’un système de taille réduite facilitant la prise
en compte d’une modification à la frontière du volume étudié. Les résultats obtenus par
cette méthode pour le cas particulier d’un écoulement en convection naturelle dans une
cavité carrée où un faible gradient de température est imposé entre deux parois
verticales, sont satisfaisants.
Références
Gresho P. M., Upson C. D., Lee R. L., Finite element simulations of thermally induced convection
in enclosed cavity, Lawrence Livermore Laboratory Report UCID-18602 (1980).
Loredo A., Calcul numérique de la quasi-inverse d’une matrice réelle symétrique semi - définie
positive, C. R. Acad. Sci. Paris, t. 321, série I (1995), 247-252.
Peureux B., La réanalyse en thermique, Thèse de doctorat de l’université de Bourgogne (1995).
Peureux B., Durastanti J. F., Martin B., Loredo A., La réanalyse de l’équilibre .thermique, Revue
Générale de Thermique (1995), 679-689.
Rabearivelo P. M., Contribution à la résolution des équations de Navier Stokes par une méthode
de réanalyse, Thèse de doctorat de l’université de Paris XII (2005).
Vahl Davis G., Natural convection in square cavity. A benchmark solution, International journal
for numerical methods in fluids, 3 (1983), 249-264.
Winters K. H., A numerical study of natural convection in a square cavity, United Kingdom
Atomic Energy Authority, AERE-R9747 (1980).