Aymeric Vié
email: aymeric.vie@em2c.ecp.fr

Stéphane Jay
email: stéphane.jay@ifpen.fr

Bénédicte Cuenot
email: benedicte.cuenot@cerfacs.fr

Marc Massot
email: marc.massot@em2c.ecp.fr

Accounting for Polydispersion in the Eulerian Large Eddy Simulation of the Two-Phase Flow in an Aeronautical-type Burner

Keywords: Large Eddy Simulation, Two-phase flow, Polydisperse sprays, Eulerian Mesoscopic, Multifluid, Aeronautical burner

Introduction

Considering available computational resources, Large Eddy Simulation (LES) becomes more and more attractive for industrial applications. For single-phase reacting flows, LES has already shown its great potential for aeronautical combustors [START_REF] Constantinescu | A numerical method for large eddy simulation in complex geometries[END_REF][START_REF] Lartigue | Compressible large-eddy simulation of turbulent combustion in complex geometry on unstructured meshes[END_REF] or automotive engines [START_REF] Laget | LES Calculations of a four cylinder engine[END_REF].

For two-phase reacting flows, simulating the whole physics from injection to turbulent dispersion is still a challenge as the simulation of the primary atomization requires accurate models and high computational resources [START_REF] Ménard | Coupling level set/VOF/ghost fluid methods: Validation and application to 3D simulation of the primary break-up of a liquid jet[END_REF][START_REF] Desjardins | An accurate conservative level set /ghost fluid method for simulating turbulent atomization[END_REF]19,[START_REF]Primary atomization of a liquid jet in crossflow[END_REF]. Furthermore most methods used to simulate the primary atomization are not adapted to subsequent turbulent spray. Fortunately, the primary atomization occurs in a limited zone near the injector and, in first approximation, can be accounted for by representative boundary conditions [START_REF] Martinez | A model for the injection boundary conditions in the context of 3D simulation of diesel spray : methodology and validation[END_REF][START_REF] Sanjosé | Fuel Injection Model for Euler-Euler and Euler-Lagrange Large-Eddy Simulations of an evaporating spray inside an aeronautical combustor[END_REF].

Assuming point-particles by comparison to the smallest flow resolved scales (estimated by the mesh size), the disperse phase statistics can be described by the number density function (NDF), which is solution of the Williams-Boltzmann Equation (WBE, [START_REF] Williams | Spray Combustion and atomization[END_REF]), describing its evolution in phase space. Solving the WBE directly is completely unreachable, due to the high number of dimensions of the phase space. Instead Lagrangian methods, like the Direct Simulation Monte Carlo [START_REF] Bird | Molecular gas dynamics and the direct simulation of gas flows[END_REF] method, build statistics by tracking stochastic particles, or like the Discrete Particle Simulation (DPS, [START_REF] Riley | Diffusion experiments with numerically integrated isotropic turbulence[END_REF]), consist in tracking particles of one individual realization of the disperse phase. Numerous studies have used this approach with success [START_REF] Apte | Large-eddy simulation of evaporating spray in a coaxial combustor[END_REF][START_REF] Patel | Simulation of spray-turbulence-flame interactions in a lean direct injection combustor[END_REF][START_REF] Jaegle | Eulerian and Lagrangian spray simulations of an aeronautical multipoint injector[END_REF], and shown its great potential (ease of modeling, computational time with small number of particles), but also its main drawbacks (projection on the Eulerian grid to couple with the gas phase, computational time with high number of particles, statistical convergence). Finally, Eulerian approaches solve the moments of the NDF (such as mean number density or volume fraction). This implies a loss of information, as moments are integrals of the NDF, but results in Eulerian equations to which classical methods for the gas phase may be directly applied. The CPU cost then depends only on the mesh size. In such methods, two main issues arise, concerning the modeling of the effect of turbulence on the spray, and the polydispersion in size of the spray. Those two issues are of primary importance for industrial applications, as they are responsible for the fuel mass fraction distribution in the combustion chamber, controlling ignition and combustion in the full device.

Turbulent spray dispersion models in the Eulerian formulation depend on the Stokes number range. For small Stokes number (St << 1), equilibrium approaches can be used [START_REF] Ferry | Equilibrium Eulerian approach for predicting the thermal field of a dispersion of small particles[END_REF], which only consider the diffusive effect of turbulence on the zeroth order moment for low Stokes numbers. For moderate Stokes number (St ≈ 1), the disperse phase dynamics become more complex: strong segregation effects occur and the disperse phase motion exhibits correlated and uncorrelated components [START_REF] Simonin | On the spatial distribution of heavy particle velocities in turbulent flow: from continuous field to particulate chaos[END_REF][START_REF] Reeks | On the probability density function equations for particle dispersion in a uniform shear flow[END_REF]. To describe such physics, Février et al. [START_REF] Février | Partitioning of particle velocities in gas-solid turbulent flow into a continuous field and a spatially uncorrelated random distribution: theoretical formalism and numerical study[END_REF] have developed the Mesoscopic Eulerian Formalism (MEF), which considers one droplet size only, but introduces a velocity dispersion around the mean velocity. Kaufmann et al. [START_REF] Kaufmann | Comparison between Lagrangian and mesoscopic Eulerian modelling approaches for inertial particles suspended in decaying isotropic turbulence[END_REF] and Moreau et al. [START_REF] Moreau | Development of gas-particle eulereuler LES approach: a priori analysis of particle sub-grid models in homogeneous isotropic turbulence[END_REF] addressed the problem of accurate closures for the uncorrelated stresses imposed by the velocity dispersion in the LES context. The resulting model was applied to automotive engine injection [START_REF] Martinez | Large eddy simulation of fuel sprays using the Eulerian mesoscopic approach[END_REF], aeronautical jets in crossflow [START_REF] Jaegle | Large eddy simulation of evaporating sprays in complex geometries using Eulerian and Lagrangian methods[END_REF], or aeronautical configurations such as the one considered here [START_REF] Sanjosé | LES of two-phase reacting flows[END_REF][START_REF] Sanjosé | Évaluation de la méthode Euler-Euler pour la simulation numérique aux grandes échelles des chambres carburant liquide[END_REF]. More recently, Zaichik et al. [START_REF] Zaichik | An Eulerian approach for large eddy simulation of particle transport in turbulent flows[END_REF] proposed a LES formalism based on the work of Reeks [START_REF] Reeks | On a kinetic equation for the transport of particles in turbulent flows[END_REF] and Zaichik et al. [67], taking into account the effect of subgrid turbulence directly at the kinetic level. This approach is still to be applied in more complex configurations.

Taking into account size polydispersion adds another complexity to Eulerian models. A first method consists in sampling the NDF by a set of monodisperse Eulerian phases [START_REF] Laurent | Multi-fluid modeling of laminar poly-dispersed spray flames: origin, assumptions and comparison of sectional and sampling methods[END_REF], as is done in a Lagrangian stochastic approach. This method is relatively simple when Eulerian phases do not interact, i.e. when no collision occur. However this approach suffers from intrinsic poor statistical convergence, since only a very limited number of samples can be used. Another method is based on the quadrature approximation [START_REF] Mcgraw | Description of aerosol dynamics by the quadrature method of moments[END_REF][START_REF] Fox | A quadrature-based third-order moment method for dilute gasparticle flow[END_REF][START_REF] Yuan | Conditional Quadrature Method of Moments for Kinetic Equations[END_REF] of the moments. Gaussian quadratures are usually used, leading to simple closures for the evolution in real and phase space. Although good results have been obtained for coalescing sprays [START_REF] Fox | Numerical simulation of spray coalescence in an Eulerian framework: Direct quadrature method of moments and multi-fluid method[END_REF], even with turbulence [START_REF] Belt | Quadrature method of moments for the pdf modeling of droplet coalescence in turbulent two-phase flow[END_REF], Fox et al. [START_REF] Fox | Numerical simulation of spray coalescence in an Eulerian framework: Direct quadrature method of moments and multi-fluid method[END_REF] demonstrate the difficulties of such methods to reproduce evaporation for continuous distributions, as they cannot recover point wise values of the NDF. Today's most promising formalism is the Multifluid Approach (MA, [START_REF] Laurent | Multi-fluid modeling of laminar poly-dispersed spray flames: origin, assumptions and comparison of sectional and sampling methods[END_REF]). Dividing the size space in sections [START_REF] Greenberg | On the origin of spray sectional conservation equations[END_REF], this formalism reproduces the evaporation of polydisperse sprays, and captures primary aspects of the size/velocity correlations in complex flows. Validation of this formalism was done through 0D and 1D simple test cases with one or two size moments per section [START_REF] Laurent | Numerical analysis of Eulerian multi-fluid models in the context of kinetic formulations for dilute evaporating sprays[END_REF][START_REF] Dufour | A second-order multi-fluid model for evaporationg sprays[END_REF], showing the robust mathematical basis.

Previous work on a real configuration using both the Lagrangian approach [START_REF] Senoner | Simulation aux grandes échelles de l'écoulement diphasique dans un bruleur aéronautique par une approche Euler-Lagrange[END_REF] and the MEF [START_REF] Lamarque | Schémas numériques et conditions limites pour la simulation aux grandes echelles de la combustion diphasique dans les foyers d'hélicoptre[END_REF][START_REF] Roux | Simulation aux grandes échelles d'un statoréacteur[END_REF][START_REF] Sanjosé | Évaluation de la méthode Euler-Euler pour la simulation numérique aux grandes échelles des chambres carburant liquide[END_REF] has shown that the size polydispersion does not greatly influence the mean velocity, but has significant effects on fluctuations and on the distribution of kerosene vapor. A good agreement was obtained between Lagrangian monodisperse and MEF results, showing the importance and good modeling of velocity dispersion. Introducing polydispersity in Lagrangian simulation however resulted in a significant deviation from monodisperse results.

Starting from these conclusions, the goal of this work is to introduce the size polydispersion in the Eulerian framework, by coupling MEF and MA approaches, to take into account the size dependent effect of turbulence on velocity dispersion of the spray dynamics, through size-conditioned closures. This new formalism is called the Multifluid Mesoscopic Eulerian Formalism (MMEF), and is implemented in the AVBP code to simulate the MERCATO configuration [START_REF] Garcia-Rosa | Phénomènes d'allumage d'un foyer de turbomachine en conditions de haute altitude[END_REF].

This paper is organized as follows. First the MMEF is introduced. Emphasis is made on the differences between strategies to solve the WBE. In section 3, the aeronautical configuration MERCATO is presented with the associated operating point and numerical methods. In section 4, a 0D evaporation case and a 2D vortex case are proposed to analyze the importance of polydispersity in the MERCATO conditions on similar reduced configurations. Finally, results obtained with the Lagrangian approach, the MEF and the MMEF are compared to experiment, with regards to mean and RMS velocities, as well as size-conditioned statistics.

Modeling approach

The modeling approach considered here aims at the description of a spray, i.e. a dispersed liquid phase resulting from the disintegration of a liquid core, and is restricted by the following assumptions:

droplets are spherical and smaller than Kolmogorov length scale, -the spray is dilute, -the liquid-gas density ratio is high.

These assumptions lead to a point-particle description without collisions and with drag force as the only external force acting on droplets (effect of gravity is neglected).

In this framework, the goal is to predict the NDF f p of the droplets, driven by the Williams-Boltzmann Equation (WBE):

∂f p ∂t + ∂c p,j f p ∂x j + 1 m p ∂F p,j f p ∂c p,j + ∂E S f p ∂S p + ∂E T f p ∂T p = 0 ( 1
)
where t, x j , c p,j , m p , S p , and T p are respectively the time, position, velocity, mass, surface and temperature of the droplet. The quantities F p,j , E S , E T are respectively the drag force and the mass and energy transfer due to evaporation:

F p,j = c p,j -u g@p,j τ p (2)
E S = 4 √ π ρ l S 3/2 p ṁp (3)
E T = φ c p m p C p,l . (4)
where u g@p,j , τ p , ṁp and φ c p are respectively the gas velocity at droplet position, the droplet relaxation time, the mass rate due to evaporation, and the heat transfer rate due to conduction. The evaporation model is derived for an isolated spherical droplet at the thermodynamic equilibrium with a surrounding quiescent mixture of ideal gas. Within the assumption of infinitely fast thermal conductivity in the liquid, the evaporating rate ṁp and the heat transfer rate φ c p are [START_REF] Sirignano | Fluid dynamics and transport of droplets and sprays[END_REF][START_REF] Kuo | Principles of Combustion[END_REF]:

ṁp = - p πS p Sh[ρD f ] ln(1 + B M ) (5)
φ c p = -φ c g + ṁp L v (T p ) = p πS p Nuλ(T ∞ -T p ) 1 + B T B T + ṁp L v (T p ) (6)
where L v (T p ) is the latent heat of vaporization at temperature T p , B M and B T stand for the Spalding numbers of mass and temperature respectively, accounting for the gas mixture modification in the gas layer surrouding the droplet, following the 1/3 rd -2/3 rd rule [START_REF] Hubbard | Droplet evaporation: effects of transient and variable properties[END_REF]. The Nusselt Nu and Sherwood Sh numbers are modified following [START_REF] Frössling | Über die Verdunstung fallender Tropfen (On the evaporation of falling drops)[END_REF] to account for the effect of non-vanishing slip velocity between droplet and carrier phase.

In the LES context, the Favre (.) and Reynolds (.) filtering procedures are introduced:

mΦ c f Φ (x) = Z ∆ f m Φ f Φ (x)G ∆ f (x -x ′ )dx ′ (7) fΦ (x) = Z ∆ f f Φ (x)G ∆ f (x -x ′ )dx ′ (8)

Lagrangian approach

A first way to solve Eq. 1 is to build statistics from a sample of individual particles. The Lagrangian equations associated to Eq. 1 are [START_REF] Crowe | Multiphase Flows with Droplets and Particles[END_REF]:

dx i dt = c p,i (9)
dc p,i dt = F p,i (10)
dS p dt = E S ( 11
)
dT p dt = E T (12)
LES filtering of Eqs. 9-12 leads to unclosed terms, due to the non-linear dependence of the right hand side on the gas phase. These terms can be accounted for like in [START_REF] Patel | Simulation of spray-turbulence-flame interactions in a lean direct injection combustor[END_REF][START_REF] Pozorski | Filtered particle tracking in isotropic turbulence and stochastic modeling of subgrid-scale dispersion[END_REF], but they are often neglected as in [START_REF] Apte | Large-eddy simulation of evaporating spray in a coaxial combustor[END_REF][START_REF] Senoner | Simulation aux grandes échelles de l'écoulement diphasique dans un bruleur aéronautique par une approche Euler-Lagrange[END_REF], so that all the source terms are written from filtered quantities. This type of simulation solves only one realisation of the disperse phase, and thus suffers from a slow statistical convergence. To overcome this problem, Eulerian methods can be envisaged, as they solve the NDF directly through its moments.

Mesoscopic Eulerian Formalism

The MEF developed by Février et al. [START_REF] Février | Partitioning of particle velocities in gas-solid turbulent flow into a continuous field and a spatially uncorrelated random distribution: theoretical formalism and numerical study[END_REF] is designed to solve turbulent two phase flows with one-way coupling, i.e. the gas phase is not influenced by the disperse phase. This formalism is an extension of the Bi-fluid model [START_REF] Drew | Theory of multicomponent fluids[END_REF], and accounts for the effect of turbulence on the disperse phase by solving the Random Uncorrelated Motion (RUM), which plays a similar role than pressure and viscosity in gas dynamics. More details about the derivation and the main assumptions are available in [START_REF] Février | Partitioning of particle velocities in gas-solid turbulent flow into a continuous field and a spatially uncorrelated random distribution: theoretical formalism and numerical study[END_REF].

Instead of solving for Eq. 1 on the NDF f p , one can solve for the moments, by integrating f p over velocities, size and temperature:

n l = Z f p dc p dS p dT p (13)
m l = Z m p f p dc p dS p dT p ( 14
)
u l,i = 1 n l Z c p,i f p dc p dS p dT p ( 15
)
δΘ l = 1 2n l Z (c p,i -U l,i ) 2 f p dc p dS p dT p (16)
H l = 1 n l Z C p,l T p f p dc p dS p dT p ( 17
)
where n l , m l , U l,i , δΘ l , and H l are respectively the number density, mass, velocity, Random Uncorrelated Energy (RUE), and enthalpy of the liquid phase at position x and time t. The RUE is the key-point of MEF, as it allows to capture velocity dispersion of the spray for moderate Stokes number, and represents most of the modeling effort. The MEF is monodisperse, i.e. there is only one droplet surface S l at a position x and a time t:

S l (x, t) = " m l n l 6 √ π ρ l « 2/3 (18)
Introducing the moments in the WBE and applying the LES filter to the resulting equations lead to LES filtered equations [START_REF] Moreau | Development of gas-particle eulereuler LES approach: a priori analysis of particle sub-grid models in homogeneous isotropic turbulence[END_REF]:

∂ ∂t n l + ∂ ∂x m n l b u l,m = 0 (19) ∂ ∂t m l + ∂ ∂x m m l b u l,m = E(1) (20)
∂ ∂t m l b u l,i + ∂ ∂x m m l b u l,i b u l,m = E(b u l,i ) + F d - ∂ ∂x m m l ( c δR l,im + τ l,im ) (21) ∂ ∂t m l δ b Θ l + ∂ ∂x m m l δ b Θ l b u l,m = E(δ b Θ l ) + W d + Π + Π SGS - 1 2 ∂ ∂x m m l ( c δS l,iim ) (22)
∂ ∂t m l b H l + ∂ ∂x m m l b H l b u l,m = E( b H l ) + b φ l (23)
where E(Φ) = ṁl Φ are the evaporation source terms, F d is the drag force, δR l is the uncorrelated stress tensor, τ l is the subgrid stress tensor, δS l is the RUE diffusion tensor, W d is the RUE source term due to drag force, Π and Π SGS are the RUE production due to resolved and subgrid motions, and φ l is the heat transfer.

Closures

The drag force F d and its RUM contribution W d are:

F d = -m l b u l,i -b u g,i τ p (24)
W d = -2 m l τ p δ b Θ l ( 25
)
where the relaxation time τ p is determined by the droplet surface S l . The second order tensor c δR l,im is decomposed into isotropic and deviatoric parts:

c δR l,ij = c δR l,ii 3 δ ij + c δR * l,ij = P RU M δ ij + c δR * l,ij (26)
where

P RU M = 2/3 c δΘ l .
The deviatoric part is closed using a viscous assumption [START_REF] Kaufmann | Comparison between Lagrangian and mesoscopic Eulerian modelling approaches for inertial particles suspended in decaying isotropic turbulence[END_REF]:

c δR * l,ij = -ν RU M " ∂ b u l,i ∂x j + ∂ b u l,j ∂x i - 2 3
∂ b u l,k ∂x k δ ij « ( 27
)
where ν RU M = τ p c δΘ l /3. The third order tensor c δS l,iij is closed using gradient diffusion law:

c δS p,iij = -κ RU M ∂ c δΘ l ∂x j ( 28
)
where κ RU M = 10/27τ p δ b Θ l . Finally, the RUE production Π is modeled as:

Π = -m l c δR l,im ∂ b u l,i ∂x m (29)
Following [START_REF] Moreau | Development of gas-particle eulereuler LES approach: a priori analysis of particle sub-grid models in homogeneous isotropic turbulence[END_REF], the subgrid tensor T l,ij is closed using a Smagorinsky model for the deviatoric part and a Yoshizawa model for the isotropic part:

τ l,ij = -C 2 s ∆ 2 f | Ŝl | " Ŝl,ij - δ ij 3 Ŝl,kk « + C 2 l ∆ 2 f | Ŝl | 2 δ ij ( 30
)
where b S l,ij = ∂ b u l,i ∂x j + ∂ b u l,j
∂x i , C l = 0.116 and C s = 0.16 are constants, and ∆ f is the grid size. Subgrid diffusion term is neglected in Eq. 22, whereas the subgrid production term considers that the subgrid scale motion is totally transfered to the RUM:

Π SGS = -m l T l,im ∂ b u l,i ∂x m (31)
The evaporation rate ṁl and heat transfer φ l are calculated using the droplet diameter, and taking the gas properties at infinity as the local gas properties:

ṁl = - p πS l Sh[ρD f ] ln(1 + B M ) (32)
φ l = p πS l Nuλ(T g -T l ) 1 + B T B T + ṁl L v (T l ) ( 33
)
where T l = H l /m l C p,l is the liquid phase temperature. Filtered evaporation terms are calculated as E( b Φ) = ṁl b Φ.

Multifluid Mesoscopic Eulerian Formalism

To extend the MEF to polydisperse sprays, it is proposed to couple the MEF to the MA of Laurent et al. [START_REF] Laurent | Multi-fluid modeling of laminar poly-dispersed spray flames: origin, assumptions and comparison of sectional and sampling methods[END_REF]. The derivation of this new formalism called Multifluid Mesoscopic Eulerian Formalism (MMEF) is further detailed in [START_REF] Vié | Simulation aux grandes échelles d'écoulements diphasiques turbulents phase liquide dispersée[END_REF] and is presented in the following. The main idea is to divide the size space in

N intervals Ω k = [S (k) p , S (k+1) p
] of fixed length, called sections [START_REF] Greenberg | On the origin of spray sectional conservation equations[END_REF]. Each section is calculated with its own set of equations, and interacts with the others sections through droplet size changes only, for example due to evaporation. Moments are then now defined for each section k:

n (k) l = Z Ω k f p dc p dS p dT p (34)
m (k) l = Z Ω k m p f p dc p dS p dT p (35)
u (k) l,i = 1 m (k) l Z Ω k m p c p,i f p dc p dS p dT p (36)
δΘ

(k) l = 1 2m (k) l Z Ω k m p (c p,i -U (k) l,i )(c p,i -U (k) l,i )f p dc p dS p dT p (37)
H (k) l = 1 m (k) l Z Ω k m p C p,l T p f p dc p dS p dT p ( 38
)
where superscript . (k) qualifies quantities integrated over the section [S

(k) p , S (k+1) p
]. The resulting LES filtered equations for the MMEF are similar to MEF equations:

∂ ∂t m (k) l + ∂ ∂x m m (k) l b u (k) l,m = b E (k) (1) (39)
∂ ∂t m (k) l b u (k) l,i + ∂ ∂x m m (k) l b u (k) l,i b u (k) l,m = b E (k) (u (k) l,i ) + F (k) d,i - ∂ ∂x m (m (k) l c δR (k) l,im + τ (k) l,im ) (40)
∂ ∂t m (k) l c δΘ (k) l + ∂ ∂x m m (k) l c δΘ (k) l b u (k) l,m = b E (k) (δΘ (k) l ) + W (k) d + Π (k) + Π (k) SGS - 1 2 ∂ ∂x m m (k) l ( c δS (k) l,iim ) (41)
∂ ∂t m (k) l b H (k) l + ∂ ∂x m m (k) l b H (k) l b u (k) l,m = b E (k) (H (k) l ) + b φ (k) l (42)
with the same definitions for source terms as in 2.2, but restricted to each section k.

Closures

All closures are written as for the MEF approach, applying the definitions to the sectional moments quantities. Only the evaporation terms take a different form, as they imply sectional interactions. These terms are closed using the first order multifluid approach [START_REF] Laurent | Multi-fluid modeling of laminar poly-dispersed spray flames: origin, assumptions and comparison of sectional and sampling methods[END_REF]:

E (k) (Φ (k) ) = E (k+1) 1 m (k+1) l Φ (k+1) -E (k) 1 m (k) l Φ (k) -E (k) 2 m (k) l Φ (k) (43)
where

E (k) 1 = 2E (k) S S (k) /(S (k+1) 2 -S (k) 2
) corresponds to the droplet flux from section k to section k -1 due to surface decrease and

E (k) 2 = 3E (k)
S /(S (k+1) +S (k) ) corresponds to the droplet flux from section k to the gas phase due to droplet mass loss, and where E (k) S is the mass transfer rate defined in Eq. 1. For the RUE equation, an additional term is needed to account for the modification of the mean velocity induced by the evaporative fluxes:

E (k) ( c δΘ (k) l ) = E (k+1) 1 m (k+1) l c δΘ (k+1) l -E (k) 1 m (k) l c δΘ (k) l -E (k) 2 m (k) l c δΘ (k) l +E (k+1) 1 m (k+1) l (b u (k+1) l,i -b u (k) l,i ) 2 2 ( 44
)
The MA accounts for the effect of polydispersion on evaporation and drag force. The MMEF goes further in the modeling, because it considers the effect of polydispersion on the turbulent motion by the resolved scales through the sizeconditioned random uncorrelated energy. As this motion is directly related to the Stokes number, a more precise description of the overall dynamics is expected.

Equations for the gas phase

Knowing now the coupling terms between the liquid and gas phases, the LES equations of the gas phase are introduced. Indeed, all models considered here are in a two-way coupling framework in terms of drag and evaporation. As stated in [START_REF] Février | Partitioning of particle velocities in gas-solid turbulent flow into a continuous field and a spatially uncorrelated random distribution: theoretical formalism and numerical study[END_REF], a rigorous derivation of the moment equations would have to consider the effect of two-way coupling at the kinetic level. However, in the case of small Stokes number relatively to the smallest resolved scales, Février et al. [START_REF] Février | Partitioning of particle velocities in gas-solid turbulent flow into a continuous field and a spatially uncorrelated random distribution: theoretical formalism and numerical study[END_REF] suggest that MEF stays valid for a two-way coupling. Gas phase equations are then similar to pure gaseous flows, with additional source terms resulting from exchanges with the liquid:

∂ ∂t ρ g + ∂ ∂x m ρ g b u g,m = Φ ρ g (45) ∂ ∂t ρ g b u g,i + ∂ ∂x m ρ g b u g,i b u g,m = Φ ρ g u g,i - ∂ b P ∂x i + ∂ ∂x m " ρ g ν g b D g,im + τ g,ij " (46)
∂ ∂t ρ g b E g + ∂ ∂x m (ρ g b E g + b P )b u g,m = Φ ρ g E g + ∂ ∂x m " ρ g ν g b D g,im b u g,i + Q g,ij " (47)
where ρ g , u g,i , E g , P and ν g are the density, velocity, total energy, pressure and kinetic viscosity of the gas phase, b

D g,ij = ∂ b u l,i ∂x j + ∂ b u l,i ∂x j -2 3 ∂ b u l,k
∂x k δ ij is the strain tensor, τ g,ij is the subgrid scale stress tensor modeled with the WALE approach [START_REF] Nicoud | Subgrid-scale stress modelling based on the square of the velocity gradient[END_REF], and Q g,ij is the filtered heat flux. Φ ρ g , Φ ρ g u g,i and Φ ρ g E g are the coupling terms between liquid and gas phases. In the MEF, the coupling terms are:

Φ ρ g = ṁl ( 48
)
Φ ρ g u g,i = ṁl b u l,i -F d (49)
Φ ρ g E g = -F d,i u l,i + ṁl " b H l + 1 2 b u 2 l,i « -φ l (50)
whereas in the MMEF, the coupling terms are the sum of the contributions of each section:

Φ ρ g = N X k=1 E (k) 2 m (k) l (51)
Φ ρ g u g,i = N X k=1 h -F (k) d,i + E (k) 2 m (k) l b u (k) l,i i (52)
Φ ρ g E g = N X k=1 » -F (k) d,i b u l,i + E (k) 2 m (k) l " b H (k) l + 1 2 b u (k)2 l,i « -φ (k) l - (53)
In the Lagrangian approach, the coupling terms are obtained by simple projection on the gas flow grid [START_REF] Senoner | Simulation aux grandes échelles de l'écoulement diphasique dans un bruleur aéronautique par une approche Euler-Lagrange[END_REF]:

Φ ρ g = - N p X p=1 ṁp δ(x -x p ) ( 54
)
Φ ρ g u g,i = - N p X p=1 [ ṁp u p,i + F p,i ] δ(x -x p ) (55)
Φ ρ g E g = - N p X p=1 » F p,i u p,i + 1 2 ṁp u 2 p,i -φ c p - δ(x -x p ) ( 56
)
where N p is the number of Lagrangian particles.

3 The aeronautical-type configuration and numerical methods

The configuration

The MERCATO test rig is an experimental swirl combustor (Fig. 1) operated at ONERA (Fauga), for which measurements are available [START_REF] Garcia-Rosa | Phénomènes d'allumage d'un foyer de turbomachine en conditions de haute altitude[END_REF]. Air is injected inside the plenum and enters the combustion chamber through the swirler. The liquid injection is located at the end of the swirler stage, and at the center of the swirled gaseous jet. The studied case is a two-phase evaporating flow without combustion.

The injection system is a pressure-swirl atomizer. It generates a hollow cone spray with a half angle of 40 • and orthoradial motion, which produces a polydisperse cloud of droplets. The droplet size distribution at 13mm downstream the injector is shown in Fig. 3. The injected fuel is kerosene, modeled by a a surrogate [START_REF] Sanjosé | Fuel Injection Model for Euler-Euler and Euler-Lagrange Large-Eddy Simulations of an evaporating spray inside an aeronautical combustor[END_REF] in the simulations.

The liquid flow rate is 1 g/s at a temperature of 300K, and the gas flow rate is 15 g/s at a temperature of 463K.

Numerical methods

Gas phase

The AVBP code solves compressible Navier-Stokes equation for reactive flows on unstructured grids in a cell-vertex formulation. To perform Large Eddy Simula-Fig. 1: Scheme of the MERCATO test rig. tions, a 3 rd order in time and space finite element scheme is used [START_REF] Rudgyard | Development of high-order Taylor-Galerkin schemes for unsteady calculations[END_REF]. To handle the oscillatory behavior of such a scheme, artificial viscosity [START_REF] Von Neumann | A method for the numerical calculation of hydrodynamic shocks[END_REF] is added where the solution exhibits strong gradients. Non-reflecting NSCBC boundary conditions [START_REF] Poinsot | Boundary conditions for direct simulations of compressible reacting flows[END_REF] and wall law [START_REF] Nicoud | Subgrid-scale stress modelling based on the square of the velocity gradient[END_REF] associated to the LES Smagorinsky model are used.

Lagrangian disperse phase

The Lagrangian equations are solved in AVBP using a first order Euler explicit time integration. Considering that the computation is driven by the acoustic time step, a first order in time scheme is sufficient. The interpolation of gas properties at particle positions is made by a linear least square approximation.

Eulerian disperse phase

The Eulerian liquid phase solver of AVBP is similar to the gas phase solver in terms of equations, but it behaves like a highly compressible flow, so that vacuum and strong gradients occur, difficult to handle numerically. A classical way to handle such numerical problems is to use artificial viscosity, associated to sensors that detect the difficulties. This methodology has been applied to the simulation of a decaying Homogeneous Isotropic Turbulence [START_REF] Vié | Validation of the Eulerian mesoscopic approach in particle-charged homogeneous isotropic decaying turbulence in the scope of large eddy simulation of fuel sprays[END_REF], and has proven to be very accurate. However more adapted schemes are needed in order to improve the robustness of the solver, like kinetic schemes [START_REF] De Chaisemartin | Modèles eulériens et simulation de la dispersion turbulente de brouillards qui s'évaporent[END_REF], but their derivation in unstructured cell vertex 3-dimensional formulation is not straightforward and requires further developments.

Boundary conditions for the disperse phase

Injection is a keypoint in the simulation of two-phase flow burners. Here the injection system is modeled by the FIMUR model [START_REF] Sanjosé | Fuel Injection Model for Euler-Euler and Euler-Lagrange Large-Eddy Simulations of an evaporating spray inside an aeronautical combustor[END_REF], which using autosimilarity assumptions. Fig. 2 shows a sketch of this model [START_REF] Sanjosé | Fuel Injection Model for Euler-Euler and Euler-Lagrange Large-Eddy Simulations of an evaporating spray inside an aeronautical combustor[END_REF]. The model gives the following Fig. 2: Definition of parameters and sketch of the FIMUR methodology [START_REF] Sanjosé | Fuel Injection Model for Euler-Euler and Euler-Lagrange Large-Eddy Simulations of an evaporating spray inside an aeronautical combustor[END_REF]. injection condition in radial coordinates at fixed x position:

u i l,r (θ, r) = ṁl ρ l A p v u u t 1 - R 0 S R i S ! 2 r R i S ( 57
)
u i l,θ (θ, r) = ṁl ρ l A p R 0 S R i S r R i S ( 58
)
u i l,x (θ, r) = ṁl πR 2 0 A u (59)
m i l (θ, r) = πR 2 0 I α A u exp -(r -R i S ) 2 σ 2 ! ( 60
)
where R 0 S , R i S are the middle position of the liquid sheet at the injector exit and the liquid spray at the numerical injection location, R 0 the radius of the injection hole, A p the injection surface, A u a parameter which accounts for contraction effects on the discharge orifice as well as air entrainement effects along the axial distance and σ the variance of the Gaussian profile of liquid mass. This model was initially designed for monodisperse injection, and its extension to multifluid injection is made by distributing the liquid volume fraction and droplet number density over all diameters, following the size distribution of Fig. 3. First the number density f (p) num of each section is calculated by integrating the experimental size distribution f exp over each interval:

f (k) num = Z S (k+1) S (k) f exp (S p )dS p (61)
Then the mass density of each section is obtained as:

f (k) mass = f (p) num ρ l 24 √ π (S (k+1)2 -S (k)2 ) ( √ S (k+1) - √ S (k) ) (62)
The total mass predicted by the monodisperse FIMUR model is then distributed using f

(k) mass : m (k),i l = f (k) mass m i l ( 63
)
The velocity of each section is equal to the monodisperse velocity, in order to keep the same momentum for multifluid injection and monodisperse injection:

u (k),i l,r = u i l,r , u (k),i l,θ = u i l,θ , u (k),i l,x = u i l,x (64)
X u (k),i l,x m (p) l = u i l,x X m (p) l = u l,x m l (65)
This assumption allows to compare monodisperse and polydisperse results with an equivalent injection condition in terms of velocity. To obtain a more accurate injection condition, it would be desirable to adapt the FIMUR model to a polydisperse spray. This will be a future work on this model.

Importance of the polydispersion for the MERCATO configuration on reduced test configurations

In this section, the importance of the polydispersion in the MERCATO conditions is investigated on two reduced configurations: a 0D evaporation case and a 2D frozen vortex. The former is expected to show the influence of the polydispersion on the evaporation time and the mean properties of the liquid phase, and the latter is expected to show its influence on the spatial dispersion induced by the drag-forced size-conditionned dynamics. The parameters of the two cases are set to represent the conditions in the MERCATO test rig.

0D Evaporation

The initial temperature is 300K for the droplets and 463K for the gas phase. For the MEF simulations, the initial diameter is 55µm. For the MMEF, the NDF is the one of Fig. 3. Because of transient heating, no analytical solution is available.

Figure 4 shows the time evolution of the total mass and mean temperature for MEF and MMEF. First, one can notice the important difference on the evaporation time for the two formalisms. Where the MEF exhibits an evaporation time of 50ms, the MMEF predicts an evaporation time of 400ms. This difference appears also on the heating time (time to reach the maximum temperature), for which the MEF predicts 25ms whereas the MMEF predicts 150ms. Those two effects are directly related to the polydispersion, as the MEF cannot capture the large evaporation and heating times of big droplets.

2D Vortex

It is proposed here to analyse the spatial dispersion effects of a 2D vortex having the same characteristics than the swirl motion in MERCATO at 13 mm. The test case is a two-dimensional frozen vortex, with a spatially homogeneous initial droplet distribution having the same velocity as the gas phase. The radial velocity is zero and the tangential velocity of the gas phase corresponds to an isotropic vortex:

U g,θ (r) = Γ r r 2 v exp - r 2 2r 2 v ! (66)
where Γ = 0.6638 is the vortex strength, and r v = 0.012 its radius. Figure 5 shows that Eq. 66 reproduces well the velocity profile of MERCATO. As the edge of the vortex has a non-zero radial velocity, Eq. 66 is used only for r < 0.02. For r ≥ 0.02, a linear interpolation of the true MERCATO profile is used. The domain of [-0.12, 0.12] 2 m is discretized with 800 2 quad cells. Ten sections are used for the MMEF, with the same droplet number in each section. The computation is run for 5ms, which is the characteristic time needed by the injected spray to impact the wall in the MERCATO configuration.

On Fig. 6, the total droplet number density obtained at 5 ms is shown for MEF and MMEF. The high vorticity of the central zone ejects the droplets with a force depending on their sizes. With the MEF approach, an accumulation zone is generated, which corresponds to the accumulation of droplets in the low vorticity zone, latter convected due to the non-zero velocity at the edge of the vortex. With the MMEF approach, one can see several accumulation zones, each corresponding to one section, with its proper Stokes number and so its proper dynamics. Figure 7 shows the number density and mean diameter radial profiles for both MEF and MMEF at 5 ms. While the droplet number density is widely distributed spatially in MMEF, it is concentrated around 0.038m in MEF. Due do the weaker influence of the gas phase, bigger droplets keep their initial velocity out of the vortex, and so their initial radial deviation. On the contrary, small droplets are in equilibrium with the gas phase, and reach a small velocity out of the vortex. This impacts the mean diameter profile, with small diameters near the vortex center and bigger diameters far from it. Two difficulties are shown in Fig. 7. The first one is linked to the mean diameter at the center of the vortex, where the absence of droplets should lead to a zero mean diameter. However the numerical methods impose to keep a non-zero droplet number density everywhere, even in empty zones where it is set to a value several orders of magnitude lower than the mean value, leading to a non-zero mean diameter in these empty zones. The second difficulty concerns the discretization in sections, resulting in discrete accumulations corresponding to each section trajectories. However in practical cases, it is expected that turbulent mixing and species diffusion will considerably reduce this effect.

Results on the MERCATO configuration

Here are compared the Lagrangian results of [START_REF] Senoner | Simulation aux grandes échelles de l'écoulement diphasique dans un bruleur aéronautique par une approche Euler-Lagrange[END_REF], the MEF results of [START_REF] Sanjosé | Fuel Injection Model for Euler-Euler and Euler-Lagrange Large-Eddy Simulations of an evaporating spray inside an aeronautical combustor[END_REF] and the present MMEF results. For the MMEF approach, the mesh is composed of 1, 351, 767 tetrahedrons (Fig. 8), while 3, 527, 637 tetrahedrons were used for the MEF approach and 3, 934, 364 tetrahedrons were used for the Lagrangian. The main difference between these meshes concerns the region far from the injection (after 56mm downstream the injection), all meshes having the same level of refinement in the region close to the injection location. First the gas phase is validated on a pure gaseous computation (no disperse phase), verifying that the conditioning flow field for the dispersed phase is well resolved. Then liquid phase results are analyzed, qualitatively in terms of liquid volume fraction and gas mass fraction, and quantitatively in terms of mean and RMS velocity at 3 positions in the chamber. Size-conditioned statistics are all also investigated at 8 locations in the spray. Fig. 8: Cut plane of the mesh of the combustion chamber Fig. 9: Cut plane of the mean velocity field for the gas phase with isolines of zero axial velocity.

Gas phase

Mean quantities for the gas phase are obtained after averaging over 800 ms. Figure 9 shows a cut plane of the mean axial velocity, with isolines of zero axial velocity. The flow exhibits two recirculation zones : the Central Toroidal Recirculation Zone (CTRZ) and the Corner Recirculation Zone (CRZ). One can also notice the high shear stress zone induced by the gas injection in the combustion chamber, which impacts on the wall, and follows it until the exit of the combustion chamber.

Figures 101112131415show mean and RMS velocites for axial, radial and orthoradial velocity components. Globally, numerical results are in good agreement with experiment for the first three measurement positions (6, 26 and 56 mm), where liquid phae measurements were made, whereas the RMS and velocities are underestimated far from the injection location (86mm and 116mm), due to the coarser mesh. The present results are considered sufficiently accurate to allow the spray modeling study. A more detailed characterization of the gas phase is available in [START_REF] Sanjosé | Évaluation de la méthode Euler-Euler pour la simulation numérique aux grandes échelles des chambres carburant liquide[END_REF] and [START_REF] Senoner | Simulation aux grandes échelles de l'écoulement diphasique dans un bruleur aéronautique par une approche Euler-Lagrange[END_REF]. Fig. 15: Orthoradial RMS velocity profiles for the gas phase at 6, 26, 56, 86 and 116mm locations (Fig. 8): experiment (square) and simulation (full line) 5.2 Liquid phase

Mean flow

The results of EL [START_REF] Senoner | Simulation aux grandes échelles de l'écoulement diphasique dans un bruleur aéronautique par une approche Euler-Lagrange[END_REF], MEF [START_REF] Sanjosé | Évaluation de la méthode Euler-Euler pour la simulation numérique aux grandes échelles des chambres carburant liquide[END_REF][START_REF] Sanjosé | Fuel Injection Model for Euler-Euler and Euler-Lagrange Large-Eddy Simulations of an evaporating spray inside an aeronautical combustor[END_REF] and MMEF are confronted to experiment. Figure 16 shows the mean liquid volume fraction for the three simulations. The three approaches exhibit similar distributions. The liquid spray impacts the wall at nearly the same position. Droplets are ejected by the CTRZ, leading to a zone empty of droplets in the center of the chamber. As expected, Eulerian cases appear more diffusive than the Lagrangian case. The major difference is linked to polydispersion: in EL and MMEF, large droplets persist in the whole chamber, while in MEF (monodisperse), the mean droplet lifetime is too short to allow them to reach the downstream zone.

Figures 171819show mean velocity profiles at three locations for the three simulations and the experiment. All results are globally in good agreement with experiment. Some discrepancies are observed in the central zone for the axial velocity and at the spray boundary for the radial velocity. The two polydisperse approaches (MMEF and EL) misestimate the axial velocity in the central zone, whereas the MEF reproduces it accurately. On the contrary, MEF predicts non-zero values near the spray boundary at z = 26mm, where MMEF and EL show values close to zero, in agreement with experiment. However discrepancies are small and all these results confirm that the size distribution has a weak effect on the mean velocity, as mentioned in [START_REF] Senoner | Simulation aux grandes échelles de l'écoulement diphasique dans un bruleur aéronautique par une approche Euler-Lagrange[END_REF].

Figures [START_REF] Greenberg | On the origin of spray sectional conservation equations[END_REF][START_REF] Hubbard | Droplet evaporation: effects of transient and variable properties[END_REF][START_REF] Jaegle | Large eddy simulation of evaporating sprays in complex geometries using Eulerian and Lagrangian methods[END_REF] shows RMS velocity profiles at the same locations. Here polydispersion has a more important effect. All models are similar and good close to the injection, but MEF underestimates the RMS velocities at further locations, while the two polydisperse approaches (MMEF and EL) give better predictions even downstream. Better results are obtained by the MMEF compared to EL, except for the orthoradial component in the last plane. This is possibly due to the way the RMS velocity is estimated. Indeed the RMS velocity is composed of the fluctuating part of the resolved field, the subgrid velocity and the random uncorrelated velocity. The latter is evaluated assuming isotropy, which is not true in shear flows, as shown in [START_REF] Masi | The mesoscopic Eulerian approach for evaporating droplets interacting with turbulent flows[END_REF]. A better estimation of RMS velocities in MEF and MMEF would be achieved by defining the orientation of the stress tensor, but this is not investigated here. Table 1: Probe location for size-conditioned statistics.

In Fig. 24, MMEF and EL size distributions are compared to experimental results. At 26mm and 56mm, EL and MMEF are in good agreement with experiment. At 6mm, both approaches reproduce the size distribution at the boundary of the spray, but not on the axis and in the spray, where the distribution is shifted towards big sizes, with worse results for the MMEF. This may be due to the absence of secondary break-up model, as shown by Senoner et al. [START_REF] Senoner | Simulation aux grandes échelles de l'écoulement diphasique dans un bruleur aéronautique par une approche Euler-Lagrange[END_REF]. One can also notice that MMEF overestimates the probability of big droplets, which may be due to the absence of secondary break-up, which is suppose to break these big droplets, but can also be related to the phase space diffusion occuring in such a first order Multifluid approach, and resulting in overestimated mean sizes [START_REF] Laurent | Numerical analysis of Eulerian multi-fluid models in the context of kinetic formulations for dilute evaporating sprays[END_REF].

Figures 252627show size-conditioned velocities. The axial component is relatively well captured by both approaches. Some discrepancy is observed for the velocities of small droplets at 6mm, on the axis, where both simulations predict a high negative velocity, whereas the experiment show a positive velocity. Considering that small droplets have a strong tendency to follow the gas phase, and that the gas phase velocity is about -15m/s at this location, the experimental results are surprising. A possible explanation of such a behavior is that the liquid injection strongly modifies the gas flow close to the injector, due to a dense and coherent liquid core, shifting the CTRZ downstream and thus leading to positive velocities at the first measurement location. As we inject a dispersed phase directly, the impact will be limited compared to a dense structure, and a modification of the FIMUR model may be needed in this perspective. For radial and orthoradial components, EL and MMEF are in good agreement with experiment. Largest differences are found at 6mm on the axis, where MMEF predicts a quasi-zero velocity for those two components, while EL shows high velocities of small droplets. These two results differ from the experimental results, and may also be explained by a lack of the injection strategy. It is interesting to notice that experiment and EL computations suffer from a poor statistical convergence on the axis where the droplet flux is low, whereas the statistical description of both eulerian formalisms (MEF and MMEF) achieves converged statistics.

The axial velocity conditioned by size (Fig. 25) also gives informations about the spatial dispersion of droplets. Basically, velocities on the axis (A1-A3) shows that droplets smaller than 30µm are sensitive to the CTRZ, as they always have a negative axial velocity, contrary to droplets bigger than 30µm, which have a positive axial velocity. Furthermore, the axial velocities at the point W3 show that droplets smaller than 30µm are also captured by the CRZ.

Gaseous kerosene mass fraction

On Fig. 28, instantaneous gaseous kerosene mass fraction for MEF, MMEF and EL are compared. No experimental results are available but the three simulations show a great influence of the polydispersity on the kerosene distribution. In particular, MEF predicts a high mass fraction close to the injection, whereas the two polydisperse approaches lead to a higher mass fraction downstream the combustion chamber. This is again due to the shorter life time of the mean size droplets in the monodisperse approach, leading to the complete disappearance of the liquid phase at 60mm downstream the injection location. These results in conjunction with the liquid phase distribution will lead to a different behavior in terms of ignition and combustion [START_REF] Versaevel | Combustion laminaire diphasique: Etude théorique et expérimentale[END_REF][START_REF] Boileau | Simulations aux grandes échelles de l'allumage diphasique des foyers aéronautiques[END_REF].

Conclusions

Accounting for polydispersion in Eulerian-Eulerian simulation of turbulent twophase flows has been proposed by combining the mesoscopic Eulerian formalism [START_REF] Février | Partitioning of particle velocities in gas-solid turbulent flow into a continuous field and a spatially uncorrelated random distribution: theoretical formalism and numerical study[END_REF] and the Multifluid Approach [START_REF] Laurent | Multi-fluid modeling of laminar poly-dispersed spray flames: origin, assumptions and comparison of sectional and sampling methods[END_REF]. The resulting Multifluid Mesoscopic Eulerian Formalism is able to capture size-conditioned statistics on drag, turbulence and evaporation. This formalism has been integrated in the AVBP solver and applied to the simulation of an aeronautical test rig operated by ONERA and called MERCATO [START_REF] Garcia-Rosa | Phénomènes d'allumage d'un foyer de turbomachine en conditions de haute altitude[END_REF]. The gaseous flow is first validated and results are in good agreement with the experiment for the three first measurement planes.

Results on the two-phase flow are compared to experiment and to two other modeling approaches: the monodisperse mesoscopic Eulerian formalism and a polydisperse Eulerian-Lagrangian approach.

Mean liquid velocities are well captured by all approaches, confirming the weak influence of polydispersity on these quantities, as stated in [START_REF] Senoner | Simulation aux grandes échelles de l'écoulement diphasique dans un bruleur aéronautique par une approche Euler-Lagrange[END_REF]. RMS velocities are however significantly impacted by polydispersity, as shown by the good agreement obtained with the two polydisperse approaches (EL and MMEF), where the MEF tends to underestimate RMS velocities.

Results in terms of size-conditioned statistics show a good agreement of both approaches with experiment, even if some discrepancies appear close to the injection. The resulting field of gaseous kerosene mass fraction confirms the great influence of polydispersity, showing a totally different repartition between polydisperse and monodisperse approaches. This is critical for the prediction of ignition and combustion, and comparison with experiment is crucially needed here.

Finally, even if results obtained by EL and MMEF are similar, it is important to note that the computational time required by the 10 sections MMEF simulation is about six times higher than the EL simulation with 285000 particles (MEF and EL simulations require approximatively the same resources for this configuration [START_REF] Sanjosé | Évaluation de la méthode Euler-Euler pour la simulation numérique aux grandes échelles des chambres carburant liquide[END_REF]). This overcost could be reduced by using moment methods such as high order multifluid approaches [START_REF] Laurent | Numerical analysis of Eulerian multi-fluid models in the context of kinetic formulations for dilute evaporating sprays[END_REF][START_REF] Dufour | A second-order multi-fluid model for evaporationg sprays[END_REF]. An interesting alternative is to use moments methods like the Coupled Size-Velocity Moments method (CSVM) of Vié et al. [START_REF] Vié | Size-velocity correlations in high order moment methods for polydisperse evaporating sprays: modelling and numerical issues[END_REF]. This model is an extension of the Eulerian Multi-Size Moment method (EMSM) of Kah et al. [START_REF] Kah | Taking into account polydispersity in the framework of a coupled Euler-Lagrange approach for the modeling of liquid fuel injection in internal combustion engines[END_REF][START_REF] Massot | A robust moment method for evaluation of the disappearance rate of evaporating sprays[END_REF], which captures the evaporation process more accurately than the 10 sections multifluid approach with only 4 size moments. It enables to capture size-velocity correlations, which permits to reconstruct the velocity for each droplet size, and thus is of primary importance in real flows like the MERCATO configuration. Its adaptation to turbulent regimes and implementation in AVBP will be interesting for complex industrial simulations. lence, Vol. 10 (4), pp. 1-21, 2009 67. Zaichik, L., and Alipchenkov, V.M. , Statistical models for predicting particle dispersion and preferential concentration in turbulent flows, International Journal of Heat and Fluid Flow, Vol. 26, pp. 416-430, 2005

Fig. 3 :

 3 Fig. 3: MERCATO: Number Density Function at 13mm downstream the injector location (histogram) and Multifluid representation with 10 sections (solid line).

Fig. 4 :

 4 Fig. 4: Evaporation test case in MERCATO conditions: Evolution of total mass (upper left), droplet number density (upper right), liquid temperature (lower left) and mean diameter (lower right) for MMEF (full line) and MEF (dashed line).

Fig. 5 :

 5 Fig. 5: Gas velocity profile for the 2D vortex (solid line) and in MERCATO test case (13mm downstream the injection location)(squares).

Fig. 6 :

 6 Fig. 6: Droplet Number Density at t = 5ms for MMEF (left) and MEF(right).

Fig. 7 :

 7 Fig. 7: Radial distribution of droplet number density for MEF (top left) and MMEF (top right), and arithmetic mean diameter (bottom) for MEF (dashed line) and MMEF (solid line).

Fig. 10 :Fig. 11 :Fig. 12 :Fig. 13 :Fig. 14 :

 1011121314 Fig.10: Axial mean velocity profiles for the gas phase at 6, 26, 56, 86 and 116mm locations (Fig.8): experiment (square) and simulation (full line)

Fig. 16 :

 16 Fig. 16: Mean liquid volume fraction for MMEF (top), EL (middle) and MEF (bottom).

Fig. 17 :

 17 Fig.17: Axial mean velocity profiles for the liquid phase at 6, 26 and 56mm locations (Fig.8): experiment (square), MMEF (full line), MEF (dashed line)and EL (dot-dashed line).

Fig. 18 :

 18 Fig. 18: Radial mean velocity profiles for the liquid phase at 6, 26 and 56mm locations (Fig. 8): experiment (square), MMEF (full line), MEF (dashed line) and EL (dot-dashed line).

Fig. 19 :

 19 Fig. 19: Orthoradial mean velocity profiles for the liquid phase at 6, 26 and 56mm locations (Fig. 8): experiment (square), MMEF (full line), MEF (dashed line) and EL (dot-dashed line).

Fig. 20 :

 20 Fig. 20: Axial RMS velocity profiles for the liquid phase at 6, 26 and 56mm locations (Fig. 8): experiment (square), MMEF (full line), MEF (dashed line) and EL (dot-dashed line).

Fig. 21 :

 21 Fig. 21: Radial RMS velocity profiles for the liquid phase at 6, 26 and 56mm locations (Fig. 8): experiment (square), MMEF (full line), MEF (dashed line) and EL (dot-dashed line).

Fig. 22 :

 22 Fig. 22: Orthoradial RMS velocity profiles for the liquid phase at 6, 26 and 56mm locations (Fig. 8): experiment (square), MMEF (full line), MEF (dashed line) and EL (dot-dashed line).

Fig. 24 :

 24 Fig. 24: Size distributions: experiment (histograms), MMEF (full line) and EL (dashed line).

Fig. 25 :

 25 Fig. 25: Size-conditioned axial velocity: experiments (histograms), MMEF (full line) and EL (dashed line).

Fig. 26 :

 26 Fig. 26: Size-conditioned radial velocity: experiments (histograms), MMEF (full line) and EL (dashed line).

Fig. 27 :

 27 Fig. 27: Size-conditioned orthoradial velocity: experiments (histograms), MMEF (full line) and EL (dashed line).

Fig. 28 :

 28 Fig. 28: Instantaneous kerosene mass fraction for MMEF (top), EL (middle, image from [56], black lines represent iso-lines of mass transfer) and MEF (bottom).

Acknowledgments

We would like to thanks Marlène Sanjosé and J.M. Senoner for providing their numerical set-up and results for the MERCATO configuration, and C. Angelberger and S. De Chaisemartin for several fruitful discussions. We also acknowledge the computational resources of IFP énergies nouvelles and CERFACS where simulations were performed.

Size-conditioned statistics

Size-conditioned statistics are investigated here for MMEF and EL at select locations, and compared to experiment. Table . 1 and Fig. 23 show the probe positions: three are on the axis (A1-A3), three in the spray (S1-S3), one at the spray boundary (B1) and one close to the wall (W3).


P Delage

N Sultan

Y J Cui

Permeability changes in Boom clay with temperature

Keywords: Boom clay, temperature, permeability, constant head test, free water

In the framework of research into radioactive waste disposal, this paper presents some tests carried out to investigate the effects of temperature on the changes in permeability of Boom between 20 and 90°C. Constant head permeability tests were carried out in a high pressure isotropic compression tests at various temperature and isotropic stresses. The results show that the changes in permeability of Boom clay with temperature are only due to the changes in viscosity of free water with temperaure. This demonstrates on the one hand that the water involved in transfer at various temperatures is free water and on the other hand that there is apparently only few changes in the status of adsorbed water with respect to temperature between 20 and 90°C.

Introduction

Relatively few data are available about the changes in permeability in clays with temperature. This problem is of significant importance when considering the feasibility of high level radioactive waste disposal in clays at great depth. It is well known that the changes in permeability in inactive porous media are directly related to the changes in viscosity of water with temperature. This is less evident in clay and particularly plastic clays, in which significant clay-water interaction exists. It is known that clay-water interactions are sensitive to temperature changes and the temperature dependence of permeability in these clays hence requires further investigation. Permeability tests at various temperatures were carried on Boom clay samples extracted from the HADES facility in Mol (Belgium). Constant head tests were carried out at various temperatures between 20 and 100°C. by applying a constant pressure gradient across a standard triaxial sample by using the pressure-volume controllers. Tests were carried out in a high pressure temperature controlled isotropic cell in which confining and back pressures were controlled by using digital pressure-volume controllers.

2.Materials and methods

Tests were performed on samples, extracted from a depth of 223 m in the Boom clay deposit [START_REF] Decleer | Relationships between chemical, physical and mineralogical characteristics of the Rupelian Boom clay, Belgium[END_REF], in the underground facilities of SCK-CEN in Mol (Belgium). Boom clay is a stiff clay, with a plasticity index of about 50%, a natural porosity around 40% and a water content between 24 and 30%. Previous work on the thermal behaviour of Boom clay includes [START_REF] Baldi | Thermal volume changes of the mineral-water system in low-porosity clay soils[END_REF], De Bruyn and Thimus (1995), [START_REF] Belanteur | Etude des comportements mécanique, thermo-mécanique et hydro-mécanique des argiles gonflantes et non gonflantes fortement compactées[END_REF], [START_REF] Delage | On the thermal consolidation of Boom clay[END_REF], [START_REF] Cui | A thermomechanical model for saturated clays[END_REF], [START_REF] Sultan | Temperature effects on the volume change behaviour of Boom clay[END_REF][START_REF] Cui | Investigating the time dependent behaviour of Boom clay under thermo-mechanical loading[END_REF]. Tests were carried out in an isotropic compression cell (see Figure 1) designed to support high pressures and high temperatures (up to 60 MPa and 100°C respectively). A heating coil was placed on the outer wall of the cell, and the temperature was controlled by a thermocouple thermometer placed inside the cell, in the confining water. The precision of the temperature regulation system was ±0.05°C. The isotropic stress was applied by a high pressure volume controller (PVC, GDS brand, 60 MPa) and the back pressure (up to 2 MPa) by a standard PVC. The main advantage of this system is its ability to monitor volume changes while applying pressures. Due to the thermal volume changes of water under non isothermal conditions, it was decided to monitor the sample volume changes from the changes in volume of the confining water. As discussed in [START_REF] Delage | On the thermal consolidation of Boom clay[END_REF], this technique is much more satisfactory at high stresses than in the range of stresses used in standard soil mechanics.

Experimental results

In order to ensure significant flow rates, a sample smaller than the standard triaxial sample was used. The standard diameter (φ = 37.7 mm) was kept but the height was reduced down to h = 22.95 mm. The sample was saturated under a low all round pressure (110 kPa) with a back pressure equal to 40 kPa. A Skemton B coefficient of 0.97 was attained after 24 hours. Afterwards, the sample was consolidated up to 2.5 MPa. The first permeability measurement was carried out at T = 20°C by increasing the back pressure up to 1 MPa at the bottom of the sample, while putting the top porous stone in contact with atmospheric pressure. It results in a high value of the hydraulic gradient (i ≈ 4000), as compared to that of variable head permeameters (i ≈ 100). The high gradient was found necessary in order to create a measurable flow in the dense plastic clay, and to get a satisfactory precision in the measurement of flow rate and, hence, of the permeability.

In Figure 2, the volume of water injected at the base is given as a function of time. In spite of the high gradient applied, the figure shows that 10 hours were necessary to achieve permanent flow and 15 hours to obtain a satisfactory determination of the slope of the curve. This section corresponds to a constant flow from which the permeability can be calculated. At 20°C, a permeability value of 2.5 10 -12 m/s was obtained with a porosity value equal to 39%. Possible problems related to the high value of the hydraulic gradient were examined by comparing this result with that of a variable head test, run in a rigid oedometer cell (i ≈ 100). A permeability value of 3.5 10 -12 m/s for a porosity of 44.1% was obtained, a result consistent with the former result. The two points are plotted in Figure 3. The test carried out for the determination of the permeability at various temperatures is presented in Figure 3, in a diagram giving the porosity versus the logarithm of the permeability. Constant head permeability tests were successively performed under an allround pressure of 2.5 MPa at 20°C, 60°C, 70°C, 80°C and 90°C. Heating phases were very progressive, at a heating rate of 0.1°C/20 minutes, in order to fulfil drained conditions (Cui et al. 2000). At 90°C, the sample was isothermally loaded from 2.5 up to 4 MPa, and permeability tests were carried out in a cooling phase at 90°C, 80°C, 70°C, 60°C and 30°C. At 30°C, subsequent isothermal loading from 4 MPa up to 6 MPa was carried out, and permeability tests were performed at 30°C and 60°C. One observes in Figure 3 that heating from 20°C up to 90°C results in a contraction of the sample (the porosity n decreasing from 39 % down to 37.2 %) and in an increase in permeability from 2.5 10 -12 m/s to 6.2 10 -12 m/s. The loading sequence at 4 MPa at a temperature of 90°C decreases the permeability down to 4.4 10 -12 m/s. During the cooling phase under 4 MPa, the permeability decreases from 4.4 10 -12 m/s down to 2.4 10-12 m/s and the porosity remains almost constant (34 % to 33.9 %). Loading up to 6 MPa under a temperature of 30°C decreases the porosity down to 31.9%. Further heating at 60°C under 6 MPa gives a final porosity of 30.6%. The changes in permeability presented in Figure 3 are due to the coupled effect of changes in temperature and porosity. In order to separate these effects, the intrinsic permeability values (K) were calculated according to :

( ) T ) T ( K k w µ γ = (1)
where µ (T) is the water viscosity and γ w (T) is the unit weight of water. The following relation, valid for free water, was derived from experimental values reported by [START_REF] Hillel | Fundamentals of soil physics[END_REF]:

( ) ( ) ) s . Pa ( 00239138 . 0 T ln 00046575 . 0 T + - = µ (2)
In the range of temperatures considered, there is no significant change in the value of γ w (T).

Observation of the results presented in Figure 4 in a semi-logarithmic diagram shows that the relationship between the intrinsic permeability K and the porosity n is unique and linear. It appears reasonably independent of temperature. The data obtained from the variable head test performed in the rigid ring oedometer also correctly matches with the results. Consequently, the intrinsic permeability of a sample loaded at a given temperature is only dependent on its porosity, independently on the thermo-mechanical path previously followed in a (p';T) plane. In other words, volume changes created by stress and/or temperature changes have the same effect on the intrinsic permeability of a sample. Data of Figure 4 are also consistent with the interpretation of [START_REF] Habibagahi | Temperature effect and the concept of effective void ratio[END_REF] and confirm that the water put in movement during a permeability test is the free water that circulates through the channels limited by the water tightly adsorbed on the clay minerals. Similar conclusions regarding the intrinsic permeability were drawn by [START_REF] Morin | The effects of high pressure and high temperature on some physical properties of ocean sediments[END_REF] from permeability tests carried out at various temperatures (between 22 and 220°C) on 4 soils of different plasticity indices and higher void ratios (I p between 52 and 179, and e between 1 and 9). They also showed that, at a given void ratio, a more plastic soil is more impervious, due to a greater number of layers of water molecules adsorbed on the clay minerals, which do not participate to the flow, as commented by [START_REF] Habibagahi | Temperature effect and the concept of effective void ratio[END_REF].

Conclusion

Constant head permeability tests carried out on Boom clay samples at various temperatures between 20° and 90°C and under various isotropic stresses demonstrated that the permeability changes observed were only linked to changes in the viscosity of free water with temperature. As a consequence, no temperature influence was observed on the changes in intrinsic permeability with temperature. This means firstly that the water molecules involved in the transfers at various temperature are free molecules and secondly there is no significant changes in the status of adsorbed water with temperature between 20 and 90°C.

Figure 1 .

 1 Figure 1. Schematic diagram of the triaxial cell.

Figure 2 .

 2 Figure 2. Volume of water injected during a permeability test

Figure 3 .

 3 Figure 3. Permeability tests performed on a sample of Boom clay at various temperatures andstresses.

Figure 4 .

 4 Figure 4. Results of the permeability tests, in terms of intrinsic permeability.


C Crobach

E J Harmanus

Kuijper

Marjolein P M Hensgens

Abraham Goorhuis

Caroline M J Van Kinschot

Monique J T Crobach

Celine Harmanus

Ed J Kuijper
email: e.j.kuijper@lumc.nl

Clostridium difficile infection in an endemic setting in the Netherlands

Keywords: Clostridium difficile infection, risk factor, endemic, fluoroquinolones, mortality 40-

infection in an endemic setting in the Netherlands

Introduction

Since 2002, outbreaks caused by Clostridium difficile infection (CDI) have been reported in Canada, the USA and Europe, associated with the emergence of a new hypervirulent type. This type has been characterized as North American pulsed-field type 1, restriction-endonuclease analysis group type BI, toxinotype III and PCR-ribotype 027 (type 027) [START_REF] Joseph | First isolation of Clostridium difficile PCR ribotype 027, toxinotype III in Belgium[END_REF][START_REF] Kuijper | Emergence of Clostridium difficile-associated disease in North America and Europe[END_REF][START_REF] Loo | A predominantly clonal multi-institutional outbreak of Clostridium difficile-associated diarrhea with high morbidity and mortality[END_REF][START_REF] Mcdonald | An epidemic, toxin gene-variant strain of Clostridium difficile[END_REF][START_REF] Muto | A large outbreak of Clostridium difficileassociated disease with an unexpected proportion of deaths and colectomies at a teaching hospital following increased fluoroquinolone use[END_REF]. During outbreaks in the USA and Canada, the reported incidences of CDI varied between 155 and 225 per 10,000 hospital admissions [START_REF] Loo | A predominantly clonal multi-institutional outbreak of Clostridium difficile-associated diarrhea with high morbidity and mortality[END_REF][START_REF] Kazakova | A hospital outbreak of diarrhea due to an emerging epidemic strain of Clostridium difficile[END_REF]. Peak incidences of CDI due to type 027 during outbreaks in the Netherlands were remarkably lower, around 50 per 10,000 hospital admissions [START_REF] Kuijper | Clostridium difficile ribotype 027, toxinotype III, the Netherlands[END_REF][START_REF] Paltansing | Characteristics and incidence of Clostridium difficile-associated disease in The Netherlands, 2005[END_REF].

Most recent studies on risk factors of CDI focussed on outbreaks, whereas less is known about CDI in settings with a low incidence. Well described risk factors for CDI in outbreak situations are prior use of antibiotics, increased disease severity, and, in case of outbreaks caused by type 027, advanced age and prior use of fluoroquinolones [START_REF] Goorhuis | Spread and epidemiology of Clostridium difficile polymerase chain reaction ribotype 027/toxinotype III in The Netherlands[END_REF][START_REF] Dubberke | Clostridium difficile--associated disease in a setting of endemicity: identification of novel risk factors[END_REF][START_REF] Hubert | A portrait of the geographic dissemination of the Clostridium difficile North American pulsed-field type 1 strain and the epidemiology of C. difficile-associated disease in Quebec[END_REF].

The aim of our study was to identify risk factors for CDI in a true endemic setting.

A second aim was to establish risk factors specific for CDI, in comparison with factors for diarrhoea in general. To answer these questions, we performed a prospective casecontrol study at the Leiden University Medical Center during a period of 34 months.

Methods

Patients.

From July 2006 through April 2009, all hospitalized patients with CDI were included in the study. Tests for CDI were performed daily upon request and on all unformed faecal samples from patients admitted for two days or more, regardless the physicians' request. For each hospitalized CDI patient, two controls were included, matched for ward at which CDI was diagnosed and time of admission. The controls included one control patient without diarrhoea (control patient) and one control patient with diarrhoea and a negative C. difficile toxin test (non-CDI patient). Controls were consecutive patients on the alphabetical ward list.

Definitions.

Definitions as proposed by the European and American Centres of Disease

Control were used [START_REF] Kuijper | Emergence of Clostridium difficile-associated disease in North America and Europe[END_REF][START_REF] Mcdonald | Recommendations for surveillance of Clostridium difficile-associated disease[END_REF]. Diarrhoea was defined as ≥3 unformed stools per 24 hours. CDI was defined as the presence of diarrhoea in combination with a positive toxin test for C. difficile. A community association was defined as development of CDI outside the hospital or within 48 hours after admission, without a history of admission in the previous three months. We defined diarrhoea as severe, when it occurred with one or more of the following: bloody stools, hypovolemia, fever (T>38.0 0 C) and leucocytosis (>12.0x10 9 /l), hypo-albuminemia (<20 g/l), pseudomembranous colitis. A complicated course of CDI was defined as: admission to the intensive care unit (ICU), a surgical intervention in association with CDI, or death within one month. Mortality was considered contributable to CDI when a patient died during admission, partly due to the consequences of CDI.

Isolation and characterization of Clostridium difficile.

C. difficile toxins in stools were detected by VIDAS C. Diff. toxin A during the first 12 months of the study and VIDAS toxin A/B assay during the ensuing 22 months (BioMérieux, France). Each positive sample was cultured. Available isolates were identified as C. difficile using a PCR to detect the presence of gluD and were PCRribotyped as previously described [START_REF] Paltansing | Characteristics and incidence of Clostridium difficile-associated disease in The Netherlands, 2005[END_REF][START_REF] Bidet | Comparison of PCR-ribotyping, arbitrarily primed PCR, and pulsed-field gel electrophoresis for typing Clostridium difficile[END_REF].

Data collection.

Approval was obtained from the Medical Review Ethics Committee to collect demographical and clinical patient data. Information was collected on patients' age, sex, co-morbidity, ward of acquisition, disease severity, clinical course and mortality.

Furthermore, data were collected on surgery, invasive procedures, admissions, use of antibiotics and other medications in the 3 months prior to CDI. We gathered this information through consultation of the physician in charge, as well as by using patient records and the hospital electronic medical information system. The period of 3 months prior to CDI was determined by calculating backwards from a reference date. For CDI and non-CDI patients, this reference date was defined as the day on which the diarrhoea started. The reference date for control patients was determined by adding the hospitalized period of the matched CDI patient (time between admission and start of diarrhoea) to the admission date of the control patient. Co-morbidity was assessed by both the Charlson co-morbidity index and the ICD-10 classification in ten disease groups; mentioned in table 1 [START_REF] Charlson | A new method of classifying prognostic comorbidity in longitudinal studies: development and validation[END_REF].

Statistical analysis.

Continuous data were compared between groups using the T-test. The Pearson's-chi-square test and the Fisher's exact test were used for the analysis of proportions. Factors that were associated in univariate analysis (UVA) with a p-value <0.10, as well as putative risk factors from earlier studies, were analyzed in a multivariable model. Here, associations were always adjusted for age, sex, ward and Charlson co-morbidity index. To evaluate the effect of medications and interventions on (CDI) diarrhoea, we performed additional adjustments for co-medication and other interventions. When comparing non-CDI patients with control patients, we also corrected for the time between admission and the reference date. Relative risks were estimated as odds ratios (OR) and presented with a 95% confidence interval (95% CI). Statistical significance was reached with a 2-sided p-value <0.05; trends were defined by a p-value <0.10. All analyses were performed using the SPSS for Windows software package, version 17.0.

Results

During the 34 month study period, 93 patients were diagnosed with CDI. The incidence varied from 0 to 43 per 10,000 hospital admissions with an average of 17.5.

During this period, no outbreaks were observed. CDI was community associated in four patients (4.3%). Most patients (n=30; 32.3%) were hospitalized at the department of internal medicine, followed by the general surgery ward (n=15; 16.1%). Eighty-nine CDI patients were positive on both toxin testing and culture (95.7%). Isolates from 63 (67.7%) patients were available for PCR-ribotyping: type 014 was the most frequently found type (n=10; 15.9%), followed by types 078 (n=8; 12.7%) and 015 (n=5; 7.9%). Type 027 was not present. Three patients with CDI had a co-infection with an enterovirus, norovirus, and Cryptosporidium, respectively.

The 93 CDI patients were compared to 76 control patients and 64 patients with non-CDI diarrhoea. Of all patients, physicians responded and records were available, however, in some cases (the exact number is depicted in the subscript of table 1) no information about use of nasogastric intubation, surgery or endoscopy was noted.

In the group of non-CDI patients, two patients were diagnosed with a rotavirus and Giardia Lamblia, respectively. Among the other 62 patients no causal agent was found. CDI patients had a median age of 56 years; non-CDI diarrhoea and control patients had a median age of 50 years. Of the CDI patients, 60% were male, compared to 50% and 54% of the non-CDI and control patients, respectively. The time span between admission and start of diarrhoea did not significantly differ between CDI and non-CDI patients.

Characteristics and risk factors.

We present baseline characteristics and risk factors for CDI and non-CDI diarrhoea in tables 1 and 2. The use of antibiotics as a risk factor for CDI and non-CDI is depicted in table 3 Compared to control patients, patients with non-CDI diarrhoea more frequently used third generation cephalosporins but less frequently used first generation cephalosporins.

Interventions and admissions.

Patients with CDI, compared to control patients, were more frequently admitted in the previous 3 months, either at the hospital or ICU department. They also more frequently had a nasogastric intubation or an endoscopy (significant only in UVA).

Patients with non-CDI diarrhoea more frequently had a nasogastric intubation (significant only in UVA), and were more frequently admitted to the ICU in the previous 3 months (trend in MVA).

Clinical course.

Severe diarrhoea was present among 51 hospitalized patients with CDI (58.6%) and 25 patients with non-CDI diarrhoea (39.7%) (OR 2.22; 95% CI 1.14-4.30). No significant differences between CDI and non-CDI diarrhoeal patients were found regarding the frequency of fever (55.6% resp. 43.3%), bloody stools (12.2% resp. 12.9%) or abdominal pain (54.5% resp. 48.2%). CDI patients did however have a higher white blood cell count (≥15 x 10 9 /l: 49.9% resp. 30.0%, OR 2.28; 95% CI 1.13-4.59).

Most patients with CDI were treated with metronidazole (n=57; 63.3%), two patients (2.2%) were treated with vancomycin and in 27 patients (30.0%) no specific CDI

Discussion

In this 34 months prospective case control study, risk factors for CDI were studied in an endemic setting with a low incidence rate of CDI. The inclusion of a control group of patients with diarrhoea, tested negative for CDI, enabled us to discriminate between risk factors for CDI and for diarrhoea in general.

Common risk factors for CDI outbreaks, such as age above 65 years and a high comorbidity index, were recognized as trends in our study. This may be due to the fact that these risk factors are of less importance in endemic settings, resulting in a lack of power to discern these risk factors. Other well known risk factors for CDI, such as the use of second generation cephalosporins and previous (ICU) admission were also found in this endemic situation [START_REF] Loo | A predominantly clonal multi-institutional outbreak of Clostridium difficile-associated diarrhea with high morbidity and mortality[END_REF][START_REF] Dubberke | Clostridium difficile--associated disease in a setting of endemicity: identification of novel risk factors[END_REF][START_REF] Nelson | Epidemic Clostridium difficileassociated diarrhea: role of second-and third-generation cephalosporins[END_REF]. Conversely, the use of fluoroquinolones or PPIs was not a risk factor for CDI. Furthermore, the previous use of third generation cephalosporins was a risk factor for diarrhoea in general.

The CDI incidence in our hospital was lower than that described in other studies in endemic situations, but comparable to the incidence of 18 per 10,000 hospital admissions found in other Dutch hospitals [START_REF] Hensgens | Decrease of hypervirulent Clostridium difficile PCR ribotype 027 in the Netherlands[END_REF]. Recently, a retrospective study analyzing risk factors for CDI in an endemic setting in USA reported an incidence rate of CDI of 106 per 10,000 hospital admissions, which is a factor 5 higher than what we found in this study [START_REF] Dubberke | Clostridium difficile--associated disease in a setting of endemicity: identification of novel risk factors[END_REF]. There seems to be a considerable difference, per hospital and per country, in the application of the definition of endemic CDI. Therefore, reported rates of endemic CDI may merely reflect a baseline incidence.

In outbreak situations, the previous use of fluoroquinolones has been recognized as an important risk factor for CDI [START_REF] Goorhuis | Spread and epidemiology of Clostridium difficile polymerase chain reaction ribotype 027/toxinotype III in The Netherlands[END_REF][START_REF] Hubert | A portrait of the geographic dissemination of the Clostridium difficile North American pulsed-field type 1 strain and the epidemiology of C. difficile-associated disease in Quebec[END_REF][START_REF] Bourgault | In vitro susceptibility of Clostridium difficile clinical isolates from a multi-institutional outbreak in Southern Quebec, Canada[END_REF][START_REF] Pepin | Emergence of fluoroquinolones as the predominant risk factor for Clostridium difficile-associated diarrhea: a cohort study during an epidemic in Quebec[END_REF]. This association may be due to disruption of the gut flora by newer fluoroquinolones or the high fluoroquinolone resistance found among hypervirulent type 027 strains [START_REF] Weiss | Clostridium difficile and fluoroquinolones: is there a link?[END_REF]. Although fluoroquinolones (mainly ciprofloxacin) were frequently prescribed in this study, we found no association with CDI. An explanation could be that we did not encounter type 027 in our hospital.

The most frequently found PCR-ribotypes in our study (types 014, 078 and 015) are commonly found in the Netherlands and Europe and are more susceptible to fluoroquinolones than type 027 [START_REF] Hubert | A portrait of the geographic dissemination of the Clostridium difficile North American pulsed-field type 1 strain and the epidemiology of C. difficile-associated disease in Quebec[END_REF].

The use of vancomycin was previously recognized as a risk factor for endemic CDI [START_REF] Dubberke | Clostridium difficile--associated disease in a setting of endemicity: identification of novel risk factors[END_REF]. Instead, in this study, the association between vancomycin and CDI was strongly confounded by concomitant use of second and especially third generation cephalosporins (the combination is part of the in-house empirical sepsis therapy) and was not a risk factor for CDI.

PPIs raise the gastric pH, which is associated with enhanced bacterial colonization of first part of the GI-tract. Studies to the use of PPIs in association with CDI revealed conflicting conclusions [START_REF] Dial | Use of gastric acid-suppressive agents and the risk of community-acquired Clostridium difficile-associated disease[END_REF][START_REF] Lowe | Proton pump inhibitors and hospitalization for Clostridium difficile-associated disease: a populationbased study[END_REF]. In our study, we found no association of the use of PPIs with CDI. It should be noted that half of the non-CDI and control patients also used PPIs.

Earlier studies have found high contamination and colonization rates with C. difficile spores in the hospital environment, among hospitalized patients and asymptomatic carriers [START_REF] Dubberke | Prevalence of Clostridium difficile environmental contamination and strain variability in multiple health care facilities[END_REF][START_REF] Mcfarland | Nosocomial acquisition of Clostridium difficile infection[END_REF]. A high colonization pressure on a ward (exposure in time to multiple colonized or infected patients) is associated with an increased risk of CDI [START_REF] Dubberke | Clostridium difficile--associated disease in a setting of endemicity: identification of novel risk factors[END_REF]. To insure that CDI and control patients were exposed to a similar colonization pressure, we selected control patients from the same ward as CDI patients using the same time period between admission and reference date [START_REF] Goorhuis | Novel risk factors for Clostridium difficileassociated disease in a setting of endemicity?[END_REF].

We observed contributable and overall mortality rates of 3.2% and 7.5% after 30days follow-up, respectively. These rates are in between the overall 30-day mortality rate of 4.7%, found in an endemic setting in Canada, and a reported rate of 20% after 60 days in a U.S. study [START_REF] Dubberke | Attributable outcomes of endemic Clostridium difficile-associated disease in nonsurgical patients[END_REF][START_REF] Pepin | Clostridium difficile-associated diarrhea in a region of Quebec from 1991 to 2003: a changing pattern of disease severity[END_REF]. These mortality rates are much lower than those reported during outbreaks caused by the type 027 strain [START_REF] Loo | A predominantly clonal multi-institutional outbreak of Clostridium difficile-associated diarrhea with high morbidity and mortality[END_REF][START_REF] Hubert | A portrait of the geographic dissemination of the Clostridium difficile North American pulsed-field type 1 strain and the epidemiology of C. difficile-associated disease in Quebec[END_REF][START_REF] Pepin | Clostridium difficile-associated diarrhea in a region of Quebec from 1991 to 2003: a changing pattern of disease severity[END_REF][START_REF] Pepin | Mortality attributable to nosocomial Clostridium difficile-associated disease during an epidemic caused by a hypervirulent strain in Quebec[END_REF]. In the Netherlands, a complicated course due to type 027 was described in 12,5%, with an attributable mortality rate of 6.3% [START_REF] Goorhuis | Spread and epidemiology of Clostridium difficile polymerase chain reaction ribotype 027/toxinotype III in The Netherlands[END_REF].

Our study has several limitations. First, we used the presence of toxins in faeces as a screening test for CDI which is in agreement with the European recommendations [START_REF] Crobach | European Society of Clinical Microbiology and Infectious Diseases (ESCMID): data review and recommendations for diagnosing Clostridium difficile-infection (CDI)[END_REF]. An alternative standard for diagnosing CDI is the detection of C. difficile in faeces by toxinogenic culture or PCR. Application of this definition could have resulted in a different case and non-CDI control group. However, none of the patients with non-CDI diarrhoea developed CDI at a later moment during admission, which was in accordance with the high negative predictive value of our toxin test. Second, although the endemic incidence found in our study is comparable to that in other Dutch hospitals, it is lower than incidence rates reported in other studies in endemic situations, which can imply that our findings may not be applicable to endemic situations in other countries [START_REF] Paltansing | Characteristics and incidence of Clostridium difficile-associated disease in The Netherlands, 2005[END_REF][START_REF] Pepin | Clostridium difficile-associated diarrhea in a region of Quebec from 1991 to 2003: a changing pattern of disease severity[END_REF][START_REF] Alfa | Survey of incidence of Clostridium difficile infection in Canadian hospitals and diagnostic approaches[END_REF].

In conclusion, in this endemic setting, some risk factors for CDI were similar to those found in outbreak situations, but some risk factors that have earlier been ascribed to CDI were in this study not specific for CDI, but for diarrhoea in general. Use of fluoroquinolones and PPIs did not influence the risk of endemic CDI. CDI patients were more severely ill than non-CDI diarrhoeal patients, as illustrated by a higher leukocyte count and the relatively high 30-and 60-day mortality. Because CDI is the most important cause of nosocomial diarrhoea, more studies are needed to determine the long term outcome associated with C. difficile infections.

Figure 1 :

 1 Figure 1: Survival curve of patients with CDI, non-CDI diarrhoea and control patients, in a period of 60 days after the reference date.

  . All following results reached statistical significance in multivariable analysis (MVA), unless otherwise stated. Age. Patients with CDI were older than control patients (age > 65 years in 35.5% vs. 23.7%; trend in MVA). . Both CDI and non-CDI diarrhoeal patients had a higher Charlson co-morbidity index (index of 3-4 or >5) than control patients (not significant). CDI patients were more likely to have haematological malignancies, diseases of the urogenital tract or other diseases (all trends in MVA). The category 'other diseases' comprised organ transplants in 69.7%. Use of medications. Compared to control patients, patients with CDI more frequently used antibiotics, specifically second and third generation cephalosporins. CDI patients also more frequently used penicillin and vancomycin (all significant only in UVA). Furthermore, CDI patients used proton pump inhibitors (PPIs) more frequently (significant only in UVA). The use of antacids (17.2% vs. 18.4%; OR 0.68; 95% CI 0.26-1.79) or the combined use of PPIs and antacids (74.2% vs. 59.2%; OR 0.75; 95% CI 0.29-1.95) was not significantly more frequent in patients with CDI in MVA (data not shown in the table).

	Comorbidity

Table 1 : Baseline characteristics of patients with CDI, patients with non-CDI diarrhea and control patients.

 1

			CDI	Non-CDI	Control
			patients	patients	patients
		Risk factors	(N=93) •	(N=64) ••	(N=76) •••
			N	(%)	N	(%)	N	(%)
	Age > 65 years	33 (35.5)	18 (28.1)	18 (23.7)
	Male sex	56 (60.2)	32 (50.0)	41 (53.9)
	Charlson co-morbidity index
		0	14 (15.1)	12 (18.8)	19 (25.0)
		1-2	38 (40.9)	26 (40.6)	32 (42.1)
		3-4	28 (30.1)	15 (23.4)	18 (23.7)
		5+	13 (14.0)	11 (17.2)	7	(9.2)
	Any underlying disease	90 (96.8)	61 (95.3)	70 (92.1)
		Malignancy	24 (26.1)	18 (28.1)	21 (27.6)
		Solid tumor	10 (10.9)	5	(7.8)	11 (14.5)
		Hematologic malignancy	15 (16.1)	13 (20.3)	10 (13.2)
		Endocrine diseases	26 (28.0)	16 (25.0)	20 (26.3)
		Respiratory tract diseases	14 (15.1)	9	(14.1)	8	(10.5)
		Gastro-intesinal tract diseases	36 (38.7)	16 (25.0)	21 (27.6)
		Cardiovascular tract diseases	42 (45.2)	27 (42.2)	30 (39.5)
		Urogenital tract diseases	42 (45.2)	21 (32.8)	24 (31.6)
		Nervous system diseases	6	(6.5)	4	(6.2)	6	(7.9)
		Infectious diseases	13 (14.3)	6	(9.4)	7	(9.2)
		Muscular / conn. tissue diseases	10 (10.8)	4	(6.2)	7	(9.2)
		Other diseases	36 (39.1)	24 (37.5)	22 (28.9)
	Any antibiotic	87 (93.5)	48 (75.0)	51 (68.0)
	Proton pump inhibitors	64 (68.8)	36 (56.2)	38 (50.0)
	NSAIDs	11 (11.8)	3	(4.7)	7	(9.2)
	Immunosuppressive agents	54 (58.8)	38 (59.4)	34 (44.7)
	Cytostatic agents	21 (22.6)	13 (20.3)	11 (14.5)
	Nasogastric tube	39 (44.3)	29 (45.3)	20 (28.2)
	Abdominal surgery	35 (37.6)	24 (37.5)	20 (28.6)
	Endoscopy	28 (31.5)	16 (25.0)	10 (13.2)
	Previous admission	68 (74.7)	19 (30.2)	30 (41.7)
	Previous admission to ICU	26 (28.0)	12 (18.8)	5	(6.6)
	•	N between 88 and 93
	•• N between 62 and 64
	••• N between 71 and 76
	ICU: intensive care unit
	NSAIDs: non-steroidal anti-inflammatory drugs

Table 2 : Crude and adjusted odds ratios for development of CDI and non-CDI diarrhea.

 2

			CDI vs. Control	Non-CDI vs. Control
		Risk factors	Crude odds ratio	Adjusted odds	Crude odds ratio	Adjusted odds
			(95% C.I.)	ratio (95% C.I.)	(95% C.I.)	ratio (95% C.I.)
		Age > 65 years	1.77 (0.90-3.49) ^1.82 (0.92-3.62) ^1.26 (0.59-2.69)	1.17 (0.54-2.55)
		Male sex	1.29 (0.70-2.39)	1.30 (0.70-2.43)	0.85 (0.44-1.67)	0.88 (0.45-1.72)
		Charlson co-morbidity index
		0	Reference	Reference	Reference	Reference
		1-2	1.61 (0.70-3.72)	1.78 (0.73-4.37)	1.29 (0.53-3.13)	1.35 (0.50-3.66)
		3-4	2.11 (0.85-5.24)	2.42 (0.87-6.73) ^1.32 (0.49-3.57)	1.32 (0.41-4.30)
		5+	2.52 (0.80-7.95)	2.57 (0.76-8.65)	2.49 (0.76-8.19)	3.10 (0.82-11.7)	Âny
		underlying disease	2.57 (0.62-10.7)	2.45 (0.58-10.4)	1.74 (0.42-7.27)	2.10 (0.46-9.56)
		Hematologic malignancy	1.27 (0.54-3.01)	2.33 (0.86-6.23) ^1.68 (0.68-4.15)	2.19 (0.70-6.88)
		Urogenital tract diseases	1.78 (0.95-3.36) ^1.97 (0.97-4.02) ^1.06 (0.52-2.16)	0.99 (0.42-2.34)
		Other diseases	1.58 (0.83-3.02)	1.47 (0.72-3.00) ^1.47 (0.73-2.99)	1.41 (0.65-3.07)
		Any antibiotic	6.82 (2.62-17.8) ^^5.41 (1.79-16.3) ^^1.41 (0.67-2.98)	0.99 (0.40-2.42)
		Proton pump inhibitors	2.21 (1.18-4.14) ^^1.14 (0.51-2.58)	1.29 (0.66-2.51)	1.01 (0.46-2.22)
		NSAIDs	1.32 (0.49-3.60)	0.86 (0.27-2.73)	0.49 (0.12-1.96)	0.34 (0.07-1.57)
		Immunosuppressive agents	1.71 (0.93-3.15) ^1.39 (0.64-3.06)	1.81 (0.92-3.54) ^1.44 (0.64-3.22)
		Cytostatic agents	1.72 (0.77-3.85)	1.61 (0.61-4.24)	1.51 (0.62-3.64)	1.64 (0.58-4.63)
		Nasogastric tube	2.03 (1.04-3.95) ^^1.50 (0.66-3.43)	2.11 (1.04-4.31) ^^1.77 (0.70-4.50)
		Abdominal surgery	1.51 (0.77-2.94)	1.17 (0.56-2.45)	1.50 (0.73-3.10)	1.28 (0.57-2.84)
		Endoscopy	3.03 (1.36-6.75) ^^2.64 (1.00-6.96) ^2.20 (0.92-5.27) ^2.63 (0.90-7.64)	Previous
		admission	4.14 (2.13-8.05) ^^4.49 (2.23-9.01) ^^0.61 (0.30-1.23)	0.55 (0.26-1.17)
		Previous admission to ICU	5.51 (2.00-15.2) ^^5.47 (1.95-15.3) ^^3.28 (1.09-9.87) ^^2.64 (0.83-8.37)	Trend
		detected (p<0.10)
		^^Significant difference (p<0.05)
		ICU: intensive care unit
	349	NSAIDs: non-steroidal anti-inflammatory drugs

Table 3 : The use of antibiotics expressed in Defined Daily Doses (DDDs) in patients with CDI and non-CDI diarrhea and control patients.

 3

	Antibiotics	CDI patients (N=93) •	Non-CDI (N=64) • patients	Control (N=76) • patients	CDI vs. Control	Non-CDI vs. Control
		N	(%)	N

treatment was initiated. The 30-day and 60-day mortality rates are depicted in figure1.At one month follow-up, a complicated course was observed in 9 CDI patients (10.3%), comprising two colectomies, four ICU admissions due to CDI and seven deaths (7.5%). CDI contributed directly to three of these deaths, but was not the primary cause. One non-CDI patient (1.6%) and none of the control patients died at one month follow-up. No significant association were detected between the severity of the diarrhoea, treatment or outcome.

Transparency declarations

This study was funded by a grant from ZonMw (grant number 4726)