Experimental study of stable imbibition displacements in a model open fracture . I . Local avalanche dynamics

Xavier Clotet,1,2 Jordi Ortı́n,1,* and Stéphane Santucci2,† 1Departament ECM, Fac. de Fı́sica, Universitat de Barcelona, C. Martı́ i Franqués 1, 08028 Barcelona, Catalonia, Spain 2Laboratoire de physique, CNRS UMR 5672, École Normale Supérieure de Lyon, 46 Allée d’Italie, 69364 Lyon Cedex 07, France (Received 5 September 2015; revised manuscript received 30 November 2015; published 28 January 2016)


I. INTRODUCTION
Fluid flow through disordered media occurs in geological, agricultural, and industrial processes of capital importance.It often involves the presence of an interface separating two different phases, when a fluid originally residing in the disordered medium is displaced by an invading fluid.In secondary oil recovery, for instance, the front dynamics between the two immiscible phases (oil and water) determines the efficiency of the recovery process.Printing, coating, filtering of chemicals and contaminants, impregnation, and soil irrigation are other examples of two-phase fluid displacements in disordered media [1][2][3][4][5].The invasion process depends first on the relative ability of the two fluids (displacing and displaced) to wet the walls of the disordered medium.When the invading fluid wets preferentially the medium, the displacement-favored by capillary forces-is called imbibition.Second, the displacement may be either stable or unstable, depending on the relative viscosity of the fluids involved.The displacement is stable when the displacing fluid is more viscous than the displaced fluid.Front disturbances, due to capillary pressure fluctuations at the pore scale, cannot grow because the viscous pressure gradient on the side of the displacing fluid is larger than on the side of the displaced fluid [6,7].
Stable-imbibition displacements in disordered media have received a lot of attention from a fundamental perspective, since they represent an interesting realization of the more general problem of the motion of a slowly driven front in a random potential.In this context, the competition between stabilizing and destabilizing (random) forces on different time and length scales leads to the unbounded growth of correlations and to the emergence of scale-invariant asymptotic properties [8,9].In stable-imbibition displacements, specifically, the stabilizing forces are the surface tension along the front and the change in viscous pressure drop across it.In contrast, capillary pressure fluctuations and permeability variations play the role of a quenched disorder, responsible for the front distortions.In the absence of gravitational forces (e.g., for horizontal displacements) the main features of these flows are controlled by the capillary number, a ratio of viscous to surface tension forces.The capillary number is defined as Ca = μv/σ , where μ is the dynamic viscosity of the invading fluid (much larger than the dynamic viscosity of the receding fluid), v is the average velocity of the invading front, and σ is the interfacial tension between the two fluids.Large-scale front fluctuations are damped preferentially by the fluid viscosity, while small-scale fluctuations are damped by surface tension.The condition of mass continuity implies that these two stabilizing mechanisms cross over at a characteristic length scale c = √ κ/Ca, where κ is the permeability of the disordered medium.This length scale was predicted theoretically and numerically to determine also the lateral extent of correlations of the front distortions [10,11].We verified this prediction experimentally for constant flowrate displacements at different v, and we could show, moreover, that this correlation length c plays a crucial role in controlling the dynamics of the invasion process [12].Indeed, the front motion takes place by avalanches, localized front displacements in which the front velocity is very high and spatially correlated, with a maximum lateral extent set by c .The limit v → 0 corresponds to a critical depinning transition, characterized by diverging lateral correlations c ∼ 1/ √ v and an infinite susceptibility of the system to front distortions in the thermodynamic limit.Avalanches in stable-imbibition motions have been studied experimentally [12][13][14][15][16][17], theoretically [10,[18][19][20], and numerically with phase-field simulations [18][19][20].
We present here a large ensemble of new experimental results on the dynamics of stable imbibition displacements.In the present experiments, a two-dimensional model of an air-filled disordered medium is invaded by a viscous silicone oil at constant flow rate.Displacements are horizontal, so gravitational effects are absent.The experiments are performed in the statistically stationary state of saturated front roughness, in which the spatial fluctuations of the front height h(x,t) do not depend on t in average.In this paper (Part I) we focus our attention on the spatiotemporal behavior of the instantaneous velocities v(x,t) ≡ ∂h/∂t at a local scale, defined by our spatial resolution.A corresponding study of the dynamics through increasing length scales, up to the lateral system size L, is reported in the accompanying paper (Part II) [21].
In comparison to previous works we have expanded significantly the space of experimental parameters: We have used five silicone oils of very different dynamic viscosity and increased the range of driving velocities to span more than one decade in the lowest accessible range.In this way we have been able to investigate imbibition displacements in a very wide range of capillary numbers, 6 × 10 −5 < Ca < 2 × 10 −3 .We have also made an important effort in improving the statistics: For each experimental condition, defined by (v,μ), we have performed 15 to 20 experiments with different disorder realizations.Improved statistics has led to very robust scalings and to accurate values of the various critical exponents that characterize the local spatiotemporal dynamics in this kind of displacements.Accurate values of complete sets of critical exponents are most valuable to test new and existing theoretical approaches of imbibition front dynamics.The main effort in the present work has been to use invading fluids of various dynamic viscosities.Indeed, by changing systematically the dynamic viscosity of the invading fluid-while keeping the oil-air surface tension practically unchanged-we have been able to show the distinct effect of the flow rate v and the viscosity μ on the invading process.Specifically, we have verified the scaling c ∼ 1/ √ Ca independently for v and μ, but also we have demonstrated that the temporal correlations of local front velocities-and, consequently, the duration of local avalanches-do not scale identically with v and μ and therefore are not controlled only by Ca.
This article is organized as follows.The experimental setup and the experimental methods are described in Sec.II.The front dynamics is analyzed at the (local) scale of our spatial resolution in Sec.III, and the corresponding avalanche behavior is quantified in Sec.IV.Finally, in Sec.V, we discuss the results and draw the conclusions of this work.glass plate.The topography consists of copper patches of size 0.4 × 0.4 mm 2 and height b = 0.06 mm.These patches are randomly distributed in a square grid, occupying 35% of the surface.Adjacent patches form islands of disorder.Their lateral sizes are exponentially distributed, with a characteristic length d = 0.6 mm [12].

Our
The fluid is injected through the medium at constant flow rate by means of a syringe pump connected to the inlet, displacing the air initially present.Lateral edges of the HS cell are sealed in order to allow fluid flow in the y direction only.The imposed flow rates explored produce mean velocities in the range 0.04 < v < 0.6 mm/s.The fluids used are silicone oils of dynamic viscosities μ = 10,50,100,169,350 cP, density ρ 1000 kg/m 3 , and surface tension around σ = 21 mN/m in all cases.We thus explore a wide range of capillary numbers Ca = μv/σ between 6 × 10 −5 and 2 × 10 −3 and nominal correlations lengths along the front 2.6 < c < 15.3 mm.All experimental conditions explored and the corresponding values of Ca and c are quoted in Table I.For each given experimental condition (v,μ) we have performed 15 to 20 TABLE I. Experimental conditions: μ is the dynamic viscosity of the fluids, v the imposed mean velocity, Ca their corresponding capillary numbers, c the nominal correlation length along the front, and c /L the fraction of the lateral system size that is correlated.experiments with different disorder realizations.This leads to very large data sets, required for an accurate statistical analysis.The oil-air front advancement is recorded from the top using a high-speed, high-resolution camera.The spatial resolution is r = 0.106 mm/pixel.The image acquisition rate is adjusted (up to 200 frames/s) to ensure that the front locally advances one pixel at most between consecutive images.The region of interest (ROI) spans L = 136 mm in the transverse direction x (smaller than the cell width to avoid boundary effects) and 25 to 45 mm in the propagation direction y.
We base our study on the analysis of the local waiting time fluctuations along the front during its propagation.This method was originally developed to study the dynamics of interfacial cracks [22] and successfully applied to characterize the burstlike dynamics observed in forced-flow imbibition [12][13][14].It consists on computing the amount of time wt(x,y = h(x,t)) spent by the front at each position (x,y) of the ROI (image recorded).The front lines h(x,t) are obtained by applying an edge-tracking algorithm on the recorded frames.The map of local velocities is then computed as v(x,y) = r/wt(x,y), where r is the spatial resolution.Each local velocity determined in this way is actually the mean value of the local front velocity at the pixel resolution.

III. LOCAL VELOCITIES
We perform initially a qualitative analysis of the properties of the local velocities.Figure 2 displays local velocity maps for three experiments performed at the same imposed mean velocity v = 0.13 mm/s but different viscosities μ = 10, 50, 350 cP.In all cases regions of high velocities (above v) alternate with regions of low velocities.We clearly observe that these velocity fluctuations around v become larger at smaller viscosities.In particular, regions of local velocities v(x,y) v become bigger, and pinned regions with v(x,y) v clearly appear.Besides, low-viscosity fronts are rougher than highviscosity fronts.Reducing v instead of μ produces equivalent effects, suggesting that Ca controls indeed the fluctuations of v(x,y).

A. Statistical distributions of v(x, y)
Figure 3 shows the distributions of local velocities rescaled by their mean value P (v/ v ) for different μ and v.We choose this quantity in order to compare experiments performed at very different v-covering more than one decade.First, local velocities are widely distributed, spanning a large range of more than three decades.In all cases, the distributions show a steep decay at large values of v, close to (v/ v ) −3. 5 .In experiments performed at large v or μ, the distributions have a maximum slightly below the mean value.However, for lower v and μ the distributions get wider: The probability of observing extreme velocities, both much smaller and much larger than the average velocity, increase.This was already noticed qualitatively in Fig. 2. The distributions develop a power law around v well fitted by (v/ v ) −1. 25 for almost two decades in the experiments of lowest viscosity and velocity.
This systematic evolution of the distributions P (v/ v ) with μ and v is analyzed in terms of the capillary number Ca ∼ μv in Fig. 4. The distributions shown correspond to four different values of Ca, including the two extreme values studied Ca = 6.3 × 10 −5 and 2.2 × 10 −3 .Remarkably, for different experimental conditions but similar Ca, the distributions P (v/ v ) collapse on the same curve.This provides strong evidence that the capillary number controls the statistical behavior of the local velocities of the front.
We can also study the statistical distributions of the normalized fluctuations of the local velocities, defined as Y = (v − v )/σ , where σ stands here for the standard deviation of P (v).The probability density functions (pdf) for all experiments are shown in Fig. 5. Again, all of them show a peak around the mean value and an abrupt decay at large v, as previously observed in Fig. 3. Interestingly, now all data collapse for Y > 1, that is, when fluctuations are one σ larger than the mean velocity of the experiment.We will make use of this result to analyze local avalanches in Sec.IV.

Spatial correlations along the front
Lateral spatial correlations can be studied by computing the autocorrelation function of the local velocities along the front: v( ) is the average velocity within the window of size , • • • is the average over all windows ∈ [0,L − ] along the interface, and • • • is the ensemble average over interfaces recorded at different times.
Figure 6 shows the systematic evolution of C v with imposed flow rate for a given viscosity (left panel) and with viscosity for a given mean velocity (right panel).We observe that the autocorrelation of the local velocities decays with towards a minimum below zero-i.e., maximum anticorrelation-at * , and then increases again towards zero, until the signal becomes fully uncorrelated at large .The evolution of the length scale * is plotted as a function of the nominal correlation length c in Fig. 7 (solid symbols).The length scale min at which the signal becomes fully uncorrelated is also displayed in the figure (open symbols).For all experimental conditions (μ,v) explored * and min are proportional to c .The values of * are systematically smaller than c while the values of min are systematically larger than c .This result confirms that the capillary number Ca = μv/σ controls the lateral spatial correlations of the interface, through the characteristic length scale c ∼ 1/ √ μv.
Equivalent information can be extracted from the analysis of the statistical width W v ( ) of the fluctuations of v over windows of size along the front.W v ( ) increases with until it saturates above a crossover length scale * as observed in Fig. 8.The transition between the scaling and the saturated regimes is smooth, making more difficult to extract * from this plot than to obtain the characteristic length * from C v .However, we observe that the crossover length scale * evolves similarly with μ and v, and the results are compatible with * ∼ c .Nevertheless, we observe that μ and v play an opposite role on the amplitude of the velocity fluctuations, W v ( ): Fluctuations are larger the smaller μ but also the larger v.

Spatial correlations in the direction of propagation
We study also the correlation of v(x,y) in the direction of propagation of the front.We do not consider the temporal information of the position of the front explicitly.Instead, we analyze directly the correlation of the velocities of the interface in y, and compare it with the underlying quenched disorder.Figure 9 displays the autocorrelation functions C v ( y) for different oil viscosities (solid lines).Again, in all cases the correlation decays with y until it reaches a minimum negative value (maximum anticorrelation) and then it increases and We can compare these results to the underlying disorder in the cell.We consider the disorder landscape as a matrix with elements d(x,y) equal to 1 on the disorder patches and 0 otherwise.Then the autocorrelation function of d(x,y) in the direction of propagation of the front (dashed line in Fig. 9) can be compared to the autocorrelation function C v of the fronts (solid lines in the same figure).In tune with the nominal lateral size of single disorder patches (0.4 mm) and with the characteristic lateral size of the clusters of disorder patches ( d = 0.6 mm), the disorder is highly correlated up to 0.3 mm; the autocorrelation of d(x,y) decays abruptly beyond this value, and it becomes zero at y d = 0.54 ± 0.05 mm.This value is compared to the experimental values y * and y min in the figure inset.A strong influence of the disorder patches on the local dynamics of the interface in the propagation direction y is clearly observed: The maximum anticorrelation is achieved at y * = y d , indicating that the largest fluctuations of the local velocity occur within the characteristic size of the disorder patches d .However, the relaxation of those large fluctuations is controlled by the capillary number.

Temporal correlations
We now analyze the temporal correlation of local velocities by computing the autocorrelation of the signal v(x,t) in time and taking the ensemble average over different positions x along the front: v(t) is the average velocity within the temporal window t, • • • t is the average over all the temporal windows of v(x,t) such that t ∈ [0,t max − t], and • • • is the ensemble average over velocity signals analyzed at different positions x and different realizations of the experiment.C v ( t) presents a similar qualitative behavior (not shown here) than the spatial autocorrelation function displayed in Fig. 6: C v decays as t increases, until it reaches a minimum value below zeromaximum anticorrelation-at t * .For larger time lags, local velocities become uncorrelated in time and C v → 0.
Figure 10 displays the evolution of t * for all experimental conditions explored as a function of the capillary number Ca ∼ vμ.Although the range of imposed velocities explored (v max /v min = 15) is smaller than the range of viscosities (μ max /μ min = 35), t * varies much more with v than with μ.Those results make clear that velocity and viscosity do not play the same role in controlling t * -unlike what happens with * , as discussed before.The dashed line goes through experiments performed at different v but the same μ and is proportional to a power law v −1.4 .Experiments at the same v but different μ correspond to the dotted line, which is FIG.10.Maximum anticorrelation time t * as a function of the capillary number Ca. Dashed and dotted lines are guides to the eye proportional to v −1.4 and μ −0.4 , respectively.Open symbols correspond to experiments at v = 0.13 mm/s.proportional to μ −0.4 .These scaling relations are discussed in Sec.V in connection with the statistical properties of the local avalanches studied in the next section.

IV. LOCAL AVALANCHES
Our purpose in this section is to study the spatiotemporal dynamics of imbibition fronts in terms of local avalanches.We define local avalanches as spatially connected clusters of velocities v(x,y) higher than a given velocity threshold v c .The threshold used is v c = v + cσ v , where v is the mean value of the local velocities in the experiment [i.e., the temporal and spatial average of v(x,t)], σ v is its standard deviation, and c is an arbitrary constant called the clip level.
Local avalanches are characterized by their size A, given by the surface area of the cluster of high velocities; their duration D, given by the time interval elapsed from the moment that a first interface enters the cluster of high velocities to the moment that a last interface leaves the cluster; and their extension, L x and L y , in the corresponding directions x and y.These definitions are represented graphically in Fig. 11.

A. Nontrivial role of the disorder
Thanks to the precise control of the disorder in our setup, and the high spatial resolution of the acquired images, we are able to compare directly disorder patches with local bursts of velocity.Figure 12 displays the clipped velocity map corresponding to an experiment with v = 0.13 mm/s and μ = 50 cP, on top of the disorder landscape explored by the front.Velocities above v c = v + 2σ v are shown in white, while velocities below are not displayed.Hence the white clusters observed in the figure correspond to local avalanches as defined earlier.Disorder patches-copper islands that reduce the gap spacing of the cell-are represented in black and disorder-free regions in orange (dark gray).Since the disorder matrix does not coincide with the matrix of pixels of the image, some pixels are only partially filled by disorder patches.The pixel color code (gray scale) ranges from orange (dark gray) (no disorder patch) to black (disorder).The bottom panel shows the whole velocity map v(x,y) of the experiment, and the top panel on the left shows a close-up of a single avalanche.We observe that the geometry of the local disorder bears no direct correlation with the clusters of high velocities.To quantify this effect we measure the amount of disorder swept by the front during a local avalanche, i.e., the amount of copper patches within a cluster of fast motion.The top panel on the right of Figure 12 shows the areas of disorder swept by local avalanches versus the areas of the corresponding avalanches.The most probable values of these magnitudes are proportional to each other with a proportionality coefficient 0.35; this means that the amount of disorder below a local avalanche is 35% in average.This value corresponds precisely to the nominal filling fraction of the disorder patches in our model.We conclude that avalanches occupy any region of the cell with no particular preference for the underlying disorder patches.

B. Statistical characterization of local avalanches
In this section we analyze the statistical properties of sizes, durations, and shapes of local avalanches.We first study the of the threshold v c on the statistical distributions.In earlier experimental work on imbibition avalanches the threshold for the local velocities was defined as ṽc = c v [12].However, we have shown in Fig. 3 that the distributions of v/ v for different experiments do not collapse, meaning that the underlying statistics of v/ v differs for each experimental condition (v,μ).If we instead consider the normalized variable Y = (v − v )/σ v , the distributions of Y practically collapse for Y > 1 for all the experimental conditions explored (Fig. 5).The choice v c = v + cσ v , with c > 1, is thus more appropriate to consider local avalanches of statistically equivalent local velocities.It allows to compare experiments performed at very different experimental conditions (μ,v).
The statistical distributions of the various properties of the local avalanches-size, duration, and lateral extents-can be approximated by the general form P x (x) = a x x −m x e −x/ξ x , a power law with an exponential cutoff.
First, for given experimental conditions μ and v, we verified that the power-law exponent, α A for sizes and α D for durations, is independent of the clip level within the range 1 c 3 explored.The clip level, however, does play a systematic role in the cutoff of the distributions.As the clip level c increases, the sizes and durations of the high-velocity clusters decrease, as well as the power-law range of the distributions.Given that c does not modify the scaling exponent, we have selected c = 2 in all subsequent analysis in order to compare experiments performed with different (μ,v).
Second, for fixed clip level, the cutoffs of the power-law distributions depend on the controlling parameters of the experiment.This dependence will be studied in detail later.When distributions with different cutoffs are available, a good procedure to obtain a reliable value of the power-law exponent of the distributions is to look for their best collapse in the appropriate reduced variables.This method is described in Refs.[23,24].It has been applied here to obtain the power-law exponents of the distributions of A, D, L x , and L y .For some particular experimental conditions, however, the probability distributions display only the exponential decay.The corresponding data have been discarded from the data set used to extract the corresponding power-law exponent.

Distributions of sizes and durations
The size of a disorder patch introduces a natural lower cutoff of avalanche sizes.For this reason we consider as avalanches only those high-velocity clusters that span at least one disorder patch.This implies that the minimum avalanche size is A min = 0.4 × 0.4 mm 2 .Figure 13 shows the distributions of sizes for various experimental conditions.P (A) follows a power law with a decaying cutoff.The power-law exponent does not depend on the mean velocity of the front nor on the viscosity of the fluid, as observed on the left and the right panels of the figure, respectively.In all cases the power-law regime is represented by the same exponent, This exponent is compatible with the value reported earlier by our group for one single oil (μ = 50 cP) [12].However, it differs from the value α A = 1.54 ± 0.05 obtained in phasefield simulations [19,20].The scaling range of the distributions and thus the cutoffs, however, do depend on v and μ.The evolution of ξ A with v is quite clear and systematic.The slower the fluid is pushed through the medium, the bigger the avalanches can be.The evolution of ξ A with viscosity is less clear.In general we observe that the smaller μ the bigger the avalanches can be.
The distributions of durations of local avalanches are displayed in Fig. 14 for the same experimental conditions of the previous figure.Again, the distributions show a power-law regime with an exponential cutoff.The power-law exponent does not depend on the experimental conditions: In this case the exponent differs from the value obtained in earlier experiments where the statistics of durations was poorer [12].The result from phase-field simulations again differs considerably, α D = 1.62 ± 0.06 [19].
The effect of the controlling parameters (v,μ) on the avalanche durations is manifested in the evolution of the cutoffs, which measure the maximum duration of the bursts of fast motion.ξ D evolves systematically with the flow rate v, as clearly observed in the left panel of Fig. 14.The faster the front, the shorter the avalanches.On the contrary, the evolution of ξ D with viscosity is not obvious.

Morphology of avalanches
The shape of local avalanches is not isotropic.Highvelocity clusters are larger in the lateral direction x than in the direction of front advancement y, as can be observed in perpendicular and parallel to the advancement of the front, are displayed in Fig. 15 for various experimental conditions.Both P (L x ) and P (L y ) are power-law distributed with exponential cutoffs.The power-law exponents α L x = 0.83 ± 0.10 and α L y = 1.08 ± 0.15 do not depend sensibly on the experimental conditions.The exponent of P (L x ) is in very good agreement with earlier experimental results with only one viscosity [12].Contrary to the power-law exponents, the range in which these distributions follow a power law depends on both v and μ.The evolution of the cutoffs ξ L x and ξ L y with experimental conditions will be discussed at the end of the section.
In order to quantify the anisotropic shape of the local clusters of high velocities we analyze the joint distribution P (L x ,L y ).We expect a scaling relation L y ∼ L H x between the most probable values, where the exponent H is a measure of the anisotropy.The joint distributions P (L x ,L y ) for all experiments studied are shown in Fig. 16.The left panel displays the average of the 2D histograms of the distributions for each experiment.We observe a scaling range where L y ∼ L H x .At small values of L x the distributions clearly bend.This is a consequence of windowing, resulting from the fact that we consider only avalanches larger than the size of a disorder patch A d = 0.4 × 0.4 mm 2 .Not taking into account avalanches of A < A d implies that very narrow avalanches-small L x -are under-represented and, therefore, must be disregarded.We choose to consider as lower cutoff L x = 1.2 mm, i.e., two times the characteristic length of the medium heterogeneities ( d = 0.6 mm).Bending of the joint distribution is also observed at large values of L x .We choose as upper cutoff L x = 30 mm, corresponding to 2 c for the largest c explored.Data within these two limits are represented using full symbols in the figure, while the excluded points are shown in open symbols.To obtain a more reliable value of H , we consider the joint distribution of L x and L y rescaled by their cutoffs ξ L x and ξ L y .The right panel of Fig. 16 shows the 2D histogram corresponding to all the avalanches within the range of L x chosen.Fitting those data gives the exponent: This result is in agreement with the value reported earlier for the same kind of imbibition displacements using one single oil viscosity [12].The uncertainty now is smaller thanks to the improved statistics.

Maximum sizes and durations of avalanches
The cutoffs ξ A , ξ D , ξ L x , and ξ L y are representative values of the maximum size, duration, and extent of local avalanches in the x and y directions.
On experiments with one single oil, Santucci and coworkers observed that ξ L x ∼ 1/ √ v [12].We extend this result now by including the effect of viscosity.The dependence of ξ L x on capillary number is shown in the top left panel of Fig. 17.Given that Ca ∼ μv, this representation of the data allows us to compare the dependence on both μ and v on the same plot.We observe that the cutoffs diverge as Ca → 0 and that the evolution is compatible with ξ L x ∼ 1/ √ Ca.In the bottom left panel of the figure we compare these cutoffs (solid symbols) with the distance * of maximum anticorrelation of the front (open symbols) obtained in Fig. 6. ξ L and * are plotted versus the nominal value of c .Both quantities are approximately proportional to c , reflecting that the lateral extent of local avalanches is indeed bounded by the correlation length c .However, the actual values are slightly larger (ξ L x ) or smaller ( * ) than the nominal c , as a consequence of the dependence on the clip level of the former or on the particular criterion chosen to define the latter.At any rate, the result reinforces the notion that μ and v play the same role in controlling the lateral extent of the avalanches.
It is worth noting that for the largest v and μ (largest Ca) the smallest avalanches may be under-represented in Fig. 15.In these cases the correlation length can be as small as 2.6 mm.Given that the maximum lateral extent of avalanches is bounded by c , we expect a large number of tiny avalanches.Because of the constraints arising from the disorder patches, that impose a minimum size A d = 0.16 mm 2 , and data windowing, we only consider events whose lateral extent is L x > 1.2 mm.Hence the expected tiny avalanches are disregarded.Taking into account that the distributions are normalized, the occurrence of large events is overestimated and the distributions are distorted.This applies to the distributions of sizes, durations, and lateral extents, whose cutoffs are presumably overestimated for the largest capillary numbers.
Having established that ξ L x ∼ c ∼ 1/ √ Ca, Santucci and coworkers proposed several other scaling relations between the cutoffs of the distributions of different properties of the local avalanches [12].We check these scaling relations again on our data, thereby extending their validity to a larger range of velocities and also to different viscosities.
From the anisotropic shape of the avalanches, L y ∼ L H x , we infer that the cutoffs of the corresponding distributions will also scale as ξ L y ∼ ξ H L x , and thus ξ L y ∼ H c ∼ Ca −H/2 ∼ Ca −0.41 .The cutoffs ξ L y versus Ca are shown on the top right panel of Fig. 17.The data are compatible with the expected scaling.ξ L y can be compared also to the correlation length in the direction of propagation y min obtained in Fig. 9.These quantities are plotted versus H c in the bottom right panel of Fig. 17. ξ L y and y min are proportional to H c , showing that the local dynamics is strongly controlled by the avalanches.
The size of local avalanches can be approximated by the product of their lateral dimensions, A ∼ L x L y .Its maximum value, therefore, is expected to scale as ξ A ∼ ξ L x ξ Ly ∼ ξ 1+H L x ∼ Ca (1+H )/2 ∼ Ca −0.91 .The evolution of ξ A with Ca is shown in the left panel of Fig. 18.The result is consistent with the scaling proposed.
Finally, we consider the maximum durations of the local avalanches.The average velocity of the front, v, is a lower bound for the local velocity in an avalanche, and ξ L y is an upper bound of the corresponding displacement in the direction of propagation.Hence, the maximum duration can be estimated as clusters (sizes, durations, and lateral extensions) are power-law distributed, with cutoffs that diverge as Ca → 0. Such behavior is a signature of the proximity to a depinning transition, where avalanches are expected to be scale invariant at the critical point [11].Close to the critical point avalanches are also expected to be scale invariant, but only up to the cutoff scale imposed by a finite correlation length c .In this framework, the following scaling behavior is expected [12,20]: where γ av = 1+H z av and z av is the avalanche dynamic exponent.The top and middle panels of Fig. 20 show the joint distributions P (L x ,D) and P (A,D) in double-logarithmic scale.In both cases the distributions have been rescaled by their cutoffs to include all experiments.We observe a clear scaling of A/ξ A with D/ξ d (top panel).The exponent obtained from the fit is γ av = 1.8 ± 0.1.This value agrees within error bars with the experimental result obtained in Ref. [12], γ av = 1.6 ± 0.2.However, it differs significantly from the exponent γ av = 1.28 ± 0.04 obtained in Ref. [19] from phase-field simulations.The scaling of L x /ξ L x with D/ξ D gives an exponent z av = 1.10 ± 0.12 (middle panel of Fig. 20).The value is in very good agreement also with previously reported results in Ref. [12], z av = 1.1 ± 0.1.A value z av 1 corresponds to a ballistic lateral propagation of the avalanches L x ∼ D. Again, the measured value differs from the one obtained in phase-field simulations, z av = 1.59 ± 0.15, halfway between a ballistic and a diffusive propagation of avalanches [19].Using P (A), P (D), and P (L x ) can be obtained [19]: Unfortunately, it is difficult to validate these scaling relations, since the measured exponents α L x , α A , α D are close to 1.

V. DISCUSSION AND CONCLUSIONS
We have analyzed the spatiotemporal dynamics of constant flow-rate imbibition fronts in a model open fracture at the local scale of our spatial resolution [25].A wide range of experimental conditions has been explored, using different oils with dynamic viscosities μ = 10 to 350 cP and imposed flow rates that produced mean velocities of the front ranging from v = 0.036 to 0.55 mm/s.Corresponding capillary numbers ranged from 6 × 10 −5 to 2 × 10 −3 , and the lateral correlation length along the front took values 3 < c < 15 mm.The large number of experiments for each set of parameters (μ,v) with different disorder realizations has led to large data sets, allowing an accurate statistical analysis.
This new set of experiments has led us disclose the distinct effect of the imposed flow-rate v and the dynamic viscosity μ of the invading fluid on the imbibition process.In particular, we have shown that the burstlike dynamics of the imbibition front is not controlled only by the capillary number Ca.
Local velocities of the front fluctuate around the imposed mean value due to the heterogeneities of the medium and are widely distributed.We have shown that the range of amplitudes of the velocities and the specific shape of the distributions are controlled by the capillary number Ca.In particular, the local front velocities span a wider range as the critical pinningdepinning transition is approached by decreasing Ca.
In Sec.III B we have studied the large spatial and temporal correlations of the local velocities.In spite of being at the origin of front avalanches, local velocity correlations had not been studied experimentally before.We have shown that c ∼ 1/ √ Ca for both μ and v, so these two parameters play the same role in controlling the spatial correlations.The dependence of the temporal correlations, however, is much stronger on the mean imposed velocity than on the viscosity of the fluid.This is consistent with the behavior of the maximum duration of the local avalanches.The latter scales as H c /v ∼ μ −H/2 v −(1+H/2) , which corresponds indeed to the observed dependence on μ and v of the temporal correlations.
We have observed (Fig. 9) that the largest fluctuations of v(x,y) in the direction of fluid invasion are found within d , the characteristic linear size of disorder islands.Interestingly, this is consistent with our recent observation that the temporal range of intermittency of the spatially averaged velocity of the front is determined by d [13].
We have found that the statistical properties of local velocity bursts-sizes, durations, and lateral extents-are power-law distributed, with an exponential cutoff that diverges as Ca → 0, a limit that corresponds to a pinning-depinning transition of stable imbibition displacements.In the regime studied here, where front fluctuations are induced by capillary disorder, the power-law exponents are independent of Ca.The validity of this result, established earlier in Ref. [12], has been extended to a larger range of v but more importantly to different μ.The measured power-law exponents are α A = 1.09 (8), α D = 1.03 (10), α L x = 0.83 (10), and α L y = 1.08 (15).The anisotropy of the velocity clusters, L y ∼ L H x , is characterized by an exponent H = 0.82 (2).These values fulfill the expected scaling relations between exponents, within error bars, and validate our previous experimental results [12] with less uncertainty.This is important for comparing with theoretical predictions.However, they differ significantly from the values obtained in phase-field simulations [19,20].The reason for the discrepancy might lie on the disorder properties.In the experiment, single disorder patches have a finite linear extension (0.4 mm, about four times our images resolution).In phase-field simulations, in contrast, the disorder is pointlike, i.e., of the same size as the underlying discrete lattice.
We have verified that local bursts of activity do not correlate trivially with islands of disorder patches.Indeed, the disorder patches swept during a local avalanche occupy 35% of the avalanche area, a value coincident with the nominal area fraction of disorder.This result also reflects the fact that the front dynamics on a single site is highly dependent on the dynamics at all other sites.Non-locality makes forcedflow imbibition intrinsically different (and substantially more complex) from most other dynamic problems involving slowly driven fronts in disordered landscapes.
To conclude, slow and stable imbibition displacements in open fractures take place with strongly correlated local velocities in both the lateral direction along the advancing front and in the direction of propagation.Spatiotemporal correlations lead to localized velocity bursts, or avalanches, whose scale-invariant properties arise from the proximity to a critical depinning transition at Ca = 0.The present work provides an exhaustive investigation of local velocity correlations and avalanche properties of stable imbibition displacements in a laboratory model of an open fracture.Our results shed light on the physics that govern this kind of displacements.They also provide accurate quantitative results that will be useful to contrast future theoretical predictions.
FIG. 1. Sketch of the model open fracture.The liquid is driven at constant flow rate from the inlet(I) displacing the previously resident air.The model open fracture consists of two horizontal, parallel glass plates (G) and a disorder fiber-glass plate (FG).Squared copper patches (Cu) of lateral sizes 0.4 × 0.4 mm 2 are randomly distributed occupying 35% of the surface of the FG plate.The region of interest (ROI) defines the lateral system size L = 136 mm.

FIG. 2 .
FIG. 2. Velocity maps of three experiments performed at the same v = 0.13 mm/s but different μ = 10, 50, 350 cP (top, middle, bottom).The color code (gray scale) bar for v(x,y) is the same in all cases.The bottom-right white line in each panel gives the nominal correlation length c = 15.3, 6.7, 2.6 mm for each experimental condition.

FIG. 9 .
FIG. 9. Autocorrelation of local velocities in the propagation direction y for experiments with different oil viscosities.The spatial correlation of the disorder is also shown (dashed line).Inset: y * (open symbols) and y min (full symbols) vs Ca.The dotted line corresponds to y d (the length scale of decorrelation of the disorder).The dashed line is a guide to the eye proportional to Ca −0.4 .

FIG. 11 .
FIG. 11.Definition of local avalanches.Bottom: Velocity map after clipping.Local velocities v < v c are all plotted in blue (dark gray).The first and last interfaces in contact with the avalanche highlighted are shown.Top left: Close-up of the selected avalanche.Top middle: Definition of area, A, and lateral sizes, L x and L y , of the avalanche.Top right: Close-up of the same avalanche on the activity map v(x,t).The avalanche duration is D. Experimental conditions: v = 0.13 mm/s, μ = 50 cP.

FIG. 16 .
FIG. 16.Left: Joint distributions of P (L x ,L y ).Results for all experiments are shown, but different velocities are not distinguished.Full symbols represent data used to obtain H . Vertical dotted lines show the lower and upper cutoffs of the range of L x considered.Right: Rescaled joint distributions of the data within the scaling range.The power-law exponent is obtained from a linear fit of the log-log data.The color code (gray scale) corresponds to the density of events and increases from dark blue (dark gray, outer region) to red (dark gray, inner region).

FIG. 17 .
FIG. 17. Top panels: Evolution of the cutoffs of P (L x ) (left) and P (L y ) (right) with capillary number in double-logarithmic scale.Dashed straight lines are guides to the eye with the exponents shown in the figure.Bottom left: Cutoffs ξ Lx and maximum lateral anticorrelation lengths * vs the nominal correlation length c .Bottom right: Cutoffs ξ Ly and decorrelation lengths in the direction of propagation y min vs H c , where H = 0.82 is the anisotropy exponent.Dashed straight lines are to the eye.The dotted line has slope correspond to different viscosities.
FIG. 18. Evolution of the cutoffs of P (A) (left) and P (D) (right) with the capillary number in log-log scale.Dashed and dotted straight lines are guides to the eye with the exponents shown in the figure.Different symbols (either open or solid) correspond to different viscosities.The values for experiments performed at v = 0.13 mm/s have been distinguished by using open symbols.

FIG. 19 .
FIG. 19.Cutoffs of P (D) (solid symbols) and maximum anticorrelation time t * (open symbols) vs H c /v, an estimation of the maximum duration of local avalanches.The dashed line is a guide to the eye.Different symbols correspond to different viscosities.For the same viscosity, slower experiments have larger H c /v.