Benoît Camenen

Damien Pham Van Bang

Modelling the settling of suspended sediments for concentrations close to the gelling concentration

Keywords: sedimentation, consolidation, settling velocity, permeability, high concentration, cohesive sediments, gelling concentration

This paper deals with the sedimentation of highly concentrated sediment suspensions (cohesive as well as non-cohesive) and the beginning of the consolidation of cohesive sediments. Based on a comparison of existing empirical formulas and experimental data, the particle Reynolds number was shown to be of importance for the behaviour of particularly non-cohesive sediments. In addition it plays a role in determining whether one or two interfaces develop during the sedimentation phase. In the case of cohesive sediments, the estimation of the gelling concentration, although difficult, seems to be fundamental. Some suggestions on the estimation of the permeability coefficient and total settling function are then given in order to improve the modelling of the sedimentation and consolidation behaviour for concentrations close to the gelling concentration.

Introduction

The interest in understanding physical characteristics of cohesive sediment has increased significantly during the last 30 years. Cohesive sediments play an important role in river and estuary engineering because of their capability to bind pollutants. Another important issue is the transport, sedimentation and consolidation of fine sediments in reservoir, navigation channels or harbour basins.

Most of the research has been carried out on the settling of low-concentration mud suspensions [START_REF] Van Leussen | Aggregation of particles, settling velocity of mud flocs -a review[END_REF]Dyer & Manning, 1999, among others) or on the consolidation regime [START_REF] Gibson | The theory of onedimensional consolidation of saturated clays. i. finite nonlinear consolidation of thin homogeneous layers[END_REF]Been & Sills, 1981 among others). Few authors [START_REF] Winterwerp | On the dynamics of high-concentrated mud suspensions[END_REF][START_REF] Camenen | Settling velocity of sediments at high concentrations[END_REF] have studied hindered settling of cohesive sediments where large concentrations of sediments do affect the settling velocity. [START_REF] Dankers | Hindered settling of mud flocs: theory and validation[END_REF] attempted to link the hindered regime with the consolidation regime intro-ducing a total settling velocity or settling function including permeability effects. It seems however that many unknowns still exist about the behaviour of a mud suspension at the onset of its gelling concentration (defined as the concentration for which the flocs form a space-filling network, [START_REF] Winterwerp | Introduction to the physics of cohesive sediment in the marine environment[END_REF]. In estuaries, wave induced turbulence can strongly affect the hindered regime [START_REF] Gratiot | On the determination of the settling flux of cohesive sediments in a turbulent fluid[END_REF], which can create highly concentrated fluid mud layer with concentrations close to the gelling concentration. On the other hand, the mechanics of lutocline (interface where a sharp gradient in sediment concentration exists and generally maintained thanks to turbulence) and fluid mud are also strongly affected by the hindered settling velocity of the sediments [START_REF] Wolanski | Mixing across a lutocline[END_REF]. [START_REF] Mehta | Review notes on cohesive sediment erosion[END_REF] also showed that erosion processes are deeply influenced by the state of the mud (fluid, non-consolidated, and consolidated). It is thus fundamental to better understand the behaviour of cohesive sediments close to the gelling concentration.

A study on the effects of the choice of hindered settling formulas for the sedimentation regime is presented in this paper using the Kynch theory (1952). Based on experiments of sedimentation in a quiescent fluid (where turbulence is negligible), the ability of the formulas to induce one or two interfaces depending on the particle Reynolds number is emphasized. The second part of this paper establishes a permeability equation and a total settling function that can be used in an advection-diffusion model.

1DV equations for a sedimentation-consolidation model

Equation for sedimentation and consolidation

Assuming that solid particles are of the same size, shape and density (a1), that both the solid particles and the fluid of the suspension are incompressible (a2), that the flow is one-dimensional (a3), and that the settling velocity of particles in a dispersion and the dissipation coefficient can be determined by the local particle concentration only (a4), a general 1D equation for the sedimentation and consolidation regimes may be [START_REF] Toorman | Sedimentation and self-weight consolidation: general unifying theory[END_REF][START_REF] Toorman | Sedimentation and self-weight consolidation: constitutive equations and numerical modelling[END_REF] :

∂c ∂t + ∂ ∂z [W sk (c)c] + ∂ ∂z D(c) ∂c ∂z = 0 (1)
where W sk (c) is the total settling velocity function including permeability effects and D(c) is the total dissipation coefficient, and c the volumetric concentration of matter. The total dissipation coefficient is the sum of the molecular diffusion effects (D m ), consolidation diffusion effects (D c ), and eddy diffusivity or turbulence effects (D t ). The molecular diffusion D m is generally negligible and the eddy diffusivity D t may be estimated using a kǫ turbulence model. However, for an experiment with no hydrodynamic constraints (waves or current) and considering the flow is one-dimensional, the turbulence effects are negligible. A sedimentation consolidation experiment is thus represented by three main regimes [START_REF] Toorman | A hindered settling model for the prediction of settling and consolidation of cohesive sediment[END_REF][START_REF] Dankers | Hindered settling of mud flocs: theory and validation[END_REF]; cf. Fig. 1) :

1. hindered regime, where the settling velocity is mainly driven by the concentration of particles ( the total pressure minus the neutral pressure of water in pores or effective stress σ ′ = 0) ; 2. permeability regime, where the settling velocity function is mainly driven by the permeability. This mode physically represents compression and expulsion of pore water (σ ′ ≈ 0) ; 3. effective stress regime, where particle deformation causes further compression (σ ′ > 0; diffusion term due to D c , cf. Eq. 1).

The first regime corresponds to the sedimentation regime, whereas the two last regimes correspond to the consolidation regime. However, the permeability regime may be modelled in a similar way as the sedimentation regime, i.e. with D = 0. The following study will focus on the estimation of the total settling velocity function. [START_REF] Kynch | A theory on sedimentation[END_REF] developed a theory on sedimentation based on the batch experiment. This experiment consists of an initially well-mixed suspension of particles with a concentration uniformly initialized to c 0 in a settling column subject to gravity only. This suspension will separate into three different phases (one approaching the maximum packing particle concentration c = c max at the bottom of the container, one with a well-mixed suspension c = c 0 in the middle of the container at the initial average concentration, and clear water at the top of the container). Using the same assumptions as for Eq. 1 and adding that solid particles are all small with respect to the settling column, the mass-balance of the suspended sediments (Eq. 1) simplifies into the simple wave equation in z-direction, known as the Kynch equation:

The Kynch model for the sedimentation of particles

∂c ∂t + ∂Φ ∂z = 0 (2)
where Φ = W sh c is the sediment flux or the Kynch batch flux density function, and W sh , the hindered settling velocity, which is assumed to be a function of c and W s0 (settling velocity of a single particle) only.

Sedimentation and batch theory

The Kynch sedimentation experiment is a simple and interesting test to quantify and calibrate the hindered settling velocity formulas, including the beginning of the permeability regime for cohesive sediments. As larger uncertainties are observed for cohesive sediments because of their stochastic characteristics (density, size and shape of the population of sediments), the hindered settling velocity formulas were first tested against non-cohesive data. The hindered settling formulas studied in this paper were previously presented and compared by [START_REF] Camenen | Settling velocity of sediments at high concentrations[END_REF]. Equations are reported in appendix A.

Non-cohesive sediments

In case of non-cohesive sediments, the Kynch equation (Eq. 2) is fully valid. Assuming W sh is a function of concentration only, [START_REF] Kynch | A theory on sedimentation[END_REF] showed that Eq. 2 may be rewritten as follows:

∂c ∂t + W s0 F(c) ∂c ∂z = 0 (3)
where f (c) = W sh /W s0 corresponds to any semi-empirical equation for the hindered effects on the settling velocity, and

F(c) = ∂[c f (c)] ∂c (4)
Eq. 3 may be solved by integrating along the characteristic lines dz/dt = W s0 F(c) (method of characteristics). Concentration gradients increase when the characteristic lines converge [START_REF] Winterwerp | Introduction to the physics of cohesive sediment in the marine environment[END_REF]. Moreover, F , and therefore ∂Φ/∂c, may have a minimum at c = c cr . For c < c cr , dF/dc < 0, and two interfaces (upper interface corresponding to the settling front and lower interface corresponding to the bed front) will develop (cf. Fig. 2a). For c > c cr , dF/dc > 0, and only the upper interface develops (cf. Fig. 2b).

The choice of the hindered settling formulas may significantly influence the results for the Kynch theory, even if the formulas yield similar prediction for the settling velocity. Indeed, if two hindered functions f a and f b yield same result for

c = c 1 ( f a (c 1 ) ≈ f b (c 1 )), F a (c 1 ) = {∂[c f a (c)]/∂c} c=c 1 may differ significantly from F b (c 1 ) = {∂[c f b (c)]/∂c} c=c 1 .
It results from the variability of the semi-empirical expressions for the hindered functions (see appendix A). Thus, the prediction of c cr , which separate the two regimes or of c m for which the maximum flux density Φ = cW s is reached may vary significantly with the choice of the hindered settling formula as c cr and c m correspond to the 1st and 2nd derivatives of the function F with respect to c, respectively.

The particle Reynolds number R e * = W s d/ν (ν is the kinematic viscosity of water), which characterizes the flow regime around a settling particle, is a measure of the relative weight of the inertial forces and frictional forces to the total drag forces. For a single particle, the drag coefficient decreases with the particle Reynolds number. And for large particle Reynolds numbers, the velocity of a particle may thus be less affected by the wake of the surrounding particles. Based on experimental results, [START_REF] Richardson | Sedimentation and fluidisation: Part I[END_REF] suggested the formula f = (1c) n where the index n decreases with the particle Reynolds number from 4.8 (R e * ≪ 1) to 2.2 (R e * ≫ 1). With this formula, the maximum of the function c f increases from 0.06 to 0.14 for large R e * ; and, as a consequence, knowing that d(c m f (c m ))/dc = 0 and f (c max ) = 0, the curve steepens when c > c m , and the coordinate of the inflection point (c = c cr ) also increases with the particle Reynolds number (see also in Fig. 4). Using [START_REF] Rowe | A convenient empirical equation for estimation of the Richardson-Zaki exponent[END_REF] [START_REF] Oliver | The sedimentation of suspensions of closely-sized spherical particles[END_REF][START_REF] Shannon | Batch and continuous thickening: Prediction of batch settling behavior from initial rate data with results for rigid spheres[END_REF][START_REF] Baldock | Setting velocity of sediments at high concentrations[END_REF], the flux density function (c f = Φ/W s0 ) is plotted in Fig. 3 versus the relative concentration c/c max with the particle Reynolds number emphasized. The general behaviour of the Richardson & Zaki formula is confirmed by the data. A hindered settling formula should thus take into account the effect of the particle Reynolds number. Due to the scatter in the data and relatively small number of data points, it is however difficult to precisely identify the position of the critical concentration c cr (inflection point). Only an analytical expression of f (or a curve fitting of the data) can be used to estimate c cr . Expressions for c f and F are reported in the literature (cf. Camenen, 2008 for details). The [START_REF] Richardson | Sedimentation and fluidisation: Part I[END_REF]), Souslby (1997, p.135-136), [START_REF] Winterwerp | On the dynamics of high-concentrated mud suspensions[END_REF] and [START_REF] Dankers | Hindered settling of mud flocs: theory and validation[END_REF] formulas, as well as the modified Richardson & Zaki formula and Camenen formula based on mixture theory [START_REF] Camenen | Settling velocity of sediments at high concentrations[END_REF] are plotted as a function of concentration in Fig. 4 (see Appendix for the description of the formulas). Following [START_REF] Camenen | Settling velocity of sediments at high concentrations[END_REF], the modified Richardson & Zaki formula will be used as a reference : The flux density function c f using several formulas is plotted against the volumetric concentration in Figs. 4a and4c for two different particle Reynolds numbers (R e * = 0.1 and 100, respectively). It is observed that the formulas yield quite different results for the prediction of the hindered settling velocity. Consequently, the formulas also yield large differences in the prediction of c m , and especially c cr . In Table 1, the results (c m , c m f (c m ), and c cr ) are presented for two different particle Reynolds numbers (R e * = 0.1 and 100, respectively) depending on the formula used, together with an estimation based on polynomial fits over experimental data. An estimation of the error made on the estimation of c m , c m f (c m ), and c cr is also provided. The [START_REF] Winterwerp | On the dynamics of high-concentrated mud suspensions[END_REF] and [START_REF] Dankers | Hindered settling of mud flocs: theory and validation[END_REF] formulas (obtained for cohesive sediments) seem inaccurate for non-cohesive sediments as they do not take into account the effect of the particle Reynolds number, even though it clearly influences the results. Even for small Reynolds number, the Winterwerp formula fails to predict c m or c cr (see Tab. The Souslby formula (1997, p.135-136), modified Richardson & Zaki formula and Camenen formula based on mixture theory [START_REF] Camenen | Settling velocity of sediments at high concentrations[END_REF] yield very similar results, all implying that c m and c cr do vary with the particle Reynolds number. They show, however, different behaviour for large concentrations close to c max : the Soulsby formula tends to overestimate the flux as it yields Φ > 0 when c = c max ; for relatively high particle Reynolds number (R e * > 100, cf. Figs. 4d and5), the modified Richardson & Zaki yields no value for c cr , whereas for relatively low particle Reynolds number, it yields two values for c cr (as observed by [START_REF] Shannon | Batch and continuous thickening: Prediction of batch settling behavior from initial rate data with results for rigid spheres[END_REF]; the Camenen formula yields a function F = 0 when c = c max (c f tangential to the line y(c) = 0), which seems to disagree with experimental data.

W s W s0 = (1 -c) n-1 1 - c c max c max (7) (a) (b) 0 0.2 0.4 0.6 0.8 1 0 0.05 0.1 0.15 c / c max c f(c) R&Z R&Zm Winterwerp D&W Soulsby Camenen 0 0.2 0.4 0.6 0.8 1 -0.5 0 0.5 1 c / c max F(c) R&Z R&Zm Winterwerp D&W Soulsby Camenen (c) (d) 0
R e * = 0.1 R e * = 100 c m /c max Φ(c m ) c cr /c max c m /c max Φ(c m ) c cr /c
In Fig. 5, c m and c cr were estimated using the Richardson & Zaki, modified Richardson & Zaki, Soulsby, and Camenen formulas as well as experimental data. Results based on experimental data were obtained using polynomial fits over data points (cf. Fig. 3). For small R e * , c m is less than 0.3c max , and increases to 0.4-0.5c max at higher R e * . The four formulas are in good agreement with experimental estimations. A slight dispersion occurs depending on the formula for large particle Reynolds numbers. For small R e * , c cr is less than 0.55c max , and increases to 0.9c max at higher R e * . The modified Richardson & Zaki formula predicts a second value for c cr (called c cr,2 ) for low R e * which decreases from 0.95 (R e * ≪ 1) to 0.85 at a critical particle Reynolds number ℜ * ,cr ≈ 50. Above this critical value, no solution exist for c cr , which means that the modified Richardson & Zaki formula will always result in the occurrence of two interfaces in a settling suspension, whatever the initial concentration. Again, the four formulas are in good agreement with experimental estimations when R e * < 200. For R e * > 200, it is difficult to conclude as there are large uncertainties on the experimental estimation; c cr was actually not always observed depending on the polynomial fit.

Cohesive sediments

Batch theory is not exactly applicable to cohesive suspensions, because they correspond to a population of sediments with varying size, density and shape. Moreover, as a floc is made of sediment matter and water, the volumetric concentration c does not properly describe the floc volumes. Assuming that this population may be represented by a single particle with fixed characteristics, [START_REF] Kranenburg | Hindered settling and consolidation of mud -analytical results[END_REF] applied Eq. 2 to the volumetric concentration of flocs φ (φ = (ρ sρ)/(ρ fρ)c , where ρ, ρ s , and ρ f are the water, sediment, and flocs densities, respectively), and rewrote the equation as follows:

∂φ ∂t + W s0 F(φ) ∂φ ∂z = 0 (8) F(φ) = ∂[φ f (φ)] ∂φ (9)
where f (φ) = W sh /W s0 .

Eq. 8 may be solved by integrating along the characteristic lines dz/dt = W s F(φ) as for the non-cohesive sediments case. In the same way as for the non-cohesive sediments, the function F may have a minimum at φ = φ cr . For φ < φ cr , dF/dφ < 0, two interfaces will also develop (cf. Fig. 6a). For φ > φ cr , dF/dφ > 0, only the upper interface is visible (cf. Fig. 6b). One main difference with sand is that cohesive particles are deformable. It explains why φ max > c max (cf. Camenen, 2008). First, water is ejected from pore space only (0.6 φ < φ max ). When the gelling concentration is reached (φ = φ max ≈ 0.85 when c = c gel ) and consolidation begins, the volumetric concentration of flocs ceases to increase; the floc density ρ f increases, and so, the volumetric concentration of matter increases (c gel < c bed,1 < c bed,2 < c bed,3 < c max ).

In Fig. 7, flux density functions φ f and F are plotted versus φ/φ max using the formulas considered earlier, and using a typical particle Reynolds number for cohesive sediments : R e * = 0.1. For all of the following tests, a constant gelling concentration C gel = 100 g/l (or c gel = 0.038) was assumed, however, this value may also vary with the particle Reynolds number. The maximum volumetric concentration of flocs φ max was set to 0.85. Following [START_REF] Camenen | Settling velocity of sediments at high concentrations[END_REF], the reference formula was the modified Richardson & Zaki. For cohesive sediments, assuming the size and density of the floc to be constant, and n > 2, Eq. 7 may be modified as : This differs slightly from the original suggestion (Eq. 11 in [START_REF] Camenen | Settling velocity of sediments at high concentrations[END_REF], which had a factor of (1c) n-2 instead of (1c) n/2 ; Because of buoyancy effects, one should have (1c) x with x ≥ 1 whereas n -2 could be less than 1). In Table 2, the results (φ m , φ m f (φ m ), and φ cr ) are presented for R e * = 0.1 depending on the formula used, together with an estimation based on polynomial fits over experimental data. An estimation of the error made on the estimation of φ m , φ m f (φ m ), and φ cr is also provided. The Soulsby and the Richardson & Zaki equations yield similar results as both equations include the effect of the hindered settling velocity using a coefficient (1φ) n with n ≈ 4.7 for small R e * . These minimize F at φ cr = 2/(n + 1) (cf. Eq. 6), i.e. φ cr ≈ 0.4φ max for fine sediments, which appears to be too low. Both equations underestimate Φ-values as they underestimate the settling velocity [START_REF] Camenen | Settling velocity of sediments at high concentrations[END_REF]. The Winterwerp, modified Richardson & Zaki and Camenen relationships behave similarly for the function f at low volumetric concentrations (φ/φ max < 0.5). As the Winterwerp formula does not take into account the possibility that φ max < 1 (the assumption φ max = 1 is stated), it yields a larger flux for large concentrations. However, a different calibration of C gel for the Winterwerp formula would produce estimates similar to Richardson & Zaki and Camenen relationships. On the other hand, it does not yield any minimum for the function F, and thus will always predict the existence of two interfaces in a settling suspension. Dankers & Winterwerp proposed a modification that improved this behaviour and derived a minimum for the function F. However, this modification significantly decreased the magnitude of the density function.

W s W s0 = (1 -c) n/2 (1 -φ) n/2-1 1 - φ φ max φ max (10) (a) (b) φ 0 φ cr < c bed
φ m /φ max Φ(φ m ) φ cr /
To confirm these results, experimental data obtained with cohesive sediments were analysed [START_REF] Thorn | Physical processes of siltation in tidal channels[END_REF][START_REF] Ross | Vertical structure of estuarine fine sediment suspensions[END_REF][START_REF] Wolanski | Mixing across a lutocline[END_REF][START_REF] Wolanski | The role of turbulence in settling of mudflocs[END_REF][START_REF] Mory | CBS layers in a diffusive turbulence grid oscillation experiment[END_REF][START_REF] Gratiot | Floculation des vases de gironde: mesure experimentale du flux de chute en colonne[END_REF][START_REF] Dankers | A preliminary study on the hindered settling of kaolinite flocs[END_REF]. All of these studies used natural mud with the exception of the [START_REF] Wolanski | Mixing across a lutocline[END_REF] and [START_REF] Dankers | A preliminary study on the hindered settling of kaolinite flocs[END_REF] data sets (kaolinite) and one data set from [START_REF] Mory | CBS layers in a diffusive turbulence grid oscillation experiment[END_REF] where the mud was pretreated. In the case of cohesive sediments, large uncertainties occur as the size and density of the floc cannot be measured accurately and because these values represent a statistical measure of a population. In dilute suspensions, the size (and the density) of the floc varies with the concentration (due to flocculation). It is assumed that the floc characteristics (size and density) did not change during the sedimentation experiments, and were identical to the flocs at the maximum concentration, below which hindered effects were not observed (C ≈5 to 10 g/l). The assumption of a constant mean floc size may be justified since the flocculation effects (increase of the mean size of the flocs with concentration) may be balanced by the concentration effects (decrease of the mean size of the flocs with concentration due to floc breakup). The density of the floc is then calculated using the estimated gelling concentration, i.e. ρ f = ρ + (ρ sρ)c gel /φ max (with φ max = 0.85). This latter was evaluated to maximize the fit to experimental data points (C gel varying from from 30g/l to 90g/l depending on the data set), with the exception of the Dankers et al. data set, where the gelling concentration (C gel ≈ 85 g/l) was estimated experimentally.

In Fig. 8, the flux density function (Φ/W s0 = φ f ) is plotted versus the relative concentration φ/φ max using the data collected for cohesive sediments, with particle Reynolds number R e * emphasized. The results show a lot of scatter due to experimental uncertainties (including estimation of gelling concentration) and the variability in floc characteristics, as they were assumed constants. The data corresponding to φ > φ max should correspond to the consolidation regime. As data points are based on volumetric concentration of matter c, physically unrealistic values of φ are observed in Fig. 8 when c > c gel as φ should be lower than φ max . When c > c gel , φ = φ max but the floc density increases as water is expelled from the floc. It is difficult to determine the influence of the particle Reynolds number from these experimental results, because of the scatter in the data and also because they cover a narrow range of particle Reynolds numbers (2 × 10 -2 < R e * < 2). Using polynomial fits over data sets, φ m and φ cr were approximated (see also Tab. 2 and Fig. 9) and a rough estimation of φ m and φ cr may be suggested : 0.2 < φ m /φ max < 0.4 and 0.6 < φ cr /φ max < 0. [START_REF] Dankers | A preliminary study on the hindered settling of kaolinite flocs[END_REF] found from observations in their settling column (where C gel = 85 g/l and R e * = 3 × 10 -2 ) and by comparing the rising bed in the concentration time series, that the settling behaviour changed for 0.53 < φ cr < 0.55, which means φ cr /φ max ≈ 0.68. These values are consistent with results of the modified Richardson & Zaki and Camenen equations (φ cr ≈ 0.6-0.7φ max , cf. Tab. 2).

In Fig. 9, φ m and φ cr were estimated using the Richardson & Zaki, modified Richardson & Zaki, Soulsby, and Camenen formulas as well as experimental data. Results based on experimental data were obtained using polynomial fits over data points when possible (cf. Fig. 8). The φ m and φ cr values vary significantly among empirical formulas. The sen-

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 φ / φ max Φ / W s0 R e* ≤ 5 10 -2 5 10 -2 < R e* ≤ 0.2 0.2 < R e* ≤ 1 R e* > 1
Eq. DW Eq. RZm Eq. Camenen

Permeability regime

As proposed by [START_REF] Toorman | Sedimentation and self-weight consolidation: general unifying theory[END_REF][START_REF] Toorman | Sedimentation and self-weight consolidation: constitutive equations and numerical modelling[END_REF] and also [START_REF] Winterwerp | On the dynamics of high-concentrated mud suspensions[END_REF], the initial phase of the consolidation regime (ejection of pore water) may be included in the 1D equation for sedimentation through the total settling velocity function W sk (cf. Eqs. 1 and 2). The difficulty lies in determining the permeability of mud.

Estimation of the permeability

Many relationships have been proposed to estimate permeability. Permeability is often a function of the void ratio e vr = 1/c-1. Among others, Townsend & McVay (1990) used a power law relationship (cf. Eq. 11) whereas [START_REF] Bartholomeeusen | Compound shock waves and creep behaviour in sediment beds[END_REF] used an exponential function (cf. Eq. 12) :

k = A kp e vr B kp (11) k = A ke exp(B ke e vr ) (12)
Both formulations, however, require the calibration of two coefficients. It should also be noted that these relations (Eqs. 11 and 12) are mainly empirical and derived from data where effective stress are expressed in kPa or even MPa (far from the permeability regime defined in section 2.1). [START_REF] Merckelbach | Consolidation and strength evolution of soft mud layers[END_REF] and [START_REF] Merckelbach | Equations for effective stress and permeability of soft mud Űsand mixtures[END_REF] suggested a fractal approach to estimate permeability, and their formulation includes the effect of sand content :

k = A k f c mud 1 -c sand -2 3-n f (13)
where c mud and c sand are the volumetric concentrations for mud and sand, respectively (c mud = (1p sand )c and c sand = p sand c where p sand is the sand content (in %) and c is then the total volumetric concentration of matter), and n f ≤ 3 is the fractal dimension. When sand content is negligible, c mud /(1c sand ) ≈ c. The degree to which permeability depends on concentration may be indicated by n f . The coefficient A k f (A kp or A ke ) remains to be determined. [START_REF] Bartholomeeusen | Sidere: numerical prediction of large-strain consolidation[END_REF] made a prediction exercise (named "Sidere"), using a batch of sediments from the river Scheldt (Antwerpen, Belgium). The grain size distribution were given approximately by the parameters d 10 , d 50 , d 90 , equal to 6, 70, and 210 µm, respectively, with a specific gravity equal to s = 2.72. Different experiments were performed with an initial density 1300 < ρ mix,init < 1550kg/m 3 greater than the structural density, which corresponds approximately to a gelling concentration c gel ≈ 0.1-0.15. Several scientists (Bartholomuseeusen, Carrier, Lin & Penumadu, Masada & Chan, Merckelbach, Van Kesteren, Winterwerp, and Znidarcic;cf. Bartholomeeusen et al., 2002 for details) calibrated their own model using the same experimental data (cf. Fig. 10). Though large differences were observed for the calibration of the coefficients "A" and "B" depending on the specific formulations used (Eq. 11, 12 or 13, see also [START_REF] Bartholomeeusen | Sidere: numerical prediction of large-strain consolidation[END_REF], similar results were observed when 0.2 < c < 0.5, i.e. c gel < c < 4c gel , which corresponds to the limit of validity of the present model (consolidation included as a sedimentation effect). However, one can observe in Fig. 10 that all these formulas yield a large variability in the results when c ≈ c gel .

To fit Eqs. 11,12,and 13 to the data requires the estimation of two coefficients with significant risks of disconnection with the settling characteristics of the sediment. One way to estimate the permeability coefficients may be use the hindered settling characteristics of the sediment since the permeability function reads W k = (s -1)kc (where s = ρ s /ρ is the relative density of matter). Eq. 13 with p sand = 0% gives k = A k c B k (with B k = -2/(3-n f )). Fixing B k , it is possible to find a concentration c = c gel /χ (with χ a fixed parameter to be calculated, χ > 1) such that the slope of the hindered velocity function is the same as the permeability function, i.e. d{log(W sh for c > c gel /χ, we obtain :

)}/dc = d{log(W k )}/dlog(c) = B k + 1. The parameter A k is determined from the relationship W k (c = c gel /χ) = (s -1)k(c = c gel /χ)c gel /χ. Then,
k = W sh (c = c gel /χ) (s -1)c gel /χ c c gel /χ B k (14)

Total settling function

As hindered settling functions yield W s = 0 when c ≥ c gel and permeability functions (Eqs. 11,12, and 13) cannot be linked to the settling velocity when c ≤ c gel , [START_REF] Winterwerp | On the dynamics of high-concentrated mud suspensions[END_REF] suggested defining a fitting function to obtain an unique equation for W sk :

W sk = W sh + (s -1)kc 1 + ξ(s -1)kc (15)
with ξ ≈ 10 4 -10 5 being a heuristic parameter to obtain a smooth transition between the descriptions for the hindered settling regime and the permeability regime (cf. Fig. 11). Using Eq. 14, a second method to compute the total settling function may be suggested : Estimates of the total settling velocity function W sk are presented in Fig. 11 and are compared to experimental data (a) from [START_REF] Wolanski | Mixing across a lutocline[END_REF] and (b) from [START_REF] Dankers | Hindered settling of mud flocs: theory and validation[END_REF]. For both cases, Eq. 10 has been plotted versus data (after an estimation of W s0 and c gel ) as well as W k based on Eq. 13 (with p sand = 0%) and fitted with the experimental data (fit of the coefficients A k and B k ). The total settling function W sk was then plotted using Eqs. 15 and 16. For Eq. 16, a second equation for the permeability function was estimated by using Eq. 14 with the same slope B k as for W k .

W sk =              W sh if c ≤ c gel χ W k = W sh (c = c gel /χ) χc c gel B k +1 if c > c gel χ (
If W k > W sh ∀c (cf. Fig. 11b), both equations yield unsatisfactory results. The Winterwerp method tends to overestimate results for the hindered settling regime. On the other hand, the present method tend to underestimate results for the permeability regime. On the other hand, if the W k function intersects the W sh function (cf. Fig. 11a), both equations yield satisfactory results. Eq. 15 does however induce a discontinuity in the slope dW sk /dc at c = c gel which seems to be unrealistic. Indeed, in reality, because of the variability of the cohesive particles in the suspension, the consolidation regime may start locally before the concentration reaches its gelling point (which is defined globally). As shown in Fig. 11, data points present a smooth transition between the hindered regime and the permeability regime , which is not reproduced by Eq. 15. Moreover, Eq. 15 needs to fit three parameters (A k , B k , and ξ) whereas Eq. 16 needs to fit only one parameter (B k ).

It is important to realize how the use of the total settling velocity function influences sedimentation dynamics. In Fig. 12, the flux density functions F (calculated using total settling velocity function W sk based on Eq. 15 or on Eq. 16 fitted with the experimental data of Dankers & Winterwerp, see also Fig. 11) have been plotted versus the relative concentration φ/φ max . It appears that Eq. 16 induces a slight decrease of the value for φ cr compared to the original formula for the hindered settling (Eq. 10). On the other hand, Eq. 15 yields a discontinuity in the function F for φ = φ max = φ gel with a change in the sign of dF/dφ. Following the study of [START_REF] Concha | A century of research in sedimentation and thickening[END_REF], it yields an interface for φ = φ max = φ gel , which does not seem to be physically realistic.

Conclusions

The effect of the hindered settling formula on predicting sedimentation has been presented. For non-cohesive sediments, the particle Reynolds number strongly affected estimates of the critical concentrations c m , where the flux reaches a maximum, and c cr , above which only one interface will be observed during a batch experiment (if c 0 > c cr ). For cohesive sediments, as the particle Reynolds number is generally smaller than 1, the effect of R e * was not as significant. One main difficulty remains in determining the gelling concentration. The modified Richardson & Zaki formula suggested by [START_REF] Camenen | Settling velocity of sediments at high concentrations[END_REF] appeared to yield the best overall results among the studied formulas compared to the experimental data, i.e. the best representation of the density function for both non-cohesive (c m , Φ(c m ), and c cr ) and cohesive (φ m , Φ(φ m ), and φ cr ) sediments. However, the valida-tion with cohesive experimental data remains difficult because of the uncertainties and the inherent variability of cohesive sediments.

Several propositions for the estimation of permeability were also discussed for concentrations two to three times larger than the gelling concentration. A fundamental issue, which is the link between the hindered regime and the permeability regime, was discussed. It appeared preferable to obtain a relationship for W sk that remains smooth at c = c gel to avoid prediction of a physically unrealistic interface at this concentration. A solution was suggested using a permeability function that is tangential to the hindered settling function.

Many uncertainties remain in the description of cohesive sediments, hindered settling velocities and consolidation of cohesive sediments : i.e. the influence of organic content and flocculation, the interaction between particles of different sizes (as muds often have a large particle size distribution), and sediment history for consolidation. The model proposed in this paper suggests some improvements compared to the existing literature. Some efforts are still needed to better understand the dynamics of suspensions close to the gelling concentration, and to be able to represent them in a relatively simple way hat could be applicable in engineering models.

Figure 1 :

 1 Figure 1: Schematic view of the three different prevailing processes in a sedimentation consolidation experiment from the flux curve (Φ = cW sk ) as a function of concentration c.

Figure 2 :

 2 Figure 2: Schematic view of the development of two (a) or one (b) interfaces in the Kynch sedimentation test with non-cohesive sediments depending on the initial concentration c 0 (the lines with an arrow correspond to the iso-concentration lines with c 0 < c 1 < c 2 < c max ).

Figure 3 :

 3 Figure 3: Flux density function (Φ/W s0 = c f ) as a function of the relative concentration c/c max for a range of particle Reynolds numbers R e * .

Figure 4 :

 4 Figure 4: Variation of the flux density functions Φ = c f (a and c) and F (b and d) with the relative concentration c/c max using the different studied formulas for two different particle Reynolds numbers R e * = 0.1 (a and b) and R e * = 100 (c and d).

  1). The Dankers & Winterwerp formula yields better results when R e * = 0.1 although it slightly overestimates c cr . Moreover, the Winterwerp and Dankers & Winterwerp formulas, as well as the Richardson & Zaki and Soulsby formulas, yield a flux density function Φ > 0 when c = c max , whereas it should be zero.

Figure 5 :

 5 Figure 5: Concentrations c m and c cr as a function of the particle Reynolds number using the Richardson & Zaki (RZ), modified Richardson & Zaki (RZm), Soulsby, and Camenen formulas as well as experimental results (cf. Fig. 3).

Figure 6 :

 6 Figure 6: Schematic view of the development of one (b) or two (a) interfaces in the Kynch sedimentationtest with cohesive sediments depending on the initial concentration φ 0 (the lines with an arrow correspond to the iso-concentration lines with φ 0 < φ 1 < φ 2 < φ max , and c gel < c bed,1 < c bed,2 < c bed,3 < c max ).

Figure 7 :

 7 Figure 7: Variation of the flux density functions φ f (a) and F (b) with the relative concentration φ/φ max using various formulas and assuming R e * = 0.1 and a constant gelling concentration C gel = 100 g/l for all the tests.

  8, which confirms the results from the modified Richardson & Zaki, Camenen, and Dankers & Winterwerp formulas. Using a particle Reynolds number R e * = 0.2 (median values of the collected data) and C gel = 100g/l, the curves obtained from the modified Richardson & Zaki, Camenen and Dankers & Winterwerp equations are plotted together with the data in Fig. 8. As observed before, the Dankers & Winterwerp formula seems to underestimate values for flux density. The Winterwerp and Dankers & Winterwerp formulas also yield a constant value for φ m and φ cr , no matter what the particle Reynolds number is (see Tab. 2).

Figure 8 :Figure 9 :

 89 Figure 8: Flux density function (Φ/W s0 = φ f ) as a function of the relative concentration φ/φ max using the collected data for cohesive sediments with the particle Reynolds number R e * emphasized.

Figure 10 :

 10 Figure10: Estimation of the permeability k as a function of the concentration c using various empirical formulas calibrated by several scientists for the "Sidere" exercice ((a) refers to Eq. 11; (b) to Eq. 12 and (c) to Eq. 13).

Figure 11 :

 11 Figure 11: Estimation of the total settling velocity function W sk (based on Eq. 15 or on Eq. 16) compared to experimental data (a) from Wolanski et al. (1989; C gel = 70g/l, A k = 7 × 10 -9 , B k = 3.4, and ξ = 1 × 10 4 ) and (b) from Dankers & Winterwerp (2007; C gel = 85g/l, A k = 4 × 10 -9 , B k = 3.5, and ξ = 3 × 10 4 ).

Figure 12 :

 12 Figure 12: Variation of the flux density function F with the relative concentration φ/φ max . The flux density function is calculated using total settling velocity function W sk based on Eq. 15 or on Eq. 16 fitted with the experimental data of Dankers & Winterwerp, see also Fig.11b)

  empirical fit for n, solutions may be obtained with the Richardson & Zaki equation as functions of the particle Reynolds number R e * :

	c m =	1 n + 1	=	3/4 5.7 + 0.586 R e * 1 + 0.175 R e * 3/4	(5)
	c cr =	2 n + 1	=	3/4 5.7 + 0.586 R e * 2 + 0.35 R e * 3/4	(6)
	c m and c cr are increasing functions of R e * (0.18 < c m < 0.30 and 0.36 < c cr < 0.60).
	Using a large data set from batch experiments

Table 2 :

 2 Estimation of the φ m , Φ(phi m ) and phi cr values for cohesive sediments with R e * = 0.1 based on experimental data and using hindered settling formulas.

	φ max

Author-produced version of the article published in Continental ShelfResearch, vol. 31,. The original publication is available at http://www.sciencedirect.com/ -doi:10.1016/j.csr.2010.07.003

Author-produced version of article published in Continental ShelfResearch, vol. 31,. The original publication is available at http://www.sciencedirect.com/ -doi:10.1016/j.csr.2010.07.003

Acknowledgement

The manuscript improved considerably as a result of comments by two unknown reviewers (#1 and #3) and by the guest-editor (R. Verney).

Appendices

A. Settling functions

The following settling functions were used in the paper. They are written for cohesive sediments; in case of non cohesive sediments, φ = c.

• [START_REF] Richardson | Sedimentation and fluidisation: Part I[END_REF] : f = (1φ) n ;

• Souslby (1997) : f = [10.36 2 + 1.049(1φ) 4.7 d * 3 ] 1/2 -10.36

[10.36 2 + 1.049d * 3 ] 1/2 -10.36

where d * = [(s -1)g/ν] 1/3 d and d grain diameter;

where A, B, and N are coefficients for the settling formula, which are function of the grain shape and roundness (see [START_REF] Camenen | Simple and general formula for the settling velocity of particles[END_REF]; for natural sands, A = 24.6, B = 0.96, and N = 1.53; for flocs, A = 27, B = 2.1, and N = 1.2) and subscript mix yields that the kinetic viscosity and density were calculated for a mixture (water + sediments).

B. Table of notation

The following symbols are used in this paper : A, B

[-] : coefficients c

[-] : volumetric concentration of matter c bed [-] : volumetric concentration of matter in the bed for a batch experiment c mud , c sand [-] : volumetric concentration of matter for mud and sand, respectively

[m] : sediment grain diameter for which k% of the grain by mass is finer


Alexis Michaud

Christine Mathieu

Cindy Ho

Pictographs and the language of Naxi rituals

Introduction

In the field of Sino-Tibetan studies, there are few languages with a long-standing written tradition. Of these, the Naxi pictographic tradition encapsulates unique information about the Naxi and their language and holds special promise for research.

The language of the Naxi rituals raises a range of issues, such as: How old are the characters of the Naxi script? What is the origin of the seemingly strange words and turns of phrase found in the rituals? The philological study of the Naxi tradition is greatly complicated by the fact that Naxi books were passed from one generation to the next as mnemonic summaries of the rituals rather than complete transcriptions, unlike Tibetan or Chinese texts. The absence of standardization of Naxi texts allowed the Naxi priests some freedom when copying books; this resulted in great diversification.

The field of linguistics can contribute some evidence and provide some tools to address these complex topics. The approach adopted here consists in looking at the Naxi facts in the light of a comparison between several dialects of the Naxi language and, beyond Naxi proper, a comparison with other Sino-Tibetan languages that are closely related to Naxi.

Naxi is the largest and best-described language of the Naish branch of Sino-Tibetan languages. Two other languages of this branch, spoken in areas that do not have pictographs, will be mentioned: Na and Laze, spoken in Yongning and Muli, respectively (these two language names are printed in capital letters on the map). The analysis of detailed examples reveals how this comparative work sheds light on seemingly obscure aspects of the language used in the Naxi rituals.

Naxi books and the Naxi language

The language of the Naxi rituals is close to the ordinary language, what is known as colloquial Naxi. This is due to the fact that, strictly speaking, there was no clergy among the Naxi: the practitioners of the Naxi religion, the Naxi priests, did not constitute a separate social class. i Monks at a Tibetan Buddhist monastery learned written Tibetan and performed rituals in Tibetan, irrespective of their native language; by contrast, Naxi priests did not need to study another language in order to learn the rituals, even though they did need to learn by rote some proper names, such as those of mythological creatures, and some unusual terms and turns of phrase. This means that research into the language of the rituals can build on the research being done on colloquial Naxi and vice versa.

Anthropologists and linguists have a common interest in transcribing large amounts of texts. The transcription of oral renderings of Naxi texts has been a concern of researchers ever since these texts came under scholarly scrutiny at the turn of the twentieth century. As a result, philological and linguistic work progressed hand in hand, sometimes as a coordinated effort, such as when the linguist Chang K'un (张 琨) worked out the sound system of Naxi to devise a transcription system for recording the pronunciation of Naxi characters. ii The dictionaries by Fang Guoyu and He Zhiwu and by Joseph Rock constitute important lexicographic references as well as importance resources to understand the Naxi writing systems: the pictographic system-called /to˧mbɑ˩/ in Naxi (hereafter to-mba) and Dongba (东 巴) in Chinese-and the syllabic system. iii To this day, annotated editions of Naxi ritual texts actually constitute the largest available translated corpora in Naxi. In addition to collections such as that by Fu Maoji, which predate the official definition of a transcription for so-called "Standard Naxi" in the 1980s, there now exists a massive, one-hundred-volume edition of annotated and translated pictographic manuscripts. iv However, anthropologists and linguists diverge in the amount of attention that they pay to fine details in pronunciation. To the linguist, these details provide crucial hints about the history of the language and the rituals. But first, some essential notions about the Naxi language should be mentioned. Like all of the world's languages, Naxi has various dialects. The city of Lijiang can be considered the center of the Naxi-speaking area; in Chinese scholarship, its dialect is called "Standard Naxi." As one moves away from Lijiang, there are differences from place to place in the pronunciation of the vowels, consonants, and tones; there are also some differences in the lexicon and in the grammar. The point of view adopted in Chinese scholarship consists of using the phrase "Naxi dialect" for all the language varieties closely related to Naxi-although quite a few of these language varieties are not called Naxi by the speakers themselves and are actually so distant from it that they are not intelligible to Naxi speakers from Lijiang. Chinese linguists distinguish two dialect areas: Western and Eastern. v The Western area is relatively homogeneous linguistically; throughout this area, the name of the language (autonym) is "Naxi" (the exact transcription is /nɑ˩hi˧/). (In this article, the phrases "Naxi language" and "Naxi dialects" will be restricted to this area, following the practice of the speakers themselves.) This area largely corresponds with the territory controlled by the feudal lords of Lijiang: the plain of Lijiang had semi-independent rulers as early as the tenth century and until the eighteenth century. vi The relative homogeneity of Naxi dialects is due to historical causes. The centralization of power in Lijiang created a degree of linguistic convergence, and the enforcement of conscription presumably played a role in leveling dialect differences: it is generally observed that nonstandard dialects tend to be stigmatized or derided-typically, the dialect of a major city is more prestigious than that of smaller cities and villages-and that conscripts are under strong social pressure to adopt the most prestigious dialect. As for the languages referred to in Chinese scholarship as "Eastern Naxi," they are so diverse that it appears more adequate to refer to them as distinct languages, rather than as dialects of a single language-keeping in mind that there is no clear boundary between a language and a dialect. These languages, referred to in Chinese scholarship as "Naxi dialects," will be referred to here as "Naish languages." vii It is of paramount importance to linguistic research to possess data on various Naish languages, not on one variety only. Minor differences in pronunciation within the Naxi area proper are best studied through comparison with other Naish languages. Suppose, for instance, that a to-mba priest from the village of Nda-le (to the south of Lijiang; see map) pronounces "horn" as /kʰɑ˧/ whereas a priest from Old Town of Lijiang pronounces this word as /kʰo˧/, with a different vowel. This difference in pronunciation is not reflected in the writing; it is solely due to a difference in the Naxi dialects of the priests, who, when reciting rituals, pronounce such common words as they would in everyday conversation. It is tempting for researchers whose interest lies in the Naxi tradition to adopt a standard transcription, overlook dialect differences, and convert all the transcribed data into Standard Naxi. A standard transcription system was indeed developed in the 1950s, selecting the dialect of the city of Lijiang as the norm. viii International Phonetic Alphabet symbols for the sounds of Naxi were selected, and a Romanized script was created. Although this transcription system is not actually used among the Naxi, it has served as a standard for researchers. One consequence is that a transcriber working with a priest whose pronunciation differed from the standard did not actually transcribe that dialect on its own terms: instead, the transcriber converted the sounds into the standard transcription. This choice was made for the sake of convenience: standardization allows the relatively rapid training of large teams of transcribers. Requiring each transcriber to work out the entire phonetic system of the dialect of each speaker would have been unrealistic since this requires solid training in linguistics and a painstaking analysis of the language variety under study.

For linguists, however, the fine details whereby dialects differ from one another are precious clues to the history of the Naxi language. When a dialectal difference is observed, all the words that belong in the same phonetic set are examined systematically. For example, the word "horn" is presented in Table 1 The first two examples, "horn" and "to kill," have the same vowel in Standard Naxi.

However, the fact that they have different vowels in Nda-le Naxi shows that they have two distinct historical origins. It is a basic principle of historical linguistics that sounds are not expected to change at random: the presence of a vowel /ɑ/ in the word for "horn" in Nda-le Naxi is not due to chance. The presence of other examples exhibiting the same sound correspondences, such as "noise" and "hole," puts to rest any suspicion that the vowel correspondence o:ɑ for "horn" is due to some unaccountable anomaly.

"Horn" and "to kill" must therefore be assumed to have had different pronunciations at a time before Standard Naxi and Nda-le Naxi split apart (due to historical events such as migrations). If only one variety of Naxi were documented, these pieces of historical information would be irretrievably lost. As more dialects of Naxi are described, linguistic investigation can proceed further. In the case of the data in Table 1, crucial evidence comes from another Naish language: the language spoken in Yongning- Evidence of the age of to-mba characters How old are to-mba characters? This is not a simple question. According to Fang Guoyu and He Zhiwu, the Naxi pictographs developed more than one thousand years ago; Guo Dalie and He Zhiwu hypothesize that major texts of the written Naxi tradition date from the Song dynasty, 960-1279. On the other hand, A. Jackson provocatively argues that this tradition only flourished much later, in the 18 th and 19 th centuries. xi Some hints to the answer could in theory be gleaned from phonetic components within the characters.

These phonetic components were developed to remedy the inherent limitations of pictographic writing systems: not many objects, persons, and actions can be represented in a simple pictograph. Pictographic systems therefore tend to be enriched through the combination of two or more basic pictographs, typically using one as a clue to meaning and one as a clue to pronunciation. This is especially well studied in the case of the Chinese writing system. xii To take an example in the Naxi script, "village" (pronounced /mbe˧/) is represented by the character for "house," which is a highly stylized representation of a house (its roof and timber work and its two side walls), to which the symbol for "snow" is added, drawn inside the "house" character. The "house" pictograph is used for its semantic value: "house" and "village" belong in the same semantic field. The symbol for "snow" serves as a phonetic cue: "snow" is pronounced /mbe˧/ and is homophonous with "village." Also, a pictograph is sometimes used simply for its phonetic value, without any visual signal that it is being used for a pronunciation cue and not for its original meaning. For instance, the characters for "monkey" and "life" are both pronounced /y˩/; the pictograph for "monkey" can be used to stand for "life." xiii This is of special interest in assessing the age of to-mba characters. In the history of the Chinese writing system, simple characters were likewise used for their phonetic value; adding a semantic clue created a new character, and the resulting sets of characters are extremely useful to understand the history of the Chinese language. xiv Since the time when the characters were created, the Chinese language underwent considerable phonetic changes, but the characters remained basically unchanged. The phonetic components are therefore unreliable as a cue to the present-day pronunciation of the characters-hence the staggering complexity of Chinese characters for present-day Chinese language learners, who must learn characters whose initial phonetic motivations have long disappeared. For linguists, however, Chinese characters are of great interest because they provide evidence as to which words had similar pronunciations in ancient Chinese. Applying the same method to the Naxi pictographs, characters which share the same phonetic element call for a closer examination: they could potentially provide hints on sets of homophones or quasi-homophones from when the Naxi pictographs developed. Returning to the example above, does the use of the pictograph for "monkey" to mean "life" in some Naxi texts demonstrate that these words were already pronounced in the same way (or in a very similar way) centuries ago?

The answer to this question is somewhat disappointing: almost all the pictographs used for their phonetic value are transparent, that is, they are used for words that have the same consonants and vowels in the present state of the language (tone is simply overlooked). xv This leads to the hypothesis that the Naxi writing system underwent a continuous historical evolution: unlike the Chinese writing system, which remained largely unchanged since its official codification, Naxi pictographs were not standardized, and priests modified them where they found it necessary. When a priest copied a character in which the phonetic element did not coincide neatly with the actual pronunciation of the intended word in his own speech, practical considerations would encourage him to substitute another, more adequate phonetic element. These decisions were not affected by the forces of tradition since the script was not standardized.

This leads to a question for future research: how did the Naxi script evolve from the earliest Naxi manuscripts? The study of the evolution of the characters could reveal changes in the sound system of Naxi and could provide insights into the Naxi language spoken several centuries ago. An ideal basis for addressing this question would consist in studying a large sample of manuscripts, from the oldest (perhaps from the Ming period, 1368-1644) to twentieth-century manuscripts. Using manuscripts of the same rituals would facilitate this type of philological research.

The origin of "strange" words and turns of phrase in to-mba rituals Some words are specific to Naxi rituals: they do not exist in colloquial Naxi and are unintelligible to people who have not been trained in the practice of rituals. Table 2 provides some examples. Linguistic tools can be extremely useful to interpret these words. For example, the word /le˩/ for "moon" has equivalents (cognates) in two other Naish languages: Yongning Na and Laze, a Naish language spoken in the county of Muli, to the northeast of Lijiang Municipality. xvi In Yongning Na, "moon" is pronounced /ɬi˧mi˧/; in Laze, it is pronounced /ɬie˧mie˧/. The second syllable in these words can be analyzed as a suffix, like the /mi˧/ in the colloquial Naxi words for "sun" and "moon." Leaving aside this suffix, and only considering the consonants and vowels of the main syllable, the correspondence between the language of Naxi ritual, Na, and Laze is /le:ɬi:ɬie/. The vowel correspondence /e:i:ie/ between Naxi, Na, and Laze is very common, as is the consonant correspondence /l:ɬ:ɬ/. The word /le˩/ in ritual Naxi cannot be a recent borrowing from Yongning Na or Laze, or the vowel in Naxi would be the same as in the donor language. The presence of the vowel /e/ reveals that this is an archaic word. xvii

The conclusion to be drawn from this comparison is that the word for "moon" used in Naxi rituals, /le˩/, is a word of some antiquity, which has been replaced in colloquial

Naxi by the term /he˧me˧/. The word /le˩/ has disappeared from colloquial Naxi whereas it has been preserved in the language of rituals. At the time when the word /he˧me˧/ was coined, it was probably perceived as more familiar or vulgar than the older word, and the priests therefore continued to use /le˩/, even after this word vanished from colloquial Naxi.

As for the word "chicken," /hjɤ˩˥/, it has long been analyzed as a Tibetan borrowing, xviii though the exact conditions of borrowing need to be worked out (the classical Tibetan word for "chicken" is bya; it was brought into Naxi via a Tibetan dialect which remains unidentified). Quite a few notions of the to-mba religion are borrowed from Tibetan Buddhism along with numerous proper names. It is somewhat surprising that Tibetan In addition to individual words (typically nouns) in the rituals that differ from the everyday language, there also exist special turns of phrase (morphosyntax). A literal translation could be: "A white egg, a green egg, a black egg, a yellow egg, a red egg: there appeared five kinds (of eggs), but there was no hatcher (to hatch them)."

The word /hi˧/ means "person" or "human being" in colloquial Naxi; in this passage, it has a grammatical function, referring to the agent hatching the eggs. It means "that which hatches eggs" rather than "a person who hatches eggs." This point is made clear by the Chinese translation, 没有什么 东 西来孵化它, which means "there was no thing to come and hatch them" and not "there was no one to come and hatch them."

In present-day Naxi, this construction does not make sense: one would need to use the possessive /gɤ˧/ to convey this meaning. However, looking at languages closely related to Naxi, this strange construction becomes clear. The nominalizer is actually pronounced in the same way as the word for "man" in the two Naish languages mentioned earlier, Yongning Na and Laze. These two languages do not appear to have been in contact with each other in the past few centuries; the presence of the same nominalizer in both languages is therefore probably due to a retention from Proto-Naish, the reconstructed common ancestor of Naxi, Na, and Laze, and not to a borrowing from one language into another. In this light, Naxi can be hypothesized to have had the same type of grammatical construction and to have lost it in the course of history; its survival in some to-mba rituals can safely be interpreted as an archaic feature. Thus, linguistic evidence leads to an unambiguous conclusion about the origin of this peculiar linguistic trait: it is not borrowed from another language, and it is archaic. Such analyses open new perspectives for research into the historical depth of Naxi rituals by providing tools for establishing their relative chronology.

Concluding perspectives

Studying the Naxi books and the Naxi rituals with increasing philological precision is a formidable challenge; at stake is the interpretation of Naxi culture and the deconstruction of widespread clichés. The Naxi pictographic texts are commonly considered in present China as a homogeneous, self-contained body of traditional lore.

This point of view is reflected in the phrase "Dongba culture" (东 巴文化), xxiv which identifies the pictographic writing system as the defining characteristic of "traditional Naxi culture." xxv This presents the culture of the Naxi as an idealized heritage from the distant past that has somehow survived into the present age. This is in keeping with the expectation that ethnic minorities should be steeped in tradition, preserving qualities of authenticity and peculiarity. xxvi In contrast to this conception of "Dongba culture" as fixed and timeless, philology and linguistics approach texts from quite the opposite perspective. There are numerous parallels between the philological-linguistic approach and the ethnohistorical approach. To the ethnohistorian, there is not one fixed Naxi identity but an ever-continuing process of construction and reinterpretation of identities in response to changing historical circumstances and opportunities. xxvii Likewise, to the linguist, there is not one Naxi language but a complex landscape of interacting languages and dialects-each dialect lending itself, in turn, to an analysis of a variety of speaking styles (idiolects). Philology and linguistics can therefore be strong allies of history in its arduous task of analysis and reconstruction. Learning the origins of new

  referred to by different authors as Yongning Na, Mosuo, or Moso. ix In Yongning Na, the words that have a vowel correspondence o:ɑ between Standard Naxi and Nda-le Naxi all have an uvular initial consonant: an initial pronounced further back in the mouth than the velar /k/. It is plausible that these uvular consonants are of some antiquity in the Naish languages. The crucial evidence comes from a highly conservative language of the Sino-Tibetan family, which is distantly related to the Naish languages: the words that have an uvular consonant in Yongning Na also have an uvular initial in Rgyalrong, a language spoken in Sichuan, China. x The interpretation is that, in the course of the historical evolution which led to Standard Naxi, the uvular consonants merged with the velar consonants without leaving any traces; in Nda-le Naxi, they also merged with the velars, but before they disappeared, they influenced the pronunciation of the following vowel in such a way that the vowel *a following uvulars and the vowel *a following velars had different evolutions. (Note that the asterisk, *, is used to indicate reconstructed forms, i.e. hypotheses about past forms, as opposed to forms found in the Naish languages as we observe them today.) Therefore, dialectal differences are essential for research; in this respect, the official creation of a standard transcription for Naxi was a mixed blessing. It would be a valuable resource for linguists to have transcriptions of rituals as chanted by priests from the widest possible range of places within the Naxi area. Language documentation projects should ideally include the collection of the rituals that remain alive in Naish language areas without a written tradition, in particular the counties of Muli (木里), Yanyuan (盐 源), Yanbian (盐边 ) and Mianning (冕宁). If substantial numbers of transcriptions can be produced before these traditions fade from memory, this will yield precious historical insights about the evolution of the Naxi language. Conversely, linguistic research may provide useful insights into Naxi rituals.

  borrowings are mostly restricted to proper names and religious notions given the influence of Tibetan culture in this area, which borders on Tibet. This may be partly accounted for by the limitations of Naxi priests' command of Tibetan; however, another factor may be the way in which Tibetan influence was looked upon by the rulers of the Naxi area and by the Naxi priests. Language is very important in religion, as abundantly documented in the history of Christianity and the use of Latin. It is clear that Buddhist monks trained in Tibetan rituals actively resisted the local language and local cultural practices, which were perceived as threats to orthodoxy. A monk in the village of La-bai (拉伯) who translated Tibetan rituals into the vernacular-using the Tibetan alphabet to transcribe his native language (a Naish language close to Yongning Na)-was persecuted by the monastery authorities with such violence that he was driven to madness. xix While the pressure of orthodoxy was probably less intense within the to-mba religion, concerns of political orthodoxy no doubt played a role in the evolution of tomba rituals. The to-mba religion has sometimes been considered an emanation of folklore and indigenous culture-an age-old primitive religion directly addressing the forces of nature, which survived into the twentieth century as part of popular culture. xx This is only part of the story: anthropological and historical research into the Naxi religion reveals the ties between the to-mba rituals and the political projects of the Naxi state of Lijiang from the late fourteenth to the mid-eighteenth century. The to-mba religion appears to have been used as part of a strategy of territorial expansion, playing, in effect, the role of an official cult. xxi Time and again in the history of the Naxi, the influence of Tibet was a sensitive issue. The creation of family names based on Chinese characters (such as Mu [木] for the ruling clan and He [和] for commoners) is one significant episode in the sinicization of Naxi culture.xxii The pressure to adopt Chinese family names makes itself felt to this day in areas where names of Tibetan origins are still used, for example in Yongning, where political tensions lead some local people to insist (in the face of linguistic evidence to the contrary) that their names are purely local and to deny their Tibetan origin. The observation by linguists that the proportion of Tibetan words in Naxi rituals is relatively low may be interpreted as evidence of resilience to Tibetan influence, which is worth investigating further from a historicalanthropological point of view.

Table 1

 1 as part of the entire phonetic set to which it belongs.

	Standard Naxi	Nda-le Naxi	Meaning	Correspondence
	(Lijiang Old			between vowels
	Town)
	kʰo˧	k ʰɑ˧	'horn'	o:ɑ
	kʰo˥	k ʰo˥	'to kill (animal)' o:o
	kʰo˧	k ʰɑ˧	'noise'	o:ɑ
	kʰo˧	k ʰɑ˧	'hole'	o:ɑ
	lɑ˩jɤ˥kʰo˧	lɑ˩jɤ˧kʰɑ˥	'sleeve'	o:ɑ
	o˩	ɑ˩	'turquoise'	o:ɑ
	o˧	ɑ˧	'trunk (of tree)' o:ɑ
	ʂo˧ɭo˧	s ɑ˧lɑ˧	'bone'	o:ɑ
	ho˩	ho˩	'rib'	o:o

. Some correspondences between two Naxi dialects: Lijiang Old Town (defined as the standard variety of Naxi in China) and Nda-le (in the Lijiang plain).

Table 2 .

 2 Some words that differ in the language of the rituals and in colloquial Naxi.

	Meaning	Naxi	Learned	Page number in	Page number in
		(colloquial)	reading of to-	Fang Guoyu	Rock 1963-
			mba character	and He Zhiwu	1972
				1995
	moon	he˧me˧	le˩	91	223
	sun	ŋi˧mi˧	bi˧	91	22
	forest	çi˩	bi˧	144	22
	child	zy˥ly˩	by˧	-	26
	chicken	ae˩	hjɤ˩˥	167	191


Alexandre Lebouvier

François Cauneau

Laurent Fulcheri
email: laurent.fulcheri@mines-paristech.fr

2D Axisymmetric Coupled Computational Fluid Dynamics-Kinetics Modeling of a Nonthermal Arc Plasma Torch for Diesel Fuel Reforming

come   L'archive ouverte pluridisciplinaire

Introduction

For many years, many advanced coupled CFD-kinetics models have been developed, mainly for combustion applications (O/C > 2). For instance, the development of homogeneous charge compression ignition (HCCI) engines needs a lot of knowledge on real kinetic reaction happening in the cylinder. These models are highly time-consuming and essentially depend on the size of the kinetic mechanism. For the reforming application (O/C ~ 1), recent CFDkinetics models have been developed to study partial oxidation, steam or autothermal reforming of methane [1][2][3][4] , autothermal reforming of biogas 5 or hexadecane 6 with very simple kinetic mechanism. CFD-kinetics models have been also developed for different reforming reactors: solid oxide fuel cell (SOFC) 7 , membrane reformer 8 , fluidized bed reactor 9 or microreactor 2,4 . Most CFD models developed for plasma torch study are non-reactive models 10,11 . A recent study by Serbin et al. has focused on a CFD-kinetics model of plasma coal gasification 12 .

The literature is very poor concerning the modeling of low current -high voltage reforming plasma torches under reactive conditions 13 while experimental plasma-assisted reforming has been widely studied [14][15][16][17] . Benilov et al. has developed a low current plasma discharge model in atmospheric pressure to study the non-equilibrium effects 18 . Most CFD models developed for plasma-assisted reforming processes were mainly dedicated to the flow inside the torch [START_REF] Bromberg | CFD modeling of Plasmatron Methane Reformer[END_REF][START_REF] Chirokov | Proceedings of the 15th International symposium on plasma chemistry[END_REF] . Chemical reactions were not taken into account. In the first paper, the plasma is not considered and the study has been focused on the methane/air mixing at the location of the plasma. In the second one, the plasma is considered as a source term of power and the reverse vortex flow influence on the plasma is investigated. Most of the other plasma-assisted reforming numerical models are based on chemical kinetics modeling with CHEMKIN package [START_REF] Kee | Sandia National Laboratories Report[END_REF] , in 0D or 1D. A comprehensive review of these models is given in Ref. 13.

For several years, intensive researches have been dedicated at CEP, both experimentally and numerically, on reforming processes for two different applications (i) onboard production of hydrogen for fuel cell powering and (ii) NO x trap regeneration, using a nonthermal arc discharge torch. Let us point out that these cases correspond to very different reacting conditions. In the first case, diesel fuel reacts with air while in the second case it reacts with diesel engine exhaust gas which is much less oxidative than air. Experimentally, we have studied the plasma-assisted reforming of gasoline [START_REF] Rollier | [END_REF] , ethanol and E85 23 for H 2 fuel cell feeding, and the plasma-assisted exhaust gas fuel reforming of diesel fuel for a NO x trap regeneration application 24 . Indeed, using a plasma to produce onboard reductants, such as CO and H 2 , is an alternative to the classical reforming catalysts. During a regenerating period, the syngas is injected to the NO x trap catalyst to convert the NO x into N 2 . A comprehensive description of the experimental bench and the application are given in Ref. 24.

Thermodynamics model 25 , multistage kinetics model 26 and non-reactive MHD model have been devoted to plasma-assisted syngas production. These models suffer from drawbacks. The thermodynamics and multistage kinetics models consider the perfect mix of reactants and do not consider the thermal losses. The MHD model considers only air as plasma gas because, so far, it is not possible to incorporate chemical kinetic mechanism into MHD models due to the high complexity of involved mechanisms and the lack of basic kinetic data.

In order to move forward a completely integrated model, we present in this paper the simulation of a reforming plasma torch under reactive conditions. This model allows getting significant information about the mixing zone between cold and hot gas.

The plasma-assisted diesel fuel reforming model has been studied for two different applications: (i) onboard H 2 production for fuel cell feeding and (ii) NO x trap regeneration application. Both applications have been investigated under adiabatic and nonadiabatic conditions. The results of the nonadiabatic model have been compared to the experimental results in terms of syngas production. The experimental setup and results have been detailed in Ref. 24.

Mathematical Model

Assumptions. The two-dimensional axisymmetric model studied is based on the following main assumptions:

-The plasma is modeled as a homogeneous high-density power source under thermodynamics equilibrium and thus the non-equilibrium phenomena are not taken into account.

-The kinetic mechanism does not consider the excited or charged species.

-Gravitational effects are not taken into account.

-The gas is treated as semi-compressible. The mass density depends on temperature.

Governing Equations. The stationary conservation equations can be written in the generalized form as:

S v ) ( ) (     , ( 1
)
Where ρ is the fluid mass density, v the velocity vector, Γ Φ the diffusion coefficient, S Φ the source term and Φ represents the scalar variable that must be solved in the various conservation equations. Φ, Γ Φ and S Φ are given for each conservation equation in Table 1. u, v, w, μ l , μ t , P are x-, y-and z-velocity, laminar viscosity, turbulent viscosity and pressure, respectively. h, λ, C p , Pr are enthalpy, thermal conductivity, specific heat and Prandtl number, respectively. X i is the mass fraction. k, P k , ε, C 1 and C 2 * are turbulent kinetic energy, production rate of turbulent kinetic energy, dissipation rate, 1.42 (model constant) and model non-constant coefficient, respectively. The viscosity and the specific heat depend on temperature. The thermodynamic and transport properties have been taken from the Chemkin thermodynamic [START_REF] Kee | The Chemkin Thermodynamic Data Base[END_REF] and transport databases [START_REF] Kee | The Chemkin Transport Database[END_REF] .

The equations have been solved using the commercial CFD code Fluent v. 12.1. The SIMPLE algorithm and Pressure-based solver have been used. The ReNormalization Group (RNG) k-ε turbulence scheme has been chosen. The RNG k-ε model is a classical turbulence scheme used for swirl-dominated combustion flows. The RNG k-ε model considers several scales of turbulence, contrary to the standard the standard k-ε model which considers only one scale of turbulence. It allows guaranteeing that turbulent mixing is computed at the finest possible case regarding the mechanical conditions. As we are in the case of a rotating flow, we considered using the RNG formalism of the k-ε model as it is known to give better results in this case [START_REF] Yakhot | [END_REF] . Both schemes are not different, as they both use a closure on the same two variables. We did not compared the quality result for the two schemes, as there is no evidence in the literature showing results should be degraded in our case.

The turbulence-chemistry interaction model is based on the Eddy Dissipation Concept (EDC).

It assumes that reactions occur in the small turbulent scales, called the fine scales. We know that the EDC model may be influenced by the mixing conditions: mixing time and the fine scale fraction at the cell level. We know that specific works show the influence of this effect 30 , but our goal is to evaluate, at this step, the limits and needs of classical computing schemes.

The radiations emitted by the plasma have not been taken into account because, in a developed MHD model, assuming an optically thin plasma, we have shown that the arc radiations are negligible in the energy equation compared to Joule heating and convection (about 4 orders of magnitude).

Computational grid and boundary conditions. The torch geometry is presented in Figure 1 and corresponds to a real experimental set-up, shown in Figure 2. It is composed of two separated zones: the plasma zone and the postdischarge zone with 70 mm / 4 mm and 500 mm / 11 mm length / inner radius respectively with a divergent between both zones. Let us point out that the plasma zone is the part where the arc plasma really takes place whereas the postdischarge zone is an active or passive zone, located downstream the plasma zone where, reactions ignited in the plasma zone continue to take place. The whole reactor domain has been meshed to obtain significant information about the entry zone, the reactive plasma zone and the mixing zone. The postdischarge has also been entirely meshed to obtain the syngas production reached at the exit in order to compare them to experiments. The grid mesh is exclusively composed of hexa-elements and contained 13 456 cells. The mesh is structured and exclusively composed of rectangles elements. In the nozzle, the mesh cells size is 0.2 mm x 0.2 mm. The number of cells is 20 in the radial direction in the whole torch. In the postdischarge zone, the mesh is released in the axial direction. The mesh is not refined at the wall.

A sensibility study has been made on the mesh density and showed that our mesh was refined enough. The plasma source, modeled by a source term of power, is implemented by means of User Defined Function (UDF) in a limited zone of 1 mm / 70 mm radius / length respectively in the center of the nozzle.

The boundary conditions are detailed in Table 2. The atmospheric pressure is imposed at the outlet. Experimentally, the diesel fuel is vaporized externally by a plug flow system and then mix with the exhaust gas before entering the torch. The mix temperature measured by a thermocouple before entering the torch is about 573 K, and is used as the inlet parameter of the model. A swirl injection is applied at the inlet. The swirl allows the wall-stabilization of the plasma, the rotation of the arc root to avoid fast destruction of the electrode, and allows covering a higher reaction volume. The swirl intensity has not been evaluated.

Experimentally, the flow is injected tangentially to the anode. That corresponds to an infinite swirl number. In the system, we injected a small component of the axial velocity in order to run the simulation and stay close to the experimental parameters. After the bottleneck, the swirl intensity (tangential velocity over axial velocity) varies between 1 and 2. Thus the bottleneck leads to a strong decrease in the swirl intensity, and the boundary conditions before the bottleneck will not slightly influence the results. The injected species for both cases studied are detailed in Table 3. The first case corresponds to the partial oxidation reaction, which is highly exothermic:

-1 2 2 16 7 kJ.mol 586 - = H ; H 8 + CO 7 7O + H C , (2) 31
The second operating condition corresponds to a real diesel engine exhaust gas for a low engine load. The high content of H 2 O and CO 2 in the plasma gas (up to 5 %) lead to energetically unfavorable conditions compared to partial oxidation commonly studied. Indeed, in addition to PO x , the steam and dry reforming of n-heptane can occur and are significantly endothermic reactions. These reactions are represented by the two equations hereafter:

- 33 which is very similar to conventional diesel fuel. This choice has been detailed in Ref. 24.

The plasma-assisted reforming process is far from combustion processes due to low amount of oxygen, O/C ~ 1 and O/C > 2 for reforming and combustion, respectively. Furthermore, due to the presence of H 2 O and CO 2 , the conditions of exhaust gas fuel reforming are very far from those of partial oxidation. However, the relevance of a combustion kinetic mechanism has been proven by Benilov and Naidis for low current arc discharge 34 .

To avoid large computational times, a reduced mechanism of n-heptane, developed for combustion, has been chosen. However, reduced mechanisms are valid on a narrowest conditions range with respect with detailed mechanism. In a first step, several mechanisms 35- 40 have been first tested numerically for both PO x and exhaust gas fuel reforming conditions :

Results and discussion

First, the effect of the inlet composition has been investigated in terms of temperature, syngas production and turbulence for an adiabatic reforming torch. Then, the effects of thermal losses on kinetics evolution and temperature have been studied and compared with experimental results. N 2 is never shown for convenience reasons.

H 2 fuel cell feeding application.

A first inlet composition, corresponding to the PO x condition, has been studied in order to produce syngas which can be used in fuel cells, after a treatment unit. The input power is 1000 W. The thermal losses are not considered.

Figure 3 displays the H 2 molar fraction in the fluid domain in which we can observe a delay in H 2 production in the nozzle because of the significant temperature gradient in the radial direction. In this zone, the torch cannot be considered as a plug flow reactor. Once the flow is well mixed, this delay in conversion is not very important due to lower temperature on the axis. At the outlet, the temperature and species gradients in the radial direction are relatively low and the H 2 concentration is homogeneous in the postdischarge zone. So the results represented along the axis allow studying the phenomena occurring in the plasma zone and are representative of the postdischarge sections. Figure 4 represents the evolution of main species and temperature along the reactor axis and Figure 5 the minority species present in the first 150 mm of the torch. We can distinguish 4 distinct zones. The first zone corresponds to the heating of the reactants (C 7 H 16 and O 2 ). In the very beginning of the plasma zone, a strong mixing is observed due to the geometry of the plasma torch. In this heating zone, the radical densities do not vary much until an axial distance of 37 mm. Then, C 7 H 16 attains a sufficient temperature to be subjected to pyrolysis. The thermal decomposition of n-heptane leads to a growth of many species and radicals such as CH 3 , H, OH, O which are known as playing a significant role in the reforming reactions.

The second zone begins in the middle of the plasma zone located at 44 mm from the input.

Thanks to the high concentration of radicals, the reforming reactions take place in the plasma zone. In the beginning of PO x reaction, while most of species densities decrease, we can see a bump of CH The third zone, just downstream the plasma zone, is a mixing zone. In the PO x condition, the mixing zone corresponds to the divergent zone. The hot gas coming from the plasma zone and the cold gas which pass around the plasma zone are mixed. The strong mixing downstream the plasma zone is highlighted in Figure 7 The reactions follow the same scheme with a jump of these precursors just before the PO x ignition and then a decrease of these precursors to the benefits of H 2 and CO.

The velocity profiles are shown in Figure 8. The profile near the wall has been plotted at a distance of 0.7 mm from the wall in the plasma zone (y = 3.3 mm). However, in the postdischarge zone, it corresponds to a distance of 7.7 mm from the wall. The gas is highly accelerated till 60 m/s in the injection zone because of the bottleneck. In the plasma zone, the high temperature increase leads to a very low mass density, and thus a high volumetric flow rate and velocity up to 110 m/s. Consequently, the species residence time is very low in the plasma zone. The total average residence time is comprised between 23 ms for the part of the gas going through the plasma to 30 ms for the gas passing next to the walls. In the postdischarge zone, except in the boundary layer, the velocity is homogeneous in radial direction.

NO x trap regeneration application. Figure 9 shows the evolution of main species and temperature along the axis for the NO x trap regeneration condition. The first zone, corresponding to the heating of the gas, is longer due to the high fractions of H 2 O and CO 2 in the inflow and a higher mass flow rate. As a result, the second zone is shorter and leads to a lower temperature peak (4200 K). A part of the plasma energy is absorbed by CO 2 and H 2 O, which are already oxidized species, and have a high heat capacity. As a consequence, C 7 H 16 and O 2 are not totally consumed in the second zone and a high fraction of these species appear in the mixing zone. Downstream the mixing zone, the speed of reforming reactions is lower than in the PO x condition. This phenomenon, already observed in Ref. 24, is due to the low concentration of oxygen in the gas leading to lower temperature.

We can distinguish two different zones in the postdischarge zone. In the first part, the PO x reaction continues occurring, raising the temperature, and therefore accelerating the reaction.

All n-heptane is consumed before O 2 .

Influence of thermal losses and comparison with experimental data.

The postdischarge walls are mainly composed of ceramics mullite C530 which has a density of 2400 kg/m 3 , a specific heat of 753 J.kg -1 .K -1 and a thermal conductivity of 2 W.m -1 .K -1 .

The external temperature is set to 300 K. Convection is neglected. For the PO x condition, the experimental data used for the comparison are taken from not published results. For the second case, the experimental data are taken from Ref. 24. In the latter case, the injected power has been increased to 1200 W because 1000 W are not sufficient to ignite the reforming reactions.

Figure 10 shows the evolution of main species and temperature for the PO x condition. In order to compare both models, with and without thermal losses, this figure displays also the temperature, H 2 and H 2 O fractions for the adiabatic and nonadiabatic reactor. As the thermal losses are only implemented in the postdischarge zone, we do not observe any difference in the plasma zone. The temperature decreases quickly from 2200 to 1100 K downstream the mixing zone. In the end of the postdischarge zone, we observe a slight increase of H 2 O molar fraction leading to a slight decrease of H 2 . The thermal losses lead to a too low temperature to ignite the steam reforming reaction and lead a lower H 2 and CO production.

The results for the NO x trap regeneration condition are shown in Figure 11. In the plasma zone, as the input power applied is higher in the nonadiabatic model, the temperature is higher and the reforming reactions occur quicker. H 2 O and CO 2 reach the same values but the thermal losses lead to a decrease of syngas in the second part of the postdischarge zone.

Consequently, a good insulation is needed to improve performances of the system.

The Table 4 shows the comparison of both cases with experiments. The model predicts relevant syngas fraction rates in both cases. In the PO x condition, the difference in CO 2 production comes mainly from gas nonhomogeneity in the experimental reactor. The phenomenon is weaker in the exhaust gas reforming condition where a significant fraction of CO 2 is injected.

Even if the trends are similar for the nonadiabatic case, we observe a significant shift between the model and the experimental data for the temperature. The experimentally measured temperatures given by thermocouples are several hundreds kelvins lower than the predicted temperatures obtained by the model. One explanation is certainly linked to the fact that the thermocouples do not directly measure the gas temperature but a temperature resulting from the thermal equilibrium of the thermocouples between the convective gas temperature and the radiation exchanges between the thermocouple and the walls which are much colder than both the gas flow and the thermocouple.

Conclusions

We have successfully implemented a coupled CFD-kinetics model in order to study plasmaassisted reforming of diesel fuel in partial oxidation condition and a poor oxidative condition corresponding to two different applications. The plasma torch used experimentally has been modeled. First, the study has focused on the influence of thermal losses in the reforming processes for both cases. Then, the numerical results have been compared with experimental results. Among the most important results to point out, one can quote:

The model shows very good correlation with the experimental data in terms of H 2 , CO and CO 2 outlet fractions.

Four zones can be extracted for each case: a reactant heating zone, a plasma zone, a mixing zone and a postdischarge zone.

The main precursors of the reforming reactions are very reactive species: atomic hydrogen (H), atomic oxygen (O) and hydroxyl radical (OH).

O 2 and the temperature are the key points of the plasma-assisted reforming process.

The higher the inlet oxygen rate, the higher the temperature, the quicker the reforming reactions, the higher the syngas production.

The thermal losses in the postdischarge lead to a more difficult ignition of the reforming reactions and lead to the decrease of syngas production in the postdischarge zone. A good insulation of the postdischarge zone is needed to raise the performances.

The perspectives of the model are the implementation of a three-dimensional model able to better reproduce the swirl flow and the interactions with kinetics. A more detailed kinetic mechanism may be implemented depending on the growth of numerical resources.

  rich condition, atmospheric pressure and low temperature. The kinetic model used the Senkin module of the Chemkin II package which simulates a plug flow reactor (PFR). The model has been run for a long numerical time (~10 8 s) corresponding to the steady state. The results are checked against a thermodynamic model based T&TWinner database41 , with the same inlet molar fractions. The output data of the thermodynamic model are the species molar fractions as a function of temperature. For the outlet temperature of the Chemkin kinetic model, we compare the outlet molar fractions of the main species with the molar fractions of these species for the same temperature from the thermodynamic model.Unfortunately, most of these mechanisms turned to be not adapted to the reforming conditions. The kinetic outlet molar fractions are very far from the thermodynamic model because of the low outlet temperature. The only kinetic mechanism, with less than 50 species and 50 reactions, which gave coherent results in terms of outlet molar fractions for exhaust gas fuel reforming conditions, is presented on Appendix A at the end of the paper. It is taken from Ref. 42 and contains 18 species and 41 reactions. No charged or excited species are considered. The kinetic mechanism does not include NO chemistry. From an energy balance point of view, nitrogen-based species will not intervene. The NOx, NH 3 and CN compound production are of course of interest regarding environmental impacts of gas emissions. Studies made with a GC-MS apparatus have shown that cyanides were in trace amounts but no ammonia was observed.

  by the turbulent intensity. It leads to a significant drop of the temperature and a drop of H 2 and CO molar fractions. We can also observe a jump of H 2 O production (+6.2%) and a slight increase of the CO 2 production (+1%). The steam production comes mainly from thermal decomposition of species such as CH 2 O and from the reaction 40 around the arc where the temperature is around 2200 K. The radicals H, O and OH which were present only close to the axis because of high temperature needed are mixed together with other species which exist only around the plasma. As a consequence, we can observe a bump of most the minority species. As these species have a short life time, they decrease in the end of the mixing zone.The fourth zone, downstream the mixing zone, is the postdischarge zone. The reactions ignited in the plasma zone continue to take place. Along the postdischarge zone, CO and H 2 slightly increase, H 2 O and the temperature slightly decrease, and CO 2 remains constant. The reversed reaction of reaction 3 results in a slight H 2 increase and H 2 O decrease. Indeed, the pre-exponential factor (A) of reaction 3 is very low and the temperature is not sufficient to increase the reaction speed. As a consequence, this endothermic reaction has a long characteristic time which explains the observed trends. At the outlet, the H 2 and CO molar fractions reach 19% and 19.8%, respectively. These results can be compared to the results of the 1D multistage model from Ref. 26 which use a n-octane detailed mechanism. The main PO x reaction precursors are similar: H, O, OH and CH 3 .

C 7 H

 7 16 is successively cracked into smaller hydrocarbons as C 7 H 15 , C 3 H 6 and CH 3 . It results in a slight increase of the temperature due to the exothermic property of oxidation. In the middle of the postdischarge zone, we get a runaway of the PO x reaction which leads to the high increase of the temperature and syngas production. From an axial distance of 0.4 m, we can observe a slight consumption of CO 2 and H 2 O and the temperature decreases. It shows that the dry and steam reforming slightly occurs in this part. A sub-mechanism of the dry reforming is modeled by the reversed reaction of reaction 15 which produces CO and OH. It is known that OH radical has a significant role in the reforming reaction ignition.

Figure 1 :

 1 Figure 1: Scheme of the computational grid (unit: mm).

Figure 2 .

 2 Figure 2. Photograph of the experimental set-up.

Figure 3 .

 3 Figure 3. H 2 molar fraction in the fluid domain. PO x condition. O/C =1. P plasma = 1000 W.

Figure 4 .

 4 Figure 4. Evolution of the temperature and main species along the reactor axis. P plasma = 1000

Figure 5 .

 5 Figure 5. Evolution of the minority species in the beginning of the torch axis.

Figure 6 .

 6 Figure 6. Temperature field in the plasma zone and mixing zone. PO x condition.

Figure 7 .

 7 Figure 7. Turbulent intensity in the plasma zone and mixing zone. PO x condition.

Figure 8 .

 8 Figure 8. Velocity magnitude on the axis and on a line near the cathode wall (y = 3.3 mm).

Figure 9 .

 9 Figure 9. Evolution of the temperature and main species along the reactor axis for the NO x

Figure 10 .

 10 Figure 10. Evolution of the temperature and main species along the nonadiabatic reactor axis

Figure 11 .

 11 Figure 11. Evolution of the temperature and main species along the nonadiabatic reactor axis

  3 and CH 3 HCO which species oxidation helps the PO x reaction to ignite by increasing quickly the temperature. The reaction 38 is ignited by the H specie produced and leads to a chain reaction dominated by reaction 1, reaction 38 and the reversed reaction 3. As a consequence, we can also observe a high increase of very reactive radicals H (up to 10 %), O and OH. A sudden decrease of C 7 H 16 and O 2 , and a significant increase of the temperature are observed. The significant increase of the temperature is due to the very exothermic properties of the partial oxidation reaction. The temperature in the domain is displayed in Figure6. The decrease of C 7 H 16 and O 2 implies a sudden increase of H 2 and CO by radical recombination until the end of the plasma zone but also a small amount of H 2 O. The bump of H 2 O (+2.2%) is mainly due to the reactions involving OH.

Table 1 .

 1 Φ, Γ Φ and S Φ for each conservation equation.

	Conservation equations Φ	Γ Φ					S Φ
	Mass	1	0						0
	Momentum	u, v & w μ l + μ t					P
	Energy	h					t	0
			C	p	t Pr	,	h
	Mass fraction	X i				t		0
			C		Pr
				p		t
	Kinetic energy	k			t			P k -ρε
			l		t Pr	,	k
	Dissipation rate	ε	l		t Pr t ,			k	( C 1	P k	C	2	*	)

Table 3 :

 3 Inlet composition for the two cases studied.

		Case 1	Case 2
	Conditions
		PO x	NO x trap regeneration
	Total reactants mass flow (g/s)	0.75	1.04
	C 7 H 16 (% mol )	5.66	3.85
	O 2 (% mol )	19.81	13.11
	N 2 (% mol )	74.53	74.12
	CO 2 (% mol )	0	4.16
	H 2 O (% mol )	0	4.76

Acknowledgments

The authors gratefully acknowledge Renault SAS for their financial support.