Ken Tore Tallakstad

Henning Arendt Knudsen

Thomas Ramstad

Grunde Løvoll

Jørgen Knut

Renaud Måløy

Eirik Grude Toussaint

Flekkøy

Knut Jørgen Måløy

Renaud Toussaint

Eirik Grude Flekkøy

Steady-state Two-phase Flow in Porous Media: Statistics and Transport Properties

Keywords: numbers: 47.56.+r, 47.55.dd, 47.55.Ca, 89.75.Fb

come   L'archive ouverte pluridisciplinaire

Different types of immiscible multi-phase fluid flow in porous media play an important role in many natural and commercial processes [START_REF] Bear | Dynamics of Fluids in Porous Media[END_REF][START_REF] Dullien | Porous Media Fluid Transport and Pore Structure[END_REF][START_REF] Sahimi | Flow and Transport in Porous Media and Fractured Rock[END_REF]. The complex fluid patterns observed in such processes have been extensively studied and modeled over the last decades, see [START_REF] Bear | Dynamics of Fluids in Porous Media[END_REF][START_REF] Dullien | Porous Media Fluid Transport and Pore Structure[END_REF][START_REF] Sahimi | Flow and Transport in Porous Media and Fractured Rock[END_REF][START_REF] Lenormand | [END_REF][5] and references therein.

The vast majority of work up to now has focused on invasion processes; either pure drainage or pure imbibition. These inherently transient processes give different displacement patterns and are classified into capillary fingering [6], viscous fingering [5,[7][8][9][10][11][12], and stable front displacement [13,14]. These are non-stationary processes, and to understand them in a broader context there is a need to understand the stationary case which has received far less attention: steady-state flow, which is in equilibrium in the sense that average flow properties and distribution functions are invariant in time. This stationary system is in statistical equilibrium although it is a dissipative process; an external energy input balances the internal energy loss to maintain the equilibrium.

With some notable exceptions there is to our knowledge very little pore scale experimental data available for such problems [15]. The Payatakes group did pore scale steady-state experiments using network models etched in glass [16], and later theoretical modeling predicting the non-linearity of such flows [17]. In addition some numerical work has focused on a steady-state regime; a pore scale Lattice-Boltzmann study by Rothman et al. [18] and network simulations at larger scales by .

In this Letter we experimentally demonstrate that an equilibrium flow situation results after simultaneous injection of two fluids into a porous medium. This allows for the combination of mean-field approximations of local quantities and energy dissipation considerations. As a result we analytically obtain the highly non-trivial steady-state pressure-flowrate relationship. Furthermore, in contrast to transient flows, steady-state mass conservation gives a normalization condition based on flowrate rather than saturation. From this we derive a scaling law of the cluster size distributions of nonwetting fluid.

FIG. 1: Depicted inside the model frame (L×W =85×42 cm 2 ) is the initial transient stage of an experiment. There are 15 independent inlet holes with tubes and syringes attached. Every second syringe injects the wetting phase, the others inject the nonwetting one. An outlet channel with four exit holes allows the fluid mix to leave the system. Three SensorTechnics 26PC0100G6G flow-through pressure sensors (indicated by the solid rectangles) are attached alongside the model at the positions y = 0, L/2, L. Additionally one pressure sensor is attached to the first "air"-tube from the bottom.

Our system is shown in Fig. 1. The horizontal porous model consists of a mono-layer of glass beads of diameter a=1 mm, which are randomly spread between two transparent contact papers [7,13]. The model dimensions are L × W =85 × 42 cm 2 , with thickness a and volume V = aLW . The porosity and absolute permeability are measured to be φ = 0.63 and κ 0 = (1.95 ± 0.1)× 10 -5 cm 2 respectively. The wetting fluid used is a 85% -15% by weight glycerol-water solution dyed with 0.1% Negrosine, and has a viscosity µ w = 0.11 Pa • s at room temperature. Air is used as the nonwetting fluid with viscosity µ nw =1.9 × 10 -5 Pa • s, giving a viscosity ratio M =µ nw /µ w ∼10 -4 . The surface tension is measured to be γ = 6.4 × 10 -2 N • m -1 .

The tuning parameter in the experiments is the total flow-rate, i.e. the sum of the flow-rate of the wetting and nonwetting fluid, and can during steady-state be written as

Q tot = Q w + Q nw = (8 + 7)Q 0 ,
where Q 0 is the flowrate from every single syringe.

Gray scale images of the flow structure are taken at regular intervals with a Pixelink Industrial Vision PL-A781 digital camera. An image contains 3000×2208 pixels, corresponding to a spatial resolution of ∼0.19 mm per pixel (27 pixels in a pore of size 1 mm 2 ). All analysis is done on the basis of black and white thresholded images [14], and the measured pressure signals.

The porous model is initially saturated with the wetting phase. An experiment is started by injecting the fluid pair from every other inlet hole. The initial structure consists of bubbles or clusters of air distributed over various sizes, but always much smaller than the system size. The clusters are embedded in a background field of percolating wetting fluid. Usually the smallest air clusters are immobile and trapped, whereas larger clusters are mobile and propagate in the porous medium. However, trapped clusters can be mobilized when they coalesce with larger migrating clusters and migrating clusters can be fragmented and thereby trapped.

We divide an experiment into two regimes. A transient regime where the mix of nonwetting clusters and wetting fluid gradually fills up the model, as seen in Fig. 1. During this time the measured average pressure difference between y = 0 and y = L, ∆P L increases. This is due to the presence of more and more air clusters trapped in the system, effectively lowering the relative permeability for the viscous wetting fluid. At some characteristic time, shortly after both phases are produced at the outlet, ∆P L starts to fluctuate around a constant value. This marks the start of the steady-state (or statistically stationary) regime. The whole model now contains a homogeneous mix of the two phases, transported through the model without "long time" flow parameter changes.

Through six experiments we have studied how the measured steady-state pressure difference ∆P L varies with the Capillary number Ca, defined as

Ca = µ w Q w a 2 γκ 0 A , (1)
where A = W a is the cross-sectional area. This is shown in Fig. 2, for a span in the Ca-number of roughly two decades. It is evident that the pressure is consistent with a power law in the Ca-number, ∆P L ∝ Ca β , where the exponent is found to be β = 0.54 ± 0.08. This is a nontrivial result, and we will return to the discussion shortly. A general trend in the experiments, in passing from high to low Ca-numbers, is that the size or area of the largest air clusters increases. This means that the geometry of the clusters depends on the steady-state pressure gradient. To quantify this, we have found the normalized probability distributions of cluster extension in the xand y-direction, P (l x ) and P (l y ) respectively (see Fig. 3). We define the extension lengths l x and l y as the sides of the smallest rectangle (bounding box) that can contain a cluster. For clarity, l y lays parallel to the average flow direction, whereas l x lays transverse to the average flow direction.

Analysis shows that for a cluster of a given area s, the extension lengths have well defined means l x and l y increasing monotonically with s [23]. The corresponding standard deviations are small and proportional to these means (relatively 20%) for all l x and l y values [23]. Furthermore we find that l x = l y up to a characteristic length scale in the system l * . Above l * , l y > l x , as shown in the upper right inset of Fig. 3. The exact same behaviour is seen in the distributions P (l x ) and P (l y ). Fig. 3 shows, for the same Ca-numbers used previously, a collapse of the P (l x ) and P (l y ) distributions by the rescaling l * φ P (l i ) vs. l i /l * . The scaling exponent φ = 2.8 ± 0.3 is taken as the value that gives the best collapse. Apart from the expected crossover when the extension lengths reaches the pore length scale, a, the collapse is very good. The above results reveal important information of our system; particularly that there is only one length scale l * , dependent on the pressure difference ∆P L , that controls the steady-state displacement structure.

In the following we shall give a simple and minimal scaling theory for the purpose of predicting the exponents β and φ from Fig. 2 and3 respectively. Consider a nonwetting cluster in the porous medium surrounded by flowing viscous wetting fluid. The cluster perimeter is made up of several menisci standing in different pores. A single meniscus at a particular position has a surface pressure given by the pressure difference of the nonwetting and wetting fluid on each side. The nonwetting fluid pressure is assumed constant inside the cluster due to its low viscosity, whereas the wetting fluid pressure is position dependent, decreasing in the y-direction. The neighbouring pores inside and outside the meniscus can be either imbibed or drained respectively if the surface pressure exceeds one of the capillary pressure thresholds for imbibition or drainage. The imbibition and drainage threshold pressures depend on pore geometry, and are thus distributions due to the randomness of the porous medium [24]. If all menisci along the perimeter are in mechanical equilibrium, the cluster is immobile. This is typically the case for clusters with only a small l y extension. However, for large enough l y , the viscous pressure drop on the wetting side of the perimeter is sufficient for a migration step to take place. Migration is the process of drainage in one pore and imbibition in another along the perimeter. The result of several migration steps is that the cluster moves, and perhaps also changes shape. The onset of migration will depend on the difference between the mean threshold pressure for drainage and imbibition [24]; a pressure we denote Pt . Furthermore, Pt predicts a characteristic length of extension l * for cluster mobility

|∇P | l * = ∆P L l * L = Pt , (2)
where we make the mean field assumption that the pressure gradient |∇P | is constant. The scaling of the last equality in Eq. ( 2) is verified experimentally, as shows the lower left inset in Fig. 3. Note that l * also determines when clusters become unstable against breakup, since the mechanism of cluster mobilization is the same as that of cluster fragmentation. This is important because it links the single crossover length that collapses both P (l i ) distributions in Fig. 3 to the extension of mobilized clusters.

Particularly, it means that the characteristic l x -extension of mobilized clusters is l * . Equilibrium conditions require the total dissipation in the system to be balanced by the work rate done through the external pressure drop: Q tot ∆P L = D f . In obtaining D f , we assume that the main contribution to dissipation is in the volume of the wetting fluid, and that the dissipation in the nonwetting fluid is negligible. Visual observation indicates that most of the wetting fluid is restricted to flow through narrow channels at some typical spacing. The motion and configuration of the nonwetting clusters seem to show that the channel width is of the order of a pore size ∼ a, and that the permeability in between channels is made relatively low by the presence of lowly mobile nonwetting clusters. Motivated by these observations we define a dissipative wetting fluid volume.

V dis = LA dis = La 2 W l * = aV l * , (3)
where l * is taken as the spacing between channels, making W/l * the number of channels through the system. This simplification of channel flow of the wetting fluid is a strong assumption but supporting numerical simulations also show that the dissipative volume is constrained to a small fraction of the total volume, changing with the flowrate.

Since the overall interface area between the wetting and nonwetting phase is fluctuating around a constant value in steady-state, changes in the potential energy stored in the interfaces do not contribute to the average dissipation, and we are justified in writing

Q tot ∆P L = D f = - V dis dV u|∇P | = µ w κ 0 V dis dV u 2 , (4)
where we have applied Darcy's law locally, in the dissipative part of the wetting fluid. Taking the local Darcy velocity u = (Q w /A)(V /V dis ) as a constant, and using Eqs. ( 3), ( 1) and ( 2) respectively we obtain

∆P L = 8γV l * 15a 3 A Ca ⇒ |∇P | 2 = 8γ Pt 15a 3 Ca , (5)
i.e. |∇P | ∝ √ Ca, consistent with the exponent β in Fig. 2. An alternative interpretation of this result is that the wetting fluid experiences an effective permeability, assuming a Darcy law κ eff (Ca) = µ w Q w /A|∇P |, due to the flow of air:

κ 2 eff = 15γ 8a Pt κ 2 0 Ca . (6)
We turn now to the distributions of cluster extension lengths, and the found exponent φ. From the collapses in Fig. 3 it is seen that the distributions can be written

P i (l) = l -φ h i (l/l * ) , (7)
where the cutoff functions h i (x) are dominating. Note that this scaling form should only be expected to hold for l above the lower cutoff scale ∼ a. To obtain φ we use the fact that the nonwetting flowrate is an imposed quantity, and must in steady-state be equal to the accumulated flow of all mobile clusters on the average. The contribution of a single cluster of extent l to the total nonwetting flowrate is alU (l), where U (l) is the average centre of mass velocity. Further, the average number of clusters of extension (l, l + dl) that intersects any given cross section A, is given as dlP i (l)N l/L. Hence we can write

Q nw = aN L ∞ a dl l 2 U (l)P i (l) , (8)
where N is the total number of clusters. This number is measured and found to depend only weakly, at most logarithmic, on Q nw . For simplicity, N will be treated as a constant in the following.

To obtain U (l) we make the general assumption that it is linear in Q nw and has some functional dependence on l/l * U

(l) = Q nw A f (l/l * ) . (9)
As a first order approximation f (x) would be a step function, since clusters of size l < l * usually are immobile. By Eqs. ( 7) and ( 9) and the substitution x = l/l * , Eq. ( 8) can be written

1 = l * (3-φ) aN V ∞ a/l * dxx 2-φ f (x)h i (x) . (10)
Since the right hand side of Eq. ( 10) must be independent of l * , we obtain φ = 3 consistent with the experimental value in Fig. 3.

In conclusion, experiments have been done on twophase flow in a porous medium under steady-state conditions. In contrast to invasion processes and other inherently transient phenomena of two-phase flow, steadystate is in a statistical sense an equilibrium situation. Whereas the description of transient behaviour is a whole range of separate loosely attached cases, depending on flow parameters, the description and formalism for steady-state should be more integrated and universal. Our work explores a part of its parameter space and we find a robust power law behaviour: pressure increases as ∆P L ∝ Ca 0.5 , alternatively for permeability κ eff ∝ Ca 0.5 .

The power law is valid over roughly two decades, but there should be cutoffs for large and small flow rates.

For high flow rate, the cutoff cluster size will approach the pore size, and from that point the permeability must reach a plateau. The same is the case for flow rates low enough that the largest clusters are limited by the system size. These limits were not realized experimentally, but in numerical work these cutoffs have been seen [20,22].

The scaling behaviour of the system was explained by theoretical arguments, relying on a high viscosity ratio. Numerical work with a lower viscosity ratio (steady-state but somewhat different boundary conditions) indicates a lower exponent: ∆P L ∝ Ca β ′ , where β ′ < 0.5 [20]. We conjecture that the presented theory is a limiting case, thus suitable as a starting point for further theoretical developments, aiming at incorporating the more complex case where the two fluids have more similar viscosities.

The work was supported by NFR through PETRO-MAKS, a CNRS PICS, an ANR ECOUPREF project and a regional REALISE program. A special thanks to Alex Hansen for useful comments.

FIG. 2 :

 2 FIG.2: Mean value of the pressure difference ∆PL during steady-state as a function of Ca. The fluctuations in ∆PL are of the order of 1 kPa, i.e. very small compared to the mean values. A power law dependence is found, with exponent β = 0.54 ± 0.08.

FIG. 3 :

 3 FIG. 3: Air cluster extension length distributions P (ly) (filled markers) and P (lx) (empty markers), collapsed by the rescaling l * φ P (li) vs. li/l * , where i ∈ {x, y} and φ = 2.8. The dashed vertical line at li = l * indicates the start of the different cutoff behaviour in the two directions. Lower left inset shows (solid line) ∆PL ∝ 1/l * with the corresponding experimental values (squares). Upper right inset shows ∆l/ ly vs.lx ly /l * 2 , where ∆l = ly -lx .


Nadjime Pindra

Veronique Lazarus

Jean-Baptiste B Leblond

Jean Schmittbuhl

Renaud Toussaint

Deformation of the crack front during propagation in some disordered medium: theoretical and experimental studies

à 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.

In a perfectly ideal homogeneous material, an initial straight crack front remains straight during propagation. But in an heterogeneous disordered materials, it becomes rough. The aim of the present paper is to derive an analytical description of the evolution of this roughness and to compare it to experimental results. The assumption of quasi-static brittle crack propagation will be done. Among the experimental works, one may cite on the one hand, the pioneer work of Daguier et al. [START_REF] Daguier | Roughness of a crack front pinned by microstructural obstacles[END_REF] in which the crack front is obtained postmortem, the crack surface being marked by ink and on the other hand, the works of Delaplace, Maloy and Schmittbuhl [START_REF] Schmittbuhl | Direct observation of a self-affine crack propagation[END_REF][START_REF] Delaplace | High resolution description of a crack front in a heterogeneous plexiglas block[END_REF] in transparent plexiglas in which the crack front can be observed in situ during its evolution. They deal mainly with the universal self-affine character of the crack front. The roughness exponent ξ was measured between 0.5 and 0.6. Here, we have again used the experimental framework of [START_REF] Schmittbuhl | Direct observation of a self-affine crack propagation[END_REF][START_REF] Delaplace | High resolution description of a crack front in a heterogeneous plexiglas block[END_REF] to measure the time evolution of the fluctuations in addition to its roughness. All the theoretical studies of the statistical properties of the crack front performed in quasi-static, use Bueckner [START_REF] Bueckner | Weight functions and fundamental fields for the pennyshaped and the half-plane crack in three-space[END_REF]-Rice [START_REF] Rice | First-order variation in elastic fields due to variation in location of a planar crack front[END_REF] weight function theory, also called line elastic models, to evaluate the stress intensity factors along the perturbed crack front. Among them, one may distinguish two groups depending on the type of the advance law used. The first ones [START_REF] Schmittbuhl | Interfacial crack pinning: effect of nonlocal interactions[END_REF][START_REF] Schmittbuhl | Interfacial crack front wandering: influence of quenched noise correlations[END_REF][START_REF] Pindra | Slight in-plane perturbation of a system of two coplanar parallel tensile slit-cracks[END_REF] deal with crack advance governed by brittle fracture Irwin's criterion with a slightly heterogeneous toughness. This criterion is a threshold type one: the crack propagates only if the stress intensity factor becomes equal to the local toughness. In particular, by a first order analysis the roughness exponent was derived and found to be ξ = 0.37 or ξ = 0.5 depending on the papers. This apparent discrepancy will be considered further. The second group deals with crack advance governed by Paris' law (fatigue or subcritical fracture). It is a time dependent type criterion: the rate of crack advance is proportional to a power law of the stress intensity factor. Lazarus, Leblond and coauthors have performed the study of a tensile tunnel-crack [START_REF] Favier | Statistics of the deformation of the front of a tunnel-crack propagating in some inhomogeneous medium[END_REF] and of a tensile semi-infinite interfacial crack [START_REF] Pindra | Slight in-plane perturbation of a system of two coplanar parallel tensile slit-cracks[END_REF]. Contrary to the case of threshold advance law, their first order study in crack advance was not sufficient to obtain the crack front roughness. However Adda-Bedia and Katzav [START_REF] Katzav | Roughness of tensile crack fronts in heterogenous materials[END_REF] performed the second order study for a semi-infinite crack and obtained ξ = 0.5.

Here, the work of Pindra, Lazarus and Leblond [START_REF] Pindra | Slight in-plane perturbation of a system of two coplanar parallel tensile slit-cracks[END_REF] is applied to experiments made with the same framework as Delaplace, Maloy and Schmittbuhl [START_REF] Schmittbuhl | Direct observation of a self-affine crack propagation[END_REF][START_REF] Delaplace | High resolution description of a crack front in a heterogeneous plexiglas block[END_REF]. For Irwin's advance law, using Bueckner-Rice formulation for a semi-infinite crack subjected to line loading on its faces, the evolution of the variance and power spectrum, so as the roughness exponent are derived and compared to previous theoretical results of Schmittbuhl, Vilotte and coauthors [START_REF] Schmittbuhl | Interfacial crack pinning: effect of nonlocal interactions[END_REF][START_REF] Schmittbuhl | Interfacial crack front wandering: influence of quenched noise correlations[END_REF]. Then, comparison with experiments are performed, not solely on the roughness exponent as in previous papers but also on the evolution of the crack front amplitude. The experimental framework is the one still used by Delaplace, Maloy and Schmittbuhl [START_REF] Schmittbuhl | Direct observation of a self-affine crack propagation[END_REF][START_REF] Delaplace | High resolution description of a crack front in a heterogeneous plexiglas block[END_REF]. Two PMMA transparent plates are sandblasted and annealed together by increasing temperature to obtain one block with a weak plane where the crack will propagate. The larger plate is clamped. A normal displacement is applied to the smaller one that creates a stable propagating crack in the weak plane of the specimens (fig. 1). The crack front is observed with a microscope at several stages of the propagation.

Theoretical determination of the roughness 2.1 Evolution of the perturbation of the crack fronts

Since [START_REF] Pindra | Slight in-plane perturbation of a system of two coplanar parallel tensile slit-cracks[END_REF] is a general paper (an interfacial crack, several type of loadings, fatigue and fragile fracture are considered), and hence is quiet complicated, we give here the reasoning in the particular case of an homogeneous medium (ε = 0), the particular line loading of the experiments, brittle fracture, although most of the results can be directly (but not painlessly) derived from the paper. Let us consider a semi-infinite plane crack in an infinite homogeneous linear elastic medium. Line tractions ±P e y are applied on the crack faces at a distance a of the crack front (fig. 3). We suppose that for each equilibrium position of the crack front, Irwin's criterion is satisfied, that is:

K(M) = K c (M) (1)
at each point M of the crack front. Here K(M) is the stress intensity factor (SIF) and K c (M) the toughness, both at point M. The SIF along the crack front are then :

K = P 2 πa (2)
so that the loading corresponding to position a of the front verifies :

P = K c πa 2 (3)
Notice that the propagation is stable, in the sense that the loading has to be increased for the crack to advance (a increases). Now introduce some small fluctuations of the toughness :

K c (z, x) = K c (1 + κ(z, x)), |κ| ≪ 1 (4)
It produces small fluctuations δa(z, a) and δK(z, a) of the crack front position a(z, a) and of the SIF K(z, a) (see fig. 3) so that :

a(z, a) ≡ a + δa(z, a), |δa(z, a)| ≪ a K(z, a) ≡ K(a) + δK(z), |δK(z, a)| ≪ K(a) (5)
where a and K(a) denote their mean values. Expanding Irwin's criterion (1) to first order and identifying terms of order 0 and 1, one gets :

       K(a) = K c δK(z, a) K(a) = κ(z, a) (6)
Now we use Bueckner [START_REF] Bueckner | Weight functions and fundamental fields for the pennyshaped and the half-plane crack in three-space[END_REF]-Rice [START_REF] Rice | First-order variation in elastic fields due to variation in location of a planar crack front[END_REF] formalism, expanded to line tractions as in [START_REF] Pindra | Slight in-plane perturbation of a system of two coplanar parallel tensile slit-cracks[END_REF], to express the perturbed SIF to first order in δa :

δK(z, a) K(a) = - 1 2a δa(z, a) + 1 2π PV ∞ -∞ δa(z ′ , a) -δa(z, a) (z ′ -z) 2 dz ′ (7)
Inserting equation ( 7) into (6.2) and taking the Fourier transform of the equation, one gets:

δâ(k, a) = - 2aκ(k, a) 1 + |ka| (8)

Statistical study of the deformation of the front

We consider a large number of random possible realizations of the heterogeneous medium and the crack geometry. Statistical invariance of the functions κ(z, x) in the direction z and x, and δa(z, a) in the direction z being assumed, the two point auto-correlation functions of these functions depend only on the relative position of the points considered:

E [κ(z 1 , a)κ(z 2 , a)] = K(z 2 -z 1 ) (9) E [δa(z 1 , a)δa(z 2 , a)] = A(z 2 -z 1 , a) (10)
where E[X] denotes the mathematical expectation of any quantity X. The functions K(z) and A(z, a) can be identified with the average values of κ(z ′ , a)κ(z ′ + z, a) and δa(z ′ , a)δ(z ′ + z, a) over the crack front, provided an ergodic hypothesis is made. The function A(z, a) and its z -Fourier transform A(k, a) (the spectral density of the perturbation) provide statistical information about the geometry of the crack front.

Using the property that

A(k, a) = 1 2π ∞ -∞ E [δâ(k 1 , a)δâ(-k, a)] dk 1 (11)
and equation ( 8), one obtains the power spectrum of the perturbation of the crack front :

A(k, a) = 4a 2 (1 + |ka|) 2 K(k) (12)
For a ≫ a 0 this yields :

     A(k, a) = 4 k 2 K(k), for k = 0 A(0, a) = 4 K(0)a 2 (13)
Taking the inverse Fourier transform of equation ( 12), one gets for a ≫ a 0 :

A(z, a) ∼ 4 π a K(0) (14)
The squared fluctuation σ(z, a) ≡ E [(δa(z, a)δa(0, a)) 2 ] 1/2 can also be derived. For a ≫ a 0 and z ≫ b, b being the correlation length of κ (b = 0 for white noise):

σ(z, a) 2 ∼ 4 K(0)z (15)

Roughness exponent

Suppose that the toughness fluctuation are given by a white noise: K(k) = K 0 = Cst.. Equations ( 13) or (15) gives the value of the so-called wandering or roughness exponent :

ξ = 0.5 (16)
However, numerical simulations performed in [START_REF] Schmittbuhl | Interfacial crack front wandering: influence of quenched noise correlations[END_REF] using also Bueckner-Rice formula [START_REF] Rice | First-order variation in elastic fields due to variation in location of a planar crack front[END_REF] yield two different values: ξ = 0.37 as in [START_REF] Schmittbuhl | Interfacial crack pinning: effect of nonlocal interactions[END_REF] and ξ = 0.5, the first value being obtained if in their simulations δa max > b and the second if δa max < b. The difference hence is due to the numerical parameters used. May be that for δa max > b, the results are in fact illicit since Bueckner-Rice formula, that supposes the perturbation and all is derivatives small, is then no more valid for ∂δa ∂z ∼ δamax b > 1. Moreover, physically it may be more pertinent to suppose that δa max < b when the crack front is trapped, as in the experiments by fluctuations of thoughness (roughness smaller than the typical size of the zones of higher thoughness). Whatever, further investigations are under study to clarify this point. In figure 4, the experimental power spectrum k → A(k, a) evolution is given for several crack positions a. Several remarks can be done:

Comparison experiments and theoretical results

Power spectrum evolution

1. For all the fronts considered, A(k, a) is independent of a. This is in agreement with the behavior for large a obtained in eq. ( 13).

2. The roughness exponent as in previous measures is near ξ = 0.5 in agreement with our theoretical prediction. In figure 5, the experimental evolution of the variance is given. One notices that:

Variance evolution

1. The values for several front oscillates around the mean value during propagation. This corresponds to fluctuations of each position of the front around a steady state situation, due to the heterogeneities.

2. The steady state situation observed is also in agreement with the long time behavior of the variance found in equation (15).

3. The mean value varies linearly with z as predicted by equation (15). Hence one retrieves a roughness exponent of ξ = 0.5.

Conclusion

We considered a semi-infinite crack embedded in an infinite elastic medium subjected to traction line loading on its faces. Using Bueckner-Rice weight function theory for this geometry (elastic line model), the power spectrum and the variance of the crack fluctuations arising from small toughness fluctuations has been derived analytically. In particular, a roughness exponent of ξ = 0.5 has been found in contradiction with some previous numerical works but in agreement with the results found by [START_REF] Schmittbuhl | Interfacial crack front wandering: influence of quenched noise correlations[END_REF] in the case of small roughness amplitude toward the correlation length of the thoughness fluctuations.

In spite of the approximations made in the theoretical modelization (finite geometry replaced by infinite one, uniform displacement loading along the line replaced by uniform force loading), one recovers at least qualitatively the experimental behavior in particular for the steady state situation and for the value of the roughness exponent.

Figure 1 :

 1 Figure 1: Two PMMA plates sintered together: a thick clamped one (21 cm long, 11 cm wide, and 1 cm thick) and a thin one (23 cm × 9 cm × 0.6 cm) subjected to a normal displacement.

Figure 2 :

 2 Figure 2: Experimental picture of the crack front for one equilibrium position. In blue the extracted crack front

Figure 3 :

 3 Figure 3: Semi-infinite plane crack with a slightly perturbed crack front and line tractions at a distance a of the crack front

Figure 4 :

 4 Figure 4: Power spectrum. Marker lines: for several crack fronts with a between 35 mm and 100 mm. Continuous line: A(k) ∝ k -1-2ξ , ξ = 0.5.

Figure 5 :

 5 Figure 5: Variance. Marker lines: for several crack front with a between 36 mm and 58 mm. Continuous line: the mean of all these values


S Santucci

M Grob

Renaud Toussaint

J Schmittbuhl

Knut A Hansen

Jørgen Måløy

K J Måløy

Crackling dynamics during the failure of heterogeneous material: Optical and acoustic tracking of slow interfacial crack growth

Introduction

The failure of heterogeneous materials has been extensively studied in the recent past (see the recent reviews [1,2] and references therein). Those various studies have been motivated by the numerous practical benefits to many engineering fields, but also by the fundamental questions brought forward in particular in statistical physics [1,2]. Indeed, they belong clearly to a more general framework concerning the dynamics of elastic disordered systems, concerning many different physical situations, such as crack fronts [3,4], magnetic domain walls in disordered ferromagnets [5], or wetting contact lines on rough substrates [6], where elasticity and disorder are competing to shape the interface. A fracture front propagating in a 3d heterogeneous medium has itself a 3-d shape with different in-plane and out-of-plane roughnesses, denoted respectively ζ || and ζ ┴ . In order to simplify this problem and shed some light on the interactions between the crack front and material heterogeneities, a simplification to a two dimensional configuration -an interfacial crack, where the crack front is constrained to propagate within a weak heterogeneous plane of a 3d block -has been proposed both experimentally [7,8] and theoretically [4].

In this work, we will focus on this simpler configuration, studying both the structure and the dynamics of those planar cracks. Indeed, using a transparent material and both a high resolution fast camera and two microphones, the detailed complex crack dynamics can be fully captured. Indeed we can follow the propagation of the crack front along a weak heterogeneous plane of a transparent block with a high precision both in time and space. First, we will show that the fracture front displays a complex, jerky dynamics also called -Crackling Noise‖: the fracture dynamics is governed by irregular avalanches with very large size and velocity fluctuations. We characterized the intermittent dynamics observed -the local pinnings and depinnings of the crack front-by measuring both the acoustic emission and the local waiting time fluctuations along the crack front during its propagation. Then, we will focus on the morphology of the interfacial crack fronts -at rest, in order to obtain a high resolution in the description of those interfaces -by studying in details the scaling properties of the front. 2

2.

A complex dynamics: "Crackling Noise"

Experimental set-up and method

We describe here experiments where two Plexiglas plates are annealed together to create a single block with a weak interface, as shown on Fig. 1. (for more details see [7][8][9][10]). The Plexiglas plates have the following dimensions: the main plate is 16 cm long, 14 cm large and 1 cm thick; the other ones are a set of 17 cm long, 0.4 cm thin and 1 cm wide bands. Before annealing, both plates are sandblasted on one side with 200 µm large glass beads. Sandblasting introduces a random topography which induces local toughness fluctuations during the annealing procedure. Another consequence of the blasting process is that the transparency of the Plexiglas plates is lost: light scatters because of the microstructures created. However, after annealing, the newly formed block recovers its transparency, since the contrast of the refraction index along the interface disappears, allowing a direct observation of the crack interface when opening the block in mode 1. Then, we impose a normal displacement at a constant speed (between 1 μm/s up to 100 μm/s) to the lower plate, resulting in stable crack propagation in mode I. A sketch of the experimental set-up is shown in Fig. 2. Using a high speed camera (Photron Ultima) mounted on top of a microscope, up to 12288 consecutive images of the crack advance can be recorded with a spatial resolution of 1024 x 512 pixels (1 pixel ≈ 10 x 10 μm2) at an acquisition rate of 1000 frames per second (fps). It is important to notice that the width of the image acquired by the fast camera is adjusted to the width of the PMMA band, and therefore we can capture the full dynamics of the front. Moreover using two contact microphones (Wide Band acoustic sensors in the ultrasonic range 100kHz -1MHz), placed on top of the thick Plexiglas plate, on the sides of the band, we record simultaneously at an acquisition rate of 1MHz the acoustic emissions during the propagation of the crack.

A complex intermittent dynamics: "Crackling Noise"

On Fig. 2, we show a typical image recorded by the high resolution fast camera during an experiment with an average crack speed ‹V› ≈ 200 µm/s. Moreover, we superimpose on the same figure 6000 fronts recorded by the fast camera at a frame rate of 1000 fps. The spatially random toughness along the weak interface generates a rough crack line in pinning the crack front, and triggers a rich burst activity on a wide range of length scales. The various regions of gray levels (the darker parts show the longer waiting time) reveal this complex dynamics.

We characterized the intermittent dynamics observed, i.e. the local pinnings and depinnings of the crack front by measuring the local waiting time fluctuations along the interface during its propagation [10]. This analysis allowed us to define local velocity at each point of the front during its propagation, and finally avalanches or quakes, as clusters of large local velocities.

On Fig. 3, we show the spatio-temporal map of the local velocity of the front, suggesting once again an intermittent dynamics on a wide range of temporal and length scales. These (optically detected) quakes can be ranked in catalogs with analogous characteristics to seismicity catalogs (epicenter location, time of occurrence and rupture area of each quake). Interestingly, in a recent work [START_REF] Grob | Måløy to be[END_REF], we have shown that the spatial and temporal distributions of the events for our experimental quake catalogs share similar statistical features to those observed for large scale, real seismicity data, suggesting a tight link between this laboratory experimental model and natural events. We also show on Fig 3 . an acoustic signal that was acquired continuously (during the same experiment) at a rate of 1 MHz using a wide-band acoustic sensor in the ultrasonic range 100kHz-1MHz. We clearly observe a qualitative correlation between the acoustic activity and the propagation of the front. Work is in progress in order to build catalogs of acoustic events, pursuing the statistical analysis of both acoustic and optical bursts in order to investigate in detail the correlation between the acoustic emissions and the propagation of the front. In particular, we are currently measuring the acoustic energy released for each optical events and comparing their statistical features.

An ideal experimental system

We have shown here that our experimental set-up allows to fully capture the complex crackling dynamics of a fracture front propagating along a weak heterogeneous interface. We used a double experimental technique developing a simultaneous acoustic and optical tracking of the failure of the block. This experimental study appears relevant both in Earth sciences and statistical physics. Indeed, it is an ideal system for testing the various theoretical description proposed in the literature [12][13][14] in order to clarify the role of the disorder and in particular the elastic interactions between the heterogeneities within the material which give rise to such a complex intermittent behavior.

After studying the crackling dynamics of interfacial crack fronts, let us examine now the morphology of those interfaces.

A complex morphology -Different self-affine behavior at small and large scales-

Motivation

Since the pioneering work of Mandelbrot et al [15] revealing the self-affine character of fracture surfaces of metals, numerous studies have been devoted to the morphological properties of crack surfaces (see the reviews [1,2]). In particular, the roughness exponent characterizing this self-affinity was shown to be very robust and further on conjectured to be Universal [1,16,17]. This apparent Universality led to interpret the fracture of heterogeneous materials within the framework of critical phenomena [1,18]-following the work of Ertas and Kardar for a very different system i.e. the motion of vortex lines in superconductors [19]. In this approach, the failure of heterogeneous material appears as a depinning transition of an elastic line through a disordered landscape, and the crack surface corresponds to the trace left by the front during its propagation.

Although qualitatively appealing, this description failed to reproduce quantitatively the scaling exponents experimentally observed and therefore triggered an upsurge of studies [1]. Despite the plethoric literature concerning the morphology of rough crack surfaces, a clear understanding of the role of material heterogeneities leading to those self-affine long range correlations is still lacking, with up to now an apparent disagreement between the scaling exponents obtained both theoretically and numerically and those measured experimentally [1,2,4,7,8,20,21,22]. In particular, the morphological properties of in-plane crackswhich were proposed to simplify the original 3d problem remain a very controversial issue. Indeed, the in-plane roughness exponent measured experimentally till now ζ exp ~0.6 [7,8,23] is significantly larger than the one predicted either theoretically via renormalization group calculations or numerical simulations ζ th ~0.39 [4,20,21,22]. A percolation-like approach recently proposed [12,13] to account for the local damage ahead of the crack front was able for the first time to reproduce the scaling behavior of those interfacial fronts. However, this description has been tremendously debated [24,25].

Therefore, it appears of crucial importance to clarify this controversy by obtaining a clear and detailed experimental description of the scaling behavior of those planar cracks.

Experimental set-up and method

We will focus here on the morphology of the interfacial crack fronts at rest in order to obtain a high resolution description of the height fluctuations of the interface. The experimental set-up is identical to the one described in the previous section. We just changed the dimensions of the Plexiglas plates (32 cm x 14 cm x 1 cm and 34 cm x 12 cm x 0.4 cm), and the loading procedure: we impose a normal displacement to the lower plate, resulting in a stable crack propagation in mode I; when the front has advanced of a few millimeters (typically 5 mm), the plate is unloaded in order to arrest the front. Then, we take high resolution pictures of the crack interface using a Nikon camera (3871 x 2592 pixels) mounted on a microscope, placed itself on a translation stage that can move in the direction parallel and perpendicular to the propagation plane (x,y). During each loading stop, the table is translated along the front x-axis, and neighboring pictures are taken. The assemblage of up to 15 high resolution pictures of an interfacial crack front at rest provides a large front description with around 16 000 data points and a pixel size of a = 1 µm. This procedure was repeated 20 times in order to obtain 20 different and independent fracture fronts. The total length of each analyzed crack was around 16 mm, with a typical example shown on Fig. 1. Finally, in order to modify the heterogeneities within the weak plane of the block, we prepared different samples using glass beads of different diameters during the sand-blasting procedure, respectively Ø 1 ~ 50 µm and Ø 2 ~ 200 µm.

A complex self-affine behavior

We aim at characterizing the scaling properties of the height fluctuations of the interfacial crack fronts Δh(δ) = h(x+δ) -h(δ). On Fig. 4, we plot the standard deviation σ of the height fluctuations Δh(δ) of the crack interfaces as a function of the scale δ for samples roughened with glass beads of diameter Ø 2 ~ 200 µm during the blasting procedure. At small scales, below δ * ~ 100 µm, we observe a self-affine scaling behavior σ(δ<δ * ) ~ δ ζ-with a roughness exponent around ζ -~ 0.63 which is in agreement with previous observations [7,8] and the percolation-like approach recently proposed [12,13] to account for the local damage ahead of the crack front. However, interestingly, at scales larger than δ * ~ 100 µm, we observe a cross-over to another scaling regime σ(δ > δ * ) ~ δ ζ+ with a smaller roughness exponent ζ + ~ 0.37. It is important to emphasize that this is the first time that such a scaling regime is reported for in-plane crack fronts. Moreover, this value appears in agreement with the prediction of a depinning model [4,20,21,22].

In order to analyze in more details this scaling behavior at large scales, we define height fluctuations at a given scale as Δh(δ,x)=h(x+ δ)-h(x), and propose to examine the statistical distribution of the height fluctuations P[Δh(δ)] on that range of scales δ > δ * ~100µm . It is important to remark that our experimental data set is large enough to perform such a direct reliable measurement. In Fig. 2, we show the distributions of the height fluctuations for logarithmically increasing length scales δ. On that range of scales, above δ * , δ>δ * ~100µm, we clearly observe that the distributions P[Δh(δ)] follow a Gaussian statistics. The self-affine character of the height fluctuations of the fronts appears through the scaling behavior of the standard deviation σ of those distributions P[Δh(δ)], σ ~ δ ζ+ with a unique roughness exponent ζ + ~ 0.37.

Effect of the disorder and the material microstructure

Let us now examine and discuss the crossover length scale δ* separating the large and small scale behavior of the height fluctuations. In order to check the effect of disorder and material microstructure on the scaling behavior of the interfacial crack fronts, we modified the heterogeneities within the weak plane of the block, preparing different samples using glass beads of different diameters Ø 1 ~ 50 µm and Ø 2 ~ 200 µm for the blasting procedure. We show on Fig. 4. the scaling behavior of the interfacial crack fronts for those different samples. We observe mainly the same features as on previous figure with the two different scaling behaviors separated by a characteristic size respectively δ * Ø1 ~ 35 µm, and δ * Ø2 ~ 100 µm. It is important to note that the scaling regime at small scales up to the characteristic sizes δ * Ø1 and δ * Ø2 can be associated to the material heterogeneities. Indeed, when using larger glass beads for roughening the Plexiglas plates, the amplitude of the height fluctuations of the fronts increases (vertical shift of the blue curves) as well as the scaling range at small scales, providing a roughness exponent ζ -~ 0.6 up to the scale δ * Ø2 ~ 100 µm. It is important to notice that those results appear by analogy, in agreement with recent experimental observations on 3-d fracture surfaces of phase separated glasses [26]. However, more work is needed, in order to obtain a characterization of the disorder produced by the indentation impacts during the blasting process and finally to be able to compare this typical size δ * (separating the two different scaling regimes) with the typical scale of the heterogeneities. In particular, it appears important to study the influence of the crack growth velocity, the mechanical properties of the plates on this typical size δ * and its eventual link to the damage process zone [27] .

Conclusion

Finally, we can conclude by insisting on the fact that our high resolution experimental set-up allows to obtain a detailed description of the failure of heterogeneous materials. Indeed, we can fully capture the crack dynamics by observing directly the interface with a high resolution fast camera, and recording simultaneously the acoustic emissions during the front propagation along a weak heterogeneous plane of a transparent Plexiglas block. The role of heterogeneities during crack propagation is of central importance since they induce local pinnings of the crack front and subsequently trigger a very complex burst activity, characterized by avalanches, spanning over a broad range of temporal and length scales. This intermittent dynamics also called -crackling noise‖ is also suggested from the acoustic emission accompanying the failure process. Another consequence of this phenomenology is the roughness of fracture surfaces left by the crack. Indeed, cracks in heterogeneous media exhibit a selfaffine morphology, with long range correlations. In particular, we have shown here that the interfacial crack fronts follow different scaling regime at small and large scales. For scales larger than a typical size δ * , that could be interpreted as a characteristic size of disorder, we observed a self-affine regime with a roughness exponent ζ + = 0.35 ± 0.05 much smaller than the value measured up to now and observed here at small scales ζ -= 0.6 ± 0.03. Those different scaling regimes can be understood within the framework of recent statistical models [28].

Figure 1 :

 1 Sample preparation procedure. The Plexiglas plates are sand-blasted and sintered together to form new single blocks. (a) The blasting procedure creates roughness fluctuations on the plate surface that control the local toughness during the annealing process. It also induces a loss of transparency by unpolishing the surfaces. Transparency is however recovered thanks to the sintering process. (b) Sintering process: in white, the PMMA sample; in gray, the two aluminum plates; in black, the solid frame with the screws to apply a torque of = 10 N.m. The whole block is put in an oven at 205 o C during 30 minutes.

Figure 2 :

 2 Figure 2: Sketch of the experimental setup (top left: side view; top right: view from above). The thick plate of the sample (in white) is clamped into a rigid aluminum frame (in black). A cylindrical rod imposes a displacement of the thin Plexiglas band at a constant speed. The crack propagates in the y direction and is observed from above using a high speed camera mounted on a microscope. Moreover, using two microphones, placed on top of the thick Plexiglas plate on the sides of the band, we can record simultaneously the acoustic emissions during the propagation of the crack. The figure shows a typical image recorded by the high resolution fast camera during an experiment with an average crack front speed ‹V› ≈ 200 µm/s. We superimpose 6000 fronts recorded at an acquisition rate of 1000 fps.

Figure 3 :

 3 Figure 3: (a) Spatio-temporal of the local velocity of the crack front during its propagation within the weak heterogeneous plane of the Plexiglas block at an average speed of ‹V› ≈ 200 µm/s. (b) Acoustic emissions continuously recorded by one contact microphone, output signal in volts (Arbitrary Units) (Wide Band acoustic sensor 100kHz -1MHz) at 1MHz during the crack propagation. (c) The average crack front speed shows an intermittent dynamics.

Figure 1 :

 1 Figure 1 : Assembling up to 15 high resolution pictures of an interfacial crack front at rest provides a large front description with around 16 000 data points and a pixel size of 1 µm.

Figure 2 :

 2 Figure 2 : Scaling behavior of interfacial crack fronts at rest. We plot the standard deviation σ of the height fluctuations Δh(δ) of the crack interfaces as a function of the scale δ. The glass beads used for preparing the samples -blasting and roughening of the Plexiglas plates -had a diameter of Ø 2 ~ 200 µm.

Figure 3 : 2 ) 1 / 2 .

 3212 Figure 3 : Statistical distributions of the height fluctuations Δh(δ)-‹Δh(δ)›, for logarithmically increasing length scales δ. We observe on that range of scales (δ ≥ 100 µm) that the statistical distributions P[Δh(δ)-‹Δh(δ)›] follow a Gaussian shape, with a width increasing with the range of scales studied δ. The self-affine behavior of those distributions enters through the scaling behavior of the standard deviation σ of the height fluctuations Δh(δ)-‹Δh(δ)›, σ ~ δ ζ+ , with a roughness exponent ζ + ~ 0.37, as we can observe on the inset top-right. The Gaussian statistics is confirmed in the figure top-left where we observe a collapse of all the various distributions when plotting on a semilog scale P[Δh(δ)-‹Δh(δ)›].(2πσ 2 ) 1/2 as a function of [Δh(δ)-‹Δh(δ)›]/(2πσ 2 ) 1/2 . The blue dotted line represents the curve y=exp(-x 2 ) and fits perfectly the experimental rescaled distributions.

Figure 4 :

 4 Figure 4 : Effect of disorder on the self-affine scaling of interfacial crack fronts. We plot the standard deviation σ of the height fluctuations Δh(δ) of the crack interfaces as a function of the scale δ. We analyzed two different sets of data, prepared with glass beads of different diameters respectively Ø 1 ~50µm and Ø 1 ~200µm.


Jan Ludvig Vinningland

Eirik Øistein Johnsen

Grude Flekkøy

Renaud Toussaint

Knut Jørgen Måløy

Eirik G Flekkøy

Granular Rayleigh-Taylor instability

A granular instability driven by gravity is studied experimentally and numerically. The instability arises as grains fall in a closed Hele-Shaw cell where a layer of dense granular material is positioned above a layer of air. The initially flat front defined by the grains subsequently develops into a pattern of falling granular fingers separated by rising bubbles of air. A transient coarsening of the front is observed right from the start by a finger merging process. The coarsening is later stabilized by new fingers growing from the center of the rising bubbles. The structures are quantified by means of Fourier analysis and quantitative agreement between experiment and computation is shown. This analysis also reveals scale invariance of the flow structures under overall change of spatial scale.

Introduction

Improved understanding of granular flows would be of essential importance to a range of industrial applications, to the study of geological pattern forming processes, and, in general, to the theoretical description of disordered media.

As grains become smaller the effect of the interstitial fluid becomes more important. The result is a combination of dry granular dynamics and the hydrodynamics of the fluid. These systems give rise to a variety of exotic and most often poorly understood phenomena such as fluidization [START_REF]Physics of Dry Granular Media[END_REF] and bubble instabilities [START_REF] Gendron | Bubble propagation in a pipe filled with sand[END_REF], quicksand and jet formation [START_REF] Lohse | Creating a dry variety of quicksand[END_REF], and sandwich structures in systems where different particle types segregate [START_REF] Zeilstra | Simulation study of air-induced segregation of equal-sized bronze and glass particles[END_REF]. While the study of dry granular media has been extensive over the past decades, the exploration of fluid-granular systems has been of more limited scope.

A wide range of granular instabilities where various structures form along fluid-grain interfaces have been reported the last few years [START_REF] Johnsen | Pattern formation during air injection into granular materials confined in a circular hele-shaw cell[END_REF][START_REF] Gendron | Bubble propagation in a pipe filled with sand[END_REF]. Notably, the patterns formed by grains falling in a highly viscous liquid were investigated experimentally and theoretically by Völtz et al. [START_REF] Völtz | Rayleigh-taylor instability in a sedimenting suspension[END_REF]. The instability reported by Völtz shares its main qualitative characteristics with the classical Rayleigh-Taylor instability, i.e. a single dominating wavelength growing right from the start. In contrast, the instability discussed here arises along an air-grain interface where falling granular fingers and rising air-bubbles form [START_REF] Vinningland | Granular rayleigh-taylor instability: Experiments and simulations[END_REF][START_REF] Vinningland | Experiments and simulations of a gravitational granular flow instability[END_REF]. A coarsening of the finger-bubble pattern is observed right from the start of the instability. The coarsening subsequently stabilizes due to the continuous formation of new fingers at the centers of the rising and growing bubbles.

Our experimental and numerical data compare favorably both qualitatively and quantitatively, despite the simplifications of the numerical model: Only two spatial dimensions, zero interparticle friction, no gas inertia, and a resolution of the gas-grid of 2.5 grain diameters. The finger-bubble structures further exhibit scale invariance under change of particle size, a feature which is supported both by observations and theoretical considerations.

The experiment is presented in Sec. 2, and the numerical model is outlined in Sec. 3. The analysis and the consistency check are presented in Section 4, together with the scale invariance observed when the grain diameter increases. The conclusion is given in Section 5.

Experiment

The experimental setup, illustrated in Fig. 1, consists of a closed Hele-Shaw cell mounted on a hinged bar which enables the cell to rotate approximately 130 degrees from a lower to an upper vertical position (from A to B in Fig. 1). The Hele-Shaw cell is made of a 1 mm thick silicone frame sandwiched between two 8 mm thick glass plates. The internal dimensions of the cell are 56 mm × 86 mm × 1 mm, and it is filled with polystyrene beads and air at atmospheric pressure. The cell is rotated manually and it takes about 0.2 seconds to bring the cell to an upright, vertical position. The off-center pivot of the cell causes the rotation to slow down the falling motion of the grains due to centrifugal forces. A bar is mounted at position B in order to have some control of the final vertical position of the cell. The evolution of the instability is recorded by a high speed digital camera (Photron Fastcam-APX 120K) taking images with a resolution of 512x512 pixels at a rate of 500 frames per second.

Monodisperse polystyrene beads of 140 µm in diameter (Microbeads ® Dynoseeds TS 140-51) are used in the experiments. The filling of the cell is performed with one glass plate lying down horizontally with the silicone frame adhered on top. Small portions of beads are carefully deposited on the plate and leveled with the frame before the upper plate is attached and fasted with clamps. The cell is flipped a few times after closure to allow the grains to form a random loose packed configuration before the cell is rotated.

The humidity in the lab is important in order to keep the electrostatic and cohesive properties of the grains at a suitable level to prevent the grains from clustering or sticking to the glass plates. During the filling of the cell and throughout the experiment the humidity was kept constant at about 30%.

A series of images from the experiment is shown in Fig. 2. Due to experimental imperfections the initial air-grain interface in Fig. 2(a) is not a well-defined, horizontal front. This is mainly caused by three factors: (1) Os-nodes grain Fig. 3. Illustrating the smoothing function: A single grain i and its four neighbouring grid-nodes are shown. The four subdivisions of the cell-volume (i.e. the area between the stapled lines) correspond to the fraction of vi (or mi) that contributes to the continuous u(x, y) (or ρ(x, y)) fields at the opposite nodes. The gray-scale of the subdivisions indicates which node they belong to.

cillations in the cell caused by the abrupt stop as it reaches the vertical (and hits the bar), ( 2) the falling of grains during the rotation, and (3) the granular Boycott effect during rotation (i.e. grains and air pass each others in layers parallel to the plates) caused by the fact that up to 7 layers of beads may fit between the plates. However, all these experimental effects are non-persistent and in Fig. 2(b) the first finger-bubble structures are clearly visible. As the instability evolves a coarsening of the structures is observed. However, new fingers emerging from the center of the rising bubbles (see Fig. 2(c)) prevent the size of the structures to grow indefinitely and gives rise to a characteristic bubble-size. The coarsening is driven by fingers merging and forming an inverted Y as is indicated in Fig. 2(c).

Simulation

The numerical model, first introduced by McNamara et al. [START_REF] Mcnamara | Grains and gas flow: Molecular dynamics with hydrodynamic interactions[END_REF][START_REF] Anghel | Erratum: Grains and gas flow: Molecular dynamics with hydrodynamic interactions[END_REF], combines a continuum description of the air with a discrete description of the granular phase. The granular phase is a collection of rigid spheres that constitute a deformable porous media described by coarse grained solid fraction ρ(x, y) and granular velocity u(x, y) fields on the grid where (x, y) are the two dimensional space coordinates. The continuous fields ρ(x, y) and u(x, y) are obtained from the positions and velocities of the individual grains by means of a linear smoothing function: To calculate u(x, y) a fraction of v i is distributed to each of the four neighbouring grid-nodes of grain i. This is illustrated in Fig. 3.

The continuum gas phase is described solely by its pressure P (x, y). The inertia of the gas, and hence its velocity field, is not considered. This is justified for small particle Reynolds numbers as is the case for our system. The pressure is governed by the equation [START_REF] Mcnamara | Grains and gas flow: Molecular dynamics with hydrodynamic interactions[END_REF][START_REF] Anghel | Erratum: Grains and gas flow: Molecular dynamics with hydrodynamic interactions[END_REF] where φ(x, y) = 1 -ρ(x, y) is the porosity, κ(x, y) the permeability, u(x, y) the granular velocity field, and µ the gas viscosity. This equation is derived from the conservation of air and grain mass, using Darcy's law [START_REF] Darcy | Les Fontaines Publiques de la Ville de Dijon[END_REF] to obtain the pressure drop over a volume with permeability κ. The Carman-Kozeny relation [START_REF] Carman | Fluid flow through granular beds[END_REF] is assumed for the permeability, and the isothermal ideal gas law is assumed for the compressible gas. See Refs. [START_REF] Mcnamara | Grains and gas flow: Molecular dynamics with hydrodynamic interactions[END_REF][START_REF] Anghel | Erratum: Grains and gas flow: Molecular dynamics with hydrodynamic interactions[END_REF] for a detailed derivation. The empirical Carman-Kozeny relation is given by

φ ∂P ∂t + u • ∇P = ∇ • P κ(φ) µ ∇P -P ∇ • u , (1)
κ(ρ, d) = d 2 180 (1 -ρ) 3 ρ 2 , ( 2
)
where ρ is the local solid fraction, and d is the diameter of the grains, and 1/180 is an empirical constant valid for a packing of beads. The grains are governed by Newton's second law

m dv dt = mg + F I - V ∇P ρ , (3)
where m, v, and V are respectively the mass, velocity, and volume of the grain. Contact dynamics [START_REF] Radjai | Force distributions in dense two-dimensional granular systems[END_REF] is used to calculate the interparticle force F I which keeps the grains from overlapping. The dynamics of the grains are simplified by neglecting particle-particle and particle-wall friction. A lower cutoff is imposed on the solid fraction because the Carman-Kozeny relation is not valid as the solid fraction drops below 0.25 [START_REF] Zick | Stokes-flow through periodic arrays of spheres[END_REF]. This cutoff causes the permeability of the most dilute regions of the system to be slightly lower than the true permeability. The effect is a slight overestimation of the pressure forces acting on the grains in the dilute regions. Figure 4 shows a series of snapshots from a simulation with the same grain diameter, grain density, and cell dimensions (except for the cell thickness since the simulation is 2D) as in the experiment shown in Fig. 2. A qualitative comparison of the image-series in Figs. 2 and4 renders the simulation and experiment consistent in many respects: The sizes of the bubbles and the fingers are comparable, and the dynamical processes of finger merging and finger nucleation are observed in both cases. Nevertheless, some discrepancies are observed, particularly at the start and toward the end of the instability. The initial differences are mainly attributed to imperfections of the experimental setup: In the experiments the instability start from a homogeneous sheet of grains rather than from a sharp air-grain interface as in the simulations. Further differences are observed toward the end of the instability: In the numerical snapshots bigger bubbles reach the upper surface before the smaller bubbles (see e.g. Fig. 4(f)). In the experiment, however, the row of bubbles stays more or less horizontally aligned from start to finish. We believe this is an effect of the zero grain-wall and grain-grain friction used in the simulation.

σ k = k 2 S(k) S(k) -k 2 k k (1/cm) S(k) (b)

Results

To study the coarsening of the observed structures quantitatively we use the Fourier spectrum of the solid fraction field to obtain a characteristic average wave number k by the following procedure. The discrete Fourier transform and the power spectrum of each horizontal line of ρ(x, y) is calculated. The averaged power spectrum, S(k), is then obtained by averaging over all these horizontal power spectra. An average wave number is defined

as k = k S(k) • k/ k S(k)
, where 1/k is the wavelength. Likewise, we define the squared standard deviation

σ 2 k = k S(k) • k 2 / k S(k) -k 2 .
For the experimental data the image pixel values are used to estimate the solid fraction. Figure 5(a) shows ρ(x, y) for t = 0.09 s and three power spectra obtained at different vertical positions (indicated by black lines). Figure 5(b) shows the averaged power spectrum S(k) and mean wave number k of the solid fraction field in Fig. 5(a).

Figure 6 shows the temporal evolution of k and σ k (inset) for the numerical and experimental data. An additional set of experimental data is added to the plot. The numerical curve starts out with a significantly higher wave number than the experimental curves. This behaviour is however expected given the difference in initial conditions. After about 0.2 seconds the experimental and numerical data coincide and continue to show a similar coarsening behavior. The fluctuations of k and σ k are associated with the continuous nucleation of new fingers and merging of existing fingers.

We further investigate the behavior of the system as the overall spatial scale is changed: Keeping all length ratios and the particle number fixed, the size of the system will scale according to the diameter d of the grains. We measure the characteristic inverse length scale k as d is changed and observe a scale invariance of the evolution. A series of seven simulations are performed where d varies from 70 µm to 490 µm in steps of 70 µm. The dimension of the numerical cell confining grains of 70 µm in diameter is 28 mm x 34 mm. To compare, a series of experiments using polystyrene beads of 80, 140, 230, and 570 µm in diameter, confined in Hele-Shaw cells that scale proportionally with d in all directions, are performed.

In these simulations we have introduced the larger density of glass, rather than polystyrene in order to minimize the numerical artifacts associated with the solid fraction cutoff in the permeability. As the inertia of the grains decreases the overestimated pressure forces in the dilute regions of the system will act even stronger on the falling granular fingers.

Data-collapse plots of the rescaled mean wave number d k are shown in Figs. 7(a) (simulation) and 7(b) (experiment). These plots indicate that the characteristic size of the structures is invariant when size is measured in units of d; the number of grains that spans the width of the bubbles is the same for a wide range of grain sizes.

Theoretically, the scale invariance of the product d k may be interpreted as follows: Compared to the other terms of Eqs. ( 1) and ( 3) the m dv/dt, F I and P ∇ • u terms may be shown to be small [START_REF] Vinningland | Size invariance in granular flows: simulation and experiment[END_REF]. For that reason, these equations exhibit an approximate invariance under system size scaling. If we take δP to be the pressure deviation from the background pressure, express the velocity of grain i as v i = δv i + u 0 and the locally averaged granular velocity as u = δu + u 0 , where u 0 is the constant sedimentation velocity of a close packed system, this scaling may be expressed as x → λx, δP → λδP , u 0 → λ 2 u 0 , δu → λδu and κ → λ 2 κ, where λ is a scale factor. The structure formation of the system is governed by δu and, since this velocity scales the same way with λ as the length scales themselves, the evolution of any structure measured in units of d will be scale invariant. In particular this is true for the structures measured by the length 1/ k , and so d k is scale invariant. However, the invariance deteriorates both when particle size is increased, and when it is decreased. In the first case, the relative effect of granular inertia is increased, in the second, the relative effect of the P ∇ • u term is increased.

The convergence of the numerical data-collapse in Fig. 7(a) is quite good. The deviation of the 70 µm curve for small t is probably explained by the increase in the relative importance of the P ∇ • u-term. The divergences of the 350, 420, and 490 µm curves for greater t in the same plot arise because the bubbles in the coarser packings disappear before they reach the surface due to the increase of u 0 with λ 2 [START_REF] Vinningland | Size invariance in granular flows: simulation and experiment[END_REF]. The experimental data in Fig. 7(b) have a wider distribution but collapses satisfactorily given the standard deviation error bars. The experimental data are obtained by averaging over three experiments for each diameter d. The standard deviation is calculated over a time window of 0.3 seconds. The accuracy of the experiments is at its lowest during the initial coarsening of the structures, but as the mean wave number stabilizes around 0.2 seconds the accuracy improves. Nevertheless, the data points are, with a few exceptions, within a distance of one standard deviation from one another. The loss of precision for small times is most likely caused by the inaccuracy involved with the manual rotation.

Conclusion

In conclusion, we have presented experimental and numerical results of a gravity-driven granular flow instability which is significantly different from its classical hydrodynamic analog. The simulations reproduce the characteristic shape and size of the experimentally observed structures and provide fine patterns in the early phase of the process that are not resolved experimentally. Data-collapse plots of the mean wave number k indicate that the flow and the resulting structures are invariant when measured on a scale proportional to the grain diameter d for a range of diameters that spans from 70 µm to 570 µm.

Fig. 1 .

 1 Fig. 1. Front and side view of the experimental setup. Two cell positions are superimposed in the side view to illustrate the manual rotation from position A to B.

Fig. 2 .

 2 Fig. 2. Experimental images of the granular Rayleigh-Taylor instability where a layer of monodisperse polystyrene beads of 140 µm in diameter (black) displaces air (white) in a Hele-Shaw cell of dimensions 56 mm × 86 mm × 1 mm.

Fig. 4 .

 4 Fig. 4. Snapshots from a simulation where the input parameters are determined by the experiment in Fig. 2. The grain density is 1.05 g/cm 3 , grain diameter is 140 µm, and the cell dimensions are 56 mm × 68 mm.

Fig. 5 .

 5 Fig. 5. (a) Power spectra superposed on the solid fraction field from which they are obtained. (b) The averaged power spectrum of the solid fraction field shown in (a) with the mean wave number k indicated.

Fig. 6 .

 6 Fig.6. Mean wave number k and standard deviation σ k (inset) for two experiments and one simulation, all using polystyrene beads of 140 µm in diameter.

Fig. 7 .

 7 Fig. 7. Data-collapse plot of d k for a series of (a) simulations using glass beads, and (b) experiments using polystyrene beads. The grain diameters d are given in the legend box. The inset shows the evolution of the unscaled k .


Keywords: Christian Salim, Data Reduction based energy-efficient approaches for secure priority-based managed wireless video sensor networks, Networking and Internet Architecture [

Data Reduction based energy-efficient approaches for secure priority-based managed wireless video sensor networks

Christian