3.2 Numerical aerodynamics
3.2.1 Meshing Techniques
Grid generation is a crucial problem for the computation of complex aircraft
configurations using a body fitted structured code. Furthermore, due to the data
management of structured grids, the local refinements around the geometry and
special flow regions (boundary layers, stagnation lines, wakes) tend to spread
through the whole domain even in zones where gradients are expected to be weak.
This can lead to very large grids, especially for complex geometries. Two approaches can reduce this drawback:
-
i) conservative non coincident interface boundary condition;
- ii) automatic mesh refinement.
CERFACS develops and maintains those techniques to help Airbus in reducing
simulation turn-around times and in improving the flow accuracy. Those two meshing techniques have been implemented and validated in the elsA
software.
Conservative non coincident interface boundary condition (M. Montagnac)
The purpose is to develop an efficient way to simplify the grid generation and
to deal with complex configurations using moderately sized grids.
In that technique also called patched grid technique, two domains must have a
common interface, or adjacent interface, but grid points of both interfaces do
not have to be at the same location or coincident. Grid lines through the
interface may be not continuous.
This approach prevents mesh points from spreading from a block to others.
Figure 3.2: Meshes with non coincident interface boundary
conditions for a RAE 2822 configuration (left) and mach number isolines
(right).
A more restrictive technique than the previous one is the conservative
near-matching interface boundary condition. Two domains must also have an
adjacent interface but grid points of both interfaces must match in a periodic
way.
Figure 3.3: Meshes with non coincident interface boundary conditions for a nozzle configuration (left) and
mach number isolines (right).
In the elsA software, those two methods are implemented differently for CPU
performance and memory requirement issues and they are both fully operational in
sequential and parallel modes.
Those techniques are compatible with all spatial numerical schemes, all
time marching methods, multigrid acceleration convergence techniques and
all turbulence models. Important efforts have been carried out in order to
maintain good CPU performances on both scalar and vector computers as well as in
a parallel mode.
Figure 3.4: Two sliding O-meshes of a pitching NACA0012 airfoil (left) and
the lift coefficient versus the angle of attack (right).
Figure 3.5: Sliding meshes of an axisymmetric four-block configuration for a laminar Taylor-Couette flow (left) and
time-independent Taylor vortices (right).
Fig. 3.2 shows the 2D turbulent transonic flow around the RAE
2822 airfoil computed with the RANS k-w turbulence model at
 = 0.73
and 2.79 degree angle of attack. The Euler flow across a nozzle
configuration is shown on Fig. 3.3.
The non coincident interface boundary condition is the core of the sliding mesh
option. This functionality is important in turbo machine activities, in
rotor-stator interactions or in aerodynamics around advanced high-speed
propellers for aeroelastic analysis. The actuator disk boundary condition often
used to model a propeller can then be replaced by the mesh of a propeller itself.
The first test case is a 2D transonic Euler flow around a pitching NACA0012
airfoil at
= 0.796
Fig. 3.4 shows two O-meshes
with the inner one oscillating in a harmonic motion with an amplitude of 1.01
degree and a frequency of 0.202. It also shows the limit cycle in the
evolution of the Cz aerodynamic coefficient that appears after a short
transitional period.
The second test case is a laminar Taylor-Couette configuration at Re=80 with
an internal moving cylinder and an outer immobile cylinder. The configuration is
decomposed in four meshes with two domains for each cylinder. One domain of each
is shown on Fig. 3.5 as well as the four time-independent Taylor
vortices splitted by the sliding boundaries.
Figure 3.6: AS28 wing, level 2, mesh and mach number contours.
Figure 3.7: Hierarchy of 3 grids for a RAE 2822 configuration (left) and
mach number contours (right).
Adaptive Mesh Refinement (M. Montagnac, J.-C. Jouhaud)
Another approach to handle complex configurations with a reasonable number of
grid points is to enrich a basic coarse mesh by means of a hierarchical
structure of grids, that is to say to use an Adaptive Mesh Refinement (AMR)
strategy
([Benkenida, 2002].
Areas of interest (shocks, boundary layers,
wake vortices,...) are dectected using physical sensors and
are refined to build new blocks. Then,
multigrid techniques can be applied on these grids.
The mesh adaptation is currently performed outside the flow solver with the
Airbus GAME software. Extension of this method to turbulent flows has been introduced in the
elsA software with the development of the local multigrid technique
([Jouhaud, 2002].
The AS28 wing (Fig. 3.6) in a transonic Euler flow at
= 0.80
and 2.2 degree angle of attack and the RAE 2822 airfoil
case (Fig. 3.7) in a turbulent transonic configuration with the
RANS k-w turbulence model at
= 0.73
and 2.79 degree angle of
attack, have demonstrated the efficiency of that AMR technique [1].
[1] J.-C. Jouhaud, M. Montagnac and L. Tourette. Multigrid adaptive mesh refinement for structured meshes for 3D aerodynamic design. Submitted in International Journal of Numerical Methods in Fluids, 2003.
3.2.2 Aerodynamic Data Models
The so-called ``Données Aérodynamiques'' (DA) characterize the global aerodynamic behavior of an aircraft in all its flight envelope.
These data are essential to qualify Performances (flight dynamics, aircraft dimensioning (sploilers, high-lift devices, jacks ...)), Flight Qualities (low and high speeds) and Loads (particularly for airframe integration) of aircrafts.
The current competing environment imposes to supply these data as soon as possible in the design cycle, as often as possible, with an increased precision and for a controlled cost.
The improvement in precision requires the use of simulation tools by engineers implied in the production of DA.
The exploration of a rather broad flight envelope led to a significant volume of calculations whereas the effort related to wind tunnel or flight testing tends to be reduced.
To preserve acceptable production costs of DA, it is necessary to reduce the restoration time of CFD calculations while improving their precision.
One of the aims of this project is consequently to enrich the "Model" engineers toolbox by adding efficient and robust complex CFD methods.
These aspects of speed and precision of numerical simulations are addressed during the preliminary work (first phase) of this program between AIRBUS FRANCE, ONERA.
CERFACS, within the framework of its commitment in the elsA software development, has joined this step while bringing in its experience and its know-how.
CERFACS's commitment in the first phase of this program (from january 2003 to december 2004) relates to generic aspects, primarily dedicated to the improvement of the elsA software
[Champagneux, 2003;
Champagneux, 2003a].
Three working directions are considered :
-
Algorithmic performances (improvement of the FAS (Full Approximation Scheme) multigrid algorithm within the framework of wall functions usage)
- Precision (verification and validation of wall functions, low speed preconditioning, extension of the AMR strategy and patch grids connectivities to wall functions)
- Automation and Pre/Post processing (complete polar calculations, convergence analysis)
Among these directions, the following tasks are already completed or still under progress.
Figure 3.8: AS28G aerodynamic coefficients convergence
Anisotropic multigrid (J.-F. Boussuge)
For steady CFD problems, multigrid algorithms are some of the fastest and most
adaptable ways of solving a system of partial differential equations on a fixed grid.
elsA uses the Full Aproximation Storage multigrid scheme which has been
developped for a full coarsening strategy (this means that the mesh is coarsened isotropically).
It implies, that for a given multiblock configuration, the global multigrid level which
drives the convergence is fully determined by the dimension of the most constraining boundary conditions.
Moreover, for anisotropic problems as including notably external viscous flows, the mesh used to
accurately resolve the physics presents an anisotropic refinement. An isotropic coarsening
of such a mesh can deteriorate the convergence rate of multigrid.
To overcome the problems quoted above, multigrid methods are much more efficient
if coarsening is done anisotropically. This permits to increase the global multigrid level
for configurations with different multigrid block levels and to reduce the cell aspect ratio
near the wall and the associated numerical stiffness for anisotropic refinement.
CERFACS has extended the standard FAS to anisotropic formulation and this development
has been successfully validated on two industrial test cases given by AIRBUS
(see Fig. 3.8 for the AS28G configuration).
Wall functions (N. Denève, G. Puigt)
As part of the effort for time restoration's reduction of high Reynolds number RANS simulations, wall functions appear as an unavoidable way to assess. Indeed, such wall functions both participate to the problem's size reduction since their use implies a lower resolution effort and also to the problem's stiffness reduction because of the replacement of the no-slip condition at walls. Moreover, their use sometimes alleviates the grid dependency of turbulence models.
Wall functions are based on the statement that, provided that certain assumptions are verified among which the flow is attached, the internal zone's structure of the boundary layer is independant of the flow in its external zone. A self-similarity velocity profile can be exhibited for the logarithmic region and the linear viscous sub-layer.
The implementation of a specific kind of wall functions
[Champagneux, 2003]
developped in the framework of this project is based on the approximation that the distance between the real wall and the location where boundary conditions should be applied is neglected. This approach differs from the initial wall functions of the elsA software based on an apriori agglomeration of near-wall cells of a low Reynolds grid that include the linear and logarithmic zones. Advantages of this method consist of grid generation and reusability of all available algorithms [1].
To evaluate this implementation, some computations have been realized on the RAE2822 test-case. Fig. 3.9 presents a comparison between the pressure coefficient obtained with a low-Reynolds approach (use of a no-slip boundary condition and of an adapted mesh : y+
1., 89010 nodes) and the one obtained with the implemented wall functions approach (on a coarser mesh : 55890 nodes and y+
[30,130] ). The results are shown to match correctly.
[1] S. Champagneux, J.-F. Boussuge, N. Denève, G. Hanss. Contribution au DTP Modèles de Données Aérodynamiques - Rapport d'avancement à un an. Contract Report, CR/CFD/04/15, 2004.
Figure 3.9: RAE2822 test-case : iso-Mach number field (left) and pressure coefficient (right)
Error estimations (G. Hanss)
The quality of a numerical solution with respect to the mesh on which it has been obtained can be assessed by an a posteriori error estimation. This post-processing analysis give us a mean to quantify the incertainty of the solution, and could provide us with a criteria for a local mesh adaptation process. An error estimate is expected to give information about the accuracy of the method, therefore it should estimate the absolute error level as well as the distribution of the error throughout the domain.
Two estimators have been implemented in elsA: the first one (Richardson extrapolation) is based on the analysis of the numerical solution in terms of the Taylor series expansion. An approximation of the leading term in the truncation error is deduced from solutions on two meshes with different cell sizes. Hence this method naturally couples with the use of multigrid acceleration techniques, where solutions on grids with different cell size are already available. One drawbacks is the fact that Richardson extrapolation assumes that spacing in each direction scales with h. This is valid only for geometrically similar control volumes. If the control volumes of two meshes are not similar in shape, the extrapolation accuracy decreases.
The second estimator implemented in elsA is a Residual error estimate. This estimator measures the disprepancy between the model equations to be solved and the discretized ones.
The following example presents the typical result on a 3D configuration. Fig. 3.10 presents isosurfaces corresponding to 80% of the maximum error throughout the field for the ONERA M6 Wing in a transonic, turbulent flow
( = 0.84, angle of attack  = 3.06º
As expected, the maximum errors are found in the lambda shock, leading and trailing edges and the wing tip.
Both estimators have proven a capability to predict the region of maximum error either by using the L2 and
norms or by using the field extraction. They can also highlight problems in the flow like bad boundary conditions or mesh irregularities. These capabilities could be use for an a posteriori analysis of the simulation but also as criteria in a mesh adaptation process as proposed in further studies.
Figure 3.10: Error estimation : isosurfaces corresponding to 80% of the maximum error
Automatic polar curve (G. Hanss)
The evaluation of polar data (Pressure vs Incidence for example) still constitutes a significant investment in terms of resources and know-how. Indeed, because of its non-linear nature, it is necessary to optimize the number of calculations. Generally, one evaluates three points of polar in the linear zone and approximately five in the non-linear one. Today, the total cost of a polar is of the order of the number of points to evaluate because the user computes several independant steady simulations for each flight point.
To minimize this cost, an automatic sequence of stationary calculations has been implemented. Each simulation is initialized with the solution of the previous flight point. The user is offered the possibility to choose the number of flight points and the condition of these points (angle of attack, slip angle, mach number, ...).
The solution proposed here presents several advantages:
-
each simulation of a flight point is initialized with the results from the nearest flight point, which should enable the simulation to converge more quickly,
- since it is fully implemented in the core source and not in the Python interface, it is possible to fully use the parallelism already implemented in elsA,
- on the output side, elsA provides a polar curve with each flight point, and a polar curve caracterizing the evolution during the steady calculation (Fig. 3.11), thus allowing the user "to reconsider" flight points which did not converge perfectly.
Figure 3.11: Automatic polar : RAE Profile
Figure 3.12: Comparison of the first harmonic of Cp given by elsA unsteady simulation and experimental data
3.2.3 Fluid-structure interaction (J. Delbove, G. Chevalier)
The field of aeroelasticity encompasses the interaction of aerodynamic, elastic, and inertia
forces acting on a flight vehicle. It had become recognized as an important
part of the aircraft design process. CFD emerged as a
practical technique for numerically solving the partial differential equations of fluid
flow in the compressible flow regimes.
The first area of interest is the static fluid-structure interaction.
The shapes of the wings of an aircraft in steady flight condition are
deformed by the constraint of aerodynamic loads.
CERFACS, in collaboration with Airbus, has developed an algorithm which,
for a given flexibility matrix representing the wings structure and a set of
flight conditions, computes the deformed shape of the wings and the fluid flow.
It implies that a robust mesh deformation algorithm must be available in elsA.
A first algorithm in this type has been successfully developed (see
Fig. 3.13 for an example of mesh deformation).
Figure 3.13: Shape deformation of a structural mode
The second area of interest is unsteady simulation. Efficient numerical
methods have been developed in order to enable unsteady simulation with elsA.
Direct industrial application is flutter prediction. Flutter is a destructive
fluid-structure interaction due to a transfer of energy from the fluid to the aircraft structure.
It is characterized by a growing oscillation which can lead to the destruction of the
aircraft. It can be predicted either with the use of unsteady aerodynamic loads
provided by unsteady simulations (see Fig. 3.12 for an example of the results
given by an unsteady simulation), or with a direct temporal fluid structure simulation.
The increasing accuracy and efficiency of high performance computing favor the
simulation of the temporal response of the nonlinear aeroelastic system.
3.2.4 Codes coupling (F. Loercher, J. Delbove)
In the framework of internal collaborations and in-house software's development and valorisation, the Group has conducted some coupling experiments with the PALM software (see the "Climate and Modelling and Global" chapter). Two types of simple coupling have been achieved : aeroelasticity with elsA and a simple in-house Computational Structural Mechanics (CSM) code and coupling elsA with aerodynamic analysis tools for on-the-fly post-processing.
These assessments have shown that PALM is relevant, easy to use and efficient for these applications.
Figure 3.14: Mach number iso-contours for the MARCO nozzle configuration
3.2.5 Treatment of axisymetric configurations (S. Champagneux, M. Fontaine)
A rigorous and effective two-dimensional method for the specific treatment of axisymetric configurations has been implemented in the elsA software
[Champagneux, 2003a].
This method offers an alternative to the three-dimensional one that was initialy available and previously used.
The main drawback of the original approach was the CPU penalty due to the mandatory 3D numerical treatement despite the intrinsically 2D configuration. Moreover, the need for specific boundary conditions and the definition of an arbitrary azimuthal angle are case-dependent procedures which have limited scientific interest.
A finite volume consistent approach relying on the addition of an extra source term and a metric modification has been implemented and verified for several steady test cases (incompressible laminar and turbulent Poiseuille flows, laminar hypersonic blunt body flow).
Finally, the MARCO configuration (a ventilated nozzle from SNECMA-moteurs, Fig. 3.14) where flow conditions lead to an under-expanded jet is computed with this approach. A speed-up factor of three is obtained concerning the CPU time for convergence and obtained results agree with those provided by SNECMA-moteurs.
3.2.6 Prediction of dynamic derivatives of full aircraft including large winglets (D. Margerit, S. Champagneux, G. Chevalier)
CERFACS also studied the influence of large winglets on aerodynamics characteristics of the wing, and more precisely the use of advanced CFD computation to determine such characteristics through the use of the elsA and Arbitrary Lagrangian Eulerian (ALE) techniques.
-effect
(incidence) and pitching effect (q-effect) computations have been validated on the M6 wing geometry. Values of quasi-static derivatives
 and 
are in good agreement with previous results. Full multi-grid, grid sequencing and solution continuation are used in order to accelerate convergence to steady state solutions.
Invicid (Euler) computations of
and
were performed for a coarse grid of the A340 MSR1 geometry (see Fig. 3.15). The process of transforming the Euler A340 grid to a NS (viscous) one is also done. This work will go on with NS computations and large/standard winglets comparisons.
Figure 3.15: A340 MSR1 configuration.
|
|
|
