3. Simulations et résultats

Mise en place des simulations

Dans ce cas de validation, nous allons faire l'hypothèse que la température n'a aucune influence sur l'hydrodynamique, c'est à dire que la température est considérée comme scalaire passif. L'écoulement est incompressible et les vitesses considérées en entrée sont telles que la transition vers un état turbulent ne sera jamais observée le long de la plaque.

Dans ces conditions, nous avons choisi d'utiliser le solveur incompressible icoFoam et d'y intégrer l'équation de la chaleur (avec terme de convection), le nouveau solveur porte ainsi le nom de my_icoFoam.

$$ \frac {\partial T} {\partial t} + Div(\phi T) = \frac {\lambda} {\rho c_p} \Delta T $$

Cette modification du solveur a déjà été étudiée : un tutoriel en anglais est disponible à cette adresse. L'idée pour nous était de voir comment il est possible d'intégrer simplement de nouvelles variables et équations à un modèle déjà existant. Ceci nous permet d'évaluer le potentiel d'OpenFOAM en terme de flexibilité, un des atouts majeurs de tout logiciel "open source".

De plus, avant toute simulation, il faut veiller à avoir un nombre de Courant inférieur à 1 :

$$Co=\frac {U  \Delta t} {\Delta x} < 1$$

$U$ : vitesse dans la direction $x$

$\Delta t$ : intervalle de temps à choisir pour avoir un nombre de Courant au plus égal à 1

$\Delta x$ : distance entre deux mailles

Dans notre cas, le pas de temps à choisir est au plus :

$$\Delta t=\frac {Co  \Delta x} {U} =\frac{1 \times  1.10^{-3}} {50.10^{-2}}=2  ms$$

 

Visualisation des résultats

La procédure pour effectuer la simulation et visualiser les résultats se fait par le biais du terminal en trois parties principales :

0. Se placer dans le répertoire de travail contenant le cas de calcul

         

Contenu du cas de calcul "Validation" avant simulation

 

1. Mailler le domaine en appelant l'outil blockMesh

 

2. Résoudre les équations mises en jeu en appelant le solveur approprié, ici my_icoFoam

3. Afficher les résultats dans le logiciel Paraview en appelant l'outil de post-traitement paraFoam

Champ de température sur la plaque avec l'iso-contour $T(\delta_T)=293,55  K$ avec une vitesse d'entrée $U_0=50  cm/s$ en régime établi (t > 0,4 s)

 

  

Champ de vitesse sur la plaque avec l'iso-contour $U(\delta)=0,495  m/s$ avec une vitesse d'entrée $U_0=50  cm/s$ en régime établi (t > 0,4 s)

 

 

Visualisation des couches limites

Afin de comparer les épaisseurs des couches limites thermique et dynamique, on peut les visualiser sous Paraview en utilisant la fonction isosurface en spécifiant les valeurs correspondantes en vitesse et température.

En exportant les résultats sur Matlab, il est alors possible de superposer les profils numérique et théorique :

 

Tracé des couches limites thermique et dynamique sur la plaque plane d'après les iso-contours OpenFOAM

 

On remarque d'une part que la couche limite thermique se développe plus vite que la couche limite dynamique. Ceci s'explique par le fait que le nombre de Prandtl est inférieur à 1.

D'autre part, on constate un certain écart avec la théorie : la couche limite thermique est sur-estimée tandis que la couche limite dynamique est quant à elle sous-estimée. Deux raisons peuvent expliquer cet écart :

  • Vitesse et température loin de la plaque sont mal évaluées (par exemple, la vitesse loin de la plaque n'est pas exactement égale à la vitesse imposée en entrée car la plaque modifie l'écoulement)
  • Les paramètres $\nu$ et $\alpha$ ne sont pas estimés à la bonne température

 

 

Evolution du Nusselt local

Afin de calculer le Nusselt local le long de la plaque, il nous faut connaître l'évolution du gradient de température sur la plaque, perpendiculairement à celle-ci.

Ceci a été réalisé en ajoutant un nouveau champ dans le fichier icoFoam_gradT.C et calculé seulement le long de la plaque. Le code à rajouter doit s'effectuer sur le fichier précédent placé dans le dossier nommé icoFoam_gradT localisé dans les dossiers applications puis solvers dont le contenu est le suivant :

Les lignes de code à ajouter sont les suivantes :

Capture du fichier texte "icoFoam_gradT.C"

 

La nouvelle variable ainsi créée apparaît désormais dans Paraview, nous permettant de connaître l'évolution de la dérivée normale à la paroi de la température. En adimensionnalisant cette grandeur (on la multiplie par $ \frac{x}{\Delta T} $), on obtient le profil du Nusselt local.

La figure suivante montre les résultats obtenus pour différentes vitesses en entrée :

vitesses nombre de Reynolds
$U_0=5  mm/s$ $33,44$
$U_0=5  cm/s$ $334,4$
$U_0=50  cm/s$ $3344$
$U_0=5  m/s$ $33440$

 

Le tracé du Nusselt montre d'excellents résultats vis-à-vis de la théorie, même dans le cas d'une vitesse d'entrée de 5 m/s qui donne un nombre de Reynolds s'approchant de la valeur critique de transition laminaire/turbulente. On peut cependant observer quelques différences au niveau du bord d'attaque de la paroi : cela s'explique par la résolution du maillage qui n'est pas assez élevée dans cette zone où les gradients sont les plus forts.