Coefficients de pression sur un cube - Comparaison entre approches CFD avec Code_Saturne et résultats expérimentaux

Résumé

Cet article traite de la modélisation CFD des contraintes en pression générées par le vent sur des infrastructures. Le travail présenté a été réalisé avec le logiciel open source code_saturne[1], avec lequel plusieurs modèles de turbulence ont été utilisés. L’objectif est de mettre en valeur les forces et faiblesses de chaque modèle testé grâce à une comparaison avec des données expérimentales, ceci afin de pouvoir sélectionner un bon compromis entre précision et vitesse de calcul lors d’études industrielles.

Mots-clés : CFD, wind engineering, code_saturne

Introduction

Actuellement, dans le monde du génie civil, l’impact du vent sur des infrastructures est principalement évalué via des méthodologies prescriptives, comme Eurocode[2]. Cependant, il est parfois nécessaire de disposer de résultats plus précis, et la CFD est l’un des moyens pour y parvenir. Toutefois, pour modéliser une situation donnée, de nombreuses approches sont possibles. En particulier, des compromis s’avèrent nécessaires entre la précision requise et le coût de la simulation, que ce soit en termes de temps ou de puissance de calcul. Ces compromis concernent particulièrement deux points :

  • La modélisation de la turbulence. Historiquement, les modèles de type RANS avec une approche basée sur la notion de viscosité turbulente, comme les modèles k-ε ou k-ω, sont souvent utilisés grâce à leur faible coût numérique et à leur relative facilité d’utilisation. Au fil du temps, de nombreuses variantes de ces modèles ont été développées, chacune étant adaptée à des situations données.
    Dans la famille RANS, des modèles du second ordre ont également été développés, nommés « Reynolds Stress Models » (RSM). Ces derniers sont plus complets, car ils résolvent les équations de transport du tenseur de Reynolds.
    Enfin, depuis quelques temps, la croissance des ressources informatiques permet progressivement l’utilisation de la LES. Toutefois, pour des cas industriels de grande taille, cela reste la plupart du temps très coûteux et globalement irréaliste.

  • La technique de maillage choisie, dépendant des schémas numériques disponibles et du modèle de turbulence choisi. Les maillages semi-automatiques avec génération de couches limites ont fait d’énormes progrès durant les deux dernières décennies. Toutefois, dans la plupart des logiciels de CFD, certaines modélisations (notamment la LES) nécessitent encore l’utilisation d’un maillage de « bonne qualité » avec une structure hexaédrique.

Les avantages et inconvénients des différentes approches sont décrits dans de nombreuses publications. Dans cette étude, nous proposons d’obtenir un aperçu de différentes possibilités proposées par le logiciel open source code_saturne[1], développé par EDF R&D. L’objectif de l’analyse est de comparer les résultats obtenus pour un certain nombre de modèles de turbulence avec deux jeux de données expérimentales concernant les coefficients de pression sur un cube soumis au vent. Trois types de modèles sont considérés :

  • Les modèles basés sur une approche « haut-Reynolds »
  • Les modèles basés sur une approche « bas-Reynolds »
  • Les modèles de type LES

La différence de résultats entre les modélisations CFD et les données expérimentales, ainsi qu’entre les différents modèles de CFD eux-mêmes, permettra d’identifier les meilleures approches en termes de précision, mais également en termes de coûts numériques. En particulier, nous nous attarderons sur la pertinence d’une approche de type « haut-Reynolds », qui est la moins contraignante concernant la taille et la qualité du maillage à utiliser.

Données expérimentales utilisées

Deux jeux de données expérimentales sont considérés dans cette étude :

  • Les données d’un cube sur site réel
  • Les données d’un cube de taille réduite utilisé lors de tests en soufflerie

Ces deux ensembles de données fournissent des résultats assez comparables, ce qui permet de donner du crédit aux résultats obtenus (cf. paragraphe de résultats). Les sections suivantes proposent un descriptif succinct de ces deux expérimentations.

Cube de Silsoe sur site réel

La configuration expérimentale sur site réelle est celle du « Silsoe Experimental Building » (SEB) présentée sur la Figure 1. Il s’agit d’un cube de dimension [6m×6m×6m], situé à l’institut de recherche de Silsoe (SRI), et spécifiquement construit pour réaliser des expérimentations en taille réelle. Ce cube a été construit sur un terrain plat avec de l’herbe rase, ce qui permet d’avoir une idée assez précise de la longueur de rugosité (environ 0.01m). Par ailleurs, le profil de couche limite peut être facilement modélisé dans l’optique d’une approche numérique.

Des capteurs en pression sont situés sur certaines faces du cube, selon des lignes à mi-hauteur ou mi-profondeur (carrés gris sur la Figure 1). Davantage de détails sur les capteurs utilisés et les mesures réalisées sont disponibles dans  Richards and Hoxey, 2008[5]. Au final, les données expérimentales disponibles contiennent les propriétés moyennes, fluctuantes et spectrales de la pression (cf. Richards et al., 2007[3], Richards and Hoxey, 2008[5]).

Figure 1 - Aperçu du cube de Silsoe. Image tirée de Richards et al. [3]

Cube réduit pour soufflerie

Le second jeu de données expérimentales utilisé dans cet article vient de l’étude d’un cube de dimension [0.2m×0.2m×0.2m] en soufflerie par Irtaza et al., 2013[4]. Ce cube, montré sur la Figure 2, possède donc une échelle 1:30 par rapport au cube de Silsoe.

experimental scaled cube
Figure 2 - Aperçu du cube à échelle réduite pour soufflerie. Image tirée de Irtaza et al.[4]

La vitesse moyenne à la hauteur de l’arête supérieure du cube a été fixée à 4.81m/s pour les expérimentations. Cette valeur permet d’avoir un nombre de Reynolds compris entre 0.72×105 et 1.09×105, une plage de valeur similaire à celle obtenue lors de l’expérience menée sur site réel. La localisation des capteurs de pression sur la face supérieure du cube est indiquée sur la figure 3. Tous les détails concernant les expérimentations réalisées et l’assimilation des données sont disponibles dans Irtaza et al., 2013[4]. La configuration sélectionnée pour nos modélisations est celle d’un vent de face avec une incidence nulle.

Figure 3 - Position des capteurs de pression sur la face supérieure du cube. Image tirée de Irtaza et al.[4]

Profiles expérimentaux obtenus

Le profil de vitesse et l’intensité turbulente du cube à taille réelle ont été mesurés par Richards et al., 2007[3] et sont représentés respectivement sur la Figure 4 et Figure 5 pour un vent de face avec un angle d’incidence nul. Le profil de vitesse est similaire à un profil logarithmique, tandis que l’intensité turbulente longitudinale affiche une valeur d’environ 18%.

Des données similaires sont disponibles pour le cube à échelle réduite dans Irtaza et al., 2013[4], données également indiquées sur la Figure 4 et la Figure 5. Le profil de vitesse moyenne, tout comme l’intensité turbulente, sont assez similaires à ceux obtenus lors de l’expérimentation à taille réelle.

experimental results - velocity ratio
Figure 4 - Profil de vitesse moyenne pour les expérimentations à taille réelle et à taille réduite. Les données de l'expérimentation à taille réelle sont tirées de Richards et al. [3], tandis que celles de l'expérimentation en soufflerie sont tirées de Irtaza et al.[4]
experimental results - turbulent intensity
Figure 5 - Profil d'intensité turbulente pour les expérimentations à taille réelle et à taille réduite. Les données de l'expérimentation à taille réelle sont tirées de Richards et al. [3], tandis que celles de l'expérimentation en soufflerie sont tirées de Irtaza et al.[4]

Modélisation CFD

Le logiciel code_saturne

Les simulations ont été réalisées avec code_saturne, un logiciel open source développé par EDF R&D. Comme écrit sur le site officiel de code_saturne, « It solves the Navier-Stokes equations with a finite volume approach for 2D, 2D-axisymmetric and 3D flows, steady or unsteady, laminar or turbulent, incompressible or weakly dilatable, isothermal or not, with scalar transport ». De plus amples informations sont disponibles dans la publication de Frédéric Archambeau et al., 2004[1], ou encore sur le site https://www.code-saturne.org/cms/web/. C’est la version 7.0.5 du logiciel qui a été utilisée durant cette étude.

code_saturne logo

Description du modèle

Géométrie et conditions aux limites

Le cube modélisé est le cube à échelle réduite, de dimension [H×H×H] avec H=0.2m. Le domaine de calcul utilisé est présenté sur les Figures 6, 7 et 8 :

  • Le vent en entrée de domaine est selon la direction X, sans angle d’attaque.
  • Le domaine possède les dimensions suivantes : 29H selon la direction principale de l’écoulement X, 13H selon la direction latérale Y et 5H selon la direction verticale Z.
  • Le cube est situé 6H en aval de l’entrée du domaine, 22H en amont de la sortie, et il est centré dans la direction Y avec 6H de chaque côté du cube.

Avec de telles dimensions, l’écoulement autour de l’obstacle n’a pas d’effet sur les flux entrant ou sortant au niveau des conditions aux limites, ce qui est un prérequis indispensable dans de nombreuses applications de CFD.

Figure 6 - Aperçu du domaine de calcul
config simulation - lateral view
Figure 7 - Vue latérale du domaine de calcul
config simulation - top view
Figure 8 - Vue de dessus du domaine de calcul

La condition à la limite en entrée est configurée selon une loi logarithmique pour la vitesse :

U(z) = \cfrac{U_{\tau}}{\kappa} ln \Big(\cfrac{z+z_0}{z_0}\Big)

où Uτ = 0.648m/s est la vitesse de friction, κ=0.41 est la constante de Von Karman et Z0 =0.01m est la longueur de rugosité. Le profil obtenu, montré sur la Figure 9, assure une valeur U(z=H)=4.81m/s, comme pour l’expérimentation en soufflerie précédemment décrite.

D’après la Figure 5, l’intensité turbulente relevée lors des deux expérimentations décroît pour z/H < 2, avec une valeur d’environ 18% en z/H=1. Un profil simplifié,  présenté sur la Figure 10, est utilisé dans nos simulations, avec une décroissance linéaire de  I(z/H=0)=20% à I(z/H=2)=16%, et une intensité turbulente constante à 16% pour z/H > 2.

Figure 9 - Profil de vitesse imposé en entrée
Figure 10 - Profil d'intensité turbulente imposé en entrée

Les paramètres de turbulence à spécifier pour la condition à la limite d’entrée dépendent du modèle de turbulence choisi (cf. section suivante), toutefois les relations utilisées sont toujours les mêmes :

  • L’énergie cinétique turbulente est définie comme suit :

k(z) = \cfrac{3}{2} \big(U_{avg} I(z)\big)^2

où Uavg est la vitesse moyenne en entrée, et I(z) est l’intensité turbulente.

  • Le taux de dissipation est défini selon l’équation suivante :

ε(z) = C_{\mu}^{3/4}\cfrac{k^{3/2}}{l_s}

où Cμ=0.09 et ls est l’échelle de longueur caractéristique de la turbulence.

  • Le taux de dissipation spécifique est relié aux précédentes définitions de la manière suivante :

ω(z) = \cfrac{ε}{C_{\mu}k}

  • Dans le cas d’un modèle de type RSM, le tenseur de Reynolds est défini comme suit :

R_{11}=\cfrac{2k}{3}; R_{22}=\cfrac{2k}{3}; R_{33}=\cfrac{2k}{3}; R_{12}=0; R_{13}=0; R_{23}=0

  • Concernant la LES, les paramètres par défaut de code_saturne sont conservés.

Paramètres physiques et numériques

Les hypothèses suivantes sont adoptées pour les simulations :

  • L’air est considéré incompressible avec ρ=1.2kg/m3.
  • L’étude concerne avant tout les effets dynamiques dus au vent . En conséquence, la gravité n’est pas prise en compte dans les simulations.
  • Aucune résolution thermique n’est considérée.
  • Pour toutes les simulations de type RANS, une approche quasi-stationnaire est adoptée, avec un pas de temps local et adaptatif, permettant de garantir un nombre CFL inférieur à 1.
  • Pour les simulations LES, une approche instationnaire est privilégiée, avec un pas de temps constant de 5×10-4s.
  • Un schéma centré du second ordre est utilisé pour la résolution spatiale de la vitesse, tandis qu’un schéma upwind du premier ordre est utilisé pour les variables liée à la turbulence dans les modélisations RANS.

Caractéristiques des maillages

Deux maillages sont utilisés, le choix dépendant du modèle de turbulence sélectionné (voir la section suivante). Dans tous les cas, ils sont non structurés mais à dominante hexaédrique, avec une taille de maille allant de 40mm à 2.5mm, le principal raffinement se situant autour du cube et pour un ratio de hauteur z/H < 2, comme montré sur la Figure 11. Les maillages ont été générés via des outils semi-automatiques avec un algorithme de type Octree. La différence entre eux concerne la résolution de la sous-couche visqueuse :

  • Le premier maillage « HR_Mesh » est configuré pour ne pas résoudre la sous-couche visqueuse.
  • Le second, « LR_Mesh« , inclut la définition d’une couche limite, avec une hauteur de 0.5mm pour la première cellule dans le sens normal à la paroi.

Pour valider la pertinence de l’utilisation de ces maillages, les tests habituels de convergence ont par ailleurs été réalisés.

Realized mesh of the domain
FIgure 11 - Aperçu du maillage autour du cube

Modèles de turbulence sélectionnés

De nombreux modèles de turbulence sont disponibles sur code_saturne (cf. le guide théorique[6] et le guide pour utilisateur[7]). Six d’entre eux sont retenus pour cette étude :

  • Modèles RANS du premier ordre :  k-ε linear production, k-ω SST et v2f BL v2k
  • Modèles RANS du second ordre (RSM) : SSG et EBRSM
  • Modèle LES : WALE

Selon leurs caractéristiques, ils ont été séparées en trois catégories, décrites ci-dessous.

Approche "Haut Reynolds"

Dans cette approche, la sous-couche visqueuse n’est pas résolue : une loi de paroi de type « scalable log law » d’ordre 2 est utilisée sur les faces du cube et au niveau du sol, avec le maillage « HR_Mesh ». De cette façon, les valeurs de y+ se situent essentiellement autour de 30. Cette stratégie est utilisée pour les deux modèles de turbulence suivants :

  • k-ε linear production
  • RSM-SSG

Approche "Bas Reynolds"

Dans l’approche « Bas Reynolds », la sous-couche visqueuse est résolue. Le maillage « LR_Mesh » est utilisé, sans loi de paroi. Ceci permet d’obtenir des valeurs de y+ autour de 1. Cette stratégie est utilisée pour les quatre modèles de turbulence suivants :

  • k-ε linear production
  • k-ω SST
  • v2f BL v2k
  • RSM-EBRSM

Il est à noter que cette approche n’est théoriquement pas recommandée avec le modèle k-ε linear production. Ceci permettra néanmoins de comparer les résultats des deux approches avec ce dernier.

Approche LES

Pour la simulation LES, le maillage « LR_Mesh » est utilisé avec une loi logarithmique d’ordre 1. Cela permet de maintenir les valeurs de y+ proches de 1, ainsi que Δx+ et Δz+ < 20.

  • Le modèle LES retenu est le suivant : LES-WALE

Résultats

Convergence et variables post-traitées

Pour s’assurer de la bonne convergence des simulations, des sondes ont été placées dans certaines zones d’intérêt :

  • Proche de la section d’entrée, afin de contrôler les profils injectés.
  • Autour du cube, afin de vérifier la stabilité des résultats.
  • Proche de la section de sortie, afin de contrôler d’éventuels retours d’écoulement ou instabilités numériques venant de cette condition à la limite.

Toutes les simulations réalisées ont correctement convergé vers des résultats stables. Une fois la convergence atteinte, des valeurs moyennes ont été construites durant la fin des simulations, comme illustré sur les Figures 12 et 13 pour la simulation LES. L’analyse des résultats présentée dans les sections suivantes ne concerne que ces valeurs moyennées.

sensors location
Figure 12 - Sondes sur la face au vent pour le contrôle de la convergence
Pressure over time for all sensors
Figure 13 - Evolution de la pression sur les sondes de la face au vent pour la simulation LES, et illustration de la plage temporelle retenue pour réaliser les valeurs moyennes analysées par la suite

Enfin, les pressions calculées sont exprimées sous une forme adimensionnelle, suivant la définition du coefficient de pression Cp :

C_p = \cfrac{p_{mean}-p_0}{\cfrac{1}{2} \rho U_h^2}

où pmean est la pression moyennée temporellement, p0 est la moyenne spatiale de la pression au niveau de la section d’entrée, et Uh = 4.81m/s est la vitesse en entrée à la hauteur H=0.2m.

Discussion

Comme il existe une bonne correspondance entre les résultats à échelle réelle de Richards et al.[3] et ceux à échelle réduite d’Irtaza et al.[4], l’analyse est dans un premier temps basée sur la comparaison avec ces résultats expérimentaux. La Figure 14 and la Figure 15 présentent respectivement les coefficients de pression sur deux chemins selon la mi-largeur et la mi-profondeur du cube pour toutes les simulations réalisées, ainsi que la comparaison avec les résultats des tests expérimentaux mentionnés.

Figure 14 - Coefficient de pression au niveau de la mi-largeur du cube

En se basant sur la Figure 14, il est possible de constater que dans la direction principale de l’écoulement :

  • Le modèle k-ε linear production donne les tendances les moins bonnes, que ce soit avec une approche « haut Reynolds » ou « bas Reynolds ». Il surestime les coefficients de pression sur l’essentiel de la face au vent, avec des valeurs de Cp deux fois plus importantes que dans les tests expérimentaux sur la partie supérieure de cette face (distance adimensionnée entre 0.8 et 1). A l’inverse, les coefficients de pression sur la face supérieure ont des valeurs négatives moins prononcées que dans  Richards et al.[3] et Irtaza et al.[4], en dehors de la présence d’un important minimum sur la partie avant de la face (distance adimensionnée entre 1 et 1.2).

  • Les autres modèles fournissent des résultats relativement proche des données expérimentales sur la face au vent (distance adimensionnée entre 0 et 1), mais seule la simulation LES reste proche de ces données sur la face supérieure (distance adimensionnée entre 1 et 2). Sur cette dernière, on peut toutefois noter des différences entre les simulations k-ε/v2f d’une part, avec des valeurs de Cp environ 50% plus importantes que dans les données expérimentales, et les simulations k-ω SST/SSG/EBRSM d’autre part, avec des valeurs et tendances plus proches desdites données.

  • Enfin, sur la face sous le vent (distance adimensionnée entre 2 et 3), toutes les simulations fournissent des résultats assez proches, avec des valeurs de coefficient de pression demeurant un peu plus élevées que celles obtenues dans les tests expérimentaux.
Figure 15 - Coefficient de pression au niveau de la mi-profondeur du cube

En se basant sur la Figure 15, il est possible de constater que sur un profil tracé à mi-profondeur du cube :

  • Seule la simulation LES fournit des résultats relativement proches des données expérimentales, bien qu’il subsiste une différence de l’ordre de 10 à 20%.

  • L’ensemble des autres modèles testés prévoient des coefficients de pression plus faibles relativement parlant, avec par ailleurs des différences substantielles de comportement en comparaison des résultats expérimentaux. C’est particulièrement visible au niveau des deux arêtes de la face supérieure du cube (distance adimensionnée = 1 et distance adimensionnée = 2).

  • Comme observé précédemment avec la Figure 14, les modèles k-ε linear production et v2f fournissent les résultats les plus éloignés des données expérimentales.

Les Figures 16 à 18 présentent maintenant les cartographies et contours des coefficients de pression obtenus pour toutes les simulations réalisées sur trois des faces du cube (respectivement la face au vent, la face supérieure et l’une des deux faces latérales). En assumant le constat que la simulation LES est celle qui se rapproche le plus des données expérimentales, c’est cette dernière qui sert de base pour la comparaison des résultats de l’ensemble des simulations.

Selon ces figures, il est ainsi possible de constater les choses suivantes :

  • Sur la face au vent, le modèle k-ε linear (avec les deux approches) semble être hors jeu, avec des valeurs de coefficient de pression beaucoup plus élevées qu’avec tous les autres modèles. Les résultats des autres modèles sont par ailleurs relativement proches, notamment concernant la localisation du point de stagnation de l’écoulement, bien que les valeurs de Cp soient globalement plus basses dans la simulation LES.

  • Sur la face supérieure, en accord avec les observations tirées des Figures 14 et 15, aucun des modèles RANS ne propose de résultats proches de ceux de la LES. Si l’on se base sur les valeurs très négatives des coefficients de pression observées proche de l’arête face au vent pour l’ensemble de ces modèles, on peut émettre l’hypothèse que la différence de comportement est potentiellement causée par une surestimation de la boucle de recirculation juste derrière cette arête, ce qui engendre un manque de dynamisme sur la partie restante de la face supérieure. Ceci est particulièrement vrai pour les modèles k-ε LP et v2f models, les plus proches de la LES étant les modèles  RSM-EBRSM et k-ω SST. Cette hypothèse semble confirmée par la Figure 19, qui présente les lignes de courant selon une coupe à mi-largeur du cube pour trois cas (les simulations LES, k-ε LP et k-ω SST).

  • Sur la face latérale, les observations sont assez similaires à celles de la face supérieure : les modèles RSM-EBRSM et k-ω SST sont les modèles RANS qui se rapprochent le plus des résultats de la simulation LES, bien que les différences demeurent notables.
Figure 16 - Cartes et contours des coefficients de pression obtenus sur la face au vent pour l'ensemble des simulations
Figure 17 - Cartes et contours des coefficients de pression obtenus sur la face supérieure pour l'ensemble des simulations
Figure 18 - Cartes et contours des coefficients de pression obtenus sur l'une des faces latérale pour l'ensemble des simulations
Figure 19 - Lignes de champ sur une coupe à mi-largeur du cube pour les simulations LES, k-ε linear production et k-ω SST

Conclusion & perspectives

Dans cette étude, nous avons pu tester plusieurs modèles de turbulence implémentés dans Code_saturne sur une configuration avec deux jeux de données expérimentales disponibles. Cela a permis de mettre en avant les avantages d’une approche LES pour ce genre d’étude, étant donné qu’aucun des modèles RANS testés n’a pu prévoir avec précision les coefficients de pression obtenus expérimentalement sur cette configuration basique. En particulier, les résultats observés ont montré que les modèles k-ε linear production et v2f ne devraient pas être utilisés pour ce type d’application. Au regard des résultats des autres modèles, il apparaît que le modèle k-ω SST avec une approche « bas Reynolds » pourrait être un bon compromis, étant donné que les modèles RSM sont plus sujets à de potentielles instabilités numériques dans des configurations complexes.

Pour compléter cette approche, plusieurs points doivent toutefois être approfondis, comme l’impact d’un maillage totalement non structuré (tétra-dominant avec présence d’une couche limite prismatique, par exemple), et la stabilité des simulations en cas de configuration plus complexe. Enfin, si cette étude s’est avant tout concentrée sur les coefficients de pression moyens pour discriminer les différents modèles testés, il faut garder à l’esprit que l’un des apports souhaitables de la CFD dans le génie civil serait surtout d’avoir des informations sur les fluctuations instationnaires (coefficients de pression de pointe).

Nomenclature

CFD : Computational Fluid Dynamics
LES : Large Eddy Simulation
RANS : Reynolds Averaged Navier Stokes
RSM : Reynolds Stress Models
SEB : Silsoe Experimental Building
SRI : Silsoe Research Institute

Références

[1] Frédéric Archambeau, Namane Méchitoua and Marc Sakiz,  ‘Code_Saturne: a Finite Volume Code for the Computation of Turbulent Incompressible Flows’, International Journal on Finite Volumes, volume 1, 2004.

[2] Eurocode 1: Actions on structures — Part 1-4: General actions — Wind actions, NF EN 1991-1-4, 2005.

[3] Richards P.J., Hoxey R.P., Connell B.D., Lander D.P. 2007. ‘Wind-tunnel modelling of Silsoe Cube’, Journal of Wind Engineer-
ing and Industrial Aerodynamics , pp. 1384-1399.

[4] H. Irtaza, R.G. Beale, M.H.R. Godley, A. Jameel, ‘Comparison of wind pressure measurements on Silsoe experimental building fro full-scale observation, wind-tunnel experiments and various CFD techniques’, International Journal of Engineering, Science and Technology, Vol. 5, No. 1, 2013, pp. 28-41.

[5] Richards P.J., Hoxey R.P. 2008. ‘Wind loads on the roof of a 6 m cube’, Journal of Wind Engineering and Industrial Aerodynam-
ics , pp. 984-993, 96.

[6] code_saturne 7.0 Theory Guide, 2021.

[7] code_saturne 7.0 practical user’s guide, 2021.

Laisser un commentaire

Votre adresse e-mail ne sera pas publiée. Les champs obligatoires sont indiqués avec *