N . R A K O T O M A L A L A J M A R T I N

L T A L O N A N

D D S A L I N

Viscous lock-exchange in rectangular channels

In a viscous lock-exchange gravity current, which describes the reciprocal exchange of two fluids of different densities in a horizontal channel, the front between two Newtonian fluids spreads as the square root of time. The resulting diffusion coefficient reflects the competition between the buoyancy driving effect and the viscous damping, and depends on the geometry of the channel. This lock-exchange diffusion coefficient has already been computed for a porous medium, a 2D Stokes flow between two parallel horizontal boundaries separated by a vertical height, H, and, recently, for a cylindrical tube. In the present paper, we calculate it, analytically, for a rectangular channel (horizontal thickness b, vertical height H) of any aspect ratio (H/b) and compare our results with experiments in horizontal rectangular channels for a wide range of aspect ratios (1/10 -10). We also discuss the 2D Stokes-Darcy model for flows in Hele-Shaw cells and show that it leads to a rather good approximation, when an appropriate Brinkman correction is used.

Introduction

The lock-exchange configuration refers to the release, under gravity, of the interface between two fluids of different densities, confined in the section of a horizontal channel.

This physical process has prompted renewed interest, as a part of the carbon dioxide sequestration issues [START_REF] Neufeld | Modelling carbon dioxide sequestration in layered strata[END_REF]). The top of Fig. 1 shows the initial lock-exchange situation of a so-called full-depth release. The two fluids, initially separated by a vertical barrier (the lock gate), fill the whole section of the tank. When the gate is withdrawn (bottom of Fig. 1), buoyancy drives the denser fluid along the bottom wall, while the lighter one flows in the opposite direction at the top of the channel. The so-called lock-exchange results in the elongation of the interface between the two fluids along the horizontal direction. Different regimes have been reported for the velocity and shape of the elongating interface. The slumping phase refers to the initial regime where inertia dominates over viscous forces, which typically applies for the case of salted and fresh water in a tank. In this regime, [START_REF] Benjamin | Gravity currents and related phenomena[END_REF], and more recently [START_REF] Shin | Gravity currents produced by lock exchange[END_REF] showed that in the presence of a small density contrast (i.e. in the Boussinesq approximation ∆ρ << ρ), the two opposite currents traveled at the same constant velocity. When the Boussinesq approximation does not apply, [START_REF] Lowe | The non-Boussinesq lock-exchange problem. part 1. theory and experiments[END_REF], [START_REF] Birman | The non-Boussinesq lock-exchange problem. Part 2. High-resolution simulations[END_REF], [START_REF] Cantero | On the front velocity of gravity currents[END_REF] and [START_REF] Bonometti | Wall effects in non-Boussinesq density currents[END_REF] showed that the two opposite fronts did travel at constant, but with different velocities. However this interface elongation, proportional to the time, is slowed down at later stages, in the viscous phase, where dissipation prevails over inertia. In the latter regime, the interface elongates as t α , where the exponent α, smaller than unity, may take different values depending on the geometry and the confinement of the flow [START_REF] Didden | The viscous spreading of plane and axisymmetric gravity currents[END_REF]; [START_REF] Huppert | The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface[END_REF]; [START_REF] Gratton | Self-similar viscous gravity currents: phase-plane formalism[END_REF]; [START_REF] Cantero | On the front velocity of gravity currents[END_REF]; [START_REF] Takagi | The effect of confining boundaries on viscous gravity currents[END_REF]; [START_REF] Hallez | A numerical investigation of horizontal viscous gravity currents[END_REF]). In porous media, [START_REF] Bear | Dynamics of Fluids in Porous Media[END_REF] and [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF] predicted an interface spreading proportionally to the square root of time that [START_REF] Séon | Transient buoyancy-driven front dynamics in nearly horizontal tubes[END_REF] observed in a horizontal cylindrical tube. Such a spreading can be quantified with a diffusion coefficient, which reflects the balance between the buoyancy driving and the viscous damping. This coefficient, which depends on the nature and the geometry of the flow, has been computed for a porous medium by [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF], for a 2D Stokes flow between two parallel horizontal boundaries separated by a vertical height, H, by [START_REF] Hinch | Private communication[END_REF] and [START_REF] Taghavi | Buoyancy-dominated displacement flows in near-horizontal channels: the viscous limit[END_REF], and for a cylindrical tube by [START_REF] Séon | Transient buoyancy-driven front dynamics in nearly horizontal tubes[END_REF]. However, to our knowledge, such a diffusion coefficient has not been derived for a rectangular channel (horizontal thickness b, vertical height H, Fig. 1), for which one expects to recover the porous medium regime for b ≪ H, and, possibly, the Stokes flow regime for b ≫ H. In order to gather the limiting cases in the same paper, we first recall the results for porous media and 2D Stokes flows, together with the tube case, for the sake of comparison. Then we compute, for a rectangular channel of aspect ratio, H/b, the dependence of the interface h(x, t) and the corresponding viscous lock-exchange diffusion coefficient. We also test the so-called Stokes-Darcy 2D model to this lock-exchange configuration. Finally, we test and validate our theoretical results with experiments in horizontal rectangular channels for a wide range of aspect ratios (1/10 -10).

Lock-exchange in different geometries

Let us first recall the basic hypotheses on the viscous gravity currents, common to the different geometries, used for instance by [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF] or [START_REF] Hinch | Private communication[END_REF]).

As sketched in Fig. 1, the interface between the two fluids is assumed independent on the y direction. Its distance from the bottom boundary of the vessel is denoted h(x, t).

This interface can be a pseudo-interface between two miscible fluids for which molecular diffusion can be neglected or between two immiscible fluids, provided that the interfacial tension can be neglected. The flow is assumed to be quasi-parallel to the horizontal x axis. This is a key hypothesis. Neglecting accordingly the vertical component of the fluid velocity implies that the vertical pressure gradient follows the hydrostatic variations:

∂p/∂z = -ρg. This hypothesis is violated at short times, immediately after the opening of the gate, but should become valid at later stages, as soon as the interface has slumped over a distance larger than H, thus ensuring a small enough local slope ∂h/∂x. Then, the pressures, P + and P -, in the lower layer, 0 < z < h(x, t), and in the upper one, h(x, t) < z < H, respectively, write

P + = p(x, t) -ρ + gz ; P -= P + + ∆ρg(z -h(x, t)) (2.1)
where p(x, t) denotes the pressure at the lower wall, z = 0. The difference between the horizontal pressure gradients in the two fluids is therefore linked to the interface slope by:

∂P + ∂x - ∂P - ∂x = ∆ρg ∂h ∂x (2.2)
The time evolution of the interface, h(x, t), is governed by the mass conservation of each fluid (see Fig. 1 for notations). For instance, for the heavier bottom layer we have:

∂h ∂t + ∂q ∂x = 0 (2.3)
where q(x, t) = q + (x, t) is the horizontal flux (m 2 /s) of the denser fluid at the location

x: q(x, t) = q + (x, t) = h 0 1 b u x (x, y, z, t)dydz (2.4)
with u x (x, y, z, t) the x-velocity component and b the spanwise length (the y integration is along this spanwise length). Moreover, in our configuration of uniform section along the horizontal axis x, u x (x, y, z, t) must also satisfy the no net flux condition:

q + (x, t) + q -(x, t) = H 0 1 b u x (x, y, z, t)dydz = 0 (2.5)
We will see in the following that, in the viscous regime of interest, the horizontal velocity component u x , solution of either a Darcy or a Stokes equation, is proportional to the pressure gradient in each fluid layer. Such solutions, combined with eq. (2.2), eq. (2.4) and eq. (2.5) allow then to eliminate the pressure gradients and to derive an expression of the flux q, of the form:

q = -Df ( h H ) ∂h ∂x (2.6)
where D writes

D = τ ∆ρg η (2.7)
and where the constant τ (scaling with a volume) and the function f depend on the geometry and the flow equation and η is the dynamic viscosity. Using the expression (2.6) for the flux, eq. ( 2.3) admits a self-similar solution h(ζ) = H ψ(ζ) with the similarity variable ζ = x/ √ Dt, which obeys:

-ζ dψ dζ = 2 d dζ (f (ψ) dψ dζ ) (2.8)
This equation may alternatively be rewritten, in terms of ζ(ψ):

ζ( dζ dψ ) 2 -2 f d 2 ζ dψ 2 + 2 df dψ ( dζ dψ ) = 0 (2.9)
The solution of the above equations can be found analytically or numerically, depending on the complexity of the normalized flux function f (ψ). In the following, we will first recall the case of porous media, treated by [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF] and the 2D Stokes flow, addressed by Hinch (2007) (unpublished) and [START_REF] Taghavi | Buoyancy-dominated displacement flows in near-horizontal channels: the viscous limit[END_REF]. We note that the latter paper included the effects of the rheological properties of the fluids. However, in order to focus on the geometrical aspects, we will assume in the following that both fluids are Newtonian and have the same viscosity.

Lock-exchange in porous media

For a homogeneous layer of porous medium of permeability κ (see for instance [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF]), the flow in each fluid is given by Darcy's law which relates the velocity in each phase to the local pressure gradient :

u x± = - κ η ∂P ± ∂x (2.10)
At a given location x, the velocity is then uniform in each layer, and the no net flux condition (eq. (2.5)) simply writes: hu x+ +(H -h)u x-= 0. The latter equation, combined with eq. (2.10) and eq. (2.2), leads to eq. (2.6), and thus (combined with eq. (2.3)) to eq.

(2.8) with a diffusion coefficient and a flux function:

D P M = κ H ∆ρg η (2.11) f P M (ψ) = ψ(1 -ψ) (2.12)
The solution of eq. (2.8), in the similarity variable ζ = x/ √ D P M t, is then a linear profile [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF]):

ψ = h(x, t)/H = (1 + ζ)/2 (2.13)
The so-obtained front profile in homogeneous porous media is displayed in Using the no net flux condition (eq. (2.5)) and eq. ( 2.2), we obtain eq. (2.6), which enables to rewrite eq. ( 2.3), in the form of eq. (2.8) or eq. (2.9), with:

D = H 3 ∆ρg 3η (2.16) f (ψ) = ψ 3 (1 -ψ) 3 (2.17)
Note that the polynomial development of the solution of eq. ( 2.9) around ψ = 0 gives:

ζ = -ζ 0 + 2ψ 3 /(3 ζ 0 ). Thus, the location, ζ(ψ = 0) = -ζ 0 , of the leading edge of
the interface is indeed constant in the similarity variable. Moreover, the development shows that the slope of the interface is vertical at the bottom wall (ψ = 0). This is also the case at the upper wall (ψ = 1), as the problem is symmetric with respect to the centre of the cell. We note that in the presence of such a vertical slope, our (horizontal) quasi-parallel flow assumption falls locally, but it is still valid upstream and downstream, where the slope of the interface remains small. The solution of eq. (2.9) can be found numerically using a shooting method similar to the one used by [START_REF] Hinch | Private communication[END_REF] and [START_REF] Taghavi | Buoyancy-dominated displacement flows in near-horizontal channels: the viscous limit[END_REF]. It was computed starting the integration of eq. (2.9) from (ψ = 0.5, ζ = 0) and matching the asymptotic development in the vicinity of ψ = 0.

From the so-obtained solution, one can deduce the spreading diffusion coefficient between the leading edge (h = 0, -ζ 0 = -0.1607) and the trailing edge

(h = H, ζ 0 ) of the front, from [0.5 (x(h = 0) -x(h = H))] 2 = D t, which gives D 2D = D ζ 2 0 , so that: D 2D = 0.0086 H 3 ∆ρg η (2.18)
This result is in agreement with the one found by [START_REF] Hinch | Private communication[END_REF]. [START_REF] Taghavi | Buoyancy-dominated displacement flows in near-horizontal channels: the viscous limit[END_REF] provide five ψ(ζ) plots in their Fig. 9, corresponding to different viscosity ratios and including our case. From that figure we may obtain a value of their similarity variable, η 0 ∼ 0.09, which is consistent with our finding ζ 0 = 0.1607, when taking into account their definition of the similarity variable, ζ = η √ 3.

For completeness, D 2D may be compared to the result for a cylindrical tube of diameter d [START_REF] Séon | Transient buoyancy-driven front dynamics in nearly horizontal tubes[END_REF]):

D T = 0.0054 d 3 ∆ρg η (2.19)
The above expression is indeed very close to the 2D result, with d playing the role of H.

Lock-exchange in a rectangular cross-section channel

This article aims to extend the computation of the lock-exchange diffusion coefficient to rectangular cells of arbitrary cross-sections H × b (see Fig. 1). In the following, the crosssection aspect ratio will be denoted Γ = H/b. As in the previous section, a quasi-parallel flow approximation is assumed (i.e. small interface slope) which leads to eq. ( 2.2) for the pressure gradient. We will also assume the invariance of the interface location along the gap direction y. This requires that the deformation of the interface, induced by the flow profile along the direction y, relaxes much more quickly than the gradient along x. The flow in each fluid obeys a 3D Stokes equation:

η ∇ 2 u x± (x, y, z) = ∂P ± ∂x (2.20)
In order to solve this equation, we follow the series decomposition in Fourier modes of the velocity field used by [START_REF] Gondret | Viscous parallel flows in finite aspect ratio Hele-Shaw cell: Analytical and numerical results[END_REF]. This paper addressed the issue of the parallel flow of two fluids of different viscosities in a rectangular cell. This issue is very closed to ours, as it requires to solve the Poisson equation (eq. (2.20)), but with different viscosities and the same pressure gradient for both fluids in [START_REF] Gondret | Viscous parallel flows in finite aspect ratio Hele-Shaw cell: Analytical and numerical results[END_REF]. The method used was to split the velocity into two terms,

u x± (x, y, z) = u * x± (x, y) + u * * x± (x, y, z). Here, the first term, u * x± (x, y) = b 2 8η ∂P± ∂x [1 -( 2y b ) 2
] is the Poiseuille-like unperturbed velocity far away from the interface. The second term satisfies the Laplace equation,

∇ 2 u * *
x± (x, y, z) = 0 and vanishes far away from the interface. Its expression in terms of a sum of Fourier modes leads to a velocity profile of the form: Combining the so-obtained expressions for the velocity with eq. (2.2) and the no net flux condition (eq. (2.5)), one obtains the horizontal flux of the heavy fluid (eq. (2.6)):

u x± (x, y, z) = b 2 8η ∂P ± ∂x 1 -( 2y b ) 2 + ∞ n=1 32(-1) n (a ±n e (2n-1) π(z-H/2) b + b ±n e -(2n-1) π(z-H/2) b ) (2 n -1) 3 π 3 cos[(2n -
q = -D HS f Γ (h/H) ∂ h ∂ x (2.22) with D HS = b 2 H∆ρg 12η (2.23) f Γ (ψ) = ψ + α Γ (ψ) 1 -γ Γ (1 -ψ -α Γ (ψ) -γ Γ ) -δ Γ (ψ) (2.24) α Γ (ψ) = 1 Γ ∞ n=1 96 1 + e Γ(1-ψ) (2 n-1) π 1 -e Γ ψ (2 n-1) π 1 + e Γ (2 n-1) π (2 n -1) 5 π 5 (2.25) δ Γ (ψ) = 1 Γ ∞ n=1 48 1 -e Γ(1-ψ) (2 n-1) π 2 1 -e Γ ψ (2 n-1) π 2 -1 + e 2 Γ (2 n-1) π (2 n -1) 5 π 5 (2.26) γ Γ = 1 Γ ∞ n=1 192 tanh ( Γ (2 n-1) π 2 ) (2 n -1) 5 π 5 (2.
= -ζ 0 + 8 Γ 2 ψ 3 /(3 ζ 0 )
in the vicinity of the boundary, ψ = 0. We solve this equation using the shooting method previously described and using Mathematica Software. The solutions h(ζ) are plotted in Fig. 2 for different values of the cell aspect ratio Γ. We notice that, in contrast with Darcy predictions (straight line in Fig. 2), but similarly to the case of the 2D Stokes flow, the profiles, h(ζ), exhibit vertical slopes at the edges of the cell. We note also that such vertical slopes were observed in the experiments by [START_REF] Séon | Transient buoyancy-driven front dynamics in nearly horizontal tubes[END_REF] and [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF]. When comparing their experiments in a Hele-Shaw cell with Darcy predictions, the latter authors reported that "Some discrepancies develop near the leading edge of the current as a result of the increasing importance of the bottom friction at the nose" (Fig. 2 of [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF]). This mismatch will be addressed in the next section. According to the so-obtained profiles, stationary in the similarity variable ζ, the leading and trailing edges of the front spread as the square root of time, and a lock-exchange diffusion coefficient, dependent on the cell aspect ratio, can be defined: obtained for the 2D homogeneous porous medium case (eq. (2.14)). We note that the latter case, which corresponds to a Hele-Shaw cell of infinite aspect ratio, overestimates lock-exchange diffusion coefficient, by a relative amount of about 30%, for aspect ratios as large as Γ = 10 -20.

D R = b 2 H∆ρg 12η F ( H b ) = D HS F ( H b ) (2.

2D

- 12η b 2 U x± (x, z) + β η∇ 2 U x± (x, z) = ∂P ± ∂x (3.1)
The first term on the left hand side of eq. (3.1) and the pressure gradient correspond to the Darcy's law (eq. (2.10)) with a permeability κ = b 2 12 for the Hele-Shaw cell as mentioned above (eq. (2.14)). The second term on the left hand side of eq. (3.1) is the Brinkman correction to the Darcy equation (see [START_REF] Brinkman | A calculation of the viscous forces exerted by a flowing fluid on a dense swarm of particles[END_REF]), which involves an effective viscosity, βη. This effective viscosity may be taken equal to the one of the fluid (β = 1) for the sake of simplicity (or to enable the matching with a 2D Stokes regime at Γ → 0). However, [START_REF] Zeng | On the Brinkman correction in unidirectional Hele-Shaw flows[END_REF] showed that in the Hele-Shaw cell regime (at large Γ), the effective viscosity was slightly higher, with β = 12/π 2 ≃ 1. 2.2), we obtain the horizontal flux (eq. (2.4)):

q(x) = -D HS f SD (h/H) ∂ h ∂ x (3.2)
where D HS was already given in eq. (2.14) and the reduced flux function is equal to:

f SD (ψ) = 1 4d(d -tanh d) 2 + 4 d 2 (1 -ψ) ψ -d 3 cosh(2 d) + cosh(2 d (1 -2 ψ)) sinh(2 d) +4 d (1 -ψ) cosh(2 d (1 -ψ)) + ψ cosh(2 d ψ) sinh(2 d) -2 cosh(d (1 -2 ψ)) cosh(d) (3.3) where d = H 2 4κβ = 3 β Γ (3.4)
and κ = b 2 /12. A comparison of the full 3D calculations for a rectangular channel of aspect ratio Γ = H/b with this 2D approximation can be performed on the flux functions, f Γ (ψ) (eq. (2.24)) and f SD (ψ) (eq. ( 3.3)). These two flux functions are close to each other, within a few per cents. In order to address the comparison in the range of interest for the Hele-Shaw assumption, i.e. Γ ≫ 1, let us analyze the limit Γ → ∞ (d → ∞), which gives

f SD,Γ→∞ ≃ ψ(1 -ψ) -( 3 4 -ψ(1 -ψ)) β 3 b H + O ( b H ) 2 (3.5)
for the Stokes-Darcy flux and

f Γ,Γ→∞ ≃ ψ(1 -ψ) -( 3 4 -ψ(1 -ψ)) 186 Zeta(5) π 5 b H + O ( b H ) 2 (3.6)
for the full 3D rectangular cell flux (with Zeta(5) = ∞ 1 n -5 = 1.03693, the value of the Riemann-Zeta function). We note that the leading term of both series corresponds to the expected porous media Darcy limit (eq. (2.12)) with a permeability κ = b 2 /12. However, the next order term (O(b/H)) is not the same, unless one chooses for the factor β,

β = 3 ( 186 Zeta(5) π 5 ) 2 ≃ 1.192 (3.7)
which is very close to the value 12/π 2 ≃ 1.215 found by [START_REF] Zeng | On the Brinkman correction in unidirectional Hele-Shaw flows[END_REF] and to the value 6/5 proposed by [START_REF] Ruyer-Quil | Inertial corrections to the Darcy law in a Hele-Shaw cell[END_REF]. The lock-exchange diffusion coefficient has been computed, with the same procedure as above, by integrating eq. (2.9), using data over the whole range of aspect ratios, Γ = 1 -100. The inset of Fig. 4 gives the percentage of error for the two values of β. We note that these data were obtained by the difference between values of accuracy of the order of a few 10 -3 , which results in the small dispersion observed in the inset of Fig. 4. As expected, the results at large Γ (in the Hele-Shaw regime) are closer to the 3D full problem for β = 12/π 2 than for β = 1. We point out that, whereas the Brinkman term does bring a significant correction, the exact value of the Brinkman viscosity factor β is however not crucial: For instance, for Γ = 10, we obtained a diffusion coefficient 3.5% smaller than the 3D value for β = 1 and 0.1% larger for β = 12/π 2 , to be compared to the 30% of error if the cell was assumed to be of infinite aspect ratio (Hele-shaw limit) as in [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF]. In conclusion of this comparison, we have shown that the 2D Stokes-Darcy model for lock-exchange in a rectangular cell captures quite accurately the effect of the finiteness of the cross-section aspect ratio. By using the correct β value, the error in the model is smaller than 5% for aspect ratios larger than Γ > 1.

Experiments

In this section, we will present experimental measurements of the diffusion coefficient in Hele-Shaw cells of different aspect ratios and we will compare them with our computed values.

We used borosilicate rectangular cells of height H and thickness b and typical length 30 cm (Fig. 1). The rectangular cross-sections of the cells were (in mm edge is fuzzy. This can be attributed to the stick condition at the upper wall: The dark dense fluid does stay at the walls for a long time, in particular in the corners of the crosssection. The same phenomenon takes place at the bottom of the cell, but the presence of transparent light fluid has little effect on the turbidity of the heavy dark fluid, and is therefore not noticeable on the pictures. It is worth noting that the shape of the leading edge evolves from an edge at large aspect ratios Γ to a more and more step-like shape as Γ decreases. It should be noticed that for small aspect ratios, although it is rather difficult to take pictures, a top view of the cell reveals a mild spanwise dependency of the interface, but we do not observe the spanwise lobe-and-cleft instability reported by [START_REF] Simpson | Effects of the lower boundary on the head of a gravity current[END_REF]. For each experiment, the locations of the leading and trailing edges of the front were measured in time. with a typical accuracy of 20%. Fig. 7 (top) displays the so-obtained normalized lockexchange diffusion coefficient as a function of Γ. One can see that the agreement with the 3D calculations over the two decades of our measurements is rather good. We note that for the large aspect ratio limit of the experiments (up to Γ = 10), the Hele-Shaw cell limit is not reached, and would underestimate, by 30%, the lock-exchange diffusion coefficient.

This result thus confirms that for such aspect ratios, one should either compute the full 3D Stokes equation or use the Stokes-Darcy model to obtain the correct behaviour.

We also note that our calculation still holds for aspect ratios as small as Γ = 0.1. This result is quite unexpected since for such aspect ratios, some spanwise dependency of the profile was observed, and the hypothesis of the interface surface, h(x, y, z), invariant in the y direction is certainly broken. The bottom of Fig. 7 displays the superimposition of the theoretical and the experimental interfaces between the fluids, for an aspect ratio Γ = 4. The agreement between the two is rather good, thus validating our model. Such an agreement is rather surprising as our small slope assumption is violated at the edges of the gravity current. This agreement, already emphasized by [START_REF] Huppert | The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface[END_REF] and [START_REF] Séon | Transient buoyancy-driven front dynamics in nearly horizontal tubes[END_REF], is likely to be common to viscosity dominated gravity current without surface tension.

Conclusion

The viscous lock-exchange diffusion coefficient reflects the competition between the buoyancy driving effect and the viscous damping, and depends on the geometry of the channel. We give the backbone to calculate this coefficient in different configurations: We recall its computation for a porous medium already found by [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF], This analysis provides a cross-over between the 2D Stokes flow between two parallel horizontal boundaries separated by a vertical height, H, and the Hele-Shaw cell limit (applying for H/b > 100). Moreover, the shape of our profiles allows to account for the discrepancy observed at the nose of the gravity current in the experiments by [START_REF] Huppert | Gravity-driven flows in porous layers[END_REF]. The agreement, obtained despite the failure of the lubrication assumption at the edges of the current, should deserve however further theoretical investigation. Our calculations of the diffusion coefficient and of the shape of the profile have also been convincingly compared to new experiments carried out in cells of various aspect ratios

(1/10 -10). We have also calculated the lock-exchange diffusion coefficient for the same rectangular cells, using the 2D Stokes-Darcy model. This model is shown to apply to aspect ratios H/b > 1, provided that the appropriate Brinkman correction is used. Such a 2D model may be useful to describe gravity currents with a finite volume of release, with fluids of different viscosities, or in heterogeneous vertical fractures.

Figure 1 .

 1 Figure 1. Sketch of the rectangular cell of height H in the gravity direction and width b. Top: Initial configuration where the heavy fluid of density, ρ, is separated from the light fluid of density, ρ -∆ρ, by a vertical gate (plane yz). Bottom: After removing the gate, buoyancy differences cause the denser fluid to flow in one direction along the bottom of vessel, while the lighter one flows in the opposite direction at the top of the vessel. The goal is to determine the space and time dependencies of the pseudo-interface, h(x, t).

Figure 2 .

 2 Figure 2. Lock-exchange interface ψ(ζ) between the two fluids for rectangular cells of different aspect ratios Γ = H/b. The straight line of slope 0.5 corresponds to the Darcy porous media limit (eq. (2.13)). The other curves correspond to the case of a rectangular cross-section channel (subsection 2.3). From the bottom left to right, Γ = 10, 4, 2, 1, 0.5 and 0.2.

Figure 3 .

 3 Figure 3. Log-log plot of the normalized lock-exchange diffusion coefficient, DR/DHS, versus the aspect ratio Γ = H/b, for rectangular cells. DR is normalized by the diffusion coefficient DHS = b 2 H∆ρg/(12 η) (eq. (2.14), obtained for the Hele-Shaw cell limit, also corresponding to a 2D porous medium. Accordingly, the latter regime corresponds to the dashed straight line (DR/DHS = 1) in this representation. The solid straight line corresponds to the 2D Stokes limit (eq. (2.18)), leading to a slope 2 in this representation.

  Fig. 3 displays a log-log plot of the normalized rectangular cell lock-exchange diffusion coefficient, D R /D HS = F (H/b), versus the aspect ratio Γ = H/b. At small aspect ratios, Γ < 1, the diffusion coefficient falls on top of the full line of slope 2, which corresponds to the 2D Stokes flow between boundaries distant of H (b → ∞, eq. (2.18)). At large aspect ratios, Γ → 100, the diffusion coefficient approaches the dashed line, D R /D HS = 1,

  Stokes-Darcy model for lock-exchange in a rectangular cross-section channel The above-mentioned failures of the 2D Darcy model at finite aspect ratios may come from the velocity slip condition at the bottom and top edges of the cell (z = 0 and z = H, respectively). This non physical condition is indeed required by the use of Darcy equation for the flow, which neglects the momentum diffusion in the presence of velocity gradients, in the plane of the cell (zx plane). The momentum diffusion may however be taken into account in 2D, through the so-called Stokes-Darcy equation (see Bizon et al. (1997); Ruyer-Quil (2001); Martin et al. (2002a); Zeng et al. (2003)), which is similar to the Darcy-Brinkman equation used in porous media (see Brinkman (1947)). This 2D model enables to handle discontinuities such as cell edges, gap heterogeneities and fluid interfaces (Ruyer-Quil (2001); Martin et al. (2002a); Zeng et al. (2003); Talon et al. (2003)) and was successfully applied in the study of Rayleigh-Taylor instability (Martin et al. (2002a); Fernandez et al. (2002); Graf et al. (2002)), of dispersion in heterogeneous fractures (Talon et al. (2003)) and of chemical reaction fronts (Martin et al. (2002b)).Although our present case of interest can be handled with 3D Stokes calculations, it is of interest to test the applicability of the 2D Stokes-Darcy model to the case of deep and narrow cells. Indeed, such a 2D model, once validated, could be a useful tool to address the issue of more complicated cases, such as gravity currents in the presence of viscosity contrasts, or in fractures with aperture heterogeneities.In this model, the flow in the rectangular cell (Fig.1) is assumed to be parallel to the plates ( u(x, y, z) = (u x (x, y, z), 0, u z (x, y, z)) with a Poiseuille parabolic profile across the gap (the key assumption). Using the Stokes equation with this y dependency, the gap-averaged fluid velocity U (x, z) 2 u(x, y, z)dy, follows a Stokes-Darcy (SD) equation which reads here for the horizontal component of the velocity:

Figure 4 .

 4 Figure 4. Log-log plot of the normalized lock-exchange diffusion coefficient, versus the aspect ratio Γ = H/b, for rectangular cells. The data points, obtained using the 2D Stokes-Darcy model (eq. (4)) with Brinkman viscosity factors β = 1 (open squares) and β = 12/π 2 (open triangles) are compared to the full 3D results already shown in Fig. 3 and depicted here as a solid line. The inset gives the relative difference between the 2D Stokes-Darcy model (β = 1, open squares and β = 12/π 2 , open triangles) and the 3D calculations. Note that these data were obtained by the difference between values of accuracy of the order of a few 10 -3 , which results in the small dispersion observed in the figure.

  215. At a given location x, integrating eq. (3.1) leads to the two velocity profiles matching the no slip boundary conditions at the bottom and the top boundaries (u x+ (x, 0) = u x-(x, H) = 0) and the continuity of the velocity (u x-(x, h) = u x+ (x, h)) and of the shear stress (β η ∂u x-(x, h)/∂z = β η ∂u x+ (x, h)/∂z) at the interface. Using the no net flux condition, 0 u x+ dz + H h u x-dz = 0, and eq. (

f

  SD (ψ), from ψ = 0.5, and matching the asymptotics, ζ = -ζ 0 + 8 d 2 ψ 3 /(9 ζ 0 ) in the vicinity of the boundary, ψ = 0. The so-obtained lock-exchange diffusion coefficients, calculated for two different values of β (β = 1 and β = 12/π 2 ) are compared to the 3D calculations in Fig. 4. The data for both values of β are indeed very close to the 3D

Figure 6 .

 6 Figure 6. Plot of the square of the spreading distance of the leading edge of the front versus time for five cells: (•) H = 3 mm, b = 30 mm;( ) H = 6 mm, b = 2 mm; ( * ) H = 6 mm, b = 4 mm; (•) H = b = 6 mm; (+) H = 30 mm, b = 3 mm. The solid line is a linear fit to the data, the slope of which gives the diffusion coefficient plotted in Fig. 7 (top).

Figure 7 .

 7 Figure 7. Top: Log-log plot of the normalized measured lock-exchange diffusion coefficient, Dexp/DHS (open squares), versus the aspect ratio, Γ = H/b, of the cell. The solid line corresponds to the 3D results (same as in Fig. 4). Bottom : Superimposition of the measured pseudo-interface between the two fluids (grey fuzzy line) and of the theoretical profile (full dark line) corresponding to Fig. 2. The aspect ratio is Γ = 4.

  and compute it for a 2D Stokes flow between two parallel horizontal boundaries by a vertical height, H. This result is in agreement withHinch (2007) (unpublished) and in reasonable agreement with recent computations by[START_REF] Taghavi | Buoyancy-dominated displacement flows in near-horizontal channels: the viscous limit[END_REF]. Using a quasi-parallel flow assumption, we have calculated the pseudo-interface profile between the two fluids and the diffusion coefficient of viscous lock-exchange gravity currents for a rectangular channel (horizontal thickness b, vertical height H) of any aspect ratio (H/b).

  HS t, which obeys eq. (2.8) or eq. (2.9). As previously, it is easier to compute the solution ζ(ψ) of eq. (2.9) subject to the corresponding asymptotics, ζ

		27)
	Eq. (2.22) admits a self-similar solution, h(ζ) = H ψ(ζ), with the similarity variable
	ζ = x/	√ D

Acknowledgement

This work was partly supported by CNES (No 793/CNES/00/8368), ESA (No AO-99-083), by Réseaux de Thématiques de Recherches Avancées "Triangle de la physique", by the Initial Training Network (ITN) "Multiflow" and by French Research National Agency (ANR) through the "Captage et Stockage du CO 2 " program (projet CO-LINER No ANR-08-PCO2-XXX). All these sources of support are gratefully acknowledged.


Mischa Dohler
email: mischa.dohler@orange-ftgroup.com

Dominique Barthel
email: dominique.barthel@orange-ftgroup.com

Florence Maraninchi
email: florence.maraninchi@imag.fr.

Laurent Mounier
email: laurent.mounier@imag.fr.

Stéphane Aubert
email: stephane.aubert@coronis.com.

Christophe Dugas
email: christophe.dugas@coronis.com.

A Buhrig

Franck Paugnat
email: franck.paugnat@imag.fr.

Marc Renaudin
email: marc.renaudin@imag.fr.

Andrzej Duda
email: andrzej.duda@imag.fr.

M Heusse

F Valois

The ARESA Project: Facilitating Research, Development and Commercialization of WSNs

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

I. INTRODUCTION

The ARESA project is a French national project with the aim to facilitate research, developments and commercialization of wireless sensor networks (WSNs) and embedded systems. It is comprised of a healthy mix of industrial and academic partners, all of which have a wide range of complementing expertise. The project background and aims, as well as related projects and the partners' involvements are exposed in the following.

A. Background

Sensor networks have been researched and deployed for decades already; their wireless extension, however, has witnessed a tremendous upsurge in recent years. This is mainly attributed to the unprecedented operating conditions of WSNs, i.e. a potentially enormous amount of sensor nodes reliably operating under stringent energy constraints.

WSNs allow for an untethered sensing of the environment. It is anticipated that within a few years, sensors will be deployed in a variety of scenarios, ranging from environmental monitoring to health care, from the public to the private sector, etc. They will be battery-driven and deployed in great numbers in an ad hoc fashion, requiring communication protocols and embedded system components to run in an utmost energy efficient manner.

Prior to large-scale deployment, however, a gamut of problems has still to be solved which relates to various issues, such as the extraction of application scenarios, design of suitable software and hardware architectures, development of communication and organization protocols, validation and first steps of prototyping, until the actual commercialization.

B. Related Projects

Numerous past, current and emerging projects worldwide focus on a subset of above mentioned problems, so as to facilitate a successful deployment of WSNs. A non-exhaustive but related set of international projects along with the projects' main goalsis summarized below.

E-Sense: The aim of the E-Sense project is to capture ambient intelligence through WSNs by means of interaction between body sensor networks, object sensor networks and environmental sensor networks.

CRUISE: It is a European Network of Excellence (NoE), which deals with a wide range of scenarios and applications of WSNs.

Bridges: The Sustainable Bridges project aims at developing a reliable and cost-effective solution for detecting structural defects in bridges. The technology also allows for the estimation of the remaining lifetime of the bridges.

SDIC: It stands for Smart Dust Inventory Control, which allows the exact tracking of your products -from the packaging to the truck carrying the palettes.

Car Control: The Cortex's car control project automatically finds you the optimum route to the destination, where various sensor readings are included in the calculation; for instance, traffic, weather, constructions, jam reporting, etc.

GlacsWeb: This project monitors the drifting behavior of glaciers by aggregating pressure, movement and temperature data collected by sensors on top and within the glacier.

SMWF: It stands for Smart Mesh Weather Forecasting and is deployed in the Yosemite National Park. Various factors are measured, analyzed and predicted over the years, relating, for instance, to snowmelt.

MyHeart: Sensors measure variables impacting upon cardio-vascular diseases, thereby facilitating an early diagnosis. The system extends beyond the measurements of body internal variables, such as blood sugar level, pressure, etc; lifestyle factors ranging from choice of clothes, environment, etc., are also considered to draw preventive conclusions.

WiSeNts: It is an EU IST FP6 project with the aims to develop a new research domain, integrating the broad context of embedded systems with ubiquitous computing and wireless sensor networks in support of Cooperating Objects. They are integrating existing research in the field and related fields and developing a technology roadmap to drive the vision forward.

These projects are usually furnished by a large and knowledgeable partnership. The interests and backgrounds of participating partners, however, are usually fairly correlated and homogenous.

C. ARESA's Project Aims

The ARESA project had been assembled keeping in mind that designing a highly-efficient WSN is a crosscommunity exercise. Therefore, in contrast to abovementioned projects, the ARESA project is constituted of partners with a fairly complementary expertise with the aim to:

• explore new event-driven and asynchronous software and hardware architectures, tailored to extremely low power consumptions; • propose new communication and organization protocols, which are optimized in terms of energy consumption and robustness; • study new network structures which facilitate autoconfiguration and auto-organization;

• find new application protocols that are designed for data fusion and aggregation; • provide tools of modeling and validation, which also take into account the physical environment and the interaction thereof with the wireless sensor nodes; and • validate the developed concepts, protocols and mechanisms by means of a testbed. It is the aim of the ARESA consortium to propose an integrated industrial and applicative solution in the emerging area of wireless sensor networks. We wish to knowledgably influence the design of potentially tobe-standardized communication protocols for energyconstrained wireless sensor networks, which, we believe, form the basis of ambient data processing and communication systems.

D. Project Partners

To achieve these goals, a healthy mix of industrial and academic partners has been assembled with a strong track record in sensor and embedded communication systems. The industrial partners are comprised of France Telecom R&D, which has a leading expertise in ambient systems and wireless sensor network design, and CORO-NIS Systems, which are a European leader in sensor networks with invaluable real-world roll-out experiences. The academic partners, on the other hand, are constituted of the research centers LSR/INPG, CITI/INSA Lyon, VERIMAG/UJF and TIMA/INPG.

The France Telecom Group constitutes one of the biggest integrated operators worldwide. It offers a variety of services to its clientele, including mobile and fixed telephony, wired and wireless internet, as well as integrated home and business solutions. Through wireless sensor networks (WSNs), it hopes to offer more complete services by creating and facilitating ambient environments, which interface with incumbent and emerging services. For this reason, France Telecom has strong R&D activities in the area of WSNs -corroborated by the leadership of the ARESA project. The expertise of France Telecom is on the inputs of commercially viable sensor scenarios, design of low-power physical (PHY) and medium access control (MAC) layers, routing protocols, as well as cross-layer optimized communication mechanisms.

Coronis provides solutions for ultra-low-power (ULP) and long-range wireless applications. The company is the creator of Wavenis , a technology core for its radio frequency (RF) transceiver and wireless communication protocol. With Wavenis, Coronis offers a complete line of wireless platforms for original equipment manufacturers (OEMs), meeting technical, operational, and cost requirements of ultra-low-powered wireless mesh sensor networks. Major markets for Coronis technology are remote utility meter monitoring, home comfort, alarms for protecting people and property, home healthcare, industrial automation, centralized building management, access control, temperature monitoring, as well as longrange UHF RFID applications for the identification, tracking, and locating people and objects. Over 1 million devices made with Wavenis will be deployed by the end of 2006. The expertise of Coronis hence lays on development and commercialization, as well as realworld protocol design.

LSR/INPG has worked for years on issues related to MAC and routing for WSNs. In particular, their work relates to routing protocols in dense energy-constrained networks. The laboratory has also strong competences on wireless local area networks and ad hoc networks. Their main implication within ARESA is to design suitable protocols, which are tailored to extracted application scenarios and implementable into the developed hardware.

CITI/INSA, whilst also having a strong expertise in MAC, their implication in ARESA is on autoorganization and auto-configuration mechanisms for energy-constrained WSNs. Their research includes theoretical studies and modeling, where developed solutions will have to fit neatly with developed routing protocols.

VERIMAG/UJF is specialized in methods and tools for the development of safe and efficient embedded systems, ranging from embedded control to communication protocols. Their prime role in ARESA is to develop and analyze a formal model of a sensor network's behavior. This model should be accurate enough to estimate the energy consumption, and it should include the hardware, the application software, the protocols, and the physical environment. Such a modeling approach uses formal languages for reactive systems that have proven both efficiently executable and formally analyzable in a wide variety of embedded contexts.

TIMA/INPG is specialized in the design and the integration of autonomous wireless sensor systems. It has strong competencies in the methodology of designing hardware systems, in the field of the energy harvesting and in the field of systems very low fuel consumption. Their implications within ARESA are mainly related to optimizing the power consumption of asynchronous and event-driven hardware architectures, i.e. the algorithms developed and verified by the other partners will be realized and optimized for their hardware.

The complementary expertise of the involved partners gives ARESA the much needed synergy to facilitate application-oriented research, development and commercialization.

II. MOST CHALLENGING RESEARCH PROBLEMS

The research challenges listed in this section are the result of Coronis' experiences obtained during the rollout of the world-wide first commercial wireless sensor meter-reading service in the South of France. Almost one Million nodes are now operational, giving Coronis the knowledge on pertinent deployment problems. Resulting research challenges and associated problems are listed below.

A. Design Drivers

To accelerate roll-outs, the WSN ought to be decentralized, i.e. have no central point of control prior to deployment. This behavior is well studied within ad hoc networks and derived conclusions and insights can hence be utilized. Also, the information flow within WSNs is usually highly directed, i.e. from a large set of sensors towards a few data-collecting processing units. This behavior is contrasting a typical ad hoc network data flow, however, is well studied within traditional and multihop cellular communication systems. Finally, just as in any incumbent system, WSNs need to exhibit robustness (reliable delivery of data), integrity (correct delivery of data), and confidentiality (secure delivery of data); among several other requirements.

In contrast to known and well understood systems, however, a WSN bears also some fundamental design differences. These are summarized below.

Number of Nodes:

The number of nodes involved is very large, where current rollout examples include a few thousands; however, roll-out expectations are in the range of a few hundred thousand nodes communicating simultaneously. This is also atypical any wireless system today.

Energy: WSNs are nowadays battery powered and, because changing batteries in a few thousand nodes on a regular basis is clearly impractical, they are required to have a long lifetime and are hence considered to be highly constrained in energy. This is in contrast to any deployed wireless system, where nodes are usually either powered by the mains or rechargeable on a regular basis.

Applications: The gamut of applications is vast, hence requiring very different solutions to be developed for different applications. This problem is further enhanced due to the stringent energy constraints, requiring subtle solutions to be developed for different requirements.

This means that, unlike incumbent systems, our WSNs need to be:

• highly scalable (protocols ought to work at arbitrary number of nodes);

• highly energy efficient (at all layers and functionalities); and • highly application tailored (efficient for a given task). To achieve this, we shall subsequently look at some fundamental trade-offs and research challenges.

B. Impact of Scaling Laws

The first question posed by the consortium was on the issue of scalability of WSNs, i.e. when a network has to be considered large. To exemplify this problem, we have presupposed systems with and without internal conflicts. For instance, two systems without conflicts are our circle of true friends (comprised by a small number of elements) and the soldiers of an ant colony (comprised by a small number of elements). On the other hand, two systems with conflicts, frictions and competition are, for example, a few children left on their own (comprised by a small number of elements) and state without government (comprised by a small number of elements). As such, 'large' is hence not about size. It is rather about managing existing and emerging conflicts, and hence the amount of overheads needed to facilitate (fair) communication.

A management of conflicts is facilitated by a centralized entity, as for instance a base station in a cellular system. WSNs, however, do not have centralized entities prior to their deployment. This leads to scalability problems, as quantified below.

Kumar & Gupta's Throughput Scaling Law. This milestone contribution [START_REF] Gupta | Towards an Information Theory of Large Networks: An Achievable Rate Region[END_REF] quantifies the theoretically achievable network capacity assuming that everybody talks with everybody. They have determined that, assuming random deployment of N nodes in a large network, the throughput scales proportionally to 1/ √ N log N . The throughput hence decreases rapidly with an increasing number of nodes. In other words, no matter what we try, we cannot design a scalable protocol for large networks and hence topologies different from pure ad hoc have to be invoked.

Odlyzko & Tilly's Value Scaling Law. Mainly economically driven, various efforts in the past have been dedicated to establishing the value of a network in dependency of the number of its elements N . Sarnoff's Law quantifies the value of a broadcast network to be proportional to N Since a WSN is not truly of broadcast nature, nor do sensors form all possible communities, nor does every sensor communicate with every other sensor, another value scaling law is required to quantify the network's behavior. To this end, Odlyzko and Tilly have proposed a value scaling which is proportional to N log N [START_REF] Briscoe | Metcalfe's Law is Wrong[END_REF]. Their argumentation bases on Zipf's Law [START_REF] Zipf | Some determinants of the circulation of information[END_REF], which states that if one orders a large collection of entities by size or popularity, the entity ranked k-th, will be about 1/k of the first one. The added value of a node to talk to the remaining nodes is hence 1 + 1/2 + 1/3 + • • • + 1/(N -1) ∝ log N ; the value of the total network with N members hence scales with N log N . Among others, this law has been found to describe accurately the merging and partitioning of companies of unequal size [START_REF] Briscoe | Metcalfe's Law is Wrong[END_REF].

Using Odlyzko and Tilly's value scaling law, we wish to determine a clustering of a large WSN such that its value is maximized. With N nodes in the total network and C clusters, we have M = N/C nodes per cluster. Assuming that the value of the nodes within a cluster as well as the cluster heads obeys Zipf's Law, the value per cluster is M log M and the value of the clustered network is C log C • M log M . We shall normalize this value by the average number of links needed to maintain all nodes and clusters at any time, which is C log C • +M log M . The relative network value for different cluster sizes is depicted in Figure 1 with N = 10, 000 nodes. Clearly, clustering increases the normalized network value. For instance, if we assume a WSN with 10, 000 nodes and a directed information transmission among nodes to the cluster heads and among cluster heads towards the information sink, an optimal cluster number is 100 with about 100 nodes per cluster. For N = 1000, the optimal cluster number would about 12, etc.

Practical Hierarchical Scaling Law. This simple law quantifies the network throughput for practical systems with a given topology and gives insights into the design of data pipes between nodes. Among many possible topologies, we assume that each node communicates only with its respective cluster head and all cluster heads communicate among each other.

We hence assume a 2-tier hierarchy with N total nodes, C clusters and M = N/C nodes per cluster. As before, this 2-tier hierarchy requires two communication phases. In the first phase, all nodes communicate with their respective cluster-heads, and in the second phase, all cluster-heads communicate among each other. For subsequent analysis, we first assume all data pipes to have equal rates and then extend this to unequal pipes.

For equal data pipes, in the first phase, there shall be M time slots to transmit c • N bits, where c is a constant assumed to be one. In the second phase, there are hence M • C • (C -1) time slots to transmit these N bits to every cluster head. The efficiency is hence N/(M • C • (C -1)); remember that no new information is injected in the second phase.

For unequal data pipes, let us assume the clusterheads' pipes to be r times stronger than the data pipes between nodes towards the cluster heads. Therefore, in the first phase, there are again M time slots to transmit N bits; and in the second phase there are now M • C • (C -1)/r time slots to transmit these N bits. The efficiency is hence N/(M • C • (C -1)/r).

The relative network throughput for different cluster sizes is depicted in Figure 2 with N = 10, 000 nodes. Clearly, clustering increases the normalized network throughput only if the data pipes among the cluster heads are stronger or data aggregation [START_REF] Chen | A Meta-Data-Based Data Aggregation Scheme in Clustering Wireless Sensor Networks[END_REF] is performed to decimate the information shared among cluster heads. For instance, if we assume a WSN with 10, 000 nodes, then an optimal cluster number is 12 assuming the cluster heads' data pipes to be 1000 times stronger. If not all cluster heads communicated, as in the previous example, then the optimal cluster number would be larger.

Above quantification of throughput and value of a large WSN hence stipulate the use of clustered approaches. This is corroborated by real-world roll-outs, all of which use hierarchical and/or clustered network topologies with stronger data pipes between cluster heads. For example, the currently functioning meter reading application of Coronis uses a hierarchical approach [START_REF] Dugas | Configuring and managing a large-scale monitoring network: solving real world challenges for ultra-lowpowered and long-range wireless mesh networks[END_REF] (see also Figure 3) and so does Intel's WSN. The challenge is hence to design self-organizing protocols for the intra and inter cluster communication topologies. Parts of this shall briefly be discussed subsequently. Normalised Network Throughput sensors @ 1 kbps & cluster-heads @ 1 kbps sensors @ 1 kbps & cluster-heads @ 100 kbps (Zigbee) sensors @ 1 kbps & cluster-heads @ 1 Mbps (Bluetooth) Fig. 2.

Normalized network value with different clustering approaches.

C. MAC & Routing Protocols

Once deployed, wireless sensor networks are supposed to convey information from sensor nodes to sinks. This task should be carried out for the longest time possible given a fixed amount of initial energy at the sensors. Therefore, forwarding information in the network must be done according to the remaining energy of the nodes along the chosen paths. Dynamic routing, in which paths depend on the state of the network links, or adaptive routing that operates on finer-grain metrics, need to be adapted so that they not only reach a satisfactory solution, but also they do so at a controlled cost in terms of transmissions.

It is not also not apparent why routing needs to be dynamic in wireless sensor networks. In fact, routing can be computed at the network deployment, but then, any failure of nodes close to the sink, which are inherently more loaded than others, may cause the entire sets of sensors to become isolated. Early failure detection and reverting to backup paths is an immediate improvement and we believe that an autonomous routing protocol can provide a satisfactory solution.

Routing in wireless sensor networks may also benefit from new addressing schemes based on virtual coordinates. If a set of sensors is organized in such a way that their addresses correspond to coordinates in a virtual space, routing is straightforward and does not require any overhead nor signaling traffic for constructing routing tables. Information forwarding just consists of sending packets towards the coordinates of a given destination. Obviously such an approach needs to take into account disappearing nodes, but signaling traffic is in this case only needed to repair or reconstruct the addressing space after some node failures.

Moreover, such an approach can support organizing nodes in a logical way so that a high level querying of nodes will be made possible. For instance, we can explore a content-based approach for forwarding information within a wireless sensor network in which routing is done based on some attributes and not on addresses.

Furthermore, the routing problem is closely related to the way the MAC layer operates. For instance, the relative cost of broadcast transmissions compared to unicast transmissions influences the main design choices for the routing protocol. Conversely, lowering transmission costs is more or less crucial depending on the routing involved.

D. Self-Organization & Self-Healing

Because wireless sensors networks suffer from a random deployment with the least possible human intervention, self-'mechanisms' ought to be provided, such as self-configuration, self-organization, self-healing, etc. Self-configuration allows autonomous nodes configuration in order to set up communication. We commence by assuming that each node owns a unique identity; thus we are not focused on self-configuration [START_REF] Zhou | Prophet address allocation for large scale MANETs[END_REF].

The self-organization paradigm deals with an emergent behavior coming from local nodes interactions [START_REF] Prehofer | Self-organization in communication networks: principales and design paradigm[END_REF] and leading to a logical view of the network. The goal of such logical view is to structure the whole network in order to provide more efficient communication protocols (routing or data dissemination). In our point of view, an emergent behavior can be a connected structure as a virtual backbone or non connected one as clusters [START_REF] Lu | Low-Energy Self-Organization Scheme for Wireless Ad Hoc Sensor Network[END_REF]. The literature has provided intensive works in the area of distributed algorithms in order to structure the network: Connected Dominating Set, Maximum Independent Set, Local Minimum Spanning Tree, Relative Neighborhood Graph, Max-Min clusters, etc. Nevertheless, part of these works are focused only on virtual topology construction without consideration for a maintenance procedure in order to maintain this structure connected. Because of power-energy saving mechanisms or network dynamics the topology evolve. Thus, an event-driven maintenance algorithm should be provided, which out to be the goal of self-healing.

Another goal of self-organization is to build a framework to enable more efficient routing protocols, data aggregation mechanism, etc. Taking into account the applications coming from Coronis and France Telecom R&D, we deem important to specify a self-* architecture in order to build autonomously a logical structure with dedicated applications, such as data dissemination and data aggregation with self-healing mechanisms in order to be able to deal with network dynamics.

III. PROTOTYPING AND DEVELOPMENT

A strength of the ARESA consortium is that the researched protocols are directly fed into development efforts, as described below.

A. Virtual Prototyping and Formal Analysis

Design Issues for Sensor Networks. A sensor network may be considered as a whole, as a new kind of computer system dedicated to one particular application. It is an embedded system, reacting to the stimuli of some physical environment. It is also subject to the usual constraints of embedded system design: resources are scarce, and it is very difficult, if not impossible, to modify a sensor network's behavior once it has been deployed. Moreover, the sensors are usually powered by a battery that cannot be recharged. They should therefore have the lowest consumption possible to maximize the network lifetime.

Several issues may be considered for the design of sensor networks, among which:

• How to program the network? For example, knowing that the ultimate goal of the network is to detect a fire, how to design the application code for each node? • How to take security constraints into account? Indeed, a sensor network is probably very easy to attack. • How to perform energy-aware design? In ARESA, we are mainly interested in the last point. The problem is difficult because all the elements of a sensor network have an influence on energy consumption: the hardware of a node, the sensors, the mediumaccess-control and routing protocols, the application itself, the initial self-organization phase, and even the physical environment that stimulates the sensors (see, for instance [START_REF] Samper | On the Importance of Modeling the Environment when Analyzing Sensor Networks[END_REF], where we showed that a precise modeling of the physical environment is compulsory for a realistic estimation of the energy consumption).

The design of an energy-"optimal" solution is probably out of reach because of all the interacting criteria, and the complexity of some of the elements that have to be taken into account. If there is no way to compute an optimal solution, then the only available method seems to build complete solutions and then to evaluate them. However, since a sensor network includes dedicated hardware, it may be long and costly to build a complete solution before evaluating it.

Virtual Prototyping Approaches. For all these reasons, the usual approach is to build a virtual prototype of a sensor network, and then to perform simulations or mathematical analyzes in order to evaluate the energy consumption. This is the approach taken by people who design new protocols, and show their benefits using a network simulator. In all these approaches, a lot of abstractions are necessary, in order to build manageable models of very large systems (thousands of nodes). For instance, the energy consumption may be evaluated by counting packets, and associating a worst-case estimated energy with the transmission of one individual packet.

However, there seems to be a wide agreement on the fact that traditional network simulators like NS are not sufficient for ad-hoc sensor networks. In particular, they cannot be used to describe the hardware in an accurate way, which seems compulsory for power analyzes. A lot of approaches have been proposed for simulating ad-hoc sensor networks in a both accurate and efficient way. None of these approaches is formalized. Libraries have been developed for some reusable elements of the models, like the protocols, but it is still hard to obtain an accurate and efficient simulator while preserving the faithfulness of the model. Moreover, those simulators do not help in modeling the environment of the network, i.e., the physical phenomena that have some influence on the sensors. Finally, as soon as the power analysis needs an accurate simulation of the hardware of a node, the problem becomes the same as simulating efficiently a large piece of hardware. Simulating 1000 nodes at the Register-Transfer-Level (RTL) is probably hopeless. People in the hardware design domain have tackled this problem by defining new levels of abstraction (like the so-called "transaction-level modeling" [START_REF] Ghenassia | Transaction Level Modeling With SystemC: TLM Concepts And Applications for Embedded Systems[END_REF]) that are both accurate enough for a first approximate timing or power analysis, and fast to simulate. Developing such an approach for ad-hoc sensor networks requires a clear understanding of the abstractions than can be made on their behavior, while retaining their main power characteristics.

Formal Virtual Prototypes. In ARESA, we address the following challenge: we aim at developing accurate prototypes of sensor networks, that can also be formally analyzed. What does it mean? Consider a property related to the lifetime of a sensor network, for instance: "is this possible to spend more than the energy E in a time period less than T ?". This property can be tested on a lot of simulations, and we may find a scenario that indeed spends more than E in less than T . If we do not, however, it does not mean that there exists no such scenario. It may well be the case that we did not try a particularly complex scenario that exhausts all the nodes within a very short time. The aim of formal validation is to give exhaustive answers to such questions, for instance by exploring a mathematical structure that represents all the possible behaviors of the system. However, automatic formal validation usually faces two major problems: first, the mathematical structure to be explored may be huge; this is called state-explosion, and this is mainly due to the fact that representing the behavior of a parallel system amounts to considering all the states in the Cartesian product of the states of the parallel activities. Even if a node has only two states of interest, a 1000-nodes sensor network has 2 1000 potential states. Second, even if the size of the model is manageable, it may be the case that the properties we want to assess are undecidable. In other words, there exists no algorithm able to answer the question. It happens as soon as one includes general computations on numbers in the model, and it is likely to be the case when modeling the energy consumption of a sensor network.

To address these two problems, one needs either human interaction, or aggressive abstractions. In ARESA, we favor automatic techniques and tools, hence we are interested in the development of abstract models. What does it mean for a model to be abstract? First, since the exact behavior of a sensor network depends on the physical environment, no model or prototype can be exact with respect to energy consumption; all models are bound to make some hypothesis on the physical environment. For evaluating energy consumption, we should develop worst-case models. For instance, the radio can be modeled by considering several states depending on the emission power. A worst-case consumption is then associated with each state.

For the virtual prototype to be trusted, all the worstcase abstractions have to be well-understood and wellidentified in a global model. This is why we think that a formal model should be organized into components, in such a way that replacing a component C in a global model M by a more abstract version C yields a new global model M which is indeed more abstract than M . This abstraction preservation property is essential when playing with various abstractions of the individual components. We should also be able to prove the property: C is more abstract than C, i.e., that the power estimations computed with C are always greater than that of C.

Prototyping Challenges. To summarize, the challenges we exposed related to the development of virtual prototypes of sensor networks are as follows:

• The prototype should include all the elements that determine energy consumption: node hardware, protocols, application code, and physical environment. • The prototype should be executable, so as to allow efficient simulations • The prototype should be described in a formallydefined language, so as to be analyzable by automatic tools

• The prototype should be made of well-defined components, so as to allow modular abstractions as defined above. We will use techniques and tools developed at VER-IMAG, and which have proven useful in other domains of embedded system design, ranging from embedded control to communication protocols. In particular, the modeling and verification of various embedded reactive systems has led to the design of a formal modeling language for physical environments, which is well adapted to sensor networks. We plan to develop several forms of formal models, depending on the properties to be proved. For instance, the details of a node are best described in some synchronous formalism, while the way several nodes behave together is best described with some asynchronous formalism. We will also use the new application domain of sensor networks to make our tools evolve.

B. Hardware Implementation

As mentioned in the previous section, it is important to tackle the power consumption issue at hardware level in order to increase the WSNs' lifetime. In traditional synchronous design, computations are sequenced by a clock that is connected to every digital hardware blocks which perform calculations at each rising edge, even if no valid data are present. Although some techniques, such as clock gating, exist to reduce activity of such circuits, lots of unneeded signal transitions remain.

In asynchronous QDI technology [START_REF] Martin | The Limitations to Delay-Insensitivity in Asynchronous Circuits[END_REF], the parts of the circuit that are performing an operation have an activity. The rest of the circuit consumes very little energy (only static leakage) and are immediately wokenup when an event occurs on its inputs. Synchronization between blocks is performed by requests and acknowledgements.TIMA's asynchronous circuits are delay insensitive (functional correctness independent on gates and wires delays). Therefore, they are very robust to voltage changes and are particularly suited for dynamic voltage scaling (DVS) and can support very low voltages since signal rise and fall times does not affect the correctness.

The DVS technique is an efficient way for reducing power consumption. As the energy consumed is proportional to the square of the voltage, reducing a little the supply voltage decreases a little the speed of the hardware (including the microcontroller) and allows important energy savings at the same time [START_REF] Pering | Dynamic Voltage Scaling and the Design of a Low-Power Microprocessor System[END_REF], [START_REF] Rios-Arambula | On the use of Feedback Systems to Dynamically Control the Supply Voltage of Low-Power Circuits[END_REF]. As a result, we have experienced that for an asynchronous AES crypto-processor (130nm), reducing the computation speed by a factor 2 reduces energy consumption by a factor 9 [START_REF] Bouesse | A Clockless low-voltage AES crypto-processor[END_REF].

To exploit the benefits of asynchronous circuits and to have very low software overheads, TinyOS 2.0 and an appropriate scheduler are considered to apply a power consumption policy that should fit the application requirements at best.

The goal of our work is to efficiently combine TinyOS and asynchronous circuits by, first, supporting hardware tasks within the programming model of TinyOS at quasi no cost in terms of energy and timing overheads; and, second, defining an appropriate scheduler and regulation mechanism that will apply a power consumption/speed policy that should fit the application requirements at best.

C. Product Development

Coronis' real-world meter-reading application [START_REF] Dugas | Configuring and managing a large-scale monitoring network: solving real world challenges for ultra-lowpowered and long-range wireless mesh networks[END_REF] requires solutions with least possible power consumption. To this end, Coronis' devices are equipment by a single lithium battery running on 2.7Ah. Because of leakage, the passivation effect, and internal resistance, the actual usage, however, is limited to only 1.6Ah over 10 years. This yields an average consumption of 18µA over 10 years.

Further elaborating on the hardware design approach mentioned in the previous section, higher energy savings can be achieved if the synchronization procedures in the meter-reading products are optimized. To this end, it is important to remember that when a sensor network has to operate for several years, device synchronization cannot be maintained continuously. If synchronization is lost in a synchronized network, devices try to re-establish synchronization by entering into high power consumption mode, thereby shortening the network's lifespan. One solution is for devices to synchronize only when necessary, with communication between devices backed down immediately afterwards. Synchronization can then be performed at predetermined times depending on the application; for instance, one second latency for meterreading. Application devices hence toggle sequentially between receive mode and standby mode for most of their lives.

It is also possible to extend range by reducing the bandwidth, while keeping the average operating consumption at an ultra low 15µA [START_REF] Dugas | Configuring and managing a large-scale monitoring network: solving real world challenges for ultra-lowpowered and long-range wireless mesh networks[END_REF]. New algorithms, tailored to these needs and requirements, are hence key to automatically configuring and managing nodes. As corroborated by currently running applications and briefly elaborated upon in the next section, routing and remote access functionality can be maintained without sacrificing neither technical goals nor key marketing requirements.

IV. BEING THE FIRST IN THE MARKET

Being the first in the market meant to identify potential markets first, where we have distinguished between:

• governmentally support market; • mass market; • specialized market; and • niche market. Governmental support is usually granted for applications which facilitate the health and security of its citizens. Money and deployment costs are usually not of prime interest. The technological solutions are often required to be robust, reliable and secure, but not necessarily optimum and efficient. For instance, a WSN deployed in border surveillance must be secure and reliable, but may have access to a wide frequency spectrum.

Mass markets are constituted by a non-negligible part of the population. That generally includes the human population, but -in the case of WSNsmay very well extend to the entire flora and fauna. Technological solutions need to be cheap and to a certain degree secure and reliable. The deployment and usage cycle of mass market products rarely exceeds half a dozen of years.

Specialized markets are application and need tailored markets, where a special solution might be applicable to a specific problem. Specialized markets infer high revenues, but not necessarily size. For example, monitoring the growth of grapes is a specialized market, however, with potentially large revenues. The technological requirements are usually very diverse and vary from application to application.

Niche markets address needs which are very specialized and small in size, but -if WSNs are deployed in these markets -life can be made easier and some savings achieved. For instance, the measurement of wind speed for households which use private wind turbines to generated energy. The technological requirements here are also very diverse and vary from application to application.

In Table I, a classification of these markets against typical application domains is given. It is clearly visible that the majority of applications will be in the specialized market, i.e. requiring application-tailored solutions. With this in mind, the market requirements of one of Coronis' commercial products are discussed from a technical point of view.

By operating in the 2.4GHz, 868-870MHz, and 902-928MHz bands (like the recent ZigBee proposal), wireless technologies continue to allow critical equipment to use the less crowded and more regulated 868MHz and 915MHz bands when participating in regular Bluetooth piconets. Coronis' approach is to connect With these high level technical requirements, Coronis were the first to provide a commercially viable meterreading solution. The toughest bit about evaluating the feasibility of this new technology has been to see how it actually performed in the field. To this end, 25,000 water meter-reading nodes were installed in the South of France, which constitutes a fully operational automated network capable of gathering water consumption data via a wireless mesh network and transmitting it back to the home office. It is the first large-scale remote-controllable, fixed wireless network of its type in Europe.

With this technology, not only cost-effective solutions for the utility company have been provided, but remote monitoring actually offers many advantages for water customers as well, such as accurate billing and increased meter security. The current customers can even check their usage via Internet and even receive alerts via mobile or pager regarding abnormal consumption patterns such as those caused by leaks. The topology of this commercial application is depicted in Figure 3 Novel applications are currently being envisaged, ranging from environmental monitoring to health care applications. For these applications, research from ARESA will feed directly into the optimization of the data flows, handling of node failure and optimizing the network's lifetime. The proposed solutions are then verified by means of formal methods, so as to ensure 100% reliability under all conceivable operating conditions. These algorithms are then optimized in terms of hardware and energy consumption, which finally facilitates the development of commercially viable products and shortens the crucial time-to-market.

V. CONCLUDING REMARKS

The aim of this paper was to expose some crucial issues related to the research, development and commercialization of wireless sensor networks. The coherent insights given here stem from the diverse background of the partners involved in the ARESA project. We have looked at the impact of clustering on the maintenance of a large-scale wireless sensor network, where node numbers are often in the thousands. Using some known scaling laws, we have determined that scalable protocols for flat sensor network topologies cannot exist and that clustering or hierarchical approaches ought to be used instead. We have also identified some approaches which may be useful in determining an optimum cluster size. Thereafter, we have identified crucial research problems at medium access control and routing levels, as well as related to auto-organization and self-healing mechanisms. The concept of virtual prototyping has then been introduced -a tool which proves to be very useful in evaluating the reliability of designed protocols. Furthermore, various hardware issues have been identified which require highest attention if energy consumption was to be minimized. Finally, some market and commercial insights have been given and it has been identified that research and development need to be closely coupled so as to shorten time-to-market for novel wireless sensor networks technologies. We hope that these insights are of significance for emerging and future real-world installations, such as data collection, remote control solutions, wireless telemetry, automatic monitoring, metering solutions and smart environments such as homes, hospitals, and buildings of all kinds.

  [START_REF] Bilby | The General: David Sarnoff and the Rise of the Communications Industry[END_REF]. Reed's Law claims that with N members you can form communities in 2 N possible ways; the value hence scales with 2 N [3]. Matcalfe's Law, unjustifiably blamed for many dot-com crashes, claims that N members can have N (N -1) communication connections; the value of the network hence scales with N 2 .

  Fig. 1.Normalized network value with different clustering approaches.

  C [logarithmic]

Fig. 3 .

 3 Fig.3. An optimal topology is automatically configured for each local group in a large-scale ultra-low power network.

TABLE I TAXONOMY

 I OF APPLICATION DOMAIN VERSUS MARKET TYPE.

	govern. mass special. niche
	home, office
	control, automation
	logistics, transport
	environm. monitoring
	healthcare
	security
	tourism, leisure
	education, training
	entertainment
	ultra-low power devices to the outside world without
	defining another standard in the crowded 2.4GHz band.
	The real challenge is to devise a solution that keeps
	networks manageable, while offering both longer range
	and far lower power consumption. Gaussian frequency
	shift keying modulation, frequency-hopping, hop rate,
	and timing in the commercial solution are identical to
	those for Bluetooth, but the data rate has been reduced.
	Furthermore, as discussed above, synchronization has
	been adapted to match the requirements of stringent
	energy constraints.