# Simulation numérique bidimensionnelle d'un écoulement turbulent autour de deux cavités

Khouloud Hallek<fnref n="1" />, Mourad Bouterra<fnref n="1" />*, Afif El Cafsi<fnref n="1" />, Ali Belghith<fnref n="1" /> et Patrick Le Quere<fnref n="2" />

(reçu le 29 Septembre 2007 – accepté le 25 Décembre 2007)

**Résumé** - Le but de ce travail est d'étudier l'interaction d'une couche limite avec deux cavités et de caractériser la structure dynamique de l'écoulement. Pour le cas d'une seule cavité large, la structure assure l'existence d'une grosse recirculation où il y a la naissance de deux tourbillons; une petite au pied de la première marche et l'autre plus grande qui couvre toute la cavité avant de s'échapper après recollement. Ce phénomène a été bloqué, dans le deuxième cas, par l'existence de l'obstacle et la création d'une nouvelle cavité. L'étude des profils des vitesses, dans les différentes sections du domaine, montre que le profil d'entrée n'est plus laminaire.

**Abstract** - The aim of this work is to study the interaction of a boundary layer with two cavities and characterize the dynamical structure of the flow. In the case of a single cavity large, the structure ensures the existence of a large recirculation where there is the birth of two vortices: a little at the foot of the first step and the second largest, which covers the cavity before to escape after recollement. This has been blocked, in the second case, by the existence of the barrier and the creation of a new cavity. The study of the profiles speeds, in different sections of the area shows that the profile of entry is more laminar.

**Mots clés:** Simulation – Ecoulement turbulent – Cavité – Obstacle - Profil de vitesse.

## 1. INTRODUCTION

En milieu urbain, la dynamique spatio-temporelle des mouvements de fluide y est très complexe, il s'est produit des mouvements tourbillonnaires, d'échelle très variable. La turbulence qui s'est produite peut avoir plusieurs causes. En effet, La présence des obstacles, leur géométrie et leur orientation sont autant des facteurs qui modifient la dynamique et le niveau de la turbulence. En effet, la présence de bâtiments au voisinage de la canopée urbaine, modifie considérablement le niveau de la turbulence et par conséquent la structure de l'écoulement à l'intérieur de la rue canyon. Les travaux de Baranger [1], montrent l'effet de l'orientation du vent sur la dynamique du mouvement au dessus d'une rue canyon. Oke [2], dans une étude expérimentale en soufflerie portant sur deux blocs parallèles, a mis en évidence trois régimes d'écoulement à l'intérieur du canyon pour différentes géométries. A l'échelle urbaine, l'interaction d'une couche limite avec une cavité constitue un phénomène physique que l'on rencontre fréquemment. Cette interaction engendre des variations des paramètres physiques qui conditionnent la nature de l'écoulement. Ces écoulements sont turbulents, instationnaires et sont à l'origine de la présence des mouvements tourbillonnaires d'échelles spatiales variables. A l'échelle locale, une cavité est une configuration illustrant une rue canyon formée par un ensemble de bâtiment et de rue.

L'objectif essentiel de cette étude est la caractérisation des structures tourbillonnaires prenant naissance lors de l'interaction entre une couche limite avec une cavité, puis avec deux cavités en fonction du nombre de Reynolds.

<footnote n="1">¹ Faculté des Sciences de Tunis, Campus Universitaire 1060 Tunis, Tunisie</footnote>

<footnote n="2">² Laboratoire d'Informatique pour la Mécanique et les Sciences de l'Ingénieur 'LIMSI', Bâtiment 508, B.P. 133, 91403 Orsay Cedex, France</footnote>

<footnote>* mourad.bouterra@fst.rnu.tn</footnote>
## 2. CONFIGURATION GEOMETRIQUE ET CONDITIONS AUX LIMITES

Le domaine d'étude est rectangulaire et caractérisé par son rapport de forme $R = \frac{L}{H} = 4$ où L et H sont respectivement la longueur et la hauteur du domaine. Quant à la première cavité, elle est de section carrée et caractérisée par le rapport de forme $r = \frac{1}{h} = 4$. La seconde cavité est caractérisée par le rapport de forme $r_1 = \frac{l_1}{h} = 1$ (Fig. 1). A l'entrée du domaine, nous imposons un profil de type couche limite laminaire de la forme:

$$U(z) = 1.92 (z/10)^{0.25} \quad \text{pour } 0.25 \le z \le 0.75$$

$$U(z) = 1 \quad \text{pour } 0.75 \le z \le 1$$

<figure><img src="image_1.png" /><figcaption>Fig. 1: Domaine d'étude</figcaption></figure>

## 3. FORMULATION MATHEMATIQUE

L'écoulement turbulent est simulé en résolvant les équations de Navier Stokes. Le fluide est supposé obéir à l'approximation de Boussinesq avec le choix des grandeurs de référence suivantes: $L_{ref} = H$ pour la longueur, $V_{ref} = U_0$ (vitesse maximale à l'entrée du domaine) pour la vitesse et $t = \frac{H}{U_0}$ pour le temps, les équations adimensionnelles s'écrivent:

$$\frac{\partial u}{\partial x} + \frac{\partial w}{\partial z} = 0$$

$$\frac{\partial u}{\partial t} + \frac{\partial(uu)}{\partial x} + \frac{\partial(wu)}{\partial z} = -\frac{\partial p}{\partial x} + \frac{\partial}{\partial x}\left(\frac{1}{\text{Re}}\left(1+\frac{v_t}{v}\right)\frac{\partial u}{\partial x}\right) + \frac{\partial}{\partial z}\left(\frac{1}{\text{Re}}\left(1+\frac{v_t}{v}\right)\frac{\partial u}{\partial z}\right) + \frac{\partial}{\partial x}\left(\frac{v_t}{v}\frac{\partial u}{\partial x}\right) + \frac{\partial}{\partial z}\left(\frac{v_t}{v}\frac{\partial w}{\partial x}\right)$$

$$\frac{\partial w}{\partial t} + \frac{\partial(uw)}{\partial x} + \frac{\partial(ww)}{\partial z} = -\frac{\partial p}{\partial z} + \frac{\partial}{\partial x}\left(\frac{1}{\text{Re}}\left(1+\frac{v_t}{v}\right)\frac{\partial w}{\partial x}\right) + \frac{\partial}{\partial z}\left(\frac{1}{\text{Re}}\left(1+\frac{v_t}{v}\right)\frac{\partial w}{\partial z}\right) + \frac{\partial}{\partial x}\left(\frac{v_t}{v}\frac{\partial u}{\partial z}\right) + \frac{\partial}{\partial z}\left(\frac{v_t}{v}\frac{\partial w}{\partial z}\right)$$

où u et w sont les deux composantes longitudinale et transversale de la vitesse et p est la pression motrice divisée par le carré de la vitesse de référence.

$Re = \frac{U_0 \cdot H}{v}$ et $v_t$ sont respectivement le nombre de Reynolds et la viscosité turbulente de la sous maille.

Pour calculer $v_t$, on utilise un modèle d'échelles mixtes proposé par Ta Phuoc [3] et Sagaut [4] basé sur la viscosité de la sous maille: $v_t = v_s v_d^{1-\alpha}$, où $v_s = c_s \Delta^2 |\bar{s}|$ et $v_d = c_d \Delta \sqrt{q_c}$, $v_s$ et $v_d$ sont respectivement les viscosités turbulentes du modèle de Smagorinsky et du modèle dynamique.
$c_s$ et $c_d$ sont les paramètres du modèle. $\Delta$ est la largeur du filtre et $|\bar{s}| = (\bar{s}_{i,j} \bar{s}_{i,j})^{1/2}$ est le module du tenseur des contraintes, $\bar{s}_{i,j} = \frac{1}{2} \left( \frac{\partial \bar{u}_i}{\partial x_j} + \frac{\partial \bar{u}_j}{\partial x_i} \right)$ et $q_c = \frac{1}{2} (\bar{u}_i)^{\prime} (\bar{u}_i)^{\prime}$, $q_c$ est l'énergie turbulente. On retrouve le modèle de Smagorinsky pour $\alpha = 0$ et le modèle dynamique pour $\alpha = 1$.

A la suite d'une étude comparative entre les différents modèles de sous maille, on a retenu le modèle dynamique mixte ($\alpha = 0.5$) pour ses résultats très proches de ceux de la simulation directe. Il est à remarquer que les solutions retenues pour les comparaisons sont obtenues lorsque le régime asymptotique est atteint.

## 4. METHODE DE RESOLUTION

Le schéma temporel adopté pour résoudre le système d'équations précédent est de type différences finies d'ordre 2. Ce schéma, appelé Euler d'ordre 2 Adams-Bashford, consiste à traiter explicitement les termes non linéaires (termes convectifs) et implicitement les termes diffusifs. Les équations de conservation de la quantité de mouvement sont discrétisées suivant la méthode de volumes finis à maillage décalé. Chaque inconnue est définie sur une grille qui lui est propre.

Les termes non linéaires ont été approchés par un schéma Quick (Quadratic Upstream Interpolation for Convection Kinematic). Le traitement des équations de mouvement associées à la condition d'incompressibilité est effectué par la méthode de projection. La méthode aux directions alternées (ADI) est utilisée pour résoudre la partie prédiction. La partie projection est résolue par une procédure multigrille linéaire.

## 5. RESULTATS ET DISCUSSION

La résolution des équations de Navier Stokes et de l'équation de l'énergie pour un écoulement turbulent bidimensionnel instationnaire et incompressible a été faite en utilisant un maillage de 66 × 66 points, des nombres de Reynolds croissants et un pas de temps $\Delta t = 5 \cdot 10^{-3}$ s.

La figure 2 (a, b, c, d et e) représente le champ de vorticité pour le cas d'une seule cavité, pour des nombres de Reynolds croissants ($Re=10^4$, $Re=5.10^4$, $Re=10^5$, $Re=5.10^5$ et $Re=10^6$).

<figure><img src="image_2.png" /></figure>
<figure><img src="image_3.png" /><figcaption>Fig. 2: Isovaleurs du champ de vorticité pour différents nombres de Reynolds</figcaption></figure>

<caption>Fig. 3: Isovaleurs du champ de vorticité pour différents nombres de Reynolds</caption>

On constate que le décollement de la couche cisaillée a provoqué l'apparition d'un tourbillon au pied de la première marche. Ce dernier devient de plus en plus intense tout en augmentant le nombre de Reynolds. Il est intéressant de remarquer la tendance à l'organisation en structures allongées tout au long de la cavité avec l'existence d'une petite structure face à l'écoulement formé sous l'effet du recollement du fluide avec la marche. La petite structure disparaît au fur et à mesure qu'on augmente le nombre de Reynolds. Ceci peut être du à l'augmentation de la vitesse de l'écoulement qui fait élargie la couche cisaillée.

La figure 3 (a, b, c, d et e) représente le champ de vorticité dans le deuxième cas. On constate que l'écoulement est bloqué par l'obstacle et que le fluide ne pénètre pas dans la première cavité. La structure de sillage dans la seconde cavité révèle l'existence d'une allée tourbillonnaire composée des structures contra-rotatives.

<figure><img src="image_4.png" /></figure>
<figure><img src="image_5.png" /><figcaption>Fig. 4: Lignes de courant pour différents nombres de Reynolds</figcaption></figure>

Sur la figure 4 (a, b, c, d et e), sont représentées les lignes de courant de l'écoulement, pour les différents nombres de Reynolds. On remarque que la structure tourbillonnaire, qui s'établit dans le canyon, est de plus en plus complexe au fur et à mesure qu'on augmente le nombre de Reynolds. En effet, on observe deux structures tourbillonnaires, l'une derrière la paroi sous l'écoulement dont l'extension horizontale est la plus grande, l'autre proche de la paroi face à l'écoulement.

La figure 4 montre la présence d'une couche cisaillée, qui se décolle à partir du bord amont. Cette couche s'élargit au fur et à mesure qu'elle s'approche de la paroi avale. On constate l'existence d'un point de rattachement de cette couche avec le fond de la cavité. Cette zone appelée zone de rattachement correspond à une adhérence de l'écoulement à la paroi. Ensuite, on observe l'impact de cette couche avec le coin aval qui va donner naissance, d'une part, à de petites structures localisées dans les coins inférieurs de la cavité (zone de recirculation secondaire), d'autre part, aux petits chapelets de tourbillons caractérisants le phénomène d'échappement tourbillonnaire.

Ces structures tourbillonnaires évoluent suivant un processus bien établi: l'une est éjectée dans la couche limite rattachée dans la partie avale de la cavité: c'est le phénomène d'échappement tourbillonnaire. Lorsque un premier tourbillon est formé à l'angle amont, il croît tout en étant convecté par le champ moyen de vitesse jusqu'à atteindre l'angle aval. Lorsqu'il impacte le coin aval, il se scinde en deux parties: l'une est éjectée dans la couche limite rattachée dans la partie avale de la cavité: c'est le phénomène d'échappement tourbillonnaire. L'autre partie est aspirée vers la cavité et vient alimenter la zone de recirculation déjà existante.

Sur la figure 5, on représente la variation du nombre de Strouhal $St_{ec} = \frac{f_{ec} h}{u_c}$, (où $u_c$ désigne la vitesse moyenne de convection calculée, $h$ est la hauteur de la cavité et $f_{ec}$ est la fréquence d'échappement tourbillonnaire), pour une position de mesure au milieu de la cavité large et pour les deux configurations. Les résultats montrent que, pour les deux cas, le nombre de Strouhal croît au fur et à mesure que le nombre de Reynolds augmente et les deux courbes ont la même allure. La présence de l'obstacle dans la cavité a un grand effet sur la variation du nombre de Strouhal et par suite sur la fréquence tourbillonnaire à l'intérieur de la cavité et ceci est confirmé par l'enregistrement d'un écart important entre les deux courbes.

<figure><img src="image_6.png" /><figcaption>Fig. 5: Variation du nombre de Strouhal en fonction du nombre de Reynolds en un point, au milieu de la large cavité</figcaption></figure>
Sur la figure 6, on présente les profils de vitesse longitudinale pour deux valeurs du nombre de Reynolds: $Re=10^4$ et $Re=10^6$, représentant sa plus faible et sa plus grande valeur. Sur la figure 6 ($a_1-a_2$), on remarque que les profils enregistrés sont presque les mêmes pour les deux configurations. Ces profils deviennent de plus en plus différents pour les deux configurations à partir de la figure 6 ($b_1-b_2$). Cette différence de profils peut être expliquée par l'existence de l'obstacle dans la cavité. Cette différence est bien observée sur la figure 6 ($c_1-c_2$) et pour la section de mesure prise au milieu de la deuxième cavité (Fig. 6 ($d_1-d_2$)). Pour la section à la sortie du domaine (Fig. 6 ($e_1-e_2$)), ces profils deviennent identiques. Ces résultats montrent que la perturbation due à la présence de l'obstacle s'atténue à la sortie de la cavité.

<figure><img src="image_7.png" /><figcaption>Fig. 6: Profils de vitesse longitudinale pour les deux configurations pour $Re=10^4$ et $Re=10^6$</figcaption></figure>
## 6. CONCLUSION

Cette étude nous a permis de décrire la structure de l'écoulement dans une cavité large, d'étudier l'effet d'un obstacle placé dans cette cavité et d'analyser l'allure des structures tourbillonnaires ainsi que le changement de ce phénomène en fonction du nombre de Reynolds.

## NOMENCLATURE

L : Longueur du domaine, m

H : Hauteur du domaine, m

p : Pression modifiée

Pr : Nombre de Prandtl

R : Rapport de forme du domaine

r : Rapport de forme

Re : Nombre de Reynolds

Ra : Nombre de Rayleigh

Ri : Nombre de Richardson

U₀ : Vitesse maximale, m².s⁻¹

t : Temps, s

Pr_t : Nombre de Prandtl turbulent

x : Abscisse, m

z : Ordonnée, m

Lettres Grecques

v : Viscosité cinématique, m².s⁻¹

v_t : Viscosité turbulente, m².s⁻¹

ρ₀ : Masse volumique de l'air

## REFERENCES

[1] P. Baranger et al., 'Improved Turbulence Models Based on Large Eddy Simulation of Homogeneous Incompressible Turbulent Flows', Report, TF19, Stanford University, 1983.

[2] T.R. Oke, 'Street Design and Urban Canopy Layer Climate', Energy and Building, Vol. 11, N°1-3, pp. 103–113, 1988.

[3] L. Ta. Phuoc, 'Modèles de Sous Maille Appliqués aux Ecoulements Stationnaires et Décollés', Journée Thématique DRET: Aérodynamique Stationnaire Turbulente, Aspects Numériques et Expérimentaux, 1994.

[4] P. Sagaut, 'Introduction à la Simulation à Grandes Echelles pour les Ecoulements de Fluide Incompressible', Ouvrage: Mathématiques et Applications, 1998.