Je souhaite tout d'abord remercier ma directrice de thèse, Bérengère Dubrulle. Au cours de ce travail, elle a su me guider et me transmettre son énergie inépuisable ainsi que ses innombrables idées qui m'ont permis d'avancer à grands pas. J'ai véritablement apprecié l'attention particulière qu'elle m'a portée ainsi que sa pédagogie dans ses explications sur les aspects mathématiques et physiques de la turbulence.

Je remercie également mon directeur de thèse, Michel Stanislas, pour m'avoir encadré dans les meilleures conditions possibles. Qu'il trouve dans ces quelques lignes ma reconnaissance et ma gratitude. J'ai vraiment apprécié son aide précieuse dans l'organisation du travail de l'équipe, ses nombreuses relectures de ce manuscrit ainsi que son reconfort dans les moments difficiles.

Ensuite, je tiens à remercier Jean-Philippe Laval, avec qui j'ai eu de nombreuses discussions à propos de l'approche LES, de la turbulence, de différentes problèmes informatiques qui m'ont pris énormément de temps... Sa présence permanente au cours de ces années m'a permis d'avancer dans le développement du modèle LES dans de bonnes conditions.

J'exprime ma reconnaissance envers les rapporteurs de thèse, Pierre Sagaut et Sébastien Galtier. Leur travail de lecture et leur participation dans le jury m'ont aidé à approfondir mes connaissances et d'avoir une vision plus large du problème étudié. Je remercie également le président du jury, Daniel Buisine, pour sa participation à la soutenance et pour la discussion intéressante que nous avons eue.

J'apprecie énormément le travail et l'apport des connaissances de Matthieu Marquillie qui m'a formé sur différents aspects du calcul scientifique. Il m'a souvent épaté en répondant précisemment à toutes les questions auxquelles je me posais.

Je remercie tous les membres de l'équipe de Michel Stanislas avec qui j'ai passé des bons moments au sein du LML et à Lille, une ville qui restera dans mon coeur.

Y up c$ &P 6" F § ¦¥ '£ 2 ¦¥ ¢ ¦¥ '0 ¥ ( )$ & 0 ¦2 ¦R ¤¥ 93 w ¡ ¤2 ¦ r( %0 6 Qp '4 &© ¦ XB ¦0 ¦¡ &s 82 ¦R e$ & g 72 ¦ R 0 @$ 1p '£ ¦x b5 a 93 ¤0 @$ &£ ¦$ #5 2 £ ¦$ 0 @$ #p c£ ¦x b5 uP i¡ ¤£ 3 &2 ¦£ ¦¥ '£ 63 #$ #5 ) w$ 1£ ¦2 6" F¥ v c( %0 © f

Summary

The LES-Langevin model was studied on turbulent channel flows. The approach is based on the dynamics of the subgrid velocity scales previously studied in the frame of Rapid Distortion Theory. The subgrid stress tensor is modeled through a combination of a turbulent force and eddy-viscosity. The turbulent force equation is derived from the subgrid scales velocity dynamics with the hypotheses of the Rapid Distortion Theory. The complex nonlinear terms containing the subgrid scale pressure were modeled by a stochastic forcing. The well-known eddy-viscosity closure was chosen to improve the resolved scales dissipation subgrid Reynolds stress.

The advantage of the model is that the dynamics of turbulent force is prescribed by a dynamical equation derived from the Navier-Stokes equations. This allows a possibility to include all the important physical effects of the subgrid scales like the backward subgrid scale energy transfer. The direct modeling of the gradient of the subgrid scale tensor allows a reduction of the computational time compared to direct numerical simulation, which is the global objective of LES. The work demonstrates a need for further studies of the equations of the turbulent force in order to better estimate the subgrid scale action.

iii CHAPTER 0. SUMMARY

Résumé

Le modèle LES-Langevin a été étudié pour les écoulements turbulents de canal. Cette approche se base sur la dynamique des échelles sous-mailles qui d'abord a été étudié dans le cadre de la théorie de la distorsion rapide. Le tenseur de contrainte des échelles sous-mailles (τ ij ) a été modélisé par la combinaison d'une force turbulente et d'une viscosité turbulente. L'équation de la force turbulente est dérivée de la dynamique des échelles sous-mailles avec les hypothèses de la théorie de la distorsion rapide. Les termes non-linéaires contenant les échelles sous-mailles de la pression ont été modélisés par une force stochastique. La viscosité turbulente permet d'améliorer la dissipation des échelles résolues qui est produite par le tenseur de Reynolds sous-maille.

L'avantage du modèle est la modélisation de la force turbulente par une équation dynamique dérivée des équations de Navier-Stokes. Cela donne la possibilité d'inclure tous les effets importants des échelles sous-mailles représentés par les équations de Navier-Stokes. Par exemple, le transfert inverse d'énergie des échelles sous-mailles vers les échelles résolues. La modélisation directe du gradient du tenseur de contrainte des échelles sousmailles (τ ij ) réduit le temps de calcul par rapport à une simulation directe, ce qui est l'objective principal d'une simulation des grandes échelles (LES). Le travail a montré un besoin pour une étude plus approfondie de l'équation de la force turbulente afin de mieux reproduire l'action des échelles sous-mailles. v

Résumé étendu

Le travail de thèse porte sur le développement d'un nouveau modèle LES basé sur la théorie de la distorsion rapide. Le mémoire est composé de cinq chapitres incluant l'introduction et la conclusion.

Le premier chapitre concerne la définition de l'approche LES. Une revue de plusieurs fermetures sous-mailles est proposée. Les contraintes physiques et mathématiques qui s'appliquent sur la modélisation sous-maille sont également abordées et discutées. La non-commutation entre le filtrage et la dérivation spatiale et les conditions de symmetrie de l'écoulement sont en particulièrement évoquées. L'importance du transfert d'énergie due aux échelles sous-mailles est discutée. Ensuite, les principaux modèles LES sont représentés dans l'ordre de leur complexité mathématique. Cette partie débute par la présentation du modèle classique de Smagorinsky, son amélioration par l'introduction d'une fonction de Van Driest, les modèles basés sur cette fermeture (Modèle dynamique, "Mixed model", "Planar-averaged scale dependent dynamic model", "Lagrangian scale dependent and independent models"), puis le modèle de similarité et les modèles utilisants la séparation implicite et explicite des échelles de la vitesse ("Approximate deconvolution model", "Implicit LES", "Subgrid scale estimation model", "Variational multiscale approach" et "Resolvable subfilter scale model"). Les performances de ces différents types de modèles sont discutées sur la base des résultats decrits dans la littérature.

Le deuxième chapitre concerne le développement du modèle LES-Langevin. Comme toutes les approches LES, le modèle est basé sur la séparation explicite des échelles résolues (ū) et non-résolues (u ) de la vitesse: u = ū + u . L'approche consiste à dériver une équation dynamique pour le gradient du tenseur C ij = ūi u j + ūj u i qui fait partie du tenseur sous-maille τ ij = C ij + R ij , où R ij = u j u i . Le gradient du tenseur C ij est applé la force turbulente l: l i = ∂ j C ij et la dynamique de ce vecteur est la pierre angulaire du modèle LES-Langevin. La solution exacte des équations obtenues pour la force turbulente semble être très complexe, donc les hypothèses de la théorie de la distorsion rapide:

|∂ t ū| |∂ t u | |∂ t ∇p |
sont appliquées pour simplifier les équations obtenues (p étant la pression sous-mailles). La dérivation de l'équation pour la force turbulente est basée sur le modèle RDT-Langevin de la turbulence, développé par [START_REF] Laval | Langevin models of turbulence: Renormalization group, distant interaction algorithms or rapid distorsion theory?[END_REF]. Cette approche consiste à représenter la dynamique des petites échelles u par une équation approximée. Le développement de l'équation pour la force turbulente l utilise des approximations phénoménologiques pour simplifier les équations à résoudre numériquement:

• Les termes du modèle RDT-Langevin contenants la viscosité turbulente ν t sont vii CHAPTER 0. R ÉSUM É ÉTENDU regroupés dans un terme diffusif avec une nouvelle viscosité turbulente ν t .

• L'hypothèse sur le rapport des échelles de temps caractéristiques du gradient de la pression non-résolue et la vitesse non-résolue (|∂ t u | |∂ t ∇ • p |) permet de modéliser les termes de la pression non-résolue croisés avec la vitesse résolue par une force stochastique ξ. Les propriétées générales de cette force sont: < ξ i >= 0, < ξ i (t, x)ξ j (0, 0) >= T ij (t, x), où les propriétées spatio-temporelles de la fonction T ij dépendent de la physique de l'écoulement. Les études a priori de la force stochastique ξ sont fait sur le cas d'un écoulement de canal à Re τ = 600. L'influence du maillage et du filtre est mise en évidence. La forme spatiale du vecteur ξ et ses caractéristiques temporelles sont prises en compte par la modélisation.

• Le terme de friction -l/τ f est introduit dans l'équation du vecteur l pour satisfaire l'asymptotique du spectre moyen de la vitesse aux temps longues. Pour la turbulence isotrope, l'intégrale de collision a été estimé comme τ f = 27/22 < | S-1 | > [START_REF] Laval | A LES-Langevin model for turbulence[END_REF].

L'hypothèse de la rapidité d'évolution des échelles non-résolues de la vitesse u par rapport aux échelles résolues ū (|∂ t ū| |∂ t u |) permet de modéliser le terme ∂ j R ij par un terme diffusif avec une viscosité turbulente µ t .

Une analyse a priori des propriétés spatio-temporelles de la force turbulente pour l'écoulement de canal à Re τ = 600 a été réalisée. Il est montré que la valeur de la force turbulente est sensible pour toutes les échelles de l'écoulement. Cela montre que le transfert d'énergie reproduit par la force turbulente dépend de toutes les échelles. Les maxima des composantes du vecteur l sont situés dans la partie haute de la zone tampon (40 < y + < 60) et dans la partie basse de la zone log (60 < y + < 100). Cela montre que, même si le maximum de l'énergie cinétique est atteint au milieu de la zone tampon (13 < y + < 17), les variations du vecteur l continuent de croître avec la distance à la paroi. La forme spectrale du vecteur l est completement défini par le filtre. La valeur maximale du spectre se situe près de l'échelle de coupure, ce qui montre bien l'importance des échelles au voisinage de coupure dans le transfert d'énergie [START_REF] Kraichnan | Eddy viscosity in two and three dimensions[END_REF]. La propriété inattendue de ces spectres est que les échelles situées après la séparation sont dominantes énergétiquement. Ce comportement semble être indépendant de la distance à la paroi.

Le troisième chapitre présente une description du code numérique ainsi que sa validation pour la turbulence de canal. Le code permet de calculer l'écoulement dans un canal avec un profil donné sur la paroi basse. Les discrétisations spatiales et temporelles sont discutées, ainsi que la méthode de parallélisation utilisant les librairies MPI ("Message Passing Interface"). Les filtres explicites utilisés au cours des tests numériques, leurs propriétés et leur compatibilité avec le schéma numérique sont expliqués. La validation du code est présentée pour des résultats de simulation de canal à Re τ = 400 avec un modèle de Smagorinsky.

Le quatrième chapitre regroupe les résultats obtenus par la modélisation LES-Langevin de l'écoulement de canal à Re τ = 600. Dans un premier temps, la DNS de référence est comparée avec une autre DNS obtenue par Alamo et al., 2003. Il est remarqué que le profil moyen de la vitesse et le tenseur de Reynolds de ces simulations sont presque identiques sauf au centre du canal. Les raisons possibles de ce désaccord sont discutées.

Dans un second temps, les propriétées de la paramétrisation du modèle de référence de Smagorinsky avec une fonction de Van Driest sont étudiées. Il a été montré que même le modèle de Smagorinsky standard produit des résultats relativement corrects pour l'écoulement considéré. Néanmoins, sur les spectres d'énergie au voisinage de la paroi, une dissipation excessive des échelles résolues a été retrouvée.

Ensuite, la modélisation du transfert inverse d'énergie par une viscosité turbulente négative a été étudiée. Le travail théorique de [START_REF] Dubrulle | Low-viscosity lattice gases[END_REF] a montré que, pour les nombres de Reynolds suffisamment grands, l'écoulement cisaillé produit naturellement une viscosité turbulente négative. Cette viscosité est à l'origine des instabilités des grandes échelles, donc autorise le transfert d'énergie des petites échelles vers les grandes. Ces résultats sont consistants avec les études de Härtel et al. sur l'écoulement de canal. Härtel et al. ont montré, par des tests a priori à partir des données de DNS, que le transfert d'énergie produit par la partie fluctuante du tenseur sous-maille peut être négatif dans la zone tampon (5 < y + < 60). Son intensité dépend de la taille du filtre. Les tests numériques ont montrés que cette méthode permet d'augmenter l'énergie cinétique de l'écoulement près de la paroi. Néanmoins, le modèle appliqué de la viscosité turbulente ne permet pas de prendre en compte le caractère local de l'écoulement, ce qui rend difficile l'obtention de résultats quantitatifs satisfaisants.

L'influence des paramètres du modèle a été étudiée. Il a été montré que l'intensité du vecteur l est sensible à la variation du temps de friction τ f , des viscosités turbulentes µ t et ν t et du temps de corrélation de la force stochastique τ stf . Par rapport au temps de corrélation τ stf de la force stochastiqe, la force turbulent l s'est montrée plus sensible près de la paroi qu'au centre du canal. Les changements plus importants de la force turbulente produits par la variation du terme de friction et de la viscosité turbulente ν t ont été retrouvés au centre du canal. La viscosité turbulente µ t agit par la vitesse résolue sur le vecteur l de la même manière à toutes les distances de la paroi. La paramétrisation de l'équation de la force turbulente n'est pas considérablement affectée par la dissipation des échelles résolues bien que l'amplitude du vecteur l est comparable au terme non-linéaire de la vitesse. Ce résultat nécessite une meilleure analyse du mécanisme local de transfert d'énergie.

Tous les résultats du modèle LES-Langevin ont été comparés avec le modèle de Smagorinsky, des résultats de DNS (simulation directe des équations de Navier-Stokes) et de DNS mal résolues (sur le maillage de la LES). L'influence de la viscosité turbulente µ t sur les statistiques de l'écoulement est plus importante que l'influence de la force turbulente l. L'énergie cinétique moyenne produit par le modèle LES-Langevin est meilleure que celle obtenue par le modèle de Smagorinsky. Par contre, la dissipation excessive des petites échelles résolues au voisinage du maximum de l'énergie cinétique (13 < y + < 17) est persistante quelque soit la paramétrisation de l'équation de la force turbulente.

L'influence du maillage et du filtre sur la modélisation a été également étudiée. Il a été montré qu'un filtre discret dissipe considérablement les échelles résolues. Cela rend difficile la discussion sur les performances du modèle et sur le transfert d'énergie. En même temps, l'application d'un tel filtre est utile pour la modélisation de l'écoulement dans une géometrie complexe. Le maillage trop grossier dans la direction transversale peut induire l'apparition de structures non-physiques ("superstreaks") et la surestimation du tenseur de Reynolds dans la zone de l'écoulement non-homogène. Pour une résolution suffisamment fine dans la direction transversale (∆ + z ≈ 20), l'effet de superstreaks est négligable. Cela montre la consistance du modèle LES-Langevin par rapport à la limite de la résolution fine de DNS.

Le cinquième chapitre est la conclusion. Premièrement, la performance du modèle LES-Langevin appliqué à un écoulement de canal a été comparée avec les performances des modèles LES connus présentés dans le chapitre 1. Ensuite, la modélisation appliquée par le modèle LES-Langevin a été discutée du point de vue de ces propriétés locales et moyennes. Il a été remarqué que le cas particulier de la turbulence de canal nécéssite de prendre en compte le caractère locale des structures cohérentes. Ces structures sont très importantes dans le bilan d'énergie. En plus, elles ont des propriétés différentes par rapport à la turbulence isotrope où le modèle LES-Langevin donne des meilleurs résultats [START_REF] Laval | A LES-Langevin model for turbulence[END_REF]. De même, l'analyse des modèles LES existants (ex: VMS et SGEM) montre qu'une modèlisation efficace du transfert d'énergie due aux échelles sousmailles est souvent liée à une formulation directe du tenseur sous-maille. Les problèmes rencontrés par l'approche LES-Langevin sont nouveaux pour la méthodologie LES. Les développements futurs du modèle LES-Langevin nécessitent la vérification des hypothèses utilisées dans l'approche actuelle pour mieux comprendre les mécanismes de transfert d'énergie qui doivent être reproduits par la force turbulente.

Les perspectives du travail sont nombreux. Tout d'abord, un travail théorique et une analyse de données de DNS et d'expériences doivent permettre de raffiner les approximations du modèle LES-Langevin et prendre en compte les caractéristiques des structures cohérentes de l'écoulement. Ce travail théorique est interessant pour plusieurs modèles LES de la turbulence de paroi. La puissance de l'approche LES-Langevin se manifeste par la modélisation directe par un vecteur du gradient du tenseur C ij , ce qui ouvre la possibilité de réduire considérablement le temps des calculs en utilisant un maillage assez grossier pour la résolution des équations du modèle.

x

Introduction

Introduction to turbulence modeling

Turbulence is an important and unsolved problem of basic physics. Most of the flows which can be found in nature are turbulent (atmosphere, oceans, rivers, etc.). One can find also many turbulent flows in engineering (in the turbines, different engines, pipe networks, etc.), aerodynamics, transports, ships, submarines, chemical reactions, etc. Practical interest of turbulence study can be drag reduction, control of heat transfer in engineering flows, control of chemical reactions by the variation of turbulent diffusion . . .

From the theoretical point of view, turbulence is a out of equilibrium system, where the number of degrees of freedom is very large. It contains a large number of space and time scales that is difficult to reduce. The small perturbations in the initial turbulent state can produce an uncontrolled deviation in the time evolution of the largest most energetic structures. So, the state of a turbulent flow at a given moment is very sensitive to the history of most of the flow scales. This makes it difficult to reduce the degrees of freedom and to propose a simple description of the turbulent dynamics.

Measurement techniques in the study of the turbulence are complicated and are not yet able to return a complete information about the fluid motion. As mentioned by Ferziger et al. [START_REF] Ferziger | Computational methods for fluid dynamics[END_REF], experimental techniques allow the measurement of time-averaged drag and heat transfer but the pressure fluctuations measurement is rather difficult. Particle Image Velocimetry (PIV) techniques supplies the instantaneous velocity. Nevertheless, it can not provide all information about 3D velocity gradients. Another problem of turbulence measurement is the robustness of the registered data. The highly complicated techniques need a careful use and can considerably contaminate the velocity or temperature measurement. Sometimes, the primary data can not be directly used to analyze or to compute the turbulent statistics because of the measurement noise. A preprocessing of raw data has then to be done, which is also a possible source of the errors. This complicates the theoretical analysis of experimental data. The consequence to the modeling of the turbulent flow is the need for hypothesis which compensate the lack of experimental information. Nevertheless, only experiments can provide the information about the turbulent flows at high Reynolds number and the ultimate reference for the theoretical study and numerical modeling.

Experimental studies of turbulence can be efficiently completed by numerical simulations. In fact, the incompressible Navier-Stokes equations are based on the conservation hypothesis of mass and momentum and can successfully describe the turbulent unsteady motion. The problem is that these equations are very difficult to solve analytically except in some special cases, where the equations can be linearized or considerably simplified. As turbulence is essentially a 3D system, the Navier-Stokes equations have to be solved in 3D space and time. The unstable behavior of turbulent motion does not allow general linearized solutions. Robust results on turbulence with complex boundary conditions can be obtained only from numerical solutions resolving all the turbulent scales from the largest comparable to the physical domain to the smallest (dissipative) ones. Such computation needs large resources, to solve the nonlinear terms of the Navier-Stokes equations. So, the set of turbulent flows which can be described by Direct Numerical Simulation (DNS) is limited to flows with relatively small interval of physically important scales. The level of instability growth in turbulent flows is proportional to the interval of the excited scales and characterized by the Reynolds number.

Nowadays, DNS of the turbulent channel flow can reach Re τ = 2000 for a Reynolds number based on the half height of the channel and the friction velocity (see for example J. C. del Alamo et al. [START_REF] Del Alamo | Scaling of the energy spectra of turbulent channels[END_REF]), which is not sufficient to reproduce the turbulent flows in a wind tunnel. Nevertheless, DNS provides a complete time-resolved information about all the velocity scales.

The domain of DNS applications is restricted by the need for large computer resources. This is the price to pay for a complete information about turbulent scales provided by DNS. Nevertheless, in real industrial applications the domain of interest can be limited by the calculation of mean parameters like the mean shear stress, mean fluctuating velocity RMS, in a statistically steady flow. These quantities do not necessary depend on all turbulent scales, so the DNS provides too much information. In practical computational applications, the Reynolds Average Navier-Stokes (RANS) method was developed. The approach is based on the time-averaging of the Navier-Stokes equations and the separation of flow variables into mean and fluctuating parts. Averaged Navier-Stokes equations are solved on a coarse grid and the fluctuating terms are approximated using subgrid models (the main concepts can be found in the book of Ferziger et al. [START_REF] Ferziger | Computational methods for fluid dynamics[END_REF]). These methods are suitable for engineering purpose, like an optimization of a geometry or an investigation of the turbulent properties of a flow in a pipe network, external aerodynamic flows, etc. RANS simulations are successful in computation of the mean statistics of the flow, which does not require a very powerful computation even for high Reynolds number turbulence.

The compromise between the complete resolution of turbulent scales by DNS and only the mean fields resolution by RANS is Large Eddy Simulation (LES). This approach resolves the spatially filtered Navier-Stokes equations with modeling of the subfilter scales. The length of the spatial filter can vary in different applications but it has to leave a sufficient number of the large scales containing the most energic part of the flow. LES provides an instantaneous information about the fluid motion and is less expensive in computing resources than DNS. This approach has a potential of particular interest for industry. Nevertheless, the LES models are usually less stable especially at high Reynolds number and complex boundary conditions. That is why RANS modeling remains more popular for practical applications.

The fields of investigation for LES are numerous: the estimation and correction of the numerical error, the near wall treatment of turbulence, the modeling of complex physical processes in turbulent flow, the study of the universality of the model parameters, the influence of the resolution and scale separation in turbulent dynamics, etc. The common properties of the LES models are often tested on well-known subgrid scale models like a Smagorinsky eddy-viscosity model. At the same time, new models are proposed in CONTENTS order to reduce some of the problems and to refine the hypotheses of LES models. For example, the implicit LES approach is looking for the unification of the truncation errors and the modeling of the subgrid scales. Other examples, can be found in the numerous improvement of the Smagorinsky model in order to satisfy the near-wall turbulent flow properties or to propose a model more suitable with complex boundary conditions.

The present work is devoted to the development of a new LES-Langevin model. The motivation of the work is to develop a LES model based on the study of the subgrid scale dynamics. The main issue is to model the energy transfer between the resolved and non-resolved scales. This is in of particular interest for wall turbulence models.

Chapter 1

State of the art of Large Eddy Simulation (LES)

1.1 General properties 1.1.

The LES idea

The development of the LES strategy is based on the necessity to model numerically unsteady flows with limited computer resources. Due to these limited resources, all the scales of a turbulent flux cannot be computed at present for a Reynolds number representative of practical applications. Consequently the strategy is to compute only part of the turbulence spectrum directly. The large scale dynamics must be reproduced through the time integration of the momentum equations. The small scales action on the large scales is modeled. It is in this modeling of small scales dynamics that LES saves computation resources. Without a closure model one should resolve all the physical scales down to the dissipation ones. The number of degrees of freedom in a turbulent flow can be estimated as the ratio between the largest (integral) and the smallest (Kolmogorov, or dissipative) scales. The Kolmogorov'41 theory divides all the velocity scales into three parts: -injection scales l 0 correspond to the largest scales, generally comparable to those of the flow; -inertial range scales cascade the energy to the smaller ones; -dissipation (Kolmogorov) scales η, correspond to a range where the energy is dissipated through molecular interaction. The ratio between the largest (l 0 ) and smallest (η) turbulent scales indicates the number of degree of freedom in the flow. From l 0 and a turbulent velocity scale u 0 , a dimensional reasoning allows to estimate the mean turbulent energy dissipation rate ε [START_REF] Frisch | Turbulence: the legacy of A. N. Kolmogorov[END_REF]:

ε ∼ u 3 0 l 0 .
This dissipation rate is then used by Kolmogorov to estimate the dissipation scale η as:

η ∼ ν 3 ε 1/4
, where ν is the kinematic viscosity. Consequently, the number of degrees of freedom increases as a power of the turbulence Reynolds number Re = l 0 u 0 /ν [START_REF] Frisch | Turbulence: the legacy of A. N. Kolmogorov[END_REF]:

N ∼ l 0 η 3 = Re 9/4 . (1.1)
To give an example, on an aircraft flying at cruise speed, the free stream velocity is of the order of 250 m/s. The turbulence velocity scale u 0 is about 10% of the mean velocity, so u 0 ∼ 25 m/s. The integral scale in the wing boundary layer is comparable to this boundary layer thickness, which is about l 0 ∼ 1 cm on a wing of about 2 m cord c. The Kolmogorov scale is then η ∼ 10 -3 m. Considering a computational domain which is (10c) 3 and a spatial resolution of order η, the number of grid points needed to resolve all the scales would be 10 13 . This is far beyond the capacity of any present computer.

The modeling of the small scale dynamics permits to reduce the number of resolved degrees of freedom and so the number of numerical grid nodes. For example, in the above case, if one can model all the scales below 100η the number of grid points goes down to 10 7 , which is more tractable.

Of course, the choice of the cutoff scale, separating the resolved (large) scales from the modeled (small) ones is critical. To do it properly, the separation between l 0 and η must be large enough for an inertial range to exist in between. Equation (1.1) tells us that this is true if Re is sufficiently large. In that case, the turbulent energy spectrum in the inertial range was shown by Kolmogorov to be of the form E(k) = C K 2/3 k -5/3 . In that range, the turbulent structure becomes independent of the large scales of order l 0 and the energy correlation decreases rapidly (as the power 5/3). This opens the possibility of modeling the small scales in a approximative way.

Energy scaling E l ∼ (εl) 2/3 indicates the possibility of the small scale modeling because of their relative small contribution to the total energy bulk. When the separation between l 0 and η decreases, the energy of the modeled scales vanishes in comparison to the energy of the resolved scales. In the limit of maximum scale separation, the modeling becomes less efficient and all turbulence must be resolved. This corresponds to Direct Numerical Simulation (DNS) which is the only solution in that case. If the separation is large enough, the question of the position of the cutoff comes into play. If the cutoff is too close from η, the modeling is simple but the benefit in terms of number of grid points is low. When the cutoff is moved to l 0 , more and more energy is contained in the subgrid scales. The contribution of the modeling scales becomes important and can significantly affect the resolved scales. Any errors on the non-resolved scales modeling will affect directly the large scale dynamics.

There are two ways to look at the small scales modeling. The simplest approach is to resolve the Navier-Stokes (NS) equations on a coarse grid (coarser than η). The truncation error on this coarse mesh can generate non-physical instabilities in the computation. In order to reduce the effect of the truncation error, one needs a regularization procedure representing the subgrid scale action. In this approach, the model cannot vanish in the limit of completely resolved grid. This is due to the fact that there is no a priori criterion for the necessity of the subgrid closure.

The other route to a LES closure is to filter the NS equations and model the subgrid

GENERAL PROPERTIES

scales action on the resolved scales. In this case, in the limit of a fine mesh, we can ask the model to produce exactly the same subfilter action as by direct simulation. This can be a criterion to estimate the quality of the proposed closure. The closure has to vanish in the limit of the well-resolved DNS, when the filter size approaches the grid mesh size.

The distinction between the "truncated" and "filtered" definition of the nonresolved scale models is not an inherent property of the LES models. A lot of closures can be formalized in both representation, but sometimes the results vary significantly as a function of scale separation method or numerical integration algorithm. The slow dissipative numerical scheme without spurious dissipation of the filter permit clarify the subgrid model influence on the resolved scale dynamics.

Subgrid scale (SGS) stress tensor

The "filtered" LES approach uses a filter procedure to separate velocity scales into large and small scales parts. The large resolved scales are obtained by the integration of the truncated NS equations. The small subfilter scales contribution to the large scale velocity dynamics is modeled through the subgrid scale stress tensor. The model performances can be estimated from the following physical and numerical behavior:

• Stability. The turbulent flow simulation provided by the truncated NS equations can be numerically unstable. The unstable character depends on the physics of the flow (boundary conditions, Reynolds number, viscous and compressible properties of fluid, etc.), numerical scheme (spectral, finite difference, energy conservative, etc.) and filter form.

• Large scale velocity statistics (a posteriori estimation). As the dissipation scales and a part of inertial scales are truncated in the LES equations, the energy transfer mechanism below scale separation is completely dependent on the SGS tensor action. Also the fluctuating, temporal and symmetry properties of the large scales are influenced by the subgrid scale tensor action.

• In the limit of a filter characteristic length smaller than the Kolmogorov scale, the SGS tensor should go to zero. The fully resolved DNS simulation is the ideal case of the turbulent flow modeling, so SGS stress tensor should not contribute to the energy scale dynamics.

• A priori space and time characteristics of SGS stress tensor. The intensities and anisotropy of its components should corresponds to experimental or DNS data as function of the boundary conditions and the fluid flow characteristics. Nevertheless, good a priori features does not guarantee sufficient resolved scale dynamics. It is important to correctly reproduce the energy transfer mechanism between the scales.

• Filtering action properties. The explicit filtering could significantly affects the energy transfer and change the SGS stress tensor contribution to the model. The non-commutation of the filter with the time and spatial derivations increases the errors of the model. These errors also interact with the closure model and can depend on the geometry of the flow.

• Numerical errors, provided by the spatial differential operator, filtering procedures, time integration mechanism, mesh resolution, aliasing effect, etc. Modeling of the turbulent flow depends on the practical realization of the approach. The numerical scheme choice correlates with the boundary conditions of the equations.

Why subgrid scale terms modeling?

As it is impossible to simulate the turbulent flow at the coarse resolution without modeling of the SGS stress tensor? Or, can the filtering (additional explicit or implicitly imposed by the mesh) provide the good resolved scale velocity characteristics without SGS stress tensor? Corresponding simulations show the spurious increase in the smallest resolved scale energy [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF]. This is not surprising because of the need, for the truncated NS equations, of an energy cascade mechanism near and behind the scale separation cutoff. It could be supplied by the non-linear interaction mechanism of the resolved and subgrid scales.

Some turbulence characteristics, as for example, skin friction or mean velocity profile, can be well produced by a coarse DNS [START_REF] Meyers | Evaluation of Smagorinsky variants in large-eddy simulations of wall-resolved plane channel flows[END_REF] even if the turbulent kinetic energy transfer between the resolved and subgrid scales are wrong. Nevertheless, these results depend strongly on the particular choice of a Reynolds number, a scale separation and a integration algorithm.

The dynamic equations

We consider the flow of a Newtonian incompressible fluid. The space domain of the flow Ω is infinite or limited by the boundary Γ. The conditions of momentum and mass conservation in Ω leads to the system of Navier-Stokes equations:

(∇ • u) = 0, ∂ t u + (u • ∇)u = -∇p + ν∆u + f , (1.2)
where u = u(x, t) is the fluid velocity, p = p(x, t) is the pressure, ν = const is the kinematic viscosity and f = f (x, t) is a volume force. The density of the liquid is supposed constant ρ = const = 1. The nonlinear term (u • ∇)u represents advection, the first term at r.h.s. represents the pressure influence, the volume force f represents all possible external or internal forces: electric, magnetic, gravity, etc. This system of equations needs initial and boundary conditions:

u(x, t) = u Γ (t), where x ∈ Γ u(x, t) = u 0 (x)
, where t = t 0 .

In the present work, we neglect the volume force f and it will be omitted from now. One can eliminate the pressure from equation (1.2) by application of the divergence operator and using the incompressibility condition:

∆p = ∇ • ((u • ∇)u)
with the Neumann boundary conditions calculated from the u Γ and Eq. (1.2). So, the initial velocity u 0 (x) and the boundary conditions u Γ (t) are sufficient to find solution of

1.1. GENERAL PROPERTIES Eq. (1.2).
Equations (1.2) describe both laminar and turbulent flows. The dominance of the nonlinear terms (u • ∇)u is responsible for the growth of instabilities. When the diffusion terms ν∆u is dominant, the solution stays laminar. Dimensional arguments helps to estimate the ratio between advection and diffusion. If the characteristic velocity and length scales are defined by u 0 and l 0 , the ratio of convective to viscous terms of Eq. (1.2) becomes:

| (u • ∇)u | ν∆u ∼ u 0 l 0 ν = Re,
where Re is Reynolds number. When the Reynolds number is greater than a critical value Re c , the flow becomes turbulent. For lower Reynolds number it remains laminar. In the experiment of Reynolds (1883) with a circular pipe (circular Poiseuille flow), the critical value was found to be Re c ∼ 2000 for the Reynolds number based on the diameter of the pipe and the mean velocity across the pipe [START_REF] Lesieur | Turbulence in fluids[END_REF]. For a plane Couette flow, the Reynolds number based on the distance between the walls and on the entrainment velocity of the moving wall has the critical value Re c = 1000 [START_REF] Lesieur | Turbulence in fluids[END_REF].

The Reynolds number appears as a coefficient of the diffusion terms in the nondimensional form of the NS equations (1.2):

(∇ • u) = 0, ∂ t u + (u • ∇)u = -∇p + 1 Re ∆u, (1.3)
when the coordinates and functions are non-dimensionalized by the characteristic velocity u 0 and length l 0 : u -→ u 0 u, x -→ l 0 x, t -→ l 0 u 0 t, p -→ u 2 0 p. A large Reynolds number tends the diffusion to zero and the advection becomes dominant.

There is no unique definition of characteristic velocity and scale u 0 and l 0 . The choice depends on the flow under study: isotropic turbulence, boundary layer, channel flow, shear flow with adverse pressure gradient, etc. For example, we note the Reynolds based on the Taylor microscale λ [START_REF] Frisch | Turbulence: the legacy of A. N. Kolmogorov[END_REF]:

λ 2 ≡ u 2 (∇ × u) 2 , Re = λu 0 ν ,
where • is an ensemble average. An other important Reynolds number is the one based on the shear velocity u τ (for constant density):

u 2 τ ≡ τ w ρ = ν∇ ⊥ u ,
where subscripts "⊥" and " " mean the normal and parallel directions to the boundary Γ. The Taylor microscale is particularly meaningful in isotropic turbulence, and the corresponding velocity scale is generally the fluctuating RMS velocity: u 0 = (uu ) 2 . The shear stress velocity is significant for the wall turbulence like channel flow or boundary layer. Then, the characteristic length l 0 can be the boundary layer thickness or the channel height.

In the example of turbulence around an aircraft, the chord and the aircraft wing are reasonable characteristic parameters of the flow.

The analytical solution of system (1.3) is an open mathematical problem except for some simple laminar flows like Couette or Poiseuille flows. The numerical solution of this system is possible and called DNS (Direct Numerical Simulation). Its performance depends on the computer power and rapidly increases with the Reynolds number. The LES approach approximates the system (1.3) by the large scale velocity equations.

Filtering

A filtering operator G defines the scale separation of an arbitrary space-time function f (x, t):

f = G * f, f = f -f , (1.4)
where f (x, t) and f (x, t) are the large and small scale parts of f (x, t).

The LES equations are derived for resolved velocity ū by the filtering of NS equations. The filter must satisfy some properties [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF]:

• To be sufficiently derivable: G(x) ∈ C n . This is important not to produce additional error when computing ∂ n i ū.

• To be constant conservative: G * c = c where c is a constant. Consequently, the filter must be normalized:

G * 1 = 1.
• To be a linear operator:

G * (f 1 + f 2 ) = G * f 1 + G * f 2 .
• In the limit of small G-function support ∆, the convolution filter should became a unit operator. For convolution filter this means: lim

∆→0 G(x, ∆) = δ(x).
The solution of truncated equations should reproduce the large scale velocity of the complete system. The most used filter is the space convolution filter:

[G * f ](x) ≡ Ω G(x, y)f (y) dy, (1.5)
where the integration is done on the whole Ω domain. The description of more exotic filters like time or differential filters can be found in [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF].

The performances of LES model depends on the filter properties. The filters can be classified by the criteria [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF][START_REF] Vreman | The filtering analog of the variational multiscale method in large-eddy simulation[END_REF]:

• Commutation with the spatial derivatives:

G ∂ ∂x i = ∂ ∂x i G.
• Property of scale separation. A Filter is a projector if G 2 ≡ G * G = G. This filter supply a perfect scale separation, which is interesting to estimate the proper closure contribution to the modeled turbulence.

• Symmetry of the flow. Inhomogeneous filter has a form and characteristic length which depends on the space: ∆ G = ∆ G (x). Homogeneous filter has a constant characteristic length ∆ G . The homogeneous filter is naturally applied in the homogeneous directions of turbulence.

• The filter is conservative, if it does not affect the integral value of an arbitrary nonfiltered field f :

Ω [G * f ](x) dx = Ω f (x) dx.
The filter of consistent LES closure must respect this property.

• The filter is statistically consistent f ensemble = f ensemble . It is true if the volume integral is the appropriate representation of the ensemble average and the ergodic hypothesis holds.

• No singularity production. The smooth character of the filters could be formalized in not increasing the global maximum, neither decreasing the global minimum of the variable: ∀f : min

Ω (f ) ≤ G * f ≤ max Ω (f ).
• The second smoothing property

Ω [Gf ] 2 dx ≤ Ω f 2 dx
signifies not increasing the L 2 norm of the signal by filtering.

All LES models have at least one filter which corresponds to the discretization of equations (1.3) on a coarse mesh. This discretization acts as an implicit filter because the LES grid (or, filter) size is much larger than Kolmogorov scale. In DNS this filtering has no physical effect because all the flow scales are represented on the mesh. Now, an explicit filter will contaminate the modeled flow if it does not verify the above mentioned criteria. The best choice of filter depends on the physics of the flow (its symmetry and boundary conditions) and the numerical algorithm. The ideal filter should respect all the criteria described previously which is usually not possible in practice. So, the problem is to find a more appropriate filter and to understand its influence on the modeled flow.

Leonard and double decomposition of the SGS stress tensor

In a first step, we can look at the general LES equations in the case of one scale separation. This separation can be a special operator acting on all scales of the turbulent field or a projection of the momentum equations to the coarse LES mesh. Assuming the filter-derivative commutation, we apply the filter (1.4) to the continuity and momentum equations (1.3):

∂ j ūi = 0, ∂ t ūi + ∂ j N ij = -∂ i p + 1 Re ∂ 2 j ūi -∂ j τ ij .
The SGS stress tensor τ ij represents the terms which cannot be explicitly calculated in the filtered system:

τ ij = u i u j -N ij .
The nonlinear term N ij can be defined as [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF]: N ij = ūi ūj or N ij = ūi ūj . The Leonard decomposition (here "LCR") of nonlinear terms does not take into account the reduction of scales below the effective filter size:

u i u j = ūj ūi + τ ij , where τ ij ≡ u i u j -ūi ūj = L ij + C ij + R ij .
The Cross-Stress Tensor C ij :

C ij ≡ ūi u j + u i ūj corresponds to the nonlocal resolved-subgrid scale velocity interaction. The Subgrid Reynolds Stress Tensor R ij : R ij ≡ u i u j is responsible for the local subgrid velocity interaction. Leonard Tensor L ij :

L ij ≡ ūi ūj -ūi ūj
represents the aliasing modes situated below the filter characteristic scale. An other possible definition of the SGS stress tensor is called the double decomposition [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF], for which N ij = ūi ūj :

u i u j = ūj ūi + τ ij , where τ ij ≡ u i u j -ūi ūj = C ij + R ij .
The relative intensity of C ij , R ij and L ij tensors depends on the scale separation and the length of the spectral range. The main interest of LES is to separate the scales in the inertial part of the spectra and to integrate the equations for the most energetic scales. Thus, in practice, C ij is larger than R ij and is more important to model. Nevertheless, the relative intensity of C ij and R ij depends on the ratio between the integral, scale separation and dissipation lengths and boundary conditions.

If the mesh sampling is the only filter, both nonlinear terms formulations are identical: ūi ūj = ūi ūj . In this case both the Leonard and the double decomposition of the nonlinear terms gives the same results.

It can be noted that the Leonard stress can be exactly computed if the grid size is two times smaller than the explicit filtering scale.

Explicit filtering formulation

In general, the filter operator is different from the mesh projection. The LES models can use both scale separation by an explicit filter and by the mesh sampling. In this case, the filter length ∆ LES must be larger than the cell size ∆ mesh . This divides the spectral 1.1. GENERAL PROPERTIES range into three intervals, which introduce more hierarchical organization of the velocity scale interactions.

The definition of the LES model in the case of an explicit scale separation, can be found in [START_REF] Domaradzki | Direct modelling of subgrid scales of turbulence in large eddy simulations[END_REF][START_REF] Domaradzki | Large eddy simulations using the subgridscale estimation model and truncated Navier-Stokes dynamics[END_REF]. The total velocity can be represented as a sum of grid-represented (u L i ) and subgrid (u S i ) components:

u i = u L i + u S i .
The grid-represented velocity is also named truncated velocity. The filtering and meshsampling applied to the NS equations (1.3) gives:

∂ t ūL i + ∂ j (u i u j ) L = -∂ i pL + 1 Re ∂ 2 j ūL i .
The filtered and truncated product of the velocities can be rewritten as:

(u i u j ) L = (ū L i ūL j ) L + τ sim ij + τ phy ij ,
where

τ sim ij ≡ (u L i u L j ) L -(ū L i ūL j ), τ phy ij ≡ [u L i u S j + u S i u L j + u S i u S j ] L . (1.6)
The subgrid scale stress is expressed as the sum of two terms:

τ ij = τ sim ij + τ phy ij .
The first one has the form of a similarity model (see below). The second one accounts for the physics of the nonlinear interactions between the mesh-represented and subgrid scale components.

The explicit filtering formulation has some specificities in comparison to the mesh implicit filter approach. First of all, it allows to separate the resolved scale (ū i ) range from the subgrid scales (u S i ) by an intermediate scale gap (u L iūi , subfilter scales). The intensity of the distant scale interaction decreases with the extension of the scale separation interval. Thus, the direct interactions between the largest and smallest scales can be neglected for a sufficiently large intermediate scales range. The energy is transfered to the subgrid scales by a two-step process through the intermediate scales.

Secondary, the explicit filtering can eliminate instabilities in the subfilter mesh-represented scales. In addition to that, if the filter is not a projection it has a dissipation action, which depends on its form. These properties allow to damp the possible instabilities of the SGS closure.

Physical and mathematical constraints

Physical and mathematical constraints describe some a priori analytical properties of the subgrid closure, as was discussed by C. Fureby et al. in [START_REF] Fureby | Mathematical and physical constraints on large-eddy simulations[END_REF] and S. Ghosal in [START_REF]Mathematical and physical constraints on large-eddy simulation of turbulence[END_REF]. This does not mean that each model which respects these constraints correctly predicts the energy transfer between the resolved and the subgrid scales. With the development of LES, a number of the necessary conditions of the SGS closure construction have been derived with respect to the flow symmetries and a priori features of the SGS stress tensor τ ij , which depends also on the filter properties. Once again, these constraints do not ensure satisfactory predictions. The physical constraints can be of two types:

• first: to force the model to have the same invariance properties (such as reference frame invariance) as the NS equations (1.3)

• second: to ensure the realizability condition, which means that the modeled SGS tensor has to respect the same mathematical properties as the tensor itself.

The mathematical constraints of LES approach are linked mainly to the commutativity of the derivatives and filters. The commutation depends on the boundary conditions of the simulated flow and the filter form.

Symmetry constraints

The principle of symmetry means the invariance of the physical process with respect to changes of coordinates. Let us consider an infinite space domain. The filtering procedure of the LES equations and the subgrid scale closure should not violate the basic symmetries of the NS equations:

• translational symmetry;

• rotational symmetry;

• parity invariance;

• Galilean invariance;

• scale invariance.

The most common filtering procedure in LES is the convolution, Eq. (1.5). In this section we will see, what conditions the kernel function must satisfy to fulfill the symmetries of the NS equations. The translational symmetry describes insensitivity of the velocity field to a coordinate translation by a constant vector a

x = x + a, u i (x ) = u i (x).
According to the filtered velocity definition, we have

ū i (x ) = G(x , y )u i (y )dy = G(x+a, y+a)u i (y+a)dy = G(x+a, y+a)u i (y)dy.
The chosen filter will respect the translation symmetry, if G(x + a, y + a) = G(x, y).

Consequently, the kernel function should depend only on the separation between the two points:

G(x, y) = G(x -y).
The rotational symmetry corresponds to the invariance with respect to any rotation of the coordinate system:

x i = A ij x j .

GENERAL PROPERTIES

The norm conservation of the rotation fixes the operator's modulus value |det(A)| = 1. The filtered velocity space is defined from convolution integral:

ū i (x ) = G(x , y )u i (y )dy = G(Ax, Ay) | det(A) | u i (Ay)dy.
The conservation of rotational symmetry by the LES equations restraints the arbitrary kernel function to functions which do not depend on the coordinates rotation:

G(Ax, Ay) = G(x, y).
Coupled with the translation symmetry condition, this leads to kernel functions depending only on the modulus of the separation between the two points:

G(x, y) = G(| x -y |).
The parity invariance means insensitivity to the left-right change of the coordinate system. The coordinate change is described by the matrix A ij with a negative value of the determinant:

x i = A ij x j with | det(A) |= -1.
This constraint is automatically satisfied by filter kernel functions which respect the translation and rotational symmetry conditions:

G(x, y) = G(| x -y |).
It must be noted that the SGS tensor is not a tensor (does not respect the rule of coordinate transformation:

τ ij = dx i dx k dx j dx l τ kl = A ik A jl τ kl )
unless the filter has a spherical symmetry. Nevertheless, as it was remarked by Ghosal [START_REF]Mathematical and physical constraints on large-eddy simulation of turbulence[END_REF], if the filter length is sufficiently large, the particular form of the filter stops to play an important role for the filtered values. An analogy maybe shown with a gas thermodynamic state in a vessel which does not depend on the geometry of the vessel, if the size is sufficiently large compared to the mean free path.

Galilean invariance is with respect to uniform translation of the coordinates:

x i = u 0 i t + x i .
For time-independent kernel functions, the Galilean invariance of the filtered velocity field is equivalent to translational invariance. So, it is satisfied by the previous assumption.

The scale invariance ensures the universality of the fluid dynamic equations with respect to physical units. The equations describing same physical phenomenon in different conditions should have the same non-dimensional form. The applicability of the filtering procedure to the non-dimensional equations imposes the filter kernel function, and its variables to be dimensionless:

G(x, y) = (1/∆ 3 )G 0 (x/∆, y/∆).
The subgrid models have also to satisfy this property.

In the presence of a boundary, certain symmetries of NS equations (1.3) can be broken and the filter becomes dependent on the coordinates:

G(x, y) = (1/∆ 3 )G 0 ( x -y ∆ , x ∆ ).
For example, in a channel flow, the filter width depends on the distance to the wall and satisfies all the above mentioned symmetries in the streamwise/spanwise plane.

Realizability Conditions

If the SGS model has the same a priori properties as the SGS tensor τ ij itself, the realizability conditions are satisfied. Precisely, the realizability conditions are formalized as a statement of positiveness of τ ij eigen values [START_REF] Fureby | Mathematical and physical constraints on large-eddy simulations[END_REF][START_REF]Mathematical and physical constraints on large-eddy simulation of turbulence[END_REF].

In the reference frame of the eigen vectors, the SGS tensor is diagonal:

λ 1 = τ 11 = u 2 1 -ū2 1 ≥ 0, λ 2 = τ 22 = u 2 2 -ū2 2 ≥ 0, λ 3 = τ 33 = u 2 3 -ū2 3 ≥ 0. ( 1
a i = √ ω i and b i = √ ω i f i .
Following [START_REF]Mathematical and physical constraints on large-eddy simulation of turbulence[END_REF] the conditions (1.7) are obtained in the limit

| x i -x i-1 |-→ 0 of the expression : i ω i f i 2 ≤ i ω i i ω i f 2 i = i ω i f 2 i ,
where the normalized filter condition imposes i ω i = 1.

Mathematical constraints

The commutation of the filter and derivation operators is characterized by the commutation error:

C[f ] ≡ df dx - d f dx (1.8)
In an homogeneous flow, the filter width is a constant (∆ = const). By application of the definition (1.5) and using the derivation property of the kernel function

d dx G x -y ∆ = - d dy G x -y ∆ ,
using the Stokes theorem and supposing the vanishing of the surface integrals (which is the case, for example, of a periodic boundary conditions), one can show that the commutation error is zero. So, there is no additional terms in the LES equations produced by the commutation errors.

In the presence of the boundaries, the filter width changes as a function of the distance to the wall (∆ = ∆(x)). The computational mesh must be refined near the boundary to take into account a sufficient number of flow degrees of freedom. Then the commutation error (1.8) becomes non-zero and should be taken into account. The influence of the noncommutation errors on the LES modeling were studied by Ghosal et al. [START_REF] Ghosal | The Basic equations for the large eddy simulation of turbulent flow in complex geometry[END_REF] and Vasilyev et al. [START_REF] Vasilyev | A general class of commutative filters for LES in complex geometries[END_REF]. A complete description can also be found in [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF].

GENERAL PROPERTIES

Let us consider a 1d problem. The convolution integral in a domain [α, β] with a regular grid is defined as: ] is defined as δ(x) = ∆/f (x). For the function φ defined as φ(ξ) ≡ ψ(F (ξ)), the definition (1.9) is:

ψ(ξ) = 1 ∆ β α G( ξ -η ∆ )ψ(η)dη, (1.9
φ(x) = 1 ∆ b a G( f (x) -f (y) ∆ , f (x))φ(y)f (y)dy.
If the boundary integral vanishes, the commutation error (1.8) becomes:

C[ψ] = 1 ∆ b a G f (x) -f (y) ∆ ψ (y)f (y) 1 - f (x) f (y) dy.
This expression can be represented as the power series of ∆ n :

C[ψ] = ∞ n=1 c n ∆ n ,
where the first two coefficients are:

c 1 = f f 2 ∞ -∞ ζG(ζ)dζ and c 2 = 2f f ψ + f f ψ -3f 2 ψ 2f 4 ∞ -∞ ζ 2 G(ζ)dζ.
Here, ζ defined as:

f (y) = f (x) + ∆ζ
supplies the representation of the "y" coordinate as a series of ∆ n with ζ-dependent coefficients:

y = ∞ n=1 y n ∆ n = x + ∆ζ f - ∆ 2 f 2f 3 ζ 2 + • • • .
For the spectral representation of ψ: ψ = ψk e ıkx , with the assumption that ∆ k∆ one gets:

C[ψ] = α f f 3 ∆ 2 ψ + (k∆) 4 , (1.10)
where α = ∞ -∞ ζ 2 G(ζ) dζ.
In the case of a symmetrical function G, the first non-zero term of the series (1.10) corresponds to ∆ 2 . If the commutation error is non zero, it can be included in the LES equations with the help of the derivation operator D:

∂ i ψ = D i ψ.
It was found by Ghosal et al. [START_REF] Ghosal | The Basic equations for the large eddy simulation of turbulent flow in complex geometry[END_REF] that D can be expressed as a series of ∆ n and derivations

∂ ij••• : D i = ∂ i -α∆ 2 Γ ijk ∂ 2 jk + (∆ 3 ) (1.11)
where Γ ijk are the integral coefficients depending on the convolution function G. The change ∂ i -→ D i modifies the LES equations (1.3) and takes into account the commutation error up to the order of the approximation (1.11):

D i ūi = 0, ∂ t ūi + D j (ū i ūj ) = -D i p -D j τ ij + Re -1 D k D k ūi .
An other computation of the commutation error was proposed by Vasilyev et al. [START_REF] Vasilyev | A general class of commutative filters for LES in complex geometries[END_REF]. The computation is based on the same mapping of the uniform grid to the physical domain with a non-uniform mesh size. The development of the convolution function G on the Taylor series of ζ n with the change of the integration variable in (1.9) from ξ to ζ gives:

C[ψ] = ∞ n=1 A n M n (ξ)∆ n + ∞ n=0 B n dM n dξ (ξ)∆ n , (1.12)
where the moments M n (ξ) is defined by

M n (ξ) = ξ-α ∆ ξ-β ∆ ζ n G(ζ, ξ)dξ
and A n and B n are, in general, non-zero coefficients. The radius of convergence of the series (1.12) can be considered as infinite [START_REF] Vasilyev | A general class of commutative filters for LES in complex geometries[END_REF]. The commutation properties are dependent on the mapping function F (ξ) and the moments M n (ξ).

If the filter satisfies the following properties:

M 0 (ξ) = 1 for ξ ∈ [α, β]; M n (ξ) = 0 for ξ ∈ [α, β], and k = 1, • • • , N -1; M n (ξ) = 0 for k ≥ N. (1.13)
the commutation error is of order ∆ N . It is possible to construct discrete filters which satisfy conditions (1.13). For homogeneous (periodic) direction, the commutation error (1.12) is zero.

Energy transfer

A detailed analysis of the energy transfer in turbulence and of its links with the LES approach can be found in [START_REF] Frisch | Turbulence: the legacy of A. N. Kolmogorov[END_REF] and [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF]. The modeling of subgrid scales influences the energy transfer process in the resolved scale flow. Let us consider a turbulent flow forced at the largest scales by an external force f . The energy equation can be obtained by multiplying the momentum equations in Eq. (1.3) by ūi :

∂ t q r = -ū i ∂ j N ij + 1 Re ūi ∆ū i -ūi ∂ j τ ij + fi ūi , (1.14)
where N ij and τ ij are defined through the double or Leonard decomposition of the filtered nonlinear velocity terms, q r is the resolved-scales energy and fi ūi represents the source of energy. Equation (1.14) supposes the vanishing of the surface integral:

Γ pū i n i dx,
where n is a normal vector to the domain Ω. The energy cascade is described by the triple and subgrid terms: ūi ∂ j N ij and ūi ∂ j τ ij , respectively. The dissipative scales are truncated in the LES approach, although the dissipative term contributes to the energy balance.

In the energy conservative system, ∂ t Ω q r dx = 0 and the energy budget equation has the form:

Ω ūi ∂ j N ij - 1 Re ūi ∆ū i + ūi ∂ j τ ij + fi ūi dx = 0, (1.15)
where the diffusion term goes to zero for very large Reynolds numbers.

Let us suppose that the domain Ω is infinite. The spectral form of the energy budget equation (1.15) can be obtained from the LES momentum equation in the spectral space:

∂ t ūi (k) = -F T (∂ j N ij ) -ik i p(k) - 1 Re k 2 ūi (k) -F T (∂ j τ ij ) + fi (k),
where "FT" means "Fourier transform" and ūi (k) is a velocity in the Fourier space. Multiplying this equation and the same equation for the complex conjugate by ū * (k) and ū(k) respectively and adding them, one gets the spectral form of Eq. (1.14):

(∂ t + 1 Re k 2 )q(k) = T (k) + S f (k) ≡ T r (k) + T sgs (k) + S f (k), (1.16)
where the energy transfer function of resolved scales is:

T r (k) ≡ - 1 2 (ū * i (k)F T (∂ j N ij )) + 1 2 k i (p(k)ū * i (k)),
the energy transfer function related to the subgrid scale motion is:

T sgs (k) = 1 2 (ū * i (k)F T (∂ j τ ij ))
and the source term of the energy is:

S f (k) = 1 2 (ū * i (k).
The diffusion term moved to the r.h.s., gives a strictly negative contribution to the energy balance:

- R 3 k 2 q(k)d 3 k ≤ 0.
Thus, the total energy transfer is zero, which corresponds to the energy conservation in the non-viscous fluid:

Ω T (x)d 3 x = R 3 T (k)d 3 k = 0, (1.17)
where the energy transfer in the physical space is defined as:

T (x) ≡ -ū i ∂ j N ij -ūi ∂ j τ ij .
So, the energy conservation writes:

R 3 k 2 q(k)d 3 k = R 3 S f (k)d 3 k
Let us consider turbulence without forcing (f = 0) at very large Reynolds number (Re 1), with a large inertial range. Most of the energy transfer is located in this range with a low dissipation, as dissipation occurs mostly at the Kolmogorov scale. Thus, the diffusion term can be neglected for a large range of wavenumbers. For isotropic turbulence the energy and energy transfer functions depend only on the modulus of the wave vector k : k = |k| and equation (1.16) can be rewritten as:

∂ t q(k) = T (k) ≡ T r (k) + T sgs (k) ≡ T (k), (1.18)
where the functions q(k), T r (k) and T sgs (k) are obtained through integration over the sphere S(k) of radius k:

q(k) = S(k) q(k)dS(k), T r (k) = S(k) T r (k)dS(k), T sgs (k) = S(k) T r (k)dS(k).
The total energy transfer function T (k) can then be interpreted as the energy variation at the scale 2π/k during the time dt. In an energy conservative system ∂ t ∞ 0 q(k)dk = 0, we can always write:

kc 0 T (k) dk = ∞ kc T (k) dk.
The l.h.s. is the energy received through the cascade at the intermediate scale ∼ π/k c and transfered to the smaller scales without accumulation or dissipation at k c . The r.h.s. represents the received energy by the scales below k c . This energy transfer is hold only for the inertial range of turbulence scales, if there is no source of energy and dissipation does not exist. The dissipation becomes important for the largest wavenumbers k ∼ k max , where part of the energy cannot "go back" from the smaller to the larger scales.

The wavenumber k of the energy transfer function is the vectoral sum of three wavenumbers of the velocity fields. As velocity scales with wavenumber (u l ∼ (εl) 1/3 ), important interaction features can be derived from the topology of interacting wavenumbers, as was studied by Kraichnan [START_REF] Kraichnan | Eddy viscosity in two and three dimensions[END_REF]. The mechanism of triadic interactions describes the internal structure of the energy transfer. This interaction when all scales have approximately the same absolute value is known as local. The interaction of distant characteristic length is called non-local. In the LES only the local large scale interactions are resolved directly. The large-small nonlocal interactions and small-small local interactions are modeled by the SGS stress tensor. The small scale interactions act more like a small intermittent force which dissipates the energy through the destruction of the large scale structures. The non-local large-small scale interactions are stronger and their phase is defined by the variation of the rapid small scales. Thus, these interactions contribute mostly to the forward energy transfer, but sometimes they contribute to the backward energy transfer.

The LES approach separates the energy transfer between the resolved scales T r = ūi ∂ j N ij and the subgrid scales T sgs = ūi ∂ j τ ij . The sum of the two gives the total transfer function T (k) = T r (k) + T sgs (k). The properties of T sgs for bounded flow were studied by Härtel et al. [START_REF] Härtel | Galilean invariance and filtering dependence of nearwall grid-scale/subgrid-scale interactions in large-eddy simulations[END_REF][START_REF] Härtel | Subgrid-scale energy transfer in the near-wall region of turbulent flows[END_REF], Domaradzki et al. [START_REF] Domaradzki | Energy transfer in numerically simulated wall-bounded turbulent flows[END_REF] and others. It was shown that the total resolved scale dissipation due to the subgrid scales is positive at all wall distances:

T sgs (h) := R 2 T sgs (k, h) d 2 k > 0 (k = (k 1 , k 2
) is a vector in the plane parallel to the wall, h is a distance to the wall). If a LES model does not produce a correct level of integral dissipation, this leads to a non-physical increase of the resolved scale energy caused by coarse resolution (bottleneck effect, [START_REF] Falkovich | Influence of dissipation on kolmogorov spectra of wave turbulence[END_REF]). In the same time, the backward energy transfer is produced mostly in the buffer layer (see sec.2.3) from the fluctuating part of T sgs . For the shear flow, this effect is aligned with the mean rate of strain

(S ij = 1/2(∂ i u j + ∂ j u i )).
The intensity of backscatter is maximal for the scale separation of large scales and vanish for the small ones. A performant LES model should allow these complex properties of the forward and backward energy transfer.

The SGS stress decomposition into Cross Stress C ij and Subgrid Stress R ij could have also phenomenological explication. Let us represent the LES truncation as the limiting case of the zero spectral gap between the resolved and subgrid scales. For the non-zero separation interval we can imagine the fluid as a two-component matter. The first component is the heavy structures giving the most important contribution to the kinetic energy of the flow. The second component is the small scale rapid motion acting like a noise. The statistics of the small scales are more universal and we can see their action on the energetic scales as a friction. In this representation, the Subgrid stress R ij describes the noisy influence of the "thermal bath". Its modeling by an eddy-viscosity approach seems to be reasonable. On the other side, the Cross Stress C ij is more energetic and represents directly the subgrid-resolved scales interactions. So, the small intensity variations of the subgrid velocity generates an important increase of the cross resolvedsubgrid velocity terms acting like a volume forcing on the large structures trajectory. The spontaneous change of the subgrid velocity direction and intensity leads to a dissipation of the large scale kinetic energy of the system through a stochastic forcing action. So, the Cross Stress tensor action can generate backscatter as well as forward energy transfer.

Universality of small scale motion

In real turbulent flows the injection of energy into the turbulent motion is done at large scales. This scales are often comparable to the mean flow size. It can be the size of the obstacle provoking the instability of the flow passing around, the submarine size in water turbulence, etc. So, the properties of integral scales of turbulence are often imposed by the boundary conditions. Through the energy-cascade process the system "forgets" about these initial conditions. For sufficiently large Reynolds number the inertial interval is established in the intermediate wave number range. The universality of this spectral interval shows that its features are defined only by the internal turbulence interactions and independent of the boundary conditions. Without additional external forcing at the smaller scales, the universality is conserved through all the inertial scale range down to the dissipation scales. This dissipation is supposed to be universal for Newtonian flows. So, only the large integral scales of turbulence are external-force dependent. This allows to hope that a good understanding of the internal turbulence mechanism could provide the necessary information for the subgrid model construction.

Therefore, most of the LES models are based on the inertial range effective cutoff. The modeling of the subgrid scales should not be different for sufficiently large and smooth boundary conditions. The situation changes if we try to model the scales before the beginning of the inertial scale range, as in the case of a Very Large Eddy Simulation methods.

Some remarkable LES models

The central question of LES is how to model the Subgrid Stress Tensor τ ij in order to reproduce at best the physics of the large scale motions, minimizing at the same time the computer resources. The large number of existing SGS closures [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF] indicates that, at the moment, there is not a unique LES model working for a large variety of physical flows and boundary conditions. The development of different LES methodology is motivated by the different physical and mathematical aspects which are of different importance for various flows. The continuous increase in available computer resources, stimulates the modeling activity to take into account more and more complicate features of the fluid motion. Nevertheless, the efficiency of the complicate models is a matter of discussion. Sometimes an old and simple model gives a better prediction of the velocity statistics and energy transfer than a new approach based on a more intelligent background. Also, often two models can each be better on one specific flow property. For example, a good prediction of the Reynolds Stresses does not guarantee the correct spectral behavior. A chronological presentation of some remarkable models will naturally introduce the physical basis of the models and their problem to reproduce the large scale velocity dynamics.

The Smagorinsky model

For sufficiently large Reynolds numbers (Re ∼ 10 4 ) the velocity spectra can be divided in three successive parts: large scale motion, inertial range and viscous scales. The most energetic scales are comparable to the characteristic scales of the mean flow and cannot be universal for all fluid motions. The filtering procedure applied to the total velocity should not significantly perturb large scale dynamics. On the other hand, the inertial range is characterized by the internal turbulent energy transfer and is less dependent on 1.2. SOME REMARKABLE LES MODELS the mean flow. These scales are more universal and it is reasonable to try to estimate their action by Reynolds Stress modeling.

To produce the forward energy cascade the eddy-viscosity concept is introduced. The effective kinematic turbulent viscosity ν t (x) can provide the resolved scale energy diffusion. The Smagorinsky model [START_REF] Lilly | On the application of eddy-viscosity concept in the inertial sub-range of turbulence[END_REF][START_REF] Smagorinsky | General circulation experiments with the primitive equations[END_REF] constructs the additional dissipative terms in order to regularize the truncated Navier-Stokes equations:

τ ij - 1 3 δ ij τ kk = -2ν t (x) Sij , ν t = (c s ∆) 2 | S|,
where ∆ is a characteristic space filter length, Sij is the filtered strain rate tensor:

S ij = 1 2 (∂ i u j + ∂ j u i ) (1.19)
with the norm |S| = (2S ij S ij ) 1/2 and c s is the Smagorinsky constant which must be optimized to dissipate the resolved scale energy. The Smagorinsky model gives good results for isotropic turbulence. In this case the boundary conditions are periodic and there are no interaction of the subgrid scales with walls. By the definition, the Smagorinsky closure models both the Cross Stress and Subgrid Stress tensors using the nonlinear resolved strain-rate terms. As the subgridresolved scale interaction and subgrid-subgrid scale interaction are less energetic than the nonlinear resolved terms, the modeling constant adjusts the energy transfer level. We can also see the c s ∆ term as an effective filter size. It has been noted by Sagaut [67] that LES models have to vanish at the limit of DNS, which is not the case of Smagorinsky model. So, this model is not a priori consistent.

The optimal value of the constant depends on the filter size and the boundary conditions. An analytical prediction of this constant as a function of the turbulent spectra, filter size and boundary conditions can be found in Hughes et al. [START_REF] Hughes | The multiscale formulation of large eddy simulation: Decay of homogeneous isotropic turbulence[END_REF]. For the turbulent channel flow the constant was estimated as c s = 0.1, [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF].

Smagorinsky model with Van Driest Damping

After isotropic turbulence, the next step in complexity is the plane channel flow which was studied to investigate the ability of the Smagorinsky model to predict wall turbulence. In fact, a problem arises near solid boundaries because the model over-estimates the asymptotic value of the turbulent viscosity. In order to minimize its action (or change the characteristic length of the additional dissipation terms) different types of wall dumping were introduced [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF]. Improved asymptotics were obtained with Van Driest damping f vd , which reduces the filter size of the model (∆):

∆ → ∆f vd ≡ ∆(1 -exp(-y + /A + )).
(1.20)

Here y + = y u τ /ν is the standard non-dimensional distance from the wall, and

A + = 25
is the Van Driest constant.

Even with this improvement, the main drawback of the eddy viscosity models is that they correlate poorly with the real stress tensor.

Dynamic model

The previous ad hoc regularization of the near-wall behavior of the Smagorinsky model can be replaced by a computation of the Smagorinsky constant as a function of resolvedscale velocity. The forward energy transfer should be less important near the boundaries. The local energy dissipation provided by the subgrid scales is defined as:

sgs = τ ij Sij (1.21)
in the case of a divergence free model (i.e. model providing ∂ i (∂ j τ ij ) = 0). On average, the dissipation must have a positive value ( Ω (x) dx > 0) that corresponds to the integral energy transfer from the large to the small scales at any time. Nevertheless, it is possible to have locally a backward energy transfer from the small to the large scales. The time evolution of this backscatter process is intermittent. The subgrid scale dissipation of the Smagorinsky model sgs = ν t (x)| S| 2 has the same sign as the turbulent viscosity or the Smagorinsky constant. In the previous consideration the constant was always positive, such as only forward energy transfer was produced. If a new procedure allows locally inverse energy transfer, the square of the Smagorinsky constant (c 2 s ) should become a new function (C) of space and time C ≡ c 2 s = C(x, t) which can have locally a negative value. This model was first proposed by Germano et al. [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF]. Here we use the Lilly version of the Dynamic model [START_REF] Lilly | A proposed modification of the Germano subgrid-scale closure method[END_REF]. The mathematical formulation is based on the explicit filtering of the grid-represented scales f = G * f L which is additional to the grid filter f L = G L * f . The width of the filter G is assumed to be larger than that of the grid filter G L . Let's apply G * G L to the NS equations (1.3) (we suppose the commutativity of the filters and derivation operators):

∂ t ūL i + ∂ j (ū L i ūL j ) = -∂ j pL -∂ j T ij + 1 R e ∆ū L i , (1.22)
where the subgrid-scale stress tensor is now

T ij = u i u j L -ūL i ūL j .
The relations between the filtered tensor T ij and subgrid stress tensor

τ ij = (u i u j ) L -u L i u L j is L ij ≡ T ij -τij = u L i u L j -ūL i ūL j . (1.23)
in which the modeling of the divergence free parts of the subgrid τ ij and subfilter T ij tensors are based on the standard Smagorinsky model:

τ ij -1 3 δ ij τ kk = 2C∆ 2 |S L |S L ij , T ij -1 3 δ ij T kk = 2C ∆2 | SL | SL ij .
(1.24)

The elements of the L ij tensor are the resolved components of the Stress Tensor associated

SOME REMARKABLE LES MODELS

with scales of motion between the test scale and the grid scale. Using the Smagorinsky closure, we can approximate the divergence free part of L ij by:

L ij - 1 3 δ ij L kk = 2CM ij ,
where

M ij = ∆2 | SL | SL ij -∆ 2 |S L |S L ij . (1.25)
Let us define the modeling error as:

e ij = L ij - 1 3 δ ij L kk -2CM ij (1.26)
and impose the minimization of e ij with respect to the C constant: ∂ C e 2 ij = 0. This leads to the Germano identity for the dynamic constant:

C = 1 2 (L ij M ij /M 2 ij ). (1.27)
This model shows a better near-wall behavior than the Smagorinsky model with Van Driest damping. The value of the dynamic constant far from the boundary is closed to the standard Smagorinsky constant for isotropic turbulence. Near the wall it provides a good asymptotic behavior in comparison to DNS data. The model does not need an additional ad hoc procedure like the Van Driest damping. The stability of the Dynamic model depends on the practical realization. If it is possible by the symmetry of the flow motion, the constant should be average in the homogeneous directions in order to not provide too much negative dissipation effect. But even in the plane channel flow, where the dynamic constant depends only on the wall-normal direction and time (C = C(y, t)), it is necessary to clip all negative values to zero. So, the theoretical prediction of backscatter stays unachieved.

The additional explicit filtering has no direct influence on the simulated fields and statistics. It is used only for the dynamic constant estimation. So, the commutation and symmetry problem which could arise for the explicit filtering procedure mostly are limited to the implicit numerical grid filtering.

One should also notice that it is not possible to model effectively with a single universal constant the variety of phenomena of turbulent flows.

Another important assumption was done about the scale-invariance of the Germano identity. In fact, it is not clear that for the modeling of the subgrid stress τ ij and subfilter stress T ij (Eq. (1.24)) one can use the same universal function C(x, t). If the physics imposes the importance of the certain modes (ex. the integral scales becoming comparable to the wall distance near the boundaries, or the characteristic filter scale ∆ approaches transitional scales) the Germano identity should depends on the filter size support ∆.

Similarity model

Another well-known model is the Similarity model introduced firstly by Bardina [START_REF] Bardina | Improved subgrid scale models for large eddy simulation[END_REF]. This model explicitly separates the resolved scales into large resolved scales and small resolved or large subfilter scales with the help of an explicit filter. The hypothesis supposes the similar behavior of the SGS tensor specific elements with the subfilter ones. It was successfully verified in experiments on the jet, the plane wake and rapid straining turbulence (see [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF], p. 202). This hypothesis is statistically consistent with the energy cascade mechanism. Let us consider three consecutive spectral intervals. Energy transmission from the lowest to the fastest scales is done through the intermediate scales. So, the dynamical mechanisms of the nonlinear interaction are similar for distinct intervals from the inertial range.

The similarity model estimates the total velocity by the grid-represented velocity and uses the following definition of the SGS stress:

τ ij = τ sim ij = u L i u L j -ūL i ūL j .
(1.28)

It was shown that the filter choice can be significant for simulated velocity features. In order to adjust the dissipation by the subgrid terms, a corrective constant can be added. Also, the model allows the inverse energy cascade. In fact, for a non-divergent model (like Self-Similarity) the most convenient detection of backscatter due to the subgrid scale dynamics corresponds to negative values of sgs (x) ≡ u L i ∂ i τ ij [START_REF] Vreman | The adjoint filter operator in large-eddy simulation of turbuleny flow[END_REF]. As we can see, the Bardina closure (1.28) does not fix the sign of the energy transfer.

The implicit mesh-projection filter may be changed to an explicit one. It does not change the closure formulation.

Mixed model

An important characteristic of the subgrid terms is their a priori behavior. One can remark that eddy-viscosity models have the divergence free form in opposition to the SGS tensor definition which is divergence-free only in the case of the filter-derivative commutation. One can expect to have a good a priori estimation of subgrid scale stress as well as good a posteriori tests. The comparison of the Smagorinsky SGS tensor components and the truncated DNS data shows a significant disagreement for the channel flow [START_REF] Piomelli | Model consistency in large eddy simulation of turbulent channel flow[END_REF]. One can expect that this disagreement correlates with the too dissipative character of the model in the wall regions.

The similarity model has better a priori statistics and spectra of the subgrid scale stress, than the Smagorinsky closure. Nevertheless, the model overestimates the backward energy cascade and this makes it instable. The total value of energy transfer = Ω (x) dx has no definitive sign. So, the stability of the model depends on the filtering and boundary conditions, which is not satisfactory. In order to regularize this problem, an additional Smagorinsky viscosity term was added:

τ ij - 1 3 δ ij τ kk = 2ν t S L ij + C B (L ij - 1 3 δ ij L kk ),
where

L ij = u L i u L j -ūL i ūL j 1.2

. SOME REMARKABLE LES MODELS

This is the definition of the mixed model. C B is a modeling constant.

The Smagorinsky part of the mixed model provides only dissipation by subgrid scales. At the same time, the Bardina term can produce negative and positive energy transfer. Physically, the nonlinear subgrid scale interaction can be responsible only for the dissipation mechanisms. The nonlocal interaction of resolved and subgrid velocities are responsible for both forward and backward energy transfer. This shows that the cumulative dissipative effect of the Cross Stress and Subgrid Reynolds tensors is modeled by the cumulative dissipative action of the Smagorinsky and similarity terms. On the other side, the inverse energy cascade can be due only to the nonlocal interactions and its modeling is provided only by the similarity term. This can fix the two modeling constants of the mixed model.

The Planar-Averaged Scale Dependent Dynamic model

The refinement of the dynamic model was done by the Porté-Agel et al. [START_REF] Porté-Agel | A scale-dependent dynamic model for large-eddy simulation: application to a neutral atmospheric noundary layer[END_REF] for atmospheric boundary-layer flow simulation. Let us consider the support of the additional explicit filter as a function of the implicit filter support and a coefficient: ∆ = α∆. The parameter α is constant for a regular grid and could vary with the refinement of the mesh in the regions with increasing local Reynolds number. The simple dynamic model does not account for the possible scale dependence of the dynamic constant c s . To correct this behavior one can generalize the relation (1.24) to the case of two different dynamic constants:

τij - 1 3 δ ij τkk = 2c 2 s,∆ ∆ 2 |S L |S L ij , (1.29)
T ij - 1 3 δ ij T kk = 2c 2 s,α∆ (α∆) 2 | SL | SL ij , (1.30)
and express the M ij tensor (defined by Eq. (1.25)) as:

M ij = 2∆ 2 [α 2 β| SL | SL ij -|S L |S L ij ],
where β ≡ c 2 s,α∆ /c 2 s,∆ can supply information about the scale dependence of the dynamic constant.

The model is based on the assumption of that a power-law behavior describes the scale dependence of the coefficient: c s,∆ = c s,α∆ ( ∆ α∆ ) φ , so β = α 2φ . The Germano identity, written between scales ∆ and α 2 ∆ yields:

C ∆ = c 2 s,∆ = Q ij P ij P ij P ij , (1.31)
where Q ij is the resolved stress tensor between ∆ and ∆ ≡ α 2 ∆. Q ij and P ij are given by:

Q ij = u L i u L j -u L i u L j , P ij = 2∆ 2 ( |S L |S L ij -α 4 β 2 | S L | S L ij ).
(1.32)

The system of two Smagorinsky identity equations (1.27) and (1.31) yields a fifth-order polynomial equation in β. In non-homogeneous flows one must find the β-roots at each points of the numerical grid. The number of algebraic equations to solve can be reduced by the averaging of the tensor products L ij M ij and M ij M ij (the same for Q ij P ij and P ij P ij ) in the flow homogeneous directions. In the case of planar averaging one obtains a planar-averaged scale depended model. The averaging serves to stabilize the constant C ∆ which could diverge locally due to the small values of M ij M ij and P ij P ij .

Lagrangian Scale Dependent and Independent models

In order to generalize the Germano-like models for complex geometries, Meneveau et al. [START_REF] Meneveau | A Lagrangian dynamic subgrid-scale model of turbulence[END_REF] proposed to perform the integrations in time to compute the model constant. The corresponding model is called Lagrangian-averaged scale-invariant model (LASI). The characteristic time of subfilter and subgrid scale interactions was approximated independently of the space. The scale-dependent version of the model was applied to the boundary layer flow with and without rough surface [START_REF] Bou-Zeid | A scale-dependent Lagrangian dynamic model for large eddy simulation of complex turbulent flows[END_REF]. Let us define the parameter α for two explicit filters in comparison to the implicit mesh size as α = 2 and α = 4. The minimization of the weighted time average of the local error contraction e 2 ij over pathlines yields the Germano identity coefficients C 2∆ ≡ c 2 s,2∆ and C 4∆ ≡ c 2 s,4∆ . The weighted time average can be written as:

E = t -∞ e ij (z(t ), t )e ij (z(t ), t )W (t -t ) dt, (1.33)
where z(t ) are the previous positions of the fluid elements and W (τ ) is a relaxation function decreasing when τ goes up. The minimization of E with respect to c s,2∆ and c s,4∆ provides the system of integral equations:

       ∂E ∂C 2∆ = t -∞ 2e ij ∂e ij ∂C 2∆ W (t -t ) dt = 0 ∂E ∂C 4∆ = t -∞ 2e ij ∂e ij ∂C 4∆ W (t -t ) dt = 0
with the solutions:

     C 2∆ = J LM J M M J LM = t -∞ L ij M ij (z(t ), t )W (t -t ) dt , J M M = t -∞ M ij M ij (z(t ), t )W (t -t ) dt . and      C 4∆ = J LP J LP J LP = t -∞ L ij P ij (z(t ), t )W (t -t ) dt , J P P = t -∞ P ij P ij (z(t ), t )W (t -t ) dt .
The weighted function is approximated by an exponential:

W (t -t ) = (1/T )exp(-(t -t )/T )
with the time scale T .

The relaxation transport equations thus obtained for J LM and J M M are:

DJ LM Dt = ∂J LM ∂t + u L • ∇J LM = 1 T 2∆ (L ij M ij -J LM ), DJ M M Dt = ∂J M M ∂t + u L • ∇J M M = 1 T 2∆ (M ij M ij -J M M ).
Using first-order numerical approximations in space and time, these equations can be discretized. The resultant formulation to update from time step "n" to time step "n+1" at the position x is:

J i+1 M M (x) = 2∆ (M ij M ij ) n+1 (x) + (1 -2∆ )J n M M (x -u Ln ∆t), (1.34) J i+1 LM (x) = H{ 2∆ (L ij M ij ) n+1 (x) + (1 -2∆ )J n LM (x -u Ln ∆t)}, (1.35)
where

2∆ = ∆t/T n 2∆ 1 + ∆t/T n 2∆ , T n 2∆ = 1.5∆(J n LM J n M M ) -1/8 (1.36)
and H is a ramp function:

H(x) = x, if x ≥ 0 10 -32 otherwise,
where ∆t is the time discretization. For C 4∆ Germano identity the same equations can be written with 2∆ → 4∆,

L ij → Q ij , M ij → P ij , 2∆ → 4∆ .
In the work of Bou-Zeid et al. [START_REF] Bou-Zeid | A scale-dependent Lagrangian dynamic model for large eddy simulation of complex turbulent flows[END_REF], in contrary to the Plane-averaged scale-dependent (not Lagrangian) model, the assumption of the scale invariance was supposed

C 4∆ /C 2∆ = C 2∆ /C ∆ . The coefficient β was defined as β = C 4∆ /C 2∆ which fixes C ∆ = C 2∆ /β.
Numerical experiments have shown the unstable behavior of the dynamic coefficient for small Germano's identity ratio (β → 0). This limit violates the local viscosity stability conditions [START_REF] Bou-Zeid | A scale-dependent Lagrangian dynamic model for large eddy simulation of complex turbulent flows[END_REF]. The clipping limit was chosen considerably lower than the physically expected limitation: β ≥ 0.125. The other limit of large dynamic constant ratio β → ∞ has not shown any instability problem. Finally, the model coefficient to be used was calculated for each space points and time steps:

C = C 2∆ max(β, 0.125) = J LM /J M M max( J QM J M M J M M J LM , 0.125) .

An Approximate Deconvolution Model (ADM)

The Approximate deconvolution model [START_REF] Domaradzki | Direct modelling of subgrid scales of turbulence in large eddy simulations[END_REF][START_REF] Stolz | An approximate deconvolution model for large-eddy simulation with application to incompressible wall-bounded flows[END_REF] is based on a defiltering procedure. It reconstructs the truncated velocity field u L i from the filtered velocity ūL i :

u L i = G -1 ūL i .
The ADM approximation of the total velocity field u i :

ũL i = Q N * ūL i = N ν=0 (I -G) ν * ūL i .
The series converges if (I -G) < 1. The leading-order error term of the deconvolution is function of the filter kernel form:

δu L i = ũL i -u L i = C j (-1) (N +1) ∆ r(N +1) ∂ r(N +1) u L i ∂x r(N +1) j .
Here the coefficient C j = C j (∆) depends on the filter kernel which is of order r.

The modeling of τ phy ij (Eq. 1.6) can be done by the traditional eddy-viscosity approach or by a mathematical regularization procedure. The second choice provides the numerical dissipative mechanism to preserve the accuracy of large-scale dynamics. In the ADM model, the relaxation term was added to regularize the τ phy ij term [START_REF] Domaradzki | Direct modelling of subgrid scales of turbulence in large eddy simulations[END_REF]:

Φ u = ∂ j τ phy ij = χ u (I -G 2 ) * ūi .
(1.37)

The parameter χ u is an inverse characteristic dissipation time of the filtered solution ūi to some filtered solution G 2 * ūi where the G 2 filter has a larger support that the filter G.

In [START_REF] Stolz | An approximate deconvolution model for large-eddy simulation with application to incompressible wall-bounded flows[END_REF] G 2 was constructed using Q N and G operators:

G 2 = Q N * G.
The regularizing term was calculated dynamically through the iteration procedure (on χ u and χ u0 ) in order to conserve the kinetic energy level:

χ u = χ u0 F 2 (x, t + τ )| χu=0 -F 2 (x, t) F 2 (x, t + τ )| χu=0 -F 2 (x, t + τ )| χu=χ u0 ,
where

F 2 (x, t) = Φ u (x + ∆x, t) -Φ u (x, t) 2
∆x =h is a 3D local second-order structure function, calculated on the neighbors of x point.

Subgrid scale estimation model (SGEM)

Another existing SGS modeling strategy is based on an estimation of the SGS velocity u (x, t) and the reconstruction of the full velocity u = ū+u . The SGS tensor is computed

SOME REMARKABLE LES MODELS

by the definition τ ij = τ sim ij . The SGS velocity can be represented as a synthetic field or through an evolution equation. One of the existing SGS closure is the Subgrid scale estimation model developed by Domaradzki et al. [START_REF] Domaradzki | The subgrid-scale estimation model in the physical space representation[END_REF][START_REF] Domaradzki | A subgrid-scale model based on the estimation of unresolved scale turbulence[END_REF] for the spectral-and real-space representation of the channel flow. . Necessary correction of the model for high Reynolds number flow was developed for the isotropic turbulence [START_REF] Domaradzki | The subgrid-scale estimation model for high Reynolds number turbulence[END_REF]. The model was successively applied to the Rayleigh-Bénard convection [START_REF] Kimmel | Large eddy simulations of Rayleigh-Bénard convection using subgrid scale estimation model[END_REF] and compressible turbulent flow [START_REF] Dubois | The subgrid-scale estimation model applied to large eddy simulations of compressible turbulence[END_REF]. The summary of the hypothesis and a clear presentation of the model can be found in [START_REF] Domaradzki | Large eddy simulations using the subgridscale estimation model and truncated Navier-Stokes dynamics[END_REF].

The approach reconstructs the total grid-represented velocity scales from the filtered scales. The large scales of the approximated velocity are calculated through the deconvolution procedure (ũ i = G -1 ūi ). The deconvolved velocity ũ0

i should have the same value on the coarse mesh as the resolved velocity ū(N) i :

ũ0 i (x n ) = ū(N) i (x n ). (1.38)
In the case of the Simpson's rule for spatial filtering, the equation (1.38) takes the form:

1 6 [ũ 0 i (x n-1 ) + 4ũ 0 i (x n ) + ũ0 i (x n+1 )] = ū(N) i (x n ).
Consequently, the deconvolved velocity is sampled to the finer grid. The subfilter meshrepresented (perturbation) velocity scales is supposed proportional to the small scale nonlinear interaction of the resolved velocity:

ũ i = θN i , N i = N 0 i -N 0 i .
(1.39)

The time scale θ is supposed to be an eddy turnover time which depends locally on the space coordinates. It is calculated from the assumption of the local proportionality of SGS energy to the energy of the smallest resolved scales:

θ = R (ũ 0 i -ũ0 i ) 2 N 2 i .
The proportionality constant R was estimated by Domaradzki et al. [START_REF] Domaradzki | A subgrid-scale model based on the estimation of unresolved scale turbulence[END_REF] R ≈ 0.5 in order to preserve the Kolmogorov inertial range spectra

E(k) ∼ k -5/3 .
In this procedure it is supposed that only a single nonlinear action determines the growth-rate of the subgrid scales. In order to better reproduce the coupling between the resolved and subgrid scales a SGS stress time-estimation procedure was added. It consists in integration of coarsely resolved NS equations on the finer grid. The SGS stress computation is performed according to definition:

τ ij = ũi ũj -ũi ũj .
(1.40)

The finer grid is twice as fine as the rough grid. The total model computation consists in two parts. The first one is an estimation step of the subfilter velocity; the second one is a parallel time resolution of the LES equations with the truncated NS equations:

ū(N) i (t) Estimation Model ↓ ũ(2N) (t) -→ τ ij (t) ū(N) i (t) NS ↓ ↓ LES ũ(2N) (t + ∆t) -→ τ ij (t+∆t) ū(N) i (t + ∆t) NS ↓ ↓ LES . . . . . . ũ(2N) (t + T -∆t) -→ τ ij (t+T -∆t) ū(N) i (t + T -∆t) NS ↓ LES ū(N) i (t + T ) ↓ Repeat Procedure, (1.41)
where T is the time between two initializations.

The SGS stress contribution is evaluated through the NS nonlinear mechanism of large and small scales interaction. The LES velocity solution is "corrected" at each time step. On the contrary, the backward velocity contribution is estimated at the same time. The re-initialization procedure allows to regularize the excessive energy of the velocity scales closed to the cutoff scale. As their action on the LES-evaluated velocity filed is not direct, this does not yield any divergence problems. The re-initialization is repeated at each time T , which is a fraction of the large eddy turnover time.

The model contains three modeling parameters: the constant R, the number of time steps N f in the forced truncated NS equations used to couple the nonlinear interactions between resolved and subfilter scales and the period T of the truncated NS equations initialization. The authors [START_REF] Domaradzki | Large eddy simulations using the subgridscale estimation model and truncated Navier-Stokes dynamics[END_REF] state the low sensitivity of the model to variation of these parameters R, N f and T by a factor 2.

The total model evaluation is approximately 8 times cheaper than the simple LES model. This is mostly due to the truncated NS equations evaluated on the finer grid containing 2 3 times more points. A simplified version of the model was also proposed by the authors in order to reduce excessive computations:

ū(N) i inver. → ũ0 i (t) → ũ i (t) → [ũ i (t) = ũ0 i (t)+ũ i (t)] N-S → ũi (t+T ) → [ũ i (t+T )-ũ i (t)] → ū(N) i (t+T ) (1.
42) Here the SGS-estimated procedure is initialized at each time step. In order to conserve the zero mean subfilter scale velocity the authors use regularization of the deconvolved velocity field [START_REF] Domaradzki | Large eddy simulations using the subgridscale estimation model and truncated Navier-Stokes dynamics[END_REF].

The advantage of this model is the variation of the filter width in order to adapt the solution to the numerical grid. Also, the calculated SGS stress has exactly the same invariance properties as the real SGS stress.

The model does not assume local isotropy of the subgrid velocity. It is advantageous in the case of VLES (Very Large Eddy Simulation) where only the largest scales can be resolved and the effective filter cutoff can be situated before the inertial-interval range.

The problem of the commutation of the filter and derivation is not crucial for the total model closure. The SGS stress is computed from an approximate DNS field and one can filter the derivated velocity rather than derivate the filtered field.

Implicit LES (ILES) in the case of the adaptive local deconvolution model (ALDM)

The description of the model can be found in Hickel et al. [START_REF] Hickel | An adaptive local deconvolution method for implicit LES[END_REF]. The implicit closure uses the truncation error of the NS equations. This approach does not need an additional explicit modeling. As the truncation error depends strongly on the numerical scheme, the implicit LES model is naturally consistent with the discretization method. Let us collect the pressure and quadratic velocity terms of the NS equations in a nonlinear term:

∂ j N ij ≡ ∂ j (u i u j ) + ∂ i p,
and impose this terms in general non-dimensional LES equations:

∂ t ūL i + G * ∂ j N L ij - 1 Re ∆ū L i = -∂ j τ ij = -∂ j τ sim ij -∂ j τ phy ij .
The ILES attempts to change the filtered nonlinear term G * ∂ j N L ij and SGS stress by the combination of new approximation velocity ũL i and the approximated filter-derivative term G * ∂j :

ζ ≡ G * ∂j Ñ L ij ≡ G * ∂ j N L ij + ∂ j τ ij .
The Modified NS equations use this approach:

∂ t ūL i + G * ∂j Ñ L ij - 1 Re ∆ū L i = 0.
The Adaptive local deconvolution method applied to the ILES consists in three steps:

• Local deconvolution procedure: ũL i = G-1 * ūL ; • Computation of the numerical flux function F L ij ;
• Numerical integration and approximation of G * ∂ i by G * ∂i .

The optimization of the filter and numerical scheme parameters is done by the minimization of the non-dimensional numerical viscosity ν + num in comparison with the reference one.

Let us write the kinetic energy equation in the Fourier space:

∂ t E(k) -T (k) + 2νk 2 E(k) = -u * (k)[∂ j τ ij ](k) = -num (k),
where num (k) is the numerical dissipation provided by the numerical error of the energy transfer; E(k) is the kinetic energy and T (k) is the energy transfer function. The subgrid scale dissipation through the eddy-viscosity hypothesis provides:

sgs (k) = 2ν sgs (k)k 2 E(k).
By analogy, one can define the numerical dissipation:

ν num (k) = num (k) 2k 2 E(k)
.

We can eliminate the velocity dependence of the k-vector orientation by the integration of the numerical viscosity over a sphere of radius k = |k|:

ν num (k) = 1 2k 2 E(k) |k| dku * (k)[∂ j τ ij ](k).
The scaling of the viscosity is done in a standard way:

ν + num (k + ) = ν num k k C k C E(k C )
.

The minimized functional for the optimization procedure compares the difference between the numerical viscosity and a reference eddy-viscosity (ex. Chollet viscosity:

ν + (k + ) Chollet = 0.441C -3/2 K (1 + 34.467 e -3.03/k + )): Cost function(parameters) = 1 k C -1 ( ν + num -ν + Chollet ) 2 .

Variational multiscale approach (VMS)

The method is based on the variational scale separation concept [START_REF] Collis | Monitoring unresolved scales in multiscale turbulence modeling[END_REF][START_REF] Hughes | Large eddy simulation and the variational multiscale method[END_REF][START_REF] Hughes | Large eddy simualtion of turbulent channel flows by variational multiscale method[END_REF]. Initially the model was formulated in a weak variational form. We introduce the model as it was done by Collis [START_REF] Collis | Monitoring unresolved scales in multiscale turbulence modeling[END_REF]. Let us define the computational domain Q ≡ Ω×]0, T [ with the boundary P ≡ Γ×]0, T [, where the time interval is ]0, T [ and the spatial domain Ω is an open, connected, bounded subset of R d , d = 3. The state vector U ≡ {u, p} t is defined on the closure of the space-time domain Q and in the function space V. To construct the variational form of the equations we define a second function space W, represented as a set of the test functions W ≡ {w, r} .

The model description needs an inner-product definition:

(f , g) Q ≡ Q f • g dQ, (f , g) P ≡ P f • g dP.
Let us re-write the system of the NS equations with the incompressible condition (1.3)

SOME REMARKABLE LES MODELS

in a generalized form:

N S(U) ≡ ∂ t u + ∇ • (u ⊗ u) + ∇p - 1 Re ∆u ∇ • u = f ψ ≡ F.
Here f is the body force, ψ denotes the volume mass source.

For the weak (or variational) formulation of the NS equations:

B(W, U) ≡ (W, N S(U)) Q = (w, F) Q we define B(W, U) as B(W, U) = (w, ∂ t u) Q -(∇w, u ⊗ u) Q -(∇ • w, p) Q +(∇ s w, 2 1 Re ∇ s u) Q + (r, ∇ • u) Q +(w, n • (u ⊗ u)) P + (w, pn) P -(w, 2 1 Re ∇ s u • n) P ,
where (∇ s u) ij ≡ S ij ≡ (u i,j + u j,i )/2 is a symmetric part of the velocity gradient tensor (or strain rate tensor) and n is the unit outward-pointing normal vector on the boundary Γ. The integration by parts is applied to get the surface-dependent terms and change the spatial derivation of the function space V to the spatial derivation of the second functional space W.

The variational operator linearized about the field U is defined as:

B (W, U, U ) ≡ (w, ∂ t u ) Q -(∇w, u ⊗ u + u ⊗ u ) Q -(∇ • w, p ) Q + (∇ s w, 2 1 Re ∇ s u ) Q +(r, ∇ • u ) Q + (w, n • (u ⊗ u + u ⊗ u )) P +(w, p n) P -(w, 2 1 Re ∇ s u • n) P ,
where U is a linear perturbations and B (W, U, U ) is linear with respect to the first and the third arguments and affine with respect to the second argument.

The space of solutions is partitioned as:

V = V + Ṽ + V.
Let's separate the total velocity on the large Ū and small Ũ resolved scales and the unresolved scales Û: U = Ū + Ũ + Û and W = W + W + Ŵ.

The projection of the NS equations on the three subset ( V, Ṽ, V) of the function space

V is:

B( W, Ū + Ũ + Û) = ( W, F) Q B( W, Ū + Ũ + Û) = ( W, F) Q B( Ŵ, Ū + Ũ + Û) = ( Ŵ, F) Q
The Subgrid Stress:

R( w, u ) ≡ (∇ w, u ⊗ u ) Q -( w, n • (u ⊗ u )) P
and Cross Stress:

C( w, ū, u ) ≡ (∇ w, ū ⊗ u + u ⊗ ū) Q -( w, n • (ū ⊗ u + u ⊗ ū)) P .
are defined by analogy with the strong LES formulation. Large resolved scale projection is:

B( W, Ū) + B ( W, Ū, Ũ) -R( w, ũ) = ( W, F) Q -B ( W, Ū, Û) + R( w, û) + C( w, ũ, û),
where the second line terms are not represented for the resolved scale projection and have to be modeled. Small resolved scale projection is:

B ( W, Ū, Ũ) -R( w, ũ) = -[B( W, Ū) -( W, F) Q ] -B ( W, Ū, Û) + R( w, û) + C( w, ũ, û),
where the second line terms must also be modeled like in the previous case. Subgrid scale projection is:

B ( W, Ū + Ũ, Û) -R( ŵ, û) = -[B( Ŵ, Ū + Ũ) -( Ŵ, F) Q ]
Two assumptions were done for the VMS modeling:

• It seems reasonable to neglect the direct subgrid scale influence on the large resolved scales. The main energy transfer properties can be modeled by the large-small resolved scales interaction;

• The dissipation of the small resolved scales is conditioned by the subgrid scales and can be modeled in the traditional Smagorinsky eddy-viscosity closure.

These hypotheses lead to the VMS modeling equations:

B( W, Ūh ) + B ( W, Ūh , Ũh ) = R( w, ũh ) + ( W, F) Q B ( W, Ūh , Ũh ) -R( w, ũh ) = -(∇ s w, 2ν t ∇ s ũh ) Q -[B( W, Ūh ) -( W, F) Q ] 1.2.

SOME REMARKABLE LES MODELS

Combining the two equations for the large and small resolved scale projection and defining Ū ≡ Ū + Ũ, the system takes a compact form:

B( W, Ūh ) + (∇ s w, 2ν t ∇ s ũh ) Q = ( W, F) Q ,
where the Smagorinsky constant is calculated via the large νt = (c s ∆) 2 |∇ s ūh | or small resolved scales νt = (c s ∆) 2 |∇ s ũh |. These are large-small and small-small closures.

A resolvable subfilter scale (RSFS) model

The model appears to overcome the algebraic closure difficulties in highly underresolved homogeneous turbulence [START_REF] Zhou | A resolvable subfilter-scale model specific to large-eddy simulation of under-resolved turbulence[END_REF]. The model is looking for a separate dynamical system for the subfilter scale (SFS) terms. It approximates the subfilter scale dynamics from the NS equations. The authors separate the velocity into the resolved filtered part u r , subfilter resolved part u sf and subgrid scales u sg . The total velocity is the sum of all this components: u = u r + u sf + u sg . The model resolves approximately the exact evolution equations for u sf . The main LES equations are:

∂ t u r i = -∂ j (u r i u r j ) -∂ i (p r + p s ) + ν∆u r i .
The exact SGS stress has the form:

τ r ij = (u r i u s j + u s i u r j + u s i u s j ) r
In the resolved velocity governing equations u s is replaced by u sf . The modeled subfilter scale stress has the form:

τ r ij = β(u r i u sf j + u sf i u r j + u sf i u sf j ) r ∂ t u sf i = -∂ j τ r,r;sf ij -∂ j τ r,sf ;sf ij -∂ j τ sf,sf ;sf ij -∂ j (τ r,sg;sf ij + τ sf,sg;sf ij + τ sg,sg;sf ij ) -∂ i p sf + ν∆u sf i ,
where () sf implies a band-pass-like filter with a filter function dictated by the explicit and grid filters. The τ ij terms are composed of the corresponding velocity interactions (before the semicolon sign superscript) and the filter operation (after semicolon sign superscript):

τ r,r;sf ij = (u r i u r j ) sf , τ r,sf ;sf ij = (u r i u sf j + u sf i u r j ) sf , τ sf,sf ;sf ij = (u sf i u sf j ) sf , τ r,sg;sf ij = (u r i u sg j + u sg i u r j ) sf + (u sg i u sg j ) sf .
The terms representing the interactions between the resolved or subfilter scales with the subgrid scales σ ij ≡ (τ r,sg;sf ij + τ sf,sg;sf ij + τ sg,sg;sf ij ), must be modeled. The main physical property that should be satisfied is the energy transfer between the subfilter and the resolved velocities. So, the modeling terms should reproduce the action of the nonlocal subgrid-interaction terms σ ij and the non-resolved scale approximation by the subfilter scale field. The standard Smagorinsky closure is applied for modeling:

σ ij - 1 3 tr(σ ij )I = ν t S sf,sf . (1.43)
The results show the importance of all the nonlinear interactions in the resolved velocity governing equations. When neglecting the local (τ r,r;sf ij ) or nonlocal (τ r,sf ;sf ij ) terms, poor features of the subgrid scale induced energy transfer τ ij S r ij were obtained. The parameters proposed by Zhou et al. [START_REF] Zhou | A resolvable subfilter-scale model specific to large-eddy simulation of under-resolved turbulence[END_REF] are c s = 0.16 and β = 1. The model leads to better results in the prediction of the subfilter pressure force which is under-predicted with eddy-viscosity closure models.

Comparison of the models

Energy transfer.

The main feature of the SGS modeling is to extract the energy from the resolved scales closed to the wavenumber cutoff. Many models do this through the additional subgrid scale dissipation terms. This could be a good phenomenological approach of the energy transfer due to the subgrid scales. On the other side, the backward energy transfer is impossible if the additional viscosity is strictly non-negative.

The Smagorinsky closure is the simplest dissipative models. In isotropic turbulence it provides good spectra and third-order statistics. The modeling constant in this case can be estimated as c s = 0.18 which is approved by the simulations. The model provides an integral energy dissipation rather than a dissipation of the energy local fluctuations.

For near-wall modeling, the Smagorinsky model gives an excessive dissipation of the resolved structures. It can be due to the zero-level of the backward energy transfer modeling. But there is no guarantee that the intensity of the forward scatter is correctly produced for the different cases of turbulent flow. For channel flow turbulence the modeling constant was estimated as c s = 0.074 [START_REF] Ramakrishnan | Partition selection in multiscale turbulence modeling[END_REF] or c s = 0.1 [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF].

The application of the Van Driest damping with the Smagorinsky model does not produce the backward energy transfer. But it is used to decrease the dissipative effect of the closure near the boundaries. A good asymptotic comportment of the near-wall velocity RMS can be achieved if the dissipation level of the model is smaller than the real one. Nevertheless, the damping function imposes globally the modeling results, which can be good for different applications. The drawback of this damping function is the need for its ad hoc optimization for different flows. There is no guarantee that the optimal value of the damping parameters for the channel flow can be applied to arbitrary different boundary conditions. The damping reduces the near-wall diffusion effect of the Smagorinsky closure, so the modeling constant should be larger than in an undamped case.

The Dynamic model brings the possibility to produce the backward energy transfer. Numerical experiments show that the model gives a better near-wall asymptotic behavior of the square velocity statistics for the plane channel flow than the Smagorinsky model.

On the other hand, the Dynamic model is not numerically stable. The stability is directly dependent on the filtering procedure, the boundaries, the numerical schemes, etc. To achieve a reasonable stability level, a regularizing procedure has to be implemented: the dynamical constant is clipped to zero if the Germano procedure prescribes a negative value. So, the backscatter is forbidden for the Dynamic model. Nevertheless, better wall statistics are produced without ad hoc optimization in comparison to the Smagorinsky closure.

The similarity model provides two directions of energy transfer. The applications of the model show an unstable behavior of the computation. The SGS stress generates excessive backward energy transfer, which leads to excessive turbulence intensity. The mixed model is more stable because of the extra diffusion terms.

A further development of the Smagorinsky constant optimization by the Lagrangian scale dependent model, also provides the possibility of backscatter. The model takes into account the time-dependent properties of the subgrid scale interactions. During the authorized time of interaction between the subgrid and resolved scales one can get an increase of the local fluctuation of the "constant" C. Negative values are allowed and do not make the simulation diverge.

The Deconvolution model also provides a backward energy transfer. This procedure generates a synthetic field to produce the SGS stress and is based on no phenomenological assumption. The problem of the deconvolution operation is its unstable behavior leading to a regularizing procedure. This means that the deconvolution itself is insufficient to extract the necessary level of resolved scale energy from the system. For the plane channel flow case, the main instability was detected proportional to the largest subfilter wavenumber [START_REF] Stolz | An approximate deconvolution model for large-eddy simulation with application to incompressible wall-bounded flows[END_REF]. The integral value of the forward energy transfer is provided with the help of the regularizing friction term. The resolved-subfilter scale transfer is modeled through the standard nonlinear velocity interaction. The existence of the subfilter spectral gap permits to isolate the direct resolved-subgrid scale interaction. The regularizing procedure acts mostly on the subfilter energy [START_REF] Stolz | An approximate deconvolution model for large-eddy simulation with application to incompressible wall-bounded flows[END_REF] and its resolved scale action depends on the filter properties. In the ideal case of a completely invertible filter, the regularizing term is zero ((I -G 2 ) * ūi = (I -Q N * G) * ūi = (I -I) * ūi = 0). In real simulations, the regularizing spectral interval width depends on the order of the deconvolution approximation N (ex.: N = 5 in [START_REF] Stolz | An approximate deconvolution model for large-eddy simulation with application to incompressible wall-bounded flows[END_REF]).

From the energy transfer point of view, the Subgrid estimation models are interesting by their time-interaction correction. In numerical experiments with the high Reynolds number, it was shown by Domaradzki et al. [START_REF] Domaradzki | The subgrid-scale estimation model for high Reynolds number turbulence[END_REF] that the non-correct coupling between the resolved and modeled scales leads to an incorrect energy spectra. The parallel coarse DNS simulation was shown to produce the good locally estimated energy transfer. The procedure does not need an additional friction or diffusion term. The backward energy transfer is also possible in this model through the nonlinear velocity interaction term in the coarse DNS simulation.

The Implicit adaptive local deconvolution method uses the effective energy dissipation procedure to compute the numerical scheme parameters by comparing the effective modeled viscosity action with the known turbulent viscosity closure. So, a part of the subgrid scale force is a priori estimated by the known SGS closure. It guarantees that the integral value of the modeled energy transfer should not be far from the known SGS energy transfer function. This closure could provide a backscatter if the reference model is convenient. The parameters of the weighting functions, determined for isotropic flow at infinitely large Reynolds number, give improved characteristics of the transitional flows, Taylor-Green vortex, plane channel flow at different Reynolds numbers in comparison to the Dynamic model. So, there is a certain insensitivity of the modeling parameters to the Reynolds number. In the case of a bounded flow, an excessive dissipation was found. To improve the mean shear stress value, a wall damping function was used. IALDM shows stable turbulent statistics even for extremely coarse simulation (12 3 for Re τ = 590). This allows to estimate a good forward energy drain without "energy containing" modes at the resolved scale range. The parametrization of damping functions is also insensitive to the Reynolds number variation. The negative value of the numerical flux function produces inverse energy transfer. Example of this case was shown for the isotropic turbulence at infinite Reynolds number [START_REF] Hickel | On implicit subgrid-scale modeling in wall-bounded flows[END_REF].

The VMS model demonstrates a good quality of the velocity square statistics. It uses Smagorinsky approach to close the small scale resolved velocity equations. The energy transfer of the resolved scales interval is produced by the nonlinear term of the resolved and subfilter scales. Only resolved-subgrid scale interactions are neglected. It is sufficient to produce a good small scale drain of resolved scale energy for a sufficiently wide resolved and subfilter scale interval, which size depends on the Reynolds number and boundary conditions. The subfilter-subgrid scale energy transfer is modeled by the Smagorinsky closure. The three-level scale separation allows to separate the direct interaction between the resolved and subgrid scales. So, the SGS motion cannot directly affects the energetically important scales. This coincides with a Kraichnan argument that the near cutoff local triadic interactions are more important to produce a good dissipation level at the characteristic filter (or grid) size. Finally, a sufficiently large subfilter scale range provides a correct level of energy drain of the resolved scales. The choice of large-small or smallsmall closure based on ∇ s ūh or ∇ s ũh influences the modeling constant estimation. For the isotropic turbulent flow the optimal value of the constant was calculated as [START_REF] Hughes | The multiscale formulation of large eddy simulation: Decay of homogeneous isotropic turbulence[END_REF]:

c s /c s =
1.28 for smallsmall 1.62 for largesmall, where the modeling constants c s and c s correspond to the resolved scale filter size ∆ and subfilter scale filter size ∆ . The channel flow could be simulated with the standard c s = 0.1 value of the Smagorinsky constant [START_REF] Holmen | Sensitivity of the scale partition for variational multiscale large-eddy simulation of channel flow[END_REF][START_REF] Hughes | Large eddy simualtion of turbulent channel flows by variational multiscale method[END_REF].

A Resolvable subfilter scale model such as the VMS approach neglects the direct energy transfer between the resolved and subgrid scales. The resolved-subfilter energy transfer is produced to the exact NS nonlinear terms. The subgrid scale interaction is modeled by eddy-viscosity, so for backward energy transfer from subgrid to subfilter scales the local viscosity value should be negative. The RSFS model, as the VMS approach, uses three levels of spectral separation of the velocity. The errors in the subgrid scale approximation does not harm the integral scales directly. The interactions between the resolved and subfilter scales are given by the nonlinear velocity terms, so the corresponding energy transfer (forward and backward) is provided exactly as in NS equations.

Channel flow.

The VMS model (described in section 1.2.11) was tested for the turbulent channel flow at the different Reynolds numbers (Re τ = 180, Re τ = 395, Re τ = 590 and Re τ = 720 [START_REF] Holmen | Sensitivity of the scale partition for variational multiscale large-eddy simulation of channel flow[END_REF][START_REF] Hughes | Large eddy simualtion of turbulent channel flows by variational multiscale method[END_REF][START_REF] Ramakrishnan | Partition selection in multiscale turbulence modeling[END_REF]). The mean flow and RMS statistics were correctly predicted by the model. The sensitivity to the scale separation was also tested. Ramakrishnan and Collis [START_REF] Ramakrishnan | Partition selection in multiscale turbulence modeling[END_REF] noted that the minimal near-wall resolution in a streamwise and a spanwise directions must be λ +

x ≈ 200 and λ + z ≈ 50 to get an improved mean velocity field and fluctuating velocity statistics. This resolutions were defined from the analyze of the near-wall velocity correlation (Re τ = 180, y + = 14.5). Consequently, the VMS model is more sensitive to the streamwise scale separation and the preferential filter direction is spanwise. It was remarked that both the filter and mesh cutoff should belong to the inertial range of the turbulent velocity spectrum.

The comparison to the Smagorinsky and Dynamic model was done by Ramakrishnan and Collis in [START_REF] Ramakrishnan | Partition selection in multiscale turbulence modeling[END_REF] for Re τ = 180. The Smagorinsky constant was chosen c s = 0.074 in the sequence of the work of Hughes et al. [START_REF] Hughes | The multiscale formulation of large eddy simulation: Decay of homogeneous isotropic turbulence[END_REF]. The mean velocity and fluctuating velocity RMS were better produced by the VMS and Dynamic models. The Smagorinsky model overpredicted (as it was shown in [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF][START_REF] Piomelli | Model consistency in large eddy simulation of turbulent channel flow[END_REF]) the mean flow in the logarithmic inertial layer and the center of the channel. Also, the Smagorinsky model shows an excessive resolved scale energy in comparison to the complete DNS statistics for y + > 15. The reference DNS of Moser et al. [START_REF] Moser | Direct numerical simulation of turbulent channel flow up to Re τ = 590[END_REF] was not truncated to compare to the LES. So, the excessive velocity RMS of Smagorisky model should be larger compared to the truncated DNS than compared to the total one, as in [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF]. The velocity spectra at y + ≈ 10 confirm an excessive dissipation of the resolved energy produced by the Smagorinsky model. The VMS and Dynamic models have almost the same spectra in the resolved scale interval. Behind the scale separation wavenumber the velocity spectra of Dynamic model are slightly larger than the VMS ones which are better compared to DNS.

The SGEM (described in section 1.2.9) applied to the turbulent channel flow was discussed by Domaradzki et al. in [START_REF] Domaradzki | Large eddy simulations using the subgridscale estimation model and truncated Navier-Stokes dynamics[END_REF]. The authors have compared the modeled turbulent flow at Re τ ≈ 200 and Re τ ≈ 1000 to the DNS data, Dynamic model data of Piomelli [START_REF] Piomelli | High Reynolds number calculations using the dynamic subgrid-scale stress model[END_REF] and an experimental results of Wei et al. [START_REF] Wei | Reynolds-number effects on the structure of a turbulent channel flow[END_REF]. The SGEM with a parallel truncated NS run (Eq. (1.41)) produces correct mean velocity profiles of turbulent channel flows for both Reynolds numbers. The wall distances of the maximum of the Reynolds stress components are also well predicted. For Re τ ≈ 200 the Reynolds stress maximum are overpredicted for all 3 velocity components compared to the filtered DNS results. In the center of the channel, the Reynolds stress have the improved intensities except the streamwise fluctuating velocity component which exhibits too low RMS values. The maximum of the SGS dissipation is overpredicted for this Reynolds number. The model with a simplified form of the SGS stress tensor estimation (eq. 1.42) was compared to an unfiltered DNS statistics. The mean profile was well produced with a slightly higher value for y + > 30. The maxima of the streamwise and the wall normal fluctuating velocity RMS are also in a good agreement with the unfiltered DNS results. For y + > 60 these components have smaller RMS values as well as all-range RMS of the spanwise direction fluctuating velocity. The best RMS value for the simplified SGS tensor estimation is possibly due to the finer grid. It seems important to change the wall-normal direction resolution from N y = 49 to N y = 65 to get a sufficient resolution of the near-wall structures. The comparison of the Reynolds stress to the unfiltered DNS data is approved because the filtered modes do not significantly change the fluctuating velocity RMS.

The statistics of the turbulent flow at Re τ ≈ 1000 were compared to the experimental data. Both versions of the SGEM produce correctly the mean velocity flow. The Reynolds stress was compared only for the streamwise and spanwise components. The streamwise RMS velocity has a too large maximum. The vertical resolution is extremely coarse (N y = 49), which can provide a spurious increase of the resolved kinetic energy. Unfortunately, the spectral composition of the velocity as the function of the wall distance is not analyzed in the article to verify this hypothesis. Both simulations underpredict the turbulence energy for y/δ > 0.2 (δ is half of the channel length). This is probably due to the top-hat filter which is dissipative, e.g. reduces the large scale fluctuations. This filter is applied at each time step to compute the SGS stress tensor according to its definition from the estimated velocity (Eq. (1.40)). The filtering of the total experimental data to the LES resolution should reduce the intensity of the turbulent kinetic energy and change the ratio between the Reynolds stress of the modeled and measured flows.

The mixed model (described in section 1.2.5) compared to the Smagorinsky model, DNS and experiment was presented in [START_REF] Piomelli | Model consistency in large eddy simulation of turbulent channel flow[END_REF]. The Gaussian and cutoff explicit filters and resolutions were tested for both models applied for the turbulent channel flow at Re τ = 180. The refinement of the grid allows an improvement of the mean velocity profile prediction for the mixed model at y + > 15. The Smagorinsky model is less accurate at the same resolutions and does not compare to the reference DNS results. The fluctuating velocity RMS of the Smagorinsky model attain the reference curves with the refinement of the spatial resolution. Nevertheless, in the best case, the model overpredicts the wall normal and the spanwise RMS velocity at y + > 100 and y + > 50, which is in agreement with the results of Ramakrishnan et al. [START_REF] Ramakrishnan | Partition selection in multiscale turbulence modeling[END_REF]. The mixed model produces the best streamwise RMS velocity at the finest resolution but it is not the case of the vertical and the spanwise RMS velocity.

The velocity spectra near the wall y + ≈ 13 show the excessive dissipation of the resolved scales by the Smagorinsky model. The velocity spectra modeled by the mixed model are better compared to the reference ones.

The mixed model with the Gaussian filter is much better than the Smagorinsky model with the cutoff explicit filter in the a priori prediction of the SGS Reynolds stress τ 12 .

The conclusion for the properties of the mixed model is that the SGS stress tensor has a good a priori properties and the model predicts accurately the energy spectra in the near-wall region as well as the mean velocity. The model authorizes the forward and backward energy transfer which is particularly important in the near-wall region, but it is not sufficiently dissipative in the upper regions y + > 50. Also the model depends strongly on the spatial resolution and all the RMS velocity statistics can not be improved by a grid refinement.

The Lagrangian scale dependent and independent models (described in section 1.2.7) were applied to the boundary layer flow over the homogeneous and heteroge-

COMPARISON OF THE MODELS

neous rough surfaces [START_REF] Bou-Zeid | A scale-dependent Lagrangian dynamic model for large eddy simulation of complex turbulent flows[END_REF]. The Authors mentioned that the model produces the correct level of the dissipation and backward subgrid-resolved scales energy transfer in comparison to the Smagorinsky model. The advantages of the model is that it does not need an homogeneous direction. This allows a possibility to apply the model to the flow over a rough wall.

The ADM approach (described in section 1.2.8) was applied to the turbulent channel flow by Stolz et al. [START_REF] Stolz | An approximate deconvolution model for large-eddy simulation with application to incompressible wall-bounded flows[END_REF]. Two different Reynolds numbers (Re τ = 180 and Re τ = 590) and different space resolutions for both of them were applied. The shear stress was exactly produced in the simulations with the finest mesh. The mean velocity profiles match well the reference profiles of a filtered DNS. The wall distances of the maximum of the RMS velocity are well predicted but their values are slightly larger compared to the filtered DNS. The mean spectra of velocity are rather in a good agreement too. The model get a good level of the dissipation of the smallest resolved scales which produce an excessive Reynolds stress in the coarse DNS simulation, as it is shown by Stolz et al. [START_REF] Stolz | An approximate deconvolution model for large-eddy simulation with application to incompressible wall-bounded flows[END_REF]. The shear stress produced by the coarse DNS is too large. The ADM approach reduces an "extra" of the shear stress.

The dynamical iterative procedure of χ u computation provides a correct level of the mean turbulent energy.

The ALDM (described in section 1.2.10) of the channel flow was studied by Hickel et al. [START_REF] Hickel | On implicit subgrid-scale modeling in wall-bounded flows[END_REF]. Four Reynolds number were used: Re τ = 180, Re τ = 395, Re τ = 590 and Re τ = 950. To improve the resolved-subgrid scale energy transfer near the wall, the van Driest damping was used. This allows a change of the near-wall mesh anisotropy and reduces the dissipative weight of the numerical flux function [START_REF] Hickel | On implicit subgrid-scale modeling in wall-bounded flows[END_REF]. The mean velocity profile was very well reproduced for all Reynolds numbers and all spatial resolutions. The slightly underprediction of the streamwise velocity RMS was remarked for Re τ = 180 and Re τ = 950. The energy budget as a function of the spatial resolution was tested for Re τ = 590. All Reynolds stress components were shown strongly dependent on the resolution interval: ∆ +

x ∈ [51.5; 308.9], ∆ + z ∈ [25.7; 154.5] and ∆ + y ∈ [41.5; 227.2] (in the center of the channel). The refinement of the spatial resolution leads to an increase of the velocity RMS.

The comparison to the Dynamic model was done for Re τ = 395. The model overestimates the streamwise velocity RMS and underestimates the wall normal and spanwise ones. The decrease of the subgrid scale dissipation of ALDM by the Van Driest damping allows the better Reynolds stress compared to DNS reference.

The RSFS approach was applied by A. Juneja et al. [START_REF] Juneja | Characteristics of subgrid-resolved-scale dynamics in anisotropic turbulence, with application to rough-wall boundary layers[END_REF] and Y. Zhou et al. [START_REF] Zhou | A resolvable subfilter-scale model specific to large-eddy simulation of under-resolved turbulence[END_REF] to the boundary-layer turbulence over a rough wall. The model takes into account for the resolved and subfilter acceleration and resolved-subfilter energy flux. The model have shown better resolved velocity, subgrid pressure force and alignment of the acceleration vector compared to the Smagorinsky and the similarity models. Specifically, the RSFS approach improves the flow properties when the computational grid does not contain all the energetically important turbulent scales.

Chapter 2 LRDT method

Motivation for a new model

Criticism of existing models

The analysis of the previously mentioned models suggests that their efficiency depends on the complexity of the flow. For example, the Smagorinsky model [START_REF] Smagorinsky | General circulation experiments with the primitive equations[END_REF] may be accurate enough for isotropic turbulence but needs further refinements for wall bounded flows: a wall damping function or the Dynamic Germano procedure [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF] to damp the eddyviscosity near the wall. The Smagorinsky constant also is not universal and should be adapted for an optimal modeling of each flow [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF]. Another example is the refinements of the Subgrid Estimation model of Domaradzki et al. for high Reynolds number turbulence where the simplified procedure of the subgrid scale estimation is improved by the parallel evaluation of a coarse DNS simulation [START_REF] Domaradzki | Large eddy simulations using the subgridscale estimation model and truncated Navier-Stokes dynamics[END_REF]. In most cases, the tuning of the models for a new turbulent flow is done a posteriori and not before the simulation. A model which allows an a priori estimation of its parameters would be preferential.

Strategy of LES-Langevin modeling

The LES-Langevin model considers the eddy-viscosity model as the zero-level approach of the subgrid scales modeling. These models can be classified as "deterministic" in comparison to the "stochastic" models based on synthetic fields.

In eddy-viscosity methods, the action of the small scales is parameterized via a few deterministic numbers, linked to the various components of the subgrid-scale stress tensor. These models seek to reproduce the intensification of energy transport due to the action of scales widely separated from the considered one. However, they fail to reproduce backward energy transfer (backscatter) from small to large scale, created by elongated triads in the spectral space [START_REF] Ossia | Energy backscatter in large-eddy simulations of threedimensional incompressible isotropic turbulence[END_REF]. This effect has been shown to induce a stochastic behavior in LES [START_REF] Leith | Stochastic backscatter in a subgrid-scale model -plane shear mixing layer[END_REF]. In the ideal case of isotropic turbulence, this backscatter is usually viewed as secondary, and eddy-viscosity based models are generally satisfactory. However, in more realistic situations, including near wall turbulence, it has been shown that the energy backscatter is an essential feature [START_REF] Clark | Evaluation of subgrid-scale models using an accurately simulated turbulent flow[END_REF][START_REF] Meneveau | Statistics of turbulence subgrid-scale stresses -necessary conditions and experimental results[END_REF]. Some models like mixed models based on similarity hypothesis ( [START_REF] Horiuti | A new dynamic two-parameter mixed model for large-eddy simulation[END_REF][START_REF] Liu | On the properties of similarity subgrid-scale models as deduced from measurements in a turbulent jet[END_REF]) are able to reproduce a realistic backscatter in some situations. The need for backscatter modeling also leads to the development of "stochastic" strategies, where the discarded small-scale motions are random velocity field, which mimic either a random force or synthetic velocity fields [START_REF] Domaradzki | A subgrid-scale model based on the estimation of unresolved scale turbulence[END_REF][START_REF] Dubrulle | A model for rapid stochastic distorsions of small-scale turbulence[END_REF][START_REF] Scotti | A fractal model for large eddy simulation of turbulent flow[END_REF]. In many ways, this strategy is similar to the one used to describe the dynamics of a heavy particle coupled to a thermal bath involving many degrees of freedom (the so-called Brownian motion). The decimation is here performed by substituting, in place of the bath, a deterministic friction and a stochastic force, the two terms being linked through the dissipation theorem. The initial problem is then completely described through the so-called Langevin equation.

Turbulence is typically an out-of-equilibrium system, and there is probably no hope that such a simple description will ever be possible.

However, we would like to use this analogy to motivate a new strategy for LES modeling: replace the actual dynamics of the decimated degrees of freedom by a suitable noise, via a Langevin equation. Although this strategy may seem close to recent models based on synthetic fields, we would like to point out an important philosophical difference: rather than trying to estimate the actual small-scale dynamics, we aim at trying to estimate a plausible small-scale dynamics. We believe there is no unique solution for this last option. In the sequel, we present one solution based upon Rapid Distortion Theory [START_REF] Dubrulle | A model for rapid stochastic distorsions of small-scale turbulence[END_REF]. There may exist in fact more efficient models, based for example on information theory [START_REF] Jaynes | Information theory and statistical mechanics[END_REF].

To be more specific, consider a turbulent flow, with velocity field u i (x, t) and introduce a filtering procedure so as to separate it into a resolved field u i and a subfilter field u i = u iūi . The resolved field obeys a dynamical equation obtained by filtering of the Navier-Stokes equations, which may conveniently be written as [START_REF] Wu | Turbulent force as a diffusive field with vortical sources[END_REF]:

∂ t ū + (ū • ∇)ū + l + (u • ∇)u = -∇p + ν∆ū, l = (ū • ∇)u + (u • ∇)ū. (2.1)
Here, p is the resolved pressure, ν is the viscosity and l is a turbulent force. In ideal situations, including a spectral gap between resolved and subfilter scale, the vector l is zero, and one can rigorously show that the contribution of the Subgrid Reynolds Stress term (∇ • R = (u • ∇)u ) is of "diffusive" type (providing certain symmetries which exclude first order behavior such as the AKA effect [START_REF] Dubrulle | Low-viscosity lattice gases[END_REF]). Even as one departs from this idealistic situation, experimental [START_REF] Tong | Resolvable-and subgridscale measurement in the atmospheric surface layer: Technique and issues[END_REF] and numerical [START_REF] Laval | Nonlocality and intermittency in threedimensional turbulence[END_REF] studies show that this term correlates strongly with the resolved velocity gradient, thereby allowing a deterministic treatment through an eddy-viscosity of appropriate shape. In the same time, the force vector l becomes increasingly non-negligible (it can even become dominant in 2D situations see [START_REF] Dubrulle | A model for rapid stochastic distorsions of small-scale turbulence[END_REF]). It is responsible for a backscatter type of behavior and needs to be modeled through novel "non-diffusive" and "non-deterministic" strategies. In the sequel, we will therefore focus onto the modeling of the l term, via a generalized Langevin equation

∂ t l = Al + ξ,
where A is a generalized evolution operator, and ξ is a noise. However, a clear difficulty associated with this strategy is the lack of theoretical guide (equivalent of the statistical mechanics in Brownian motion) to help us devising the "best" friction, and the "best" stochastic force. For the time being, we then choose to derive our model as close as possible to the real dynamics of Navier-Stokes by trying to derive it from the original dynamical equations, rather than from pure empirical or dimensional considerations. For this, we reformulate the RDT-Langevin model of Laval et al. [START_REF] Dubrulle | A model for rapid stochastic distorsions of small-scale turbulence[END_REF] in a way suitable for LES. There are of course limitations to this approach, pertaining the need for both a simple enough model, and for tractable analytical computations. We try to formulate them as honestly as possible by pointing out the approximation we make at the various stages of the derivation of the model.

2.2

The LES-Langevin model of turbulence

Derivation of the Langevin equation

Our derivation is based on the stochastic RDT model developed by Laval et al. [START_REF] Dubrulle | A model for rapid stochastic distorsions of small-scale turbulence[END_REF][START_REF] Laval | Langevin models of turbulence: Renormalization group, distant interaction algorithms or rapid distorsion theory?[END_REF][START_REF] Laval | Fast numerical simulations of 2d turbulence using a dynamic model for subfilter motions[END_REF]. This model is based on the observation that subfilter scales are mostly linked to the resolved scales via linear processes akin to rapid distortion. This property is substantiated by various numerical simulations, and is linked with the prominence of non-local interactions at subfilter scales [START_REF] Dubrulle | A model for rapid stochastic distorsions of small-scale turbulence[END_REF]. Using incompressibility, the small-scale dynamics in this model can be written as:

∂ t u = -l -∇p + ∇(ν + ν t )∇u -f . (2.2)
Here, p is the subfilter pressure, ν t is a turbulent viscosity describing the non-linear interactions between subfilter scales and f is a forcing stemming from the energy cascade. The latter can be shown to be dominated by resolved scales non-linearities through f i = ∂ j (ū i ūjūi ūj ). Finally, we may use the observation that subfilter scales vary over short time scale compared to resolved scales to write:

∂ t l ≈ (ū • ∇)∂ t u + (∂ t u • ∇)ū, ≈ -{(ū • ∇)(l + f )) + [(l + f ) • ∇] ū} -{(ū • ∇)∇p + [(∇p ) • ∇] ū} +visc, (2.3)
where visc gathers all the term containing ν or ν t . This equation is only an approximation, in so far as the assumption of "rapidly varying scales" becomes less and less valid as the resolved scales and subfilter velocities become closer and closer in scale space. However, it seems to capture the dominant physics of the evolution of the vector l, as will be shown later. Further, it is tempting to simplify the viscous terms of the second equation of (2.3) to try and get a closed equation for l. Indeed, these terms involve an a priori rather arbitrary turbulent viscosity and one could redefine it so that the viscous terms are simply lumped into a term ∇(ν + ν t )∇l. Finally, we note that the terms involving the subfilter pressure depend on boundary conditions and on subfilter velocities, so that they can be thought to vary over a short time scale, contrary to f , which varies over a long time scale. We therefore choose to collect all the terms involving the pressure and f into a noise term, ξ 0 + ξ, such that ξ is a Gaussian-centered noise and :

ξ 0 = -(ū • ∇)f -(f • ∇)ū, < ξ > = 0, < ξ i (t, x)ξ j (t , x ) >= T ij (t -t , x -x ), (2.4)
where T ij is the noise correlation function, to be specified. Note that since ξ comes from a pressure contribution, it only affects the non-solenoidal part of the turbulent force.

Sequentially, if we study a turbulence with periodic boundary conditions, we can apply a simple projection procedure to work only with solenoidal fields, thereby discarding the ξ term. For the boundary flow turbulence, the spatial and temporal forms of T ij were found by the a priori study. Collecting all the results, we obtain the following RDT based-model for the turbulent force l as:

∂ t l = -(ū • ∇)l -(l • ∇)ū + ∇(ν + ν t )∇l + ξ 0 + ξ (2.5)
It takes the form of a generalized Langevin equation. The first and second terms of the r.h.s. are the advection and stretching of filtered l -vector by the resolved velocity.

The third term represents a diffusion of turbulent force l. The molecular viscosity ν represents the molecular friction effect and the additional turbulent viscosity ν t mimics the turbulent diffusion supplied by the viscosity-contained terms of (2.3). The force ξ 0 is generated through the energy cascade. The correlations T ij of stochastic force ξ are physically imposed by boundary conditions via pressure terms. Before implementing this Langevin equation into a LES model, we first validate it through a priori dynamical tests and seek optimal performances by tuning the parameters. The parameters estimation are different for isotropic and channel flow. Our LES-Langevin model will follow after this validation step.

A friction term

The validation tests on decaying and forced isotropic homogeneous turbulence are performed by Laval et al. [START_REF] Laval | A LES-Langevin model for turbulence[END_REF]. The Taylor Reynolds number is approximately constant and equal to R λ = 200 for the forced DNS and varies from R λ = 260 to R λ = 26 in the decaying case. Both DNS are performed with 341 3 effective wavenumbers after desaliasing (512 3 grid points). The simulation with LES-Langevin model are performed with k c = 21 (42 3 effective wavenumbers are used for the resolved scales). Complete description of simulations can be found in Laval et al. [START_REF] Laval | A LES-Langevin model for turbulence[END_REF].

The validation of the model was performed by comparing a full DNS and a simulation at the same resolution, in which the turbulent force l is replaced by the solution of Eq. (2.5). By this way, we could explore the validity of the approximations we made in its derivation. We focus here on the energy spectra, so as to capture possible deficiency of our model regarding energy transfer between scales. The test is performed for 1 < t < 2 (all the simulations are initialized with the same field at t=1). By that time, the kinetic energy is divided by a factor of approximately 2. The energy spectra at time t = 2 are shown in Fig. 2.1. The model (RDT) and the DNS agree at the largest scales (k < 3) but they significantly differ for smaller scales (the bump at k=21 is due to the coupling between the equation for resolved scales and the equation for subfilter scales running in parallel). This can be explained through the analysis of the evolution of l (Fig. 2.2). One sees that the RDT model leads to a constant increase of the smallest modes of l in time. After a given period of time, the contribution of these unrealistic resolved scales of l influence the model of the velocity field resolved scales. An explanation of this feature can be found following a study of the RDT model (2.2) by Laval et al. [START_REF] Laval | Fast numerical simulations of 2d turbulence using a dynamic model for subfilter motions[END_REF] showing that the process of small-scale stretching by random large scales is akin to a dynamo process, with exponential increase of the small-scale energy. A way to stabilize the system is to include a friction term in the RDT equation, leading to a stationary energy spectrum depending on the friction time τ f . A Kolmogorov k -5/3 spectrum is obtained for τ f = 27/22 Ω, where Ω = (S ij S ij ) -1/2 is a typical stretching rate based on a spatial average of the large-scale velocity stress tensor S ij . Using Eq. ( 2.3), one sees that such a friction term generates an equivalent friction term in the equation for the l -vector. This remark motivates the introduction of a stabilizing friction term -l/τ f in the equation for l to try to stabilize the coupled system. Indeed one observes a significant improvement with respect to the original RDT model. We therefore adopt this procedure as our starting RDT model, from which we now build our Langevin-LES model.

Derivation of the SGS model

The derivation of the LES-Langevin model of turbulence proceeds in two steps. In the first one, we replace the term (u • ∇)u by a turbulent viscous term µ t ∆U, acting only at large scales, in the spirit of standard deterministic "eddy-viscosity" models. In a second step, we derive a suitable Langevin equation of the turbulent force l through a decimation of the number of degrees of freedom corresponding to scales beyond the aliasing limit. For this, we introduce a strong hyperviscosity to damp all components beyond a given cut-off wavenumber k m . There is a priori complete freedom for the choice of k m . Here, we note that the cascade-driven forcing f has components only up to k = 2k c . Therefore, any component of l beyond 3k c will only be generated through secondary processes (stretching) rather than through the forcing. We may then hope that any k m between 2k c and 3k c provides the dominant contribution to the stochastic term l. It was found [START_REF] Laval | A LES-Langevin model for turbulence[END_REF] that k m = 2k c is in fact sufficient to capture this dominant contribution for isotropic turbulence. Our Langevin-LES (LRDT) model is therefore finally given by:

∂ t ū + (ū • ∇)ū + l = -∇p + ∇(ν + µ t )∇ū; [0 < k < k c ] ∇ • ū = 0; [0 < k < k c ] ∂ t l = -(l/τ f ) -(ū • ∇)l -(l • ∇)ū + ∇(ν + ν t )∇l + ξ 0 + ξ; [0 < k < k m ] (2.6) ξ 0 = -(ū • ∇)f -(f • ∇)ū, < ξ >= 0, < ξ i (t, x)ξ j (t , x ) >= T ij (t -t , x -x ),
where

τ f = 27/22 (S ij S ij ) 1/2 , (2.7)
ν t and µ t will be specified later and f has components f i = ∂ j (ū i ūjūi ūj ). Looking at Eq. (2.6), one recognizes a LES model where the backscatter coming from resolved subfilter scales interaction is parameterized through a noise. The latter obeys a generalized Langevin equation, with friction made of viscosity and rapid distortion by resolved scales, and with stochastic forcing ξ 0 + ξ, generated through the energy cascade and pressure processes.

Summary of hypothesis used

To clarify the used LES-Langevin assumptions we summarize the hypothesis:

• Hyp.I The resolved scale time variation can be neglected in comparison to the small scales variation:

∂ t u i ∂ t ūi .
• Hyp.II The viscous terms of the vector l can be modeled as a diffusion function ∇(ν + ν t )∇l.

• Hyp.III The SGS pressure terms can be modeled by the stochastic force ξ(x, t) with a prescribed spatio-temporal feature.

Isotropic turbulence case

In the isotropic turbulence simulation [START_REF] Laval | A LES-Langevin model for turbulence[END_REF], the model was parametrized as:

• The projection of the NS equation to the divergent-free velocity space eliminated the pressure from the equations. So, without any assumption, ξ = 0.

• The nonlocal terms (transport and stretching) were modeled by the generalized friction term:

l i /τ = l i /τ f + u j ∂ j l i + l j ∂ j u i ,
with the correction of the friction time:

τ = 27 22 < (S ij S ij ) -1/2 > +(S ij S ij ) -1/2 .
• The Reynolds Subgrid tensor was approximated by the spectral eddy-viscosity model

µ t = C(E(k c )/k c ) 1/2
, where the standard value of constant C = 0.267 (see Lesieur [START_REF] Lesieur | Turbulence in fluids[END_REF]) was reduced to C = 0.08.

• The small scale additional diffusion was neglected ν t = 0.

• A sufficient resolution of the turbulence force equation was found to be twice as fine as that of the resolved scale equation (k m = 2k c ).

Channel flow modeling

The turbulence properties in a plane channel are function of the distance to the wall. The Canonical Boundary Layer [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF] contains the inner and the outer regions. The inner regions corresponds to the interval 0 < y < 0.2δ, where δ is a boundary layer thickness. It is characterized by the viscous length l τ = ν/u τ . The inner layer can be divided in 3 sub-domains when the mean velocity profiles can be estimated [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF]:

• Viscous sublayer: 0 < y + < 5, where u + 1 (y + ) = y + ;

• Buffer layer: 5 < y + < 30, where u + 1 (y + ) ≈ 5 ln(y + ) -3.05;

• Prandtl or logarithmic inertial layer: 30 < y + ; y/δ 1, where

u + 1 (y + ) ≈ 1 κ ln(y + ) + 5.5 ± 0.1, κ = 0.41.
The outer region consists of logarithmic and wake layers. The characteristic length in this region is the thickness δ. The mean streamwise velocity dependence on the wall distance is different for the inner and outer log-layer parts.

In the viscous sublayer region the viscosity effect is predominant [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF] and the fluid flow is stable. The flow becomes instable in the buffer layer. This leads to an organisation of the streaks which propagate away from walls and interact with the mean flow. This complex processes is highly intermittent and impacts on the energy transfer between the turbulent velocity scales. In the viscous sublayer, the energy is contained only in the mean flow scales. In the buffer layer the energy production and dissipation start to increase and achieve a maxima at: y + ≈ 15. Beyond this point, the turbulence energy transfer declines. Far from the wall the turbulent energy transfer becomes the same as in isotropic turbulence.

The modeling of the initial velocity perturbations of the sublayer and their accumulation is essential for LES. The problem is that the coarse LES resolution does not allow a necessary number of freedom to compute exactly the near-wall structures. The possible modeling of the buffer-layer flow can be done by an additional stochastic forcing which corresponds to the complex nonlinear fluid interactions. The hypothesis about the behavior of the small (or subgrid) scale flow in the near wall region can changes significantly the mean flow properties. It was shown by Nazarenko et al. [START_REF] Nazarenko | Nonlinear RDT theory of nearwall turbulence[END_REF], in the frame of RDT approach, that the mean flow properties, can changes as a function of the initial small scale bursts statistics.

The LES-Langevin approach allows to include an information about the near wall flow in the model. The new wall model construction could be based on the theoretical and/or experimental analysis of the near wall small scale dynamics. The complex nonlinear interaction of the velocity and pressure is contained in the vector ξ. This term can not be computed exactly, because it depends on the small scale pressure fluctuations. This becomes clear if we project the momentum equations to the SGS set:

(∇ • u ) = 0, ∂ t u + [(u • ∇)u] = -∇p + 1 Re ∆u ,
where the prime () means the SGS part of a function. The application of the divergence operator to the SGS momentum equations and the use of the divergence free velocity condition supplies the SGS pressure equation:

∆p = -∇ • [(u • ∇)u] ,
where SGS pressure p depends on a total velocity field, a boundary condition and a filter operator. So, the analytical prediction of the vector ξ is equivalent to the prediction of the SGS velocity, which is a complicate mathematical problem. In our approach we use the a priori estimation of the vector ξ to model its spatio-temporal behavior.

The stochastic forcing ξ interacts with the resolved scale velocity through the l -vector. This turbulent force l produces the backward energy transfer from the subgrid scales to the resolved ones as well as a part of the forward energy transfer.

Estimation of stochastic forcing

The stochastic forcing

ξ(x) = (ū • ∇)∇p + [(∇p ) • ∇]ū (2.8)
is modeled via a noise correlation function T ij Eq. (2.4). The analytical derivation of T ij is difficult and needs a further hypothesis on the subgrid scale pressure. Our approach is based on an estimation of ξ(x) from the DNS data. We approximate ξ from Eq. (2.8) by a function ζ 0 , which depends separately of each space and time dimensions:

ζ 0i (t, x, y, z) = F ti (t)F xi (x)F yi (y)F zi (z), (2.9)
where index "i" means its components of the vector ζ. The functions F ti (t), F xi (x), F yi (y), and F zi (z) model the time, streamwise, wall normal and spanwise direction dependence respectively. F yi provides the intensity of the vector ξ, which is a function of the normal direction. The intensity of ξ does not change in the homogeneous x and z directions. The functions F xi and F zi supply the spectral form of the ξ-vector in these directions. Time function F ti provides a time-dependence of ξ. Without any loss of generality, we can assume:

F 2 xi x = F 2 zi z = F 2 ti t = 1.
According to the definition Eq. (2.9), the function F ij can be calculated from the a priori tests as:

F ij (x -x , t -t ) = ξ i (x, t)ξ j (x , t)
where the mean value of the r.h.s. is averaged in all 3 dimensions of space and in time.

The functions F yi are defined from (2.8) as:

F yi (y) ≡ ξ 2 i (x, t) 2 x,z,t ,
where • x,z,t stands for the averaging of the function in the streamwise, spanwise directions and time.

The functions F xi and F zi are defined from Eq. 2.8 and through the spectral form of the functions Fxi and Fzi :

| Fxi (k x ) |≡ | ξi (k x ) | 2 y,z,t | ξi (0) | 2 y,z,t 1/2 ; | Fzi (k z ) |≡ | ξi (k z ) | 2 x,y,t | ξi (0) | 2 x,y,t 1/2
, where • y,z,t and • x,y,t stands for the averaging in the wall-normal, spanwise directions and time.

At each time step, a new set of random phases are chosen for Fxi and Fzi :

arg( Fxi ) = α x ∈ (0, 2π), random, arg( Fzi ) = α z ∈ (0, 2π), random.
Time dependent function F ti (t) is approximated as a Gaussian random function with a mean value of 0:

F ti x,y,z,t = 0, F ti ∈ (-∞; ∞) σ F ti = 1.
(2.10)

If we suppose that the typical correlation time of ξ i is larger than the time step, a time correlation of ξ can be introduced using a Markovian series [START_REF] Castronovo | Subdiffusion and superdiffusion in Lagrangian stochastic models of oceanic transport[END_REF]:

   ξ i = ζ n+1 i ≡ C 1 ζ n t + C 2 ζ 0i , C 1 = 1 -dt/τ stf , C 2 = 2dt/τ stf ,
where τ stf is the expected time correlation. The time-decorrelation of the ξ-vector components corresponds to C 1 = 0, C 2 = 1.

A priori tests

DNS data

As it was noted, the approximated subgrid pressure for the particular simulation depends on the scale separation procedure. The vector ξ should reproduce the interaction of the resolved scales with all the non-resolved (subfilter and subgrid) scales. By definition it depends on the form and size of the filter. To understand the sensitivity of the forcing ξ to the variation of the separation scale size, and the filter form, the a priori tests were done for different filters. The importance of the highest subgrid scales was also tested by sampling the initial DNS pressure and velocity fields to a coarse mesh.

The sensitivity of the l -vector to the scale separation procedure was studied for the same filters and mesh samplings as for the vector ξ.

The database was generated through a DNS of channel flow at Re τ = 600 [START_REF] Marquillie | Direct numerical simulation of a separated channel flow with a smooth profile[END_REF]. The computational domain was 2π ×2×π with a resolution of 768×257×384 grid points. The a priori tests were done from 24 fields equally separated at time interval of ∆T /T 0 = 400, where characteristic time T 0 is based on the half height of the channel and the center channel velocity: T 0 ≡ h/u 0 = 1.

The tested filters are presented in Tab. 2.1. First three filters N1, N2 and N3 are applied to the complete 3D DNS velocity and pressure. They have spectral cutoff form in streamwise and spanwise directions and cutoff form in Chebyshev wall-normal direction. In turbulence at high Reynolds number, the intensity of the smallest subfilter scales is much smaller compared to the resolved scales. To understand the influence of the smallest subfilter scales on the vectors l and ξ we also reduced the subfilter scale range (N4, N5). In this case, the vectors l and ξ were calculated in two steps. Firstly, the 3D velocity and pressure fields were sampled to the coarse resolution by the spectral and Chebyshev cutoff filter. Secondly, the corresponding turbulent force l was calculated using the same cutoff form of scale separation procedure as in N1, N2 and N3 with new filter sizes. The discrete filtering in streamwise direction is represented by the test N6. This filter can be applied in the case of non-periodic boundary conditions, which is possible in our LES code. , k y m and k z m are the highest LES mesh-represented wavenumbers of Fourier (x-and z-directions) and Chebyshev (y-direction) polynomials; k x c , k y c and k z c are filter cutoff wavenumbers; the "discrete" filter is defined in Tab. 3.15, its spectral shape is shown in Fig. 3.5.

N k x m k y m k z m k x c k y c k z c

Spatial approximation of the stochastic forcing

The a priori profile and spectra of test N5 from Tab. 2.1 is presented in Fig. 2.3. The definition Eq. (2.8) allows a zero-value of the vector ξ at the wall. The function ξ 2

x,z,t exhibits a maximum and decreases at the center of the channel. The function ξ 2

x,z,t was approximated by a polynomial function P 2n , where the order of the polynomial must be even to respect the symmetry of ξ 2 with comparison to the wall-normal direction. An example of a ξ 2 polynomial fit is given in Fig. 2.3.

The spectra of Fxi , Fzi have an exponentially decreasing tail (Fig. 2.3). They were approximated by the following exponential functions:

| Fxi | 2 = exp -bx(kx-k x0 ) 2 , k x ∈ [0; k x max ] | Fzi | 2 = exp -bz(kz-k z0 ) 2 , k z ∈ [0; k z max ],
(2.11)

where k x0 and k z0 are the approximation parameters and k x max and k z max are the maximum wavenumber. These formula are not universal, but were derived to approximate the a priori spectra.

The filter and spatial resolution dependence of the stochastic forcing

Let us now estimate the susceptibility of the stochastic forcing to different LES grid resolutions and scale separations. The mean profiles of ξ 2 i corresponding to the case from Tab. 2.1 are presented in Fig. 2.4. The high near-wall calculated values of the vector ξ are due to the Chebyshev discretization near the wall. The first points in the wall-normal direction are very close one to another (at the first point ∆ + y ≈ 1), which overestimates the derivations of the velocity and pressure in this direction. We approximate the near-wall behavior of the vector ξ by a linear function which is zero at the wall.

The maxima of the mean ξ 2 are situated in the interval y + ∈ [60; 80], which is above the position of the turbulence maximum (y + ≈ 15). If all the subfilter scales are taken in account (tests N1, N2 and N3), the filter size changes the value of mean |ξ i | by a factor of 2. For the same filter size (like in tests N2 and N4), the reduction of the subfilter scales by a factor of 2 does not significantly decrease the value of the mean |ξ i |. So, the possible a priori estimation of the vector ξ for k c = (64, 49, 32) T can be calculated from the truncated velocity and pressure (test N4, Fig. 2.

3).

The a priori spectra of the vector ξ are presented in Fig. 2.5 (at y + ≈ 17) and in Fig. 2.6 (at the center of the channel). All spectra are discontinuous at the filter scale.

An example of such a discontinuity can be found for example in the a priori study of the planetary surface layer by Dubrulle et al. [START_REF] Dubrulle | A new dynamical subgrid model for planetary surface layer. Part I: The model and a priori tests[END_REF]. The k x spectra are advected by the finite differencing. The derivations of the velocity and pressure does not control the smallest scales (see [START_REF] Ghosal | The Basic equations for the large eddy simulation of turbulent flow in complex geometry[END_REF][START_REF] Hirsch | Numerical computation of internal and external flows[END_REF]). So, the a priori computation of the vector ξ provides a "noise" in the k x spectrum. This noise contributes to the nonlinear terms (ū∇)∇p + (∇p ∇)ū and produces the constant spectral behavior of the vector ξ at largest scales. The noise contribution decreases if the a priori computation is applied to the truncated velocity and pressure field (compare N2 and N4).

In most of cases, the highest modes of the vector ξ are closed to the cutoff scale. This approves the possibility to estimate the spectral form of the vector ξ from the test N4 of Tab. 2.1.

The use of the discrete filter (test N6 of Tab. 2.1) reduces significantly the subfilter scales and then the value of the vector ξ. We can see this in Fig. 2.5(a,c,e) and Fig. 2.6(a,c,e), where the vector ξ, obtained during the tests N5 and N6 has a different shape of the k x spectrum and the same form of k z spectrum.

--is N1, • • • • • • is N2, ----is N3, -• - is N4, -• • • -is N5, ---is N6. (a) (b) (c) (d) (e) (f)
Figure 2.5: A priori analysis of the stochastic forcing (vector ξ y ): mean spectra at

y + ≈ 17 of |ξ x (k x )| 2 (a), |ξ x (k z )| 2 (b), |ξ y (k x )| 2 (c), |ξ y (k z )| 2 (d), |ξ z (k x )| 2 (e) and |ξ z (k z )| 2 (f). Tests from Tab. 2.1: ---is N1, • • • • • • is N2, ----is N3, -• -is N4, -• • • -is N5, - --is N6.

Temporal dependence of the stochastic forcing

The modeled time correlation of the vector ξ was calculated from the a priori tests by

(k x )| 2 (a), |ξ x (k z )| 2 (b), |ξ y (k x )| 2 (c), |ξ y (k z )| 2 (d), |ξ z (k x )| 2 (e) and |ξ z (k z )| 2 (f). Tests from Tab. 2.1: ---is N1, • • • • • • is N2, ----is N3, -• -is N4, -• • • -is N5, ---is N6.
a classical definition:

f (t)f (0) / f (0) 2 ,
where f is an arbitrary function, the average • is performed with different realization of the initial condition f (0). An example of the correlation function is presented in Fig. 2.7. The transversal large scale velocity is compared to the small scale velocity and to the same component of the small scale pressure gradient. We see that the small scale pressure gradient correlation time is much shorter than that of the velocity. This small scale pressure gradient reduces the ξ vector time correlation. So, the delta-correlated noise should be a sufficient modeling approach for the vector ξ. If the time-dependence of the small scale pressure is taken into account, the correlation coefficient τ stf of the Markovian approach Eq. (2.11) should be less than ∆t ≈ 400T 0 .

Dynamics of turbulent force l , a priori tests.

The properties of the turbulent force l in presence of boundary conditions are different from those of isotropic turbulence [START_REF] Laval | A LES-Langevin model for turbulence[END_REF]. To analyse the properties in a new geometry and to understand the sensitivity to the filter shape, a priori computations were done, Tab. 2.1.

The mean profile of |l| obtained by a priori tests is presented in Fig. 2.8. The non-zero values of l near the wall, as in the case of vector ξ, are due to the Chebyshev discretization which increases the wall-normal direction derivatives near the boundary. The near wall behavior of |l| is not so perturbed as in the case of the vector ξ (see Fig. 2.4). This is due to the fact that the vector l depends only on the first derivatives of the subfilter scales and does not contain the second spatial derivatives of the subfilter pressure p . Also, the smallest subfilter scales (compare tests N2 and N4) are more significant for the vector l value in the center of the channel then in the wall distance of the maximum (y + ∈ [30; 60] for l x and y + ∈ [80; 110] for l y and l z ).

The spectra comparison of the a priori tests of the vector l are presented in Fig. 2.9 (at y + ≈ 17) and Fig. 2.10 (at the center of the channel). The spectra are discontinuous

--is N1, • • • • • • is N2, ----is N3 , -• - is N4, -• • • -is N5, ---is N6.
at the filter wavelength k c as it is the case in the a priori tests of vector ξ. The level of the noise provided by the streamwise finite differencing is smaller than in the a priori tests of vector ξ. We can see this by comparing Fig. 2.9(a, c, d), Fig. 2.10(a, c, d) and Fig. 2.5(a, c, d), Fig. 2.6(a, c, d). The largest modes are close to the cutoff wavenumber behind the cutoff. This confirms the argument of Kraichnan on the particular importance of the subgrid scale energy transfer by the scales close to the cutoff [START_REF] Kraichnan | Eddy viscosity in two and three dimensions[END_REF]. The backward energy transfer from the subgrid to resolved scales are modeled by the l term in LES equations 2.6. From the spectral behavior of the vector l , we can a priori estimate that | l| is much smaller than |l|. This is rather an undesirable property of the l vector, if its influence on the resolved velocity becomes too small.

The small value of the vector l calculated for the test N6 of Tab. 2.1 is explained by the streamwise spectra. We see that the discrete filtering strongly reduces the resolved and subfilter scale velocity, so the mean value of l becomes very small compared to case N5 where the spectral cutoff scale separation is used.

Figure 2.9: A priori analysis of the turbulent force (vector l ): mean spectra at

y + ≈ 17 of |l x (k x )| 2 (a), |l x (k z )| 2 (b), |l y (k x )| 2 (c), |l y (k z )| 2 (d), |l z (k x )| 2 (e) and |l z (k z )| 2 (f). Tests from Tab. 2.1: ---is N1, • • • • • • is N2, ----is N3 , -• -is N4, -• • • -is N5, - --is N6 .

Conclusions

The a priori analysis of the l and ξ vectors has shown that a maximum is achieved

(a) (b) (c) (d) (e) (f)
Figure 2.10: A priori analysis of the turbulent force (vector l ): mean spectra in the center of the channel of |l

x (k x )| 2 (a), |l x (k z )| 2 (b), |l y (k x )| 2 (c), |l y (k z )| 2 (d), |l z (k x )| 2 (e) and |l z (k z )| 2 (f). Tests from Tab. 2.1: ---is N1, • • • • • • is N2, ----is N3 , -• -is N4, -• • • -is N5, ---is N6 .
at larger wall distance (y + max ∈ [60; 110]) than the position of the turbulence kinetic energy maximum y + ≈ 15. This is surprising because the vector l describes most of the resolved-subgrid scale energy transfer which is maximal at the peak of turbulence kinetic energy [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF]. The correction can not be supplied by the vector ξ because its maximum is also a priori shifted toward larger y + .

The spectral maximum of the vector l is close to the filter size. So, the energy transfer is particularly important for the near cutoff modes. This confirms the statement of Kraichnan [START_REF] Kraichnan | Eddy viscosity in two and three dimensions[END_REF] that near cutoff turbulent kinetic energy transfer is preponderant for the whole resolved-subgrid energy transfer. Nevertheless, a priori tests shows a large decrease of the vector l in the small resolved scales range. This is crucial for LES-Langevin modeling of the resolved-subfilter energy transfer: the largest modes of the vector l are not represented in the resolved range and the resolved range spectrum has a very small intensity! The forcing ξ has the same spectral property and cannot qualitatively change the vector l spectrum, that is intensify the energy exchange between the resolved and subgrid scales.

Chapter 3

Numerical code

Geometry and equations

The numerical code for the LES-Langevin model is based on a 3D Navier-Stokes solver for 2D channel flows [START_REF] Marquillie | Simulation numérique et étude de la stabilité d'un écoulement de couche limite décollé[END_REF]. This solver can cope with a flat or curved lower boundary (allowing to compute converging diverging channels). The upper boundary is plane. The boundary conditions are no-slip on the upper and lower walls and periodic in spanwise direction. The streamwise boundary conditions are different for plane and convergingdiverging channel flows. They will be discussed later.

The spanwise discretization is done in Fourier space. The algorithm allows to decompose the 3D problem into series of 2D problems for each Fourier coefficient of the velocity and pressure. The coupling between the 2D solutions is mandatory to compute the nonlinear terms. In the wall normal direction the solution is discretized using Chebyshev polynomials. It is one of the well-known methods to allow no-slip boundary conditions. The streamwise discretization is done by explicit finite differences in the physical space.

The code is made fully parallel, using MPI. This allows to keep the computation time reasonable, even for large size problems.

Non-dimensional form of the LES-Langevin equations

Let us define the Reynolds number based on the mean streamwise inflow velocity u 1 at the center of channel and on the channel height L y :

Re = L y u 1 /ν. (3.1)
Dividing equation (2.6) by u 1 L y gives the Reynolds numbers corresponding to the eddy viscosities µ t and ν t :

Re µ = L y u 1 /µ t , Re ν = L y u 1 /ν t . (3.2)
The non-dimensional form of the LES-Langevin equations using the definitions (3.1)

         ∂ t ū + N u + l = -∇p + 1 Re ∆ū -∇.R, ∂ t l + N l = 1 Re ∆l + ξ 0 + ξ, ∇.ū = 0, (3.3)
where N u = (ū.∇)ū and N l = (ū.∇)l + (l .∇)ū are nonlinear terms, R is the Subgrid Stress tensor. Initially the viscosity ν t from Eq. (2.6) is neglected (ν t = 0). The Subgrid Reynolds Stress was chosen of the Smagorinsky form: R ij = 1 Re µ Sij .

Fourier representation in the spanwise direction

The Fourier discretization is possible because of the homogeneity of the flow in the spanwise direction. The grid points are defined by:

z j = 2L z j N z -L z , j = 0, . . . , N z -1,
where N z is a number of points and 2L z is the box size in the spanwise direction: z ∈ [-L z , L z ]. The direct and inverse Fourier representation of u(z) defined on [-L z , L z ] and its spectral coefficients

{û k }, k ∈ [-N z /2, N z /2 -1] are given by: ûk = 1 N z Nz-1 j=0 u(z j ) exp(-ıkβz j ), u(z j ) = Nz/2-1 -Nz/2 ûk exp(ıkβz j ),
where N z /2 is the number of Fourier modes and β = 2π/(2L z ) is a wave number.

The Fast Fourier Transformation (FFT) algorithm is defined for functions on the interval [0, 2π]. The definition interval of the function is mapped from [-L z , L z ] to [0, 2π]:

ζ j = π L (z j + L) = 2πj N z , j = 0, . . . , N z -1.
This change of coordinates introduces the coefficients dζ/dz and d 2 ζ/dz 2 in the derivatives.

The LES-Langevin equations have real solutions, so there are only N z /2 independent Fourier modes.

In a Fourier formulation the gradient and Laplacian operators are defined by:

∇ = (∂ x , ∂ y , -ıkβ), ∆ = ∂ 2 x + ∂ 2 y -k 2 β 2

Change of coordinates in wall normal direction

The geometry of the flow is not Cartesian because of the presence of a curved lower wall η(x). To avoid the resolution of the equation in the curvilinear geometry, the problem is reformulated in Cartesian coordinates. In order to do this, the initial curvilinear coordinates (t, x, y, z) are transformed in Cartesian ones (t, x, y, z):

t = t, x = x, y = y -η(x), z = z.
The mapping allows to represent the curvilinear grid discretization by a Cartesian one (see Fig. 3.2). The price to pay is extra curvilinear terms in equations (3.3).

The derivation operators (gradient, Laplacian) are changed to:

∇ = ∇ + G η , ∆ = ∆ + L η , ∇ × ( ∇ × u) = ∇ × ( ∇ × u) + Q η (u),
where the additional terms are

   G η = (-∂ x η∂ y , 0), L η = -∂ 2 xx η∂ y -2∂ x η∂ 2 xy + (∂ x η) 2 ∂ 2 yy , Q η (u) = (-∂ x η∂ 2 yy v, -∂ x η∂ 2 yy u -L η v).
The unit vector normal to the curved wall is:

n = (1 + (∂ x η) 2 ) -1/2 (n + n η ),
, L x ] × [0, L y ] × [-L z , L z ] are:          ∂ t ū + N Gu + l = -∇p -G η p + 1 Re ∆ū + 1 Re L η ū -(∇ + G η ).R GSM ∇.ū = -G η .ū, ∂ t l + N Gl = 1 Re ∆l + 1 Re L η l + N Gf u + ξ, (3.5)
where the nonlinear terms N Gu , N Gl , N Gf u , strain rate SG ij , and Smagorinsky stress R GSM are :

                         N Gu = (ū.∇)ū + (ū.G η )ū, SG ij = 1 2 (∂ i ūj + ∂ j ūi + G η j ūi + G η i ūj ), (∇ + G η ).R GSM = ((∇ + G η ).c 2 s (∆ x ∆ y (y)∆ z ) 2/3 | ( S + SG ) | ( S + SG )) N Gl = (ū.∇)l + (l .∇)ū + (ū.G η )l + (l .G η )ū, N Gf u = (f .∇)ū + (ū.∇)f + (f .G η )ū + (ū.G η )f , f i = (∂ j + G η j )(ū i ūj -ūi ūj )
The pressure equation was obtained by applying divergence to the LES-Langevin equation with the help of the continuity equation :

∆p = -L η p + J(u, v, w) -(∇ + G η ). l -(∇ + G η ).((∇ + G η ).R GSM ), where J(u, v, w) = 2[∂ x u∂ y v -∂ y u∂ x v + ∂ x w∂ z u -∂ x u∂ z w +∂ y w∂ z v -∂ z w∂ y v + ∂ x η(∂ y u∂ z w -∂ z u∂ y w)],
and the velocity u = (u 1 , u 2 , u 3 ) = (u, v, w).

The transformation of equations (3.5) into the spectral space, gives:

         ∂ t ū + N Gu + l = -∇ p -G η pk + 1 Re ∆ ū + 1 Re L η ū -( ∇ + G η ). R GSM ∇. ū = -G η . ū, ∂ t l + N Gl = 1 Re ∆ l + 1 Re L η l + N Gf u + ξ, (3.6)
where 2 are the gradiant and the Laplacian operators after the Fourier Transform in spanwise direction. The anti-aliasing of the nonlinear terms calculation will be explain in section 3.1.3.

∇ = (∂ x , ∂ y , ıkβ) T , ∆ = ∂ 2 x +∂ 2 y -(kβ)

Discretization in wall-normal direction

Integration of the LES equation by a high order scheme is of particular interest in LES [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF]. Low order schemes can generate spurious numerical dissipation. This does not mean that all the high order discretization scheme are not dissipative, but the use of such scheme reduces numerical errors.

The distinction of the numerical errors influence from the LES closure contribution to the modeled flow is difficult. The numerical errors advect both the computation of the extra LES terms and the instantaneous velocity, which also depends on the numerical scheme. So, it is preferential to choose the high order discretization scheme to test the effect of LES closure.

The use of a streamwise spectral scheme seems difficult because of non-homogeneous flow properties in the case of a curved wall. A high order finite difference scheme allows to minimize the numerical dissipation and to use the non-homogeneous boundary conditions. The streamwise finite difference scheme on the uniform grid with N x points was chosen of 4 th order for Laplacian term and of 8 th order for others. It is non-centered near the inflow and outflow longitudinal coordinates and centered otherwise. The computation of a Laplacian operator is implicit in time.

In the near wall region the importance of small scales increases. The largest possible scale size near the boundary can be roughly estimated comparable to the wall distance. The velocity spectrum proves it, when the intensity of the largest scales decrease with wall distance y + . The integration of the LES equations needs a refinement of the numerical grid to capture the most energetic scales. The collocation-Chebyshev method allows a high-order derivations and refines the mesh near the wall. The description of the method can be found in [START_REF] Canuto | Spectral methods in fluids dynamics[END_REF]. The collocation points are given by

ζ = cos jπ N y , j = 0, . . . , N y ,
where the number of points is equal to N y + 1 and the polynomial degree is N y . The derivation of any function f is obtained by a recursive formula using the matrices D2 i,j

and Di,j :

df dζ (ζ i ) = Ny j=0
Di,j f (ζ j ),

d 2 f dζ 2 (ζ i ) = Ny j=0 D2 i,j f (ζ j ).

Time discretization

A second order Euler backward scheme is used for the time discretization. The diffusive terms with constant coefficients and the large scale pressure term are evaluated implicitly. The other terms are evaluated by a second order Adams-Bashforth scheme. This leads to represent the equations in the form:

     ( ∆ -3σ) ūn+1 = Re ∇ pn+1 + F n,n-1 u , ∆ pn+1 = F n,n-1 p , ( ∆ -3σ) ln+1 = F n,n-1 l , (3.7)
where

     F u = -Re N Gu -Re l -Re ( ∇ + G η ) p + L η ū -( ∇ + G η ). R GSM F p = -L η p + J -[( ∇ + G η ). l] -[( ∇ + G η ).(( ∇ + G η ). R GSM )], F l = -Re N Gl -L η l -Re N Gf u -Re ξ The superscripts "n, n -1" means [•] n,n-1 = 2[•] n -[•] n-1 .
The resolution of the 3D system is represented as N z 2D equations. The equations are almost independent, except for the nonlinear terms N Gu , J k , N Gl , N Gf u , R GSM . These terms are computed in the real space and then transfered to the spectral space in z-direction:

( ū, l, l ) BF T -→ (ū, l, l ) -→ (N Gu , J, N Gl , N Gf u , R SM ) F F T -→ ( N Gu , J, N Gl , N Gf u , R GSM ).
Here "BFT" and "FFT" mean the backward and forward Fourier transformations respectively. They are performed using Fast Fourier Transform subroutines (FFT). The aliasing is removed by the 3/2 rule, as it shown by C. Canuto et al. [START_REF] Canuto | Spectral methods in fluids dynamics[END_REF]. The modes between N z and 3/2N z are set to zero after the computation of the nonlinear terms. So, the Fourier transformation is done for 3N z /2 points, which increase the time for the FFT computation.

Boundary conditions

No-slip boundary conditions are used for the top and bottom walls of the channel. The spanwise direction conditions are periodic.

The outflow velocity conditions were taken to respect the advection:

∂ t ū + u c ∂ x ū = 0, ∂ t l + u c ∂ x l = 0,
where u c is a mean outflow streamwise velocity: u c = 1/L y Ly 0 u 1 (x out ) dy. The l-vector conditions were obtained from its definition and the outflow velocity conditions. These outflow conditions do not allow a correct calculation of the pressure terms. So, the computational domain for these pressure should stop before a critical distance from the outflow boundary. The choice of the mean convection velocity u c allows to reduce the size of this buffer zone and to reduce the computer time.

The inflow conditions are of recycling type: ū(x in , y, z) = ū(x per , y, z), l(x in , y, z) = l(x per , y, z).

The choice of the recycling section x per depends on the outflow velocity which can change the mean velocity in the channel. At x per the incompressibility conditions must be satisfied.

Diagonalization method

The diagonalization of the Laplacian operator was applied in the direction normal to the wall. The 2D problem has the form of the Helmholtz equation:

(∂ 2 x + ∂ 2 y -k 2 β 2 -σ)v = f. (3.8)
The problem can be reformulated for the matrices V and F , defined as:

V ji = v(y j , x i ), F ji = f (y j , x i ), 1 ≤ j ≤ N y , 1 ≤ i ≤ N x .
Then, equation (3.8) takes the form:

D 2 y V + V D 2 T x -(k 2 β 2 + σ)V = F, (3.9)
where D 2 y is a collocation-Chebyshev operator and D 2 T x is a transposed matrix of a streamwise differencing. It is possible to find the eigen values of the y-differencing operator in a preliminary calculation:

D 2 y = P ΛP -1 , where the matrix Λ has a diagonal form. Equation (3.9) can be now reformulated as:

Λ Ṽ + Ṽ D 2 T x -(k 2 β 2 + σ) Ṽ = F ,
where

Ṽ = P -1 V, F = P -1 F. (3.10)
The problem now contains N y 1D systems in the streamwise direction:

[D 2 x + (Λ j -k 2 β 2 -σ)] Ṽ T j = F T j , 1 ≤ j ≤ N j (3.11)
with Ṽj = (ṽ j1 , . . . , ṽjNx ) Fj = ( fj1 , . . . , fjNx ).

So, the solution of the 2D problem (3.8) is obtained in 3 steps. The first one is the preliminary calculation of the eigen values and eigen vectors of the D 2 y matrix. Then F is calculated by equation (3.10). Finally, equation (3.11) gives Vj and then V j by equation (3.10).

Projection method

A Projection method was used to solve the system (3.7). The procedure of projection permits to satisfy the incompressibility condition by the projection of the preliminary calculated velocity to the divergence-free field.

The intermediate pressure p * is calculated by the formula:

∆p * = [ Ĵ (u, v, w) -L η p] n,n-1 (3.12)
from the values of two preceding time steps "n" and "n-1". Neumann boundary conditions are used and the Laplacian is represented by the form discussed by G. E. Karniadakis et al. [START_REF] Karniadakis | High-order splitting methods for the incompressible Navier-Stokes equations[END_REF] to stabilize the solution:

∆u = ∇(∇.u) -∇ × (∇ × u).
The Neumann pressure boundary conditions take the form:

∇ p * .n = - 1 2∆t (3 ūn+1 -4 ūn + ūn-1 ) .(n + n η ) -( ū. ∇) ū + ( ū.G η ) ū + G η p + 1 Re ( ∇ × ( ∇ × ū))] n,n-1 .(n + n η ) -[ ∇ p] n,n-1 .n η ,
where n and n η are defined in Eq. (3.4). The solution of Eq. (3.12) gives the intermediate velocity solution of Eq. (3.7):

( ∆ -3σ) ū * = Re ∇ p * + f n,n-1 with ū * Γ = ūn+1 Γ ,
where Γ is boundary of domain.

Let us define the scalar φ as a pressure correction term:

∇φ = - 3 2∆t ( ūn+1 - ū * ), ∇ ūn+1 = 0.
The application of a divergence operator to Eq. (3.13) gives the iterative solution for φ:

∆φ k+1 = 3 2∆t ( ∇. ū * + G η . ū * ) -L η φ k η (3.13)
with boundary conditions:

∇φ k+1 .n + G η φ k+1 .n η = ∇φ k .n η . (3.14)
In our case, Neumann boundary conditions were used. Equation (3.13) is a fixed point problem. This elliptic problem arises because of extra mapping term L η φ k η which can not be resolved implicitly. Equation (3.14) is the condition for a zero-gradient of φ on the boundary. In the case of the plane channel, this procedure is not necessary.

Finally, the correction terms are given by the relation:

pn+1 = p * + φ ūn+1 = ū * - 2∆t 3 ( ∇φ + G η φ).

Performances of the algorithm

The performances of the algorithm are imposed by the use of parallel hardware supporting the Message Passing Interface (MPI). The MPI parallelization allows to use a distributed-memory architecture. The data exchanges is done by communicators, which performances depend on the used algorithm and on the number of processors.

The advantage of the present algorithm is the Fourier space discretization of the solution in the spanwise direction. The number of calculated modes is two times smaller (N z /2) than the total number of Fourier modes N z . Other modes can be obtained by the complex conjugation of the calculated modes. So, most of the time each of the N p processor solves N z /2/N p systems of 2D equations. The MPI communications are necessary for the calculation of the nonlinear terms.

The calculation were done at the IDRIS and CRIHAN computing centers using IBM-SP4 (∼ 5,4 Gflop/s/proc) and IBM-SP5 (∼ 7.0 Gflop/s/proc) parallel computers.

Validation for turbulent flow

Comparison to a reference DNS at Re τ = 180

For stability reasons the LML code uses a symmetry in the spanwise direction. This symmetric behavior does not restrain the generality of the solution or the performances of the algorithm. This is equivalent to the simulation of a half of the domain with the same numbers of symmetrical or anti-symmetrical modes, which reduces the necessary memory by a factor 2.

The validation of the DNS algorithm for turbulent channel flow was shown by Marquillie et al. [START_REF] Marquillie | Direct numerical simulation of a separated channel flow with a smooth profile[END_REF]. The mean velocity and Reynolds stress were compared to the data of Jimenez [START_REF] Del Alamo | Spectra of the very large anisotropic scales in turbulent channels[END_REF] at Re τ = 180 in Fig. 3

Smagorinsky model statistics, Re τ =400

The quality of the LES scheme was first tested with the Smagorinsky model on a channel flow at Re τ = 395. The Smagorinsky constant was chosen as c s = 0.1 [START_REF] Hughes | The multiscale formulation of large eddy simulation: Decay of homogeneous isotropic turbulence[END_REF]. The comparison of a mean velocity profile and Reynolds stress with Kasagi DNS data [START_REF] Iwamoto | Fully developed 2d channel flow at Re τ = 400[END_REF] are given in Fig. 3.4 (a-d).

The mean velocity profile shows rather good overlap with the DNS. Small differences are located in the center of the channel and in log-region. The mean velocity is overestimated in the buffer-layer. Streamwise and wall normal RMS velocity of the LES are overpredicted behind their maxima, which agrees with a result of Ramakrishnan et al. [START_REF] Ramakrishnan | Partition selection in multiscale turbulence modeling[END_REF].

The Smagorinsky model was well resolved (∆ + x ≈ 9.7, ∆ + z ≈ 19.4, ∆ y min ≈ 0.21 and ∆ y max ≈ 12.9). The results confirm the known property of the Smagorinsky closure [START_REF] Holmen | Sensitivity of the scale partition for variational multiscale large-eddy simulation of channel flow[END_REF][START_REF] Hughes | Large eddy simualtion of turbulent channel flows by variational multiscale method[END_REF][START_REF] Ramakrishnan | Partition selection in multiscale turbulence modeling[END_REF] which was discussed in chapter 1.

Filtering

The computation of the nonlinear terms N Gu , N Gl and N Gf u of the LES-Langevin model (Eq. (3.5)) uses a scale separation. The choice of the filter changes significantly the turbulent intensity and subgrid scale dissipation. In the spectral space, the deconvolution filter of a function φ has the form:

φk = Ĝk φ k ,
where Ĝk is a Fourier transformation of the kernel function G. The velocity structures of size r are influenced by the deconvolution filter at wavenumber k ≈ π/r.

A spectral cutoff filter is a natural choice for a Navier-Stokes solver with a spectral derivative scheme. As can be seen from Fig. 3.5 it's convolution function has zero value for the subfilter scales and is equal to unity for resolved modes. The filtering is applied in three steps: a velocity transformation to the spectral space, a truncation of the sub-filter modes and a backward Fourier transformation.

A discrete filter for a function φ in physical space is defined as [START_REF] Vasilyev | A general class of commutative filters for LES in complex geometries[END_REF]:

φj = L j l=-K j w j l φ j+l ,
where w j l are weight factors. Near the walls, the summation domain is not symmetric

K j = L j .
The discrete filter can be used in physical space for a solver with non-periodic boundary conditions. In the present code, this filter is convenient for the streamwise direction, where finite differencing is used. The spectral form of this filter depends on the weights w j l . An example of discrete filter is presented in the Tab. where the first and second lines represent the weights for non centered filtering used near the boundaries. This filter scheme satisfies further properties [START_REF] Vasilyev | A general class of commutative filters for LES in complex geometries[END_REF]:

L j l=-K j w j l = 1, L j l=-K j l m w j l = 0, m = 1, . . . , n -1, L j l=-K j (-1) l w j l = 0, (3.16)
where n = 3. The first and second conditions put the first three moments of the filter to zero. The last condition is an additional constraint of vanishing of the filter cutoff frequency Ĝ(π/∆) = 0 [START_REF] Vasilyev | A general class of commutative filters for LES in complex geometries[END_REF]. The discrete filter has not a well-defined filter size and the constraints of Eq. (3.16) improve the spectral properties of the filter and therefore the scale separation. More vanishing moment would lead to a filter with spectral property closed to cutoff filter. The filter reaches the cutoff spectral form in the limit of infinite K j and L j . This reduces a damping of the turbulent kinetic energy caused by the filter.

Filter application in the code

A spectral cutoff filter was applied in the wall normal and spanwise direction. The velocity and the pressure fields in these directions are represented by Chebyshev and Fourier polynomials respectively, then the spectral cutoff scale separation is easy to implement. In the streamwise direction, two possible filters were applied: a spectral cutoff and a discrete filter of Tab. 3.15. The commutation error of this last filter with streamwise differencing is of order O(∆ 5 ).

A comparison of the discrete filter from Tab. 3.15 with a spectral cutoff one is given in Fig. 3.5. A discrete filter size can be defined in a range corresponding to 0 < G(k) < 1, but this definition does not allow a clear definition of the smallest resolved scale. Also, discrete filtering damps the resolved velocity scales if G(k) < 1. This reduces the resolved turbulent energy and should be distinguished from a dissipation caused by subgrid scales.

The use of the spectral cutoff filter in the streamwise direction provides a clear scale separation. Unfortunately, the present solver algorithm does not provide a perfectly periodic streamwise boundary conditions, as it was noted in the Sec. 3.1.4. A 2D (x-y) Neumann problem of equation 3.11 does not satisfy exactly the periodic conditions at points x = x in and x = x per . Even a decrease of the time step dt can not improve this difference. Also, to conserve a constant flux, the inflow velocity must be corrected at each time step by a factor. During numerical tests, a difference between mean inflow at x per and mean inflow after correction was of order of 0.1%. This does not significantly perturb the simulated velocity.

Two solutions were proposed to modify the cutoff filtering in the streamwise direction (see Fig. 3.6). The first one is to symmetrize the field defined on the full domain (periodic part and buffer) in the streamwise direction and to apply the cutoff filtering to this symmetrized field. This mirror symmetrization of the function is illustrated in Fig. 3.6 (a). The symmetrization forces the boundary values to be equal. The Fourier series converges to the symmetrized function if the function is sufficiently smooth. A problem

Conclusions

The presented numerical code was shown to be suitable for LES turbulence modeling. It uses a high order differencing scheme, which minimize the spurious energy dissipation by the numerical scheme. The application of the numerical algorithm to DNS agree well with the turbulence properties obtained by an other DNS solver. The simulation using a well-known LES model is also in good agreement with other simulations found in the literature. An advantage of this algorithm is that it can be applied to the channel flow with a curved wall, which gives a possibility to study the turbulence with variable pressure gradient.

Several constraints for a resolved and subfilter scale separation are due to the type of streamwise boundary conditions used. The use of discrete filter does not allow an exact resolved and subgrid scales separation. The application of the scale separation by discrete filter at each time step leads to an extra damping of the resolved scales. An other possi- bility is to use a spectral cutoff filtering. This allows a perfect scale separation without numerical damping of the resolved scales. The problem of the boundary conditions which are not exactly periodic can be resolved by the regularizing method described in Sec.

Chapter 4 Results

In this chapter, we present the results of the model in a channel flow geometry at Re τ = 600. Given the performances of the model in isotropic turbulence, we a priori expect our model to perform well in the center of the channel and that a tuning of the parameter will be necessary to obtain good performances near the wall. As will be seen in the sequel, this natural expectation is not valid, and this sets important questions regarding our LES simulations that will be emphasized in conclusion.

One advantage of the present model is the possibility to perform LES and DNS with the same numerical scheme. This makes it possible to compare the mean and RMS velocity profiles in dimensional and non-dimensional form, and have access directly to the wall shear stress. We may then use this quantity to build a convenient scalar parameter to try and estimate the performances of a given model, through [START_REF] Meyers | Evaluation of Smagorinsky variants in large-eddy simulations of wall-resolved plane channel flows[END_REF]:

δ τ ≡ u τ,LES -u τ,DN S u τ,DN S ,
where u τ,DN S and u τ,LES are the friction velocity calculated from the DNS and LES data. The reference value of the friction velocity was calculated from the DNS results of [START_REF] Marquillie | Direct numerical simulation of a separated channel flow with a smooth profile[END_REF]: u τ,DN S = 0.0490. The model we described in chapter 3 includes several control parameter.

In the sequel, we test the influence of the different parameters with reference to DNS or to other LES models. Before entering these tests, we however wish to make an important remark regarding our DNS.

Preliminary remarks on the DNS data

The reference DNS database of turbulent channel flow at Re τ ≈ 600 we use in the present thesis was obtained by Marquillie et al. [START_REF] Marquillie | Direct numerical simulation of a separated channel flow with a smooth profile[END_REF]. However, interesting questions about this reference DNS appear when comparing its statistics with the DNS data obtained by Alamo et al. for the turbulent channel flow at Re τ ≈ 550 [START_REF] Del Alamo | Direct numerical simulation of the very large anisotropic scales in a turbulent channel[END_REF][START_REF] Del Alamo | Spectra of the very large anisotropic scales in turbulent channels[END_REF], at a resolution slightly smaller. Parameters of these two DNS simulations are summarized in Tab. 4.1. In the streamwise direction, LML DNS uses finite differences, so the number of correspond-ing Fourier modes (k x,max ) is two times smaller than the number of mesh points (nx): k x,max =nx/2. DNS of Alamo et al. uses spectral discretization in streamwise direction where k x,max =nx. In wall-normal and spanwise directions both DNS use Chebyshev and Fourier discretizations. One sees that for a box that is four times larger in streamwise and wall-normal, the DNS of Alamo uses a number of modes that is less than three times the number of modes in the DNS of LML. It is therefore less resolved. The RMS profile are defined from the RMS velocity u i rms given by:

u i rms ≡ (u i -u i ) 2 .
The profiles display similar features close to the wall in both simulations. However, large discrepancy appear as one goes towards the center of the channel. First, the LML velocity profile is smaller than the mean profile of Alamo for y + > 400, by about 0.05%. Reynolds stress components of LML DNS are also smaller than the Alamo counter part for y + > 10, by a factor that can be as large as 30%. In the sequel, we shall observe similar differences in the channel center when comparing our LES simulations (less resolved than the DNS) with the DNS reference, so it is important to understand where this difference comes from. We have no definite answer yet. One possibility is that the difference comes from differences in the computational box size. The turbulence in the center of the larger box can contain more large scale structures than in the smaller domain. This could explain why Reynolds stress obtained by Alamo et al. are larger than those of Marquillie et al.. Another possibility would be a numerical instability, due to the lack of resolution, that would produce large scale features that are absent in the LML-DNS. To check this possibility, one needs to perform additional simulation changing the box size and the resolution, but these simulations have a high numerical cost and we were not able to perform them yet. After this preliminary test showing the importance of resolution and box size, we now turn to tests of our turbulence model and its various components.

Eddy-viscosity

A first parameter of our LES-Langevin is an eddy viscosity to model a part of the SGS tensor τ ij : the Subgrid scale stress R ij . In order to test this term, we put the Cross-stress tensor C ij to zero in all the simulations described in the present section. The choice of an eddy viscosity in our model is a priori vast. Many eddy-viscosity models have been proposed in the literature [START_REF] Sagaut | Large eddy simulations for incompressible flows[END_REF]. In chapter 1 we have described some of these approaches: Smagorinsky model, Sec. 1.2, its modification by a Van-Driest damping, Sec. 1.2.2, Germano dynamic procedure, Sec. 1.2.3, etc. Here, we focus on the Smagorinsky model and its modification by a Van-Driest damping, since it is easier to implement.

Smagorinsky-Van Driest model

In the sequel, we adopt a simple Smagorinsky closure with a Van-Driest damping, parametrized by a Smagorinsky constant c s and Van-Driest damping coefficient A + . In order to capture the best eddy-viscosity model, we have varied these two coefficients. A summary of the different values is given in Tab. 4.2, along with the resulting value of δ τ . The abbreviation "CDNS" means a coarse resolved DNS on LES mesh. Note that with c s = 0, we are in the coarse DNS case, where the simulation is under-resolved and no model is used to parametrize the subgrid scales.

The coarse resolved DNS is the worst situation one can have. Any model we have should at least do better than this coarse DNS, in order to be of practical use. As one intriguing feature is that whatever the value of the viscosity or of the damping, we cannot change the value of the RMS profiles with respect to their value in the coarse-DNS case. This value on the other hand is similar to that obtained in the Jimenez-code. Since the LES and the CDNS simulation are less resolved than the DNS, we may be observing an undesirable effect of the resolution, producing a spurious large scale instability. To check this point further, we have tested an additional procedure, based on negative viscosity.

The negative viscosity model

Indeed, we know from theoretical results of Dubrulle et al. [START_REF] Dubrulle | Eddy-viscosity of parity-invariant flow[END_REF] that shear geometry naturally produces negative eddy-viscosities above a certain value of the Reynolds number. This results in large scale instabilities, allowing transfer from small to large scales of the flow. This instability may be responsible for the increase of RMS energy at the center of the channel. As a test case we have therefore implemented an additional negative eddy-viscosity effect through:

ν t (x) = γ * ν sm (x), for y + ∈ [0, y + γ ], ν t (x) = f vd * ν sm (x), otherwise
where ν sm is a classical Smagorinsky viscosity, γ < 0 is a coefficient, f vd is the Van Driest function defined by Eq. (1.20). The idea is too look if this additional negative eddy viscosity produces an enhancement of the RMS velocities at the channel center.

The negative value of the energy transfer produced by the fluctuating part of SGS stress was revealed in the work of Härtel et al. [START_REF] Härtel | Subgrid-scale energy transfer in the near-wall region of turbulent flows[END_REF]. This effect was found strong in the buffer layer with an intensity which depends on the size of the explicit filter. The same conclusions can be derived by observation of the velocity profiles: mean and RMS velocity profiles (see Fig. 4.5) are rather insensitive to the modeling of the vector l. In contrast, large differences in velocity statistics are caused by the different filtering in the wall-normal and streamwise direction. When the resolved-subfilter scale separation in wall-normal direction k y moves to the small scales, we see a significant improvement in the u + RMS. The discrete filtering in streamwise direction leads to a decrease of the turbulent kinetic energy and mean velocity.

To understand the impact of the vector ξ onto the vector l in the different models, we compare in Fig. 4.6 the mean values of |l| and |ξ| (a) and the power spectra of l and ξ (b-d). The mean values of |ξ| are larger than the mean values of |l| by several orders of magnitude. The spectral variation of l and ξ are similar in the spanwise direction. The streamwise direction spectra of ξ has smaller maxima than l because of the finite differencing, as it was discussed in chapter 2. So, tests N4 and N6 of Tab. 4.4 increase the intensity of the vector l by several orders of magnitude with respect to the model N3. Yet, almost no difference can be seen on the mean or RMS velocity profiles. This means that the vector l is probably too small to significantly influence the mean velocity We can conclude that the presence of the turbulent force with the forcing ξ derived from the a priori tests does not change significantly the profiles, at least in absence of an eddy viscosity modeling. This is contradictory with the fact that the tensor C ij model should improve the SGS eddy-viscosity modeling of R ij . The possible reason for this

is ---, l y is • • • • • • , l z is ----, ξ x is -• -, ξ y is -• • • -and ξ z is ---.
results is the unrealistic value of the viscous part of the vector l equations. In order to check this fact, we now perform the same test in presence of an eddy-viscosity ν t in Eq. 2.6. This eddy-viscosity models the nonlinear resolved -non-resolved velocity terms multiplied by the kinematic viscosity ν. Nevertheless, the l vector should be filtered to produce the subgrid scales interaction with the resolved scales. The spectral form of l x is presented in Fig. 4.12. The filtering of l i decreases the amplitude of variation of the vector l by a factor of ≈ √ 100 = 10. So, the modeled interaction of the subgrid scales with the resolved scales are smaller than the interaction of the resolved scales.

The importance of the smallest resolved velocity scales for the modeling of the vector l is obvious if we compare the mean |l| from Fig. 4.11 and the velocity spectra in Fig. 4.8, Fig. 4.9 and Fig. 4.10. The excessive dissipation of the small resolved velocity scale in the range of the turbulence maximum at y + ∼ 15 does not produce a significant decrease of the mean value of nonlinear term |N i |. For the nonlinear terms of the vector l equation (ū∇)l + (l ∇)ū this allows the reduction of the largest resolved modes composed by cross terms of the scales located near the filter separation scale. This explains the large sensitivity of the vector l at wall distance around y + ∼ 15.

These results show the influence of the modeling of the vector l on the Reynolds Stress. The mean velocity and the velocity spectra are less sensitive to l. Since on the overall, none of our model are very sensitive to l, this means that we need to increase its magnitude to observe effects on the mean and RMS profiles. First of all, we study the sensitivity of the vector l to the model eddy-viscosity. As we can see in Fig variation of the constant c s does not significantly changes the variation of the velocity nonlinear terms, but is important for the vector l. The decrease of the eddy diffusion intensifies the variations of the vector l.

Parametrization of the vector l governing equation

To increase the energy transfer produced by the vector l, we study its sensitivity to the different parameters of the LES-Langevin model. In Tab. 4.6 we present test simulations of the sensitivity to the friction coefficient τ f defined by Eq. 2.7 and to the eddy-viscosity ν t of the vector l equation. The skin friction term was modeled like in the isotropic case by Eq. 2.7. In test N2 we take the friction term and eddy-viscosity ν t equal to their isotropic modeling by Eq. 2.7 and ν t = 0. In test N3, we divide the friction term by a factor two τ f /2 and the diffusion term of the vector l equation by a factor 100.

The comparison of the mean value of |l| and of the power spectra of l of the tests of Tab.4.6 are shown in Fig. 4.13. A decrease in the friction times allows for a significant increase of the intensity of l . In contrast, a finite correlation time significantly decreases the intensity of l. The relative maximum of the amplitude of l at y + ∼ 150 (Fig. 4.13, (a)) is reduced compared with its value at the center of the channel. The variation in the vector ξ correlation time changes the intensity of the modeled vector l . The spectra of the vector l in Fig. 4.13(b-e) confirms the features of |l x | mean profiles. Despite these seemingly large differences on the intensity of l, the resulting influence onto the mean and RMS profiles is very small, as can be seen in Fig. 4.14. Cases N2 and N4 produce the same mean velocity profile that differs at y + > 200 from case N3. The Reynolds Stress of cases N2 and N4 do not differ at small wall distance. Therefore, it appears difficult to tune the parameters of the model to observe a large effect. One problem is that the intensity of l is too weak to influence significantly the result. As we just have seen, we can increase its intensity by decreasing the friction coefficient, but, due to a numerical instability, we have not succeeded to go beyond a certain friction time. More work is needed to be able to work with model with larger intensity of l .

N Model Friction

c s A + ν + ν t ν ξ nx ny nz k x k y k z Re τ δ τ 1 DNS 618 0 2 LES-L τ f 0.

Comparison of full LES-Langevin model with DNS

The analysis of the previous paragraph shows that the presence of dissipation does not change significantly the contribution to the mean and RMS velocity profiles of the turbulent force in the present range of parameters. To confirm this point, we finally present a comparison of the full LES-Langevin model with a DNS, in order to try and obtain our best model. The different models are presented in Tab. 4.7. The best shear stress prediction is achieved for c s = 0.075. Reynolds number Re τ varies about 8% in comparison to the reference DNS (N1) because of the shear-stress error.

, LES-Langevin (N2) is ----, LES-Langevin (N3) is -• -, LES-Langevin (N4) is -• • • -.
does not allow a satisfactory improvement because of the shift of the energy maximum to higher y + . Velocity spectra are shown in Fig. 4.16 (u), Fig. 4.17 (v) and Fig. 4.18 (w). In the wall-normal direction (pictures (e)), k y references the Chebyshev eigenfunctions and the corresponding spectra represent the square of Chebyshev coefficients. All the curves show strong dependence of the small resolved scale to the constant c s near the turbulence maximum. The results demonstrate that it is not possible to tune the constant c s to improve the velocity spectra at all wall distances. For example, at y + ≈ 17 the best k x spectra corresponds to c s ∈ [0.065; 0.1], but at the center of the channel (y + ≈ 620) the constant c s = 0.1 would be suitable. The large RMS values found in Fig. 4.15 are also due to the wrong modeling of k z resolved scales of ū (Fig. 4.16, (c)) and w (Fig. 4. 18,(c,d)) and to the bad model of the k y spectrum of all three velocity components.

Influence of the LES resolution

As a last test of our model, we study the influence of resolution on the modeling. Indeed, as mentioned in the chapter 1, a model that captures the most energetic structure of the flow will be better suitable to describe the dynamics of the whole model than a model that leaves out the dynamically important features. In our study, we uses two resolutions of the resolved scales: (k x,max , k y,max , k z,max ) T = (64, 65, 64) T and (k x,max , k y,max , k z,max ) T = (64, 49, 32) T . The resolved turbulent energy, produced by a LES model naturally depends on the numerical grid. For the same computational domain the finer grid supplies more turbulent energy. So, the comparison of the velocity statistics is done with the DNS velocity truncated to the present LES resolutions. The comparison of the mean velocity T is better adapted for the LES computations of turbulent channel flow at Re τ = 600. We can see from the Fig. 4.19 that all the LES resolutions located between the total DNS and (k x,max , k y,max , k z,max ) T = (64, 65, 64) T are better suitable for LES approach than the coarser mesh. There is therefore an upper limit on the resolution one can use to have a correct LES model, and this remark sets some limitation on the expected performances of LES models. The coarse resolution in spanwise direction of presented LES-Langevin simulations is not sufficient to capture the turbulent structures in the buffer layer. This leads to the appearance of spurious streak-like structures ("superstreaks") and to the overestimation of the Reynolds stress, as it can be seen for u + in the Fig.

Conclusions

The present study of the LES-Langevin model shows mitigate results. On the one hand, we have seen that variation of the parameters allows for a variation of the intensity of the turbulent force, an essential ingredient of the model that aims to model non-local energy transfer between resolved and subgrid scales. However, in the range of parameters we were able to work with, the final intensity of the turbulent force was not large enough to have a significant influence on the mean and RMS velocity profiles. Another ingredient of the model is a turbulent viscosity, aiming to reproduce local interactions between the subgrid scales. However, we have seen that the adopted method, a Smagorinsky approach, was not able to reproduce the RMS and mean profile in the channel center in a satisfactory manner. In fact, quite surprisingly, the non-dimensional values at these location appear quite insensitive to the modeling and display a difference that is comparable with the difference between two DNS with different box size and resolution. To check whether this difference could be caused by a large scale instability, we have implemented a negativeviscosity procedure to see whether we could reproduce comparable differences. We indeed observed increased RMS values, but mostly located near the channel boundary. Clearly, more work would be needed to understand this point.

Chapter 5

Conclusion and perspectives

Comparison of LES-Langevin to other LES models

In the present contribution new LES-Langevin model was presented. The approach is looking for a space-time description of the subgrid scales. The model assumes the rapid variation of the subgrid scales compared to the slower dynamics of resolved scales. This RDT method was previously applied by Dubrulle et al. [START_REF] Dubrulle | A model for rapid stochastic distorsions of small-scale turbulence[END_REF] and Laval et al. [START_REF] Laval | Langevin models of turbulence: Renormalization group, distant interaction algorithms or rapid distorsion theory?[END_REF] to develop a Langevin stochastic equation for the small turbulent scales. From the LES point of view, the description of the dynamics of the subgrid scales allows the computation of the SGS stress tensor (see Chapter 1 for the definitions of the SGS stress tensor). The real interest of the LES approach is to model a three-component vector ∂ j τ ij , which is the gradient of the SGS stress tensor. Also, the SGS stress tensor contains all the scales from the largest resolved scales to the dissipative scales. It is more advantageous to compute directly the resolved scales part of τ ij . This motivate us to develop a model directly for the gradient of the SGS stress tensor. The previous successful application of the model to isotropic turbulence by Laval et al. [START_REF] Laval | A LES-Langevin model for turbulence[END_REF] motivated the validation for a more complex flow as the channel flow.

The particular interest of the LES-Langevin approach is the opportunity to model the inverse energy transfer from the subgrid to resolved scales. This is a very important mechanism for wall turbulence because the initial instability of the flow is generated by the small scales near the boundary. The amplification of the initial instabilities near the wall and interaction of these small scales with the slow large scale flow establish the properties of the wall turbulence for particular boundaries. These processes are defined through the interaction of the resolved with subgrid velocity scales and should be modeled by SGS stress tensor τ ij . So, there is a particular motivation to correctly model the backscatter by LES models and the plane channel flow is a good simple case of the wall turbulence.

The comparison of the LES-Langevin model to the self-similarity approach (Bardina version [START_REF] Bardina | Improved subgrid scale models for large eddy simulation[END_REF], sec. 1.2.4) shows some advantages of the LES-Langevin model. The selfsimilarity approach is based on the similarity hypothesis of the velocity scales. So, the modeled subgrid scales have to be governed only by the internal turbulent mechanism of the energy transfer. This is not true if the distance to the wall is comparable to the LES filter size or if there are external forces acting on the fluid. The dynamics of the turbulent force l is governed by the 3D dynamical equation. This equations can be solved in a situation when the hypothesis of the self-similarity does not hold. This is particularly the case of turbulent flow near walls or in presence external forces.

The same criticism can be applied to the Mixed model (described in section 1.2.5), which is the combination of the Similarity and Smagorinsky models. The difference between the Similarity and Mixed models is in the eddy-viscosity term which mostly represents the interaction of subgrid scales.

The comparison to eddy-viscosity LES models shows the possible advantages of LES-Langevin method. In chapter 1 we describe some eddy-viscosity approaches: Smagorinsky model [START_REF] Lilly | On the application of eddy-viscosity concept in the inertial sub-range of turbulence[END_REF][START_REF] Smagorinsky | General circulation experiments with the primitive equations[END_REF] (section 1.2.1), Smagorinsky model with Van-Driest damping (section 1.2.2), Dynamic model [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF] (section 1.2.3), PASDD model [START_REF] Porté-Agel | A scale-dependent dynamic model for large-eddy simulation: application to a neutral atmospheric noundary layer[END_REF] (section 1.2.6) and LASI model [START_REF] Meneveau | A Lagrangian dynamic subgrid-scale model of turbulence[END_REF] (section 1.2.7). These models base the computation of the eddy-viscosity from the resolved scale strain-rate (the Smagorinsky approach) and change the model constant by a function dependent on space and time. The Van-Driest damping procedure allows for a reduction of the eddy-dissipation term by decreasing the characteristic length ∆ in the Smagorinsky formulation of the eddy viscosity as a function of the shear stress. The large shear stress increases ∆ and the modeled diffusion of the resolved scale velocity, while the small shear stress tends to reduce the modeled diffusion with ∆. Our simulations of the Smagorinsky model with Van Driest damping presented in chapter 4 serve as a reference for the LES-Langevin model. We see that the friction velocity, Reynolds stress, mean velocity profiles and energy spectra are well predicted by the Smagorinsky model. So, the turbulent subgrid scale motion in a plane channel can be well modeled by a simple eddy-viscosity approach, if the numerical grid resolution is sufficiently fine. In our case, at Re τ ∼ 600, the numerical grid has to contains all the modes up to k c = [START_REF] Piomelli | Model consistency in large eddy simulation of turbulent channel flow[END_REF][START_REF] Porté-Agel | A scale-dependent dynamic model for large-eddy simulation: application to a neutral atmospheric noundary layer[END_REF][START_REF] Piomelli | Model consistency in large eddy simulation of turbulent channel flow[END_REF]. This represents a gain of 144 in terms of memory compared to the DNS. The LES-Langevin model gives also satisfactory results with this geometry, but requires more computational resources. For the same k c the LES-Langevin computation reduces the necessary memory by a factor 24 compared to DNS.

Nevertheless, in more complex geometries such as the rough wall, the computation of the mean shear stress is impossible because of the absence of the homogeneous direction, so the correction of the Smagorinsky eddy-viscosity has to be done in a different way. The other eddy-viscosity models (Dynamic, PASDD and LASI) uses the Germano identity and similarity hypothesis to allow the dependence of the constant on space and time. PASDD model also needs at least one homogeneous direction, which restricts the set of possible geometries. The Dynamic model is one of the most-popular LES models. In this approach the constant C is computed from the minimization of the square of the modeling error of the Smagorinsky closure e ij , defined by Eq. (1.26). The model improves the computation of the near-wall eddy-viscosity [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF]. For stability reason, the constant which leads to a negative total dissipation has to be modified in order to respect positive dissipation. So the backscatter mechanism can not be modeled. LASI model also minimizes e ij * e ij averaged over path-lines. The clipping procedure of LASI model does not allow the possibility of backscatter.

The LES-Langevin approach separates the modeling of the nonlocal resolved-subgrid

COMPARISON OF LES-LANGEVIN TO OTHER LES MODELS

and local subgrid-subgrid interactions parts of the SGS stress tensor. The eddy-viscosity approach is applied to model only the subgrid-subgrid interactions. Therefore, the model does not need an homogeneous direction in the flow. The turbulent force l is able to model the backscatter and should also participate in the dissipation of the resolved scale energy. In contrast to the isotropic case [START_REF] Laval | A LES-Langevin model for turbulence[END_REF], the dissipative contribution of the vector l was not detected on the Reynolds stress of the channel flow. This is surprising because in both isotropic turbulence and channel flows, the Reynolds numbers based on the Taylor scale in the center of the channel were comparable (about Re λ ≈ 200). The modeled vector l produces the resolved scale energy and is always in competition with the eddyviscosity term which allows the dissipation of this energy. As was shown in section 4.3.2, the obtained velocity statistics are more sensitive to change of the constant c s than to change of the parameters governing the dynamics of the vector l.

The ADM model [START_REF] Domaradzki | Direct modelling of subgrid scales of turbulence in large eddy simulations[END_REF][START_REF] Stolz | An approximate deconvolution model for large-eddy simulation with application to incompressible wall-bounded flows[END_REF] (described in section 1.2.8) reconstructs the subfilter scale velocity from the resolved one by the defiltering procedure. Then the reconstructed subfilter scales are more dependent on the form of the inversed filter and its support than on the boundary conditions of the flow or the possible external forces. So, the filter defines the local 3D information about the deconvolved velocity scales. The LES-Langevin model is also dependent on the filter. But in this case, the filter rather separates the resolved from subfilter scales than defines the dynamics of the subfilter scales. The change in the modeling of the subgrid pressure and diffusion of the vector l interacts with the resolved velocity, as it is shown in section 4.3.2. The LES-Langevin allows a possibility to include the information about the physics of the subgrid scales in the turbulent energy transfer. The role of regularization factor in the LES model is mostly attributed to the friction term in the equation of the turbulent force l. This term corrects the long time behavior of the largest resolved scales of the turbulent force l. This regularizing is different from the energy conservation procedure of the ADM model, represented by Eq. 1.37. In this last model, the correction is applied locally in time (at each iteration) and directly to the modeled value of all the interactions of the mesh-represented scales with subgrid scales.

The SGEM model [START_REF] Domaradzki | Large eddy simulations using the subgridscale estimation model and truncated Navier-Stokes dynamics[END_REF] (described in the section 1.2.9) interpolates the resolved velocity scales onto a finer mesh (Eq. (1.38)) and models the subgrid velocity scales as a synthetic field (Eq. (1.39)). The SGS stress tensor is computed from these reconstructed velocity according to the definition. So, the advantage of the model is that it reproduces exactly the interactions between the resolved, subgrid and subfilter scales. Nevertheless, the subgrid and subfilter scales are produced by the synthetic field, which does not satisfy the dynamical equation for the corresponding scales. The advantage of the LES-Langevin model is that it is based on the theoretical and numerical study of the subgrid scales in the frame of Rapid Distortion Theory [START_REF] Dubrulle | A dynamic subfilterscale model for plane parrallel flows[END_REF][START_REF] Dubrulle | A model for rapid stochastic distorsions of small-scale turbulence[END_REF][START_REF] Laval | Langevin models of turbulence: Renormalization group, distant interaction algorithms or rapid distorsion theory?[END_REF][START_REF] Nazarenko | WKB theory for rapid distorsion of inhomogeneous turbulence[END_REF][START_REF] Nazarenko | Nonlinear RDT theory of nearwall turbulence[END_REF]. The time interval T between two initializations (see Eq. (1.41)) can limit the life time of the modeled SGS structures, which is very long (limited by the friction factor τ f ) in the LES-Langevin approach. Both models do not need a self-similarity hypothesis about the subgrid scale flux, and they can be applied to turbulence with complex boundaries. The explicit filtering is also essential for both models. The filter in the SGEM model is necessary to eliminate the spurious increase of the kinetic energy of the smallest resolved scales which is obvious for the resolution of the Navier-Stokes equations on a coarse grid. The LES-Langevin approach only needs a filter to separate the subfilter scales from the resolved ones.

The VMS model as well as LES-Langevin model separates the large (resolved) and small (subfilter) scales of the flow. The interaction of the resolved (large) with subgrid velocity scales is completely neglected. Only additional diffusive terms represent the subgrid velocity scales interaction with the subfilter (small) scales. The interaction between the resolved and subfilter scales is directly calculated through the nonlinear term. But in VMS, the contribution of the subfilter-subfilter interaction to the resolved velocity scales is computed exactly. Setting µ t = 0 in the LES-Langevin model leads to a similar practical realization of both LES models: the turbulent force l produces partially the same interactions of the resolved velocity scales with the subfilter velocity scales as it is computed in VMS by the exact formula (Eq. 2.1). Surprisingly, this parametrization discussed in section 4.3.1 does not bring any advantage compared to the coarse DNS at the same spatial resolution. The problem was characterized in the previous chapter by the numerous comparisons of the mean |l| profile with its a priori value. In fact, the vector l a priori computation gives the exactly resolved-subfilter interaction of velocity scales computed in VMS. Comparisons of the vector l modeled with different set of parameters and its a priori value summarize the fact that the amplitude of modeled turbulent force does not achieve its a priori value for different choice of the parameters. This seems crucial for the resolved energy cascade prediction which is well modeled by the VMS method applied to the turbulent channel flow in the work of Ramakrishnan et al. [START_REF] Ramakrishnan | Partition selection in multiscale turbulence modeling[END_REF]. For the moment, we do not have a clear explication of the difference between the LES-Langevin and VMS methods. This problem may be studied in a future work.

The RSFS model [START_REF] Zhou | A resolvable subfilter-scale model specific to large-eddy simulation of under-resolved turbulence[END_REF] (described in the section 1.2.12) is looking for an extra 3D dynamical equation for the subfilter scale velocity. As in the VMS model, all the nonlocal interactions of the resolved with subgrid scales are neglected. The interaction of the subfilter with subgrid velocity scales is modeled as "small-small" Smagorinsky closure (Eq. 1.43) used also in VMS method.

Further study of the LES-Langevin model

The transfer mechanism of the kinetic energy produced by the vector l is very important in the subgrid scale modeling. Our attempt to represent its effect in a shear flow convinces that we do not know its main properties, that should be reproduced by the model. The a priori estimation of the mean profiles, correlation time and space spectra of the vector l was used as a reference to the modeled value. This does not allow any information about the local space-time mechanism of the kinetic energy transfer, especially in wall turbulence. The velocity instabilities originating in the near wall flow can be studied only in terms of the local space and time variables. The interactions of the coherent structures with the mean and turbulent flows from the outer region have also local character. Our attempt to illustrate the local time characteristics of the subfilter terms was effectuated by the model of the subfilter pressure terms as the Markovian chains with a certain correlation time. The description of the approximations used and of the results can be found in sections 2.3.1 and 4.3.2. A single time scale was attributed to the modeled terms. We do not know how many time scales have to be taken into account

PERSPECTIVES

to model correctly the energy transfer. If the number of the important time scales is sufficiently small to describe the energy transfer, this information could be included in the LES-Langevin model.

The relative comparison of the SGEM and VMS models ensures that both of them are successful because of the modeling of the interactions between the resolved and subfilter scales. SGEM approach computes directly the tensor τ ij by its definition using subfilter and resolved scales. The most energetic part of this tensor is the cross stress C ij which is a part of the resolved-subfilter velocity scale interaction of the VMS model. The subfilter scales in both models are calculated differently, but this does not seem crucial for the statistics of the resolved scale velocity. The main difference in the modeling of the subfilter scale velocity between the models is the initialization by a synthetic field in the SGEM approach. The Smagorinsky-like closure of the subfilter scale velocity equation of VMS can be viewed as a regularizing procedure, which is realized in SGEM approach by the regular explicit filtering over the time interval T. So, further study of the mechanism of the resolved-subfilter and resolved-subgrid scales interaction could lead to important consequences for turbulence in general and for LES modeling in particular.

Also, further study of the complex non-linear terms of the resolved and non-resolved velocity would be useful. In our approach these terms are modeled by an extra diffusion of the vector l. This modeling hypothesis has to be verified locally in space and time. The fluctuations of these cubic velocity terms can be relatively large compared to the nonlinear terms of the vector l equation, which does not contain the resolved scales of the vector l.

Perspectives

The present work can be continued by a research on the interaction mechanism of the turbulence velocity scales. A better understanding of this interaction can give an information of which local mechanism has to be provided by an LES model. The actual LES-Langevin formulation can be considered as the refined description of the VMS approach because VMS takes into account a maximum of the local information about the tensor C ij and the LES-Langevin uses only mean space characteristics and a simple time parameter to approximate the modeled non-linear velocity and pressure terms, described by the stochastic forcing ξ and diffusion term ∂ i ((ν + ν t )∂ i l).

The interaction mechanism between the different velocity scales can improve the actual approach of the RDT. The LES-Langevin separates the time scales into "rapid" and "slow" time scales that can be corrected by an investigation of the characteristic time scales of C ij . The theory can also be improved by a better consideration of the different velocity structures near the boundaries.

The described mechanisms of energy transfer can help in the reconstruction of the non-resolved scales wall turbulence. This can be important for the modeling of flows in presence of complex boundary conditions.

One of the purpose of the work was to change the model of the non-resolved velocity scale to the modeling of ∂ j τ ij by a Langevin equation. This approach describes the intermittent function ∂ j τ ij by a linear dynamical equation. The reduction of the degrees of freedom of the turbulent flow is done by the representation of the non-local non-linear terms as stochastic functions. Therefore, the computation of the turbulent force through the Langevin equation can significantly reduce the computational time needed to compute the tensor C ij as it is modeled in the RDT models [START_REF] Dubrulle | A model for rapid stochastic distorsions of small-scale turbulence[END_REF][START_REF] Laval | Langevin models of turbulence: Renormalization group, distant interaction algorithms or rapid distorsion theory?[END_REF] or in the VMS or SGEM models.

Conclusions

A new LES-Langevin model was studied. The model is based on Rapid Distortion Theory and previous study of the subgrid scale velocity dynamics. The different complex terms were modeled phenomenologically from their a priori statistics. The model was studied on an example of turbulent channel flow at Re τ ≈ 600. The results for the different LES-Langevin parameterizations were compared to the DNS, coarse DNS and other LES models.

The present study shows the complexity of modeling of the near-wall flow in the frame of Langevin approach. The modeled turbulent force does not provide the correct energy transfer between the resolved and unresolved scales as it is expected for gradient of the Cross Stress tenseur. Further investigation is necessary to improve the RDT hypotheses and the phenomenological deterministic and stochastic approximations of the subfilter/subgrid velocity and pressure.

The work illustrates the sensitivity of a SGS closure to the type of turbulent flow. The LES-Langevin provided good statistics for isotropic turbulence but the modeled turbulence intensity in the center of the channel does not compare well with the DNS reference value. The a priori estimation of the performance of the model to particular boundary conditions and numerical scheme can be an issue for LES in general.

The obtained results open vast perspectives for further study of the turbulent energy transfer mechanism and application of the stochastic Langevin equation to the modeling of the subgrid scale turbulence.

  ) where ξ, η ∈ [α, β]. The mapping of [α, β] to the physical domain [a, b] with a nonregular grid is done by the function F : [α, β] -→ [a, b]. The inverse function is called f : [a, b] -→ [α, β]. The non-uniform grid spacing for x ∈ [a, b

Figure 2 . 1 :Figure 2 . 2 :

 2122 Figure 2.1: Energy spectra of the RDT models (2.5) and the equivalent DNS at t=2 for a decaying isotropic turbulence (the RDT simulation is initialized with the velocity field of the DNS at t=1.)

Figure 2 . 3 :

 23 Figure 2.3: A priori estimation and approximation of ξ x : mean profile of ξ 2 (a) and k x and k z spectra (b and c) of N1 from Tab. 2.1. The ---is a priori estimation. • • • • • • is its approximation.

Figure 2 . 4 :

 24 Figure 2.4: A priori analysis of the stochastic forcing (vector ξ y ): mean profile of |ξ x | (a), |ξ y | (b) and |ξ z | (c). Tests from Tab. 2.1: ---is N1, • • • • • • is N2, ----is N3, -•is N4, -• • • -is N5, ---is N6.

Figure 2 . 6 :

 26 Figure 2.6: A priori analysis of the stochastic forcing (vector ξ y ): mean spectra in the center of the channel of |ξ x(k x )| 2 (a), |ξ x (k z )| 2 (b), |ξ y (k x )| 2 (c), |ξ y (k z )| 2 (d), |ξ z (k x )| 2 (e) and |ξ z (k z )| 2 (f). Tests from Tab. 2.1: ---is N1, • • • • • • is N2, ----is N3, -• -is N4, -• • • -is N5, ---is N6.

Figure 2 . 7 :

 27 Figure 2.7: Time correlation f (t)f (0) / f (0) 2 at y + ≈ 17 (a) and in the center of the channel (b) of test N5 of Tab. 2.1. ---is the large scale transversal velocity, • • • • • • is the small scale transversal velocity, ----is the transversal small scale pressure gradient component.

Figure 2 . 8 :

 28 Figure 2.8: A priori analysis of the turbulent force (vector l ): mean profile of |l x | (a), |l y | (b) and |l z | (c). Tests from Tab. 2.1: ---is N1, • • • • • • is N2, ----is N3 , -•is N4, -• • • -is N5, ---is N6.

Figure 3 . 1 :

 31 Figure 3.1: Example of 3d domain discretization.

  with n = (0, 1) and n η = (-∂ x η, 0).(3.4)

Figure 3 . 2 :

 32 Figure 3.2: Example of mapping in the wall-normal direction.

  .3 (a -d). The box sizes are different for the two DNS 4π × 2 × 2π and 12π × 2 × 4π for Marquillie et al. and Jimenez respectively. The mean velocity profiles given in Fig. 3.3 (a) are in a good agreement.Also, larger RMS values were found in the center of the channel for LML data (see Fig.3.3: b, d). The difference can be explained by the difference of box sizes. As the box is much larger in the DNS of Jimenez, larger scales structures are able to develop at the center of the channel.

Figure 3 . 3 :

 33 Figure 3.3: Comparison of LML-DNS (---) and Jimenez DNS (• • •) at Re τ = 180. Mean profiles (a) and Reynolds stress (b, c, d).

Figure 3 . 4 :

 34 Figure 3.4: Comparison of mean velocity profile (a) and Reynolds stress (b, c, d) from DNS and LES-Smagorinsky at Re τ =395.

Figure 3 . 5 :

 35 Figure 3.5: Fourier transform G(k) of the symmetric minimally constrained discrete filters with 3 vanishing moments and of the cutoff spectral filter.

Figure 3 . 6 :

 36 Figure 3.6: Streamwise periodization of a non-periodic function: mirror symmetrization of the function defined on the full domain [0;L per +L buff ] (a) and linear decomposition of a function from [0;L per ] (b). Initial function: -•-•-•-, periodic function: ---, linear function: • • •.

3 . 3 . 1 .

 331 High performances of the numerical code are obtained by the use of the Fast Fourier Transforms (FFT) and the Message Passing Interface (MPI).

Figure 4 . 1 :

 41 Figure 4.1: Velocity profile (a) and velocity RMS (b, c, d) of DNS-LML ( ---) and DNS-Jimenez (• • • • • • ).

Figure 4 . 2 :

 42 Figure 4.2: Velocity profile (a) and velocity RMS (b, c, d) of Smagorinsky model with different c s compared to the DNS. All the curves are referenced in Tab. 4.2. DNS (N1) ---, CDNS (N2) • • • • • • , the Smagorinsky model with Van-Driest damping (A + = 25) and c s = 0.1 (N3) ----, c s = 0.05 (N4) -• -, c s = 0.17 (N5) -• • • -.

Figure 4 . 3 :

 43 Figure 4.3: Velocity profile (a) and velocity RMS (b, c, d) of Smagorinsky model with different A + compared to the DNS. All the curves are referenced in Tab. 4.2. DNS (N1) ---, CDNS (N2) • • • • • • , the Smagorinsky model with c s = 0.1 and Van-Driest damping with A + = 25 (N3) ----, A + = 50 (N7) -• -, A + = 15 (N8) -• • • -, A + = 5 (N9) ---.

Figure 4 . 4 :

 44 Figure 4.4: Velocity profile (a) and velocity RMS (b, c, d) of negative eddy-viscosity modeling. All the curves are referenced on Tab. 4.3. DNS-M2 (N1) ---, Smagorinsky with c s = 0.1, A + = 25 (N2) • • • • • • , LES with γ = 0.001 and y + γ = 6.7 (N3) ----, LES with γ = 0.07 and y + γ = 3.0 (N4) -• -, LES with γ = 0.001 and y + γ = 26.7 (N5) -• • • -, LES with γ = 0.005 and y + γ = 26.7 (N6) ---.

Figure 4 . 5 :

 45 Figure 4.5: Velocity profile (a) and velocity RMS (b, c, d) of LES-Langevin model. All the curves are referenced on Tab. 4.4. DNS-M1-F1 (N1) ---, CDNS-M1-F1 (N2) • • • • • • , LES-Langevin without eddy-viscosity modeling with k y = 65 (N3) ----, LES-Langevin with l i = ξ i /5 (N4) -• -, LES-Langevin without vector ξ and eddy-viscosity modeling (N5) -• • • -, LES-Langevin with l i = ξ i (N6) ---.

Figure 4 . 6 :

 46 Figure 4.6: A priori comparison of the vectors l and ξ from Tab. 2.1 (N2), which correspond to the a priori study of test N6 in Tab. 4.4. Mean profile of |l i | and |ξ i | (a), mean spectra at y + ≈ 17 (b, c) and at the center of the channel (d, e).l x is ---, l y is • • • • • • , l z is ----, ξ x is -• -, ξ y is -• • • -and ξ z is ---.

  since they are in competition to produce the resolved energy transfer by the resolved and subgrid scales. The comparison of the mean values of |l i | and |N i | are presented in Fig. 4.11. The mean |l i | values are approximately the same as the nonlinear terms |N i |.

Figure 4 . 7 :

 47 Figure 4.7: Velocity profile (a) and velocity RMS (b, c, d) of LES-Langevin and Smagorinsky models. All the curves are referenced in Tab. 4.5. DNS (N1) ---, CDNS (N2) • • • • • • , Smagorinsky with c s = 0.13 (N3) ----, LES-Langevin model with c s = 0.13 (N4) -• -, Smagorinsky with c s = 0.065 (N5) -• • • -, LES-Langevin model with c s = 0.065 (N6) ---.

Figure 4 . 8 :

 48 Figure 4.8: Streamwise velocity spectra near the turbulence maximum (y + ≈ 17, a, c) and at the center of the channel (y + ≈ 620, b, d). All the curves are referenced in Tab. 4.5. DNS (N1) ---, CDNS (N2) • • • • • • , Smagorinsky with c s = 0.13 (N3) ----, LES-Langevin model with c s = 0.13 (N4) -• -, Smagorinsky with c s = 0.065 (N5) -• • • -, LES-Langevin model with c s = 0.065 (N6) ---.

Figure 4 . 9 :

 49 Figure 4.9: Wall normal velocity spectra near the turbulence maximum (y + ≈ 17, a, c) and at the center of the channel (y + ≈ 620, b, d). All the curves are referenced in Tab. 4.5. DNS (N1) ---, CDNS (N2) • • • • • • , Smagorinsky with c s = 0.13 (N3) ----, LES-Langevin model with c s = 0.13 (N4) -• -, Smagorinsky with c s = 0.065 (N5) -• • • -, LES-Langevin model with c s = 0.065 (N6) ---.

Figure 4 . 10 :

 410 Figure 4.10: Spanwise velocity spectra near the turbulence maximum (y + ≈ 17, a, c) and at the center of the channel (y + ≈ 620, b, d). All the curves are referenced in Tab. 4.5. DNS (N1) ---, CDNS (N2) • • • • • • , Smagorinsky with c s = 0.13 (N3) ----, LES-Langevin model with c s = 0.13 (N4) -• -, Smagorinsky with c s = 0.065 (N5) -• • • -, LES-Langevin model with c s = 0.065 (N6) ---.

Figure 4 . 12 :Figure 4 . 13 :

 412413 Figure 4.12: The mean spectra of |l x | at y + ≈ 17 (a, c) and at the center of the channel, where y + ≈ 620 (b, d). All the curves are referenced in Tab. 4.5. A priori ---, LES-Langevin model with c s = 0.13 (N4) -• -, LES-Langevin model with c s = 0.065 (N6) ---.

Figure 4 . 14 :

 414 Figure 4.14: Mean velocity profile (a), Reynolds Stress (b, c, d) of the tests from Tab. 4.6. DNS (N1) is ---, LES-Langevin (N2) is ----, LES-Langevin (N3) is -• -, LES-Langevin (N4) is -• • • -.

Figure 4 . 15 :

 415 Figure 4.15: Velocity profile (a) and velocity RMS (b, c, d) of LES-Langevin model normalized by different u τ corresponding to each simulations. All the curves are referenced in Tab. 4.7. DNS (N1) ---, CDNS (N2) • • • • • • , LES-Langevin model with c s = 0.13 (N3) ----, LES-Langevin model with c s = 0.1 (N4) -• -, LES-Langevin model with c s = 0.065 (N3) -• • • -.

Figure 4 . 16 :Figure 4 . 17 :Figure 4 . 18 :

 416417418 Figure 4.16: Spectra of u near the turbulence maximum (y + ≈ 17, a, c) and at the center of the channel (y + ≈ 620, b, d). All the curves are referenced in Tab. 4.7. DNS (N1) ---, CDNS (N2) • • • • • • , LES-Langevin model with c s = 0.13 (N3) ----, LES-Langevin model with c s = 0.1 (N4) -• -, LES-Langevin model with c s = 0.065 (N3) -• • • -.

Figure 4 . 19 :

 419 Figure 4.19: Velocity profile (a) and velocity RMS (b, c, d) of DNS-LML data. Total scales: ---, (k x,max , k y,max , k z,max ) T = (64, 65, 64) T : • • • • • • , (k x,max , k y,max , k z,max ) T = (64, 49, 32) T : ----.

Figure 4 . 20 :

 420 Figure 4.20: Velocity profile (a) and velocity RMS (b, c, d) of LES-Langevin model normalized by different u τ corresponding to each simulations. All the curves are referenced in Tab. 4.8. DNS (N1) ---, CDNS (N2) • • • • • • , LES-Langevin model with c s = 0.1 (N3) ----.

  Some remarkable LES models . . . . . . . . . . . . . . . . . . . . . . . . . Similarity model . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.5 Mixed model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.6 The Planar-Averaged Scale Dependent Dynamic model . . . . . . . 1.2.7 Lagrangian Scale Dependent and Independent models . . . . . . . . 1.2.8 An Approximate Deconvolution Model (ADM) . . . . . . . . . . . . 1.2.9 Subgrid scale estimation model (SGEM) . . . . . . . . . . . . . . . 1.2.10 Implicit LES (ILES) in the case of the adaptive local deconvolution model (ALDM) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.11 Variational multiscale approach (VMS) . . . . . . . . . . . . . . . . 1.2.12 A resolvable subfilter scale (RSFS) model . . . . . . . . . . . . . . . 1.3 Comparison of the models . . . . . . . . . . . . . . . . . . . . . . . . . . . Motivation for a new model . . . . . . . . . . . . . . . . . . . . . . . . . . 2.1.1 Criticism of existing models . . . . . . . . . . . . . . . . . . . . . . 2.1.2 Strategy of LES-Langevin modeling . . . . . . . . . . . . . . . . . . 2.2 The LES-Langevin model of turbulence . . . . . . . . . . . . . . . . . . . . 2.2.1 Derivation of the Langevin equation . . . . . . . . . . . . . . . . . . 2.2.2 A friction term . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.2.3 Derivation of the SGS model . . . . . . . . . . . . . . . . . . . . . . 2.2.4 Isotropic turbulence case . . . . . . . . . . . . . . . . . . . . . . . . 2.3 Channel flow modeling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.3.1 Estimation of stochastic forcing . . . . . . . . . . . . . . . . . . . . 2.3.2 A priori tests . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Validation for turbulent flow . . . . . . . . . . . . . . . . . . . . . . . . . . Filtering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3.1 Filter application in the code . . . . . . . . . . . . . . . . . . . . . 3.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Preliminary remarks on the DNS data . . . . . . . . . . . . . . . . . . . . 4.2 Eddy-viscosity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2.1 Smagorinsky-Van Driest model . . . . . . . . . . . . . . . . . . . . 4.2.2 The negative viscosity model . . . . . . . . . . . . . . . . . . . . . . 4.3 The turbulent force . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.1 Inviscid case . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.2 Viscous models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

	CONTENTS	CONTENTS
	2 LRDT method
	2.1 4 Results
	4.1
	xi

to turbulence modeling 1 State of the art of Large Eddy Simulation (LES) 1.1 General properties . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.1 The LES idea . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.2 Leonard and double decomposition of the SGS stress tensor . . . . . 1.1.3 Explicit filtering formulation . . . . . . . . . . . . . . . . . . . . . . 1.1.4 Physical and mathematical constraints . . . . . . . . . . . . . . . . 1.1.5 Energy transfer . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 1.2.1 The Smagorinsky model . . . . . . . . . . . . . . . . . . . . . . . . 1.2.2 Smagorinsky model with Van Driest Damping . . . . . . . . . . . . 1.2.3 Dynamic model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.4 1.3.1 Energy transfer. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3.2 Channel flow. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.3.3 Dynamics of turbulent force l , a priori tests. . . . . . . . . . . . . 2.3.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 Numerical code 3.1 Geometry and equations . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1.1 Change of coordinates in wall normal direction . . . . . . . . . . . . 3.1.2 Discretization in wall-normal direction . . . . . . . . . . . . . . . . 3.1.3 Time discretization . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1.4 Boundary conditions . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1.5 Diagonalization method . . . . . . . . . . . . . . . . . . . . . . . . 3.1.6 Projection method . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1.7 Performances of the algorithm . . . . . . . . . . . . . . . . . . . . . 3.2 3.2.1 Comparison to a reference DNS at Re τ = 180 . . . . . . . . . . . . 3.2.2 Smagorinsky model statistics, Re τ =400 . . . . . . . . . . . . . . . 3.3 5 Conclusion and perspectives 5.1 Comparison of LES-Langevin to other LES models . . . . . . . . . . . . . 5.2 Further study of the LES-Langevin model . . . . . . . . . . . . . . . . . . 5.3 Perspectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xii CONTENTS CONTENTS 5.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112 Bibliography 112 xiii

  Schwarz inequality | a • b |≤ a b . If in n points x 1 , x 2 , x 3 , . . . the filter function G takes the values ω 1 , ω 2 , ω 3 , • • • ≥ 0, two n-dimensional vectors can be defined as

	.7)
	For non-negative filter functions G ≥ 0, the conditions (1.7) can be satisfied in 1d with the help of the Cauchy-

Table 2 . 1 :

 21 Examples of scale separation and probable resolutions of LES: k x m

	1 384 257 384	64	65 64
	2 384 257 384	64	49 32
	3 384 257 384	32	33 16
	4 160 129 192	64	49 32
	5 128 97 128	64	49 32
	6 128 97 128 discrete 49 32

Table 4 .

 4 1: Parameters of the DNS of reference. Corresponding figures are 4.1 and 4.19 Despite this difference, both simulations predict almost the same value of shear stress. The deviation between the shear velocities u τ is about 2%. The mean velocity profiles and fluctuating velocity RMS of Marquillie et al. and Alamo et al. are shown in Fig. 4.1.

	N	DNS	nx k x,max ny	nz	L x L y L z Re τ	u τ
	1 Marquillie 768	384	257 384 2π 2	π 618 0.0490
	2	Alamo	1024 1024 257 1024 8π 2 4π 547 0.0489

Table 4 .

 4 6: Parameters used for the tests in figures 4.13 and 4.14

			13 25	1	not correlated	256 97 64 64 49 32 590 -0.04425
	3 LES-L 4 LES-L	τ f τ f /2	0.13 25 0.13 25	1 10 -2	τ st f = 7.5 * 10 -2 256 97 64 64 49 32 584 -0.05410 not correlated 256 97 64 64 49 32 591 -0.04384

Table 4 .

 4 8: Parameters used for the tests in figure 4.20

	N Model A + c s	nx ny	nz k x k y k z Re τ	δ τ
	1	DNS	768 257 384	618	0
	2 CDNS	256 97 192 64 49 96 630 0.01954
	3 LES-L 25 0.1 256 97 192 64 49 96 601 -0.02670

T : ----.

Remerciements

¡ ¤R & %B S £ ¦£ ¦¥ U i" F¥ ¨! ¤¡ ¤R ¤¡ ¤A dB S q" F¥ 'P ia c$ 1£ Q(

can see, some of the models have |δ τ | > |δ τ (CDN S)|. They are obtained by decreasing the Van-Driest constant below a critical value of the order of A + = 25, for c s = 0.1. Also, increasing c s (c s = 0.05, c s = 0.1, c s = 0.17), one observes first a reduction of δ τ until c s = 0.1, then an increase. This means that the optimal parameter set for this eddy-viscosity model is c s = 0.1 and A + = 25, resulting in an error on the shear-stress of about 0.3% in comparison to the DNS reference value. This choice of parameter is in agreement with the choice found in the literature by different authors [START_REF] Germano | A dynamic subgrid-scale eddy viscosity model[END_REF]. We present in Fig. 4.2 (a) and Fig. 4. 3 (a) the comparison between the mean and the RMS profiles from the reference DNS and the different models from Tab. 4.2. We see that c s = 0 and c s = 0.05 produce better mean velocity profiles for y + ∈ [0, 8]. The increase of c s ≥ 0.1 contaminates the mean velocity profile in the viscous sublayer and buffer layer (y + ∈ [0, 30]). At larger wall distances (y + ≥ 30), the best mean velocity profiles are obtained with c s = 0.1 and c s = 0.17. The last value corresponds to the optimal isotropic turbulence modeling by the Smagorinsky approach [START_REF] Vreman | The adjoint filter operator in large-eddy simulation of turbuleny flow[END_REF]. The best mean velocity profile as a function of A + is obtained for A + = 25 and A + = 50, see Fig. ,c,d). All the models display an excess of turbulent energy in the center of the channel similar to the discrepancy observed between the LML-DNS and the Jimenez-DNS. Surprisingly, the best Reynolds stress is obtained for c s = 0.05. The standard Smagorinsky with c s = 0.1 and A + = 25 has too large a RMS maximum that is shifted to higher y + . Increasing c s leads to an undesirable shift of the turbulence maximum. Therefore, we cannot achieve a prefect model of both the RMS maximum and its center channel value with a single c s optimization.

N Model

On the other hand, increasing A + moves the fluctuating velocity RMS maximum to the wall. We can therefore use this degree of freedom to improve the location of the RMS maximum in the model. However, it appears that high values of A + result in a reduction of the modeled resolved-subgrid scale transfer. The results become therefore strongly dependent of the resolution in the wall-normal direction. In our case of Chebyshev collocation discretization, the wall turbulence is quasi resolved. Thus the use of large A + values can be less efficient for cruder discretization.

We see that the Samgorinsky-Van-Driest model alone is not sufficient to produce the correct level of the resolved-subgrid scales energy transfer at the right height y + . One which is impossible to achieve if the effect of subgrid scales is modeled only by negative eddy-viscosity. To compare the results of Härtel et al. [START_REF] Härtel | Subgrid-scale energy transfer in the near-wall region of turbulent flows[END_REF] and Dubrulle et al. [START_REF] Dubrulle | Eddy-viscosity of parity-invariant flow[END_REF], we need to compute the dissipation of the fluctuating velocity field, which was not represented by Härtel. Nevertheless, there is no contradiction between the statement that the negative eddy-viscosity is suitable for modeling of small rapid turbulent structures and the fact that the total energy transfer is positive.

The different tests are summarized in Tab. 4.3. They were performed on the model with the best shear stress (δ τ ≈ 0.3%), obtained with c s = 0.1 and A + = 25. One sees that a first effect of the negative turbulent viscosity is to reduce the shear stress as function of The mean velocity profiles are shown in Fig. 4.4 (a). One sees that indeed, the addition of a negative turbulent viscosity is able to change the mean velocity at the channel center. For the RMS profiles, one sees that it can produce a shift of the maximum towards the wall, as well as a slight increase of the maximum. However, it does not change significantly the RMS values at the channel center. It seems therefore that the discrepancy between the LES and the DNS RMS profiles at the channel center cannot be attributed to a negative eddy-viscosity effect, and still remains mysterious.

The turbulent force

Inviscid case

We now turn to test the other parameter of our model, namely the turbulent force given by the vector l modeling the Cross Stress tensor C ij . This term should intensify the forward energy transfer from the resolved to subgrid scales as well as the backscatter. To test its influence, we therefore first turn off the eddy-viscosity and put c s = 0. The different test performed are presented in Tab. 4.4. Model N3 computes the full evolution equation for the vector l . Models N4 and N6 test the "over-damped" approximation, where the vector l is adapted instantaneously to the stochastic force, with different friction constant. Finally, in model N5, we shut down the turbulence force, but change the filter shape. One sees that neither model N3, N4 or N6 result in a significant change in the parameter δ τ with respect to a coarse DNS. In contrast, changing the filter has a larger impact on the value of δ τ .

Viscous models

The comparison to the Smagorinsky model

As a first step, to test the influence of the vector l in dissipative models, we perform direct comparison of our LES-Langevin model with pure Smagorinsky models. The different tests are summarized in Tab. 4.5. For these values of the parameters, one sees that the difference between LES-Langevin and Smagorinsky models in the shear stress predictions for the same c s is about 20%.

Fig. 4.7 shows that the vector l slightly improves the mean velocity profile for c s = 0.13. A more important difference can be seen on the RMS velocity profiles. The presence of vector l increases the RMS values mostly in the range y + > 50. The streamwise velocity RMS of LES-Langevin model has a better value in the center region of the channel as well as the maximum value of the wall normal and spanwise velocity RMS. Also, the competition between the resolved energy production by l and its dissipation by the eddyviscosity is favorite for the last for relatively high value of the constant c s = 0.13. When c s = 0.65, the turbulent viscosity does not dissipates all the resolved turbulent energy produced by the vector l .

Velocity spectra of the tests of Tab. 4.5 are shown in Fig. 4.8, Fig. 4.9 and Fig. 4.10. We see the strong dependence of the smallest resolved scales on the eddy-viscosity constant. For c s = 0.13 the dissipation of the smallest resolved scales is widely over-predicted at y + ∼ 17 for streamwise spectra. The constant c s = 0.065 under-predicts the dissipation of the resolved scales. The figures show that the spectral feature of the velocity produced by the Smagorinsky and LES-Langevin model do not differ very much. The spectra are much more sensitive to the eddy-viscosity term than to the turbulent force l .

These results motivates us to compare the l i terms with the nonlinear N i term: