Petra Á Szilágyi

Stéphane Dorbes

Gábor Molnár

José Antonio Real

Zoltán Homonnay

Christophe Faulmann

Azzedine Bousseksou

Temperature and pressure effects on the spin state of ferric ions in the [Fe(sal 2 -trien)][Ni(dmit) 2 ] spin crossover complex

Keywords: A. Inorganic compounds, B. High pressure, C. Magnetic properties, D. Phase transitions

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Spin-crossover (SC) complexes of transition metal ions (3d 4 -3d 7 ) have attracted much attention from the fundamental point of view [START_REF] Gütlich | Spin Crossover in Transition Metal Compounds I-III[END_REF][2][3][4]. These complexes represent one of the best examples of molecular bistability. They may be switched between their low spin (LS) and high spin (HS) electronic configurations, leading to important changes in their magnetic, dielectric and optical properties. This switching may be induced by varying the temperature, the pressure, the magnetic field and by light irradiation [5].

Although for both ferrous and ferric ions the change in the total spin in octahedral SC complexes is ∆S=2, the actual change in metal-ligand bond-length (for a given donor-atom set) accompanying the spin-state change is less for the latter than for the former. This is the origin of the increased rate of inter-conversion of the 2 T 2 and 6 A 1 electronic configurations for iron(III), which has a consequence that the LIESST effect (light-induced excited spin-state trapping) is rarely observed for ferric compounds even at very low temperatures, in contrast to its ubiquitous occurrence for ferrous systems [6,7]. Some of the rare examples are reported in references [8,9]. For the same reason, Mössbauer spectra of ferric SC compounds display frequently relaxation lineshapes [2]. It is also often assumed in the literature that the smaller bond-length changes may explain that the occurrence of thermal hysteresis is rather rare for ferric complexes. Indeed, the great majority of them does not display first-order spin transition, but only a gradual spin conversion [7].

Electrical conductivity is a fundamental property of materials, which has been investigated in the domain of molecular science [START_REF] Cassoux | Chemistry of Advanced Materials: An Overview[END_REF]. Compounds of type (NR 4 )[M(dmit) 2 ] x (x>1), with nonintegral oxidation state, display interesting electron transport features. Combining SC and electrical conductivity was the aim of Dorbes et al. [START_REF] Dorbes | [END_REF] when they first synthesised the 1:1 complex [Fe(sal 2 -trien)][Ni(dmit) 2 ], precursor of the fractional oxidation state complex [Fe(sal 2 -trien)][Ni(dmit) 2 ] 3 [12]. Analogously, similar properties have been studied on the salts of the [Fe(qsal) 2 ] + and Fe(R-salEen)] + complex cations as well [13][14][15][16][17][18][19].

A c c e p t e d m a n u s c r i p t

2

Since the [Ni(dmit) 2 ] n-anion is redox active [START_REF] Cassoux | Inorganic Materials[END_REF] and the [Fe(sal 2 -trien)] + cation is a wellknown SC compound [START_REF] Tweedle | [END_REF] and displays a redox stability, combination of these ions appeared as a promising way to study the eventual coupling and synergy of these molecular features in one compound. The resultant 1:1 compound is built up of slabs of [Fe(sal 2 -trien)] + cations and layers of [Ni(dmit) 2 ] -anions [START_REF] Dorbes | [END_REF], which alternate along the c-axis and spread out onto the ab plane. In the HS state, the SC sites are constituted of severely distorted octahedral [FeN 4 O 2 ] centres, whereas in the LS state, the coordination octahedron is much more regular (Figure 1). This complex exhibits a first-order thermal phase transition with a large thermal hysteresis loop of the magnetic susceptibility, around 245 K (T ½ ↓=228 K and T ½ ↑=258 K) (Figure 2).

It should be noted here that polymorphism has been observed since the first preparations of the triclinic [Fe(sal 2 -trien)][Ni(dmit) 2 ] phase: a monoclinic phase, which remains HS whatever the temperature, is now systematically obtained. Up to now, attempts to synthesise again the triclinic phase, which appeared -by the way -stable for at least one year in ambient storage conditions, have failed despite the numerous syntheses. Such a disappearance of polymorphs is not a unique case, as reflected by the literature [22][23][24][START_REF] Mccrone | Polymorphism[END_REF].

The aim of the present work was to complete the first observations of Dorbes et al. [START_REF] Dorbes | [END_REF] by investigating in more details the effects of temperature and pressure on the spin state and associated properties of the triclinic [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex, which is one of the rare Fe III complexes displaying thermal spin transition with hysteresis. This first-order transition occurs without a change of the space group of the crystal (P ), creating thus an interesting scope for the comparison of experimental observations, especially pressure effects, with theoretical predictions [26].

Experimental

Mössbauer measurements were carried out in a flow-type liquid nitrogen cryostat using a conventional constant acceleration type Mössbauer spectrometer. γ-rays were provided by a 57 Co(Rh) source (15 mCi). The spectrum evaluations were done with the assumption of Lorentzian line shapes using the MossWinn3.0 software [27]. All isomer shifts are given relative to α-Fe at room temperature. Crystallographic unit cell parameters have been collected between 150 and 290 K in the heating and cooling modes by steps of 20 K using an Xcalibur (Oxford Instruments) X-ray diffractometer equipped with an Oxford Cryosystems cooler device (Cryojet). Unit cell parameters for every temperature were determined by using the CrysAlis CCD softwares [START_REF]version 1.171.24, Xcalibur CCD system[END_REF].

Magnetic susceptibility measurements at high pressures were carried out in a magnetic field of 2 T using a clamp-type, hardened beryllium bronze (CuBe) cell. Hydrostatic conditions are obtained by mixing the powder samples with a pressure transmitting mineral oil. Details of the pressure cell have been given in ref. [START_REF] Molnár | [END_REF].

Results and Discussion

Mössbauer spectra at selected temperatures are displayed in Figure 3. Spectrum parameters and errors for all spectra recorded are listed in Table 1. The hyperfine parameters of both doublets correspond well to the typical parameters of ferric HS (S=5/2, δ 250K =0.382(3) mm/s, ∆ 250K =0.565(4) mm/s) and LS (S=1/2, δ 80K =0.182(1) mm/s, ∆ 80K =2.742(3) mm/s) complexes. Especially the high quadrupole splitting for the LS form is very characteristic [START_REF] Greenwood | Mössbauer Spectroscopy[END_REF]. The fact that only one component appears in the spectra recorded at 80 K and 250 K proves the occurrence of a complete spin transition in the system in both ways. This finding is in agreement with the magnetic susceptibility data (χ M T 250K = 4.84 cm 3 Kmol -1 and χ M T 80K = 0.66 cm 3 Kmol -1 ), but only the Mössbauer spectra can prove the completeness of the transition unambiguously, especially in the presence of other paramagnetic species. The inter-conversion of the spin states in many instances is so rapid for ferric complexes that the separate contributions to the 57 Fe Mössbauer spectra are not resolved. This is not the case for the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex, where no line broadening due to relaxation effects is displayed in the spectra. (The small, temperature- independent asymmetry of the HS and the LS doublets is a typical signature of a texture effect in the system.) One should note also that the coexistence of both doublets (HS and LS) is clearly observed in the spectrum at 230 K, i.e. in the descending branch of the thermal hysteresis loop. The HS ferric cation has a t 2g 3 e g 2 ( 6 S) electronic ground state and a spherical symmetry. Any electric field gradient (EFG) could arise solely from charges external to the central atom. It is unlikely that there would be any low-lying excited levels, thus no significant change in the quadrupole splitting is expected. The situation, however, is somewhat more complex for LS Fe 3+ systems and usually considerable temperature dependence may be observed. The t 2g 5 ( 2 T) state corresponds to a single electron hole in an otherwise cubic triplet level. In complexes with lower symmetry than cubic, further degeneracy is removed, thereby producing a nonzero EFG and quadrupole splitting consequently. Thermal population of the low-lying excited states then causes a variation in the electron configuration with temperature, and since the distinct 3d orbitals have different contributions to the electric field gradients, there is usually a strong temperature dependence of the quadrupole splitting [START_REF] Greenwood | Mössbauer Spectroscopy[END_REF]. Therefore it is noteworthy that the quadrupole splitting values for the LS state display practically no temperature dependence between 230 K and 80 K. This phenomenon may be due to the electronic ground state being well isolated.

As the SC phenomena involve a transfer of electrons between the antibonding e g and nonbonding t 2g orbitals, the metal-to-ligand bond distances increase remarkably when going from the LS to the HS state. This results in a drastic change in the density of vibrational states, mainly in the metal-ligand skeletal vibrational modes. This leads, of course, to significantly higher vibrational entropy in the HS state [START_REF] Sorai | [END_REF]. The entropy change accompanying the SC can be represented as the sum of electronic, configurational and motional, i.e. vibrational, rotational and translational terms. The two latter can be neglected in solid state, so can be the configurational term in the absence of disorder of orientation. The electronic term contains spin and orbital contributions. As the orbital contribution may be neglected in low-symmetry complexes, the electronic term is given by the spin-only value,

A c c e p t e d m a n u s c r i p t

5 ∆S el =R•ln(S HS /S LS )=R•ln3=9.134 J/Kmol. In the case of ferrous SC compounds, many data are available for the variation of the entropy, which is typically about ∆S~30-60 J/Kmol, including the electronic contribution, ∆S el =R•ln5=13.4 J/Kmol. However, for ferric complexes only a few data are available. Nevertheless, one can expect the change of entropy accompanying the spin transition to be typically smaller than in the case of ferrous compounds, since the electronic contribution is smaller and, more significantly, the vibrational contribution should be smaller as well, since the bond lengths do not undergo such a drastic change as in the case of Fe 2+ . This was shown by Sorai et al. [32] when studying the [Fe(acpa) 2 ]PF 6 (∆S=36.2 J/Kmol) complex and by Conti et al. [33] for various [Fe(3EtO-salAPA)] 2+ complexes (∆S~30-40 J/Kmol). In order to determine the change of the vibrational entropy in the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex, DSC measurements have been carried out (Figure 4). The DSC curves exhibit singularities around the spin transition temperatures. Interestingly, two maxima are observed in the heating mode and a shoulder adjacent to the singularity representing the spin transition can also be seen in the cooling mode. From the area of the peaks, the total change of enthalpy and entropy induced by the spin transition can be determined as ∆H HL =7 kJ/mol and ∆S HL =∆H/T c =29 J/Kmol (and

T c ≈(T ½ ↓+T ½ ↑)/2=243 K).
Therein, the vibrational entropy variation is: ∆S vib =20 J/Kmol (∆S tot =R•ln3+∆S vib ), which is a rather low value, even comparing to other ferric SC systems [START_REF] Gütlich | Spin Crossover in Transition Metal Compounds I-III[END_REF]. It is well known that the spin crossover is accompanied by an important change in the volume of the unit cell. Furthermore, first-order spin transitions are often accompanied by a change of the space group as well [34]. Remarkably, in the case of [Fe(sal 2 -trien)][Ni(dmit) 2 ] both the HS and LS forms display the same triclinic P space group. The main structural modification due to the spin crossover concerns the Fe-N bonds. When the electrons go from the e g to the t 2g orbital, i.e. during the spin transition from HS to LS state, a strengthening of the metaldonor atom bonds occurs (decrease in the Fe-N bond distance), but the N-Fe-N angles are also strongly affected. The change of the unit cell volume is very important (typically 1-6% in Fe II N 6 complexes) hence relatively small external pressure (i.e. a few hundred bars) can affect the spin state of the system by stabilising the LS form. The determination of ∆V HL from the

A c c e p t e d m a n u s c r i p t

6 two structures recorded at 295 K (HS) and 180 K (LS) is not accurate enough since the thermal expansion of the lattice must also be taken into account. Therefore we have determined the unit cell parameters in a wide temperature range. The upper graph in Figure 5 shows the variation of the volume of the unit cell of the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex as a function of the temperature. From this figure, the coefficients of volumetric thermal expansion -defined as α V = 1/V 235K (dV/dT) -in the two spin states are α HS =192(±1) 10 -6 K -1 , α LS =163(±2) 10 -6 K -1 and the variation of the unit cell volume upon the spin transition is

∆V HL =15.4 Å 3 .
It is worth noting that the ∆V/V 295 K ratio is only about 0.9%, which is a rather small value. Beside the volume change, the anisotropy of the unit cell contraction is also a general feature for SC materials, since they often crystallise in low symmetry space groups [35]. As a first approximation, this anisotropy may be analysed as a function of the temperature dependence of the unit cell parameters within a wide temperature range that includes the spin crossover (Figure 5). It is remarkable that in the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex one can observe the The transition temperature (T c ) and the pressure are related to the variation of the volume of the unit cell. From the change of volume in the course of the spin transition, the effect of applied pressure on the thermal SC can be estimated, using the Clausius-Clapeyron equation from which the stabilisation of the LS state under pressure can be estimated as dT c /dp=∆V HL /∆S HL =16 K/kbar (since both ∆S HL and ∆V HL are positive values). Magnetic susceptibility measurements have been carried out in order to determine the effects of hydrostatic pressure on the [Fe(sal 2 -trien)][Ni(dmit) 2 ] SC complex. Curves recorded at various applied pressure values are displayed in Figure 6.

As expected, the hysteresis loop shifts toward higher temperatures with increasing applied pressure. The pressure shift of T c is approximately 16(±1) K/kbar, which is in excellent agreement with the value estimated above, but no comparison can be safely made because of a remarkable additional phenomenon. As shown in Figure 6, when cooling the system under 0.5 kbar applied pressure the transition stops at approximately halfway and the system remains blocked in this intermediate state instead of a complete transition into the LS state. Then, in the heating mode, the system first relaxes to the LS state before returning to the HS state. This phenomenon is most probably due to an additional structural change induced by the pressure. When raising the applied pressure (i.e. at higher temperatures), the intermediate state becomes less and less apparent. Nevertheless, when releasing the applied pressure in the cell, the pressure effects turn out to be irreversible to some extent.

It is interesting to note also that at 0.5 kbar the hysteresis width is wider (∆T p=500bar =61 K) when compared to the case of atmospheric pressure (∆T p=1atm =30 K), which is exactly the opposite of what is usually observed experimentally and expected from theoretical considerations. In fact, according to the Ising-like model developed in ref. [34] the hysteresis loop width is proportional to J/T c , where J is a phenomenological, mean-field parameter for the cooperativity and T c is the transition temperature. Since such low pressures are not expected to affect considerably the value of J, the hysteresis width will decrease with increasing T c values under pressure.

A few anomalous examples have already been reported in literature when the hysteresis width was found to increase under pressure [24,[37][38][39][40][41][42]. The most obvious explanation of this phenomenon can be a pressure induced additional structural change, since the sharp discontinuity in the phase diagram below 0.5 kbar, followed by a "regular" Clausius-Clapeyron behaviour above 0.5 kbar is very difficult to explain in any other ways. However, a few different explanations for this unusual phenomenon can also be mentioned, such as cracks of particles which is not expected in our case where the hydrostatic conditions were maintained during the experiments and which should result, in any case, in a more gradual transition and narrowing of the hysteresis [43]. The variation of thermodynamic parameters could also affect the hysteresis width. However, at 0.5 kbar pressure this effect cannot be regarded as significant [24,[44][45][46][47][48]. Let us also note that further raising the applied pressure from 0.5 to 5.7 kbar leads to a monotonous decrease of the hysteresis width (from 61 to 44 K), suggesting that no further structural changes occur above ca. 0.5 kbar. This is better seen in the (p,T) phase diagram (Figure 7), which reveals clearly a "classical" behaviour between 0.5

A c c e p t e d

m a n u s c r i p t 8 and 5.7 kbar with a higher Clausius-Clapeyron slope for the HS→LS transition (dT 1/2 ↓/dp=17(±0.5) K/kbar) when compared to the reverse LS→HS transition (dT 1/2 ↑/dp=14(±1) K/kbar). From this, one can estimate that the first-order spin transition is suppressed by the applied pressure above ca. 22 kbar (660 K), which is, of course, beyond the thermal decomposition temperature of the complex. The phase diagram confirms also the obvious discontinuous behaviour between 1.013 bar and 0.5 kbar, typical of structural phase transition.

Conclusions

Thermal and pressure effects on the [Fe(sal 2 -trien)][Ni(dmit) 2 ] spin crossover complex have been investigated. A complete thermal spin transition (S=1/2↔S=5/2) of the compound has been evidenced by Mössbauer spectroscopy in agreement with the previous magnetic susceptibility measurements. From DSC measurements, the variation of the enthalpy and entropy following the thermal spin transition has been determined as ∆H HL =7 kJ/mol and ∆S HL =29 J/Kmol. This latter value points to a relatively small change of the vibrational entropy upon the spin transition. X-ray diffraction data were collected at different temperature values and the discontinuous variation of the unit cell volume accompanying the first order spin transition was determined (∆V HL =15.4 Å 3 ), which corresponds to 0.9% volume change. Furthermore, it was found that the higher unit cell volume value in the HS state is solely due to the dilatation of the ac plane, whilst the unit cell gets contracted along the b direction. Magnetic susceptibility measurements under hydrostatic pressure revealed the expected positive shift of transition temperature (~16 K/kbar) in agreement with the calculated Clausius-Clapeyron slope (∆V HL /∆S HL ). However, this agreement must only be a coincidence as a probable structural change has been observed at low applied pressure and these pressure induced changes turned out to be partly irreversible. Further work using structural and spectroscopic methods under pressure will be necessary to elucidate the mechanism of pressure effects on the system. Let us finally note here that the investigation of light induced

A c c e p t e d m a n u s c r i p t 1

 1 Temperature and pressure effects on the spin state of ferric ions in the [Fe(sal-trien)][Ni(dmit)] spin crossover complex P.A. Petra Á. Szilágyi, Stéphane Dorbes, Gábor Molnár, J.A. José Antonio Real, Z. Zoltán Homonnay, C. Faulmann, Azzedine Bousseksou To cite this version: P.A. Petra Á. Szilágyi, Stéphane Dorbes, Gábor Molnár, J.A. José Antonio Real, Z. Zoltán Homonnay, et al.. Temperature and pressure effects on the spin state of ferric ions in the [Fe(sal-trien)][Ni(dmit)] spin crossover complex. Journal of Physics and Chemistry of Solids, 2010, 69 (11), pp.2681. 10.1016/j.jpcs.2008.06.106. hal-00612401 Temperature and Pressure Effects on the Spin State of Ferric Ions in the [Fe(sal 2 -trien)][Ni(dmit) 2 ] Spin Crossover Complex 1. Introduction

Figure 1 .

 1 Figure 1. Structure of the coordination octahedron of Fe III ions in the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex, in the HS ((S=5/2) left) and LS ((S=1/2) right) states [11].

Figure 2 .

 2 Figure 2. Magnetic susceptibility data for the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex displaying a thermal hysteresis loop [11].

A c c e p t e d m a n u s c r i p t 3 Differential

 3 Scanning Calorimetry (DSC) analysis was carried out in a NETZSCH DSC 204 instrument, under He purging gas (20 cm 3 /min) at a 10 K/min scanning rate, in both, heating and cooling modes. Temperature and enthalpy (∆H) were calibrated by using melting transition of standard materials (Hg, In, Sn).

Figure 3 .

 3 Figure 3. Mössbauer spectra of the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex at selected temperatures.

Figure 4 .

 4 Figure 4. DSC curves of the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex in the heating and cooling cycles.

3 Figure 5 .

 35 Figure 5. Variation of unit cell volume and the cell parameters in the [Fe(sal 2 -trien)][Ni(dmit) 2 ]complex as a function of temperature. Squares (dots) display the data obtained when the temperature was decreased (increased). Lines are inserted to guide the eye.

A c c e p t e d m a n u s c r i p t 7 dilatation

 7 only in the ac plane while the value of b decreases in the course of the LS→HS spin transition. It is also clear that this decrease of b is related only to the spin state change and is not a "negative thermal expansion" phenomenon. Such anisotropic change of volume may complicate pressure effects and may lead to additional structural changes under pressurewhen compared to the thermal transition.

Figure 6 .

 6 Figure 6. Magnetic susceptibility curves as a function of temperature for the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex at various applied pressure values.

Figure 7 .

 7 Figure 7. (p,T) phase diagram of the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex. Open (closed) symbols display the transition temperatures obtained in the heating (cooling) modes at different applied pressures. The two crosses show the transition temperatures after releasing the applied pressure. Lines represent linear fit on the data between 0.5 and 5.7 kbar (extrapolated to 0 kbar).

A c c e p t e d m a n u s c r i p t 10 effectsA c c e p t e d m a n u s c r i p t

 10 , inside and outside the hysteresis loop, can be an interesting further perspective of this work, even if the dark colour of the sample makes these experiments challenging.

Table 1 .

 1 Hyperfine parameters (and their statistical errors) for the [Fe(sal 2 -trien)][Ni(dmit) 2 ] complex at various temperatures and spin states. N/N tot , δ and ∆ signify the relative amount of the subspectra, the isomer shift and the quadrupole splitting, respectively.

	T (K)		HS N/N tot (%) δ (mm/s) ∆ (mm/s) N/N tot (%) δ (mm/s) ∆ (mm/s) LS
		250	100	0.382(3) 0.565(4)	-	-	-
		230	73.5	0.407(2) 0.564(4)	26.5	0.121(5)	2.70(1)
	Cooling	200 160 140	---	---	---	100 100 100	0.148(2) 2.721(4) 0.157(1) 2.732(2) 0.165(2) 2.723(4)
		120	-	-	-	100	0.175 (1) 2. 73(2)
		80	-	-	-	100	0.182(1) 2.742(3)
		250	18.7	0.37(4)	0.64(5)	81.3	0.128(2) 2.705(4)
	Heating	260 275 290 300	24.0 83.7 100 100	0.38* 0.383(2) 0.553(4) 0.59* 0.365(3) 0.548(6) 0.351(4) 0.544(6)	76.0 16.3 --	0.124(2) 2.716(5) 0.11* 2.65(4) ----
		310	100	0.354(5) 0.540(7)	-	-	-
	*These parameters were fixed during the fit.

Acknowledgments

This work was supported by the European program COST D35. P. Á. Sz. is grateful to the Embassy of France in Hungary, for her PhD grant. Lionel Rechignat and Jean-François Meunier (LCC) are acknowledged for their technical help and Marc Nardone (LNCMP, Toulouse) for the conception of the high pressure cell.


Fabien Treyssede
email: fabien.treyssede@utc.fr

Mabrouk Gwenael Gabard

Ben Tahar

Fabien Treyssède

Gwénaël Gabard

Mabrouk Ben Tahar

A mixed finite element method for acoustic wave propagation in moving fluids based on an Eulerian-Lagrangian description

Keywords: numbers: 43.20.Bi, 43.28.Py, 43.20.Mv

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

A mixed finite element method for acoustic wave propagation in moving fluids based on an Eulerian-Lagrangian description

I. INTRODUCTION

Propagation of acoustic disturbances in non-uniform flows is a subject of great interest in many practical problems, particularly in transport engineering with automotive exhaust systems, aeronautical turbofan engine inlet ducts, etc. The understanding of this phenomenon is a central feature for the prediction of noise and for designing components that efficiently attenuate sound. In practice, the shape of these components is often complex and flows are not uniform. Thus, the basic equations that describe such a problem must be able to cope with those complexities. Two kinds of formulations are mainly used: the linearized Euler's equations (LEE) and the fullpotential formulation. However as discussed further below, there exists another wave equation, which is a reformulation of the LEE.

The full-potential formulation is obtained from the LEE by assuming both flow and disturbances irrotationality. Thus, it constitutes a specific case of the general LEE. The corresponding propagation equation is scalar and written only in term of the acoustic velocity potential. This makes its resolution far easier and explains why it is much more widespread in the litterature. On the basis of the full-potential equation, some authors studied the effect of flow variation and multidirectionality upon sound propagation in ducts with variable cross-sections, using a finite element method (FEM) 1,2,3 or a boundary element method 4 . Besides, the analysis of pure propagation phenomena has naturally been extended to radiation by many authors by the use of various techniques (FEM combined with a boundary element method 5,[START_REF] Baumeister | Finite element-integral acoustic simulation of JT15D turbofan engine[END_REF]7,[START_REF] Tahar | Formulation variationnelle par équations intégrales pour le rayonnement acoustique en présence d'écoulement non-uniforme[END_REF][START_REF] Zhang | A coupled FEM/BEM formulation for acoustic radiation in a subsonic non-uniform flow[END_REF] or a wave-envelope element technique [START_REF] Astley | A finite element wave envelope formulation for acoustic radiation in moving flows[END_REF][START_REF] Roy | Improved finite element modeling of the turbofan engine inlet radiation problem[END_REF][START_REF] Eversman | Aft fan duct acoustic radiation[END_REF][START_REF] Eversman | Mapped infinite wave envelope elements for acoustic radiation in a uniformly moving medium[END_REF] , dual reciprocity boundary element method [START_REF] Lee | A dual-reciprocity method for acoustic radiation in a subsonic non-uniform flow[END_REF] , etc.).

Recently, Rienstra [START_REF] Rienstra | Sound transmission in slowly varying circular and annular lined ducts with flow[END_REF] developed an analytical model based on a multiple-scales method to study pure propagation in slowly varying cross-section ducts. Rienstra and

Eversman [START_REF] Rienstra | A numerical comparison between the multiplescales and finite-element solution for sound propagation in lined flow ducts[END_REF] made comparisons between the multiple-scales and FEM solutions that validates both models (though it highlights the limits of a multiple-scales method when reflections or conversions into other modes occur). In this paper, Rienstra's analytical model is used to give a reference solution that validates the FEM implementation of Galbrun's equation.

Because of its relative simplicity, the full-potential equation is a powerful formulation and may be sufficient to study sound propagation in flows. However, its main drawback is that it cannot take into account rotational mean flows. For instance, flow rotationality cannot be neglected when the effect of boundary layer refraction is important (see, for instance, Ref. 17,18,19) or when a mean flow swirl is present (this typically happens behind a rotor stage -see Ref. 20).

Actually, if the mean flow is rotational, the decomposition of the perturbations in terms of independent acoustic and rotational modes is no longer valid. The LEE must be directly solved. These equations represent a system of five equations and five unknowns, which can be reduced to four if the flow is supposed to be homentropic (pressure and density are then directly related). Given the complexity of the system, some attempts were made to simplify its solution. For instance, Nayfeh et al. [START_REF] Nayfeh | Transmission of sound through nonuniform circular ducts with compressible mean flows[END_REF] , Uenishi and Myers [START_REF] Uenishi | Two-dimensional acoustic field in a nonuniform duct carrying compressible flow[END_REF] developed a numerical method for solving acoustic propagation in ducts with variable cross-sections and mean sheared flow.

In the late 1970s, Abrahamson [START_REF] Abrahamson | A finite element algorithm for sound propagation in axisymmetric ducts containing compressible mean flow[END_REF] , Astley and Eversman [START_REF] Astley | A finite element for transmission in non-uniform ducts without flow: comparison with the method of weighted residuals[END_REF][START_REF] Astley | Acoustic transmission in non-uniform ducts with mean flow, part II: the finite element method[END_REF] implemented the direct LEE respectively for axisymmetric and 2D geometries, using a FEM. Good results were presented for simple cases but unfortunately, their work has not been continued. Astley and Eversman [START_REF] Astley | A finite element formulation of the eigenvalue problem in lined ducts with flow[END_REF] also outlined the existence of spurious modes when using a FEM to formulate the eigenvalue problem for a straight lined duct with shear flow.

More recently, Golubev and Atassi 27 studied a straight duct containing a mean flow with swirl and showed the coupling that occurs between acoustic and rotational modes.

Cooper and Peake 20 extended Golubev's study to slowly varying lined ducts by applying a multiple-scales method. Results showed the influence of the mean flow swirl, i.e. corotating modes are always much more damped than those in a non-swirling flow and counter-rotating modes may be amplified.

Furthermore, it must be noticed that the LEE have recently been implemented by several authors (see, for instance, Ref. 28,29,30) using numerical methods based on finite difference schemes. These references proved the operator ability to adequately propagate sound wave but, so far, numerical applications are still limited to simple geometries because of numerical difficulties and high computer resources requirement.

Through what is cited above, it can be seen that the use of the LEE may be essential and inescapable in several practical cases. However, this system is far more complex than the scalar full-potential equation, and it is thus rarely solved for general cases. In this paper, it is attempted to develop a general method based on a FEM to solve sound propagation problems with rotational mean flows. Instead of choosing the LEE, a nonstandard wave equation originally established by Galbrun 31 in the early 1930s is considered. As explained in the following section, this equation is derived from an Eulerian-Lagrangian description and is an exact reformulation of the LEE. Thus, Galbrun's equation a priori describes exactly the same physical phenomenon. Its main feature is that it constitutes a second-order linear partial differential equation written only in terms of the displacement perturbation.

Although only few works deal with this equation, it may be an interesting alternative to the LEE. Only acoustic displacement is involved (even in non-homentropic cases), which yields a gain of one to two unknowns compared to the LEE ; it also provides exact expressions of acoustic intensity and energy ; besides, boundary conditions are easily expressed because acoustic displacement (whose normal component is generally continuous at any interface between two media) appears explicitly, which avoid the somewhat difficult use of Myers' condition [START_REF] Myers | On the acoustic boundary condition in the presence of flow[END_REF] .

In 1985, Poirée [START_REF] Poiree | Les équations de l'acoustique linéaire et non linéaire dans un écoulement de fluide parfait[END_REF] detailed the Eulerian-Lagrangian description in order to derive Galbrun's equation and its extension to nonlinear problems. He also derived some general continuity conditions valid at any straight interface. Recently, Peyret and Elias [START_REF] Peyret | Finite-element method to study harmonic aeroacoustics problems[END_REF] proposed a direct displacement-based formulation of Galbrun's equation solved using a FEM. They also derived the same energy conservation law as Godin but with a different approach. Bonnet et al. [START_REF] Bonnet-Ben Dhia | Analyse mathématique de l'équation de Galbrun en écoulement uniforme[END_REF] pointed out the fact that the direct displacement-based formulation associated with Galbrun's equations does not necessarily converge with standard finite elements and proposed a method to regularize the variational formulation in the case of a uniform flow. These two last references are discussed in detail later (see Sec. III).

In this paper, a mixed variational formulation based on the pressure-displacement variables is presented in order to avoid some spurious solutions. Though the overall method is quite general, finite element discretization and numerical results are presented for the axisymmetric case. A first set of results consists in comparing FEM with semianalytical solutions for a straight duct containing a sheared flow (obtained from Pridmore-Brown equation). A second set of results concerns a more complex duct geometry with a potential flow and is compared to results obtained from a multiple-scale method (which is an adaptation of Rienstra's work for the incompressible case -see Ref.

15).

II. THEORY

Compared to the classical Eulerian description, the Eulerian-Lagrangian description is not usual. This section briefly recalls the latter before giving Galbrun's equation and the associated expression of acoustic intensity. This section does not give details about calculations but reviews can be found in Ref. 31,33,34.

A. The Eulerian-Lagrangian description

In the Eulerian description, perturbations are Eulerian and expressed with Eulerian variables. This description is implicitly used when Euler's equations are directly linearized. The Eulerian-Lagrangian description consists in considering Lagrangian perturbations of quantities expressed in terms of Eulerian variables. In order to give a physical comprehension of what Lagrangian and Eulerian perturbations mean, the perturbed (or total) and non-perturbed (or mean flow) states have to be described first.

Define the geometrical position x of a given particle in the mean flow configuration and its position y in the perturbed configuration. Then, if w L denotes the displacement perturbation of this particle, x and y are related by:

ε = + L y x w (2.1)
In the remainder of this paper, Ψ represents any physical quantity (tensor of arbitrary order) and the subscript 0 is used to distinguish mean flow quantities from their total (or perturbed) counterpart. Then, two kinds of perturbations can be defined:

( ) ( ) ( ) ( ) 0 0 , , , , E L t
( ) ( ) 2 2 T 0 0 0 0 0 2 0 d c p p dt ρ ρ -∇ ∇⋅ + ∇⋅ ∇ -∇ ⋅∇ = L L L L w w w w (2.5)
and the Lagrangian density perturbation is explicitly given by: In order to derive Galbrun's equation, perfect fluid and isentropic assumptions have implicitly been made by starting from the LEE. For Lagrangian perturbations, the isentropic assumption leads to the well-known pressure-density relationship:

0 L ρ ρ = -∇ ⋅ L w (2.
2 0 L L p c ρ = (2.7)
which yields the following explicit equation for the Lagrangian pressure:

2 0 0 L p c ρ = - ∇ ⋅ L w (2.8)
Some important remarks should now be adressed. Unlike Eulerian perturbations, Eq.

(2.7) remains valid even for inhomogeneous media and/or for non-homentropic flows. In these situations, solving the LEE would require an additional equation, that is the linearized energy equation (see, for instance, Ref. 38). In fact, for Lagrangian perturbations, Eq. (2.7) still holds because the equation of state applies for a given particle so that the thermodynamical system is closed.

Thus, Galbrun's equation obviously provides a first advantage (compared to the LEE) because this equation is expressed with the Lagrangian displacement as a unique unknown, which yields a gain of one to two variables.

Another asset is the existence (but not uniqueness) of a Lagrangian density associated with Galbrun's equation (see Ref. 36). This yields an exact energy conservation law and exact expressions for the acoustic energy and intensity. These expressions can be found in Ref. 34 and 36. In particular, the acoustic intensity is given by: ( )

0 0 L d p p t dt t ρ   ∂ ∂ = ⋅ + -⋅∇   ∂ ∂   L L L L 0 w w w i v w (2.9)
It can be noticed that the Lagrangian displacement inevitably appears in this expression, which likely explains why no exact formulation for the intensity has been found based on a pure Eulerian description.

III. NUMERICAL METHOD

In this section, Galbrun's equation is solved using a FEM. From now on, fluctuations are assumed to have an e -iωt time dependence. Besides, the mean flow is steady and the mean pressure is supposed to be constant for simplicity (dropping the last two terms of Eq. (2.5)). This assumption is not valid for aeroacoustic problems but our purpose is to solve Galbrun's equation with a FEM (as outlined by Peyret and Elias 36 , terms with p o do not present a major interest from a numerical point of view). Under these assumptions, Eq. (2.5) becomes:

( ) ( ) 2 2 0 0 0 0 0 2 0 i c ρ ω ωρ ρ ρ - - ⋅∇ + ⋅∇ ⋅∇ -∇ ∇⋅ = L L L L 0 0 0 w v w v v w w (3.1)
The first subheading of this section gives a brief review about the numerical difficulties of the direct displacement based formulation associated to (3.1). In the second subheading, a displacement-pressure based mixed formulation is proposed to overcome these difficulties. A third one deals with boundary conditions and the last subheading gives some important details about the finite element discretization of the mixed formulation.

A. Displacement based formulation

Equation (3.1) is multiplied by a trial field w * and integrated over the domain Ω.

Then, integrating by part the last two terms (which imply second order derivatives) and half the second term yields the following direct displacement based variational formulation:

( )( ) ( ) ( ) ( ) ( ) ( ) 2 * 2 * 0 0 0 * * 0 0 * 0 * * 0 0 L S c d d i d i d d d p dS dt ρ ω ρ ω ρ ω ρ ρ ρ Ω Ω Ω Ω Ω ∇⋅ ∇⋅ Ω- ⋅ Ω - ⋅ ⋅∇ Ω+ ⋅∇ ⋅ Ω - ⋅∇ ⋅ ⋅∇ Ω   + ⋅⋅+ ⋅= ∀     ∫ ∫ ∫ ∫ ∫ ∫ L L L L 0 0 L 0 0 L 0 w w w w w v w v ww v w v w w w v n n w (3.2)
It can be noticed that the above formulation is presented in such a way that the first line corresponds to the no-flow case and the second and third one to the presence of mean flow (the last line is a boundary integral).

This formulation is used in Ref. 36. However, when standard finite elements are implemented to discretize the formulation, solutions are generally corrupted, even in the no-flow case. This phenomenon is purely numerical. An example of spurious solution is given in Sec. IVA.

Bonnet & al. [START_REF] Bonnet-Ben Dhia | Analyse mathématique de l'équation de Galbrun en écoulement uniforme[END_REF] have recently proposed a regularized formulation of Galbrun's equation in the uniform flow case (with this method, some specific terms are added in the formulation (3.2)). A good convergence is obtained but limitations of the method arise for the generalization to arbitrary mean flows.

The no-flow case was first studied in the 1970s when considering vibrations of coupled fluid-structure systems (see, for instance, Ref. 39) and was proved to exhibit spurious circulation modes with non-zero frequencies. Basically, this phenomenon is due to a bad accuracy of the divergence and curl (calculated from derivatives of displacements), which in turn affects the displacement prediction itself. In order to cope with this numerical phenomenon, several methods have been proposed, such as the penalty method [START_REF] Hamdi | A displacement method for the analysis of vibrations of coupled fluid-structure systems[END_REF] , edge finite element [START_REF] Bermudez | Finite element computation of three-dimensional elastoacoustic vibrations[END_REF] , and mixed finite element method [START_REF] Bathe | A mixed displacement-based finite element formulation for acoustic fluid-structure interaction[END_REF][START_REF] Wang | Displacement/Pressure based mixed finite element formulations for acoustic fluid-structure interaction problems[END_REF] .

In the presence of mean flow, the penalty method and edge finite element method cannot be directly applied because the displacement field is generally no more irrotational. Thus, the method chosen in this paper is naturally based on a mixed finite element formulation. To conclude this subheading, one emphasizes that the overall problem in the now-flow case is typically analogous to incompressible elasticity or fluid [START_REF] Bathe | Finite Element Procedures[END_REF] and electromagnetics [START_REF] Jin | The Finite Element Method in Electromagnetics[END_REF] , and is often referred to as "locking" in the litterature.

B. Mixed formulation

A mixed variational formulation based on pressure-displacement variables is now derived. Although pressure is now an explicit unknown of Galbrun's equation (as in the LEE), the efficiency of such a formulation to prevent corruption by spurious solutions has already been demonstrated by Wang and Bathe [START_REF] Wang | Displacement/Pressure based mixed finite element formulations for acoustic fluid-structure interaction problems[END_REF] in the no-flow case when considering fluid-structure interactions. It was also extensively and rigorously analyzed by Brezzi and Fortin 45 in a general way (see also Bathe [START_REF] Bathe | Finite Element Procedures[END_REF] for a more engineering-oriented approach applied to incompressible elasticity).

In the presence of arbitrary mean flow, an analogous mixed formulation of Ref. 42 can be obtained by replacing the last term of Eq. (3.1) with the pressure gradient (use of Eq. (2.8) is made). It yields the following system:

( )

2 0 0 0 2 0 0 2 0 0 L L i p p c ρ ω ωρ ρ ρ  - - ⋅∇ + ⋅∇ ⋅∇ +∇ =   + ∇⋅ =   L L L 0 0 0 L w v w v v w w (3.3)
Multiplying both equations by trial fields w * and p * respectively, integrating over the domain Ω and then by parts gives the following mixed variational formulation:

( ) ( ) ( ) ( ) ( ) ( ) { } * * * 2 * 0 2 0 0 * * 0 0 * 0 * * * * 0 1 0 , L L S S ppd p d pd d c i d i d d d dS p dS p dt ω ρ ρ ω ρ ω ρ ρ ρ Ω Ω Ω Ω Ω Ω Ω - Ω+ ∇ ⋅ Ω+ ⋅∇ Ω- ⋅ Ω - ⋅ ⋅∇ Ω+ ⋅∇ ⋅ Ω - ⋅∇ ⋅ ⋅∇ Ω   + ⋅ ⋅ - ⋅ = ∀     ∫ ∫ ∫ ∫ ∫ ∫ ∫ ∫ ∫ L L L L 0 0 L 0 0 L L 0 w w w w w v w v ww v w v w w w v n w n w (3.4)
where S is the surface enclosing the acoustic domain Ω and n is the outward normal. In order to show the efficiency of a mixed formulation, a comparison with the displacement based formulation is given in Sec. IVA.

C. Boundary conditions

Boundary conditions associated with Galbrun's equation must be carefully applied. At walls, the first surface integral of (3.4) always vanishes because the mean flow normal velocity is zero. For rigid walls, normal displacement also equals zero and the second integral also vanishes. For an absorbing wall, the adequate boundary condition is obtained from the normal displacement fluid-wall continuity. Given the wall impedance Z, this condition is simply given by the following relationship:

1 L S p i Z ω = - ⋅ L w n (3.5)
which is simply applied by replacing the normal displacement in the second boundary integral. The fact that displacement is an explicit variable in the proposed formulation makes the impedance condition (3.5) simpler to implement than in the LEE case, that would have required the use of Myer's condition [START_REF] Myers | On the acoustic boundary condition in the presence of flow[END_REF] .

When an arbitrary mean flow is present, the displacement field may not be irrotational. Thus, it is obvious that a fixed pressure inside the fluid (on S 2 ) is not sufficient to determine a unique solution. For instance, applying pressure (which represents the displacement divergence from Eq. (2.8)) at the duct inlet and outlet means that the rotational part of the displacement is left free. Consequently, one must impose the total displacement field everywhere a boundary condition is required inside the fluid (typically at the duct inlet-outlet). This condition is explicitly given by:

2 S = L w w (3.6)
and is directly enforced at nodes as a constraint in the FEM model. This makes the first surface integrals in (3.4) vanish because w * ≡δw=0 on S 2 (i.e. forced boundary condition). Note that in many instances, it is difficult to specify the particle displacement on both surfaces labeled S 2 because only the incident component is known. A boundary condition based on a multi-modal decomposition technique overcomes this difficulty.

As a side remark, in the no-flow case, a surface with fixed pressure inside the fluid is sufficient to uniquely determine the solution because displacement is implicitly irrotational.

D. Finite element discretization

It has been proven that interpolations for displacement and pressure must be adequately chosen. The choice of a mixed formulation is not sufficient to avoid locking, and one also has to choose appropriate interpolation functions. Based on Ref. 43 and 45, this subheading gives some details about what kind of finite elements has to be used.

If a bad choice of interpolating functions is made, the element may lock and/or give some spurious pressure modes. A criterion that ensures convergence and stability of the finite element is given by the so-called inf-sup condition. Though not necessary, this condition is a strong guaranty of reliability. Details about the inf-sup condition are not given in this paper but can be found in Ref. 43 and 45. For a more practical use, a numerical test of this condition has also been developed by Chapelle and Bathe [START_REF] Chapelle | The inf-sup test[END_REF] . One example of triangular 2D-element satisfying the inf-sup condition is given on Fig. 2.

This element may be referred to as the "P 1 + -P 1 ", "4-3c" or "MINI" element [START_REF] Bathe | Finite Element Procedures[END_REF][START_REF] Brezzi | Mixed and Hybrid Finite Element Methods[END_REF] and is the one used in this paper. An example of higher order element is also given on the same figure (element "9-4c"). Note that some other 2D (or 3D) elements can be found in Ref.

43, 46.

In the no-flow case, elements that satisfy the inf-sup condition have successfully been implemented by Wang and Bathe [START_REF] Wang | Displacement/Pressure based mixed finite element formulations for acoustic fluid-structure interaction problems[END_REF] . The originality of Galbrun's equation is that the pressure-displacement relationship given by (2.8) is not altered by the presence of flow and is strictly identical to the no-flow case. This enables to directly apply the infsup condition to a Galbrun-based formulation. Thus, under the assumption that the additional operators introduced by the presence of flow (terms in v o .∇) does not alter convergence properties of elements satisfying the inf-sup condition, it is expected that the proposed mixed FEM for solving Galbrun's equation is robust.

When flow is present, it should be noticed that the mass continuity equation obtained from the LEE implies material derivatives and no equation similar to Eq. (2.8) can be obtained. The inf-sup condition cannot then be directly applied when the LEE are considered. Thus, given the current knowledge about mixed formulations, Galbrun's equation seems to be more interesting from a numerical point of view.

The general variational formulation (3.4) is now restricted to axisymmetric geometries. Introducing the cylindrical coordinates, the following θ-dependence is set: The (2D) acoustic domain Ω is chosen to be meshed with 4-3c elements. On the reference element, displacement and pressure variables are thus interpolated as follows:

( ) ( ) ( ) ( ) , , , , , , im L
( ) ( ) ( ) ( ) ( ) 2 3 123 , 1 1 , 1 L u v u v u v u v uv p u v u v p up vp  = -- + + + --   = -- + +   L 1 w www a (3.9)
where the subscripts i (i=1,2,3) denotes the node number. It can be seen for the displacement, that the standard linear interpolation is enriched by introducing a generalized variable a. The term (1-u-v)uv represents a bubble function: this polynomial is null on the three side of the triangle and thus maintains the compatibility (C 0 continuity).

Then, elements have four degrees of freedom per node, plus three internal (three for each component of w L

). However, these internal degrees of freedom can be condensed out before the element are assembled, which is attractive from a computational point of view.

After assembling and applying boundary conditions, the global discretized variational formulation yields the following algebraic system:

= r r r K w f (3.10)
K r is a ω -dependent complex band matrix, unsymmetrical when flow is present. A sparse storage is chosen. For a fixed ω, the unknown nodal displacement vector w r is finally obtained by using a LU decomposition.

IV. RESULTS

In this section, the FEM numerical method is validated with two semi-analytical models. The first model corresponds to the well-known Pridmore-Brown equation, the second to Rienstra's multiple scale approximation (in the incompressible mean flow case in order to fit with the assumption of the FEM model -see Sec. III). Both models represents the propagation of a given (m,n) mode in an infinite duct (m and n denote respectively the angular and radial mode number).

They may be considered as complementary. Pridmore-Brown equation deals with a simple straight duct but a possibly sheared mean flow. Boundary layer effects upon propagation can thus be considered. In Rienstra's multiple-scales method, the mean flow must be potential but the duct is slowly varying, which permits to study more complex geometries.

Unlike Pridmore-Brown equation, Rienstra's model constitutes an approximation. It cannot be exact because of modal reflection and scattering that may occur in a varying duct. These limitations of a multiple-scales method have been highlighted by Rienstra and Eversman [START_REF] Rienstra | A numerical comparison between the multiplescales and finite-element solution for sound propagation in lined flow ducts[END_REF] and are also demonstrated in this paper.

In the following, iso-pressure contours are given in modulus for all plots. Units are chosen to be in Pa in order not to minimize errors. The averaged intensity vector may also be given on the θ=0° plane. Mean flow velocities are defined in Mach number (M).

Propagation and flow directions are also sketched in order to explicitly show if wave propagation is upstream or downstream. Besides, typical values of ρ 0 =1.2kg.m -3 and c 0 =340m.s -1 are used.

Test cases sweep a non-dimensional frequency range up to about kR=20 and the duct geometry is always meshed with a λ/10 finite element length. This is the estimated criterion for a satisfactory convergence but a λ/8 criterion may be sufficient to give good results. These criteria are yet modulated by a 2 1 M factor due to the Doppler effect when flow is present.

A. Validation for straight ducts (Pridmore-Brown equation)

For FEM computations of this subheading, the methodology is as follows. The Lagrangian displacement is calculated from the Pridmore-Brown model and is then imposed at the duct inlet of the FEM model (in the remaining, the terms "inlet" and "outlet" are used from an acoustical point of view -see the direction of propagation upon each figure -and not from a mean flow point of view). A modal non-reflective boundary condition is preferred at the outlet, which is less constraining (in particular, the phase and amplitude are left free). Hence, the mode being enforced at the inlet/outlet boundaries, the FEM solution inside the duct is computed and compared to the semianalytical one.

The non-reflective boundary condition is inspired from Peyret and Elias [START_REF] Peyret | Finite-element method to study harmonic aeroacoustics problems[END_REF] . It simulates the one-way propagation and applies only for a given (m,n) mode. In classical acoustics (no-flow case), one usually imposes on the cross-section an anechoic termination based on the wave impedance (i.e. the ratio of the pressure and the normal acoustic velocity). Here, an additional condition is required with the displacement variable in order to determine a unique solution. The velocity and the normal displacement which appear in the boundary integral of Eq. (3.4) are thus replaced with the following expressions:

( ) 0 1 ; L nr d p dt i Z ρ ω ⋅ = ⋅ = - L L L 0 nr w v n Z w w n (4.1)
Z nr is the standard modal non-reflecting impedance and Z nr can be viewed as a modal non-reflective matrix impedance. In the Pridmore-Brown model, material derivatives are simply given by d/dt=-i(ω-v 0 k z ) because an (

) z i k z m t e θ ω --
dependence is chosen for the modal acoustic variables. This yields the following explicit expressions for the modal impedances defined by (4.1):

( ) ( ) 2 0 0 0 0 0 ; z z nr z v k i v v k Z k ω ρ ω ρ ω - = - - = nr Z I (4.2)
where v 0 is the axial mean flow velocity and k z is the outlet modal axial wave-number, which is part of the semi-analytical solution. I denotes the identity matrix.

In this subheading, an annular straight duct is considered. The inner and outer radius are respectively 0.2 and 1m. Geometry with a typical mesh is shown by Fig. 3.

The first test case, given in Fig. 3, shows the efficiency of a mixed FEM compared to a displacement based one. A comparison between the semi-analytical, displacement based FEM and mixed FEM solutions is given for the pressure modulus. The example concerns a (0,1) mode at f=250Hz propagating in a hard wall duct without flow. Solution obtained with a displacement based is totally corrupted by rotational fields.

These spurious solutions are small-scale artifact which are trapped within the model grid. Unfortunately, this numerical problem does not disappear at all when the mesh is refined. For a mixed FEM, results are clearly in good agreement with the semi-analytical solution.

The second test case (Fig. 4) concerns a (8,0) mode propagating at f=500Hz in a perfectly rigid wall duct. This mode is calculated for M=0 and M=-0.4 (upstream propagation). Agreement between Pridmore-Brown and FEM solutions is perfect for both Mach number. This test case aims at pointing out one of the effect of convection upon cut-off frequencies, say a decrease by a 2 1 M factor. Figure 4 clearly shows that without flow, the (8,0) mode is cut-off whereas for M=-0.4, it becomes cut-on. In fact, for the duct dimensions, the (8,0) mode has an exact cut-off frequency of 522.0Hz without flow, decreased to 478.5Hz at M=-0.4. This frequency is lower than 500Hz, which explains why the mode fully propagates along the duct. Examining the intensity plot for both Mach number gives directly the nature of modes: unlike at M=-0.4, the acoustic intensity vector is null at M=0.0 which proves that the mode is cut-off (energy does not propagate). For the M=-0.4 case, a comparison in terms of the pressure phase is also given by Figure 5, where a good agreement between solutions can be observed.

The third test case (Fig. 6) exhibits a (10,1) mode at f=1000Hz. Walls are lined and the impedance value is Z=2040(1-i) for both inner and outer walls. The cross-section averaged Mach number is M=-0.5 (upstream propagation). In this example, a comparison between uniform and sheared flow is given in order to demonstrate the capability of the FEM approach to take into account refraction phenomena. The shear flow is arbitrarily chosen to have a boundary layer thickness of 10% (δ=0.08 m), with the same mass flow rate that the uniform profile (see Fig. 7 for the mean flow velocity profiles). This thickness is not realistic but voluntarily exaggerated in order to illustrate refraction. Note that flows with a boundary layer are obviously rotational and cannot be considered with the full-potential propagation equation.

As seen on Fig. 6, Pridmore-Brown and FEM solutions show satisfactory agreement. In particular, modal shape and attenuation are conserved in the FEM model.

Results show a strong difference in amplitude between the uniform and sheared cases. In fact, for an upstream propagation, the mean flow velocity gradient due to the presence of the boundary layer tends to refract waves toward the center, thus yielding a weaker attenuation than in the uniform case. (This is the opposite for a downstream propagation, for which waves are refracted toward walls -see, for instance Ref. 17). Attenuation coefficient values can directly be obtained from the semi-analytical model. These coefficients are 5.6 dB.m -1 and 0.9 dB.m -1 respectively for the uniform and sheared cases. Thus, for the duct being considered, an error of about 9 dB is made at the outlet when a uniform profile is supposed.

B. Validation for varying ducts (Rienstra's model)

The test geometry taken is now varying (but the flow is restricted to be potential).

This geometry is the same as in Rienstra's papers [START_REF] Rienstra | Sound transmission in slowly varying circular and annular lined ducts with flow[END_REF][START_REF] Rienstra | A numerical comparison between the multiplescales and finite-element solution for sound propagation in lined flow ducts[END_REF] and is representative of a turbofan aircraft engine. It includes a circular-to-annular transition (a central body is thus present).

There are differences between the multiple scale method and the FEM. The FEM formulation admits the propagation of many modes (reflection and scattering are integral part of the solution). On the contrary, the multiple scale approximation lies in supposing that a single given mode is propagating in a single direction inside a duct. the acoustic variables at the duct inlet and outlet via an eigenmode expansion. For FEM calculations, the method is as follows. On the inlet plane, the complex amplitude of the appropriate mode is specified via a forced boundary condition. On the outer plane, reflected mode amplitudes are set to zero, imposing a non-reflecting boundary condition.

Inlet reflected and outlet transmitted mode amplitudes are unknown and part of the solution. In this paper, ten radial modes have been used for the expansion (given the cutoff frequencies of higher order modes, this is sufficient for the test cases presented in the following). Note that the eigenmode expansion must be done with hard walls for orthogonality (soft wall modes are not orthogonal). Thus, for test cases with lined walls, the FEM geometry has been extended up to 0.5m both at the inlet and outlet by adding pieces of straight hard wall ducts (not shown on figures).

In the FEM model, the potential mean flow is first computed, also via a FEM solving Laplace's equation and then used for solving Galbrun's equation. Results

presented in this subsection give comparisons between the solutions obtained by Rienstra's multiple-scales method (mono-modal) and the FEM model proposed in the present paper (FEM solutions of Rienstra and Eversman are not shown).

The first test case is depicted in Fig. 8, where geometry, mesh and mean flow are presented. At the inlet and outlet, the axial local Mach number values are M=0.3 and 0.25 respectively. It must be observed that the outer wall is straight, which constitutes a slight difference with Rienstra's geometry (considered in the second test case). This essentially aims at minimizing reflections due to the outer wall.

In this first set of results (Fig. 9), a downstream (1,1) mode propagating in a rigid wall duct is considered. Computations are achieved for three frequencies. At 300 and 360Hz, a very good agreement is obtained between FEM and multiple-scales solutions, which validates the FEM code. At these frequencies, the only cut-on modes with m=1 are (1,0) and (1,1). Other cut-on modes with m≠1 are not considered by the FEM model because m is a fixed parameter in the code. Consequently, the good convergence between both models demonstrates that scattering into the (1,0) mode as well as reflection are almost negligible.

However, when the frequency f=420Hz is reached, a strong difference is observed.

This disagreement is likely explained by partial scattering into the (1,2) mode. In fact, analyzing local cut-off frequencies shows that this mode is cut-on at the outlet (its local cut-off frequency is 416.8Hz). Hence, it can be deduced from the difference observed that for f=420Hz, the multiple-scales approximation fails.

The last test case (see Fig. 10) concerns a (7,0) mode at f=500Hz propagating into the exact Rienstra's geometry. The outer wall is lined, the impedance value is Z=4080(1+i). The central body is left perfectly rigid. Mesh and flow are not shown for conciseness (M=0.5 and 0.49 at the inlet and outlet respectively). The goal of this example is to outline the reflection phenomenon that limits the use of a multiple-scales method. The frequency of f=500Hz is chosen in order for the (7,0) mode to be the only cut-on mode. Calculations effectively give a maximum local cut-off frequency along the duct of 411.2Hz for this mode and a minimum cut-off frequency of 558.6Hz for the (7,1) mode. This indicates that the (7,0) mode is always cut-on and the (7,1) always cut-off.

This permits to avoid any significant scattering into other modes and thus to focus on auto-reflection of the (7,0) mode only. Computations are made for both downstream and upstream propagation.

In the downstream case, a good agreement is obtained. This shows that only few reflections are produced inside the duct for this direction of propagation. In the upstream propagation case, some differences occur (on this plot, wave is propagating from the top to the bottom). Some wiggles appear and iso-pressure contours are not totally smooth. At the acoustic outlet (bottom), it can be seen on plots that the attenuation obtained by the FEM is a little greater than the semi-analytical one, which tends to prove that reflections of the (7,0) mode on itself are not negligible. This may be attributed to the central body as well as the abrupt change of outer radius located at the acoustic inlet, both viewed as a narrowing for an upstream propagation.

Finally, it may be interesting to look at the acoustic intensity vector that has been plotted for both the upstream and downstream case. Because the lining of the wall absorbs some energy, intensity is not exactly parallel to the wall but penetrates into it. This is more visible in the upstream case, which is coherent with the fact that the downstream wave is less attenuated.

V. CONCLUSION

In this paper, a mixed FEM based on Galbrun's equation has been proposed to solve sound wave propagation in arbitrary flows. Compared to the LEE, this equation may have several advantages. From a theoretical point of view, it yields an exact expression for the acoustic intensity and simpler boundary conditions (especially in case of absorbing walls). From a numerical point of view, Galbrun's equation allows the direct application of the inf-sup condition, already encountered in the study of mixed FEM for incompressible media.

Results obtained with the proposed mixed FEM have been compared with two complementary semi-analytical models and have been found to be in very good agreement. The comparison with the Pridmore-Brown equation has proven the efficiency of the numerical method to take into account convection and refraction from a boundary layer, which cannot be considered with a full-potential formulation. Comparisons with a multiple-scales method have fully validated the FEM for complex geometry and have also confirmed limitations of a mutliple-scales approach when some significant reflection or diffraction occur.

Those results show that a mixed FEM method based on Galbrun's equation could be an interesting alternative to a finite-difference method based on the LEE, for solving aeroacoustic problems.

6 )

 6 Eq. (2.6) constitutes the mass continuity equation obtained from a Lagrangian-Eulerian description. It simply states that density fluctuations are balanced by dilatation fluctuations.

  As in Eq. (3.2), the first line of the mixed formulation (3.4) represents the no-flow operators. These operators are almost identical to those used by Wang and Bathe[START_REF] Wang | Displacement/Pressure based mixed finite element formulations for acoustic fluid-structure interaction problems[END_REF] in their mixed formulation. The only slight difference is that one has chosen to integrate by parts the second term of the second equation of system (3.3) (instead of the last term of the first equation) in order to let the normal displacement appear explicitly at the boundary.Normal displacement at walls can thus be easily imposed (the second surface integral of Eq. (3.4) simply disappear for perfectly rigid walls). Besides, for fluidstructure interactions (not considered in this paper), normal displacement continuity could also be easily imposed by replacing the fluid normal displacement with the structure normal displacement.

Figure 1

 1 Figure 1 represents a typical problem of propagation inside a duct carrying flow. Two

  mode number m is a parameter of the solution. The weighting functions are taken as:

1R

  . K. Sigman, R. K. Majjigi, and B. T. Zinn,. "Determination of turbofan inlet acoustics using finite elements," American Institute of Aeronautics and Astronautics Journal 16, 1139-1145 (1978).

2R.

  Fuller, and D. A. Bies, "The effects of flow on the performance of a reactive acoustic attenuator," Journal of Sound and Vibration 62, 73-92 (1979). 3 Cabelli, "The influence of flow on the acoustic characteristics of a duct bend for higher order modes -a numerical study," Journal of Sound and Vibration 82, 131-149 (1982).4 Z. L. Ji, Q. Ma, and Z. H. Zhang, "A boundary element scheme for evaluation of four-pole parameters of ducts and mufflers with low mach number non-uniform flow," Journal of Sound and Vibration 185, 107-117 (1995).

5S

  . J. Horowitz, R. K. Sigman, and B. T. Zinn, "An iterative finite element-integral technique for predicting sound radiation from turbofan inlets in steady flight," American Institute of Aeronautics and Astronautics Journal 24, 1256-1262 (1986).

FIGURE CAPTIONS : FIG. 1 :FIG. 6 :FIG. 7 :FIG. 9 :FIG. 10 :

 :167910 FIGURE CAPTIONS:

  Godin 34 independently obtained Galbrun's equation in a quite different way. In particular, he derived exact expressions of acoustic energy and intensity, and boundary conditions for free, rigid, and absorbing surfaces in terms of the displacement. Ben Tahar and Goy[START_REF] Tahar | Resolution of a vibroacoustic problem in the presence of a nonuniform mean flow[END_REF] developed a variational formulation to study vibroacoustic problems with mean flow.