Project topics and materials
Get Complete Project Material File(s) Now! »
RANS solver implementation for turbulent applications.
Introduction
Since the early days of numerical computation, many numerical schemes and approaches have been proposed, continuously improving the precision and efficiency of solvers. Meanwhile, mod-elization followed the same trend, evolving from panel methods in the 1970’s to potential flows, inviscid flows and RANS models in the 1990’s. Since, RANS computations have become in the industry the baseline for design and verification and are expected to remain the norm for decades, despite the emergence of LES and DNS methods [Slotnick 2014]. During this process, finite volume methods has become the most popular discretization for flow solvers for di↵erent reasons. Though, with high order needs, we see today a regain in popularity in finite elements methods (continuous and discontinuous Galerkin).
Indeed, finite element discretization was considered for a long time to be inappropriate for convection dominated problems, so that it is now considered a common knowledge (be-lief) [Chitale 2014, Pirzadeh 1996, Park 2010, Frazza 2017, Sahni 2009]. Similarly, it is a com-mon believe that unstructured grids are improper for the computation of High Reynolds bound-ary layers. Many justifications have been proposed for this, like the inaccurate computation of gradients. In the facts, structured properties in the mesh significantly eases the math involved in the development of a numerical scheme, providing privileged directions in the datas to per-form finite di↵erence gradient reconstruction. Thus, such assumptions are the base of numerous schemes. These in turn perform nominally on structured grids and rather well on unstructured isotropic grids, but degenerate rapidly in presence of anisotropy without structured orthogonal elements. Hence, most of the standard flow solvers either crash or yield poor results when boundary layers are not discretized with a structured mesh [Sahni 2009]. This is generally due to an accumulation of inaccuracies in the solver, both in the residual computation and in the resolution. For example, an implicit solver with too large approximations in the linearization will be prone to divergence and be unstable.
The situation can be summed up by noting that there is somehow a constant balance be-tween the work required in the solver and the mesh generation. Easing the development of the flow solver displaces the work requirement on the mesh generation by putting strong constraints on the quality of the mesh required. Such constraints have been one of the major limitation preventing the spread of mesh adaptation methods as it is nearly impossible to generate meshes that follows both the requirement of the flow and the solver. We can cite as a common ex-ample the treatment of the mesh near a cut trailing edge of a wing [Aubry 2009]. In this first part, we show that modern strong solvers are able to compute RANS solutions on fully un-structured grids. These solvers have in common the control and minimization of unnecessary approximations in the resolution.
The flow solver Wolf was initially developed as a finite-volume inviscid compressible Euler flow solver. It was then a tool to assess the capability of mesh adaptation strategies to produce reliable meshes for CFD computations. Its whole development has thus been driven by mesh adaptation. It has to be able to run several time without crashing and without any external intervention on unstructured grids with some isolated bad elements. Wolf is now a mixed finite volume/finite element, vertex centered flow solver that handles unstationary and stationary flows on fixed grids and ALE simulations, inviscid and viscous laminar or RANS flows. In the scope of being able to perform adjoint-based mesh adaptation, the computation of the adjoint system for inviscid flows has been also implemented.
In this first part, we present the implementation of Wolf and the di↵erent parameters that will be used in our computations. We develop the various implementations that have been made to deal with the various turbulent flows simulations presented in this thesis. In particular, we introduce the choices and modifications that have been made to make it possible to perform RANS computations with Wolf on fully unstructured meshes. We also present the implemen-tation of the turbulent adjoint and consideration of periodic flows to deal with turbomachinery applications.
In this Chapter we develop the base implementation of our flow solver Wolf and the various choices that have been made to provide an efficient and robust resolution of the equations on fully unstructured grids. We detail the primal and adjoint solver that are used throughout this thesis and their validation. We also introduce the various implementations that have been made during this thesis to deal with various applications such as turbomachinery or low-Mach applications.
Main flow solver description
We solve the RANS equations using the Spalart-Allmaras one-equation model, which is a standard and well-documented option (k ! SST model is in development but not tested yet). The spatial discretization of the governing equations is based on a vertex-centered finite element-finite volume formulation, where the finite volume cells are built on unstructured meshes. Second-order space accuracy is achieved through a piecewise-linear extrapolation based on the Monotonic Upwind Scheme for Conservation Law (MUSCL) procedure with a particular edge-based formulation. Both explicit and implicit approaches are available for the time integration. During implicit simulations, a linearized system is solved at each solver iteration using an approach derived from the Lower-Upper Symmetric Gauss-Seidel (LU-SGS).
Modeling equations
The Reynolds Averaged Navier-Stokes (RANS) system relying on the Spalart-Allmaras model is composed of the compressible Navier-Stokes equations and the standard Spalart-Allmaras equation with no trip.
Turbulence modeling
The modelization of turbulence is described in more details in Section 2.1. We use the Spalart-Allmaras [Spalart 1992] or Menter k ! SST [Menter 1993] RANS turbulence models to compute an eddy viscosity µt and eddy thermal conductivity t that are added to the laminar ones. These relie on the Boussinesq approximation to model the turbulent stress tensor as proportional to the strain rate tensor. It can be eventually amended in the case of the SA QCR (see Section 2.1.2) with a quadratic constitutive relation to express the turbulent tensor with respect to the strain rate and vorticity tensors.
Spatial discretization
The spatial discretization of the fluid equations (1.1) is based on a vertex-centered finite ele-ment/finite volume formulation on unstructured meshes. It combines a HLLC approximate Rie-mann solver [Batten 1997] to compute the convective fluxes and the Galerkin centered method to evaluate the viscous terms. Second order space accuracy is achieved through a piecewise linear extrapolation based on the Monotonic Upwind Scheme for Conservation Law (MUSCL) procedure [Leer 1972] which uses a particular edge-based formulation with upwind elements. A specific slope limiter is employed to damp or eliminate spurious oscillations that may occur in the vicinity of discontinuities [Courn`ede 2006] (see Section 1.1.2).
Dual mesh construction
Let H be a mesh of domain⌦, the vertex-centered finite volume formulation consists in as-sociating with each vertex Pi of the mesh a control volume or finite volume cell, denoted Ci. Discretized domain⌦ h (see Figure 1.3) can be written as the union of the elements or the union of the finite volume.
where NK is the number of elements and NV the number of vertices.
Note that the dual mesh (composed of cells) is built in a preprocessing step. Consequently, only a simplicial mesh is needed in the input. Several choices are possible to build finite volume cells. In this work, two methods were considered: median cells and containment cells.
Median cells. In 2D, this standard method consists in building cells bounded by segments of medians (so-called median cells), see Figure 1.1. In 3D, each tetrahedron is split into four hex-ahedra (one associated with each one of its four vertices). The eight vertices of the hexahedron associated with a point Pi are given by: (i) Mi, Mj, Mk, the middle points of the three edges incident to Pi, (ii) Gfi, Gfj, Gfk, the gravity centers of the 3 faces containing Pi, (iii) G, the gravity center of the tetrahedra, and (iv) the vertex Pi considered. The cell Ci associated with vertex Pi is the union of all hexahedra of the tetrahedra surrounding Pi.
These cells enjoy a rather good robustness to distorted meshes, but they are not well-suited to stretched quasi-structured meshes. Traditional meshes for RANS computations exhibit structured boundary meshes where median cells are evenly distorted and can create spurious noise. To deal with such meshes the second method is used.
Containment cells. This method, introduced in 2D by Barth [Barth 1994] and generalized to 3D by [Gourvitch 2004], is supposed to be well-suited to discretize accurately the flow equations on highly anisotropic quasi-structured meshes (boundary layer meshes). In 2D, each triangle is split in three quadrangles, each one joining a vertex, the middle of each edge containing this vertex and the center of the containment circle of the triangle, see Figure 1.1. The center of the containment circle is the circumcenter if it is located inside the triangle, otherwise it is the center of the smallest circle containing the triangle (which, in fact, is the middle of the largest edge when the triangle has an obtuse angle). In 3D, it consists in subdividing each tetrahedron into four hexahedra cell around each vertex. The hexahedron cell vertices associated with vertex
Pi are (i) the middle of the three edges issued from Pi, (ii) the containment circle center of the three faces containing Pi, (iii) the containment sphere center of the tetrahedron and (iv) the vertex Pi considered. The containment cell of vertex Pi is the union of all of these hexahedra cells. The containment sphere center corresponds to the sphere circumcenter if it is located inside the element, otherwise it is the center of the smallest sphere containing the tetrahedra which can be the middle of the largest edge or the center of the circumsphere containing the ”largest” face.
Comparison of median and containment cells. Figure 1.3 compares the two approaches on a 2D mesh with a quasi-structured region. The same comparison in 3D is shown in Figure 1.4. These examples illustrate how the faces of containment cells are aligned with the flow direction in the presence of a boundary layer mesh, in contrast to median cells. This is why, containment cells are supposed to exhibit a better behaviour for highly stretched quasi-structured meshes.
Containment cells success and failures: Figure 1.3 shows the perfect cells obtained with containment cells on cartesian grids but such grids are seldom used. In practice containment cells are very efficient on fully structured meshes with very high aspect ratios, such as shown in Figure 1.5, where the cells are quasi perfect. Median cells are in that case, nonetheless less stable that containment cells (require lower CFL and an exact di↵erentiation of the viscous terms), but also yield a poor discretization of the equations and are more dissipatives.
Allamaras model, with HLLC solver, Fourth order numerical dissipation and the Piperno limiter. The linear system is inverted with SGS iterations with a precision of 0.1. We start with a uniform solution at CF L = 0.1, with a geometric growth of 1.02 and a maximal CF Lmax = 500.
Both seem globally similar except that with median cells (left), the shock is more di↵used and positioned aft of the shock predicted with containment cells (right). But a closer look at the turbulent variable ⌫˜ field shows that median cells predict less turbulence near the trailing edge as can be seen in Figure 1.7 (bottom). Meanwhile, Figure 1.7 also shows that the velocity field (bottom) is smoother, which leads to less turbulence production. This is due to the artificial dissipation introduced by median cells. This in turn increase the drag computed by 3%.
Following this constatation, containment cells seemed appropriate for pseudo structured meshes with a structured boundary layer and an unstructured discretization for the outer field where containment cells were expected to behave roughly like median cells. However, contain-ment cells tend to produce spurious noise and are more sensible to low-Mach on unstructured meshes, such as shown in Figure 1.8. The resulting density fields obtained on this mesh at Mach 0.1 and ↵ = 2 , with median and containment cells are shown in Figure 1.9. Both solutions are obtained with HLLC solver with no limiter and the residual is converged by 8 orders (1 ⇥10 8). The noise introduced by containment cells is clearly visible.
Moreover, in order to recover structured cells in 3D, the structured mesh has to be generated from a specific subdivision of an hexahedral mesh [Gourvitch 2004]. Any permutation in the subdivision of the hexahedra leads to bad cells. It thus yields a strong constraint on the mesh generation in 3D. Finally, in a mesh adaptation context, median cells proved to be sufficient and as we were able to discretize boundary layers with unstructured meshes, we abandoned containment cells.
Table of contents :
Introduction
I RANS solver implementation for turbulent applications
Introduction
1 Main solver description
1.1 Main flow solver description
1.1.1 Modeling equations
1.1.2 Spatial discretization
1.1.3 Boundary Conditions
1.1.4 Rotating Flow
1.1.5 Time integration
1.2 Linear System Resolution
1.2.1 Standard linear solvers
1.2.2 SGS – Line solver
1.2.3 Multigrid solver
1.2.4 Multigrid in Wolf:
1.2.5 Residual computation and limiter freezing
1.3 Adjoint
1.3.1 Functional sensitivity to the solution
1.3.2 Jacobian computation and system resolution
1.4 Verification and validation
1.4.1 Verification
1.4.2 Validation
2 RANS solver description
2.1 Turbulence modeling
2.1.1 Spalart-Allmaras Negative Model
2.1.2 Spalart-Allmaras Quadratic Constitutive Relation (SA-QCR)
2.2 Spatial discretization
2.2.1 Di↵erentiation
2.3 Wall Functions
2.3.1 RANS solution initialization
2.3.2 RANS adjoint
3 Periodic computations for turbomachinery applications
3.1 Implementation
3.2 Periodic Adjoint
3.3 Applications
3.3.1 LS89 case
3.3.2 RO37 case
4 Solution correction for error estimation
4.1 Non-linear Corrector computation
4.1.1 1D RANS Example
4.1.2 Formal discretization
4.1.3 Finite Volume Corrector
4.1.4 Finite Element Corrector
4.2 Numerical results: manufactured solutions
4.2.1 Finite elements 1D nonlinear advection di↵usion problem
4.2.2 Finite volume 2D Navier-Stokes equations
4.3 Numerical results: applied cases
4.3.1 2D Laminar NACA0012
4.3.2 2D Turbulent NACA0012
4.3.3 3D ONERA M6 wing turbulent Case
4.3.4 3D inviscid supersonic case: Sonic Boom Prediction Workshop
4.4 Conclusion and perspectives
Conclusion
II Anisotropic mesh adaptation for RANS applications
Introduction
5 Continous mesh framework for mesh adaptation
5.1 Continous mesh
5.1.1 A simple 1D example
5.2 Extension to higher dimensions
5.2.1 Euclidian and Riemannian metric spaces
5.2.2 Unit mesh
5.2.3 Operations on metrics
5.3 The continuous mesh framework
5.3.1 Duality discrete-continous: a new formalism
5.3.2 Continuous linear interpolation
5.3.3 Summary
5.4 Multiscale mesh adaptation
5.4.1 The non-linear adaptation loop
5.4.2 Optimal control of the interpolation error and optimal meshes
5.4.3 General error control
5.4.4 Error equidistribution
5.4.5 Application to numerical solutions
5.4.6 Control of the error in Lp norm
5.4.7 Conclusion
6 RANS mesh adaptation
6.1 RANS goal-oriented mesh adaptation
6.1.1 General non-linear error estimation
6.1.2 RANS goal-oriented error estimation: Linearization
6.1.3 RANS goal-oriented error estimation: Case of Spalart-Allmaras equations
6.1.4 Second derivatives error estimation
6.1.5 Hessian computation
6.2 Mesh adaptation cycle for applied cases
6.3 Applied cases: NASA Rotor 37
6.3.1 Geometry description
6.3.2 Mesh comparison
6.4 Applied cases: M6 wing
6.4.1 Case description
6.4.2 Feature-based mesh adaptation
6.4.3 Comparison of hessian reconstructions
6.4.4 Goal-oriented mesh adaptation
6.4.5 Comparison with other codes
6.5 Applied cases: Trap wing
6.5.1 Case description
6.5.2 Results
6.6 Applied cases: High-Lift Common Research Model
6.7 Improvement of RANS mesh adaptation
6.7.1 Adaptation with Wall-Laws: Subsonic NACA0012 validation
6.7.2 Hessian boundary correction: impact on a 2D High-Lift configuration
6.7.3 Conclusion
6.8 Conclusion
7 Periodic mesh adaptation
7.1 Periodic mesh adaptation
7.1.1 Handling the metric
7.1.2 Periodicity recovery
7.2 Application to turbomachinery
7.2.1 LS89 case
7.2.2 RO37 case
Conclusion
Conclusion and perspectives
Acknowledgments
Bibliography
