C Koupper

T Poinsot

L Gicquel

F Duchaine

On the ability of characteristic boundary conditions to comply with radial equilibrium in turbomachinery simulations

Setting up outlet boundary conditions in configurations which have a strong rotating motion is a crucial question for turbomachinery simulations. This is usually done using the so called 'radial equilibrium assumption' (REA) which is used before the simulation and provides an approximate expression for the pressure profile to impose in the outlet plane. This paper shows that recent methods developed for compressible flows, based on characteristic methods, including the effects of transverse terms, can capture radial equilibrium naturally without having to impose a pre-computed pressure profile. In addition, these methods are also designed to control acoustic reflections on boundaries and the present work suggests that they could replace classical REA approximations when non reflecting boundary conditions are required at the outlet of a turbomachine simulation, for example in LES. This is demonstrated in two cases: (1) a simple annulus flow with swirl imposed at the inlet and (2) a transonic turbine vane.

Nomenclature

⃗ c

Velocity

I. Introduction

Among all the difficulties encountered in the numerical simulations of fluid dynamics, prescribing adequate Boundary Condition (BC) has always been an important point. It is well known in the Computational Fluid Dynamics (CFD) community that imposing a realistic, non-disturbing outlet boundary condition is complex, if not impossible in most cases. An ideal condition should have a weak influence on the upstream flow, preserve stability, and control acoustic waves reflection and/or dissipation. [START_REF] Yoo | Characteristic boundary conditions for direct simulations of turbulent counterflow flames[END_REF][START_REF] Yoo | Characteristic boundary conditions for simulations of compressible reacting flows with multidimensional, viscous and reaction effects[END_REF][START_REF] Poinsot | Boundary conditions for direct simulations of compressible viscous flows[END_REF][START_REF] Lodato | Three-dimensional boundary conditions for direct and large-eddy simulation of compressible viscous flows[END_REF][START_REF] Pirozzoli | Generalized characteristic relaxation boundary conditions for unsteady compressible flow simulations[END_REF] This necessity becomes particularly important for high fidelity simulations: Large Eddy Simulation (LES), and Direct Numerical Simulation (DNS). Indeed, in such approaches, all (or part) of the turbulence scales are directly resolved in a highly turbulent and unsteady context, and the boundary conditions must not create spurious reflections of acoustic waves inside the domain.

Satisfying all these conditions is difficult and in practice the most common method is simply to apply a Neumann boundary condition sufficiently far from the region of interest. The convective derivative normal to the boundary is enforced to zero (∂/∂n = 0), and the value of the variables at the boundary is extrapolated in the streamwise direction. [START_REF] Tu | Computational Fluid Dynamics: A Practical Approach[END_REF] This simple approach is well suited for almost uniform steady state flows: exit of a duct, plenum etc. However, if the flow features important non uniformities, or a strong gradient (due to a preferential fluid motion for example), then a uniform Neumann BC is not adequate. This is especially true for turbomachinery flows as high inhomogeneities are found at the stage exit: wake effect, unsteady flow bubbles etc. Moreover, the swirling motion of the mean flow created by the deviation in the vane or rotor creates a positive radial pressure gradient: the so-called radial equilibrium. Assuming that pressure is homogeneous in the outlet plane is impossible. The pressure gradient in the outlet section may be estimated (as it will be demonstrated later), and can be used as BC at the exit section of a turbine: it is the Radial Equilibrium Assumption (REA). This approach is common for Reynolds Averaged Navier Stokes (RANS) simulations of turbine flows. For example, the commercial code Fluent features a radial equilibrium option for outlet BC: the user specifies the pressure at the hub and the equations for radial equilibrium are integrated using flow properties obtained inside the domain. A radial profile of pressure is then prescribed as a Neumann BC. Negru et al. [START_REF] Negru | Computation of stress distribution in a Francis turbine runner induced by fluid flow[END_REF] use this technique to impose the pressure at the exit of a Francis turbine runner. In the same context of hydraulic turbines, Ferro et al. [START_REF] Ferro | Design and experimental validation of the inlet guide vane system of a mini hydraulic bulb-turbine[END_REF] simulated the inlet guide vanes of a hydraulic bulb turbine using the Fluent radial equilibrium option. Torresi et al., [START_REF] Torresi | Accurate numerical simulation of a high solidity Wells turbine[END_REF] and similarly Kim et al. [START_REF] Kim | Numerical investigation on the effect of blade sweep on the performance of Wells turbine[END_REF] modelled a high solidity Wells turbine used for oscillatory water column device (energy conversion), using radial equilibrium in the outlet section. Garg 11,[START_REF] Garg | Heat transfer on a film-cooled rotating blade[END_REF] studied the heat transfer of a rotating blades of various aero gas turbines, imposing in a RANS solver a radial equilibrium pressure distribution one true chord after the rotor blade trailing edge. A last example can be found for the RANS computation of the heat transfer of a low pressure vane by Solano et al. [START_REF] Solano | Aero-thermal investigation of a multi-splitter axial turbine[END_REF] where the measured pressure at hub (from experiments) was used to prescribe the radial equilibrium at the outlet section.

For all these studies, imposing an outlet pressure profile is consistent with radial equilibrium but leads to a totally reflecting outlet in terms of acoustic waves. For RANS simulations, which do not capture acoustic anyway, this is not an issue but it becomes a problem for LES. Indeed, in the context of compressible LES, Navier Stokes Characteristic Boundary Conditions (NSCBC) are commonly employed. [START_REF] Yoo | Characteristic boundary conditions for direct simulations of turbulent counterflow flames[END_REF][START_REF] Lodato | Three-dimensional boundary conditions for direct and large-eddy simulation of compressible viscous flows[END_REF][START_REF] Pirozzoli | Generalized characteristic relaxation boundary conditions for unsteady compressible flow simulations[END_REF][START_REF] Poinsot | Theoretical and numerical combustion -Third Edition[END_REF][START_REF] Baba-Ahmadi | Inlet Conditions for Large Eddy Simulation of Gas-Turbine Swirl Injectors[END_REF] At an outlet, they allow to control the influence of the incoming information (through wave amplitude) and thus both the pressure level and reflection coefficient of the BC. To perform LES of turbomachines, being able to satisfy both the REA and NSCBC condition is of primary importance. This requires to investigate how NSCBC conditions perform in flows with a strong rotation. It's shown here that the NSCBC formalism can let the physical radial pressure gradient establish naturally, so that NSCBC can be used at the outlet of turbines without any additional treatment. This represents a significant improvement over existing techniques because (1) it allows to capture the REA naturally, without a priori simplified evaluation of the pressure profile and (2) it brings the power of non reflecting BC in a field where they have been used only for a few cases.

Part II introduces the equation of the simplified radial equilibrium and its physical meaning. Then, the formalism of Navier Stokes Characteristic Boundary Condition is presented in part III. Finally, in the last part the ability of NSCBC to deal with radial equilibrium is assessed in two test cases: a simple annulus (IV.A) and an industrial high pressure turbine vane (IV.B).

II. The simplified radial equilibrium

The radial equilibrium equations were first derived in the 40-50s and published by L.H. Smith Jr. in 1966. [START_REF] Smith | The Radial-Equilibrium Equation of Turbomachinery[END_REF] He showed that the whirling motion of a fluid inside a turbomachine creates a centrifugal force that has to be balanced by a centripetal one: a positive radial pressure gradient establishes. [START_REF] Wu | Application of radial-equilibrium condition to axial-flow compressor and turbine design[END_REF] By definition it only applies to axial stations between blade rows where pseudo-equilibrium can be achieved. The study was limited to axisymmetric flows, which also implies that the discrete action of the blades is not taken into account.

The usual denomination of radial equilibrium may refer either to the Full (American) or Simplified (British) definition. In the former case, it is the complete radial momentum equation arranged in a form [which is] suitable for the determination of the flow field in a turbomachine. [START_REF] Smith | The Radial-Equilibrium Equation of Turbomachinery[END_REF] Coupled with a blade aerodynamic computation, it provides the balanced distribution of the flow properties from hub to casing between all vanes. In the simplified version (most common), the full equations are simplified assuming that the radial velocity is zero: c r = 0. In order to fulfill the latter hypothesis, the axial stations of interest must be far enough from any blade. This condition is met when the outlet section is placed sufficiently far from the vanes. Only the simplified radial equilibrium is studied here. Under these conditions, the flow is governed by the Euler equations for compressible flow with the equation of state for ideal gas. In cylindrical coordinates (r, θ, x), the momentum equation in radial direction is:

∂c r ∂t + c r ∂c r ∂r + c θ r ∂c r ∂θ - c 2 θ r + c x ∂c r ∂x = - 1 ρ ∂P ∂r (1)
Applying all the assumptions made above (c r = 0; ∂/∂θ = 0; ∂/∂t = 0) to Eq. ( 1) leads to the simplified radial equilibrium equation:

1 ρ ∂P ∂r = c 2 θ r ( 2
)
where P is the static pressure, c θ the azimuthal velocity component, and ρ the density.

III. Navier Stokes Characteristic Boundary Conditions

Non reflecting boundary conditions have been extensively studied in the past decades [START_REF] Yoo | Characteristic boundary conditions for direct simulations of turbulent counterflow flames[END_REF][START_REF] Yoo | Characteristic boundary conditions for simulations of compressible reacting flows with multidimensional, viscous and reaction effects[END_REF][START_REF] Poinsot | Boundary conditions for direct simulations of compressible viscous flows[END_REF][START_REF] Lodato | Three-dimensional boundary conditions for direct and large-eddy simulation of compressible viscous flows[END_REF][START_REF] Thompson | Time dependent boundary conditions for hyperbolic systems[END_REF][START_REF] Nicoud | Defining Wave Amplitude in Characteristic Boundary Conditions[END_REF][START_REF] Lodato | Optimal Inclusion of Transverse Effects in the NonReflecting Outflow Boundary Condition[END_REF] for LES or DNS of compressible flows. In 3D turbulent unsteady flows, wave reflections from the boundaries of the domain have to be carefully controlled. Characteristic boundary conditions have the difficult task to ensure numerical stability and minimal acoustic reflections in regions where the flow dynamic activity is high.

The characteristic form of the Navier Stokes equations written for NSCBC highlights the presence of waves crossing the boundary. [START_REF] Poinsot | Boundary conditions for direct simulations of compressible viscous flows[END_REF][START_REF] Thompson | Time dependent boundary conditions for hyperbolic systems[END_REF] The outgoing waves leaving the domain can be calculated from the interior points, while the acoustic incoming wave has to be imposed and includes all effects of the domain outside the computational box. Many formulations have been proposed since the 90's. Thomson 18 applied one-dimensional approximations of the characteristic boundary conditions for Euler equations. Poinsot and Lele 3 extended the formalism for Navier Stokes equations, by introducing the Locally One-Dimensional Inviscid (LODI) assumptions. The LODI formalism is efficient for flows aligned with the direction normal to the boundary. However, when the flow or the acoustic waves do not reach the boundary at a normal angle, the one-dimensional assumptions are too restrictive, and instabilities, reflections or non-physical distortions can appear. [START_REF] Yoo | Characteristic boundary conditions for simulations of compressible reacting flows with multidimensional, viscous and reaction effects[END_REF][START_REF] Lodato | Three-dimensional boundary conditions for direct and large-eddy simulation of compressible viscous flows[END_REF][START_REF] Lodato | Optimal Inclusion of Transverse Effects in the NonReflecting Outflow Boundary Condition[END_REF] Yoo et al. [START_REF] Yoo | Characteristic boundary conditions for direct simulations of turbulent counterflow flames[END_REF][START_REF] Yoo | Characteristic boundary conditions for simulations of compressible reacting flows with multidimensional, viscous and reaction effects[END_REF] proposed to include the transverse terms in the calculation of the incoming wave. This formalism is used here, and described below.

The three-dimensional compressible Navier-Stokes and total energy conservation equations for a monospecies non-reacting flow can be written into a characteristic form 1 in cartesian coordinates (x 1 , x 2 , x 3 ), without viscous and source terms:

∂ ∂t        c 1 c 2 c 3 ρ P        +        (L 5 -L 1 )/ρa 0 L 3 L 4 L 2 + (L 5 + L 1 )/a 2 0 L 5 + L 1        +        (T 5 -T 1 )/ρa 0 T 3 T 4 T 2 + (T 5 + T 1 )/a 2 0 T 5 + T 1        =        0 0 0 0 0        (3)
with a 0 = √ γRT the speed of sound, c the velocity, and ⃗ U = (c 1 , c 2 , c 3 , ρ, P ) the primitive variable vector. The amplitude vector of the characteristic waves ⃗ L can be developed as:

⃗ L =        L 1 L 2 L 3 L 4 L 5        =           1 2 λ 1 ( ∂P ∂x1 -ρa 0 ∂c1 ∂x1 ) λ 2 ( ∂ρ ∂x1 -1 a 2 0 ∂P ∂x1 ) λ 3 ∂c2 ∂x1 λ 4 ∂c3 ∂x1 1 2 λ 5 ( ∂P ∂x1 + ρa 0 ∂c1 ∂x1 )           (4)
The characteristic velocities λ i associated with the waves amplitudes L i are:

λ 1 = c 1 -a 0 (5) λ 2 = λ 3 = λ 4 = c 1 (6) λ 5 = c 1 + a 0 (7)
where λ 1 is the velocity of acoustic waves propagating in the negative direction (ingoing wave), while λ 5 refers to the outgoing acoustic wave. λ 2 = λ 3 = λ 4 is the convection velocity in the positive direction x 1 .

The transverse terms vector ⃗ T introduced by Yoo et al. [START_REF] Yoo | Characteristic boundary conditions for direct simulations of turbulent counterflow flames[END_REF] are:

⃗ T =        T 1 T 2 T 3 T 4 T 5        =         (⃗ c t • ∇ t P + γP ∇ t • ⃗ c t -ρa 0 ⃗ c t • ∇ t c 1 )/2 ∇ t • (ρ⃗ c t ) -(⃗ c t • ∇ t P + γP ∇ t • ⃗ c t )/a 2 0 (⃗ c t • ∇ t c 2 + 1 ρ ∂P ∂x2 ) (⃗ c t • ∇ t c 3 + 1 ρ ∂P ∂x3 ) (⃗ c t • ∇ t P + γP ∇ t • ⃗ c t + ρa 0 ⃗ c t • ∇ t c 1 )/2         (8)
with the subscript t referring to tangential directions x 2 , x 3 with the transverse terms vector ⃗ c t = (c 2 , c 3 ).

At a subsonic outflow, all waves exit the domain except the ingoing wave associated with the velocity λ 1 = c 1 -a 0 corresponding to L 1 . The value of the incoming wave L 1 can not be obtained from interior points and must account for the information propagating into the domain from the exterior. The extension of the LODI method 3 by Yoo et al. [START_REF] Yoo | Characteristic boundary conditions for direct simulations of turbulent counterflow flames[END_REF] leads to the following expression for L 1 :

L 1 = K(P -P t ) + (β -1)T 1 (9)
with K a relaxation coefficient. If the outlet pressure P is very different from the target pressure P t , the reflected waves going inwards will drive the pressure towards P t . A proposal for the relaxation coefficient was made by Rudy and Strikwerda: 21

K = σa 0 (1 -M 2 )/L x ( 10
)
where M is the maximum Mach number on the boundary, σ = 0.25 and L x a characteristic dimension in the axial direction. As K → ∞ the outlet boundary condition becomes fully reflecting, while for K = 0 the ingoing wave amplitude becomes zero: the boundary condition becomes non reflecting but the pressure can drift. Selle et al. [START_REF] Selle | The actual impedance of non-reflecting boundary conditions: implications for the computation of resonators[END_REF] suggest to use values of K corresponding to 0.1 < σ < π to avoid large reflections and possible errors on pressure evaluation. Yoo et al. [START_REF] Yoo | Characteristic boundary conditions for direct simulations of turbulent counterflow flames[END_REF] introduce a damping coefficient β ∈ [0, 1] which is applied to the transverse terms in Eq. ( 9). When β = 1, the transverse terms are ignored and the LODI formalism 3 is recovered. Different authors report that β should be equal to a reference Mach number. [START_REF] Yoo | Characteristic boundary conditions for simulations of compressible reacting flows with multidimensional, viscous and reaction effects[END_REF] Granet et al. [START_REF] Granet | Comparison of Nonreflecting Outlet Boundary Conditions for Compressible Solvers on Unstructured Grids[END_REF] show that using the mean Mach number on the outlet condition gives better results than a local value. This latter option is chosen in the following.

The pressure P required to evaluate L 1 in Eq. ( 9) can be either (i) the local pressure or (ii) the average pressure on the boundary condition. In case (i), the pressure difference (P -P t ) is calculated pointwise, i.e. the relaxation locally applies a different correction at each point on the outlet. This formalism is suitable to impose a pressure profile. In case (ii), the wave amplitude is based on the difference between the spatially averaged pressure on the patch and the target pressure P t . The relaxation action is then uniform all over the boundary, whatever is the local offset to the target. This allows local differences of pressure, and is thus compatible with the pressure distortion created by a swirling flow. The latter method is illustrated in Fig. 1 for the pressure points located along the x 1 axis of a boundary. A specific property of NSCBC is that it solves the momentum equation in the outlet section: only the axial acoustic wave L 1 is modified. All other equations including the momentum equation on the radial direction, are solved on the outlet patch and therefore can naturally capture the radial equilibrium.

Two test cases presented in the next part aim at assessing this property: a swirling flow in an annulus and a transonic nozzle guide vane. The AVBP solver is used to compute the test cases. It is a 3D code solving the unsteady compressible Navier Stokes equations using a cell-vertex finite volume approximation. Outlet NSCBC accounting for the transverse terms 1 have been implemented in AVBP by Granet et al. [START_REF] Granet | Comparison of Nonreflecting Outlet Boundary Conditions for Compressible Solvers on Unstructured Grids[END_REF] and are used here. This condition will be referred to as NSCBC 3D. 1

IV. Test cases

IV.A. Annulus

First, the ability of characteristic boundary conditions to recover the REA pressure profile is assessed in a rotating flow in a simple annulus (Fig. 2). The annulus has an aspect ratio such as L/(R 2 -R 1 ) = 4, and is meshed with 4 631 800 tetrahedral cells. This geometry is chosen because is allows to find an analytical solution (Section IV.A.1). This solution is compared to simulation in Section IV.A.2.

IV.A.1. Analytical solution

The geometry and flow properties are set so that all hypothesis for simplified radial equilibrium are fulfilled: axisymmetry, steady state, no viscosity (Euler), ideal gas. Two configurations are tested: the free vortex flow (c θ = k/r with k a constant) and a constant swirl flow (c θ = C). For a ideal gas, Eq. ( 2) with the assumed expression for the swirl component of the flow c θ (r) gives the radial pressure profile expression:

x Free vortex flow c θ = k r P (r) = α 1 exp ( - k 2 2R air T r 2 ) ( 11
)
y Solid body rotation

c θ = C P (r) = α 2 r   C 2 R air T   (12)
where α i are integration constants which can be determined using a pivot point where the local pressure is known, or via the mean pressure on the outlet plane. The latter strategy is used here for the two cases, targeting the same mean pressure P t :

P t = P = 1 π(r 2 ext -r 2 int ) ∫ ∫ P (r)rdrdθ (13)
Solving Eq. ( 13) for the constants α i gives for cases x and y:

α 1 = P t (r 2 ext -r 2 int ) [ k 2 2R air T Ei ( - k 2 2R air T r 2 ) + r 2 exp ( - k 2 2R air T r 2 )] rext rint (14) α 2 = P t (r 2 ext -r 2 int ) 2   r ( C 2 R air T +2 ) C 2 RairT + 2   rext rint (15) with the common notation [Φ(r)] b a = Φ(b) -Φ(a). Ei(x) is the exponential integral function defined as Ei(x) = ∫ x -∞ e t /tdt.
The target pressure P t = P = 100kP a appears only in the expressions of α 1 and α 2 (Eqs. 14 and 15), and is identical for both cases. It is therefore only a multiplicative constant for the pressure profile as indicated in Eqs. (11) and (12). On the other hand the shape of the pressure profile, governed by radius dependant terms, is controlled by gas properties (R air = 287J.kg -1 .K -1 , T = 300K) and by the swirl profile c θ (r). Therefore, the two pressure profiles have the same mean pressure, but the free vortex flow (case x) has a stronger pressure gradient than case y.

IV.A.2. Numerical resolution

The two swirling flows x and y are simulated for six values of the reflection coefficient K, ranging from 1 to 1000, in order to assess the influence of the level of reflectivity on the radial equilibrium pressure profile. The summary of the test cases can be found in Table 1. The compressible solver AVBP is used here, with the Lax-Wendroff 24 numerical scheme, providing second-order accuracy in time and space. In order to be fully-compliant with the REA, only Euler equations are solved here, without sub-grid scale viscosity model. Non-reported results show that for the annulus flow, the viscosity effect is very low and thus solving either Euler or Navier-Stokes equations gives the same results.

The inlet boundary condition imposes the adequate swirl profile c θ (r) according to case x or y. The axial velocity is set to c x = 5m/s for x and c x = 30m/s for y. The corresponding flow through times τ = L/c x is 64ms for case x and 11ms for case y. The maximum Mach number for case x and y are respectively 0.29 and 0.17. The time step is controlled by the CFL condition (set to 0.7), and is around 1.5 • 10 -6 s. In each case, the initial velocity field in the annulus is identical to the inlet boundary condition velocity profiles (c θ ; c x ). At the outlet, NSCBC 3D 1 is used, driving the mean pressure towards the imposed value P t . The inner and outer walls use slip conditions for velocity in order to avoid any possible near-wall effect on the velocity profile and thus on the pressure distribution. The initial pressure field is uniform P (r, θ, x) = P t , and thus not consistent with the REA, in order to test the ability of NSCBC to drive the pressure towards REA . The required time to establish the pressure field is imposed by the domain size, flow properties, and relaxation parameter K. Table 1 summarizes the simulation times required to converge the mean pressure at exit. Convergence is reached when the local pressure is ±0.1% of P (t → ∞). Increasing the reflectivity of the NSCBC condition (via K) allows to reduce the convergence time down to 0.8τ for case x and 2.9τ for case y, reached for K ≥ 50. Changing from K = 1 to K = 1000 allows to reduce the convergence time by respectively 50 and 140 times for cases x and y.

As the pressure profile depends only on the radius (purely 1-D flow), it can be plotted at the domain exit, as shown in Fig. 3: a radial equilibrium pressure profile establishes at the end of the simulation. The simulations converge towards the theoretical profile as the relaxation parameter at the outlet is increased for both cases x and y. The main action of the relaxation coefficient K is to drive the mean pressure towards the target. For example in the bottom part of Fig. 3 (case y) the pressure obtained for a low relaxation simulation (K = 1) has a mean value which is much smaller than the target value P t . The reduction of the offset between mean pressure and target when increasing the BC reflectivity level is highlighted in Fig. 4.

In this log-log plot, one can observe that the error is inversely proportional to the parameter σ over the investigated range. The error, as well as the convergence time (see Tab. 1) are slightly lower for the free vortex flow x than for the solid body rotation y. This can be explained by the fact that the transverse terms are more important in case x: 20 ≥ c θ /c x ≥ 14 than in case y: c θ /c x = 1.7. For σ → ∞ (fully reflecting BC) the mean pressure would converge to the imposed value and thus the pressure profile would exactly match the theoretical one. However, this case would also lead to a fully reflecting outlet, a property which may have to be avoided, leading to an intermediate value for K providing both a limited error on pressure and a non reflecting outlet.

+ NSCBC 3D 1 -K = 5 × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × NSCBC 3D 1 -K =
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • 0.
(r -R 1 )/(R 2 -R 1 )
Case y: constant swirl flow

Theory (REA)
Mean value NSCBC 3D 1 -K = 1 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + ++ This first test case shows that a characteristic boundary condition can cope on its own with the radial equilibrium, i.e. with a non uniform outlet pressure dictated by the swirl. Note that this is obtained without having to pre-compute a pressure profile and impose it at the outlet.

+ NSCBC 3D 1 -K = 5 × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×× × NSCBC 3D 1 -K =
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • •

IV.B. Turbine vane

The ability of NSCBC to let the radial equilibrium establish at a swirled flow outlet while controlling wave reflection is now assessed on a real 3D configuration: the uncooled transonic Nozzle Guide Vane (NGV) MT1. This vane was developed and instrumented at Qinetic on the Isentropic Light Piston Facility [START_REF] Povey | The effect of hot-streaks on HP vane surface and endwall heat transfer: An experimental and numerical study[END_REF][START_REF] Povey | A hot-streak (combustor) simulator suited to aerodynamic performance measurements[END_REF][START_REF] Povey | Developments in hot-streak simulators for turbine testing[END_REF][START_REF] Simone | Analysis on the effect of a nonuniform inlet profile on heat transfer and fluid flow in turbine stages[END_REF][START_REF] Qureshi | Effect of simulated combustor temperature nonuniformity on HP vane and end wall heat transfer: An experimental and computational investigation[END_REF] for the European project TATEF 2. A periodic sector (1/32) is considered here. Note that during the test campaign, a full turbine stage was experimented whereas our study focuses only on the NGV.

The MT1 vane was selected as it does not feature technological effects, and provides an almost uniform flow deviation. Therefore, the assumption of a constant swirl flow c θ (r) = C is valid and the analytical determination of the simplified radial equilibrium which the code should predict is straightforward.

IV.B.1. Geometry

For the numerical tests, an annulus is attached immediately after the vane trailing edge (Fig. 5). Simulations are performed on three domains of increasing length: for the first one the annulus is truncated at one axial chord length (acl), the second at 2acl and the longest one at 10acl. With short domains, the flow in the outlet section is more swirled and perturbed by the vanes wakes, thereby making the task of the outlet condition more difficult. The 10acl domain is the baseline configuration: it is assumed that the outlet BC has no influence on the near-vane flow field for this case. Note that the 1 and 2 axial chord configurations use the same mesh as the baseline one, simply split at the right location. A full tetrahedral mesh containing 2.69 million nodes (13 million cells) for the 10acl domain is created using Gambit, as shown in Fig. 5. The grid is refined in the wakes region and around the airfoil (tetrahedral edge length ≃ 250µm at vane wall). The Y + values around the vane walls are between 1 and 550, and the flow in the vicinity of the wall is modelled with adiabatic law-of-the-walls [START_REF] Schmitt | Simulation aux grandes échelles de la combustion étagée dans les turbines à gaz et son interaction stabilitépolluants-thermique[END_REF][START_REF]ANSYS Fluent User's Guide[END_REF] . Inner and outer casing walls also use adiabatic law-of-the-wall. Axiperiodicity is enforced on the two lateral sides. Uniform total temperature and pressure profiles are imposed at the inlet, without swirl, as shown in Table 2. The Mach number is around 0.12 upstream at vane leading edge.

IV.B.2. Numerical resolution

In order to compare boundary conditions, two Navier Stokes equations solvers are benchmarked on this test case. Along with the AVBP solver used in section IV.A, the commercial code Fluent v13 is tested. Following the conclusions of the test case IV.A, AVBP is used with NSCBC outlet boundary condition 3D. [START_REF] Yoo | Characteristic boundary conditions for direct simulations of turbulent counterflow flames[END_REF] Fluent uses two distinct outlet BC: either a typical uniform pressure (Neumann), or a specific radial equilibrium option. For the uniform pressure BC, a Fluent in-house weak enforcement of pressure is used: the face pressure is not directly enforced to the target value, but results from a weighted balance. [START_REF]ANSYS Fluent User's Guide[END_REF] The Fluent radial equilibrium option is based on the Neumann formalism but the simplified radial equilibrium pressure profile is calculated at each radius using the domain cell data (swirl velocity and density) and then imposed at the exit. The user specifies the pressure at hub (position of smallest radius) to solve the pressure equation Eq. 2. Note that knowing the pressure at hub is not straightforward since the pressure can differ significantly from the mean pressure on the outlet section.

Fluent solves the steady-state compressible Navier-Stokes equations (density-based approach) closed with a RANS model (k -ϵ realizable [START_REF] Shih | A new k-eddy viscosity model for high reynolds number turbulent flows[END_REF] ). An implicit formalism is used allowing for CFL numbers up to 10. Turbulence variables are defined at inlet and outlet with 5% turbulence intensity and a viscosity ratio µ t /µ = 100. Convective fluxes are computed with a Roe flux-difference splitting scheme [START_REF] Roe | Characteristic-Based Schemes for the Euler Equations[END_REF] and all equations are discretized with a second order upwind scheme.

AVBP uses a Smagorinsky subgrid scale model 34 for subgrid turbulence. The numerical scheme (two step Taylor Galerkin [START_REF] Rudgyard | Development of High-Order Taylor-Galerkin Schemes for LES[END_REF] ) is third-order accurate in time and space. The global time step is controlled by the CFL number (0.7), and is around 7.5 • 10 -8 s for all the LES. The flow through time is estimated as τ = length/c x (Table 3) where length is the axial extent of each domain and c x is the mean axial velocity in the annulus. The LES are time-averaged over roughly 7 flow through times of the 1 acl domain.

The LES of the 1acl domain is difficult because of the interaction of the wakes and the shock with the boundary condition. Strong gradients exist on the boundary. To damp these gradients and related stability issues, viscosity is artificially added 36 in the near exit region (x > 0.55acl where x = 0 is the vane trailing edge).

IV.B.3. Results: flow field

Figure 6 shows a comparison of the time-averaged flow field of the baseline LES (a) with the RANS simulation using REA option (b). The vane midspan contours of density gradient |∇ρ|/ρ (similar to Schlieren image) highlight the shock on the vane suction side, and the boundary layer created on the the suction side. Both solvers seem to predict the same position for the shock, and qualitatively a similar expansion through the turbine. Overall the density gradients are more pronounced in the entire domain for the LES results. Figure 6 (c) shows the same density gradient contours for an instantaneous solution of the 10acl LES, highlighting the dynamic mixing in the wake of the vane. The view is slightly extended compared to Fig. 6 (a) and (b) in order to show how the wake mixes out with axial distance. Overall the density gradient becomes smaller because of the mixing of the wakes, and the less refined mesh which tends to increase the level of turbulent viscosity in this region. The isentropic Mach number is used to compare time-averaged pressure profiles along the vane. It is defined as:

M is = 2 γ -1 ( ( P t,inlet P ) γ-1 γ -1 ) (16)
Profiles of M is are plotted in Fig. 7 at vane midspan, for both RANS simulations and time-averaged LES. As the MT1 vane was equipped with static pressure tapings at 10%, 50% and 90% span on both pressure and suction sides, experimental results 25 are also shown for comparison. However, note that they were obtained with the full turbine stage, i.e. including rotor-stator interactions and rotor potential effect. Therefore, in the rear part of the vane, these simulations of the isolated NGV show a noticeable difference with the experiments. Overall the agreement is good. This is also confirmed at 10% and 90% span (results not shown here). The position of the BC (1, 2 or 10 acl) does not influence the flowfield in the vane, as the three LES predict an identical expansion into the turbine. The static pressure field is now investigated for different axial locations after the vane trailing edge. The post-processing planes are shown in Fig. 8. Fig. 9 shows the static pressure distribution at three locations: 0.2acl, 0.8acl and 1acl for the time-averaged solutions. For all simulations, mixing increases with axial distance from the trailing edge, and therefore the action of the vanes slowly disappears. At x ≥ 0.8acl, this distribution is mostly 1-D on the radial direction and slightly distorted by the wakes on the pitch. The positive radial pressure gradient created by the swirl motion is visible. At x = 0.2acl, the pressure field is very similar from one simulation to another, confirming that the position of the BC does not interfere with the flow.

• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • + + + +
+ + + + + × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ □ 0 0.2 0.4 0.6 0.
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • •• • • • + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × □ □ □ □ □ □ □ □ □ □ □ □

IV.B.4. Comparison with the radial equilibrium analytical profile

Before comparing the simulations results with the REA, one must check that the hypothesis required for REA to be valid are fulfilled in this case (see part II):

-No viscosity: by definition this assumption cannot be respected. However, one can suppose that for this 1-D approach, viscous effects are low, and that Navier Stokes equations are properly approached by Euler equations.

-The flow has to be axisymmetric (∂/∂θ = 0), which is essentially true at x = 1acl, as seen in Fig. 9.

-The radial displacement should be negligible c r << ||⃗ c||, which is true in the plane 1acl as shown by the contours of c r /||⃗ c|| in Fig. 10. As the radial equilibrium is a 1-D concept, the radial velocity field has to be averaged on the azimuthal direction, as shown in Fig. 10 for four axial positions after the vane trailing edge. Even if at x = 0.2acl the streamline curvature imposes a limited radial motion, at x = 1acl the radial component is definitely negligible: below 0.2% of ||⃗ c||.

-The flow is supposed to be at steady state, and therefore LES results have to be time-averaged. The level of pressure fluctuations in the plane x = 1acl is illustrated in Fig. 11 All conditions for the application of simplified radial equilibrium being essentially fulfilled in plane x = 1acl, the radial profiles of pressure can be compared with REA. The theoretical radial equilibrium profile is easily calculated if the flow features a uniform swirl velocity, which is confirmed by the contours c θ /c θ in Fig. 12. The 2-D field is circumferentially averaged, as shown on the right side of Fig. 12 and does not vary much from x = 0.2acl to x = 1acl (maximum deviation to the mean swirl is 5%). Pressure profiles are plotted in Fig. 13 between x = 0.2acl and x = 1acl downstream the NGV trailing edge. Simulations and REA are compared at x = 1acl only (bottom right in Fig. 13). The difference between analytical radial equilibrium profile and LES results is very low with the NSCBC formalism. The same difference can be found when using a Neumann condition associated with a profile of radial equilibrium (Fluent). However, this procedure is less trivial as the outlet BC pressure has to be specified at the hub, which is usually not known a priori. Finally, in this case using a classical uniform Neumann BC leads to more pronounced differences on the pressure field, and local errors between -4 and +6%. In this case, the distortion of the pressure created by the proximity of the BC is visible: at x = 0.2acl the pressure profile tends towards a radial equilibrium one (with an offset on the mean pressure), while at x = 0.8acl and x = 1acl it is clearly distorted by the uniform value imposed at x = 1acl.

• • • • • • • • • • • • • • • • • • • • • • • • • • 0.8 acl + + + + + + + + + + + + + + + + + + + + + + + + + + 1 acl × × × × × × × × × × × × × × × × × × × × × × × × × ×
Probe A • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • Probe B • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • Midspan mean < P > t,θ
• • • • • • • • • • • • • • • • • • • • • 0.8 acl + + + + + + + + + + + + + + + + + + + + + 1 acl × × × × × × × × × × × × × × × × × × × × ×
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • + + + + + + + + + + + + + + + + + + + + + + + + + × × × × × × × × × × × × × × × × × × × × × × × × × 0.
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • + + + + + + + + + + + + + + + + + + + + + + + + + × × × × × × × × × × × × × × × × × × × × × × × × ×

V. Conclusions

Characteristic boundary conditions methods developed for LES and DNS of compressible flows in the last ten years have been tested on flows with strong swirl to see whether these methods allow a control of acoustic waves on boundaries can also cope with the so called radial equilibrium condition. Results show that characteristic boundary conditions including transverse terms [START_REF] Yoo | Characteristic boundary conditions for direct simulations of turbulent counterflow flames[END_REF][START_REF] Granet | Comparison of Nonreflecting Outlet Boundary Conditions for Compressible Solvers on Unstructured Grids[END_REF] naturally allow to establish a radial equilibrium pressure profile for swirling flows. The incoming wave amplitude has to be calculated based on the difference between the target and the mean pressures in the outlet section, to allow local differences due to pressure distortion.

On the annulus test case, results show that low values of reflectivity can be used (K ≃ 50 or σ ≃ 10 -2 ) while avoiding errors on pressure. Simulations of a transonic vane show that the NSCBC formalism lets the radial equilibrium establish despite temporal and local fluctuations due to the vanes wakes. The circumferentially-averaged profile matches the analytical one with ±0.5% accuracy. It was also shown that the position of the boundary condition (10acl, 2acl, 1acl) does not influence the flowfield in the vane, which is not the case when using a uniform Neumann BC at the vane outlet.

These results confirm that characteristic methods should be used in turbine computations where they will have two major advantages:

-Unlike the usual radial equilibrium method, no assumption is required to precompute a pressure profile to impose in the outlet section. The profile is computed naturally by the code which solves the momentum equation on the outlet patch. This suggests that characteristic methods will provide a correct outlet pressure profile even in cases where the assumption of REA does not apply and the REA approximation can not be used.

-A second advantage, which has not been discussed here, is that characteristic methods allow to control wave reflection at the outlet. This property has only limited effects for RANS simulations but is crucial for LES of compressible swirling flows where non reflecting outlet sections are mandatory.

-

  Simplified radial equilibrium can be obtained by considering a flow with the following properties: No radial velocity: c r = 0

Figure 1 .

 1 Figure 1. Illustration of the patch-averaged formalism for the calculation of the incoming wave. Solid line: patch averaged pressure P , dashed line: target pressure Pt.

Figure 2 .

 2 Figure 2. Geometry and mesh for the annulus test case.

  * * * * * * * * * * * * * * * * * * * * * * * * * * NSCBC 3D 1 -K = 50

  20 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * NSCBC 3D 1 -K = 50

Figure 3 .

 3 Figure 3. Radial profile of static pressure P for case x (top) and y (bottom): convergence of the simulations towards REA as the relax is increased.

Figure 4 .

 4 Figure 4. Reduction of the error to the target pressure with increasing relaxation parameter σ, for cases x and y.

Figure 5 .

 5 Figure 5. Geometry and mesh for the MT1 vane: midspan cutview of the 10acl, 2acl and 1acl simulations.

Figure 6 .

 6 Figure 6. Density gradient |∇ρ|/ρ in the midspan plane for the time-averaged (a) and instantaneous (c) 10acl LES; and RANS REA option (b).

Figure 7 .

 7 Figure 7. Isentropic Mach number on the airfoil surface at vane midspan. AVBP (NSCBC 3D 1 ) domain 10 acl ( * ), domain 2 acl (•), domain 1 acl (•); Fluent domain 1 acl Neumann uniform (+) and radial equilibrium option (×) and Experiments full stage 25 (□).

  For lines (c),(d), (e), the pressure fields at x = 1acl are located exactly on the boundary condition. Therefore, the Fluent Neumann uniform BC shows a uniform pressure, while the distribution is clearly radial when the Fluent radial equilibrium option is enabled. The LES of the 1acl domain (c) features a smoother pressure field than longer domains (a) and (b), because of the use of artificial viscosity in the final part of the domain. The flowfield at x = 0.2acl does not seem to be affected by this non-physical damping.

Figure 8 .

 8 Figure 8. Position of the post-processing planes downstream of the vane.

Figure 9 .

 9 Figure 9. Time-averaged non-dimensional static pressure fields P/P 1acl for the different simulations (view from downstream). The interval between two pressure isolines (black lines) is 3% of P 1acl .

  by the pressure signals of two probes located at midspan. The horizontal solid line shows the time and circumference averaged pressure at this radial position. All signals are normalized by the analytical solution given by REA at this radius. The circumferentially averaged pressure is very close to 1, which indicates that the pressure is identical to the one predicted by REA, as it will be shown later. The time-averaged value of the pressure at probe A is lower than the REA value, and a bit higher for probe B. Moreover, ±2% fluctuations to the REA solution are observed. This indicates that (i) the local value of pressure is not uniform versus radius because of the wakes effect etc. (this would not be possible when imposing a 1-D REA profile) (ii) even if relatively small, local fluctuations of pressure are not incompatible with the establishment of radial equilibrium.

Figure 10 .

 10 Figure 10. Time-averaged field of radial velocity cr ∥⃗ c∥ in plane x = 1acl (left) and corresponding circumferentiallyaveraged profile (right). Baseline LES (10acl).

Figure 11 .

 11 Figure 11. Time evolution of pressure for two probes located at midspan in the plane x = 1acl, LES of the 10acl domain

Figure 12 .

 12 Figure 12. Time-averaged field of swirl velocity c θ c θ in plane x = 1acl (left) and corresponding circumferentially-averaged profile (right). Baseline LES (10acl).

Figure 13 .

 13 Figure 13. Radial profiles of circumferentially averaged static pressure P : AVBP (NSCBC) domain 10 acl ( * ), domain 2 acl (•), domain 1 acl (•) and Fluent domain 1 acl Neumann uniform (+) and radial equilibrium option (×). Radial equilibrium profile is shown at 1acl (solid line -)

Table 1 . Test cases for the annulus test case and required time to converge the pressure field

 1

	Relaxation coefficient K		1	5	20	50	500	1000
	σ =	K (1 -M 2 ) a0 R2-R1	Case x 0.0003 0.0013 0.0050 0.0126 0.1257 0.2513 Case y 0.0002 0.0012 0.0047 0.0119 0.1186 0.2371
	Physical time to	Case x	35.6	4.7	1.0	0.8	0.8	0.8
	convergence t/τ	Case y 252.1	36.2	5.6	2.9	2.9	2.9

Table 2 . MT1 vane characteristics

 2

	Aspect ratio h/C	0.65
	Pitch to chord ratio	0.85
	Inlet angle	0°I
	nlet total pressure	460 kPa
	Inlet total temperature 444K
	Mean pressure at 1acl	269 kPa
	Fluid	Pure air

Table 3 . Parameters for the LES and RANS simulations of the MT1 vane

 3

	Domain	No. of nodes τ (ms) Solver	Approach Outlet BC	CPU hours No. of cores Wall-clock
									time (days)
	10 acl	2.69 M	3.9	AVBP	LES	NSCBC 3D 1	35k	144	13
	2 acl	2.62 M	1.3	AVBP	LES	NSCBC 3D 1	19k	144	6
	1 acl	2.59 M	0.95	AVBP	LES	NSCBC 3D 1	13k	144	4.5
	1 acl	2.59 M	0.95	Fluent RANS	Neumann uniform	0.21k	8	2.2
	1 acl	2.59 M	0.95	Fluent RANS	Neumann REA option	0.26k	4	1.3

Acknowledgments

The authors would like to thank the members of the CFD Team at CERFACS and the methods department at Turbomeca. The MT1 vane data used in this work was obtained form the project TATEF II, funded under the Brite Euram Forth Framework on Turbine Aerothermal External Flows, Contract Number BRPR-CT97-0519.