E Denimal
email: e.denimal@imperial.ac.uk

C Wong

L Salles
email: l.salles@imperial.ac.uk

L Pesaresi
email: l.pesaresi@imperial.ac.uk

ON THE EFFICIENCY OF A CONICAL UNDER-PLATFORM DAMPER FOR TURBINES

Underplatform Dampers are commonly used in aircraft engines to limit the risk of High-Cycle Fatigue of turbine blades. The latter is located in a groove between two consecutive blades. The dry friction contact interface between the damper and the blades dissipates energy and so reduces the vibration amplitudes. Two common geometries of dampers are used nowadays, namely wedge and cylindrical dampers, but their efficiency is limited when the blades have an in-phase motion (or a motion close to it), since the damper tends to have a pure rolling motion. The objective of the present study is to analyse a new damper geometry, based on a conical shape, which prevents from this pure rolling motion of the damper and ensures a high kinematic slip. The objective of this study is to demonstrate the damping efficiency of this geometry. Hence, in a first part, the kinematic slip is approximated with analytical considerations. Then, a nonlinear dynamic analysis is performed, and the damping efficiency of this new geometry is compared to the wedge and the cylindrical geometries. The results demonstrate that the conical damper has a high damping capacity and is more efficient and more robust than the two others.

INTRODUCTION

Industrial requirements in the aerospace industry for gas turbine engines are becoming more stringent in terms of mass and efficiency, so that many components of such structures have reached their structural limits [START_REF] Krack | Vibration prediction of bladed disks coupled by friction joints[END_REF]. This is particularly true for turbines blades which are highly loaded, are subjected to high thermal, centrifugal stresses and vibrational stress [START_REF] Petrov | State-of-the-art dynamic analysis for non-linear gas turbine structures[END_REF]. The latter can lead to High Cycle Fatigue (HCF) and so to the failure of the blades [START_REF] Cowles | High cycle fatigue in aircraft gas turbines-an industry perspective[END_REF]. Since the aero-engines have a large operating speed range and the modal density of bladed disk is high, it is impossible to avoid all critical resonances in the system. Hence, reducing the vibrations levels at these resonances is crucial, and solutions based on friction damping have been the most widely used over the years [START_REF] Gaul | The role of friction in mechanical joints[END_REF][START_REF] Feeny | A historical review on dry friction and stick-slip phenomena[END_REF][START_REF] Griffin | A review of friction damping of turbine blade vibration[END_REF].

Dry friction can be introduced in different locations to limit the vibration levels, such as at the shrouds, the roots and blade tips. However, the most effective solution relies on the use of Underplatform Dampers (UPD) [START_REF] Griffin | A review of friction damping of turbine blade vibration[END_REF][START_REF] Szwedowicz | Numerical and experimental damping assessment of a thin-walled friction damper in the rotating setup with high pressure turbine blades[END_REF]. They consist in a metallic device placed in the groove under the platforms of two adjacent blades. The centrifugal loading keeps them in place during operation. If the blades vibrate, then the relative displacement between the platforms and the damper creates friction at the contact interface and so, energy is dissipated, and the system is damped [START_REF] Sanliturk | Underplatform dampers for turbine blades: theoretical modeling, analysis, and comparison with experimental data[END_REF][START_REF] Sanliturk | Friction damper optimization: simulation of rainbow tests[END_REF].

The identification of the most efficient UPD geometry is still an active research topic. Nevertheless, two main geometries have largely been studied and employed until now. The first shape is based on a wedge geometry [START_REF] Panning | Asymmetrical underplatform dampers in gas turbine bladings: theory and application[END_REF][START_REF] Petrov | Advanced modeling of underplatform friction dampers for analysis of bladed disk vibration[END_REF][START_REF] Panning | Optimization of interblade friction damper design[END_REF][START_REF] Jareland | A parametric study of a cottageroof damper and comparison with experimental results[END_REF][START_REF] Firrone | Effect of static/dynamic coupling on the forced response of turbine bladed disks with underplatform dampers[END_REF][START_REF] Pesaresi | Modelling the nonlinear behaviour of an underplatform damper test rig for turbine applications[END_REF] and the second shape is based on a cylindrical geometry [START_REF] Panning | Optimization of interblade friction damper design[END_REF][START_REF] Jareland | Experimental investigation of a platform damper with curved contact areas[END_REF][START_REF] Csaba | Modelling of a microslip friction damper subjected to translation and rotation[END_REF]. Some studies present a combination of the two shapes [START_REF] Bessone | Investigation on the dynamic response of blades with asymmetric under platform dampers[END_REF][START_REF] Gola | Understanding complexities in underplatform damper mechanics[END_REF][START_REF] Zucca | Range of variability in the dynamics of semi-cylindrical friction dampers for turbine blades[END_REF] and how to optimise it [START_REF] Gastaldi | Pre-optimization of asymmetrical underplatform dampers[END_REF].

The wedge geometry is the shape that is the most widely used over the years. It is easy to manufacture, install and replace, and it seals the platforms. However, with this shape, the prediction of the contact locations and conditions is difficult, which makes the modelling complex and the behaviour of the structure difficult to predict [START_REF] Pesaresi | Modelling the nonlinear behaviour of an underplatform damper test rig for turbine applications[END_REF][START_REF] Pesaresi | An advanced underplatform damper modelling approach based on a microslip contact model[END_REF]. Moreover, each damper is loaded in a different way, which complicates even more the global dynamic of the system. The interesting feature of a cylindrical geometry is that the contact pressure is controllable and so the global behaviour of the structure is more predictable. However, both shapes have the same limitation. Indeed, when the blades have an in-phase motion, then the damper tends to have a pure rolling motion and no more frictional dissipation occurs. The dissipated energy is null or low, and the damping efficiency is limited.

The objective of the present study is to present a new damper geometry based on a conical shape. This geometry makes it impossible for the damper to have a pure rolling motion when the blades have an in-phase motion. In fact, to maintain contact with the two platforms, sliding is required so friction occurs. This friction is directly related to the kinematic motion of the damper created by the in-phase motion of the platform. Hence, the energy dissipated by friction is expected to be high as well as the damping efficiency.

The paper is organized as follows. In a first part, an approximation of the expected kinematic slip of this new damper is determined. It is based on an analytical analysis with a few assumptions. It gives quickly an approximation of the efficiency of the damper with regards to its geometrical characteristics. Then, for more accurate results, a dynamic analysis based on the Multi Harmonic Balance Method (HBM) is done, where the wedge, cylindrical and conical dampers are compared. The nonlinear receptances as well as the energy dissipated at the contact interface are considered to compare the efficiencies of the different dampers.

NEW GEOMETRY PRESENTATION AND ANALYTIC PERFORMANCES

Bladed disks are composed of a disc on which blades are fixed, and UPD are located under the platforms of adjacent blades. In this study, the focus is put on a system composed of two blades and one damper. If one considers the first bending mode of the blades, two extreme cases exist: an in-phase motion of the two blades and an out-of-phase motion. It is worth noticing that in reality, the blades have a mix of these two motions, but in the following, only these two "fundamental" conditions will be considered to study the behaviour of the damper. Considering these two cases, the kinematic of the platforms can be simplified, respectively, by a pure vertical or a pure horizontal motion of the platforms as illustrated Figure 1. For a horizontal

Neutral position

Platform horizontal motion (blades OOP) Platform vertical motion (blades IP) FIGURE 1: Illustration of movement of the platforms and the rolling motion motion of the platforms, sliding occurs between the damper and the platforms, so energy is dissipated. However, in the case of a vertical motion, the damper tends to have a pure rolling motion and so the damper efficiency is very limited.

To cope with this issue, a different geometry is presented here. It is based on a conical shape as represented Figure 2. The damper is in blue and the planes defined by the platforms are in grey. The main feature of this geometry is that the rolling motion is not possible when a vertical motion of the platforms occurs. Indeed, to keep the contact with the two platforms and to keep its axis parallel to the platforms, the damper cannot have a pure rolling motion and so some sliding must occur. The latter creates friction and so energy is dissipated. This dissipation relies on the pure kinematic motion of the platforms and so a high efficiency of the damper is expected. The study is aimed to demonstrate and illustrate this basic idea. Moreover, this shape presents the advantage to have a robust description of the contact. Indeed the latter is reduced to a line which is easier to model, locate and gives a more controllable behaviour. This section is more specifically dedicated to analytical considerations.

Axial Inclination of the Damper

When the damper is positioned between the two platforms, its axial inclination is fully determined by the cone angle and the angle between the platforms. Moreover, the contact between the damper and each platform corresponds to a line which is also a generator of the cone. An illustration is given Figure 2. To determine the cone axial inclination, a few parameters need to be defined:

γ represents the axial inclination of the damper relative to the vertical. It is the angle between the plane formed by a section of the damper and the vertical, β is the characteristic angle between the two platforms, α 1 corresponds to the projection of the platform angle on the plane formed by the large horizontal section of the damper, The different angles are represented Figure 2.

From geometrical considerations, two relations can be obtained to relate the damper inclination γ to the platform angle β and the cone angle α c :

tan(α c ) = tan(γ) • sin(α 1 ) tan(α 1 ) = tan(β ) • cos(γ) (1)
From some manipulations it comes:

sin(γ) tan(β ) = tan(α c ) (1 + (cos(γ) tan(β )) 2 ) (2)
which gives the following expression of the damper inclination γ in function of the geometrical characteristic of the problem:

γ(α c , β ) = arccos (1 -tan 2 (α c ) cot 2 (β )) tan 2 (α c ) + 1 ( 3
)
From Equation 3, the effect of the platform angle β and of the cone angle α c on the damper inclination can be determined. For example, Figure 3 gives the evolution of the damper inclination depending of the cone angle α c for different platform angles β . As expected, when the damper is cylindrical, the inclination angle is null, and so the damper horizontal. On the opposite, when the cone angle becomes equal to the platform angle, then the damper sits in a vertical position at 90 degrees with respect to the platforms.

Kinematic Slip

Once the damper inclination is determined, it is possible to get an approximation of the kinematic slip that will occur between the damper and the platforms. The latter is determined when a vertical motion of the platforms occurs, i.e. for the inphase family of modes of the blades. From this analytical approximation, it is possible to determine the influence of the different parameters on the expected damping efficiency. In order to estimate the kinematic slip, a few assumptions are necessary:

-the blades in-phase family of modes is simplified as a pure translational vertical motion between the platforms and the angle between the platforms cannot change, -the damper and the platforms are always in contact, -the damper and the platforms are perfectly rigid, and no elastic deformation is allowed, -the middle of the damper is supposed to have a pure rolling motion.

The second and the last points are strong assumptions. The fullcontact assumption between the damper and the platforms comes from the idea that from the loading and the conical shape, the damper will always stay stuck to the platforms (no lift-off). The pure rolling motion of the middle section of the damper comes from the idea that the damper tends to have a symmetric motion to have similar displacements at each extremity. It is worth pointing out that these assumptions simplify significantly the problem and are only here to get a quick estimation of the kinematic slip Assuming a vertical motion of the platforms v, the vertical component in the plane of the middle section of the damper is 4a. So, if the middle section has a pure rolling motion, then a rotation θ appears with:

v t = v • cos(γ), illustrated Figure
θ ≈ v t R m cos(α 1 ) (4)
where R m is the radius of the middle section. The smaller radius is denoted R 1 and the large radius R 2 . The angle θ is assuming to be small, which is justified by the fact that the vibration amplitudes are extremely small compared to the geometrical dimensions. At the same time, each extreme section of the damper rotates of the same angle θ . However, the distances travelled by the small and the large sections are different than the distance travelled by the middle section. The difference between the distance travelled by the middle section and the small (resp. large) section of the damper is equal to δ 1 (resp. δ 2 ):

δ 1 = v t - R 1 R m v t δ 2 = v t - R 2 R m v t (5)
Once the distance δ 1 (resp. δ 2 ) is travelled by the small (resp. large) section, the latter is not in contact with the other platform. In order to maintain contact with both platforms, it must slip the remaining distance denoted s 1 (resp. s 2 ). The different configurations of an extreme section of the damper are represented Figure 4b. The black configuration corresponds to the initial configuration, the red configuration to the configuration after travelling the distance δ and the blue configuration to the final configuration (i.e. δ and s travelled). Hence, the kinematic slip corresponds to s 1 (resp. s 2 ) in Figure 4b and is equal to:

s 1 = (1 -R 1 R m ) cos(α 1 ) v t s 2 = (1 -R 2 R m ) cos(α 1 ) v t (6)
A longitudinal kinematic slip s l also appears and is equal to s l = v s cos(α c ).

As an illustration, the absolute value of the kinematic slip at each position of the damper for different cone angles α c is given Figure 5. The length of the damper is L = 40 mm and the platform angle is β = 45 o . The kinematic slip is maximum at the extremities of damper and null at the centre as expected. Indeed, the main assumption if that the centre of the damper has a pure rolling motion, and so does not slip. The evolution is linearly dependent to the distance to the centre.

More generally, the maximum expected kinematic slip for different cone angles and for different configuration of the system can be determined. The influence of different parameters is illustrated Figure 6. In Figure 6a, the influence of the platform angle β is studied. The evolution of the maximum of the kinematic slip is similar in each case. More precisely, the maximum normalized kinematic slip increases first when the cone angle increases, then it reaches a maximum before decreasing and becoming null when the cone angle is higher than the platform angle (i.e. when the damper has a vertical position). Hence, an increase of the platform angle tends to increase the maximum of the kinematic slip. Moreover, the maximum is reached for different values of the cone angle. Indeed, if β = 36 o , the maximum is reached when α c is around 55% of the platform angle, whereas if β = 60 o the maximum is reached around 44% of the platform angle. The influence of the radius of the large section R 1 and of the damper length L are given Figures 6b and6c, respectively. The evolution of the maximum of the kinematic slip is similar than in the case presented previously. However, an increase of the radius R 1 tends to decrease the maximum of expected slip: if R 1 = 5 mm then the maximum is around 0.7, whereas if R 1 = 20 mm the maximum is around 0.3. Concerning the length of the damper, an increase of the latter tends to increase the maximum of the kinematic slip: if L = 40 mm then the maximum is around 0.4, whereas if L = 60 mm then the maximum is around 0.6. These different analytical considerations provide a coarse estimation of the expected kinematic slip and so of the damper efficiency. Moreover, they allow in a preliminary process to test and compare quickly the impact of some design parameters as the length of the damper, the cone angle or the radius of the sections.

COMPARISON WITH OTHERS DAMPERS

The analytical results given in the first part can be very useful in a first design process to identify the main characteristics of a damper. However, nonlinear dynamic simulations must be performed to get the real damping performances. The following part is devoted to the analysis of the performances of the conical damper to validate the analytical approach presented previously. The latter will be compared to a cylindrical and to a wedge damper.

Model Presentation

Blades Model The model used is the model of an underplatform test rig, used in previous works [START_REF] Pesaresi | Modelling the nonlinear behaviour of an underplatform damper test rig for turbine applications[END_REF]. It has been designed to experimentally study the effect of UPDs on blade-like structure. The Finite Element Model (FEM) of the model is represented Figure 7a and is composed of two pseudo beam-like blades that are fixed on a common base which simulates a rigid disk. The dimensions of the model are chosen to mimic the behaviour of the dynamic response of a real HPT blade. The interested reader is invited to refer to [START_REF] Pesaresi | Modelling the nonlinear behaviour of an underplatform damper test rig for turbine applications[END_REF] for more information. The mesh is composed of 54972 quadratic hexahedral elements and a stainless-steel material is used, i.e. the Young modulus is equal to 197 GPa and the density is equal to 7800 kg/m 3 . The first two modes of the structure correspond to the first bending mode of the blades, that can be in-phase or out-of-phase. The mode shapes are represented Figure 7b and Figure 7c respectively. Damper Models Three damper geometries are considered for this part of the study, they are represented Figure 7. The first one is the wedge damper studied in [START_REF] Pesaresi | Modelling the nonlinear behaviour of an underplatform damper test rig for turbine applications[END_REF] and represented Figure 7d. The second one is a conical damper of cone angle 10 • represented Figure 7e. Finally, the third one is a cylindrical damper with a radius equal to the mean radius of the conical damper, and is represented Figure 7f. The meshes are constructed so that the contact points are matching the blades mesh.

Contact Pressure In all cases, a uniform contact pressure is assumed. For a wedge damper, the contact pressure due the centrifugal loading C F is modelled by a constant value [START_REF] Petrov | Advanced modeling of underplatform friction dampers for analysis of bladed disk vibration[END_REF] and

σ 0 = 1 2 C F A(cos α + µ sin α) (7)
where C F is the radial centrifugal force, A is the contact area on each side, µ is the friction coefficient and α is the damper angle.

For the cylindrical and the conical damper, it is assumed that the total contact force remains the same. Since the contact between the damper and the platform is a line in these cases, the total contact force is divided by the number of points (the points are regularly positioned).

Contact Formulation A node-to-node contact formulation is used to solve the problem. Thus, the contact surface of the wedge damper is reduced to 146 points as in [START_REF] Pesaresi | Modelling the nonlinear behaviour of an underplatform damper test rig for turbine applications[END_REF], and the contact line for both the cylindrical damper and the conical damper is composed of 23 points (i.e. 46 in total).

The dynamic contact behaviour is described at each contact node by the 3D friction contact shown in Figure 8. This contact model is composed of two decoupled Jenkins elements [START_REF] Gaul | The role of friction in mechanical joints[END_REF] allowing for 2D in plane displacements, and a third spring is added to allow normal load variations. Hence, each element is characterized by four parameters, namely the friction coefficient µ, the initial pre-load N 0 , the tangential stiffness k t and the normal stiffness k n . The normal pre-load and the stiffnesses are normalized by unit of area, so each contact element has its own characteristics. The contact formulation is the same for the different dampers and allows four different states for each point, namely stuck, stick/slip, partial separation and full separation.

x(t) y(t) z(t) k t k n ¡ gap=-N 0 /k n k t FIGURE 8: Contact model

Harmonic Balance Method

The numerical strategy used for the prediction of the vibrations levels of the system is based on the Multi-Harmonic Balance Method (MHBM) coupled with a model reduction technique included in the in-house code FORSE (FOrced Response SuitE), which is presented in details in [START_REF] Petrov | Analytical formulation of friction interface elements for analysis of nonlinear multiharmonic vibrations of bladed disks[END_REF][START_REF] Petrov | Explicit finite element models of friction dampers in forced response analysis of bladed disks[END_REF][START_REF] Petrov | A high-accuracy model reduction for analysis of nonlinear vibrations in structures with contact interfaces[END_REF]. The friction forces from the contact interface bring non-linearities in the system. The dynamic equation of the system can be written as:

M Ẍ + C Ẋ + KX + F nl ( Ẋ, X) = P(t) (8)
where M, C and K are the mass, damping and stiffness matrices respectively. X is the vector of displacements and the dot corresponds to the derivative with respect to time. P is the vector of the external excitation and F nl are the non-linear contact forces, that depend of the relative motion of nodes in contact at the interface.

The main idea of the MHBM consists in considering that the response of each DOF is periodic and so can be decomposed on a Fourier series which is truncated at the n th harmonic:

q(t) = Q 0 + n ∑ j=1 Q c j cos m j ωt + Q s j sin m j ωt (9)
where Q 0 , Q c j and Q s j are the Fourier coefficients. By injecting Eq. 9 into Eq. 8, and projecting on each harmonic, the system 8 can be approximated by:

Z(ω)Q + Fnl (Q) -P = 0 ( 10
)
where

Q = {Q 0 , Q c 1 , Q s 1 , . . . , Q c n , Q s n }
is the vector of the Fourier coefficients, Z(ω) is the dynamic stiffness matrix of the system, Fnl and P are the nonlinear forces and the external efforts in the Fourier basis respectively. The reduction of the model is based on a FRF matrix representation of the model. This FRF matrix is determined by calculating its exact value at some frequency points and adding a second term that describes its variation in a range of frequencies [START_REF] Petrov | A high-accuracy model reduction for analysis of nonlinear vibrations in structures with contact interfaces[END_REF]. An Alternate Frequency/Time (AFT) [START_REF] Cameron | The orthogonal development of non-linear functionals in series of fourierhermite functionals[END_REF][START_REF] Salles | Dynamic analysis of a bladed disk with friction and fretting-wear in blade attachments[END_REF] algorithm is used for the determination of the nonlinear forces and the problem of Eq. 10 is solved with a Newton-Raphson algorithm. In order to compute a continuous response curve in the frequency domain, continuation is performed with regard to the pulsation ω. Here, a secant predictor together with an arc-length corrector is used. In the following the harmonics 1, 2 and 3 are used for the simulations. They have been selected based on a convergence study. For the dampers and the blades, 12 modes are used for the simulation.

Comparison of the Performances of the Different Dampers

The nonlinear receptances of the three dampers are computed for different levels of excitation from 0.096 N to 30 N. The receptances are given Figure 9 for the out-of-phase family case, and Figure 10 for the in-phase case.

The nonlinear receptances for the out-of-phase mode are given Figure 9 for the three cases. The curves with a solid line correspond to the conical damper receptances, the dashed ones to the cylindrical damper and the dotted ones to the wedge damper. The first observation is that compared to the cylindrical damper, the behaviour of the conical damper is almost identical with a small frequency shift of about 0.5 Hz. In all cases, the three dampers have a similar behaviour and the damping efficiencies are comparable.

But considering the in-phase case, the behaviours of the three dampers are completely different (see Figure 10). Indeed, as a first observation, the resonance frequencies are spread from 309 Hz to 423 Hz. As it can be seen Figure 10b, all the FRF are superimposed. It demonstrates that the damper has no effect on the dynamic of the blades, in other words the cylindrical damper has a pure rolling motion for this mode and the damping efficiency for this mode is very limited. Hence, the blades have a linear behaviour. Considering the wedge damper, it damps the response since the maximum of the receptances decreases when the excitation force increases but it also presents a strong softening behaviour and an important frequency shift (around 4 Hz of difference between the peak at F=0.096 N and the peak at F=30 N). Finally, the conical damper presents a behaviour similar to the behaviour for the out-of-phase mode, i.e. a high damping efficiency without frequency shift.

The maximum of the absolute amplitude response is also given for the three dampers in Figure 11 for the in-phase case (Figure 11a) and the out-of-phase case (Figure 11b). For the outof-phase case, the three dampers have a similar behaviour. The amplitudes increase when the excitation increases and stabilize between 0.06 mm and 0.07 mm for high excitation levels. For the in-phase case, the three dampers have different behaviours. Indeed, as seen previously, the blades have a linear behaviour if the damper is the cylindrical one, so the absolute amplitude increases linearly with the excitation, which gives very high level of amplitude at high level of excitations. For the wedge damper (yellow curve in Figure 11a), the maximum of amplitude also increases with the level of excitation but more slowly and reached a maximum of 0.35 mm for F = 30 N. But for the conical damper (blue curve), at very low level of excitations, the amplitude response increases but after 5 N, it increases slightly from 7.8E-2 mm at 9.6 N to 8.7E-2 mm at 30 N of excitation amplitude, and almost does not depend of the excitation amplitude. This demonstrates the high damping capacity of the conical damper for the in-phase case. In fact, the damper is even more efficient than the wedge damper since there is a saturation effect. To investigate more deeply the efficiency of the different dampers, the energy dissipated by friction for the different dampers and the different level of excitation are determined. They are given Table 1 for the out-of-phase case and Table 2 for As expected, the level of energy dissipated by friction for the out-of-phase case is similar for the three dampers. Indeed, the quantity of dissipated energy is almost equal for the three dampers at each level of excitation (see Table 1). However, in the case of the in-phase mode family, the energy dissipated differs a lot from one damper to another. First, the cylindrical damper does not dissipate energy (see the low levels around e-18 mJ), which confirms the pure rolling motion of the damper for this case. This analysis is confirmed by the evolution of the contact status of the contact points given Figure 12 for the different dampers for excitation levels equal to 9.6 N and 30 N. For the cylindrical and the conical dampers, the points from 1 to 23 correspond to the left contact line, and the points 24 to 46 to the right. They are ordered so that the points 1 and 24 correspond to the large section and the points 23 and 46 to the small section. For the wedge damper, the points 1 to 73 (resp. 74 to 146) correspond to the left (resp. right) contact surface. Red points correspond to cases where stick/slip occurs at the contact and green points where separation occurs. As it can be seen, in the case of the cylindrical damper (top figures), friction occurs only for the out-of-phase case (i.e. at 409 Hz). For the in-phase case, the contact points remain stuck and the damper has a pure rolling motion. Considering the conical damper, stick/slip occurs in both cases, for either the in-phase and the out-of-phase modes. Moreover, for the in-phase case, stick/slip appears at the extremities of the damper when the centre has a pure rolling motion (see the stuck points in the middle of each contact line). This confirms the assumption done in the analytical part. Moreover, no full contact separation is observed for the different amplitudes of excitation, justifying the full-contact assumption done in the analytical part. Since stick/slip occurs for the in-phase case, some energy is dissipated by friction and so the damper has, as expected, a damping effect on the system. Considering the wedge damper, the dynamic of the contact is much more complex for the in-phase case since almost all points experience contact sep-aration (see Figures 12e and12f). This illustrates the rolling motion of the damper and explain the strong softening behaviour observed Figure 10c. For more insight, the evolution of the repartition of the contact condition at the first resonance frequency for both the wedge and the conical dampers for different excitation amplitudes is given Figure 13. As it can be seen, for the conical damper (solid lines) the number of contact points in a stick/slip or partial contact state increases when the excitation amplitude increases. This evolution may explain why the maximum of the amplitude response remains almost constant in Figure 11a for the conical damper. Considering the wedge damper, when the excitation amplitude increases, all the points tend to experience contact separation. The energy dissipated by friction by the wedge damper in the in-phase case is of the same order of magnitude than the energy dissipated by the conical damper. Hence, the conical and the wedge dampers have similar performances in terms of energy dissipation. However, the results obtained here for the wedge damper are based on the strong assumption that the contact surface is in perfect contact. In reality, the location and the surface area of the contact are unknown, and it results in a drop of the resonance frequency for this mode. Hence, the modelling of this damper is limited, and the resonance frequency cannot be estimated. On the opposite, the contact condition between the conical damper and the platforms is much more robust since it corresponds to a contact line that can be easily identified. It results in a strong ability to model and capture the behaviour of the structure. Moreover, the coupling between the platforms can be controlled, and so the frequency shift caused by the damper. For these different reasons, the conical damper appears more efficient and reliable than the wedge damper.

CONCLUSION

In the present study, a new damper geometry has been evaluated. Based on a conical geometry, the problem of pure rolling motion when the blades have an in-phase motion cannot occurs. Indeed, to maintain contact with the two platforms, slipping motion must takes place, and so friction appears. This friction is due to the pure kinematic motion of the platforms, and so the conical damper has a very high damping potential.

First, from geometrical considerations and some assumptions, an analytical formulation of the expected kinematic slip is determined. It allows a fast estimation of the performances of a damper for different geometric properties in a preliminary design process. Then, the real damping efficiency of the damper is determined by using a nonlinear dynamic analysis based on the HBM. The nonlinear receptances of the conical damper for different excitation levels are compared to the receptances of a cylindrical and a wedge damper. The damping properties of the different dampers are similar in the case of an out-of-phase motion of the blades. But in the case of an in-phase motion of the blades, the behaviours of the three dampers tend to be very different. The cylindrical damper has a pure linear behaviour and has almost no damping properties. The wedge damper has a strong softening behaviour, and so a non-negligible frequency shift. On the opposite, the conical damper has a high level of damping as it was expected from the analytical results. The energy dissipated by friction for the different dampers are almost compared and confirmed the pure rolling motion of the cylindrical damper in in-phase case. The energies dissipated by the wedge and the conical damper are of the same order of magnitude. As a conclusion, the conical damper has a high damping potential for both in-phase and out-of-phase motion. The contact conditions are more robust than for the wedge damper since the contact corresponds to a contact line that can be easily located. Moreover, by increasing the cone angle, the coupling effect between the two platforms is increased, therefore, the frequency shift caused by the damper can be controlled by a simple controllable parameter. Furthermore, from simple design parameters (as the cone angle), a very fine trade-off between the damping efficiency and the wear rate can be achieved.

FIGURE 2 :

 2 FIGURE 2: Illustration of the new damper geometry in contact with the platforms and the angles involved in the geometrical problem

FIGURE 3 :

 3 FIGURE 3: Damper inclination γ versus the damper cone angle for different platform angles β : β 1 = 36 • (blue)β 2 = 45 • (orange)β 3 = 60 • (yellow)

FIGURE 4 :

 4 FIGURE 4: Kinematic slip determination and position of the damper

FIGURE 5 :

 5 FIGURE 5: Evolution of the kinematic slip at each point of the damper for different cone angleα c = 5 • (blue)α c = 10 • (orange)α c = 15 • (yellow)

FIGURE 6 :

 6 FIGURE 6: Influence of different parameters on the maximum kinematic slip with regard to the cone angle α c

  (a) Mesh of the platforms model (b) In-phase mode f IP = 269.9 Hz (c) Out-of-phase mode f OOP = 269.

FIGURE 7 :

 7 FIGURE 7: Finite element model of the blades model (a) and the first two modes of the structure (b,c) -FEM of the wedge (d), conical (e) and cylindrical (f) dampers

FIGURE 9 :

 9 FIGURE 9: Nonlinear receptances for the out-of-phase mode: ( ) conical damper -( ) cylindrical damper -( ) wedge damper

FIGURE 10 :

 10 FIGURE 10: Nonlinear receptances for the in-phase mode: (left) conical damper -(middle) cylindrical damper -(right) wedge damper

FIGURE 11 :

 11 FIGURE 11: Maximum of the amplitude response for different level of excitation for the different dampers: (blue) conical damper -(orange) cylindrical -(yellow) wedge damper

FIGURE 12 :FIGURE 13 :

 1213 FIGURE 12: Evolution of the contacts status for the cylindrical (top), conical (middle) and wedge (bottom) dampers at different frequency levels for an excitation of 9.6 N (left) and 30 N (right): (•) stuck -(•) stick/slip -(•) contact separation -() separation between contact points on the left side (bottom) and on the right side (top)

TABLE 1 :

 1 Energy dissipated by friction in mJ for the different dampers and the different excitation levels for the out-of-phase mode

	Force (N)	0.096	3.84	9.6	15	20	25	30
	Wedge	2.3e-22 1.4e-3	1.2e-2 2.6e-2 3.9e-2 5.2e-2 6.6e-2
	Cylindrical 1.3e-21 1.8e-18 1.1e-2 2.4e-2 3.5e-2 4.7e-2 5.8e-2
	Conical	8.1e-22 1.7e-18 1.1e-2 2.3e-2 3.5e-2 4.6e-2 5.8e-2

TABLE 2 :

 2 Energy dissipated by friction in mJ for the different dampers and the different excitation levels for the in-phase mode Cylindrical 6.2e-23 1.7e-19 1.6e-18 1.0e-18 4.1e-18 6.2e-18 9.5e-18

	Force (N)	0.096	3.84	9.6	15	20	25	30
	Wedge	4.7e-21 3.6e-3	3.3e-2	8.1e-2	1.4e-1	2.3e-1	3.3e-1
	Conical	3.9e-20 1.8e-2	5.9e-2	1.0e-1	1.4e-1	1.8e-1	2.2e-1

ACKNOWLEDGEMENT

The third and fourth authors are grateful to Innovate UK (Grant no. MEDY.P50254) and Rolls-Royce plc (Grant no. MEDY.P42978) for providing the financial support for this work and for giving permission to publish it. This work was part of a collaborative R&T project SILOET II P19.6 which is co-funded by Innovate UK and Rolls-Royce plc. The first and third authors thank Rolls-Royce plc and the EPSRC for the support under the Prosperity Partnership Grant "Cornerstone: Mechanical Engineering Science to Enable Aero Propulsion Futures", Grant Ref: EP/R004951/1.


Andrea Ferrero
email: andrea_ferrero@polito.it

Angelo Iollo

Francesco Larocca

Field inversion for data-augmented RANS modelling in turbomachinery flows

Keywords: Field inversion, Machine learning, Turbulence modelling, Turbomachinery

The field inversion approach is investigated for improving RANS models in turbomachinery flows

• Working conditions characterised by transition and separation are considered • Some approaches to improve the robustness of the method are proposed

• The predictive ability of the method is investigated for several working conditions on different geometries

Introduction

The recent trends in the simulation of the flow field inside aerospace propulsion systems are characterised by a growing interest towards highfidelity simulations which have become feasible thanks to a significant increase in the available computational power. This paves the way to the possibility of understanding complex physical effects which characterise turbulence and combustion phenomena in modern engines. The ability to understand and control these effects can be exploited to increasing the performance and reduce the emissions of existing propulsion systems.

However, scale-resolving simulations (like for example Direct Numerical Sim-ulations (DNS) or Large Eddy Simulations (LES)) cannot be easily integrated in the design process of industrial components. This is due to two main reasons: computational cost and difficulty to manage the results. It is clear that in the first steps of a design process several configurations must be investigated and the use of high-fidelity simulations in this phase would have a prohibitive cost. For this reason, less expensive approaches like RANS equations will be probably used for several years. As far as the management of the results is concerned, LES and DNS usually generate a huge amount of data for each simulation: in order to extract the useful information required by the design process it would be necessary to perform a complex post-process step.

For example, even the computation of the average field from unsteady DNS data is not trivial because it is not known a-priori the extension of the time window required to get statistically converged results: several examples of low frequency phenomena which make difficult to compute the average field can be found in the literature, even looking to simple test cases, and special strategies to estimate the statistical error should be used [START_REF] Russo | A fast algorithm for the estimation of statistical error in dns (or experimental) time averages[END_REF]. A review of the current state of the art for high-fidelity simulations in turbomachinery was proposed by Sandberg and Michelassi (2).

Recently, several research efforts have been devoted to the development of machine learning algorithms for all those applications in which a large amount of data must be processed. In particular, several recent works in the literature have been devoted to the use of machine learning techniques to analyse high-fidelity data from experiments or high-fidelity numerical simulations. The idea behind most of these recent works is to get the physical insight hidden in the data and use it to develop or improve low order datadriven models. An example of this philosophy is represented by the work of Xie et al. [START_REF] Xie | Data-driven filtered reduced order modeling of fluid flows[END_REF] who proposed a filtered reduced order model with a datadriven closure. Dupuis et al.(4) proposed an approach in which traditional surrogate models and machine learning are combined to improve the prediction of the flow on airfoils which work in subsonic or transonic conditions. Margheri et al. [START_REF] Margheri | Epistemic uncertainties in rans model free coefficients[END_REF] performed a study on the epistemic uncertainty of some popular RANS models and used a generalised Polynomial Chaos response surface to perform the calibration of the model coefficients in the spirit of data assimilation strategies. In [START_REF] Ferrero | Reduced order modelling for turbomachinery shape design[END_REF] the Proper Orthogonal Decomposition approach is used in a discontinuous Galerkin (DG) finite element framework [START_REF] Ferrero | Global and local pod models for the prediction of compressible flows with dg methods[END_REF] together with a domain decomposition strategy [START_REF] Bergmann | A zonal galerkin-free pod model for incompressible flows[END_REF] to learn empirical local bases which are used to reduce the simulation cost of the flow field in gas turbines. An alternative path was followed by Raissi and Karniadakis [START_REF] Raissi | Hidden physics models: Machine learning of nonlinear partial differential equations[END_REF] who proposed an approach to identify the partial differential equations which govern a set of data: they applied the algorithm to an example in which they recovered the Navier-Stokes equations used to generate the database but the same approach could be used on experimental data to recover turbulence models.

While the work of Raiss and Karniadakis [START_REF] Raissi | Hidden physics models: Machine learning of nonlinear partial differential equations[END_REF] aims at discovering the full governing model, several works focus on the improvement of existing models.

For example, Wang et al. [START_REF] Wang | Physics-informed machine learning approach for reconstructing reynolds stress modeling discrepancies based on dns data[END_REF] developed a machine learning strategy to predict the discrepancy in RANS modelled Reynolds stresses starting from DNS data. Weatheritt et al. [START_REF] Weatheritt | A novel evolutionary algorithm applied to 44 algebraic modifications of the rans stress-strain relationship[END_REF] proposed the use of Gene Expression Programming to identify new expressions for the stress-strain relationship. Promising results were obtained with this technique on high pressure turbines [START_REF] Weatheritt | Machine learning for turbulence model development using a high-fidelity hpt cascade simulation[END_REF]. Duraisamy et al. (13;[START_REF] Parish | A paradigm for data-driven predictive modeling using field inversion and machine learning[END_REF] proposed a strategy based on field inversion and machine learning which allows to improve the prediction ability of RANS models. This approach is exploited in the present work in order to improve RANS modelling for low pressure gas turbine cascades.

Machine learning techniques have been investigated also on multiphase flows (15; 16), combustion (17; 18; 19) and engine modelling (20; 21). Finally, a comprehensive review of the machine learning techniques proposed for the improvement of turbulence modelling can be found in [START_REF] Duraisamy | Turbulence modeling in the age of data[END_REF].

The paper is organised as follows. In Section 2 the original RANS model is presented. In Section 3 the methods used for the discretisation of the equations are described. In Section 4 the field inversion approach is described and it is then applied to the T106c gas turbine cascade in Section 5. The data obtained by the field inversion are analysed by means of machine learning techniques in Section 6 in order to generalise the obtained results. Finally, the improved RANS model is tested on the T106c and on the T2 cascades in Section 7.

Physical model

This work is devoted to the prediction of the compressible turbulent flow in 2D turbine cascades. The study starts from the Spalart-Allmaras (SA) model implemented for compressible equations, following the guidelines of [START_REF] Allmaras | Modifications and clarifications for the implementation of the spalart-allmaras turbulence model[END_REF]. This model is widely used in the literature for fully turbulent flows.

However, the model is not suitable for the prediction of transitional flows at low Reynolds numbers. The original model gives the possibility to impose the transition location (by means of the trip term f t1 defined in [START_REF] Allmaras | Modifications and clarifications for the implementation of the spalart-allmaras turbulence model[END_REF]) but this choice is rarely followed in the literature because in general the location of transition in not known a-priori. Furthermore, when the transition trip term f t1 is used a second term f t2 for delaying natural transition (and making the trip term f t1 effective) is also activated. Further details on the effects of the term f t2 in the prediction of the flow around the T106c cascade can be found in [START_REF] Ferrero | Rans closure by artificial neural networks[END_REF].

In the present work the SA model is used without the trip terms f t1 and f t2 . With this choice the model is expected to work fine for high Reynolds numbers but to fail in predicting transition and separation at low values of Reynolds number. This model tends indeed to produce an excessive amount of turbulent eddy viscosity on this kind of flows [START_REF] Ferrero | Rans closure by artificial neural networks[END_REF]. For this reason, it represents an optimal baseline for testing the field inversion approach and evaluating how much the original model can be improved.

The mass-averaged RANS equations are reported in the following:

∂ρ ∂t + ∇ • (ρu) = 0 (1) ∂ ∂t (ρu) + ∇ • (ρuu) = -∇p + ∇ • τ (2) ∂E ∂t + ∇ • (u(E + p)) = ∇ • (τ • u -q) (3) ∂ρν ∂t +∇•(ρuν) = ρ(P -D)+ 1 σ ∇•(ρ(ν + ν)∇ν)+ c b2 σ ρ(∇ν) 2 - 1 σ (ν + ν)∇ρ•∇ν (4)
where ρ, u, p, E, ν, ν, x and t are density, velocity, pressure, total energy per unit volume, molecular viscosity, modified eddy viscosity, spatial position and time, respectively. A fluid with constant specific heat ratio γ and constant viscosity is considered. The following equation for the energy is considered:

E = P γ -1 + 1 2 ρu • u ( 5
)
where γ is the specific heat ratio.

The viscous stress tensor τ includes both the molecular and eddy viscosity contributions and its components are given by:

τ ij = 2ρ(ν + νf v1 ) 1 2
∂u i ∂x j + ∂u j ∂x i - 1 3
∂u k ∂x k δ ij (6)
The production P and destruction D terms in Eq. 4 are computed as follows:

P = c b1 S ν D = c w1 f w ν d 2 (7)
with the following definitions:

f w = g 1 + c 6 w3 g 6 + c 6 w3 1/6 g = r + c w2 (r 6 -r) r = min ν S2 κ 2 d 2 , r lim (8)
S =      S + S if S >= -c v2 S S + S(c 2 v2 S+c v3 S) (c v3 -2c v2 )S- S if S < -c v2 S ( 9
)
where S is the vorticity magnitude and S is:

S = ν κ 2 d 2 f v2 (10)
The functions f v1 and f v2 depend on the viscosity ratio χ = ν ν :

f v1 = χ 3 χ 3 + c 3 v1 f v2 = 1 - χ 1 + χf v1 (11)
The constants σ, c b1 , c b2 , c v1 , c w1 are defined in [START_REF] Allmaras | Modifications and clarifications for the implementation of the spalart-allmaras turbulence model[END_REF].

Finally, the heat flux q is described by the Fourier's law:

q = - c p µ P r + c p ρνf v1 P r t ∇T (12)
where T , c p , P r and P r t are the temperature, the constant pressure specific heat capacity, the Prandtl number and the turbulent Prandtl number. The test cases considered in this work refer to experiments performed with air and so the following values are assumed: γ = 1.4, P r = 0.72 and P r t = 0.9.

Implicit Discontinuous Galerkin discretization

The discontinuous Galerkin (DG) scheme is used in this work for the spatial discretisation on the governing equations. This approach is characterised by a significant flexibility since it allows to easily manage high-order reconstructions on unstructured meshes. The main idea behind this kind of scheme consists in adopting an high-order polynomial reconstruction inside each element without any continuity constraint at the interface between different elements. As a result, the scheme can be easily exploited in the framework of automatic adaptive approaches, in which both the size (h-adaptivity, (25; 26; 27; 28)), the order (p-adaptivity (29; 30; 31)) or both properties (hp-adaptivity, (32; 33; 34)) can be locally adapted following some error indicators.

The computational domain Ω is discretised with a hybrid mesh which contains a structured boundary layer mesh close to the body surrounded by an unstructured mesh. The grid is generated by means of Gmsh [START_REF] Geuzaine | Gmsh: A 3-d finite element mesh generator with built-in pre-and post-processing facilities[END_REF] with the Frontal-Delaunay for Quads algorithm. The management of the unstructured grid in the parallel MPI environment is performed through the DMPlex class [START_REF] Lange | Flexible, scalable mesh and data management using petsc dmplex[END_REF] provided by the PETSc library [START_REF] Balay | Petsc users manual[END_REF].

The numerical approximation of the l-th conservative variable u l (x, t) inside each element Ω e is described by a modal basis with size N e = (k+1)(k+2) 2 with a reconstruction order k:

u l (x, t) = Ne i=1 ũli (t)φ i (x) 1 ≤ i ≤ N e ( 13
)
where ũli (t) ∈ R Ne contains the degrees of freedom inside the element for the l-th conservative variable. The basis functions φ i (x) are obtained by the modified Gram-Schmidt orthonormalisation applied to a set of monomials defined in the physical space, following the approach of Bassi et al. [START_REF] Bassi | On the flexibility of agglomeration based physical space discontinuous galerkin discretizations[END_REF]. In this work a third order accurate DG scheme is used (k = 2, N e = 6).

The spatial discretisation is completed by a projection of the governing equation on the space of the approximation functions. The resulting weak formulation consists in a set of ordinary differential equations in time. The convective terms which appear in the numerical fluxes at the interface between the elements are evaluated by means of an approximate Riemann problem solver (following (39) and ( 40)). Diffusive terms are evaluated by means of a recovery-based approach (41).

Time integration is here performed by means of the linearised implicit Euler method. Since steady problems are considered the use of a first order time integrator appears suitable since it does not influence the accuracy of the final steady solution and it has good dissipative properties which are useful to accelerate the numerical transients. The solution of the linear system which is obtained at each time step is performed in parallel by means of the GMRES algorithm with the additive Schwarz preconditioner provided by the PETSc library [START_REF] Balay | Petsc users manual[END_REF]. The GMRES algorithm is employed by setting the maximum number of iterations to 200, the dimension of the Krylov subspace to 100 and the absolute tolerance to 10 -12 . The CFL number which controls the time step size is automatically adjusted according to the evolution of the residuals following the pseudo-transient continuation strategy [START_REF] Bassi | Very high-order accurate discontinuous galerkin computation of transonic turbulent flows on aeronautical configurations[END_REF]. In particular, the CFL number is allowed to vary between 10 2 and 10 4 . During the first steps of the transient, a feedback filtering procedure ( 43) is applied to remove potential instabilities which can appear due to the large CFL number. This filtering procedure is deactivated when the residuals drop under a certain threshold and so it does not influence the steady solution.

Field inversion and machine learning in a DG framework

The field inversion approach proposed by [START_REF] Parish | A paradigm for data-driven predictive modeling using field inversion and machine learning[END_REF] requires to define a goal function G which measures the distance between the experimental data and the predicted numerical results. The procedure requires the solution of an optimisation problem in which a field β(x) is found in order to minimise the goal function G. The field β(x) is then introduced in a correction term h(β(x)) which multiplies the production term in the SA transport equation:

∂ρν ∂t +∇•(ρuν) = ρ [h(β)P -D]+ 1 σ ∇•(ρ(ν+ν)∇ν)+ c b2 σ ρ(∇ν) 2 - 1 σ (ν+ν)∇ρ•∇ν (14)
In the original works of (14; 44) the correction was chosen as h(β) = β. In this work, different choices are investigated for the function h(β), as described in the next section.

As far as the goal function G is concerned, the following choice is made:

G = w (M s -M exp s ) 2 dl + λ Ω (β -1) 2 dΩ (15)
The first term is a line integral performed on the wall of the blade and allows to evaluate the norm-2 error on the wall isentropic Mach number distribution M s , which is defined as:

M s = 2 γ -1 [(p 0 i /p w ) (γ-1)/γ -1] ( 16
)
where p w is the static pressure at wall and p 0 i is the inlet total pressure.

The second term is a surface integral on the computational domain Ω which acts as a Tikhonov regularisation [START_REF] Bishop | Training with noise is equivalent to tikhonov regularization[END_REF]: it penalises the goal function when the correction factor is far from 1. This is useful to avoid unnecessary corrections which could be introduced during the optimisation process but which are not required in the final optimal solution. The choice of the penalisation constant λ will be discussed in the next Section.

In order to solve the optimisation problem, a simple gradient descent method is applied. The field β will be described in terms of the same basis functions used for the conservative variables. Starting from the original SA model (h(β(x)) = 1) the degrees of freedom related to the field β are updated with the gradient descent method:

β = β -δ dG d β ( 17
)
where δ is the step size that in this work is chosen constant for simplicity (δ = 0.1 ).

Since the dimension of the optimisation problem is related to the total number of degrees of freedom per equation the computation of the gradient dG d β by means of numerical differentiation would be prohibitive. For this reason, an adjoint-based gradient evaluation was implemented. The gradient of the goal function G which respect to the degrees of freedom of the field β(x) is computed as:

dG d β = ∂G ∂ β + ψ T ∂R ∂ β ( 18
)
where R represents the residual of the governing equations. The first term contains only the contributions related to the penalisation integral which appears in the goal function. The adjoint variable Ψ is computed by the solution of the following linear system with the GMRES iterative solver:

∂R ∂ ũ T Ψ = - ∂G ∂ ũ T (19)
in which the jacobian matrix ∂R ∂ ũ is already available from the implicit time integrator and the term ∂G ∂ ũ contains the derivatives of the goal function with respect to the fluid dynamics degrees of freedom. This last term was computed by means of automatic differentiation with the Tapenade tool [START_REF] Hascoet | The tapenade automatic differentiation tool: Principles, model, and specification[END_REF].

Summarising, the procedure works as follows. First of all, a steady solution with the original SA model is obtained. The solution is considered steady when the residuals of all the governing equations are lower than 10 -6 . Usually, the SA equation is the one which converges with the lowest speed so when the condition is satisfied the residuals of the Eqs. 1-3 are orders of magnitudes lower (typically around 10 -8 -10 -10 ). When the steady solution is reached, the gradient dG (14; 44). This last approach is followed in the present work.

In particular, the results of the inverse problem will be exploited to train an Artificial Neural Network (ANN) which can then be used to define an augmented version of the SA model.

Field inversion on the T106c cascade

The field inversion approach is here applied to the flow around the T106c The presence of separation can be easily noticed in the experimental studies on these flows by checking the wall isentropic Mach number distribution: the separation is usually related to the presence of a plateau in the distribution.

Singh et al. [START_REF] Singh | Machine-learning-augmented predictive modeling of turbulent separated flows over airfoils[END_REF] showed that the wall pressure distribution (which is directly related to the isentropic Mach number distribution) can be effectively used in the field inversion approach for improving the prediction of separated flows. They indeed showed that the field inversion based on the wall pressure distribution can significantly improve the prediction of the Reynolds stresses in the separation region [START_REF] Singh | Machine-learning-augmented predictive modeling of turbulent separated flows over airfoils[END_REF]. For these reasons, the field inversion algorithm used in this work will use the error on the wall isentropic Mach number distribution as goal function.

First of all, a convergence study is performed on the T106c cascade with the original SA model at the highest Reynolds number (Re 2s = 2.5 • 10 5 ). Three different meshes and two reconstruction orders (1 ≤ k ≤ 2) are evaluated.

The convergence level is assessed by checking the mass averaged value of the kinetic losses in a control section located 0.465c x behind the trailing edge.

The kinetic losses are defined in the following way:

ζ = 1 - 1 -(p e /p 0 e ) (γ-1)/γ 1 -(p e /p 0 i ) (γ-1)/γ (20)
where p e , p 0 e and p 0 i are the static pressure in the control section, the total pressure in the control section and the inlet total pressure, respectively. The results of the convergence analysis are reported in Table 1 which shows the number of elements n ele , the number of degrees of freedom per equation n DOF and the predicted averaged losses. It is useful to remember that in the asymptotic range mesh refinement gives a fixed convergence order (depending on k) while order refinement gives exponential convergence.

We emphasise that the losses in the wake represent a better goal function for the convergence assessment with respect to the wall isentropic Mach number distribution because the original SA model over-predicts significantly the turbulence eddy viscosity and so it gives a wall isentropic Mach number distribution which is very similar to what would be obtained by an inviscid Euler simulation, regardless of the mesh resolution. In contrast, the wake losses are influenced by the mesh resolution in the boundary layer and in the wake region.

The mesh C reported in Tab.1 will be used for all the following simulations with a third order accurate DG scheme (k = 2). The mesh contains 40436 elements and so the total number of degrees of freedom per equation is equal to 242616. The dimensionless wall cell size is y + < 1 on the entire surface. As reported in Equation 14, the field inversion approach requires to alter the production term by the presence of the correction factor h(β). In this work, different expressions for h(β) are investigated. The most straightforward approach, which was used by Singh et al. [START_REF] Singh | Machine-learning-augmented predictive modeling of turbulent separated flows over airfoils[END_REF] for the study of wind turbine airfoils, consists in setting :

n ele n DOF ζ Mesh A, k=1
h(β) = β β ∈ R (21)
In this way the correction factor is free to assume both positive and negative values and so the correction term is very general. However, this generality comes with a price: since h(β) is not limited it can lead to unstable numerical results during the transients which must be solved in predictive simulations.

An alternative approach, experimented in this work, consists in setting

h(β) = β 2 β ∈ R (22)
In this way the correction term is not allowed to assume negative values.

This means that the generality of the approach is reduced but the robustness of the simulation is increased because the correction term cannot change the nature of the production term (it can, in the limit, set the production to zero but it cannot transform the production term into a destruction term).

A third approach, which showed the most robust results in this work, is reported in the following. The idea behind this approach is to mimic the behaviour of intermittency models in which the production term of the RANS model is reduced by a factor defined in the range [0, 1] in order to reproduce transition phenomena. Following this approach, the correction term is defined as a smooth ramp function of β:

h(β) =              0 if β ≤ 0 3β 2 -2β 3 if 0 < β < 1 1 if β ≥ 1 ( 23
)
This last approach is the least general between the three alternatives examined in this work but it is the most robust. This is due to the fact that, in the end, the correction factor h will be expressed by means of an ANN.

When the SA model augmented by the ANN correction term will be used for actual predictions, the ANN will be asked to compute the correction factor for input values which could be outside of the range explored in the training database. This is very likely to happen during the numerical transient which must be solved before getting the steady solution. However, ANNs are known for their poor extrapolation accuracy and so the use of a more general expression (like for example the one defined by Equation 21) would allow the presence of unlimited values of the correction factor. In contrast, when the correction factor is limited in the range 0 ≤ h ≤ 1 the model can behave, in the limit, as the original SA model (when h → 1) or as the laminar Navier-Stokes equations (when h → 0).

In order to understand whether the limitation introduced by Equation 23 The optimal field obtained from this first step is then used as initial field for a second optimisation in which the correction factor is limited according to Eq. 23. It is useful to emphasise that, in order to apply the correction factor defined by Eq. 23, it is not possible to start with a uniform field with β = 1.

This is due to the fact that the derivative of the smooth ramp function is null for β = 1 and so it would not be possible to update the solution since the gradient of the goal function would remain to zero according to Eq. 18

( ∂R ∂ β = ∂R ∂g ∂g ∂ β , with ∂g ∂ β = 0 for β = 1).
In order to compare the two approaches, the wall isentropic Mach number distribution is reported in Figure 2 for the original SA model and the optimised solutions related to Eq. 21 and 23. The results for the correction factor defined by Eq.22 are not reported in the plot since they overlap the other results related to 21 and 23. The Figure shows also the available experimental data which are used to drive the optimisation process. The optimal solutions show a good match with the experimental data and a significant improvement with respect to the baseline model. This test confirms that the limited correction factor defined by Eq. 23 is able to provide an optimal solution which is comparable to the results provided by the unlimited correction factor. This is due to the fact that the original SA model overestimates significantly the turbulence production in this kind of flows and so the use of a correction factor limited between 0 and 1 is sufficient to correct the model.

In this sense, the correction factor proposed in this work acts exactly as a intermittency correction in the framework of laminar-to-turbulence transition. After this analysis, the limited correction factor defined by Eq. 23 was chosen for all the following simulations.

The plots in Figure 3 shows a small separation bubble followed by reattachment. Again, this separation is missed by the original SA model.

Finally, the correction field at Re 2s = 8 • 10 4 and Re 2s = 2.5 • 10 5 for the case defined by Eq. 23 is reported in Figures 4 and5 for λ = 0 and λ = 10 -3 , re-

spectively. An analysis of the pictures shows clearly that the adjoint approach obtained an optimal solution in which the production term is deactivated in the boundary layer for the first portion of the suction side: the algorithm has recovered a laminar separation just by using the knowledge on the experimental wall isentropic Mach number distribution. As far as the influence of λ is concerned, a study with λ = 0, 10 -2 , 10 -3 , 10 -4 is performed. These values are chosen by running a preliminary simulation with λ = 0 and then evaluating the order of magnitude of the two integrals which appear in the goal function defined by Eq. 15. For all these values, the optimal wall isentropic Mach number distribution does not show significant variations. The weak influence of the parameter λ can be seen in Figures 4 and5 where the higher value of λ tends to avoid unnecessary corrections at the end of the separation region.

Machine learning on the T106c cascade

The field inversion algorithm described in the previous section is able to provide a correction field which alters the original SA model in order to match very well the experimental results for two different working conditions. In this section this result will be generalised in order to express the correction factor as a function of some physical features. In particular, several choices related to the inputs and the architecture of the ANN used to express the correction factor will be investigated.

Choice of the inputs

The choice of the input variables of the ANN is not a trivial task. In particular, it is necessary to avoid input variables which would introduce a dependency on the particular frame of reference which is used to study the problem (i.e. Galilean invariance must be satisfied). Furthermore, there should not be strong correlations between the different input variables and they should be chosen as adimensional quantities in order to get general results.

A natural choice is to identify some adimensional groups which appear in the source term of the original RANS model and use them as input for the ANN. This choice was for example carried out by Singh et al. [START_REF] Singh | Machine-learning-augmented predictive modeling of turbulent separated flows over airfoils[END_REF].

A similar approach is used in this work but particular attention is here de- 44) and originally proposed by [START_REF] Spalart | A new version of detached-eddy simulation, resistant to ambiguous grid densities[END_REF] in the framework of Detached Eddy Simulations. The terms are defined as:

f d = 1 -tanh((8r d ) 3 ) f d = 1 -tanh((r d ) 0.5 ) ( 24
)
where the quantity r d is an adimensional group obtained by combining wall distance, turbulence and molecular viscosity and velocity gradient:

r d = ν + ν d 2 κ 2 ∂u i ∂x j ∂u i ∂x j ( 25
)
where κ = 0.41 is the von Karman constant.

The plot in Figure 6 explains why in this work the term f d is used instead of f d : both terms are limited between 0 and 1 but f d allows to better describes the flow features close to wall while f d tends to compress the information and does not allow to distinguish the different structures. This qualitative analysis was confirmed by quantitative analysis which shows that an ANN with f d was able to better fit the database with respect to an equivalent ANN with f d as input.

The term log(P/(D + ) + ) represents a convenient scaling of the ratio between the production P and destruction D terms of the SA model. In the work of Singh et al. [START_REF] Singh | Machine-learning-augmented predictive modeling of turbulent separated flows over airfoils[END_REF] the ratio P/D is directly used while in this work a logarithmic scaling is used: this is due to the fact that the values assumed by this ratio are distributed in a wide range which covers several orders of magnitude and some numerical experiments confirmed that the fitting significantly improves with this scaling. Furthermore, both the numerator and the denominator of this quantity can go to zero in the presence of uniform fields or where the turbulence viscosity is zero and so the constant = 10 -5 is introduced. Some numerical tests showed that the use of logarithmic scaling improves significantly the fitting of the database with the ANN.

Finally, the adimensional gradient of the modified turbulent viscosity log(|∇ν|d/(ν+ ν) + ) is considered. This quantity was not used in [START_REF] Singh | Machine-learning-augmented predictive modeling of turbulent separated flows over airfoils[END_REF] and does not appear in the production and destruction terms. However, it appears in the cross production term (the last term of Eq. 4) and allows to identify regions with strong variations in the eddy viscosity. It is normalised with respect to the wall distance and the sum of kinematic and eddy viscosity: this means that this quantity remains well conditioned even when the eddy viscosity tends to zero since the kinematic viscosity prevents the denominator to become zero. Even for this variable the logarithmic scaling was found to be useful to improve the fitting.

Choice of the ANN architecture

After choosing the input features, it is necessary to define the architecture of the ANN. In this work, feedforward ANNs are considered. As far as the activation functions are concerned, a common choice consists in using sigmoid functions for the hidden layers and linear functions for the output layer. However, since the chosen correction factor h is limited in the range [0, 1] a sigmoid activation function is adopted also for the output layer: in this way the output of the ANN will be automatically limited in the range

[0, 1] .
Particular care should be taken in choosing the number of hidden layers n HL and the number of neurons per layer n N = 10. In particular it is necessary to find a compromise between the complexity of the network (which allows to capture the correlations hidden in the database) and its ability to perform predictions outside of the database. When the complexity of the network is increased its ability to reproduce the training database is enhanced because it has more degrees of freedom which can be adjusted to fit the data.

However, if too many degrees of freedom are introduced then the network will behave poorly during predictions: this is due to the fact that when too many degrees of freedom are used then the output of the network will show strong oscillations for the points in the parameter space which do not exactly match a training point.

In order to find a suitable network by using a general criteria the following approach is used. First of all, different architectures are considered (1 The validation set is used during the training to verify that the ANN is still able to give good predictions for points which do not belong to the training set: when the validation error tends to increase the training is arrested, even if the training error is still decreasing, in order to limit the problem of overfitting. Finally, the test set is used to monitor the behaviour of the ANN on an external set of data which do not influence the training process (neither in the mean squared error computation nor in the validation checks for the overfitting). An example of training history is reported in Figure 7a in which it can be clearly seen that when the training is stopped the training error was still decreasing but the validation error just started to grow. In Figure 8 it is possible to see the regression plots for the different data sets:

≤ n HL ≤ 2 5 ≤ n N ≤ 40)
in each plot the abscissa represents the reference value in the database while the ordinate represents the approximated value computed by the network.

Another approach for avoiding overfitting was also investigated: Bayesan regularisation [START_REF] Foresee | Gauss-newton approximation to bayesian learning[END_REF]. In Bayesan regularisation the mean square error goal function is augmented by a term which penalises large values of the weights.

However, some experiments on the problems considered in this work showed that the splitting of the database in training, validation and test sets allows to achieve a better compromise between fitting and robustness with respect to the Bayesan regularisation.

A sequence of regression plots (on the full database) for the ANN 2 × 5, 2 × 10, 2 × 20 and 2 × 40 are reported in Figure 9: as the complexity of the network is increased its ability to reproduce the database is enhanced as can be cleary seen by the fact that the points tend to assume a distribution centered along the bisector of the quadrant.

In Table 2 According to the previous analysis it would seem that the larger is the network the better is the result. This is true for the fitting of the points in the database. However, it is fundamental to investigate the behaviour of the network for points which do not coincide exactly with the points in the database. In order to do this it is possible to run some CFD simulations at used for the training is built from the optimal solution at these Reynolds number and so one could aspect that the ANN should reproduce perfectly these working conditions. However, it is important to keep in mind that the regression coefficient R is always less than 1: this means that, even if the CFD simulation is initialised with the optimal solution obtained by the adjoint approach, the correction field reproduced by the ANN will not coincide exactly which the optimal one. As a consequence, the CFD solution will evolve towards a new steady solution. This introduces a perturbation in the input features given to the ANN: if the ANN is robust the new steady solution will be close to the optimal one. However, if the ANN is poorly According to this analysis, all the predictive simulations reported in the following will be performed by using the 2 × 20 ANN.

Predictions

In the previous Section the procedure for choosing the architecture of the ANN is reported. Now, the chosen network is use to perform predictive simulation for working conditions and geometries which were not included in the database. As a first step all the simulations are performed by setting h(x) = 1, i.e. with the original SA model. Then the obtained steady solution is used to initialise a simulation in which the correction term is computed with the ANN. This approach speed ups the convergence since the ANN is not employed during the strong initial transient at the beginning of the simulation.

Furthermore, the numerical experiments showed that the robustness of the method during predictive simulations can be improved by limiting the input variables to the range used for the training. This is important because the ANN has been trained only on a few steady solutions and so during the transients which can appear in predictive simulations the input features could assume values which were not observed in the training database. In particular, if h(Y ) represents the ANN approximation of the correction factor and Y is the vector of the five input variables, the modified expression h(L(Y ))is used during predictive simulations, where the limiting function L is defined as:

L(Y i ) =              Y i if Y min i ≤ Y i ≤ Y max i Y max i if Y i > Y max i Y min i if Y i < Y min i (26)
Here Y min As far as the average exit angle is concerned, the ANN augmented SA model shows a better behaviour than the original SA model at low Reynolds numbers while the two models give similar results at higher Reynolds numbers.

It is interesting to note that the asymptotic value of the exit angle for high values of Reynolds number presents an offset between experimental and nu-

T2 cascade

The prediction ability of the ANN augmented SA model is investigated also on another geometry, the T2 cascade. The simulations are carried out with a third order accurate DG scheme on a mesh with 59453 elements, corresponding to 356718 degrees of freedom per equation. The mesh resolution at wall and in the wake region is the same used for the T106c, since both cascades are investigated at similar values of Reynolds number. The T2 airfoil was designed at the VKI for the same velocity triangles of the T106 (inlet angle α = 32.7 o ) but it is characterised by a larger pitch-to-chord ratio (1.05) and an increased diffusion rate along the rear suction side [START_REF] Babajee | Detailed numerical characterization of the separationinduced transition, including bursting, in a low-pressure turbine environment[END_REF]. Also the Zweifel number is larger (Ψ = 1.46) with respect to the T106 (Ψ = 1.24).

The isentropic exit Mach number is set to M 2s = 0.65.

In Figure 13 Finally, the models are evaluated in terms of mass averaged exit kinetic losses and angle, as reported in Figure 16. As observed for the T106c, even in this case the ANN augmented SA model outperforms the original SA model at low Reynolds numbers. It is interesting to note that the numerical results obtained in the present work presents an offset in β 2 with respect to the experimental results, offset which is not observed in the results obtained from the SST-γ-Re θ model. This could be a limitation of the SA model which is inherited by the augmented model: future work will be devoted to apply the field inversion approach to other RANS models to verify whether this limitation persists.

Conclusions

The potential of the field inversion approach was investigated for the augmentation of a RANS model used in the simulation of turbomachinery flows.

In particular the approach was applied to the original SA model and the term is proposed in order to limit its value in a finite range: this, together with the introduction of a limiting on input features, significantly improves the robustness of the approach during transients and in predictions.

The results seem promising and are substantially better than the results provided by the original model. They also appears satisfactory if compared to the results obtained by a significantly more complex four equation model (SST-γ-Re θ ). In particular, even if the goal function used for the field inversion is based only on the wall isentropic Mach number, the ANN augmented model shows improvements also in terms of average losses and exit angle in the wake.

Future work will be devoted to the application of the field inversion approach to other RANS models. Furthermore, possible alternatives to the use of an ANN will be investigated for achieving a better fitting of the database with a good level of robustness in predictions.

dβ

  is computed by the adjoint approach and the correction field is updated. This generates a transient which is solved in time up to a new steady solution. Since the perturbation introduced by the correction update is small, the transient can be easily solved by marching in time with a very large CFL number. For example, in this work the constant value CFL=5000 is used for this part of the computation. The procedure is repeated until the goal function does not show any significant improvement.The correction field h(β(x)) obtained by the inversion process can be exploited for different purposes. On one hand, it gives insight for the development of new turbulence models since it shows where and how the original model fails. On the other hand, it is possible to directly generalise the correction in order to obtain a new model which can be used for predictive simulations. For example, Duraisamy and Durbin (47) used the results of field inversion to define a transport equation for an intermittency factor, where the different terms of the transport equation are computed by means of machine learning techniques. Alternatively, it is possible to find a local closure which allows to define the correction field as a function of local physical quantities

2 . 5 •

 25 gas turbine cascade. This profile is representative of high-lift low pressure gas turbines in modern turbofan engines. The cascade was experimentally investigated at the VKI and some experimental results are available from the literature(48; 49; 50). In particular, the wall isentropic Mach number distribution, the mass averaged kinetic losses and exit angle in the wake are available for several values of the Reynolds number. The flow field is studied for an inlet angle α = 32.7 o , an isentropic exit Mach number M 2s = 0.65 and different values of the exit isentropic Reynolds number 8 • 10 4 ≤ Re 2s ≤ 10 5 . The Reynolds number Re 2s is defined by using the blade chord and the isentropic exit velocity and density. The dynamic viscosity is assumed constant. The turbulence intensity during the experiments was very low (0.9%): for this reason all the RANS simulations are performed by setting a very small value of inlet eddy viscosity (ν/ν = 0.1).Houmorziadis[START_REF] Hourmouziadis | Aerodynamic design of low pressure turbines[END_REF] showed that the Reynolds number in low pressure gas turgines of turbofan engines is the range between 10 5 -4 • 10 5 where the smaller values are observed in cruise conditions and the higher values are obtained at take-off. The high-lift profiles can show large laminar separations at low values of Reynolds number. When the Reynolds number is increased the separation transforms from an open separation to a closed separation in which there is a separation bubble followed by reattached flow. The evolution from one configuration to the other takes place in a small range of Reynolds number and so the flow is quite sensitive to the working condition.

  affects the ability of the field inversion to match the experimental data, the different definitions of h(β) are tested on the T106c cascade. In particular, the gradient based optimisation process is carried out for the T106c at Re 2s = 8 • 10 4 and Re 2s = 2.5 • 10 5 . The plot in Figure1shows the history of the goal function during the optimisation process. The results shows that after approximately 50 steps of the gradient descent algorithm a minimum is reached. This optimisation is carried out by starting from the original SA model with h = 1 in all the domain and using the unlimited correction factor defined by Eq. 21 with λ = 0.

  show the Mach field for the original SA model and optimal model at Re 2s = 8 • 10 4 and Re 2s = 2.5 • 10 5 . The optimal solution at Re 2s = 8 • 10 4 is characterised by a large open separation which is completely missed by the original SA model. The optimal solution at Re 2s = 2.5 • 10 5

Figure 1 :

 1 Figure 1: Adjoint-based optimization history for T106c at Re 2s = 8 • 10 4 with β ∈ R

Figure 2 :

 2 Figure 2: Comparison between original SA model, optimized model and experimental results in terms of M is distribution for the T106c at Re 2s = 8 • 10 4 (a) and Re 2s = 2.5 • 10 5 (b)

Figure 3 :

 3 Figure 3: Mach field for T106c with the original SA model (a,c) and with optimised model (b,d) at Re 2s = 8 • 10 4 (a,b) Re 2s = 2.5 • 10 5 (c,d)

Figure 4 :

 4 Figure 4: Correction field h(x) for T106c at Re 2s = 8 • 10 4 (a) and Re 2s = 2.5 • (b) with λ = 0

Figure 5 :

 5 Figure 5: Correction field h(x) for T106c at Re 2s = 8 • 10 4 (a) and Re 2s = 2.5 • (b) with λ = 10 -3

  voted to the robustness and the prediction ability of the model. The following five input variables are used: χ, log(τ /τ ref + ), f d , log(P/(D + ) + ) and log(|∇ν|d/(ν + ν) + ). The plots in Figure 6 show the distribution for all the inputs variables in the optimised solution at Re 2s = 8 • 10 4 . The first input, χ, simply represents the turbulent intensity. The quantity τ /τ ref is obtained by normalising the module of the stress tensor with respect to a reference stress. The reference stress is here defined as τ ref = ρ(ν+ν) 2 /d 2 which makes this input a local quantity. In contrast, Singh et al. (44) used a non local normalisation in which the stress tensor is normalised with respect to the wall stress τ w . However, such non-local terms are avoided in this work since the presence of non-local terms reduces significantly the scalability of the discretisation in a parallel environment. Furthermore, the physical meaning of using τ w for the normalisation is clear for the mesh points in the boundary layer but is not so clear for other regions, like for example the wake. Finally, a logarithmic scaling of the quantity τ /τ ref was observed to significantly improve the fitting of the database. The additive constant = 10 -5 is introduced to prevent the algorithm of the logarithm to become null. The term f d is introduced in this work as a modification of the term f d used by Singh et al. (

Figure 6 :

 6 Figure 6: Input features for the neural network: f d (a), f d (b), χ (c), log(τ /τ ref + ) (d), log(P/(D + ) + ) (e), log(|∇ν|d/(ν + ν) + ) (f)

  and the ability of the networks to fit the database is investigated. Each network is trained in Matlab by means of the Levenberg-Marquadt algorithm with a goal function based on the mean squared error. The training is performed by dividing randomly the database in 3 subsets: one for training (70% of the data), one for validation (15% of the data) and one for test (15% of the data). The training set is actually used for the computation of the mean square error and for driving the training process.

Figure 7 :

 7 Figure 7: Training history (a) and architecture (b) for 2x20 ANN

Figure 8 :

 8 Figure 8: Training, validation and test error for 2x20 ANN

Figure 9 :Figure 10 : 20 and 2 ×

 910202 Figure 9: Regression plots for different ANN architectures: 2x5 (a), 2x10 (b), 2x20 (c), 2x40 (d)

i 7 . 1 .

 71 and Y max i represent the minimum and maximum values of the i-th input feature observed in the training database. T106c cascade at different Reynolds number As a first test, the ANN aumented SA model is used to predict the flow field on the T106c at Re 2s = 1.2 • 10 5 , 1.6 • 10 5 and 2.1 • 10 5 . In this range of Reynolds number a strong variation is observed in the solution due to the transition from open to closed separation. The results related to the wall isentropic Mach number distribution are reported in Figure 11 in which they are compared with the available experimental results and the original SA model. The ANN augmented SA model performs significantly better than the original model and the predictions are quite close to the experiment. Only the solution obtained at Re 2s = 1.2 • 10 5 seems to overpredict the separation.The results reported in Figure11refer to the M is distribution used in the goal function which drove the field inversion and so it is natural to expect an improvement with respect to the original model. However, the prediction ability of the model was also investigated in terms of mass averaged kinetic losses ζ and exit angle β 2 in the wake, quantities which were not included in the goal function used for the optimisation.The average is performed in a control section located 0.465c x behind the trailing edge, where c x is the axial chord, in the same location used for the experimental measurements. The results of these tests are reported in Figure12. As far as the losses are concerned, both the original SA model and the ANN augmented SA model perform well for high Reynolds values. However, for low Reynolds numbers the original SA model misses completely the separation and so it underpredicts significantly the losses. The ANN augmented SA model shows the correct trend and is quite close to the experimental re-

Figure 11 :

 11 Figure 11: Wall isentropic Mach number distribution: predictions at Re 2s = 1.2 • 10 5 (a), Re 2s = 1.6 • 10 5 (b) and Re 2s = 2.1 • 10 5 (c)

Figure 12 :

 12 Figure 12: Average losses and exit angle for T106c cascade: comparison between original SA model, SA-ANN model and experimental results

  and 14 the Mach number field at Re 2s = 1.2 • 10 5 and 2.1 • 10 5 is reported for the original SA model and for the ANN augmented SA model. The plots show clearly the presence of a open separation at Re 2s = 1.2 • 10 5 and a closed separation at Re 2s = 2.1 • 10 5 .Finally, in Figure15the predicted wall isentropic Mach number distribution is reported as a function of the curvilinear coordinate s along the blade surface, normalised with respect to the curvilinear length of the blade (s 0 ). The ANN augmented SA model shows significant improvements with respect to the baseline SA model and gives good results also with respect to the SSTγ-Re θ results from[START_REF] Babajee | Detailed numerical characterization of the separationinduced transition, including bursting, in a low-pressure turbine environment[END_REF].

Figure 13 :

 13 Figure 13: Mach field for the T2 cascade at Re 2s = 1.2 • 10 5 with the original SA model (a) and with the ANN-SA model (b)

Figure 14 :Figure 15 :

 1415 Figure 14: Mach field for the T2 cascade at Re 2s = 2.1 • 10 5 with the original SA model (a) and with the ANN-SA model (b)

Figure 16 :

 16 Figure 16: Average losses and exit angle for T2 cascade: comparison between original SA model, SA-ANN model and experimental results

Table 1 :

 1 Mass averaged kinetic losses: convergence with grid size and reconstruction order

	11480 34440 2.39E-002
	Mesh B, k=1 21195 63585 2.27E-002
	Mesh C, k=1 40436 121308 2.24E-002
	Mesh A, k=2 11480 68880 2.25E-002
	Mesh B, k=2 21195 127170 2.24E-002
	Mesh C, k=2 40436 242616 2.24E-002

Table 2 :

 2 Regression coefficient R for several architectures of the ANN

	n HL = 1 0.799	0.831	0.865	0.897
	n HL = 2 0.822	0.890	0.918	0.953

the regression coefficient R for different ANN architectures are reported. n N = 5 n N = 10 n N = 20 n N = 40

Acknowledgements

Computational resources were provided by HPC@POLITO, a project of Academic Computing within the Department of Control and Computer Engineering at the Politecnico di Torino (http://www.hpc.polito.it).

We aknowledge the CINECA award under the ISCRA initiative, for the availability of high performance computing resources and support. The simulations were performed on the Marconi Tier-0 System for the Project FIRST.