Aerospace Science and Technology 7 (2003) 239–262 www.elsevier.com/locate/aescte
Code verification/validation with respect to experimental data banks Richard Benay ∗ , Bruno Chanetz, Jean Délery Fundamental/Experimental Aerodynamics Department, ONERA, Meudon Centre, 92190 Meudon, France Received 7 June 2002; received in revised form 12 November 2002; accepted 22 January 2003
Abstract The past 40 years have known a spectacular development of CFD capabilities. It is now possible to compute complex three-dimensional unsteady flows even at the design stage by solving the Unsteady Averaged Navier–Stokes Equations (URANS approach) and progress are made every day in still more advanced approaches such as LES and DNS. However, the confidence in CFD methods is still limited because of uncertainties in the numerical accuracy of the codes and of the inadequacy of the turbulence models they use. Thus, there is still a need for well made and well documented experiments to validate the codes and to help in their improvement. Such experiments must also fulfil quality criteria to be considered as safe enough and really useful for code validations. The article presents a discussion of the strategy to be followed to ensure the reliability and accuracy of a code by placing emphasis on the experimental aspects of code validation. The purpose is illustrated by considering recent examples of CFD validation operations based on basic – or building block – experiments. The first case considers an experiment on a purely laminar shock wave/boundary layer interaction used to assess the numerical accuracy of several codes. Other examples deal with the crucial problem of the validation of turbulence models in strongly interacting flows. The conclusion stresses the importance to constitute high quality data banks on typical flows still difficult to predict. The problem of data dissemination is also briefly addressed. 2003 Éditions scientifiques et médicales Elsevier SAS. All rights reserved. Keywords: Data basis; Supersonic flow; Transonic flow; Hypersonic flow; Separated flow; Control of interactions; Turbulence modeling
1. Introduction Methods for verifying the capability of a code to solve given equations have been the object of close examinations in the past. Identification and elimination of various types of errors and use of precision criteria, methods for convergence testing, rules for establishing grid convergence, are all required when one have to assess the quality of the numerical tool. The various stages of the general process of verification that will give to the code a confidence label permitting to use it for testing theoretical models have been summed up by Roache with references to many authors [24]. This process constitutes by itself a complex program often partially carried out but that should be completely satisfied in the ideal cases. A second step is devoted to the validation of models aimed at predicting flows that cannot be presentable as clearly identified solutions of well known mathematical problems. At this stage, comparison with experimental data is indispensable. * Corresponding author.
E-mail address:
[email protected] (R. Benay).
Formerly, validation of predictive methods was made by comparison of the computed results with some measured wall properties, essentially the pressure. In many situations, this kind of comparison was sufficient since “old” predictive methods, which were either fully empirical or based on a multi-component type approach, allowed only the prediction of the wall properties (pressure, skin-friction, heat-transfer) and the global performance of the vehicle. They could also give a gross idea of the flow organisation by predicting, for example, the size of a separated region and the location of a separation point, but this information was more or less considered as qualitative. The flow prediction landscape has completely changed with the advent of theoretical models based on the solution of the Navier–Stokes equations. It is clear that this approach is the only suitable to compute complex flows containing shock waves, centred expansion waves, separated regions, shear layers, etc. . . . the dissipative region being turbulent in nearly all the practical situations. Then, not only the wall properties are computed but also field quantities including the mean velocity and the turbulent quantities. However, in its present state the Navier–Stokes approach is still far from being free of critics, difficulties persisting both on the numerical side and in the physical
1270-9638/03/$ – see front matter 2003 Éditions scientifiques et médicales Elsevier SAS. All rights reserved. doi:10.1016/S1270-9638(03)00018-X
240
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
modelling, in particular of turbulence. There is thus a strong need to validate Navier–Stokes codes before their routine use for design purposes. Though the prediction of wall properties remains a key target for most predictive methods, since the drag and the thermal loading are the quantities of most practical interest, it rapidly became obvious that a comparison restricted to wall properties was insufficient to properly validate the most advanced predictive methods. In general, Navier–Stokes codes give a faithful and impressive picture of the flow field structure. For example, in propulsion applications, the complex organisation of the jet with its pseudo-periodic pattern of shocks and expansion waves, the separated zone forming inside an overexpanded nozzle or on the afterbody in case of jet pluming are most often well reproduced. However, a more careful analysis of the data shows that the situation is far from being entirely satisfactory. Firstly, it is observed that a fairly good prediction of the wall pressure can coexist with a poor quantitative prediction of the velocity field. Frequently, the extent of the separated region is underestimated, sometimes considerably. In addition, the turbulent quantities are poorly predicted, especially if the flow is separated. Such discrepancies render suspect the validity of the code since they are indicative of some basic deficiency, either in its numerical scheme or its turbulence model, or both. Secondly, a rather fair prediction of the flow field can be accompanied by large errors in the calculation of surface properties affecting mainly transfer coefficients such as skinfriction and heat transfer. Lastly, in certain applications, the knowledge of the outer field itself is of prime interest. This is the case for problems dealing with infrared signature where detailed description of the hot propulsive jet, with an exact localisation of the Mach discs, is essential. Pollution studies necessitate a good prediction of the jet properties to allow accurate evaluation of chemical processes and species concentration. A good representation of the flow resulting from interference between the external stream and the propulsive jet(s) also requires a faithful prediction of the flow field. This is also the case of shock intersections which can be of prime importance for the performance of air-intakes at high Mach number. The problem of code validation is still more stringent in three-dimensional applications where the Navier–Stokes approach becomes mandatory. Due to the complexity of such flows, it is clear that the consideration of the surface pressure alone is completely inadequate because this information gives a very partial view of the flow (in three-dimensional flows, it is no longer possible to infer separation from an inspection of the wall pressure distributions). In these conditions, the validation of advanced computer codes requires well documented experiments providing detailed and reliable flow field measurements. The progress accomplished in the measurement techniques over the past 40 years, mainly with the advent of laser based optical methods, has operated a true breakthrough in our capacity to in-
vestigate complex turbulent flows, containing shock waves, strong expansions, thin shear layers and recirculating regions. We therefore possess now the necessary technical means for performing such investigations. This paper deals with the global problem of code and models validation by comparison with experiments. After purely numerical tests, this step is essential for determining the degree of reliability of a code using a given model. Two parts can be distinguished: • A first part is devoted to the strategy for code verification and validation with emphasis on the experimental facet of this action. This strategy is presented in the optics of a physical approach of the problem in which one focuses on the prediction of the flow fundamental properties. The performance aspect is not considered, although it is the end product of the validation chain; this ultimate step is more in the hands of engineers than scientists. • In a second part, the verification/validation strategy is illustrated by examining three fundamental problems met in modelling of supersonic flows. They concern laminar separation at high Mach number, phenomena of internal aerodynamics, affecting supersonic air intakes and supersonic base flows. The common point of these examples is flow separation, a problem presently far to be solved with all the required accuracy.
2. The code verification/validation process 2.1. What is required from the code Before considering a verification/validation action, first of all, the aim of the calculation must be clearly stressed and identified. • If calculation is used to predict the performance of a system or a sub-system, accuracy is mandatory since then the engineer has to rely on the calculation to define the object which must satisfy prescribed specifications or to evaluate the aerodynamic performance achieved by a given object, in terms of drag, thrust, maximum range, stability, fuel consumption, manufacturing and maintenance costs, etc. This consideration applies to an airfoil, a wing, an air-intake, a nozzle, a propelled afterbody, etc. Such aerodynamic calculations are frequently coupled with structural analysis, aeroelasticity, thermal calculations, flight mechanics. • In the design of machines involving complex flows whose experimental simulation is difficult, if not impossible, a calculation showing the flow field organisation is of great help for the designer. In this case, accuracy is not essential, but reliability is crucial since one must be confident on the physical features of the computed field. This applies to flow in turbomachines where a detailed experimental description of the flow is still largely out
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
of reach. This is also true for the flow past multi-body space launchers where information on the main flow features can be useful to detect possible hazards and to take precautions to avoid them. • From a fundamental point of view, the physical understanding of complex flows has to be based on a parallel theoretical analysis whose aim is to help in the interpretation of the phenomena and to establish the descriptive approach on firm rational arguments. In this case, accuracy is not needed since, for example, theoretical analyses are nearly always derived from more or less strong simplifying assumptions rendering quantitative results questionable. This is for example the case of stability analyses, perturbation methods or asymptotic expansion techniques. These theories have greatly helped in the understanding of critical phenomena, in spite of a relatively poor quantitative agreement. Numerical calculations are also a precious tool to construct the structure of complex flows by providing information on “hidden” faces of the phenomena which cannot be reached by direct observation. This point concerns for example the topological construction of three-dimensional separated flows. Here accuracy is not essential; even reliability can be limited since in the process there is a continuous exchange between experiment and calculation – or theory – which allows a cross fertilisation of the two approaches. • In the last issue, a code is used as a tool to test a new physical model (or to introduce some improvement in an existing one). Then, numerical accuracy is mandatory since it would be vain to implement a good turbulence model (if there is one) in an inaccurate code in which high gradient regions where turbulence is at work are not correctly captured. 2.2. The different steps of the validation Considering the above points, it appears that the verification/validation procedure has to be submitted to a four step strategy or tactical actions. First step: assessment of the code numerical accuracy. Before going into a more involved operation of code validation, encompassing all the aspects of the computing action, it seems obvious to first establish the accuracy and reliability of the code by focusing on its numerical aspects. A clear assessment of this point is not a straightforward issue in the sense that the numerics involves several closely linked aspects. Assessment of the code numerics can be done first from comparison with exact analytical solutions or with well established empirical results. Most often these exact solutions are only available for laminar flows. In a second phase, verification can be made by confrontation with other codes developed by independent teams. This action implies a close co-operation between the persons involved in the procedure with a complete exchange of infor-
241
mation on the calculation methods and a common evaluation of the results. A precaution for this action to be conclusive should be to agree on exactly the same configuration and to insure spatial convergence in all the calculations. If the comparison between the codes exhibits differences which cannot be resolved after a thorough analysis of the results, then one has to refer to the verdict of an independent authority which can be experiment. A further – and hopefully conclusive – step will be to run the competing codes on a configuration for which good and reliable experimental results exist. This point is far less obvious that it would appear at first sight, since the experimental data should allow to draw clear conclusions. Thus, one should avoid to mix modelling and numerical problems by considering a simple laminar case, in a calorically perfect gas in order that the basic flow physical properties (thermodynamics constants, molecular viscosity) are well known. This point is far from being easy to achieve since it is most often very difficult to maintain an entirely laminar regime in classical aerodynamics. The geometry of the body should be simple and completely defined to avoid complex meshing problems. From this point of view, a two-dimensional case is preferable but, as it will be seen below, an axisymmetric configuration is preferable. Experimental requirements will be exposed more completely below. Second step: validation of the physics implemented in the code on elementary configurations. This is the most important point for the specialist in flow physics, the first step being only a preliminary step simply aiming at verifying the tool. In the second step, the code is used to compute what can be considered as the elementary components of an aerodynamic flow: attached boundary layer, laminar-turbulent transition on a flat plate, separation induced by an obstacle, flow past a base, shock wave/boundary layer interaction, start and development of a vortex structure, vortex breakdown, shock/shock interference or shock crossing, etc. Twodimensional – preferably axisymmetric – as well as threedimensional basic situations have to be considered. For this first validation step, the numerical results are compared with building block experiments focusing on a specific elementary phenomenon. Third step: validation on more complex sub-systems. Once the code and its physical model(s) have been validated on basic cases, a more complete configuration must be considered consisting in a sub-system of a complete vehicle, where several elementary phenomena are combined. This is the case of a profile on which one encounters laminarturbulent transition, attached boundary layers, transonic shock wave/boundary layer interaction, separation, wake development, etc. The wing constitutes a three-dimensional extension with the additional problems of the vortices emanating from the wing and control surfaces extremities. The supersonic air-intake involves shock/shock interference,
242
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
shock wave/boundary layer interactions, corner flows with vortex development. After-bodies combine supersonic jets with complex shock patterns (Mach disc formation), shock induced separation, either inside the nozzle (overexpanded nozzle) or on the fuselage (jet pluming at the exit of an underexpanded nozzle). Many other examples can be cited: propulsive nacelle, compressor/turbine cascade, helicopter rotor, etc. Fourth step: validation on the complete vehicle or object. This is the ultimate target in which the code is applied to a complete airplane, an automobile, a space launcher, a helicopter, etc. This part of the exercise is the domain of engineers and will not be developed here. 2.3. Requirements for good test cases constitution We are thus naturally led to define what are the requirements for really useful experimental results and ask the question: is the data bank safe? This is a vital issue on which we will concentrate since the points which follow are the every day concern of experimentalists. Definition of the geometry. A first condition for any experiment aiming at the verification of the numerical accuracy of a code and the validation of its physical models is to focus on a configuration whose geometry must be both representative of a typical situation, precisely defined and as simple as possible. A combination of plane surface (like a twodimensional ramp) is a good choice but one should carefully avoid special situations, like a sharp leading edge, leading to singularities and to meshing difficulties. When possible, an analytical definition of the contour should be provided. It is preferable to give the dimensions in metric units to avoid risk of confusion in the reference length used to compute a Reynolds number. When possible, a two-dimensional geometry should be adopted – even for three-dimensional problems – since it offers many advantages to visualise the phenomena and to execute measurements, in addition of the lower cost of the test-set up fabrication. Furthermore, the original set-up must frequently be modified before arriving at a fully satisfactory flow; such modifications are far easier on a two-dimensional arrangement. Boundary conditions. The boundary conditions must be well identified and accurately known. This concerns the upstream flow conditions (Mach number, velocity, pressure, density) when a uniform incoming flow exists. In transonic experiments executed in a channel type arrangement, one often considers phenomena taking place on the channel walls, the test section itself being the model. In this case, a well defined origin with a uniform state at upstream infinity does not exist. Then, the data should provide all the flow conditions in a section located sufficiently far upstream of the region of interest, including the boundary layer properties (mean velocity profile, turbulent quantities). In this context, a good
example is the evaluation of the dissipation when attempting to test a two-equation turbulence model. If LDV measurements now permit to know with a good approximation the Reynolds stress profiles in moderate supersonic flows, a method must be conceived for deducing the dissipation from these data. This point will be addressed in detail in the second part of this article, where a method is proposed to evaluate the second transported variable. In all cases the stagnation conditions (pressure, temperature) and the incoming stream thermodynamic properties must be given. Downstream boundary conditions leading to a well posed problem must be provided. If the flow leaving the zone of interest is supersonic, then no-conditions have to be imposed to perform the calculation. The question of the downstream conditions is more delicate if the configuration is such that the flow leaving the test region is subsonic. When the downstream flow is again uniform, most often a downstream pressure is given, since this quantity is easily obtained. It is far more difficult to provide the pressure field in a complete plane, as some theoreticians sometimes ask for. In transonic channel experiments where a shock is produced by the choking effect of a second throat, the best way is to provide the geometry of the second throat and, in the calculation, to impose downstream conditions insuring the choking of this throat. Perturbing effects. Side effects or uncontrolled perturbations must be avoided, except if they can be taken into account by the calculation. The side effects due to the finite span of any experimental arrangement strongly affect the flow when separation occurs. Then, the experimented flow can be very different from the assumed ideal twodimensional flow which would correspond to the infinite span condition. Confrontation of such an experiment with a planar two-dimensional calculation can be deprived of any signification and lead to entirely erroneous conclusions. Nevertheless, many experiments show a limited zone in the vicinity of the wind tunnel mid longitudinal plane that can be considered as reasonably two dimensional. Even in this case, in wind tunnels, where the interactions under scope are adjusted by downstream conditions or by choking effect due to a second throat, the phenomenon is globally displaced with respect to the ideal case of the wind tunnel of infinite span. This is due to the interacting boundary layers of the lateral walls, that modify significantly the downstream effective sections of the wind tunnel when compared to the ideal case. When the boundary layer is attached or weakly separated, as in the example chosen in the following paragraph, this effect is small and 2D calculations remain acceptable. If one desires to keep the mathematical simplicity of two space dimensions, the best is to compute an axisymmetric flow. When the goal of an experiment is to reproduce closely the behaviour of an aircraft, the coincidence of the wind tunnel’s Reynolds number with that of real flight tests is a condition often hard to satisfy. Even if this condition is
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
fulfilled, premature transition can occur on the model. The boundary layers on the walls of classical wind tunnels are generally turbulent and unsteady perturbations radiate from the wall to the whole test channel. Such perturbations lead to transition of the flow around the model at Reynolds numbers significantly lower than those where it occurs in real flying tests. The conception of quiet wind tunnels, at NASA, ONERA and Purdue University, has been undertaken to overcome this difficulty. Inversely, when a plainly turbulent regime is searched for code validations, it is preferable that the wind tunnel be naturally turbulent without the help of auxiliary means like transition strips. The transition in general is indeed far from being treated satisfactorily by the methods of calculation presently at our disposal and transitional effects are often present when turbulence is set up artificially. Experimental needs. The description of the flow must be as complete as possible in order to permit a thorough validation of the code and to provide all the information useful to understand the physics of the flow and to help in the elucidation of causes of disagreement. This concerns wall quantities, like pressure, heat transfer, skin-friction, field quantities such as stagnation pressure and temperature, mean velocity, Reynolds tensor components, density, species concentration, etc. Flow visualisations are highly desirable to give a precise idea of the flow field structure: shock waves, location of shear layers, separated region, vortices. A really complete description of the flow is never possible because the experiment cost would become incompatible with the budgets currently allocated to fundamental aerodynamics. Measurements reliability and accuracy. The experimental data must be reliable. This means that the experiment is not “polluted” by an extraneous phenomenon due to a bad definition of the test arrangement or to a bad functioning of the facility. The general regulation tends now to impose to experimentalists precise information on the uncertainty margins of their measurements. This information, although important for an in depth validation of codes, is not always as essential as claimed if the objective is to test the physical validity of a complex model. The problem is different in the case of performance determination, where quantities like lift, drag, efficiency must be known with high accuracy. The physical interpretation. Constitution of safe data banks is not restricted to the execution of hopefully good experiments in relation with code development. The experimentalist must also be a physicist able to interpret its findings and to understand the physics of the investigated flow field. This interpretation, which must be based on theoretical arguments, is essential to insure the safety of the results. It must precede any numerical exploitation.
243
3. Examples of application of the code validation program 3.1. Scope of these examples High Mach number separation is a major field of investigation since it appear both in the domain of internal and external aerodynamics, the second one being mainly concerned in the case of flows around flaps. High Mach number air-intake flows involve strong shock waves and stream confluence’s which are at the origin of shock/shock interactions, extended separations, shear layer developments, large scale fluctuations, these phenomena taking place in turbulent compressible flows (except in very high altitude flight). All these points make the prediction of such flows difficult since they involve the hardest points met in applied fluid mechanics. The present predictive capacity being still limited there is a strong need to maintain a sustained research effort both on the experimental and theoretical sides for propulsion applications. These problems are particularly critical in the case of a hypersonic air-breathing spacecraft which flies in conditions where shock phenomena are much amplified, leading to still more severe interactions between the internal and outer flows. Then, integration of the propulsion unit in the vehicle architecture is a vital issue. The conception of the air-intake must take into account two major aerodynamic interactions. Firstly, compression of the capted high Mach number flow is achieved partly by taking opportunity of the shape of the vehicle front part which allows to realise a nearly isentropic compression, partly by a succession of ramps and/or shock reflections constituting the air-intake itself. The design of the compression ramp system is a complex matter combining the search for maximum efficiency, minimum cowl drag, minimum length or weight of the system, possibility of adaptation, etc. In particular, one is confronted with the risk of separation at a ramp or a shock reflection, with the subsequent loss of efficiency, possible occurrence of instabilities and air-intake un-start. This problem is not restricted to air-intake since separation is likely to occur in several other parts of the vehicle. As it is well known, separation is most often avoided by designers since it leads to a degradation of performance and a rise of the nuisance produced by the vehicle. All these reasons make essential an accurate prediction of separated flows. A satisfactory prediction of the turbulence evolution in separated flows is still not assumed by any of the presently existing turbulence models. A typical configuration where an extended zone of recirculation exists is the base flow, which constitutes an excellent test case for validating any model. Base flows have been the subject of numerous experimental and theoretical studies since the fifties in order to understand the physics of such flows which are of prime importance for projectiles, missiles or space launchers. Due to the complexity of base flows, there is still a need to constitute reliable theoretical tools for predicting the drag of after-
244
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
bodies and also the aerothermal loads on the base region of propelled after-bodies. The examples chosen in what follows present three typical test cases, constituted at ONERA for two of them, to help in the development of safer and more accurate codes. They are aimed at illustrating most of the problems evoked in the first part. The first one is a good example of a nice verification action involving five codes. For the second cases, that is relatively new, only one home-made Navier–Stokes code has been run with the main objective to check the completeness and safeness of the test case. The third case is now a well known test case for which previous calculations have been referenced. Before addressing the description of these examples, we must point out that they are two dimensional, why? Our experiment shows us that, if we had to treat simple geometry’s of the same kind in three dimensions, with the required grid resolution for space convergence, grids of 20 to 50 millions points should be necessary, even in the case of RANS calculations. According to the more optimistic specialists’ projections on future computer’s improvements, current and extensive calculations with this kind of grids will not be possible before around 2015. Turbulence modelling requires numerous parametric calculations that will therefore not be permitted on such grids before that time. Nevertheless, we have to anticipate now the subsequent future needs for validations on very detailed three dimensional experimental data. In particular, a complete scanning of the flow will necessitate facilities permitting the non intrusive probing of lateral boundary layers. The development of experimental set-ups capable of this kind of measurements in flows of limited complexity is possible now but necessitates huge investments. 3.2. The cylinder-flare case for hypersonic laminar separation 3.2.1. Aim of the operation The aim of the present operation was to assess the ability of numerical codes using various schemes to predict separation on a 30◦ axisymmetric ramp flow at a Reynolds number sufficiently low to insure laminar flow in the whole interaction domain [9,16]. The action consisted in extensive comparisons between measured and calculated wall and flow field properties. Since the calculated wall quantities are particularly sensitive to grid constitution, great attention was paid to the quality of the grids. The results of five codes have been compared in this action: • two Navier–Stokes codes from ONERA using a finitevolume approach, • a Navier–Stokes code from the University of Rome ‘La Sapienza’ using a finite-volume approach, • a Navier–Stokes code from DLR using a finite-element approach, • a Direct Simulation Monte Carlo code (DSMC) from NASA Langley Research Centre.
Fig. 1. Geometrical definition of the cylinder-flare axisymmetric model.
3.2.2. The experimental part Wind tunnel and testing arrangement. The experiments have been executed in the low Reynolds number ONERA R5Ch wind tunnel. For the nominal stagnation conditions (pressure pst = 2.5 × 105 Pa, temperature Tst = 1050 K), the upstream Mach number is equal to M0 = 9.91 , which corresponds to an upstream flow with the following properties: • • • •
static temperature: T0 = 51 K, density: ρ0 = 0.43 × 10−3 kg m−3 , static pressure: p0 = 6.3 Pa, unit Reynolds number: Reu = 186,000 m−1 .
The model used is defined in Fig. 1. It is constituted by a hollow cylinder, with a sharp leading edge, followed by a flare terminated by a cylindrical part. The flare angle is β = 30◦ , which leads to a large separated zone. The flare is followed by a cylindrical part in order to facilitate the computations. The advantage of the present configuration is that it is able to displace the base flow sufficiently far downstream of the interesting area, in such a way that the complex phenomena occurring in the base region have no effect on the interaction region. The model has a total length of 170 mm, the reference length, based on the distance between the sharp leading edge and the beginning of the flare being equal to L = 101.7 mm. The Reynolds number calculated with L is equal to ReL = 18,375. The outer diameter of the cylinder is equal to 65 mm and its inner diameter to 45 mm. Three models have been built: • One made of steel without instrumentation was used for wall flow visualisations, outer field visualisations and measurements by the electron beam fluorescence technique (EBF). • One made of steel with a ceramic insert made of MACOR on which were installed 24 platinum films for heat flux measurements. • One made of steel with 24 pressure orifices located on four generating lines.
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
245
Techniques of investigation. The flow field was visualised by EBF. Each photograph required an exposure time of 10 s during which the electron beam scanned the domain explored at 50 Hz to create a visualisation plane. The wall pressure measurements were executed with variable reluctance differential transducers installed in the testing chamber and connected to the pressure taps by rubber tubes. The accuracy of pressure measurements was ±5%. Wall pressure measurements are transformed into a pressure coefficient defined by: Cp =
p − p0 1 2 2 γp0 M0
where the upstream conditions M0 and p0 were deduced from the incoming flow stagnation pressure measured with a Pitot probe. The increasing of wall temperature during the run does not affect strongly the validity of pressure measurements as it has been observed in the framework of an ESA Technical Request for Proposal, called “Hot Measurement Techniques” [10]. The heat fluxes were determined from the surface temperature rise during the first half second of the run. A thermometer element (platinum film) was applied to an insulating support made of MACOR ceramic. From the evolution of temperature with respect to time, heat flux can be deduced for each sensor implemented on the model wall. The calculation of the heat flux hw from the surface temperature was performed by solving the heat equation assuming that the ceramic support behaves as a thin plate with parallel faces and with the hypothesis that the temperature of the ceramic internal face is equal to the model temperature (at the start of the run). The accuracy of heat flux measurements was ±7%. The results are expressed in a dimensionless form by introducing the Stanton number defined by: St =
hw ρ0 U0 Cp (Tst − Tw )
where ρ0 and U0 are calculated from the Mach number M0 and the stagnation conditions, Tw being the wall temperature, Tst the stagnation temperature (Tst = 1050 K) and Cp the constant pressure specific heat. In the above Stanton number, the adiabatic wall temperature and the local conditions at the boundary layer edge were not used because in hypersonic separated flows these quantities are hard to know. As far as calculations are concerned, it is more realistic to consider an isothermal wall than an adiabatic one. Indeed experimentally the wall is heated during the run, but the wall temperature increase for calculations can be neglected since we consider only the beginning of the run to deduce the wall heat fluxes. The wall heat fluxes are deduced from the temperature increase during the first half second of the run. The temperature evolution is recorded at a frequency equal to 100 Hz. So in one half second, we have sufficiently measurement values to deduce heat fluxes from the wall temperature evolution for each sensor. The temperature reached by the model one half second after the beginning
Fig. 2. EBF visualisation of the Mach 10 flow past the cylinder-flare model in the R5Ch hypersonic wind tunnel.
of the run is very close to the temperature at the beginning of the run, even at the sharp leading edge. A maximum of fifteen Kelvin can be measured, which leads to a variation of 5% of the nominal temperature (290 K), which justifies the assumption. Furthermore, this value is very far from the adiabatic temperature, around 850 K, in the laminar R5Ch Mach 10 flow (Tst = 1050 K). The density measurements were performed using the electron beam X-ray density system using the two-detector method: one measuring the reference density in the freestream flow outside the boundary layer and the other to explore the boundary layer. However the total duration of the run, around one minute, is used for field measurements by Xray detection. In this case six or seven points are measured along a line located perpendicularly to the wall. Once we measure the points located closer to the wall and at the end of the run we measure the points not located at the wall proximity. Thus the wall does not influence too much the measurements. Experimental results. The photograph in Fig. 2 shows the flow past the model. We see the attached shock wave at the sharp leading edge and the separation shock wave which meet above the end of the flare. The surface flow was qualified by viscous coating visualisations to check axisymmetry and to precisely locate the separation and reattachment lines. Fig. 3(a) presents the longitudinal evolution of the pressure coefficient Cp measured for ten runs. The pressure coefficient decreases slowly on the upstream part of the cylinder till abscissa X/L = 0.7. This decrease is followed by a two step compression. A first compression is due to the separation process and a second compression, far more important, takes place during reattachment. The maximum pressure value reached at the end of the flare is more than forty times higher than the static pressure p0 . The Stanton number evolution plotted in Fig. 3(b) shows that the heat flux decreases along the cylinder in the upstream part, where the flow is governed by the viscous
246
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
Fig. 3. Measured surface pressure and heat transfer distributions on the cylinder-flare model in the R5Ch wind tunnel; (a) pressure coefficient, (b) Stanton number.
interaction effects emanating from the leading edge. Further downstream, the evolution is characteristic of a laminar interaction with a large boundary layer separation: a slight decrease on the heat flux evolution occurs at the abscissa where separation begins. In turbulent flows, heat transfer increases at separation. 3.2.3. The confronted codes The FLU3M Navier–Stokes code. The FLU3M code is a general Euler and Navier–Stokes flow solver developed at ONERA for supersonic and hypersonic 3D applications on multi-block topologies [8]. The main characteristics of the FLU3M code in its Navier–Stokes option are a cell centred formulation well suited for multi-block strategies and a MUSCL (Monotonic Upstream Scheme for Conservation Law) approach for the inviscid terms which guarantees accuracy and robustness for supersonic and hypersonic applications with severe discontinuities. A set of different upwind schemes (Van Leer, Roe, Osher, . . .) and different flux limiters (minmod, Van Albada, Koren, . . .) is available and used according to the type of application. Acceleration techniques including Van Leer implicit associated to Alternate Direction Implicit and local time stepping are used to speed up the convergence of the solution. A classical central discretisation is used for the viscous terms. For the perfect gas option, the viscosity is defined by the classical Sutherland law. When equilibrium air option is applied for treating flight cases, the Srinivasan model or the Gazaix model is selected. A complete set of boundary conditions is available, including isothermal and zero pressure gradients for wall boundary treatment. In addition, different turbulence models may be selected, which have not been employed in the present study. A detailed validation of the laminar option of the FLU3M code was performed by running computations on a hyperboloid flare axisymmetric configuration in hypersonic flow conditions. The use of such an axisymmetric geometry as a test-case allowed grid convergence studies. Comparison with the wind tunnel results has shown that for purely laminar flow conditions, there is an excellent agree-
ment between computations and experiments for blunt configurations. The present computations were performed with the Flux Difference Splitting scheme of Roe associated with the Van Albada limiter. A correction of the entropy function Ψ is applied as defined by Harten. The value of the correction factor normalised by the spectral radius of the Jacobian matrices of the inviscid fluxes is very small near the wall (δ0 = 3 × 10−5 ) while larger values (δ0 = 0.2) are used out of the boundary layer. Concerning convergence, the classical 3-orders of magnitude criterion on density residual decrease was enforced by running additional iterations after stabilisation of density residual evolution until the stabilisation of the separation extent. The NASCA Navier–Stokes code. The solver NASCA is a research code developed at ONERA to be used in parallel with fundamental experimental activity, for exploitation of experimental results and for modelling studies. It is based on a finite-volume method for solving the classical Reynolds Averaged Navier–Stokes (RANS) equations [2]. The problem investigated for this study is a steady one, thus allowing the use of Beam and Warming type time discretisation, which is first order accurate in time. With regard to the space discretisation we have used an implicit approach, whereby the numerical function f (f = FE − FV , where FE and FV are, respectively, the inviscid and viscous fluxes) is obtained by means of a quasi-linearisation approach which gives: ∂ f¯ W f¯ W n+1 = f¯ W n + ∂W with: W = W n+1 − W n (W is the vector unknown), where indices n and n + 1 indicate the values at time nδt and (n + 1)δt, and ∂ f¯/∂W is the Jacobian matrix. The convective matrix ∂FE /∂W is split into two parts corresponding to negative and positive waves, giving a first order space approximation of the convective part with good stability properties. With regard to the diffusion Jacobian matrix ∂FV /∂W , a centred discretisation scheme has been
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
used. The convective contribution of the explicit part is split following Steger and Warming [30]. The scheme is an adaptation of the scheme of Osher and Chakravarthy [22] to the case of locally non-uniform grids that can be strongly oblique [2]. This scheme, both with the centred one used for the diffusive part of the flux, insures a second order space precision for these kinds of grids. The use of Osher and Chakravarthy scheme together with the carefully constructed implicit part avoids the use of the entropy correction often necessary for obtaining a physically acceptable solution. The HIG-2XP Navier–Stokes code. In this code, developed by the University of Rome, the governing equations are discretised by using a cell centred finite-volume formulation and approximating surface and boundary integrals by means of the mean value theorem and mid-point rule [17]. The following system of ordinary differential equations is obtained for each computational cell: dW +R=0 dt where R is the residual of the equations representing the inviscid (R = (F n s)) and viscous (RV = E E,num β (F n s)) flux contributions, i.e.: V ,num β R= (FE,num + FV ,num)n s − H V
V
β
where β stands for the generic cell face, n is the outward unit normal to cell face whose length is s, V is the cell area, and H accounts for terms that arise from the finite-volume discretisation of axisymmetric geometry’s: H = [0, 0, hax,v , 0]T and hax,v = (p − σ0 ), v 1 − div u . σ0 = 2µ r 3 The time integration is performed by a three-stage (m = 3) Runge–Kutta algorithm and a local time-stepping procedure. With regard to the Euler flux contribution the scheme satisfies the entropy condition. Its properties are second order upwind-biased Total Variation Diminishing. For example, at the generic cell face (i + 1/2, j ), a characteristic decomposition in the direction normal to cell face is introduced, thus obtaining: FE,num =
1 (FE )i,j + (FE )i+1,j 2 4 1 k ψ(ai+1/2 ) gik + gi+1 − 2 k=1 k k k k − ψ ai+1/2 + γi+1/2 αi+1/2 ri+1/2
k k where ai+1/2 and ri+1/2 are, respectively, the kth eigenvalue of the normal inviscid flux Jacobian and the kth right
247
eigenvector. Roe’s averaging is used to evaluate the states of (u, v, c, H ) at cell faces. The vector αi+1/2 represents the difference of the characteristic variables, and is defined as: αi+1/2 = Li+1/2 i+1/2 W with Li+1/2 the left eigenvector matrix. The entropy function Ψ here adopted is that of Harten. At all boundaries the normal derivative of g has been set equal to zero, and conservation is enforced by imposing (for all cells along the walls) that the contribution that makes the scheme satisfy the TVD property is set to zero. The FLOW Navier–Stokes code. The solver FLOW, developed et DLR, is based on a finite-element method for solving the two-dimensional/axisymmetric, unsteady, compressible Navier–Stokes equations [27]. For a solution of these equations, a given computational domain has to be subdivided into cells, the so-called elements, where the governing equations are solved by an integral approximation. The applied numerical method uses an explicit Taylor– Galerkin algorithm in a weighted residual form, with a time integration according to the two-step version of the Lax–Wendroff scheme. The first step represents a firstorder accurate Taylor expansion in time. This predictor step provides, by its cell-centred formulation, a solution for each element in turn, taking into account convective fluxes only. The second step is then the node-centred corrector step, which balances the cell-centred fluxes surrounding a node, leading to a second-order accurate solution in time. For the prediction of steady state solutions, local time step relaxation, in combination with a diagonalised, lumped mass matrix, is used. To improve shock capturing, especially with regard to high Mach numbers, the concept of fluxcorrected transport (FCT) is applied to the Taylor–Galerkin algorithm. The principle of FCT can be characterised as a limitation of the anti-diffusive fluxes of the Lax–Wendroff scheme, with the purpose of obtaining monotonic, first-order accurate solutions in the vicinity of shocks, while preserving the second-order accuracy in smooth regions of the flow. By employing unstructured grids, especially when using triangular elements, finite-element methods permit a high level of geometric flexibility, enabling a detailed resolution of high gradient variable regions. The triangulation of the computational domain is accomplished by an automatic mesh generation algorithm following the advancing front method, where points and elements are introduced simultaneously, allowing significant changes in the local mesh structure. To further improve the resolution of boundary layers, structured subgrids can be employed at wall surfaces. The applied mesh generation scheme offers the possibility of adapting the grid to the solution. The grid adaptation is achieved by a complete regeneration of the mesh, based upon information provided by the computed solution on the present grid. By using the automatic mesh generation algorithm, the new, adapted grid allows a significant variation in element size as well as a stretching of the elements. The new
248
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
Fig. 4. The two computational grids used for Navier–Stokes calculations on the cylinder flare model. (a) Grid G1 (201 × 77), (b) grid G2 (289 × 97).
distribution of the elements and nodes is controlled by an error estimation, based on the second derivative of a certain scalar key variable, i.e., density or Mach number. By means of this remeshing, the grid is accurately fit to the flow conditions, leading to an improved solution quality in a computationally efficient manner. The DSMC code. This DSMC code has been developed by NASA-Langley Research Center [21]. The DSMC method provides a numerical capability that acknowledges the discrete nature of gas and thereby provides the capability of simulating flows across the complete flow spectrum of continuum to free molecular flow regimes. However, the computing requirements can become excessive for multidimensional continuum applications. The molecular collisions are simulated by using the variable hard sphere (VHS) molecular model proposed by Bird [5]. This model employs the simple hard sphere, angular scattering law, so that all directions are equally possible for post-collision velocity in the centre-of-mass frame of reference. However, the collision cross section is a function of the relative energy in the collision. Energy exchange between kinetic and internal modes is controlled by the Larsen–Borgnakke statistical model. For the present study, simulations are performed by using a non reacting gas model consisting of two species, while considering energy exchange between translational, rotational and vibrational modes. A rotational relaxation collision number of 5 is used for the calculations. The vibrational collision number is 50. 3.2.4. The confrontation with experiment Conditions of the calculations. As far as the three finitevolume Navier–Stokes calculations are concerned (FLU3M, NASCA and HIG-2XP) two different types of grid have been used (see Fig. 4). A preliminary grid dependency study by successive dichotomic mesh refinement in both directions has been performed by means of the FLU3M code, which shows that grid independence is nearly achieved with the medium (201 × 77) and fine (401 × 153) grids (axial × radial). So it was decided to perform the calculations with the HIG-2XP and NASCA codes on the medium grid, named G1. The characteristics of grid G1 at the wall were: • minimum X-direction spacing: 0.42855 mm, • maximum X-direction spacing: 1.42223 mm,
• minimum r-direction spacing: 0.02921 mm, • maximum r-direction spacing: 0.17689 mm. Another grid G2 (289 × 97) was also used, differing from G1 in the stretching at the wall and smoothness properties. The characteristics of the grid G2 are: • • • •
minimum X-direction spacing: 0.26506 mm, maximum X-direction spacing: 1.16750 mm, minimum r-direction spacing: 0.00693 mm, maximum r-direction spacing: 0.92280 mm.
For the FLOW results, the computational domain is discretised by using a hybrid grid consisting of approximately 61,500 elements (after grid adaptation). The structured subgrid situated on the isothermal model wall is composed of 40 layers, growing with geometric progression. For the DSMC results, the calculation is made with a four-region computational domain containing 78,100 cells, where each cell is further subdivided into four (2 × 2 subcells). The collision partners are selected within the subcell; consequently, the flow resolution is much higher than the cell resolution (however, the microscopic properties are extracted from averages within the cell). Wall properties. A comparison between the three finitevolume Navier–Stokes codes (using grid G1) and the FLOW and DSMC codes (that use their own cell distribution) reveals, as far as the prediction of the wall pressure distributions (see Fig. 5(a)) is concerned, a good agreement between the three Navier–Stokes codes with only two marginal discrepancies. The first located in the growth region of the separation bubble is essentially due to a difference in the prediction of the separation point. The second discrepancy is on the prediction of the wall pressure after the reattachment point (with differences attaining 10% on the upper corner of the ramp) immediately before the final expansion fan. It can be concluded that the Navier–Stokes codes show a common tendency to slightly over-predict the bubble length and the pressure level after the reattachment point, whereas the DSMC code highly over-predicts the peak pressure. The skin friction coefficient distributions (for which no experimental values were measured) are reported in Fig. 5(b), which again shows small differences in the prediction of the separation location computed with the three
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
249
Fig. 5. Comparison of wall quantities calculated on the cylinder flare model with grid G1. (a) Pressure coefficient, (b) skin friction coefficient, (c) Stanton number. Table 1 Separation extent (mesh G1 for Navier–Stokes codes) Abscissa (X/L) Separation Reattachment
FLU3M
NASCA
HIG-2XP
FLOW
DSMC
Experiment
0.72 1.34
0.72 1.34
0.74 1.33
0.74 1.33
0.75 1.32
0.76 ± 0.01 1.34 ± 0.015
Navier–Stokes codes. The abscissas of the separation and reattachment lines (the measured values are inferred from flow visualisations) are reported in Table 1. The DSMC results show a good prediction of the separation extent, the three finite-volume Navier–Stokes codes yielding an overprediction of about 10%. The finite-element Navier–Stokes code is in good agreement with experiment. The Stanton number distributions plotted in Fig. 5(c) confirm that the main differences between the calculations are confined to the reattachment region. It is instructive to show the results obtained with the three Navier–Stokes codes on the grid G2. The wall pressure
distributions plotted in Fig. 6(a) show that the HIG-2XP and NASCA codes give results in better agreement with experiment as far as the prediction of separation extent is concerned. However, the pressure obtained on the flare is higher than that obtained on the mesh G1 and in the experiment. The figure indicates that the HIG-2XP code is in good agreement with experimental values. The FLU3M code gives nearly the same results as those obtained on the mesh G1 as far as the point of separation is concerned, while the pressure obtained on the flare is lower (after the reattachment point) than the one obtained on the grid G1.
250
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
Fig. 6. Comparison of wall quantities calculated on the cylinder flare model with grid G2. (a) Pressure coefficient, (b) skin friction coefficient, (c) Stanton number. Table 2 Separation extent (mesh G2 for Navier–Stokes codes) Abscissa (X/L) Separation Reattachment
FLU3M
NASCA
HIG-2XP
FLOW
DSMC
Experiment
0.72 1.33
0.74 1.33
0.77 1.32
0.74 1.33
0.75 1.32
0.76 ± 0.01 1.34 ± 0.015
In the vicinity of reattachment, the experimental results present a non-smooth behaviour. It was not possible to ascertain whether it is due to experimental uncertainties or to a real physical phenomenon (that neither of the codes capture). However, it is interesting to point out that the major discrepancies between codes in the prediction of the wall pressure on the flare start in the reattachment zone. The examination of the skin friction coefficient distributions (see Fig. 6(b)) shows a good agreement between the HIG-2XP results and the experimental separation point. The locations of the separation and reattachment points are summarised in Table 2.
Concerning the Stanton number distributions (see Fig. 6(c)), the NASCA and HIG-2XP codes give results on the flare showing the best agreement with experiment. The FLU3M code confirms a weaker sensitivity to the mesh. The classical method used to check the spatial convergence by doubling the cells number without modifying the grid geometry does not provide a sufficient condition for convergence. External flow properties. Since the three finite-volume Navier–Stokes codes (FLU3M, NASCA and HIG-2XP) give nearly the same results with the grid G2, only three codes are considered: the finite-volume Navier–Stokes code
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
251
Fig. 7. Density profiles in the interaction region of the cylider flare model (X-ray EBF measurements). (a) Station X/L = 0.30, (b) station X/L = 0.60, (c) station X/L = 0.76.
NASCA, the finite-element Navier–Stokes code FLOW and the DSMC code. Three density profiles have been measured. The profile X/L = 0.3, shown in Fig. 7(a), is located upstream of the separation line. At this station, the increase of density is due to the shock generated by the leading edge. There is a good agreement between numerical and experimental results for the density peak amplitude. However, the calculated radial shock position varies with the simulation used. At this station (X/L = 0.3), leading-edge effects probably subsist, their influence being not accurately represented by the Navier–Stokes codes. This can be due to slip effects that are not taken into account and the difficulty to keep the correct leading edge location with the mesh in Navier– Stokes approach. The DSMC calculation is in excellent agreement with experiment. For the profiles at X/L = 0.6, shown in Fig. 7(b), the inverse behaviour appears on the radial shock location, the best predictions being furnished by the two Navier–Stokes codes. However, this difference on shock location needs to be confirmed before concluding on
this subject. The profiles at X/L = 0.76, shown in Fig. 7(c), present the same tendencies. 3.2.5. Concluding remarks In the present operation, the performance of different methodologies for shock wave/boundary layer interactions prediction were analysed. The performed experiment constitutes a precious entirely laminar test-case for the validation of codes (ReL = 18,375). The axisymmetric geometry guarantees the absence of side effects, while the very low Reynolds number of the experiment guarantees the absence of a transitional zone which is often associated with boundary layer reattachment. Different codes have been considered: three Navier– Stokes codes using a finite-volume approach, a Navier– Stokes code using the finite-element method, and a DSMC code using the VHS model. With regard to wall heat flux and wall pressure results, a good general agreement was found between calculations and experiment. Some discrepancies concerning the pressure results have been pointed out. A new
252
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
pressure measurements campaign will be performed at ONERA to check the current results, this being an example of experiment validation by calculation! Since comparisons with wall properties is only a step in the validation process, this experimental test-case has been completed by flow field measurements thanks to an X-ray detection technique. Measurements of gas density were successfully accomplished across shocks and boundary layers and a good agreement was found between calculations and experiment. New measurements campaigns will also include flow field measurements based on the dual lines coherent Anti-Stokes Raman scattering (DLCARS) technique for the determination of temperature, velocity and density profiles.
Fig. 8. Transonic shock wave/boundary layer interaction control. Test arrangement in the S8Ch wind tunnel.
3.3. Shock wave/boundary layer interaction control in transonic flow Aim of the operation. Shock-waves and their interaction with the boundary-layer play a major role in determining the performance of propulsion systems such as air-intakes, diffusers, turbomachines cascades, etc. One way to reduce the harmful effects of these shocks is to perform a control action in the interaction region [12,29]. In the present study, the following techniques have been considered: (1) active control in which a part of the boundary-layer is sucked off through a slot, (2) passive control, (3) hybrid control which is a combination of a passive control cavity and a suction slot (or cavity) located downstream of it. The aim of slot suction is to swallow part of the low energetic flow close to the surface before interaction of the boundary layer with the shock or during the interaction process itself. The principle of passive control consists in establishing a natural circulation between the downstream high pressure face of a shock and its upstream low pressure face. This circulation is achieved through a closed cavity, placed underneath the shock foot region, the face in contact with the outer flow being made of a perforated plate. It has been shown that, in very limited circumstances, passive control may produce a reduction of an airfoil drag, while postponing to higher incidences the limit of buffet onset. However, the gain being frequently problematic, it has been proposed to combine passive and active control to realise what is called hybrid control. Control methods are more likely to be used in airintake applications to diminish stagnation pressure losses (improved efficiency) or stabilise a shock-wave. The experimental part. These experiments were executed in the S8Ch transonic-supersonic basic research wind tunnel of the ONERA Meudon Centre. This facility is a continuous wind tunnel supplied with desiccated atmospheric air mainly dedicated to LDV measurements. The stagnation conditions were: pst0 = 96,000 ± 800 Pa and Tst0 = 300 ± 4 K. A photograph of the test set-up is shown in Fig. 8. It is constituted by a transonic channel having a test section with a maximum height of 100 mm and a span of 120 mm. The lower wall is rectilinear and equipped to receive the control
Fig. 9. Transonic shock wave/boundary layer interaction control arrangement for passive control.
devices, the upper wall being a contoured profile designed to produce a uniform supersonic flow of nominal Mach number equal to 1.4. A second throat, of adjustable cross section, is placed in the test section outlet to produce by choking effect a shock-wave whose position, and hence intensity, can be adjusted in a continuous and precise manner. It also isolates the flow field from pressure perturbations emanating from downstream ducts, reducing unwanted shock oscillations. The two side walls are equipped with high quality glass windows to allow visualisations and LDV measurements. The type of control taken as example here is passive control. For this device (see Fig. 9), a 70 mm-long passive control cavity was used, the shock being centred on it. The cavity, which extends between X = 130 mm and X = 200 mm, is covered by a perforated plate, with a 5.67%porosity and 0.3 mm-diameter holes (X origin at the nozzle throat). The flows under study was qualified by schlieren visualisations and quantified by measurements of wall pressure distributions (pressure orifices being located in the vertical median plane of the test set-up) and probing of instantaneous velocity with a two-component LDV system, in the median plane too [7]. The field quantities are given with an accuracy mainly depending on uncertainties affecting the LDV sys-
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
tem calibration and the statistical treatment of the sample of the instantaneous velocity components. For the present experiment, these uncertainties are: (1) 1% of the maximum velocity modulus (i.e., the velocity U0 of the upstream external flow, equal to 384 m/s and taken as reference) for the mean velocity components, (2) less than or equal to 8% of the maximum normal stress for the normal stress components of the Reynolds tensor, and (3) less than or equal to 10% of the maximum shear stress for the turbulent shear stress component. In the S8Ch wind tunnel, for this class of weakly separated flows, a two-dimensional zone was evidenced by probing in planes located at ±20 mm from the median plane, in the interaction region [6]. A Navier–Stokes code confrontation. After the preceding example of simulations of purely laminar flows, mandatory first step in any code validation, with the present one, we tackle the important question of the boundary conditions in modelling turbulent flows. The numerical simulations were performed with the NASCA code which solves the classical Reynolds Averaged Navier–Stokes (RANS) equations. Turbulence modelling was first carried out by means of the [k–ε] transport equation models of Chien (Ch model) [11] and Launder–Sharma (LS model) [19]. The Ch model is calibrated for boundary layers in simple geometry where the distance from the wall is simple to establish, like in the present case. The LS model is the result of a generalisation of source gradient terms appearing in the well known Jones–Launder model and of a re-optimisation of its constants. These two models were compared in the reference and passive control cases with the new [k–σ ] turbulence model, where σ represents a length scale, tied to the turbulent kinetic energy k and its dissipation ε by the defining relation: σ = Cµ k 3/2 /ε, Cµ being the constant multiplying the viscosity in the [k–ε] models [3]. The calculation domain (see Fig. 10) is a part of the experimental channel extending from a far upstream section of the divergent expanding zone (at X = 60 mm), where experimental velocity and turbulent shear stress profiles are imposed to give well defined supersonic inflow conditions, to the end of the channel, where the experimental pressure is imposed. The solution adopted here for testing turbulence models accurately consists in calculating ε or σ on the inflow boundary in such a way that experimental and modelled cross correlation of fluctuations be equal. The evaluation of k from the experiment that is necessary for this purpose necessitates to calculate the third velocity fluctuation correlation component from the two effectively measured longitudinal ones u 2 and v 2 . For this, we use the classical approximation: w 2 ≈ 12 (u 2 + v 2 ) which is questionable for the totality of the flow, so, for the validations, we have focused our attention on comparisons with cross correlation values. The second fundamental point is to impose the experiment-fitted values on the entire vertical width of the channel. The calculation domain thus includes both upper and lower walls. The tests of the turbulence models are done by comparison with
253
representative experimental field values obtained in a window surrounding the interaction zone of the lower wall. At the downstream boundary, situated upstream of the convergent of the second throat, at X = 380 mm, a numerical procedure, which is an adaptation of the characteristic extrapolation method [34], is developed for fixing the static pressure level. In this downstream attached flow, the transverse static pressure variations are neglected. A complete experimental characterisation of these small pressure variations would be difficult in this subsonic zone. A three dimensional calculation should therefore include the totality of the second throat. Passive control is simulated by prescribing the unit mass flow ρv at the wall, the conditions on the other variables remaining unchanged. The ρv value at the wall is obtained by relations expressing a direct dependence of the wall vertical velocity to the pressure difference between the cavity and the external flow. The relation used for the computations is the calibration law [23]. To establish this law, a flow model is developed for an ideal hole. This means that the hole has a cylindrical shape and the flow in it is laminar, incompressible and pipe like. The cavity pressure is fixed to its measured experimental value: 56110 ± 500 Pa. The Poll wall vertical velocity law accounts for a pure dependence of the mass flow rate to the pressure difference between the cavity and the external flow [4]. One can summarise this fact by a relation: ρv = f (p). The streamwise velocity at the wall being constant and equal to zero, the continuity equation takes the form: ∂(ρv)/∂Y = 0 which implies: ∂p/∂Y = 0. The additional hypothesis of thermal equilibrium at the wall implies the classically used relations: ∂ρ/∂Y = 0, ρu = 0, ∂(ρe)/∂Y = 0, with the simple modification: ρv = f (p). The sign of the vertical velocity at the wall, positive or negative following the case of blowing or suction, is determined by the sign of the pressure difference (pc − p) between the cavity and the outer flow. The mean value of ρv being the same at both ends of the hole, the formulae can be applied at one or the other side of the perforated plate. The (178 × 361) grid used (see Fig. 10) has been retained after a study of space convergence of the calculation. The LS model, which is significantly more sensitive to the mesh than the two other ones, has necessitated the use of this particularly fine grid. Major attention has been paid to mesh + refinement at the wall. √ The value of the reduced distance Y , ρw τw equal to (Y − Yw ) µw , τ being the laminar shear stress and subscript w designating values at the wall, is 0.57 for the upstream cell close to the wall. The X-distribution of the cell length in the X-direction is shown on the top of Fig. 10. The search for an optimal grid has consisted in successive calculations with various refinements until space convergence in the viscous zones, where the gradients due to the shock spread rapidly. The best way for illustrating our final test for space convergence is to show the coincidence between the wall skin friction distributions obtained with a (178 × 381) and a (355 × 361) grids (see Fig. 11). This complete grid optimisation has been performed in the
254
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
Fig. 10. Transonic shock wave/boundary layer interaction with passive control. Computational grid (178 × 361).
Fig. 11. Transonic shock wave/boundary layer interaction with passive control. Skin friction coefficient evolution on the lower wall – reference case (no control).
reference case, without control of the interaction [4]. The quality of the resulting (178 × 381) grid has been verified, with the same results, in the case of a controlled interaction with the single LS model. A search for a good resolution of the shock in the non viscous zones, imposing considerable grid refinements and costly calculations, was beyond the scope of this study. The Mach number contours, plotted in Fig. 12, are obtained with the [k–σ ] model, the three models giving nearly identical pictures at this level. These contours show a good prediction of the large spreading of the shock system,
which begins at the origin of the cavity. The difference in the location of the crossing (quasi normal) shock is weak. The thickening of the viscous zone due to the control device is well reproduced by the calculation. Computed wall pressure distributions (see Fig. 13) show the difficulty to simulate injection through very small holes with a continuous approximation made on a discrete mesh. An entirely rigorous calculation of this problem should have been done by meshing each hole, which is unrealistic with present computing capabilities. We will see in the following results that the apparently rough approximation
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
255
Fig. 12. Transonic shock wave/boundary layer interaction with passive control. Mach number contours: confrontation between experiment and Navier–Stokes calculation.
Fig. 13. Transonic shock wave/boundary layer interaction with passive control. Computed and experimental wall pressure distribution.
made to treat this porous wall condition is almost correct. At the beginning of the perforated plate (see Fig. 13), the peak appearing on the computed wall pressure values at X = 130 mm is a consequence of the sudden change of boundary condition between the two surrounding mesh points. The numerical approximation and the experimental resolution are not sufficient to give an account of the
true physical process. Downstream of the mid cavity, the computed pressure recovers more satisfactory levels. The location where the shock takes place is determined by the pressure levels in the subsonic downstream flow. These pressure levels are themselves function of the effective walls formed by the boundary-layer displacement thickness that delimit the equivalent perfect fluid flow. The longitudinal evolution of this thickness depends on the models. We see therefore that they will influence strongly the location of the interaction zone. A correct prediction of this location was obtained only with the [k–σ ] model. In the case of the [k–ε] calculations, it was predicted more than 20 mm upstream. Comparisons in the zone of interaction between the two types of models give, in these conditions, very large errors in the case of the [k–ε]. This does not give significant information on their true global quality. So, we have preferred to force a downstream pressure level lower than the measured one to adjust the shock position calculated by the [k–ε] models at the right place. The consequence of this choice, for a given model, is only a downstream translation of the wall pressure step of Fig. 13, without other modification. The velocity profiles upstream of interaction being practically invariant, the fields in the interacting viscous zones stay the same too. The price of this adjustment of the shock is that the [k–ε] models are submitted to a pressure step lower than the experimental
256
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
Fig. 14. Transonic shock wave/boundary layer interaction with passive control. Streamwise velocity profiles in the control region.
Fig. 15. Transonic shock wave/boundary layer interaction with passive control. Transverse velocity profiles in the control region.
one between X = 170 and 280 mm. The amplitude of this shifting of pressure can be seen in Fig. 13. Longitudinal velocity levels in the area of maximum wall pressure gradient (see Fig. 14), at the beginning of the cavity (X = 140 mm), has been submitted to the effect of the oblique thin compression fan starting from the boundary between the solid and porous walls. The crossing of this compression fan, which affects the experimental longitudinal velocity profiles, can be seen in Fig. 14 between Y = 15 mm and Y = 20 mm. The first mesh point, at which the approximated porous boundary condition is applied, corresponds to the beginning of the porous plate (X = 130 mm). For the perforated plate, the first row of holes is located at 1.2 mm after the beginning of the plate, which corresponds to the distance between each row of holes. The calculation evaluates a mean normal flux continuously distributed on the plate surface. This approximation could affect weakly the main features of the global interaction system, like the observed shift between the computed and experimental traces of the compression fan crossing the external part of the velocity profiles. The evolution of the streamwise boundarylayer velocities in this area of interaction is predicted satisfactorily by the models. An effect of passive control on the
experimentally observed interaction is the occurrence of a separated flow of small size at the level of the “λ”-shock system, above the perforated plate (see Fig. 12). No recirculating flow appears in the LS calculation, a small recirculating bubble appears between X = 151 mm and 167 mm with the Ch model. The results of the [k–σ ] calculation show a bubble between X = 145 mm and X = 194 mm. This last recirculating zone extends vertically up to Y = 0.3 mm and coincides partly with an area located between X = 170 mm and X = 200 mm, where small negative values of the longitudinal velocity have been measured at Y = 0.4 mm. In any case, the calculated reverse longitudinal velocities are located too close to the wall to be compared to experimental values. More measurement points would be needed to correctly define the reversed flow region, and thus to validate turbulence models in this separated region. As a consequence of this defective capture of the reversed flow, the response of the boundary-layer flow to local compressions is too roughly simulated by the models. This fact can be observed in Fig. 14, at X = 180 mm and X = 220 mm. An important test for the validity of injection modelling at the wall is the prediction of near-wall vertical velocity profiles. As a preliminary verification, the calculated resid-
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
257
Fig. 16. Transonic shock wave/boundary layer interaction with passive control. Cross-correlation profiles in the control region.
ual mass flow rate per unit span across the perforated plate, which must approach the null experimental value, has been expressed as a ratio of the mass flow rate deficit ρ0 U0 δ0∗ in the upstream boundary-layer. In all the cases, the value of this ratio has been found to be lower than 10−5 . More detailed information is presented in Fig. 15. Upstream of the compression system (X = 140 mm), wall injection is predicted by all the models. Data at this station are compared with those obtained without control, these last data make apparent the error bar on the vertical velocity profiles. We verify indeed that these near wall values without injection, which must tend to zero, are affected of the 1% of U0 error defined before. When injection is performed, the main global effect is the occurrence of a strong vertical gradient of normal velocity. This induces a supplementary decrease of the longitudinal momentum distribution just before interaction with the shock. In the external part of the boundarylayer, global discrepancies of the non viscous flow have a more spectacular effect on these small velocities, in particular, the location of the calculated crossing expansion fan is about 5 mm higher than the experimental one. Downstream of compression (X = 180 mm), a wall suction is found by all the models, as imposed by physics and confirmed by experiment. The computed shear stresses (calculated in Favre means: ρu v /ρ) ¯ are compared to the measured Reynolds tensor components u v . It has been established that the difference between these two quantities is small for an adiabatic, flat plate boundary layer, when the Mach number varies from 0 to 5 [20]. In the transonic cases, the local relative discrepancy is not greater than 2%. At the interaction beginning (X = 140 mm), the LS and [k–σ ] models anticipate the growth of the maximum u v level, which is not the case with the Ch model (see Fig. 16). The LS and [k–σ ] models give close results at the beginning of interaction and in the separated flow (X = 140 and 180 mm) for this controlled interaction, the Ch model strongly over-predicts the maximum level of u v . Downstream of the interaction area, the rather good agreement between the LS and [k–σ ] models is confirmed.
The fact that the interaction starting point is fixed at the beginning of the porous wall facilitate the agreement between models. The global shift observed between the levels predicted by LS and [k–σ ] models on one hand and by the Ch model on the other, could be partly explained by the modification of the downstream pressure level, which implies that the profiles, particularly in the case of the Ch model, have been submitted to lower decelerations. This effect influences backwards the area of interaction via the thickened subsonic layers. After reattachment the level of cross correlation in the Ch model passes on the first grid point from 0 to a value close to its maximum (see Fig. 16, X = 260 mm). We suspect that, in this modified boundary layer, the asymptotic behaviour of the model is not in accordance with the wall value that must be 0. Such a problem, if confirmed by other cases, should have to be analysed theoretically. This kind of study is not in the scope of this work. 3.4. A supersonic base flow Aim of the operation. The main goal of the present study was to assess the ability of the [k–σ ] two equation turbulence model to predict the mean field and some fluctuating quantities in a supersonic axisymmetric base flow [3]. The chosen test case is an experiment executed by Herrin and Dutton [18] whose results are widely accepted for testing simulations and which has served as data base for previous theoretical studies on the subject [14,15,25,31]. Assessment of the model, using the NASCA code, was done by comparisons with these data and with the results given by three other well known models. The first two are the [k–ε] model of Launder and Sharma and a version of the [k–ε] model whose various constant and source terms are defined according to renormalisation group methods, named [k–ε] RNG model [33]. This last model (which we call RNG) has been improved primitively in backward facing step flows, reason why it was used for the case of base flow [28]. The [k–ω] model [32] has also been tested under the form it takes when expressed in σ variable, to give the so called [k–σ ]ω model.
258
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
Fig. 17. Supersonic base flow at Mach = 2.45. Calculation grid and overall flow organisation (streamlines in solid lines, constant pressure levels in dashed lines).
The reason for this variable change is that the wall condition for ω is replaced by the simple and well defined condition σ = 0. A modified version of the [k–σ ], named [k–σ ]2 , obtained by modifying just one constant to test its sensitivity was also compared. The flow field physics. The experiments of Herrin and Dutton were performed in a test facility specially designed to generate axisymmetric flows [26]. In particular, profiles of mean and fluctuating velocity fields in the recirculating flow are provided by 3D LDV measurements. The upstream Mach number has been determined from LDV [18] to be 2.46 ± 1%, while the static wall pressure measured just upstream of the base corner corresponds to a uniform flow with a Mach number of 2.44. During the present calculations we have chosen an upstream Mach number M0 = 2.45. A Navier–Stokes code confrontation. The present calculations have been performed by using the NASCA research code. The near-base part of the rectangular 201×281 mesh and an example of calculated streamline pattern obtained with the [k–σ ] model are both shown in Fig. 17. The rapid variations of the fields and the turbulent mixing in the shear layer forming at the base shoulder are the first challenge for the modelling of the mean and fluctuating fields. Good precision in the prediction of the nearly constant pressure in the recirculating bubble, limited by the mixing layer and the rear stagnation point on the axis, is fundamental for base drag prediction (constant pressure levels are shown in dashed lines in Fig. 17). The recompression subsequent to the flow field realignment at the rear stagnation point is vis-
ible in Fig. 17, at the point where the recirculating flow progressively changes into a wake. The mesh has been defined by careful study of the spatial convergence in which the optimum refining at the walls was determined, especially the transverse grid distribution on the base. The criterion for spatial convergence was to stabilise the solutions obtained with different grid resolutions within a 1% margin for the most sensitive variables which were the pressure and the shear stress in the reversed flow region. The numerical study with each turbulence model was performed on this mesh with typical CPU times of 20–30 minutes on the ONERA NEC SX5 computer. The experiment of Herrin and Dutton provides us with data on the evolution of the boundary layer on a restricted part of the body beginning at a distance X = −0.079D from the base (D is the cylindrical body diameter). The experimental profiles at this location are taken as upstream boundary condition for the computations. Profiles of ε or σ are calculated on this upstream boundary from the experimental values of k and u v by the same principle and methods as those of the preceding test case. The distributions of the calculated and experimental pressure distributions on the base are shown in Fig. 18. The experimental value of the ratio pc /p∞ of the mean base pressure to the uniform upstream pressure is equal to 0.55. This level is 10% higher than data on base flows at the same Mach number extracted from earlier experimental compilations [13]. These computed base pressure values are considerably closer to experiment than that obtained with a Baldwin–Lomax calculation (see Fig. 18), in agreement with previous studies [25]. However, they contradict over-
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
Fig. 18. Supersonic base flow at Mach = 2.45. Pressure distribution on the model base.
259
predictive results obtained with other two-equation models [14]. In a recent computation using LES [15], a level of 0.52 for the mean base pressure was found. Profiles obtained in the middle of the bubble, at X = 0.63D downstream of the base, are shown in Fig. 19. The velocity profiles in the mixing layer are well predicted by the [k−σ ]ω model. The evolution of k is qualitatively similar to the one given by the other models with a globally lower level which is satisfactory on the axis but unrealistic in the mixing layer. A more correct prediction of the shearstress characterises both the [k–ε] and [k–σ ] models. The coincidence between the [k–σ ] and RNG models, on the one hand, and between the [k−σ ]2 and [k−σ ]ω models on the other hand, appears clearly on the k and u v profiles. With the two [k–ε] models, non-realistic positive values of the Reynolds stress appear in the expansion fan emanating from the base corner. A bump in turbulent kinetic energy profiles is also produced at the upper boundary of the mixing layer, the level of these non-existing fluctuations being particularly important in the LS model calculation. Prediction of the maximum shear stress levels by the LS model in the mixing layer is better but accompanied by greater levels of kinetic energy than with the RNG and [k–σ ] models. A salient fact is the smoothing by the models of the shear stress evolution during the rapid transition between the mixing layer and the
Fig. 19. Supersonic base flow at Mach = 2.45. Profiles at location X/D = 0.63.
260
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
Fig. 20. Supersonic base flow at Mach = 2.45. Profiles at location X/D = 1.26.
reversed flow. This smoothing explains the incapacity of the models to predict the nearly constant reversed flow region in the immediate proximity of the base, the defect coming from an over-prediction of the negative axial velocity on the axis. This too large radial variation of the axial mean velocity is also due to an over-prediction of the eddy-viscosity in the part of the bubble situated near the axis. The results obtained by the RNG and [k–σ ] models are satisfactory at the rear stagnation point, located at X = 1.26D (see Fig. 20), except for the k profile. The good agreement of the axial velocity profile with experiment proves that the position of this stagnation point is well predicted by the two models, this result being complemented by a correct calculation of u v . One notes that the [k−σ ]ω and [k−σ ]2 models overestimate the separated bubble length. This deficiency is also observed in the calculation of the Reynolds stresses; however the profile of k obtained with these models is still satisfactory. Pathological behaviours of the [k–ε] models are confirmed in shear-stress and kinetic energy profiles at the level of the compression and expansion fans. The LS model over-predicts significantly the turbulent quantities. A local maximum of radial velocity is generally predicted at a distances from the axis of approximately 0.45D, the experimental bump being situated below 0.4D. This can be a
consequence of an anticipation of the calculated position of the origin of the compression fan with respect to experiment, this compression indeed more affecting the transverse flow.
4. Conclusion The spectacular increase in our computing capacity during the past 40 years led to a certain despise of the experimental activity. It was anticipated that “numerical wind tunnels” will soon replace the noisy, difficult to operate, dangerous and costly real wind tunnels. This is not our purpose to enter into this polemics. It would not be wise, by reaction, to despise the computational activity which has taken a considerable place in the design and development of nearly all the industrial products (and in many other sectors too!). Because of the technical and scientific difficulties encountered in the domain and the necessity to rely on safer methods, the aerospace industry has strongly invested in the development of codes since the beginning of the CFD era. Many of the numerical schemes have been devised to improve calculation methods for fluid mechanics applications, notably in the field of aerodynamics.
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
However, the development of more and more performant codes has not completely killed the experimental activity. Rather soon, it appeared that the confidence in code predictions was limited, the cause of the shortcomings being due partly to uncertainties in the numerical handling of the equations, partly to the lack of accuracy and representativity of the physical laws – or models – implemented in the codes. This perception of the reality has motivated a renewal of interest for experiments since only the confrontation with experimental data can validate – or invalidate – a code. To validate their codes, numericians need an as complete as possible information on some representative test cases. This information constitutes what is called a data bank which must respect certain rules to be useful. Thus, the data bank must contain a precise description of the configuration, along with all the necessary flow and boundary conditions. The measurements must be considered as safe, and if possible accurate. A great accuracy is not always mandatory (it costs much money), but uncertainty margins must be given. Constitution of a data bank is not a straightforward operation. In addition of technical skill to fabricate a test set up, to operate the wind tunnel and execute the experiment, to perform the measurements, it requires a solid background in fundamental fluid mechanics. The data bank constitution is not limited to the acquisition of a vast amount of results, but must be accompanied by an in depth analysis of the flow physics. Because of the investment needed by such operations and their strategic importance for the development of predictive methods, the question of the data bank dissemination inevitably arises. It is now realised that a good data bank can be as precious as a code and cannot be freely transmitted. Even basic experiments have now an economic weight and cannot be put on the market without something in exchange. Thus, dissemination rules have to be more precisely defined according to the more or less precious nature of the data bank treasure. In addition of the permanent scientific concern about more accurate, safer and less expansive predictive methods, the problem of the constitution of valuable, safe, well identified and permanent data banks is now considered as a strategic issue and addressed seriously. In this perspective, the ONERA Fluid Mechanics and Energetic Branch has started the constitution of a data bank containing the most prominent experimental results obtained in its research wind tunnels of the Chalais-Meudon Center over the last 30 years [1]. This task will be actively pursued and the data bank contents fed with new experiments satisfying the quality criteria here above defined.
Acknowledgements The results on shock wave/boundary layer interaction control have been obtained in the framework of the EUROSHOCK I and II programmes of the European Union.
261
References [1] R. Benay, La base de données du DAFE. Mise à jour 2001 (The data bank of the Fundamental/Experimental Aerodynamics Department. Update 2001), ONERA Rapport Technique N◦ RT 3/03589 DAFE, juillet 2001. [2] R. Benay, P. Servel, Applications d’un code Navier–Stokes au calcul d’écoulements d’arrière-corps de missiles ou d’avions (Application of a Navier–Stokes code to the calculation of flow past missile and aircraft afterbodies), La Recherche Aérospatiale 6 (1995) 405–426. [3] R. Benay, P. Servel, Two-equation [k–σ ] turbulence model: Application to a supersonic base flow, AIAA J. 39 (3) (2001) 407–416. [4] R. Benay, P. Berthouze, R. Bur, Modeling of controlled shock wave/ boundary layer interaction in transonic channel flow, AIAA J. 39 (12) (2001) 2293–2301. [5] G. Bird, Molecular Gas Dynamics and the Direct Simulation of Gas Flows, Clarendon Press, Oxford, 1994. [6] R. Bur, Fundamental study on passive control of shock-wave/turbulent boundary-layer interaction in transonic flow, Ph.D. Thesis, University of Paris VI, March 1991. [7] R. Bur, B. Corbel, J. Délery, Study of passive control in a transonic shock wave/boundary layer interaction, AIAA J. 36 (3) (1998) 394– 400. [8] M. Borrel, J.-L. Montagne, J. Diet, Ph. Guillen, J. Lordon, Calculs d’écoulements supersoniques autour de missiles tactiques à l’aide d’un schéma décentré (Calculation of the supersonic flows past tactical missiles with a centered scheme), La Recherche Aérospatiale (1988) 2–20. [9] B. Chanetz, R. Benay, J.-M. Bousquet, R. Bur, T. Pot, F. Grasso, J. Moss, Experimental and numerical study of the laminar separation in hypersonic flow, Aerospace Science and Technology 3 (1998) 205– 218. [10] B. Chanetz, J.-P. Davis, M. Lefebvre, M. Leplat, J. Perraud, Hot Measurements Technique phase II: synthesis of the ONERA results, ONERA RT 1/00748 DPRS/DAFE, July 2000. [11] K.Y. Chien, Prediction of channel and boundary-layer flows with a low-Reynolds-number turbulence model, AIAA J. 20 (1) (1982) 33– 38. [12] J. Délery, Shock-wave/turbulent boundary-layer interaction and its control, Prog. Aerospace Science 22 (4) (1985) 209–280, also ONERA TP-1984-27, 1984. [13] J. Délery, M. Sirieix, Ecoulements de culot (Base flows), AGARD LS98 on Missile Aerodynamics, Ankara, Rome, Bruxelles, 5–16 March 1979, ONERA T.P. n◦ 1979-14F. [14] P.I. Espina, U. Piomelli, Validation of the NPARC code in supersonic base flows, AIAA Paper 97-0032, January 1997. [15] C. Fureby, Y. Nilsson, K. Andersson, Large eddy simulation of supersonic base flow, AIAA Paper 99-0426, 1999. [16] N. Gorchakova, B. Chanetz, T. Pot, R. Bur, J.-P. Taran, D. Pigache, J. Moss, D. Schulte, L. Kuznetsov, V. Yarygin, Progress in shock waveboundary layer interaction studies in rarefied hypersonic flows using electron-beam-excited X-ray detection, AIAA J. 40 (4) (2002) 593– 598, ONERA TP n◦ 2002-76. [17] F. Grasso, M. Marini, Analysis of hypersonic shock wave/boundary layer interaction phenomena, Computers and Fluids 25 (1996) 1–21. [18] J.L. Herrin, J.C. Dutton, Supersonic base flow experiments in the nearwake of a cylindrical afterbody, AIAA Paper 93-2924, 1993. [19] B.E. Launder, B.I. Sharma, Application of the energy dissipation model of turbulence to the calculation of flows near a spinning disc, Heat and Mass Transfer 1 (1974) 131–138. [20] G. Maise, H. McDonald, Mixing length and kinematic eddy viscosity in a compressible boundary layer, AIAA J. 6 (1) (1968) 73–80. [21] J. Moss, J. Olejniczack, Shock wave/boundary layer interactions in hypersonic low density flows, AIAA Paper 98-2668, 1998. [22] S. Osher, S. Chakravarthy, Very high order accurate TVD schemes, ICASE Report n◦ 84-44, September 1984.
262
R. Benay et al. / Aerospace Science and Technology 7 (2003) 239–262
[23] D.I.A. Poll, M. Danks, B.E. Humphreys, The aerodynamic performance of laser drilled sheets, in: Proceedings of the First European Forum on Laminar Flow Technology Hamburg, Germany, 1992, pp. 274– 277. [24] P.J. Roache, Verification of codes and calculations, AIAA J. 36 (5) (1998) 696–702. [25] J. Sahu, Numerical computations of supersonic base flow with special emphasis on turbulence modeling, AIAA Paper 92-4352, 1992. [26] J.M. Sauter, J.C Dutton, Design of an axisymmetric supersonic wind tunnel and experimental study of supersonic, power-off base flow phenomena, Dept. of Mech. and Ind. Eng., Univ. of Illinois at UrbanaChampaign, UILU ENG 89-4002, Urbana, IL, March 1989. [27] D. Schulte, A. Henckels, U. Wepler, Reduction of shock-induced boundary layer separation in hypersonic inlets using bleed, Aerospace Science and Technology 2 (4) (1998). [28] C.G.Speziale, S. Thangam, Analysis of an RNG based turbulence model for separated flows, ICASE Report n◦ 92-3, January 1992.
[29] E. Stanewsky, J. Délery, J. Fulker, W. Geissler (Eds.), Euroshock – Drag Reduction by Passive Shock Control, Notes on Numerical Fluid Mechanics, Vol. 56, Vieweg, 1997. [30] J.L. Steger, R.F. Warming, Flux vector splitting of the inviscid gasdynamics equations with applications to finite difference methods, J. Comput. Phys. 40 (2) (1981). [31] P.K. Tucker, W. Shy, A numerical analysis of supersonic flow over an axisymmetric afterbody, AIAA Paper 93-2347, 1993. [32] D.C. Wilcox, Turbulence Modeling for CFD, Chapter 4, DCW Industries, La Canada, CA, 1993. [33] V. Yakhot, S.A. Orszag, Renormalisation group analysis of turbulence. I. Basic theory, J. Scientific Comput. 1 (1) (1986) 3–51. [34] H.C. Yee, R.M. Beam, R.F. Warming, Boundary approximations for implicit schemes for one-dimensional inviscid equations of gas dynamics, AIAA J. 20 (9) (1982) 1203–1211.
