# About mesh adaptation for hybrid flow calculation

B Sauvage, F Miralles, S Wornom, Bruno Koobus, Alain Dervieux

## ► To cite this version:

B Sauvage, F Miralles, S Wornom, Bruno Koobus, Alain Dervieux. About mesh adaptation for hybrid flow calculation. ETMM14 - 14th International ERCOFTAC Symposium on Engineering Turbulence Modelling and Measurements, ERCOFTAC, Sep 2023, Barcelona, Spain. hal-04372058

HAL Id: hal-04372058

https://hal.science/hal-04372058v1

Submitted on 4 Jan 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# ABOUT MESH ADAPTATION FOR HYBRID FLOW
## CALCULATION

B. Sauvage¹,*, F. Miralles², S. Wornom², B. Koobus², A. Dervieux³

¹ Institut National de Recherche en Informatique et en Automatique, INRIA, France

² Université de Montpellier, France

³ Lemma, France and INRIA

bastien.sauvage@inria.fr

## Abstract

The purpose of our study is to propose the combination of a space-time anisotropic mesh adaptation with a novel hybrid RANS/LES modelization. We first discuss the design of an hybrid model able to address a sufficiently large class of flows around obstacles. A RANS model is enriched with an intermittency equation and associated with a LES model in some regions of the flow. Second, a strategy is proposed to adapt the mesh for a better capture of unsteady behaviors. The adaptation is metric based, space and time. The new methods are applied to flows around cylinders and airfoils.

**Keywords:** turbulence, RANS, DDES, hybrid RANS/LES, Variational Multi-scale

## 1 Introduction

We propose to address the increasing complexity of industrial flows by extending the predictivity of hybrid models by the introduction of intermittency together with making easier the research of efficient and converged space and time approximations by introducing adaptation algorithms.

## 2 Modelization:

The base ingredients of our hybrid turbulence strategies are:

* **Baseline RANS components:** two low Reynolds RANS models are used in our hybrid models, more specifically they are the $k-\epsilon$ model proposed in Goldberg et al. (1998) and the $k-R$ model recently introduced by Zhang et al. (2020). They have been chosen for their abilities to properly predict separated flows with adverse pressure gradients.

* **Transition model for boundary layer:** Let $P_k = \tau : \nabla u$ and $D_k = \rho \epsilon$ be the turbulent kinetic energy production and destruction terms. Let us define the laminar and turbulent flow regions in the boundary

layer:

$$L = \{ \mathbf{x} \in \Omega_f | Re_{\theta}(\mathbf{x}) < Re_{\theta,S}(\mathbf{x}) \}, \\ Tu = complement(L),$$

where $Re_{\theta}$ denotes the Reynolds number based on the boundary layer thickness which is defined by $Re_{\theta} = 0.664\sqrt{Re|\mathbf{x}|}$.

$$Re_{\theta,S} = 163 + \exp \left( F_{\lambda} - \frac{100F_{\lambda}}{6.91} \sqrt{\frac{2}{3}k} \right).$$

with

$$F_{\lambda} = \begin{cases} 6.91 + 2.48\lambda - 12.27\lambda^2, & \lambda > 0, \\ 6.91 + 2.48\lambda + 63.64\lambda^2, & \lambda < 0. \end{cases}$$

and

$$\lambda = 0.664^2 \frac{1}{Re} \frac{\partial u_x}{\partial x} (x - x_0)$$

where $x_0$ denotes the abscissa of the front body point. Following the work of Menter et al. (2015) and the works of Akhter et al. (2007, 2009, 2015), the present transition model is defined as follows :

$$\frac{\partial \rho k}{\partial t} + \nabla \cdot (\rho \mathbf{u} k) = \tilde{P}_k - \tilde{D}_k + \nabla \cdot [(\mu + \mu_t \sigma_k) \nabla k]$$

$$\frac{\partial \rho \epsilon}{\partial t} + \nabla \cdot (\rho \mathbf{u} \epsilon) = (c_{\epsilon_1} P_k - c_{\epsilon_2} D_k + E) T_t^{-1} + \nabla \cdot [(\mu + \mu_t \sigma_{\epsilon}) \nabla \epsilon]$$

$$\frac{\partial \rho \gamma}{\partial t} + \nabla \cdot (\rho \mathbf{u} \gamma) = c_{g_1} \gamma (1 - \gamma) \frac{P_k}{k} + \rho c_{g_2} \frac{k^2}{\epsilon} \nabla \gamma \cdot \nabla \gamma + \nabla \cdot [\sigma_{\gamma} (\mu + \mu_t) \nabla \gamma].$$

The intermittency model interacts with the turbulence model by modifying the turbulent kinetic energy equation. The new production and destruction terms are defined by :

$$\tilde{P}_k = \begin{cases} 0 & \text{if } \mathbf{x} \in L, \\ \max(\gamma, \gamma_1) P_k & \text{otherwise.} \end{cases}$$

$$\tilde{D}_k = \begin{cases} 0 & \text{if } \mathbf{x} \in L, \\ \max(\gamma, \gamma_2) D_k & \text{otherwise.} \end{cases}$$
<figure><img src="image_2.png" /><figcaption>Figure 1: Time splitting of the Transient Fixed Point mesh adaptation algorithm. Sub-intervals (in green) used for the transient process and timesteps (in red).</figcaption></figure>

where the model constants are defined by $\gamma_1 = 0$ and $\gamma_2 = 0.1$. From the above equations, one can notice that the baseline RANS model is recovered for an intermittency value $\gamma = 1$ (fully turbulent mode). A zero normal flux is also imposed on $\gamma$ at the wall.

* **DDES component:** in this work, the classical DDES approach, Spalart et al. (2006), is based either on the Spalart-Allmaras model or the $k-\epsilon$ model of Goldberg et al. (1998).

* **LES-like component:** the DVMS approach proposed in Moussaed et al. (2014) is used as the LES part of our hybrid models. In this approach, the variational multiscale (VMS) model, aiming to limit the effects of the subgrid-scale (SGS) model to the smallest resolved scales, is combined with the dynamic procedure which provides a tuning of the SGS dissipation in space and time, so that the resulting DVMS model enjoys synergistic effects.

* **Hybrid models:** Our hybrid strategies blend either a RANS or DDES model with the DVMS approach.

## 3 Mesh adaptative discretization :

The baseline discretization is an upwind, vertex centered, finite-volume approximation on tetrahedra. Diffusion terms are accounted for with the $P_1$-Galerkin approximation. The numerical diffusion involved in upwinding is made of sixth-order derivatives of primitive variables.

The capture of high Reynolds number flows with a mesh adaptation algorithm demands an algorithm able to find the boundary layer when starting from a uniform flow. This is obtained by using (1) an adequate adaptation criterion and (2) a stable and convergent fixed point for the coupling of flow and mesh. Main principles for a successful adaptation can be found in Dervieux et al. (2022). The transient fixed point (TFP) was introduced and discussed in Alauzet et al. (2002) and Alauzet et al. (2007). The simulation time frame $[0, T]$ is split into several subintervals (Figure 1):

$$[0, T] = [0 = t_0, t_1] \cup \dots \cup [t_i, t_{i+1}] \cup \dots \cup [t_{n-1}, t_n = T].$$

The TFP is extended as in Sauvage et al. (2023) by defining a space-time continuous mesh ($M, \tau, nstep$) as the knowledge of the following ingredients:

- a number $nstep$ of time intervals.

- a timestep function: $t \in [0, T] \mapsto (\tau(t)) \in [0, T]$ valid in the sense that $\int_0^T (\tau(t))^{-1} dt = nstep$.

- for every $t \in [0, T]$ a spatial metric $M(t)$, constant on each subinterval, of spatial complexity $n(t) = C_{sp}(M(t))$. The space-time complexity $C(M, \tau, nstep)$ of a space-time continuous mesh ($M, \tau, nstep$) is:

$$C(M, \tau, nstep) = \int_0^T C_{sp}(M(t))(\tau(t))^{-1} dt. \quad (1)$$

To any space-time continuous mesh ($M, \tau, nstep$) corresponds an estimate of the resulting error of a calculation working on the space-time mesh build with it. With a BDF2 time advancing scheme and a second-order spatial approximation, this estimate can be (in practice we choose $p = 4$):

$$\mathcal{E}_0(M_h, \tau) = \int_0^T \int_{\Omega} \left[ (\tau^2 \left| \frac{\partial^3 u}{\partial t^3} \right| + \Delta x H_u \Delta x) \right]^p dt dx. \quad (2)$$

The sensor $u$ is in our examples taken as the local Mach number of the computed flow. Minimizing the error $\mathcal{E}_0(M_h, \tau)$ with respect to ($M_h, \tau$) under the constraint of a fixed space-time complexity $C(M, \tau, nstep) = N_{prescribed}$ can be done analytically, allowing the building of a new space-time mesh. The process is repeated in the TFP loop.

In the approach chosen here, we restrict to LES flows for which a global shedding period $T_s$ can be estimated, and we use $T_s$ as the size of the subintervals. Due to quasi-periodicity, the different subinterval meshes are very similar, and we reduce them to an unique mesh in which quasi-steady features like boundary layers are accurately followed by the mesh, while unsteady vortices travel in a rather uniformly refined region.

## 4 Impact of the intermittency model:

The new intermittency model has been tested on the flow past a cylinder at Reynolds number around the drag crisis. At Reynolds 3900, the physical boundary layer is laminar, the intermittency model inhibits the statistical model and numerically predicts the laminar boundary layer (Figure 2). The quality of the prediction is of same level as with a pure LES calculation, e.g. by Moussaed et al. (2014). At Reynolds 380K, the physical boundary layer presents a dissymmetry between top and bottom, one layer being laminar while the other one is turbulent. This is well predicted with the help of the intermittency model, Figure 3. At Reynolds 1M the physical boundary layers are both turbulent and are again well predicted by the proposed model, Figure 4. The improvement carried in drag prediction by the intermittency equation is put in evidence in Figure 5.

The new intermittency model is now applied to the flow around a 3D wing based on the NACA0018 air-
<figure><img src="image_3.png" /><figcaption>Figure 2: Flow around a cylinder: impact of the intermit-
tency model on the drag crisis prediction, vorticity
at Rey=3900</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 3: Flow around a cylinder: impact of the intermit-
tency model on the drag crisis prediction, vorticity
at Rey=380K</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Figure 4: Flow around a cylinder: impact of the intermit-
tency model on the drag crisis prediction, vorticity
at Rey=1M</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 5: Flow around a cylinder: impact of the intermit-
tency model (in red) on the drag crisis prediction,
in contrast with the same base model, without in-
termittency (blue).</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Figure 6: Flow around NACA0018 wing, α = 0 degree.</figcaption></figure>

foil. At zero angle of attack, the slight vortex shed-
ding is captured, Figure 6. For the intermediate α =
6 degrees case, the unstable extrados boundary layer
is also predicted (Figure 7). The α = 15 degrees case
is depicted in Figure 8 and global behavior of the pre-
diction is in accordance with experiments of Nakano et
al. (2007), of Du (2016), and of Boutilier et al. (2012).

5 Mesh adaptative computations

We essentially present experiments performed with
2D test cases.

• Flow around a cylinder

We present 2D computations of a flow around a
cylinder at Reynolds number 3900 with the Spalart-
Allmaras turbulence model. Mesh adaptation options
are :

- only one time interval, therefore only one adapted
  spatial mesh.

- the Space-Time complexity $N_{st}$ is prescribed succes-
<figure><img src="image_8.png" /><figcaption>Figure 7: Flow around NACA0018 wing, $\alpha = 6$ degrees.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Figure 8: Flow around NACA0018 wing, $\alpha = 15$ degrees.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 9: Cylinder at Reynolds number 3900. Example of adapted mesh.</figcaption></figure>

sively to 10M, 20M, 100M, 200M.

An example of adapted mesh is given in Figure 9.

Space-time statistics are presented in Table 1. We ob-

<table><caption>Table 1: Space-time statistics for circular cylinder case at Reynolds number 3900. k is the number of fixed points, C<sub>sp</sub> holds for the spatial complexity (≈ number of space nodes), # Δt for the number of time steps, E<sub>sp</sub> and E<sub>time</sub> for the resulting space and time error integrals.</caption><thead><tr><th>N<sub>st</sub></th><th>k</th><th>C<sub>sp</sub></th><th># Δt</th><th>E<sub>sp</sub></th><th>E<sub>time</sub></th></tr></thead><tbody><tr><td>10M</td><td>1</td><td>27K</td><td>361</td><td>5.2 10<sup>-3</sup></td><td>3.1 10<sup>-5</sup></td></tr><tr><td>10M</td><td>5</td><td>60K</td><td>166</td><td>2.9 10<sup>-3</sup></td><td>3.1 10<sup>-3</sup></td></tr><tr><td>20M</td><td>1</td><td>44K</td><td>456</td><td>3.3 10<sup>-3</sup></td><td>3.1 10<sup>-5</sup></td></tr><tr><td>20M</td><td>5</td><td>84K</td><td>237</td><td>2. 10<sup>-3</sup></td><td>8.2 10<sup>-4</sup></td></tr><tr><td>100M</td><td>1</td><td>128K</td><td>780</td><td>1. 10<sup>-3</sup></td><td>3.1 10<sup>-5</sup></td></tr><tr><td>100M</td><td>5</td><td>173K</td><td>575</td><td>9.4 10<sup>-4</sup></td><td>8.8 10<sup>-5</sup></td></tr><tr><td>200M</td><td>1</td><td>203K</td><td>982</td><td>7.1 10<sup>-4</sup></td><td>3.1 10<sup>-5</sup></td></tr><tr><td>200M</td><td>5</td><td>225K</td><td>887</td><td>7.1 10<sup>-4</sup></td><td>4. 10<sup>-5</sup></td></tr></tbody></table>

serve that space error and time error tend to equilibrate, but yet imperfectly.

### Flow around an airfoil

The second case is the flow around a NACA0021 at Reynolds number 270K and an angle of attack 60 degrees, again with the Spalart-Allmaras turbulence model. Mesh adaptation options are:

- only one adapted spatial mesh,

- Space-time complexity is prescribed to 10M, 20M, 100M, 200M.

An example of adapted mesh is presented in Figure 10. Space-time statistics are presented in Table 2.

In Figure 11, the adaptive time step (function of physical time) produced by the space-time optimization is depicted, showing the increase of time-step, reducing computational cost while preserving the global accuracy.
<figure><img src="image_11.png" /><figcaption>Figure 10: Naca0021 at high angle of attack. Example of the adapted mesh.</figcaption></figure>

<table><caption>Table 2: Space-time statistics for NACA0021 at Reynolds number 270K.</caption><thead><tr><th>N<sub>st</sub></th><th>k</th><th>C<sub>sp</sub></th><th># Δt</th><th>E<sub>sp</sub></th><th>E<sub>time</sub></th></tr></thead><tbody><tr><td>2M</td><td>1</td><td>1766</td><td>1132</td><td>1.7 10<sup>-1</sup></td><td>2.4 10<sup>-6</sup></td></tr><tr><td>2M</td><td>10</td><td>16K</td><td>123</td><td>2.6 10<sup>-2</sup></td><td>1.3 10<sup>-2</sup></td></tr><tr><td>4M</td><td>1</td><td>2804</td><td>1426</td><td>1. 10<sup>-1</sup></td><td>2.4 10<sup>-6</sup></td></tr><tr><td>4M</td><td>10</td><td>23K</td><td>172</td><td>1.9 10<sup>-2</sup></td><td>4.6 10<sup>-2</sup></td></tr><tr><td>8M</td><td>1</td><td>4451</td><td>1797</td><td>6.8 10<sup>-2</sup></td><td>2.4 10<sup>-6</sup></td></tr><tr><td>8M</td><td>10</td><td>31K</td><td>261</td><td>1.4 10<sup>-2</sup></td><td>6. 10<sup>-3</sup></td></tr><tr><td>16M</td><td>1</td><td>7066</td><td>2264</td><td>4.3 10<sup>-2</sup></td><td>2.4 10<sup>-6</sup></td></tr><tr><td>16M</td><td>10</td><td>42K</td><td>384</td><td>1.1 10<sup>-2</sup></td><td>2.3 10<sup>-3</sup></td></tr><tr><td>32M</td><td>1</td><td>11K</td><td>2852</td><td>2.7 10<sup>-2</sup></td><td>2.4 10<sup>-6</sup></td></tr><tr><td>32M</td><td>20</td><td>57K</td><td>558</td><td>7.5 10<sup>-3</sup></td><td>1.8 10<sup>-3</sup></td></tr></tbody></table>

<figure><img src="image_12.png" /><figcaption>Figure 11: Bottom: Timestep lengths of initial flow at CFL=50, 1132 timesteps on 1766 vertices, and top: first timestep lengths proposed by the adaptation algorithm, 123 timesteps on 16K vertices.</figcaption></figure>

In Figure 12 we present the evolution (with Transient Fixed point Iterations) of the two components of the error, converging to each other.

<figure><img src="image_13.png" /><figcaption>Figure 13: Flow around NACA0018 at 15 degrees angle of attack. LES-adapted mesh with an instantaneous flow obtained with it (velocity magnitude).</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Figure 14: Flow around NACA0018 at 15 degrees angle of attack. LES-adapted mesh with an instantaneous flow obtained with it (velocity magnitude).</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>Figure 12: Evolution of theoretical space error $E_{sp}$ (top) and time error $E_{time}$ (bottom) with TFP iterations.</figcaption></figure>

The novel mesh adaptation approach has been tested on the 3D unsteady flow around a NACA0018 at an angle of attack of 8 degrees. The VMS model is applied and produces a rather unsteady flow for which the adaptation is able to both concentrate anisotropically on boundary layers and manage a rather uniform for large structures propagation, Figures 13,14.

The method is also being applied to a MRF flow (Chargy and Sauvage, 2022) around a Caradonna-Tung wing rotating above a Robin helicopter shape, Figures 15,16.

## Acknowledgments

This work was supported by the ANR NORMA
<figure><img src="image_16.png" /><figcaption>Figure 15: Flow around a Caradonna-Tung wing rotating above a Robin helicopter. LES-adapted mesh.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Figure 16: Flow around a Caradonna-Tung wing rotating above a Robin helicopter. An instantaneous flow obtained with the adapted mesh (velocity magnitude).</figcaption></figure>

project, grant ANR-19-CE40-0020-01 of the French National Research Agency. The authors gratefully acknowledge GENCI for granted access to HPC resources through IDRIS (grant 2022-A0132A05067) and CINES (grants 2021-A0102A06386 and 2020-A0092A05067).

## References

Akhter, M.N. and Funazaki, K.I., 2007, Development of Prediction Method of Boundary Layer Bypass Transition using Intermittency Transport Equation, *International Journal of Gas Turbine, Propulsion and Power Systems* 1, 30–37.

Akhter, M.N., Yamada, K. and Funazaki, K.I., 2009, Numerical Simulation of Bypass Transition by the Approach of Intermittency Transport Equation, *Journal of Fluid Science and Technology* 4, 524-535.

Akhter, M.N., M. Ali, M.K. and Funazaki, K.I., 2015, Numerical simulation of heat transfer coefficient on turbine blade using intermittency factor equation. *Procedia Engineering* 105, 495-503, The 6th BSME International Conference on Thermal Engineering.

Alauzet, F., Frey, P.J. and Mohammadi, B., 2002, Adaptation de maillages non structurés pour des problèmes installationnaires. *CRAS Ser. I*, 303:773-778. in French

Alauzet, F., Frey, P.J., George, P.L. and Mohammadi, B., 2007, 3D transient fixed point mesh adaptation for time-dependent problems: Application to CFD simulations. *J. Comp. Phys.*, 222:592–623.

Boutilier, M. S. H. and Yarusevych, S, 2012, Effects of end plates and blockage on low-reynolds-number flows over airfoils. *AIAA Journal*, 50(7):1547–1559.

Chargy, D. and Sauvage, B., 2022, A mesh adaptive

method for rotating machines, INRIA Research Report 9471, https://inria.hal.science/hal-03684715.

Dervieux, A., Alauzet, F., Loseille, A. and Koobus, B., 2022, *Mesh adaptation for Computational Fluid Dynamics*, t.1. ISTE Ltd and John Wiley & Sons, New York (ISBN: 9-781-78630-832-0), 1st edition.

Du, L., 2016, *Numerical and experimental investigations of Darrieus wind turbine start-up and operation*. PhD thesis. Durham University.

Goldberg, U., Peroomian, O. and Chakravarthy, S.,1998, A wall-distance-free $k-\epsilon$ model with Enhanced Near-Wall Treatment, *Journal of Fluids Engineering*, **120**, 457–462.

Moussaed, C., Wornom, S., Salvetti, M.V., Koobus, B. and Dervieux, A., 2014, Impact of dynamic subgrid-scale modeling in variational multiscale large-eddy simulation of bluff body flows, *Acta Mechanica*, **225**, 3309–3323.

Moussaed, C., Salvetti, M.V., Wornom S., Koobus, B. and Dervieux, A., 2014, Simulation of the flow past a circular cylinder in the supercritical regime by blending RANS and variational-multiscale LES models, *Journal of Fluids and Structures*, **47**, 114–123.

Nakano, T., Fujisawa, N., Oguma, Y., Takagi, Y., and Lee, S., 2007, Experimental study on flow and noise characteristics of NACA0018 Aifoil. *Journal of Wind Engineering and Industrial Aerodynamics*, 95(7):511–531.

Sauvage, S., Alauzet, F. and Dervieux, A., 2023, A space-time mesh Transient Fixed Point for mesh adaptation. *INRIA Research Report*.

Spalart, P., Deck, S., Shur, M., Squires, S., Strelets, M. and Travin, A., 2006, A new version of detached-eddy simulation, resistant to ambiguous grid densities. *Theory and Computational Fluid Dynamics*, 20:181–195.

Zhang, Y., Rahman, M.M. and Chen, G., 2020, Development of k-R turbulence model for wall-bounded flows, *Aerospace Science and Technology*, **98**.