Matthieu Leyko

Ignacio Duran

Stphane Moreau

Franck Nicoud

Thierry Poinsot

Simulation and modelling of the waves transmission and generation in a stator blade row in a combustion-noise framework

Keywords:

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

A significant part of the noise emitted from aero-engines is generated by combustion. The relative importance of this source is expected to increase with the next generations of engines where both fan and jet noise are reduced with Preprint submitted to Journal of Sound and Vibration June 5, 2014 such technologies as chevrons, micro-jets or highly-swept and leaned composite blades. In the seventies, Candel [START_REF] Candel | Analytical studies of some acoustic problems of jet engines[END_REF],

Marble and Candel [START_REF] Marble | Acoustic disturbances from gas nonuniformities convected through a nozzle[END_REF] and Cumpsty and Marble [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF] showed that combustion noise could be both direct or indirect.

Direct combustion noise comes from the noise emitted by the turbulent flame within the combustion chamber, and is then transmitted and reflected through the upstream (compressor) and downstream (turbine) turbo-machinery stages.

Indirect combustion noise comes from entropy fluctuations generated by the flame that produce noise when passing through the downstream turbine stages. This noise source is caused by the conversion of entropy spots into acoustic waves when they are accelerated in the blade rows. Therefore, indirect noise becomes relevant for large variations of Mach numbers. Candel [START_REF] Candel | Analytical studies of some acoustic problems of jet engines[END_REF] showed that this source could be important in aero-engines and could explain the observed excess noise. Cumpsty and Marble [START_REF] Cumpsty | Core noise from gas turbine exhausts[END_REF] also confirmed these results with an analytical model [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF]. More recently, Leyko et al. [START_REF] Leyko | Comparison of direct and indirect combustion noise mechanisms in a model combustor[END_REF] obtained scaling laws for a simplified one-dimensional (1-D) model combustor, which showed that indirect combustion noise could even be ten times larger than the direct one. However, in practice both sources must be taken into account. Indeed on the one hand, the noise generation caused by entropy fluctuations passing through the turbine stages, and the transmission of the noise generated within the combustion chamber have to be managed correctly to evaluate the noise at the engine outlet. On the other hand, the acoustic reflection by the compressor stages has to be considered to predict the impedance of the combustion chamber properly. In both cases, waves transmission and generation through several blade rows are involved, and the coupling within all stages must be considered. More specifically, the transfer functions of the turbo-machinery stages for the waves involved in the combustion process (acoustic, entropy and vortical) must be calculated.

An important number of analytical and semi-analytical models have been proposed for one-(1-D), two-(2-D) and three-dimensional (3-D) flows to deal with the propagation of acoustic and vortical waves through the turbo-machinery stages, ranging from simplest 1-D model dealing with compact nozzle to the most sophisticated ones involving 3-D blades. Tsien [START_REF] Tsien | The transfer functions of rocket nozzles[END_REF] and Crocco [START_REF] Crocco | Supercritical gaseaous discharge with high frequency oscillations[END_REF] first proposed an approach to calculate propagation of acoustic waves in 1-D nozzles. Muir [START_REF] Muir | The application of a semi-actuator disk model to sound transmission calculations in turbomachinery, part I: The single blade row[END_REF][START_REF] Muir | The application of a semi-actuator disk model to sound transmission calculations in turbomachinery, part II: Multiple blade rows[END_REF] then used the actuator-disk theory that replaces the actual 2-D configurations with an infinitely thin blade row, and classical conservation laws to establish 2-D relations between upstream and downstream flow.

Kaji and Okazaki [START_REF] Kaji | Propagation of sound waves through a blade row: I. analysis based on the semi-actuator disk theory[END_REF] first proposed a semi actuator-disk theory assuming that the spacing between the blades was infinitely thin but that the chord-length was finite. They then extended the model to a finite chord-length and a finite spacing between blades [START_REF] Kaji | Propagation of sound waves through a blade row: II. analysis based on the acceleration potential method[END_REF]. For both cases, the equations are solved numerically. Many other models are available for such 2-D problems. Recently Posson and Roger [START_REF] Posson | Parametric study of gust scattering and sound transmission through a blade row[END_REF] proposed a three-dimensional rectilinear cascade model both for generation and transmission losses based on Glegg's response [START_REF] Posson | Parametric study of gust scattering and sound transmission through a blade row[END_REF], which also accounted for cascade effects and finite chord, but neglected flow deviation (flat plates at zero angle of incidence). A generalization of this model to a more realistic annular cascade has also been presented [START_REF] Posson | Prediction of sound transmission through an annular cascade using an analytical cascade response function[END_REF].

The number of models dealing with entropy waves is smaller. In the seventies, Marble and Candel [START_REF] Marble | Acoustic disturbances from gas nonuniformities convected through a nozzle[END_REF] following the work of Tsien [START_REF] Tsien | The transfer functions of rocket nozzles[END_REF] and Crocco [START_REF] Crocco | Supercritical gaseaous discharge with high frequency oscillations[END_REF], developed an approach taking into account acoustic and entropy waves for 1-D compact nozzles. At the same time, Ffwocs-Williams and Howe [START_REF] Williams | The generation of sound by density inhomogeneities in low mach number nozzle flows[END_REF] developed a Green's function to propagate entropy waves in non-compact 1-D nozzles, but with a low-Mach number assumption. Pickett [START_REF] Pickett | Core engine noise due to temperature fluctuations convecting through turbine blade rows[END_REF] established relationships based on mass and momentum conservation to propagate entropy waves in a 2-D compact blade row.

Extending the approach of Marble and Candel, Cumpsty and Marble also proposed a model based on the linearized conservation equations of mass and total enthalpy, to deal with waves transmission and generation (involving entropy waves) through an axially compact blade row in a 2-D periodic and uniform inviscid flow. Finally Stow and Dowling [START_REF] Stow | Reflection of circumferential modes in a choked nozzle[END_REF] studied circumferential modes in choked nozzles numerically.

Bake et al. [START_REF] Bake | Experimental investigation of the entropy noise mechanism in aero-engines[END_REF] have investigated the indirect noise for a 1-D nozzle flow experimentally, in the so-called Entropy Wave Generator (EWG). Mühlbauer et al. [START_REF] Mühlbauer | Numerical investigation of the fundamental mechanism for entropy noise generation in aero-engines[END_REF], Leyko et al. [START_REF] Leyko | Numerical and analytical investigation of the indirect combustion noise in a nozzle[END_REF] and Duran et al. [START_REF] Duran | Analytical and numerical study of combustion noise through a subsonic nozzle[END_REF] have simulated this configuration for several flow conditions (subsonic and supersonic choked nozzle), and have reproduced the measured indirect and direct noise with acceptable precision. Leyko et al. [START_REF] Leyko | Numerical and analytical investigation of the indirect combustion noise in a nozzle[END_REF] have evidenced the critical effect of the outlet reflections in the supersonic choked case. They have also shown that the 1-D theory of Marble and Candel [START_REF] Marble | Acoustic disturbances from gas nonuniformities convected through a nozzle[END_REF] works very well for this nozzle and flow configuration, and it could help designing future experimental set-ups to limit the effect of downstream acoustic impedance. Duran et al. [START_REF] Duran | Analytical and numerical study of combustion noise through a subsonic nozzle[END_REF] have shown similar results for the low subsonic cases provided both the direct noise from the entropy generator and the inlet reflections are also considered. Beyond a throat Mach number of 0.3, the nozzle becomes non-compact and this new effect should be accounted for to reproduce the measurements. However, no experimental or numerical studies have been proposed yet, to investigate indirect noise for a 2-D flow within a blade row. Experimental measurements of such phenomenon remain nowadays difficult on real engines because of both the compactness of the turbo-machinery stages, and of the high operating pressures and temperatures. What is proposed in the present study instead is to look at the theoretical aspects of Cumpsty and Marble's extension [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF] of Marble and

Candel's model [START_REF] Marble | Acoustic disturbances from gas nonuniformities convected through a nozzle[END_REF], and to assess its validity by comparing with dedicated numerical results. The simplest passage of entropic waves and acoustic waves through an isolated stator vane is thus studied both numerically and analytically, extending the work of Mishra and Bodony [START_REF] Mishra | Evaluation of actuator disk theory for predicting indirect combustion noise[END_REF]. The theoretical aspects of the analytical model are first outlined in Sec. 2, whereas the propagation of entropy, vorticity (because of the two-dimensional flow) and acoustic waves through an isolated blade is addressed with a Large-Eddy-Simulation (LES) code in Sec. 3. Finally, the results are presented in Sec. 4 and conclusions about the validity of the model are drawn in Sec. 5.

Theory

This section provides a detailed description of the model of Cumpsty and Marble [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF] for the propagation of waves through a stator blade, that has been subsequently used in a complete wave-propagation model within a turboengine [START_REF] Duran | Computing combustion noise by combining Large Eddy Simulations with analytical models for the propagation of waves through turbine blades[END_REF]. This model is based on the same principles as the compact nozzle of Marble and Candel [START_REF] Marble | Acoustic disturbances from gas nonuniformities convected through a nozzle[END_REF], but it assumes a 2-D configuration to take into account the circumferential component of the turbomachinery and the flow deflection inducing vorticity fluctuations. The radial component of the flow is neglected, which is valid for large hub-to-tip ratio and perfect radial equilibrium in the turbomachine. The model assumes that the blade row is axially compact and that the blades spacing is small compared to the chord. The wave lengths of the different waves involved in the process are thus large compared with the axial length of the blade row. The steady flows on both sides of the blade row are assumed to be different, but uniform, so that the former is treated as a planar interface, as shown in Fig. 1. Assuming that the perturbations are small, balance equations can be written on both sides, the outgoing waves can be evaluated as a function of the incoming ones and of the main flow parameters. The primitive variables of the flow considered x y

w 1 θ 1 w 2 θ 2 k 1 ν 1 k 2 ν 2 Blade row
Figure 1: Description of the model approach. w is the velocity vector and k is the wave vector.

here are the pressure p, the density ρ, the relative velocity vector w and the flow angle θ. Similarly, the waves are characterized by the wave vector k and the wave angle ν.

As a first step, let's consider the propagation of the perturbations in the upstream and downstream regions where the flow is assumed inviscid and uniform. The linearized Euler equations (LEE) in a steady uniform flow are given by,

D Dt ρ ′ ρ + ∂u ′ ∂x + ∂v ′ ∂y = 0, (1)
D Dt (u ′ ) = - 1 ρ ∂p ′ ∂x , (2)
D Dt (v ′ ) = - 1 ρ ∂p ′ ∂y , (3)
D Dt s ′ c p = D Dt p ′ γp - ρ ′ ρ = 0, (4)
where the material derivative is written as

D Dt = ∂ ∂t + W • ∇ = ∂ ∂t + U ∂ ∂x + V ∂ ∂y , (5)
with U = W cos θ and V = W sin θ the axial and circumferential components of the mean velocity. The fluctuating variables u ′ and v ′ are the velocity perturbations in the x and y directions respectively. They can be related to w ′ and θ ′ through the velocity triangles

u ′ U = w ′ W -θ ′ tan θ, v ′ V = w ′ W + θ ′ tan θ . (6)
In each uniform region any wave φ is assumed harmonic and can be written in the form

w φ = A φ exp[i(ωt -k φ • x)], (7)
where A φ represents the amplitude, ω the angular frequency and k φ the wave vector associated with the wave, which can be expressed as a combination of the axial and circumferential components,

k φ • x = k x,φ x + k y,φ y, with k x,φ = k φ cos ν φ and k y,φ = k φ sin ν φ , ( 8
)
k φ is the modulus of the wave vector k φ . The latter is scaled with the mean sound speed c and the angular frequency, yielding the dimensionless wave vector, K φ = k φ c/ω. The material derivative given by Eq. ( 5) then becomes

D Dt = iω -Uik x,φ -Vik y,φ . (9)
For instance, the dimensionless entropy fluctuation s ′ /c p is defined as

w s ≡ s ′ c p = A s exp[i(ωt -k s • x)], (10)
where s ′ is the entropy fluctuation and c p the thermal capacity at constant pressure. Inserting Eq. (10) into Eq. ( 4) yields the following dispersion equation for the entropy wave

ω -W cos(θ)k s cos(ν s ) -W sin(θ)k s sin(ν s ) = 0. ( 11
)
If the dimensionless wave number K s = k s c/ω is introduced, the dispersion equation can be rearranged as follows

K s M cos(ν s -θ) = 1, ( 12
)
where M is the steady Mach number based on the steady velocity magnitude W and the steady speed of sound c. By definition, the entropy wave does not generate any pressure or velocity perturbation. Therefore, the fluctuations of primitive variables (s ′ /c p , w ′ /c, p ′ /γp and θ ′ ) generated by the entropy wave are

                                   s ′ s /c p w ′ s /c p ′ s /γp θ ′ s                                    =                                    1 0 0 0                                    w s . (13)
Similarly, the perturbation of vorticity ξ ′ ≡ ∂v ′ /∂x -∂u ′ /∂y is written as follows

w v ≡ ξ ′ ω = A ξ exp[i(ωt -k ξ • x)]. (14)
Taking the spatial derivatives of Eq. ( 2) and ( 3) leads to

∂ 2 u ′ ∂t∂y + U ∂ 2 u ′ ∂x∂y + V ∂ 2 u ′ ∂y∂y = - 1 ρ ∂ 2 p ′ ∂x∂y ∂ 2 v ′ ∂t∂x + U ∂ 2 v ′ ∂x∂x + V ∂ 2 v ′ ∂y∂x = - 1 ρ ∂ 2 p ′ ∂y∂x . (15)
Subtracting both relations in Eq. ( 15) yields

∂ ∂t ∂v ′ ∂x - ∂u ′ ∂y + U ∂ ∂x ∂v ′ ∂x - ∂u ′ ∂y + V ∂ ∂y ∂v ′ ∂x - ∂u ′ ∂y = 0. ( 16
)
Finally, the vorticity fluctuation verifies the following relationship

D Dt ξ ′ = 0. ( 17
)
By similarity with the entropy fluctuations, the dimensionless vorticity wave number K ξ = k ξ c/ω checks the following dispersion equation

K ξ M cos(ν ξ -θ) = 1. ( 18
)
Because vorticity waves are divergence free,

∂u ′ ξ ∂x + ∂v ′ ξ ∂y = 0, (19)
the following relationship is obtained

v ′ ξ = -u ′ ξ k x,ξ k y,ξ . (20)
By definition

ξ ′ = -i k x,ξ v ′ ξ + i k y,ξ u ′ ξ ( 21
)
that leads to

u ′ ξ c = -i ξ ′ ω sin ν ξ K ξ v ′ ξ c = +i ξ ′ ω cos ν ξ K ξ . ( 22
)
Combining Eq. ( 22) with Eq. ( 6) for the fluctuations of the velocity components u ′ ξ and v ′ ξ yields the corresponding

fluctuations for w ′ ξ and θ ′ ξ θ ′ ξ = (-u ′ ξ sin θ + v ′ ξ cos θ)/W w ′ ξ = u ′ ξ cos θ + v ′ ξ sin θ. ( 23
)
The fluctuations related to the vorticity wave can finally be written

                                   s ′ ξ /c p w ′ ξ /c p ′ ξ /γp θ ′ ξ                                    =                                    0 -i sin(ν ξ -θ)/K ξ 0 i cos(ν ξ -θ)/(MK ξ )                                    w v . (24)
The pressure perturbation p ′ /γp is defined as follows

w ± ≡ p ′ ± γp = A ± exp[i(ωt -k ± • x)]. (25)
As entropy and vorticity fluctuations are described by independent fluctuations, the acoustic waves generate a fluctuating field that is irrotational and isentropic. Eqs.( 1)-( 4) can be combined to give the wave equation satisfied by the acoustic pressure field

D Dt 2 p ′ γp -c 2 ∂ 2 ∂x 2 + ∂ 2 ∂y 2 p ′ γp = 0, (26)
which leads to the dispersion equation verified by the dimensionless wave vectors

K ± = ck ± /ω (1 -K ± M cos(ν ± -θ)) 2 -K 2 ± = 0. ( 27
)
Indeed, two wave vectors k + and k -satisfy the wave equation [START_REF] Quartapelle | High-order Taylor-Galerkin methods for non-linear multidimensional problems[END_REF]. The wave numbers k y,± are real since the domain is periodic in the azimuthal direction. The wave numbers k x,+ and k x,-can be complex and read1

K x,± = -M cos θ(1 -K y,± M sin θ) 1 -M 2 cos 2 θ ±        (1 -K y,± M sin θ) 2 -(1 -M 2 cos 2 θ)K 2 y,± (1 -M 2 cos 2 θ) 2        1 2 . ( 28
)
By expressing u ′ ± and v ′ ± from the momentum Eqs. ( 2) and ( 3), and noticing that ρ ′ ρ = p ′ γp when considering acoustics only, the fluctuations of primitive variables can be related to the acoustic waves

w ±                                    s ′ ± /c p w ′ ± /c p ′ ± /γp θ ′ ±                                    =                                    0 K ± cos (ν ± -θ) / [1 -K ± M cos (ν ± -θ)] 1 K ± sin (ν ± -θ) / [M (1 -K ± M cos (ν ± -θ))]                                    w ± . ( 29
)
The contribution of all waves (s, v and ±) can be summed to yield the total fluctuations of the primitive variables, since the so-called waves are independent to first-order. Equations ( 13), ( 24) and ( 29) are then used to build the transformation matrix [P p w ], which relates the fluctuations of the primitive variables to these waves:

                                   s ′ /c p w ′ /c p ′ /γp θ ′                                    = [P p w ] •                                    w s w v w + w -                                    (30)
with:

[P p w ] =                                    1 0 0 0 0 -i sin(ν ξ -θ) K ξ K + cos(ν + -θ) 1-MK + cos(ν + -θ) K -cos(ν --θ) 1-MK -cos(ν --θ) 0 0 1 1 0 i cos(ν ξ -θ) MK ξ K + sin(ν + -θ) M[1-MK + cos(ν + -θ)] K -sin(ν --θ) M[1-MK -cos(ν --θ)]                                    (31)
As a second step, the blade row itself is assumed compact. Therefore the fluctuations of entropy, mass-flow and energy are the same on both sides. For small perturbations, the conservation of the entropy between the inlet and the outlet of the blade row leads to:

s ′ c p 1 = s ′ c p 2 (32)
where [START_REF] Bake | Experimental investigation of the entropy noise mechanism in aero-engines[END_REF] relates to the inlet of the blade row and (2) to its outlet, as shown in Fig. 1. For small perturbations, the mass-flow rate conservation can be expressed as

p ′ γp 1 + 1 M 1 w ′ c 1 -θ ′ 1 tan θ 1 - s ′ c p 1 = p ′ γp 2 + 1 M 1 w ′ c 2 -θ ′ 2 tan θ 2 - s ′ c p 2 . ( 33
)
The total temperature T t (representing the total fluid energy) is based here on the velocity magnitude W. Therefore, considering that T t is conserved through a stator vane row yields for small perturbations

1 1 + 1 2 (γ -1)M 2 1 p ′ γp 1 + M 1 w ′ c 1 + 1 γ -1 s ′ c p 1 = . . . 1 1 + 1 2 (γ -1)M 2 2 p ′ γp 2 + M 1 w ′ c 2 + 1 γ -1 s ′ c p 2 . ( 34
)
The axial flow upstream and downstream of the blade row is here considered to be always subsonic (like in aeroengines), therefore there are always three waves traveling downward (acoustic, entropy and vorticity) and one traveling upstream (acoustic). Four waves are entering the blade row (three at the inlet and one at the outlet) and four are also outgoing (one at the inlet and three at the outlet). For the flow deviation, the Kutta condition proposed by Cumpsty and Marble [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF] is used

θ ′ 2 = 0 (35)
or in a more general form, as proposed by Cumpsty and Marble [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF],

θ ′ 2 = αθ ′ 1 ( 36
)
where α is a constant real number defined experimentally. Finally, Eqs. ( 32), ( 33), ( 34) and ( 36) can be written with the matrix [E e p ] that relates the primitive variables fluctuations at the inlet and the outlet of the blade row:

[E e p ] 1 •                                    s ′ /c p w ′ /c p ′ /γp θ ′                                    1 = [E e p ] 2 •                                    s ′ /c p w ′ /c p ′ /γp θ ′                                    2 (37)
with:

[E e p ] 1 =                                    1 0 0 0 -1 1/M 1 1 -tan θ 1 µ 1 /(γ -1) M 1 µ 1 µ 1 0 0 0 0 α                                    (38)
and:

[E e p ] 2 =                                    1 0 0 0 -1 1/M 2 1 -tan θ 2 µ 2 /(γ -1) M 2 µ 2 µ 2 0 0 0 0 1                                    (39)
where µ = (1

+ (γ -1)M 2 /2) -1 .
The overall problem of the stator row can then be expressed in terms of waves in a matrix-form as follows

[E e p ] 1 • [P p w ] 1 [B] 1 •                                    w s w v w + w -                                    1 = [E e p ] 2 • [P p w ] 2 [B] 2 •                                    w s w v w + w -                                    2 (40)
In Eq. ( 40), the left-hand-side term contains the wave w - 1 that is travelling upward from the blade row and is unknown, whereas the other waves are imposed. The situation is reversed in the right-hand-side term, where the wave w - 2 is

an input, and the other ones have to be found. In the case of an isolated blade row, the terms of the last column of the matrix [B] can be simply permuted from left-to-right and changed of sign yielding the matrices [A] i and [A] o .

Likewise the last terms of the vectors containing the waves are also permuted. Finally, the equation to be solved is

                                   (w s ) 2 (w v ) 2 (w + ) 2 (w -) 1                                    = [A] -1 o • [A] i                                    (w s ) 1 (w v ) 1 (w + ) 1 (w -) 2                                    (41)
The general expression of the outgoing waves is rather cumbersome. However, the present simulations shown in the next section are achieved for waves parallel to the blade row (k y = 0), so that the expression of the outgoing waves in Eq. ( 41) can be simplified. Four different cases can then be defined in order to describe the problem completely and are given in Appendix A.

Numerical simulation

The goal of the numerical simulation is to calculate the actual transfer function of a turbine stator vane impinged by different perturbations, in order to evaluate the error made when using the simplified compact model of Cumpsty and Marble [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF][START_REF] Cumpsty | Core noise from gas turbine exhausts[END_REF]. However, computing the turbulent flow field in the complete annular blade rows of an actual turbine even at a mean radius (2-D), is still out of reach. Indeed validating the model involves a very large range of characteristic times (ratio in the order of 10 7 in this case). On the one hand, wavelengths of the imposed perturbations have to be sufficiently large to be able to verify the compact blade row assumption. On the other hand, the time-step of the computation has to be sufficiently small to solve the flow correctly for the small geometrical details, such as the trailing-edge for instance. The simulations are therefore done for a single 2-D stator passage as the computational domain, with periodic boundary conditions on each side of the vane to mimic an infinite blade row. As explained by Cumpsty and Marble [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF], the axial Mach number through the blade rows is subsonic in most turbine applications.

However, this is not the case for the global Mach number which is supersonic at the outlet of the stator vanes due to the strong deviation of the mean flow. The presence of a supersonic flow, and the shock waves downstream, adds another layer of complexity for both the numerical simulation of the propagation of waves and the analytical models.

The propagation of entropy and acoustic waves through shock waves has been studied analytically [START_REF] Leyko | Numerical and analytical investigation of the indirect combustion noise in a nozzle[END_REF][START_REF] Goh | Phase prediction of the response of choked nozzles to entropy and acoustic disturbances[END_REF][START_REF] Stow | Reflection of circumferential modes in a choked nozzle[END_REF][START_REF] Duran | Solution of the quasi one-dimensional linearized Euler equations using flow invariants and the Magnus expansion[END_REF] and numerically [START_REF] Leyko | Numerical and analytical investigation of the indirect combustion noise in a nozzle[END_REF] in 1-D configurations. The stator blade considered in this section is subsonic in order to compare the results with Cumpsty and Marbel's [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF] method, which does not include the acoustic and entropy waves generated caused by the interaction of perturbations with the shock wave. Moreover, the flow at the outlet of the combustion chamber is typically turbulent. This turbulence would, if introduced in the numerical simulation, generate its own noise due to the interaction with the stator blade. As this source of noise is out of scope in the present work, no turbulence has been imposed at the inlet boundary. The 2-D flow is simulated using the unsteady compressible LES solver AVBP [START_REF] Moureau | Numerical methods for unsteady compressible multicomponent reacting flows on fixed and moving grids[END_REF]. These Large-Eddy Simulations are performed using the Smagorinsky subgrid-scale model. The laminar and turbulent Prandtl numbers (Pr and Pr t respectively) are chosen large enough to minimize thermal diffusion within the flow (Pr t = 100), and thus to avoid problems in the interpretation of the results. The numerical method used for the simulations comes from the family of Two-step Taylor-Galerkin (TTG) proposed by Quartapelle and

Selmin [START_REF] Quartapelle | High-order Taylor-Galerkin methods for non-linear multidimensional problems[END_REF]. The selected TTG scheme is TTG4A that is third-order in space and fourth-order in time, and exhibits good dispersion properties. The topology of the computational case is sketched in Fig. 2. The mesh contains 115000 triangles and the origin of the reference frame is located at the blade trailing edge. The transfer functions are evaluated

L out x = 6L y L in x = 3L y L y ℓ x x y Inlet (1)
Outlet (2) by only pulsating one kind of wave (entropy, vorticity or acoustic) at a time at the corresponding computational domain boundary (inlet or outlet), whereas all the other in-going waves are set to zero. Only planar waves with a purely axial wave vector (k y = 0) are considered. The pulsed waves are imposed in the following manner

f (t) = N k k=1 A k sin(k 2π t/τ 0 ), (42)
where N k is the number of frequencies which are pulsed (50 in the present case) and where τ 0 is a fundamental characteristic time based on a Strouhal number of 0.1 (the former is based on the corresponding combustion chamber swirler properties). The amplitude A k is the same for each frequency k, and it is in the order of 10 -3 . The temporal evolution of the signal f (t) is shown in Fig. 3. The spectral decomposition of this signal is of course made up of N k harmonic frequencies based of 1/τ 0 and of the same magnitude. As shown in Tab. 1, four kinds of simulations are performed, which starts from a statistically converged unperturbed flow. The simulations with the perturbations of entropy (S-1) and acoustic (A-1) waves at the inlet, and with the perturbation of the acoustic wave at the outlet (A-2), provide the acoustic-to-acoustic and entropy-to-acoustic transfer functions of the blade row, whereas a fifth simulation N-0 evaluates the global noise level, to verify that it is lower than the one produced by the perturbation of the inlet waves. In the third simulation (V-1), the vorticity-to-acoustic transfer function is calculated imposing the vorticity perturbation at the inlet. For simulations S-1 and A-1, the frequencies of the imposed signal are representative of what can be found at the outlet of a real combustion chamber. The same frequencies are used for the V-1 and the A-2 Subsonic case: M 1 = 0.12, θ 1 = 0 and M 2 = 0.66, simulations. For each case a transient of 6 τ 0 is discarded and 66 τ 0 are used to collect the statistics of the simulations.

θ

Results

As shown in Fig. To understand the limitations of the 2-D analytical model, the numerical simulations are used to study the com-pact assumption and its influence on the transfer functions. First, only the entropy wave is imposed with the f (t) signal at the inlet in the simulation S-1, and the other incoming waves are left untouched. Fig. 5 shows the field of the entropy perturbation s ′ /c p at different times. It can be seen that the entropy wave remains planar before the blade row (there is no injection of vortical perturbations) and that it is rather perturbed after the interaction with the blade and the downstream wake. As the flow is turbulent after the blade row (Re ≈ 10 6 ), all the post-processing of waves is performed by integration of the required quantities along the transversal direction, since the waves of interest are axially oriented (k y = 0). Then, quadratic average of the waves Fourier-transform modulus is calculated along the axial direction depicted in Fig. 2. This implies that results will present only the absolute value of the waves, as the phase is lost due to the quadratic averaging. The local 1-D waves w j x are calculated as follows:

w j x (x, k) = 1 nτ 0 nτ 0 0 1 L y L y 0 w j (x, y, t)dy exp(i2πkt/τ 0 )dt (43)
whereas the global 1-D waves w j are obtained using the following relation:

w j (k) = 1 L x L x 0 w i x (x, k) 2 dx ( 44
)
In Fig. 6 the acoustic response of the blade row to an entropy perturbation is plotted as a function of the dimensionless convective wavelength λ c /ℓ x based on the inlet quantities. The transfer function obtained from the simulation S-1 tends to the compact model when the wavelengths are large, and decreases rapidly for small wavelengths, where the assumptions are no longer valid. Figure 7 shows the acoustic response of the blade row to an acoustic perturbation plotted as a function of the dimensionless acoustic wavelength λ a /ℓ x based on the inlet quantities. The frequencies in simulation A-1 are the same as in S-1, so that the wavelengths are larger than in case S-1. The transfer function obtained from simulation A-1 is thus very close to the model value since the wavelengths are large compared with the blade axial length ℓ x . As seen in Section 2, the model of Cumpsty and Marble [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF] assumes that the entropy is conserved.

It is obvious that the turbulent mixing in the flow downstream of the blade row strongly affects the structure of the wave and that the initial planar coherence is lost for higher frequencies as it is shown in Fig. 8. The global entropy wave in Fig. 8-a, which corresponds to the axially averaged wave, is not representative in this case, since the mixing process takes place after the blade row. The evolution along the axial direction of the entropy wave is more relevant of its attenuation as observed in Fig. 8-b. Large wavelengths are practically not affected by the mixing that occurs at the scales of the blade dimension, and levels remain the same until the outlet. The intensity of small wavelengths decreases more rapidly, and noticeably these waves are already attenuated at the very beginning of the downstream flow (x/l x = 0). Since the spacing between the stator and the rotor ∆x is generally very small (∆x/l x ∼ 0.1) and tends to decrease, it is interesting to understand the reasons of such a phenomenon. The attenuation in the channel between blades can be explained by the distortion of the planar waves by the mean flow, and can be calculated from the mean flow characteristics. The distortion of the waves at the outlet of the blade row is here evaluated by a Lagrangian tracking of particles seeded at the inlet. Figure 9-a shows the path-lines of the particles computed from the mean flow (the flow is quasi-steady and not turbulent in the channel between blades), whereas This strongly supports the idea that the attenuation for small convective wavelength at the beginning of the downstream flow is caused by the distortion of initially planar waves by the mean flow. This phenomenon will of course produce non-planar waves at the outlet. A more general expression of Eq. ( 46) can be: observed without attenuation. When the wavelength is close to the dimension of the blade, times m, the transversal modes are maximum, whereas the planar mode strongly decreases. These observations suggest that, to study noncompact aspects of the waves transmission and generation, not only planar waves should be looked at but also spinning modes, even if the excitation is planar. This redistribution can have a significant influence when several blade rows are involved.

D m (k) = 1 L y L y 0 exp i2π
λ c /ℓ x (w c ) o [x = x 0 ]/(w c ) i m = 0
The model of Cumpsty and Marble [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF] also predicts the generation of vorticity waves. Indeed, both acoustic and entropy waves generate vorticity waves when propagated through the stator vane, as shown in Fig. 12. Even though there is no noise theoretically caused by vorticity in the stator vane case presently investigated (due to the θ 1 = 0) this could become relevant when coupling different blade rows. It is therefore interesting to study the vorticity wave generated by both acoustic and entropy waves propagating. In Fig. 12 the vorticity waves generated by acoustic perturbations are correctly predicted for large wavelengths, while those generated by entropy perturbations are underestimated. This is caused by the low level of vorticity induced in this case, which is of the same order of magnitude as the vorticity present in the reference simulation (N-0) and due to the trailing-edge wake, as shown in Fig. 13. The transfer functions for the vorticity wave (simulations V-1) are shown in Fig. 14 as a function of the non-dimensional convective wavelength. As seen in Appendix A, when the mean flow at the inlet of the blade row is straight (θ 1 = 0) vorticity waves generate no acoustic waves. The numerical simulations give a similar result for large wavelengths (λ c /ℓ x > 5).

Waves propagating downstream through the blade row generate both downstream and upstream propagating acoustic waves. When studying a complete turbine stage, the upstream propagating acoustic waves generated by the downstream blade rows have to be taken into account and correctly propagated. Simulation A-2 allows a comparison between the analytical model and the numerical simulations for the case of an upstream propagating acoustic wave through the blade row. Figure 15 shows the corresponding transfer functions, where it can be seen that at large wavelengths the simulation tends to the analytical value. Overall, the transfer functions for all waves obtained with the simulations are in good agreement with the compact analytical model at low frequencies. This is not the case when the wavelength is of the order of the blade length, or smaller. However perturbations within a aeronautical combustion chamber are expected to be low frequency, so that this results have to be considered in the combustion-noise framework, that is to say taking into account the spectra of entropy and acoustic waves at the combustion chamber outlet. To calculate the errors made on the noise prediction by the use of the model on a real configuration the spectra of entropy and downward propagating acoustic wave, obtained with LES on a model SAFRAN combustor, are used (see Fig. 16), where it can be seen that the entropy wave at the outlet of the combustion chamber is much larger than the acoustic wave. The approach to evaluate this spectrum is the same as the one used throughout the article (planar waves and k y = 0). The noise power (P) at the outlet is calculated The result obtained by the model is found to be only 17% different from the real case when comparing the noise power as presented in Eq. ( 49), or equivalently 0.6dB. This error rises up to 56% for the entropy (indirect noise), or a 2dB overestimate as reported in Tab. 2.

Conclusions

In the present study, the generalisation to a 2-D configuration made by Cumpsty and Marble [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF] is shown to have a major effect on the acoustic transmission of both direct and indirect combustion noise through blade rows compared

Figure 2 :

 2 Figure 2: Description of the computational case.

Figure 3 :

 3 Figure 3: Temporal evolution of the pertubation imposed on acoustic or entropy waves at the computational domain inlet.

s 1 (

 1 a) Indirect noise; M 1 = 0.12, M 2 = 0.66θ 1 = [0 → 90] Direct noise; M 1 = 0.12, M 2 = 0.66 θ 1 = [0 → 90]

Figure 4 :

 4 Figure 4: Acoustic response of the blade row at the outlet (2) to a perturbation at the inlet (1) for different flow directions. (a) Entropy perturbation; (b) Acoustic perturbation. Cumpsty and Marble (--), Marble and Candel with global Mach number (-• -), Marble and Candel with axial Mach number (---).

4 ,

 4 the generalization to a 2-D configuration (blade rows) made by Cumpsty and Marble has a major effect on the propagation of both indirect (left plot) and direct noise (right plot) compared with the 1-D model of Marble and Candel. As mentioned in Sec. 2, the 2-D model can deal with any angle of attack of the incident waves. The acoustic responses at the outlet of a blade row w + 2 to a purely axial entropy perturbation at the inlet w s 1 (Fig. 4 (a)), and to an acoustic one at the inlet w + 1 (Fig. 4 (b)) are then plotted for different inlet and outlet flow angles θ 1 and θ 2 , respectively. For the 1-D model of Marble and Candel, the Mach number based either on the total velocity magnitude or on the axial component of the velocity is used. It can be seen that the 2-D and 1-D models provide quite different results for important flow deviations that can be found in turbine stages, showing the importance of the 2-D analytical model. It should also be stressed that for θ 1 = θ 2 = 0, the results of the two theories collapse as expected.

3 Figure 5 :s 1 Figure 6 : 1 Figure 7 :Figure 8 :

 3516178 Figure 5: Instantaneous field of the entropy wave at different time steps for the case when entropy is pulsed. (a) t/τ 0 = 0.1, (b) t/τ 0 = 0.2, (c) t/τ 0 = 0.3.

Fig. 9 -Figure 9 :

 99 Figure 9: Mean flow characteristics. (a) Streamlines calculated from the mean flow; (b) Iso-temporal lines of massless particles seeded from the inlet and calculated from the mean flow.

Figure 9 -Figure 10 :

 910 Figure 9-b clearly reveals the deformation of the initially planar waves. The delay in the arrival of particles at the outlet t d (y) is plotted in Fig. 10-a versus the transversal direction y normalized with the pitch length L y (y/L y = 0 and y/L y = 1 correspond to the blades trailing-edges). The particles passing by the trailing-edges should theoretically have an infinite time of arrival, but have a finite value here because of the finite number of path-lines used to evaluate the

  k(t d (y)/τ 0 ) + m(y/l y ) dy(47) where the integer m corresponds to the m-th transversal mode of propagation of the entropy wave. This function is plotted in Fig.11for different transversal modes. When the frequency tends to zero, only planar waves (m = 0) are

Figure 11 :

 11 Figure 11: Attenuation of a convective wave (w c ) through the blade row calculated from the flow distortion. Mode m = 0 (--), positive modes (---) and negative (• • •).

1 Figure 12 :

 112 Figure 12: Vorticity response of the blade row to an entropy or acoustic perturbation at the inlet (cases S-1 and A-1). (a) Entropy; (b) Acoustic. Model (--), simulation (•).

λ c /ℓ x w v 2 Figure 13 :

 213 Figure 13: Vorticity at the outlet for the three differents cases. Case N-0 (•), case S-1 ( ), case A-1 ( ).

v 1 Figure 14 :

 114 Figure 14: Acoustic response of the blade row to a vorticity perturbation at the inlet (case V-1). (a) Upstream propagating waves; (b) Downstream propagating waves. Model (--), simulation (•).

2 Figure 15 :

 215 Figure 15: Acoustic response of the blade row to an acoustic perturbation at the outlet (case A-2). (a) Upstream propagating waves; (b) Downstream propagating waves. Model (--), simulation (•).

1 Figure 16 :

 116 Figure 16: Spectra of sources from LES. (a) Entropy source; (b) Acoustic source.

Table 1 :

 1 1 = 76 o Description of the different simulations.

	Name Pertubed wave Location	Duration
	S-1	Entropy	Inlet (1)	66 τ 0
	V-1	Vorticity	Inlet (1)	66 τ 0
	A-1	Acoustic	Inlet (1)	66 τ 0
	A-2	Acoustic	Outlet (2)	66 τ 0
	N-0	None	-	18 τ 0

A mistake is present in the expression of the acoustic wave vector in the article ofCumpsty and Marble 1977

with the 1-D model of Marble and Candel [START_REF] Marble | Acoustic disturbances from gas nonuniformities convected through a nozzle[END_REF], previously used to compare with the Entropy Wave Generator (EWG) experiment [START_REF] Bake | The entropy wave generator (ewg): A reference case on entropy noise[END_REF][START_REF] Leyko | Numerical and analytical modelling of entropy noise in a supersonic nozzle with a shock[END_REF]. The 2-D model can deal with any angle of incidence of the impinging waves, and even for incident axial ones the difference can be large due to the strong deviation of the mean flow present in typical turbine blades.

Moreover, the acoustic responses at the outlet of a blade row w + 2 to a purely axial entropy perturbation at the inlet w s 1 ,

and to an acoustic one at the inlet w + 1 are found to provide quite different results for important flow deviations as found in actual turbo-machines.

However, the 2-D analytical model still relies on the compact assumption that limits it to the low frequency range.

Numerical unsteady simulations of the response of an isolated turbine stator vane have therefore been performed in a 2-D configuration to evaluate its range of validity and estimate the error made when applying the method to non-zero frequencies. The acoustic transmission and reflection are shown to be well predicted by the model in the frequency range corresponding to acoustic wavelengths λ a /ℓ x > 10. The acoustic waves generated by entropy disturbances impinging on the blade row (responsible for indirect noise) are also well predicted at low frequencies, corresponding to convective wavelengths λ c /ℓ x > 10, but the theoretical results rapidly depart from the numerical ones for higher frequencies. Another assumption in Cumpsty and Marble's [START_REF] Cumpsty | The interaction of entropy fluctuations with turbine blade rows; a mechanism of turbojet engine noise[END_REF] model is that entropy waves are just convected through the turbine blades with no distortion. Through a proper post-processing of the steady flow the initially planar entropy waves were shown to be strongly distorted by the inter-blade non-uniform steady flow of the high-pressure distributor.

This occurs for entropy wavelengths smaller than the passage length, and therefore is significant at high frequencies.

The distortion of the entropy wave scatters the initially plane wave into higher-order modes due to the azimuthal dependency of the mean flow velocity (not considered in the analytical 2-D model), reducing the amplitude of the plane mode. The attenuation of the plane mode due to this effect has been calculated in order to predict the decrease of indirect combustion noise at large frequencies. This phenomenon leads to additional disagreement with the 2-D model and can be significant when considering the multiple blade rows of a turbine. However, entropy perturbations within an actual aeronautical combustion chamber are expected to be mostly in the low-frequency range and this effect could be small. For this reason the derived transfer functions should be considered taking into account the acoustic and entropy waves as inputs coming from measurements or LES results at the exit of a combustion chamber. The calculation of the errors made on the noise prediction at the outlet of an actual distributor vane shows that an error of 0.6 dB is done on the direct noise, and 2.0 dB on the entropy one, when using the present model with disturbances spectra from an actual LES. The entropy fluctuations obtained in the latter are in the low-frequency range, therefore the overall indirect noises calculated analytically and numerically are in the same range of magnitude. The results provided by this first-order compact approach are globally quite good and make it suitable as part of a core-noise evaluation tool.

In this work only plane waves have been considered. These waves represent a large part of the contribution to combustion noise, as high-order azimuthal modes are cut-off at low frequencies. To completely validate the 2-D analytical model and its assumptions, azimuthal modes should be considered in a future work, as well as computing the transmission of waves through several blade rows, including rotor blades.

Appendix A.

The analytical solutions of Eq. ((41)) are here provided for plane incoming waves. Considering k y = 0 and the Kutta condition given in Eq. ( 36), α = 0, Eq. ((40)) leads to following relationships:

where

For incoming wave w s 1 0 (other incoming waves are equal to zero):

For incoming wave w + 1 0 (other incoming waves are equal to zero):

(A.10)

For incoming wave w - 2 0 (other incoming waves are equal to zero):


Yu Zhong

Xudong Chen

Pingping Ding

Marc Lambert

Dominique Lesselier

Fast Calculation of Response of Scatterers in Uniaxial Laminates

The volume integral equation method, usually implemented via the method of moments (MoM), is a popular approach to many scattering problems, like when inhomogeneous bodies are embedded within planarly layered media and made to interact with given sources, as henceforth considered.

When the layered media are isotropic, the construction of the impedance matrix of the MoM, the bottleneck of the method, can be fast since one is able to accelerate the calculation of the 1-D Sommerfeld integrals [START_REF] Alparslan | Closed-form Green's functions in planar layered media for all ranges and materials[END_REF]. This is due to the fact that the dispersion relations of the isotropic media enable to calculate the integral along the azimuth direction in closed form. But, if the uniaxial media are with optical axes parallel with the interfaces, e.g., [START_REF] Burke | Eddy-current induction in a uniaxially anisotropic plate[END_REF][START_REF] Wasselynck | Microscopic and macroscopic electromagnetic and thermal modeling of carbon fiber reinforced polymer composites[END_REF], or if the media are biaxially anisotropic, the above fast methods do not apply anymore. This results in possibly high computational cost to ensure accuracy when dealing with the 2-D inverse Fourier transform (IFT). However, we might still circumvent it, with the rectilinear mesh as done in traditional practice in volume integral equation methods, achieving efficient construction of the impedance matrix. In a pioneering contribution [START_REF] Catedra | Spectral domain analysis of conducting patches of arbitrary geometry in multilayer media using the CG-FFT method[END_REF], a method that uses the continuous Fourier transform and a windowing technique is proposed to analyse metal patches within an isotropic multi-layered structure, which can be used to avoid the conventional 2-D IFT. By the generalized Poisson summation formula and the windowing technique, the relation between the discrete Fourier spectrum and the continuous Fourier spectrum of the spatial response is re-derived, which, by using it, enables to efficiently and accurately calculate the response of the multilayer on a rectilinear mesh. When the continuous Fourier spectrum of the response of the laminates to the current basis function is available, which can be stably and efficiently calculated by the new recurrence relations proposed in our previous contribution [START_REF] Zhong | Electromagnetic response of anisotropic laminates to distributed sources[END_REF], we use such a fast method to construct the impedance matrix of the MoM involving uniaxial layered media whose optical axes lie parallel to their planar interfaces. Numerical tests confirm the efficacy of the proposed method.


Ferran Argelaguet

Adaptive Navigation for Virtual Environments

Keywords: I.3.7 [Computer Graphics]: Three-Dimensional Graphics and Realism-Virtual Reality H.5.2 [Information Interfaces and Presentation]: User Interfaces-Input devices and strategies, Interaction styles, User-centered design

Navigation speed for most navigation interfaces is still determined by rate-control devices (e.g. joystick). The interface designer is in charge of adjusting the range of optimal speeds according to the scale of the environment and the desired user experience. However, this approach is not valid for complex environments (e.g. multi-scale environments). Optimal speeds might vary for each section of the environment, leading to non-desired side effects such as collisions or simulator sickness. Thereby, we propose a speedadaptation algorithm based on the spatial relationship between the user and the environment and the user's perception of motion. The computed information is used to adjust the navigation speed in order to provide an optimal navigation speed and avoid collisions. Two main benefits of our approach is firstly, the ability to adapt the navigation speed in multi-scale environments and secondly, the capacity to provide a smooth navigation experience by decreasing the jerkiness of described trajectories. The evaluation showed that our approach provides comparable performance as existing navigation techniques but it significantly decreases the jerkiness of described trajectories.

INTRODUCTION

Navigation is a key task for any virtual reality application. The user should be able to modify their viewpoint in order to explore, search or maneuver in the virtual environment [START_REF] Kulik | Building on Realism and Magic for Designing 3D Interaction Techniques[END_REF]. However, although a wide range of navigation interfaces exist (e.g. locomotion interfaces, flying interfaces), few works have addressed how the navigation speed can be adjusted in order to ensure a smooth and pleasant navigation experience. In rate-control interfaces, the range of the input device is mapped into a fixed range of acceptable navigation speeds. This range might not be optimal for a wide range of virtual environments, resulting on motions that are potentially too fast (Which can induce motion sickness [START_REF] So | Effects of Navigation Speed on Motion Sickness Caused by an Immersive Virtual Environment[END_REF]) or too slow (decreasing user engagement).

In this work, we have explored how navigation speed can be automatically adjusted in order to provide an optimal navigation speed. Our approach monitors how navigation speed is perceived by the user (optical flow) and the relationship between the user and the virtual environment (time to collision). Using this information, we adjust the navigation speed according to the user and scene demands. The proposed approach can be easily integrated in any navigation interface and based on the results obtained from the conducted evaluation, it can (1) conveniently adapt the navigation speed according to the virtual content and (2) simplify the user control enabling experts and amateurs alike to efficiently navigate in virtual environments. * e-mail: fernando.argelaguet sanz@inria.fr

PREVIOUS WORK

Speed Control

The speed control of a navigation technique is linked with the scale of the environment and the user's preferences. While the scale of the environment provides a threshold for the maximum allowed speed, users, through the navigation interface, are able to adjust the speed using a number of input commands [START_REF] Jeong | User experimentation: an evaluation of velocity control techniques in immersive virtual environments[END_REF] and speed mappings [START_REF] Anthes | Navigation models for a flexible, multi-mode VR navigation framework[END_REF]. In terms of speed, we can distinguish the linear component of the motion (tangential speed) and the rotational component (rotational speed). If no speed adaptation is available, and if the scale of the environment is known a priori, then the maximum speed can be properly adjusted by the interface designer. In addition, if we are navigating through a real scale model, the tangential speed can be bounded to human walking speed (0.6m/s to 1.4m/s) [START_REF] Marchal | Joyman: A human-scale joystick for navigating in virtual worlds[END_REF]. Moreover, to better resemble the human walking speed, when rotational speed increases, tangential speed should be decreased [START_REF] Marchal | Joyman: A human-scale joystick for navigating in virtual worlds[END_REF].

However, due to the wide range of heterogeneity of virtual environments and differences in user preferences, only few papers have addressed the adaptation of navigation speed. The first approach considering dynamic speed adaptation was the "Context Sensitive Flying Interface" from Ware and Fleet [START_REF] Ware | Context sensitive flying interface[END_REF], who proposed the usage of depth information in order to adjust the navigation speed. In their implementation, the best results were achieved when the movement speed was coupled with the minimum value of the depth map. An improved version of Ware and Fleet interface is the Cubemap navigation proposed by McCrae et al. [START_REF] Mccrae | Multiscale 3D navigation[END_REF]. The Cubemap approach considers the entire surroundings of the user by computing a depth cubemap from the camera position. This information is used to adjust the speed (following a similar approach as [START_REF] Ware | Context sensitive flying interface[END_REF]) and avoid collisions. Following a different perspective, in previous work Argelaguet et al. [START_REF] Argelaguet | Automatic speed graph generation for predefined camera paths[END_REF], we proposed a perceptual-based approach based on visual attention models. Our approach was based on optimizing the navigation speed at each point of a predefined camera path according to the desired perceived motion. No interactivity was provided.

Optical Flow and Perceived Motion

Optical flow can be defined as the velocities of 3D surface points projected onto the imaging plane. Studies on optical flow have shown that it plays a key role in how camera motion is perceived by the audience (e.g. locomotion control [START_REF] Prokop | Visual influence on human locomotion. Modulation to changes in optic flow[END_REF] and motion sickness [START_REF] So | Effects of Navigation Speed on Motion Sickness Caused by an Immersive Virtual Environment[END_REF]).

Optical flow analysis has been explored for the estimation of the perceived speed when navigating in virtual environments: flight simulators [START_REF] Kolasinski | Simulator Sickness in Virtual Environments[END_REF], animations in outdoor environments [START_REF] So | A Metric to Quantify Virtual Scene Movement for the Study of Cybersickness: Definition, Implementation, and Verification[END_REF] and predefined camera paths along arbitrary scenes [START_REF] Argelaguet | Automatic speed graph generation for predefined camera paths[END_REF]. However, none of these approaches are able to provide a real-time solution which can handle arbitrary virtual environments. While the work of Konasinski [START_REF] Kolasinski | Simulator Sickness in Virtual Environments[END_REF] focused mainly on the optical flow accounted by ground motion (flight simulators), the work from So et al. [START_REF] So | A Metric to Quantify Virtual Scene Movement for the Study of Cybersickness: Definition, Implementation, and Verification[END_REF] considered that the environment always has a similar complexity. In contrast, although the approach proposed by [START_REF] Argelaguet | Automatic speed graph generation for predefined camera paths[END_REF] handles arbitrary environments, it is not suited for real-time navigation due to its complexity and the fact that it operates within predefined camera paths.

Leveraging in current visual attention models [START_REF] Argelaguet | Automatic speed graph generation for predefined camera paths[END_REF], our approach proposes a real-time adaptation of the navigation speed. While the previous approache [START_REF] Argelaguet | Automatic speed graph generation for predefined camera paths[END_REF] relied on predefined paths, we account for interactive navigations in which the user has total control of the camera movement and can also influence the camera speed. The following subsections describe the two components of the proposed approach: the quantification of the perceived speed and the adjustment of the navigation speed in order to ensure an optimal speed.

Estimation of the User's Perceived Speed

The estimation of the perceived speed is computed from (1) the state of the system (the relative distance between the 3D environment and the user, and the current navigation speed) and (2) the quantification of the optical flow between the current and the previous frame. The heuristic employed operates in image space and only takes into account the visible part of the scene. All the computations are done using vertex and fragment shaders and stored in additional render targets (256 by 256 pixels).

Time to Collision (TTC)

The TTC map is meant to gather information about the relative distance between the user and the virtual environment. For each pixel, we compute the relative speed and estimate the time when a potential collision might occur. The computation only takes into account minimum distance between each pixel and the viewing plane (z axis). Algorithm 1 describes the pseudo-code of the fragment shader used to compute the TTC map. Similar to [START_REF] Argelaguet | Automatic speed graph generation for predefined camera paths[END_REF], the time to collision is modulated through an exponential decay.

Algorithm 1 Time to Collision map computation.

Require: P t : 3D coordinates for each pixel (Observer's referenced) Require: P t-1 : 3D coordinates for each pixel in the last frame Require: ∆t: Time between frames Require: d: Weighting factor. We used 2.5s for the experiments.

for all pixels do speed =

|(||P t,i ||-||P t-1,i ||)| t ttc = P t,i (z)/speed
T TC i = e (-ttc 2 /2d 2 ) end for return TTC Optical Flow (OF) The OF map captures the perceived motion of the scene from the user's point of view. For each pixel, we compute the amount of displacement (in screen space) between two consecutive frames. Algorithm 2 describes the pseudo-code of the fragment shader used to compute the OF map. The final result is also modulated as proposed by [START_REF] Argelaguet | Automatic speed graph generation for predefined camera paths[END_REF].

Algorithm 2 Optical Flow map computation. Require: P t : 2D normalized coordinates for each pixel Require: P t-1 : 2D normalized coordinates for the previous frame Require: ∆t: Time between frames Require: pixelSize: Screen pixel size in mm (0.27mm) Require: mapSize: Map size in pixels (256px, 256px) Require: m: Weighting factor in mm/s (100mm/s) for all pixels do

distance i = |(P t,i -P t-1,i )| • mapSize speed i = ||distance|| • pixelSize/time OF i = 1 -e (-speed 2 i /2m 2 )
end for return OF Once both maps are computed, they are combined through a weighted sum (see Equation 1), although we explored different alternatives to combine the information, the weighted sum provided the best results in terms of stability and convergence. The combination of the TTC map and the inverse of the OF map accounts for objects close to the focus-of-expansion (FoE) which have a small time to collision and the optical flow accounts for the perceived camera motion.

h = n ∑ i=0 (T TC i * (1 -OF i ) + OF i n (1)
The h estimator provides the perceived user speed (no units) according to the current speed and the environment.

Speed adaptation

The speed adaptation is based on updating the current acceleration of the camera to ensure that optimal perceived speed is achieved without abrupt speed changes. The speed adaptation algorithm (see Algorithm 3) first computes the optimal speed according to the ratio (h opt /h) between the optimal perceived speed (h opt ) and the current perceived speed (h). Then it estimates the "acceptable" acceleration to reach the desired speed. The estimation depends on the (h opt /h) ratio, if the ratio is greater than 1 (too fast) the estimation ensures that a fast deceleration will ensure a fast convergence of h into h opt . On the contrary, if the ratio is smaller than 1 (too slow) a sigmoid function is used in order to rapidly increase the speed, but limiting the increase above a certain threshold.

s opt = s • (h opt /h) if h > h opt then ∆h max = -1200 + 20 • h else ∆h max = 7 + 13 • e -h 2 /(2 * 5 2 ) end if t req = |h opt -h|/∆h max a opt = (s opt -Speed)/t req a = a + 0.05 • (a opt -a) s = s + a * ∆t
However, the approach failed to provide smooth navigations for strong camera rotations. The contribution of camera rotation to the optical flow is mainly dependent on the rotation speed. Camera rotations dramatically increase the optical flow, resulting in excessive slowdowns when turning. We decided to remove the optical flow accounted for camera rotations from the computation of the OF map and provide a fixed mapping to decrease the tangential speed when turning. We applied the approach detailed by [START_REF] Marchal | Joyman: A human-scale joystick for navigating in virtual worlds[END_REF], which limits the tangential speed according to the rotational speed.

In our implementation, due to practical reasons, we used a joystick as the input device. Forward and backward motion was determined by the x-axis of the joystick while yaw rotation was controlled through the y-axis. Algorithm 4 describes how h is modulated according to the user input, where two scale factors are computed according to the tangential ( f t ) and rotational components ( f r ). The last step was the computation of the optimal perceived speed (h opt ). In several controlled environments, using the navigation technique proposed by [START_REF] Marchal | Joyman: A human-scale joystick for navigating in virtual worlds[END_REF], we computed h opt analytically. Several runs were performed and the mean h values computed were averaged to compute h opt . Our analysis determined that h opt = 65.

Collision Avoidance

Collision avoidance is a key element to ensure that the user is able to naturally navigate in the virtual environment while avoiding obstacles. In order to provide a collision avoidance mechanism, we followed a similar approach as McCrae et al. [START_REF] Mccrae | Multiscale 3D navigation[END_REF]. The information encoded in a depth map is used to determine potential collisions and avoid them. Depth information is used to compute a repulsion force altering the user trajectory and avoiding potential collisions.

The algorithm is dependent on the collision boundary and a softness parameter σ . In our implementation the collision boundary was computed according to the current tangential speed: δ = s/2, thus adapting the boundary to the speed and also to the level of scale. The softness parameter was also dependent on the the current tangential speed (see Equation 2).

w(d i ) = e (δ -d i ) 2 (0.25•δ ) 2 (2)
w(d i ) is computed for each pixel and d i represents the distance towards the observer. The repulsion force is computed for each pixel, averaged and applied to the direction of movement at each frame. If the repulsion force is greater than a certain threshold, we consider that the collision cannot be avoided and the navigation stops.

USER EVALUATION

We conducted a formal evaluation in order to explore the benefits and caveats of our approach compared to state of the art approaches. We evaluated our technique against a fixed control law described in [START_REF] Marchal | Joyman: A human-scale joystick for navigating in virtual worlds[END_REF] and an adaptive approach based only on depth information [START_REF] Mccrae | Multiscale 3D navigation[END_REF]. The aim of the experiment was to explore if the evaluated techniques allowed participants to navigate in a multi-scale environment and the jerk analysis of the described trajectories.

Procedure. The task was to traverse a virtual environment as fast as possible but limiting the number of collisions. The synthetic virtual environment build for this experiment considered two different scene configurations, first, a homogeneous zone shaped as a maze-like environment and non-homogeneous zone populated with geometrical objects (see Figure 1). Users could to control the forward-backward movement through the X-axis of the joystick and the yaw rotation was controlled through the Y-axis. For each control law, a short training session with no time limit was provided.

Apparatus. The experiments were conducted with a 24 wide screen monitor with a resolution of 1920 × 1200px, with Logitech Wingman Extreme 3 Joystick as solely input device. The distance between the user and the monitor was approximately 50cm.

Participants. Fifteen participants (thirteen male, two female) aged from 21 to 39 (M = 28.2; SD = 4.95), participated in the experiment.

Design and hypotheses

We used a factorial design in which the independent variables were (1) the speed control (Fixed, Depth, OpticalFlow), ( 2) the scale of the environment (Big, Medium, Small) and ( 3) the section of the virtual environment (Maze, Geometry). All variables were withinsubjects. We adjusted all three navigation techniques in order to ensure that all of them provided a similar behavior for the first level of scale (Big). In order to avoid learning effects, the order of the techniques was counterbalanced using a Latin square design. For the section and the scale factor, the order was fixed due to the design of the virtual environment (see Figure 1).

For each area (two environments times three levels of scale) we recorded the task completion time, the mean value of speed, and jerk [START_REF] Cirio | Kinematic evaluation of virtual walking trajectories[END_REF], and the mean value for the x-axis of the joystick (normalized values [-1..1]). At the end of each run, users were asked to provide their subjective impressions about the difficulty of the task (Control), the perceived speed, the smoothness of the trajectory in terms of speed (Jerk), and their preferences. Users had to rank each question using a 7-Likert scale.

Our hypotheses were: [H1] Higher task completion time as the level of scale decreases for the Fixed technique.

[H2] The level of scale will not influence task completion time for the Depth and OpticalFlow techniques. [H3] Different mean values of speed for each Section. [H4] Different mean values of jerk for each Section.

[H5] Lower mean values of jerk for the OpticalFlow technique.

[H6] Lower subjective ratings for the Fixed technique.

Analysis

The statistical analysis is based in ANOVA and Bonferroni pairwise tests for post-hoc analysis. Only pairwise comparisons with a significant difference of p < 0.05 are discussed.

Task Completion Time

To analyze the time required to traverse the environment, we do not consider each section independently (Maze, Geometry) but the time to traverse each level of scale. The two-way ANOVA Technique, Scale versus Time showed a strong interaction effect (F 4,116 = 35.09; p < 0.001; η 2 p = 0.54), main effects did not provide conclusive results. Post-hoc tests provided three main results. First, all pairwise tests did not show any significant differences for the Big level of scale between techniques. This result shows that the fine tuning of all techniques ensured a similar performance for the optimal level of scale. Second, pairwise tests showed significant differences for each level of scale for the Fixed technique (Fixed/Big: M = 157.8s, Fixed/Medium: M = 141.17s, Fixed/Small: M = 192.05s). Although we reject the null hypothesis, we observe that the task completion time does not increase monotonically. Participants completed the task faster for the Medium level of scale as they were still able to maneuver. However, for the Small level of scale the speed was excessive resulting in a lot of collisions with the environment, thus we cannot support [H1]. Third, pairwise tests did not show significant differences as the level of scale decreases neither for Depth nor for OpticalFlow. This result shows that the behavior of both techniques was not dependent on the level of scale, we cannot reject the null hypothesis but this result supports [H2].

Average Speed

Due to the results obtained for the task completion time, we did not consider the data from the Fixed technique for follow up analysis. The two-way ANOVA of Technique and Section versus Speed showed a main effect on Technique (F 1,14 = 29.96; p < 0.001; η 2 = 0.22) and Environment (F 1,14 = 78.75; p < 0.001; η 2 = 0.22). Posthoc tests showed that the mean speed values were higher for the Maze section (M = 1.09m/s; SD = 5.84m/s) than the Geometry section (M = 1.01m/s; SD = 8.0m/s), rejecting the null hypothesis and supporting [H3]. Post-hoc tests also showed that the mean speed was higher for the Depth (M = 1.09m/s; SD = 7.6m/s) than for the OpticalFlow technique (M = SD = 6.27m/s).

Average Jerk

The two-way ANOVA of Technique and Section versus Jerk showed a main effect on Technique (F 1,14 = 128.91; p < 0.001; η 2 = 0.34), Scale (F 1,14 = 722.20; p < 0.001; η 2 = 0.39). Post-hoc tests showed that the mean jerk was lower for the OpticalFlow technique (M = 44.5m/s 3 ; SD = 10.101m/s 3 ) than for the Depth technique (M = 76.6m/s 3 ; SD = 29.9m/s 3 ), rejecting the null hypothesis and supporting [H5]. Post-hoc tests also showed that the environment played a significant role (Maze: M = 76.6m/s 3 ; SD = 29.9m/s 3 , Geometry: M = 76.6m/s 3 ; SD = 29.9m/s 3 ), rejecting the null hypothesis and supporting [H4]. This result shows that the proposed technique reduces the amount of jerk compared to state of the art approaches. However, the environment still plays an important role.

User Behavior

The value of the x-axis of the joystick provides information regarding the effort required for the user to adapt the navigation speed. The two-way ANOVA of Technique and Section versus the mean value of the x-axis showed a strong interaction effect (F 1,28 = 902.68; p < 0.001; η 2 = 0.32). Pairwise tests did not find significant differences for the Depth and OpticalFlow techniques (for all conditions M ∼ 0.98). In contrast, for the Fixed approach significant differences appeared for all levels of Size (Big M = 0.96, Medium M = 0.61 and Small M = 0.1). These results show how participants had to adjust the speed manually.

An interesting finding about usability was the fact that by observing participants and their comments during the experiment we found that they felt more comfortable when they used the joystick as a binary input (full forward, full stop). The requirement of constantly adjusting the speed for the Fixed technique during the Medium and Small levels of scale was sometimes disturbing, especially due to the inability to efficiently control the navigation speed. Although this effect is strongly dependent on the input device used and the level of user expertise, it shows that for controlling forward speed, a binary/discreet input might be more efficient.

Subjective Ratings

Friedman rank tests showed that the ranking provided by participants was significantly different for Control (χ 2 (2) = 13.43; p < 0.001), Speed (χ 2 (2) = 13.63; p < 0.001), Jerk (χ 2 (2) = 8.10; p < 0.05) and Preference (χ 2 (2) = 8, 13; p < 0.05). Wilcoxon tests were used to follow up the analysis with Bonferroni adjustment (α < 0.0167). No significant differences were found between Depth and OpticalFlow techniques. On the contrary, the Fixed technique was the worst rated in terms of Control and Speed.

Subjective ratings were consistent with the information gathered during informal discussions. However, participants perceived differences between both adaptive techniques although they were not able to describe them. In addition, some users stated that for the Op-ticalFlow and Depth techniques, their feeling was that the camera movement was too slow, while others perceived it as reasonable.

CONCLUSIONS

In this paper we have proposed and evaluated a speed adaptation algorithm which automatically adjusts the navigation speed based on the spatial relationship between the user and the environment, and the optical flow. In comparison with existing approaches, we provide a perceptual approach which not only considers the virtual environment but also considers how the user perceives the camera motion. The evaluation of the proposed approach has shown two main results. First, it is able to provide comparable navigation performances as state of the art navigation approaches. Second, the speed adaptation approach reduces the jerkiness of the described trajectories. Decreasing the jerkiness of the trajectory ensures smoother camera motions and the described trajectories will better resemble real locomotion trajectories [START_REF] Cirio | Kinematic evaluation of virtual walking trajectories[END_REF].

However, we believe that real strength of our approach would be delivered in immersive virtual reality systems. Although not discussed in the paper, the computation of the optical flow is dependent on the field of view of the virtual camera. The amount of optical flow is highly dependent on the field of view, as we move away from the focus of expansion, the amount of optical flow increases. This is a key feature of our approach which is not handled by existing adaptive approaches. As future work, we plan on analyzing the behavior of the proposed approach in immersive systems (e.g. CAVE-like systems and HMDs). We will explore the role of the field of view on motion perception and analyze potential relationships between the perceived speed and simulator sickness.

Algorithm 3

 3 Tangential speed update. Require: s: Current tangential speed Require: a: Current acceleration Require: ∆t: Frame time Require: h opt : Target threshold for h h = computePercievedOpticalFlow()

Algorithm 4 7 h

 47 Update of the h according to the user input Require: joy: 2DoF Joystick data [-1..1] ω = joy(y) • 1 rad /s f t = || joy(x)|| f r = e -||ω|| /2 * 0.= h /ft * fr return h

Figure 1 :

 1 Figure 1: Virtual environment used for the navigation task. The environment had two different sections in order to explore the behavior of evaluated techniques in different scene configurations. Furthermore, three different levels of scale were used (1:1, 1:2 and 1:10).


Giacomo Rodeghiero

Ping-Ping Ding

Yu Zhong

Marc Lambert

Dominique Lesselier

MUSIC imaging

Keywords: simulations, signal processing, NDE applications

come

From eddy currents to test graphite-based materials to microwaves and beyond to test glass-based composite structures, one aims to obtain images of the possibly damaged parts with robust, fast inversion algorithms. In this contribution, such algorithms are tailored to detect small (compared to the local wavelength in propagative regime or skin depth in diffusive regime) inclusions affecting the structures mentioned above. These inclusions may be voids, fluid-filled cavities (i.e., isotropic) or even uniaxial ones.

Yet this requires proper models of the layerings to compute their response due to electromagnetic sources, notably electric dipoles or magnetic coils.

Based on [START_REF] Zhong | Electromagnetic response of anisotropic laminates to distributed sources[END_REF][START_REF] Rodeghiero | An efficient interpolation for calculation of the response of composite layered material and its implementation in MUSIC imaging[END_REF], it is proposed herein a method to compute in an effective fashion the dyadic Green's functions (DGF) for such structures within the framework of contrast-source integral equations. Special care is taken when the sources are far away from the origin, yielding an oscillating spectrum of the DGF. The Multiple Signal Classification (MUSIC) imaging method [START_REF] Ammari | MUSIC-type electromagnetic imaging of a collection of small three-dimensional inclusions[END_REF], which uses such DGF, is applied to find the position of small defects. Possible application of MUSIC for structure delaminations is treated also.


Achraf Kallel

Elie Hachem

Yves Demay

Jean-François F Agassant

Stability analysis of a polymer coating process

Keywords: Stability, coating, process, polymer, modeling, tracking, capturing PACS: 47.20.Gv

teaching and research institutions in France or abroad, or from public or private research centers.

INTRODUCTION

The considered coating process is used in several industrial processes such as coating and laminating. Polymer is driven in a flat die at high temperature and the extrudate is stretched at high draw ratio (above 200) on a very short distance (around 1 mm). The resulting thin polymer film is set down on a solid substrate to form a two-layer film (see Figure 1-a). A wavelike instability is observed through periodical variations in the thickness of the polymer film (draw resonance instability). It depends on processing conditions and it is obtained at draw ratios above a critical value. Numerical simulation of this process provides a better understanding of this instability. [START_REF] Barq | Stationary and dynamic analysis of film casting process[END_REF]. In the case of fiber spinning, when the stretching distance is increased, the critical Draw ratio Dr* increases, which is explained by non-isothermal effects [START_REF] Blyler | Melt spinning and draw resonance of a poly(-methylstyrene / silicone) block copolymer[END_REF][START_REF] Gupta | Draw resonance in the non-isothermal spinning of polypropylene[END_REF][START_REF] Demay | Instabilité d'étirage et bifurcation de Hopf[END_REF]. Silagy et al. investigated the stability of the cast of film process using different constitutive equations. A membrane model has been considered and both linear stability analysis and direct numerical simulation have been employed [START_REF] Silagy | Stationnary and stability analysis of the film casting process[END_REF]. In that case, increasing the drawing distance stabilizes the process. However, in our experimental studies, it was noticed that reducing the stretching distance to the order of the die gap considerably improves the stability of the process. This finding was proved by Souli et al. in the case of fiber spinning by using a 2D model and a tracking strategy combined with the linear stability method [START_REF] Souli | Finite element study of the draw resonance instability[END_REF].

In this paper, we study the flow of polymer melt between the extrusion die and the substrate during stretching (Figure 1). Since the stretching distance is of the order of die gap, thin shell assumption is no longer valid. A 2D model in the thickness of the polymer film is considered instead. Numerical two-dimensional and time-dependent models are developed to investigate the onset of draw resonance instability. A boundary-value problem with a free surface for the Stokes equations is considered and stability of the free surface is assessed using two different numerical approaches : the first one consists in a tracking strategy combined with a linear stability tool involving computation of leading eigenvalues, while the second approach is a direct simulation method involving a Level Set capturing strategy coupled with transient stability analysis. The novelty of this work is the ability to use a stabilized finite element anisotropic meshing coupled with a Level Set method to assess interface stability in coating processes.

GOVERNING EQUATIONS

Stokes Problem

Inertia and gravitational forces are neglected, isothermal problem with a Newtonian behavior is considered. Moreover, this model does not account for surface tension at fluid interfaces because it is considered negligible comparing to drawing forces [START_REF] Fortin | Numerical simulation of coextrusion and film casting[END_REF]. Therefore, the momentum equations are given by : ⃗ ⃗

where  is the Cauchy stress tensor defined by :

( ⃗ ⃗ ) (2)
and ⃗ ⃗ * + is the velocity vector, ( ⃗ ⃗ ) is the rate of deformation tensor,  is the dynamic viscosity of the fluid and p is the pressure inside the fluid. In addition, the fluid is assumed to be incompressible, so we have :

⃗ ⃗ ⃗ ⃗ (3)

The Boundary Conditions

The applicable boundary conditions are given by :  At the wall  1 : ⃗ ⃗ ⃗ (no-slip condition)

 On  2 , a Poiseuille flow velocity field is imposed  On  3 , we have symmetry conditions  On  4 , take-up velocity is imposed On the free surface H, two different boundary conditions are considered: the first condition assumes a zero stress at the free surface ⃗ , while the second condition is a non-miscibility condition and it is given by ⃗ ⃗ ⃗ , where ⃗ ,1) is a normal vector to the free surface. In general, the draw ratio is a dimensionless number that can be defined as follows:

∫ ∫ (4)

PROBLEM SOLVING

The interface-tracking and the interface-capturing techniques are both used in the computation of flow problems with moving boundaries and interfaces [START_REF] Tezduyar | Interface-tracking and interface-capturing techniques for finite element computation of moving boundaries and interfaces[END_REF]. The Stokes problem is discretized and solved using a classical stable mixed-formulation (MINI element). All the implementation details are given in the following section.

Tracking Strategy

An interface-tracking technique requires moving meshes that "track" the interfaces and are updated as the flow evolves [START_REF] Tezduyar | Interface-tracking and interface-capturing techniques for finite element computation of moving boundaries and interfaces[END_REF]. Therefore, we consider only the half flow domain of fluid Figure 1-b). For the resolution, we use a structured mesh and apply a transformation to refine it around the die exit singularity (see Figure 2).

The velocity and the pressure fields are then computed on this structured mesh which is restricted to the polymer flow. Interface with air is determined by successive iterations of Newton-Raphson's method (by adjusting the position of the nodes on the interface) to satisfy the kinematic interface equation given by : [START_REF] Silagy | Stationnary and stability analysis of the film casting process[END_REF] where ( ) is the velocity vector at interface and h is the interface position. This strategy is relatively precise but it is unable to describe transient evolution of the interface. The linear stability method is then used to predict onset of draw resonance. Since free surface is the only time-dependent variable, the steady state interface solution undergoes a perturbation of the form :

̅ ̂ ( 6
)
where ̅ is the steady-state interface solution and ̂ is a perturbation of the interface. Therefore, the kinematic interface equation is linearized near the steady state interface solution, which gives the following linearized equation:

( ̅ ) ̂ ̂

where . Thus, the stability problem is reduced to a classical eigenvalue problem. The stability of the steady-state solution is investigated through computing the leading eigenvalue.

FIGURE 1 .

 1 FIGURE 1. (a) Schematic of the coating process, (b) Half flow domain of fluidIn the case of cast film process, Barq et al. identified the occurrence of draw resonance instability by measuring thickness and width fluctuations of the film as a function of time[START_REF] Barq | Stationary and dynamic analysis of film casting process[END_REF]. In the case of fiber spinning, when the stretching distance is increased, the critical Draw ratio Dr* increases, which is explained by non-isothermal effects[START_REF] Blyler | Melt spinning and draw resonance of a poly(-methylstyrene / silicone) block copolymer[END_REF][START_REF] Gupta | Draw resonance in the non-isothermal spinning of polypropylene[END_REF][START_REF] Demay | Instabilité d'étirage et bifurcation de Hopf[END_REF]. Silagy et al. investigated the stability of the cast of film process using different constitutive equations. A membrane model has been considered and both linear stability analysis and direct numerical simulation have been employed[START_REF] Silagy | Stationnary and stability analysis of the film casting process[END_REF]. In that case, increasing the drawing distance stabilizes the process. However, in our experimental studies, it was noticed that reducing the stretching distance to the order of the die gap considerably improves the stability of the process. This finding was proved by Souli et al. in the case of fiber spinning by using a 2D model and a tracking strategy combined with the linear stability method[START_REF] Souli | Finite element study of the draw resonance instability[END_REF].In this paper, we study the flow of polymer melt between the extrusion die and the substrate during stretching (Figure1). Since the stretching distance is of the order of die gap, thin shell assumption is no longer valid. A 2D model in the thickness of the polymer film is considered instead. Numerical two-dimensional and time-dependent models are developed to investigate the onset of draw resonance instability. A boundary-value problem with a free

FIGURE 2 .

 2 FIGURE 2. Structured mesh transformation

ACKNOWLEDGMENTS

The authors gratefully acknowledge Bostik S.A. company for experimental and financial support.

Capturing Strategy

The interface between the air and liquid is now captured by solving a convected Level-Set function as introduced in [START_REF] Ville | Convected level set method for the numerical simulation of fluid buckling[END_REF][START_REF] Hachem | Stabilised Finite Element for high Reynolds number, LES and free surface flow problems[END_REF]. The basic idea of this method is to use both the physical time and the convective time derivative in the classical Hamilton-Jacobi reinitialisation equation. Consequently, a modified level-set function is given first as follows:

where  stands for the standard distance function, and E is the truncation thickness. The level-set evolution equation is then given by:

where  is a coupling constant depending on time discretisation and spatial discretisation, typically and is the convection velocity. Finally, the authors in [START_REF] Hachem | Stabilised Finite Element for high Reynolds number, LES and free surface flow problems[END_REF] show that by setting | | and √ , a rearranged form of (9) leads to the following simple convection equation:

The finite element formulation for the level set method is based on the use of the classical SUPG (Streamline upwind Petrov-Galerkin) method. It controls the spurious oscillations in the advection dominated regime. In brief, the finite element formulation of equation [START_REF] Hachem | Stabilised Finite Element for high Reynolds number, LES and free surface flow problems[END_REF] can be written as follows: find , such that,

where and are standard test and weight finite element spaces.

The classical Galerkin terms are represented by the first two integrals whereas the element-wise summation, tuned by the stabilization parameter , represents the SUPG term needed to control the convection in the streamline direction. Once the level-set function is defined all over the domain, it can be used to easily separate both polymer and surrounding air phases. At the interface, the sharp discontinuity of the viscosity is smoothed over a transition thickness using the following expressions:

where H is a smoothed Heaviside function. Here  is a small parameter such that , known as the interface thickness, and is the mesh size in the normal direction to the interface. Combined with anisotropic mesh adaptation, it leads to an accurate and precise multiphase flow framework (see [START_REF] Hachem | Stabilised Finite Element for high Reynolds number, LES and free surface flow problems[END_REF] for details).

NUMERICAL RESULTS AND DISCUSSION

In this section, two case studies are presented and analyzed. In the first case, we check the steady-state solutions obtained by both tracking and capturing techniques for a draw ratio of 18. In the second case study, we investigate the stability of the steady-state solution using a draw ratio of 10. The leading eigenvalue in each case are then computed and compared. In both cases, a shape factor is considered (Figure 1-b).

Steady-State Solution

Figure 3-a shows a comparison for the shape of the interface at the steady-state obtained by both strategies. We can clearly see that both methods are able to converge to the same solution. This is also confirmed by comparing the velocity profiles along the x-axis (figure 3-b). Using the capturing strategy, we solve a multiphase Stokes problem that accounts for both polymer and air. Figures 4-a

Stability of the free surface

In this section, the stability of steady state solution is investigated for a draw ratio of 10 and a shape factor . In the case of tracking strategy, stability analysis is conducted by computing the leading eigenvalue as explained above. The free surface is stable only if the real part of the leading eigenvalue is negative. At those conditions, a leading eigenvalue is obtained meaning that the process is stable. In the case of capturing strategy, the transient stability analysis is performed by introducing a small perturbation to the steady-state solution. The transient response due to this perturbation is followed by the direct simulation, taking into account both domains, the polymer fluid and the surrounding air. Therefore, we can capture the evolution of the final thickness of the film as depicted in figure 5-a. The real part of the leading eigenvalue is estimated from the exponential envelope of the transient response while its imaginary part is estimated from oscillation's time period. A leading eigenvalue of is obtained for the same conditions. This finding shows that both strategies were able to give the same stability results using two different stability analysis techniques. Figure 5-b shows that, at a draw ratio of 18, the leading eigenvalue decreases when the shape factor decreases. Thus, reducing the stretching distance improves the stability of the coating process.

CONCLUSIONS

The coating process was accurately simulated using tracking and capturing techniques and its stability was investigated. It was demonstrated that reducing the stretching distance improves the process stability. This corresponds to some experimental evidence. Future work will concern the development of models accounting for a pressure differential between the two sides of the polymer film and more realistic constitutive equations.