# Dynamical scaling of imbibition in columnar geometries

###### Abstract

Recent experiments of imbibition in columnar geometries show interfacial fluctuations whose dynamic scaling is not compatible with the usual non local model governed by surface tension that results from a macroscopic description. To explore this discrepancy, we exhaustively analyze numerical integrations of a phase-field model with dichotomic columnar disorder. We find that two distinct behaviors are possible depending on the capillary contrast between both values of disorder. In a high contrast case, where interface evolution is mainly dominated by the disorder, an inherent anomalous scaling is always observed. Moreover, in agreement with experimental work, the interface motion has to be described through a local model. On the other hand, in a lower contrast case, interface is dominated by interfacial tension and can be well modeled by a non local model. We have studied both spontaneous and forced-flow imbibition situations, giving a complete set of scaling exponent in each case, as well as, a comparison to the experimental results.

###### pacs:

47.56.+r, 68.35.Ct, 05.40.-a## I Introduction

Fluid–fluid displacements in porous media is a subject of much interest in industrial processes and material characterization so as in environmental problems ranging from petroleum recovery and irrigation to retention of waste waters Sahimi-95 (); AL04 (). We restrict our analysis to the case in which the invading fluid, that wets preferentially the medium, is more viscous than the displaced resident fluid giving rise to compact rough interfaces. This process can take place either spontaneously at constant pressure, spontaneous imbibition, or by the application of an external pressure at constant injection rate, forced-flow imbibition AL04 (). The surface can be represented by a single–valued function of position and time , Barabasi-Stanley (). In many cases, when fluctuations have a thermal origin or can be reduced to that, the interfacial fluctuations are self-affine and follow the dynamic scaling of Family–Vicsek Family-Vicsek-1985 (). Then, a complete description of the dynamical growth is possible with only two independent exponents. However, when disorder is relevant enough to interfere with the geometry of the interface, a more generic scaling can apply, and one more independent exponent is necessary to reach the dynamical description. It is the so-called anomalous scaling RA00 (), and it has been observed in many different experimental and numerical situations during the last decade LE93 (); LO98 (); MA97 (); SO02b (); SO03 (); SO05 ().

Several experiments of imbibition in distinct geometries have been proposed in the last years. There are experiments that use paper as a disordered medium, Buldyrev-92-I (); Horvath-95 (); Kwon-96 (); Zik-97 (); Balankin-2000 () and others performed in Hele–Shaw cells (two parallel glass plates separated by a narrow distance) with a random distribution of glass beads as a disordered medium Rubio-89 (); Horvath-91 (); He-92 (). Using the same geometry, other methods of generating a disordered medium have been explored, including random variations in gap spacing produced by a predesigned surface relief of the bottom plate SO02b (); SO03 (); SO05 (); SO02 () or by roughened plates Geromichalos-2002 (). Focusing our attention on columnar geometries we have presented in previous works SO02b (); SO03 () experimental studies of forced-flow imbibition in a Hele–Shaw cell with a columnar quenched disorder, produced by dichotomic variations in the thickness of the bottom plate. We found that the interfacial dynamics followed an intrinsic anomalous scaling with varying exponents that were incompatible with the expected results from the usual macroscopic model, that gives characteristic exponents of a non local model. Instead, we obtained good agreement with an heuristic model of diffusively coupled columns presenting local interactions. The essential change of behavior, which pass to be dominated by non local to local interactions, is obviously due to the persistence of the columnar disorder. However the detailed physical mechanism is not clear.

In the present paper we use numerical integrations of a phase-field model with columnar disorder to explore this behavior. By imposing the columnar geometry of the interface motion at the macroscopic model, we analytically derive the heuristic local equation presented in Ref. SO02b (). In addition, we show how both local and non local behaviors are observed in the numerical model through the variation of a parameter modeling the contrast between the different capillary values present in the system. The outline of the paper is as follows. In Sec. II we introduce the phase-field model and the interfacial equations obtained in the sharp interface limit. In Sec. III we show the notions of generic scaling used to characterize the interfacial dynamics. Section IV is devoted to analyze the numerical results in the low capillary contrast case for both spontaneous and forced-flow imbibition situations. Section V deals with the case of high capillary contrast and its connection with local growth models. Finally, in Sec. VI we discuss the physical contents of phase-field and macroscopic models as well as their relevance to explain the experimental results.

## Ii The phase-field and macroscopic model

Our numerical results shall be performed using the so-called phase-field model DU99 (); EL01 (); Aurora-EPL-01 (). Such a model is based on the introduction of an order parameter which can take two limit values, representing the two phases liquid/air of the system. The phase-field dynamics is controlled by a conserved equation based on a Ginzburg-Landau formulation, , where is the chemical potential and the free energy is given by the functional . The phase-field equation then reads as

(1) |

where is a parameter that is taken constant in the liquid phase and zero in the air phase. is a potential taken as and defines two stable phases through the double well potential; the destabilizing linear term accounts the effect of a capillary force that makes the interface to advance. The effect of an inhomogeneous capillarity is added by using a dichotomic capillary noise with the values:

(2) |

We consider a columnar disorder defined by single tracks of lateral size distributed along the -direction in such a way that tracks with the high disorder value occupy 35% of the system length . This is the same kind of disorder reported in the experimental work of Ref. SO02b (). Note that, for a given disorder realization, tracks wider than are obtained when two or more unit tracks are placed adjacently.

### ii.1 The macroscopic description of imbibition

The use of phase-field models to reproduce imbibition experiments is based on the ability to get the same results than those obtained from a macroscopic model. Indeed, in the so-called sharp interface limit , a matched asymptotic expansion of the field around a kink solution of Eq. (1), , can be performed, recovering then the basic macroscopic equations for the usual pressure , velocity , interfacial curvature and the columnar capillary disorder which we called :

(3) | |||

(4) | |||

(5) |

The macroscopic variables and parameters are defined from the phase-field formulation as HM03 ():

being the first order term of the expansion on
of the chemical potential and playing the role of a interfacial tension.
The three equations of the macroscopic model are well known from
phenomenological arguments involving conservation laws. Darcy’s law, Eq. (3), arises from an averaging procedure of Navier-Stokes equations at low Reynolds number, when the geometry of the Hele-shaw cell is imposed. The Laplace equation, Eq. (4), comes from imposing incompressibility of the liquid, and the Gibbs-Thomson relation, Eq. (5), comes from a principle of minimum interfacial energy. The capillary pressure at the interface can be expressed as , where corresponds to the contact angle, and is the distance between the plates of Hele-Shaw cell.

In the experimental work reported in Ref. SO02b (), the random gap distribution is constructed by using a fiberglass substrate with a pattern of copper tracks attached to the bottom plate of a Hele-Shaw cell. Tracks have a thickness with a lateral size of and are distributed along the lateral direction without overlap. Therefore, the gap of the Hele-Shaw cell has a dichotomic variation with two possible values and . Although in that case the capillarity is a 3-D effect of the cell, we can relate our numerical parameters of Eq. (2) as and . Note that the parameter is related to the capillary contrast between both values of disorder.

### ii.2 Spontaneous and forced-flow imbibition

In our study, we will consider both situations of spontaneous and forced-flow imbibition by choosing conveniently the boundary conditions into the phase-field model LA05 (). For spontaneous imbibition an applied constant pressure is imposed at the origin of the cell . In contrast, a pressure gradient has to be imposed at the origin, , to reproduce forced-flow imbibition. Therefore, the main difference between both cases is found on the mean velocity of the interface. While in the forced-flow case, the interface evolves with the imposed constant velocity , in the spontaneous imbibition case, the averaged interfacial height follows the so-called Washburn’s law . An exact expression of such an evolution can be obtained solving the following equation:

(6) |

that comes from (3), (4), and (5) with
and . The expression for reads
where . Therefore, spontaneous imbibition has a slowing-down dynamics with a mean velocity .

### ii.3 Equations for the interface

From the macroscopic equations, Eqns. (3), (4) and (5), it is possible to obtain an equation for the moving interface by means of a Green analysis. The Green function in our case follows a Poisson equation with a unit source at , evaluated at the plane 2-D. The expression for the interface is obtained using the Green identity BA89 ():

(7) |

integrated over the volume of the liquid , being the interfacial position. Depending on the intensity of the quenched noise two cases are relevant in our study.

#### ii.3.1 Low capillary contrast. Linearized equations

We assume that the quenched noise does not impose any special geometry and the interface can be linearized around its averaged value, , with . Keeping the first order on fluctuations in Eq. (7), and imposing Darcy’s law, Eq. (3), at the interface boundary condition, one can obtain the linearized interface equation AL04 (); DU99 (); LA05 (), which in Fourier space reads as

(8) |

where we have supposed that correlations do not grow up faster in time than the mean height of the interface, . Under this limit, interface fluctuations follow the same Eq. (8) both for spontaneous and forced-flow imbibition LA05 (); Aurora-EPL-01 (); PA03 (). However, it is worth to mention here that the presence of Washburn’s law () in spontaneous imbibition give rise to dynamic crossover lengths and therefore, as it has been pointed out in Ref. Pradas07 (), a rich variety of different scaling regimes can be observed. The crossover length scale can be explicitly seen as a balance between the surface tension term and the drift term DU99 ()

(9) |

In the forced-flow case, is constant and thus the crossover length is just a static length scale separating two different regimes. However, as it has been observed in several numerical results DU99 (); LA05 (); intrinsic-superrough (), this crossover length acts as a cutoff for the interface fluctuations growth due to the interface being asymptotically flat on length scales larger than .

#### ii.3.2 High capillary contrast. Coupled channel equations

The situation becomes quite different when the capillary contrast is increased. In the experimental work of Ref. SO02b (), where forced-flow imbibition is studied for high capillary contrast, it is observed that the interface motion can be modeled through a phenomenological local equation. Our purpose now is to derive such an equation directly from the macroscopic model. To do this, we shall assume that the noise is so strong that the interface adopts the columnar geometry of the disorder. The procedure consists then on integrating Eq. (7) in a closed surface along the profile of the i-th effective channel. Here, an effective channel, or simply channel , is defined as a group of single adjacent tracks of which the majority has the same disorder value, in such a way that the interface advances as a compact surface through each channel. An example is depicted in Fig. 1, where the channel is composed of several tracks of which the majority has a high disorder value (gray tracks) and therefore, the averaged disorder of the channel is larger than the mean disorder of the whole system . Note that the surrounding channels and must have a mean disorder value . The width of the channel is defined as . A numerical example of these channels can be seen in Fig. 2.

We are considering the forced-flow imbibition case with the following boundary conditions at the top and the bottom of the channel:

(10) | |||

(11) | |||

(12) |

where the pressure at the origin has been estimated by using the relation , which comes essentially from Eq. (11), taking that is the mean capillary pressure of the whole system. It means that the pressure at the origin is changing on time in order to get a mean constant velocity for the whole interface. Using these boundary conditions, and taking the general expression for the two-dimensional Green function

we can evaluate the different boundary integrals of the top and bottom segments of Eq. (7) as

where , with being a cut off due to the finite size of the system. We have also supposed enough wide channels to ensure that , which means that we are taking the initial times, before interface gets saturated. The nomenclature used for the Green function means that it is evaluated at the interface, . Therefore, Eq. (7) can be written as

(13) |

The two first terms of the right-hand side are due to the flow between neighbor channels. We have defined the parameter . In order to get an equation for the time evolution of the interface , we define the following coupling coefficients between channels as the ratio between channel flow and height differences

The coupling variable has to be understood as a diffusion coefficient which depends on each channel . As a general case, it may be taken as a random variable. Moreover, we shall assume that does not vary in time during the initial times, before interface gets saturated.

Then, we can write Eq. (II.3.2) as an inhomogeneous diffusion equation between channels:

(14) |

where corresponds to the pressure at the interface. We are assuming that the pressure at the interface is time-dependent, which is based on the experimental results reported in Ref. SO02b (). In such experimental work, the local velocity of the interface at each channel follows an expression similar to Washburn’s law until it reaches the saturation value . In order to take into account such behavior into the equation, we consider that the pressure at the interface can be expressed as , where the pressure difference is given by the usual Gibbs-Thomson relation, Eq. (5), taking a negligible atmospheric pressure. The term is an effective kinetic term due to the local capillary forces at each channel, being an arbritrary constant, and can be explained in terms of the mass conservation. When the interface goes through a channel of a high capillary disorder (), its local velocity tends to initially increase up to a nominal value. In contrast, since we are imposing a constant velocity for the whole interface, the local velocity at the neighbor channel with a lower capillary disorder tends to decrease down to a nominal value. After reaching such nominal value in both cases, the local velocity decrease or increase asymptotically to the saturation value following the Washburn’s behavior due to the capillary forces of each channel. Therefore, we are taking that for , where is a random variable defined as , which takes the values at the channel with the highest capillary value (), and at the channel with the lowest capillary value (). In addition, we also suppose that the curvature of the interface can be approximated as a constant value , and only its sign depends on each channel as .

Then, rewriting the last term on the right-hand side of Eq. (14), we get the final expression

(15) |

where we have defined the new constants as , and .

We thus conclude that when the capillary disorder is large enough, the columnar geometry of the system leads to a local description for the interface motion.

## Iii Generic scaling laws and the characterization of the interfacial dynamics

Rough interfaces grow exhibiting power laws in both a horizontal correlation length , that accounts the range of correlation, and a vertical growth length, like the interface width . is defined as the deviation of the height as (where and mean sample and spatial average respectively), and are the so-called dynamic and growth exponents, which completely characterize the growth of self-affine processes. The saturation of the surface occurs in a saturation time , when the correlation length reaches the system length . Above this time, interfacial width scales as , being the roughness exponent which it is related to the others exponents through the scaling relation . To study local growth, we shall define local widths averaged on windows of size , . Local growth is also given as a power law for and for , where is the local roughness exponent. In the case of self-affine growth, global scaling coincides with local scaling, , and the fluctuating interface is well characterized by only two independent exponents. However, following the method shown in Ref. RA00 (), we must use at least three independent exponents in a more general case. In our study we use time variation of the local width for several window lengths: , with

(16) |

obtaining direct measures of the exponents , and . Indirect measures of and can be obtained through a collapse of the individual figures as , where now the scaling function goes as for and for . Finally, although the collapse gives an indirect measure of exponents in a very robust way it is worth to have a direct measure of at least three exponents. Thus we also use the evolution of the power spectral density that scales as , with

(17) |

Depicting at different times , we have the direct measure of the spectral exponent that coincides either with the global roughness exponent , when power spectrum do not shift in time (self-affine or superrough scalings) or with the local one , when a temporal shift is observed (intrinsic anomalous scaling).

## Iv Low capillary contrast

In this section we deal with the case of low capillary contrast in both forced-flow and spontaneous imbibition cases. It is worth to mention here that the numerical parameters used in all numerical results have been , (dimensionless units), and . Equation (1) has been integrated over a system of lateral size using a spatial grid of and a time step of . The minimum length of the track disorder has been taken as .

### iv.1 Forced-flow imbibition

We start to study a regime of low capillary pressures. The used value for the capillary contrast has been of and the mean velocity has been fixed to . The shape of the interface for a given realization is shown in Fig. 2a, case . We can see a smooth interface which is slightly correlated with the disorder. The results of the roughness analysis are shown in Fig. 3. From the local width computed at different window sizes we can obtain the global and local growth exponents, and . The best collapse of theses curves (shown at the inset of Fig. 3a is obtained tuning the values of and . The slopes of the scaling function agree with the previously calculated exponents suggesting also a value for the local roughness exponent of which is corroborated at the power spectrum calculated at different times. Since there is not any temporal shift between the lines of the power spectrum, we can assume that the interface fluctuations are described within the superrough anomalous scaling, and therefore . In addition, the spectral roughness exponent obtained from the power spectrum corresponds to the global roughness exponent, which is in agreement with the value obtained previously. These measured exponents are compatible with those obtained by the linear Eq. (8) with a constant velocity . Rescaling such a linear equation by the transformation , we have trivially that and , and assuming superroughness (), we get the remainder exponents as , , in agreement with the measured exponents. Therefore, we can conclude that this regime is well modelled by the non local and linear equation (8), taking into account only the surface tension regime.

If we increase the mean velocity of the interface, then the static crossover length, Eq. (9), decreases, obtaining that interface fluctuations saturate earlier, at the time when the correlation length reaches the crossover length, AL04 (); LA05 (); intrinsic-superrough ().

### iv.2 Spontaneous imbibition

In spontaneous imbibition, the crossover length scale, Eq. (9), becomes a dynamical scale and different regimes can be observed depending on the velocity of the interface intrinsic-superrough (). For low velocities, the initial correlation length is below the crossover length meaning that the relevant mechanism to damping the interface fluctuations is the surface tension with the characteristic dynamical exponent of . On the other hand, for higher velocities, the crossover length acts as an effective correlation length of the interface fluctuations, giving rise to the genuine exponent of . In order to study both regimes, we have controlled the initial velocity of the interface by choosing the initial height of the interface conveniently. We impose an initial height of to study a low velocity regime, whereas a higher velocity regime will be achieved by simply putting .

#### iv.2.1 Low velocity regime

Typical shapes of interfaces in spontaneous imbibition are depicted in Fig. 2b. In the low capillary regime, case of , the interface is weakly correlated with disorder. A roughness analysis shown in Fig. 4 gives the exponents and from a direct measure of the growing local width. The best data collapse on these figures provides the exponents and . Likewise, the slope of the scaling function gives , which is corroborated by the power spectrum shown in Fig. 4b , where we get a spectral roughness exponent of , without temporal shift between the curves, indicating we are dealing with a superrough scaling. Hence, the measured exponents are the same than those obtained in the forced-flow case. It was actually expected, since the relevant terms of Eq. (8) at low velocities are the same in both cases.

#### iv.2.2 High velocity regime

When the initial interface velocity is increased, the velocity-dependent term of Eq. (8) starts to be relevant and the new regime adopts the dynamics of the crossover length , getting then the expected dynamical exponent DU99 (), keeping the spatial structure with the same roughness exponent as before, . The numerical results are presented in Fig. 5, obtaining , and from direct measures, and , , from the data collapse of the local widths, which are also in agreement with the linear description of Eq. (8).

## V High capillary contrast

When the capillary contrast is increased the shapes of interfaces become sharper in both cases of imbibition, as it can be seen in Fig. 2. The effect of disorder is so strong that a kind of columnar geometry is also imposed on the interface. As we will see, two important points of the dynamics of interface fluctuations arise as effect of increasing the capillary contrast. First, interface motion seems to be described by local effects instead of the non local effects observed in the case of low capillary contrast. Second, scaling type changes from superrough to anomalous intrinsic. As in the case of low contrast, forced-flow and spontaneous imbibition cases present distinct patterns of fluctuations and they need to be studied separately.

### v.1 Forced-flow imbibition

We start to study the case of a capillary contrast (see Fig. 2a). In Fig. 6a there are plotted the local interface velocity at two different points x of the system (solid lines) and the mean velocity (dashed line). The curve above the mean velocity corresponds to a channel with a high disorder value, that is , being the disorder of the channel and the mean disorder of the whole system. The curve below the mean velocity corresponds to the next channel with a lower disorder value, that is . We can see that both profiles can be locally described by Eq. (15), , where and mean the channel with and respectively. Since the averaged velocity of the interface is low enough, the velocity can be initially negative SO02b ().

As in the previous section we are interested in calculating the different scaling exponents. The scaling analysis is shown in Fig. 7. From direct measure of the interfacial local width slopes we get , and . The best data collapse is obtained using and . The power spectrum evolution changes with respect to the low contrast case since now there is a temporal shift between the curves, indicating the presence of inherent anomalous scaling. Hence, the spectral roughness exponent must be interpreted as the local roughness exponent, which is in agreement with the slopes of the collapsed scaling function. Note that these exponents are very close to the experimental exponents reported in Ref. SO02b (), , , . As shown in this reference, the shape of the interfaces and the analysis of fluctuations are well reproduced by Eq. (15), which takes into account strong diffusive coupled channels. Therefore, we can conclude that simple numerical integrations of the phase-field model reproduce both the shape and the scaling analysis of interfaces.

#### v.1.1 From low to high capillary contrast

When we increase even more the values of the capillary contrast, interfaces become more correlated to the columnar disorder presenting quantitative changes on the scaling behavior. For instance, in Fig. 8a there is plotted the local width computed in a small window of length . It allows us to calculate the local growth exponent for different capillary contrasts. The complete set of scaling exponents is presented in Table (1). We can see that for high capillary contrasts (), interface fluctuations are always described by intrinsic anomalous scaling. In addition, the dynamics of the correlation is subdiffusive, ranging from to in the highest contrast value. We must interprete this extreme case as having fluctuations completely decoupled. It is characterized by the same local and global growth exponent , which implies a dynamic exponent of , and it can be understood as a regime where the correlation length does not grow in time anymore. Experimentally, it has also been observed in the work carried out in Ref. SO03 (). It is worth noting that there are two exponents and that remain constant.

Alternatively to the phase-field model results, this phenomenon can be reproduced by Eq. (15) taking a random diffusion coefficient LO95 (). By choosing a probability density , being a cut off function, Eq. (15) predicts a growth exponent independently of the value, a roughness exponent , and a dynamic exponent , which is in accordance with the numerical values obtained using the phase-field model (see Table 1). In this sense, we can say that the diffusion coefficient of Eq. (15) depends on the capillary contrast of the system. For very high capillary contrasts () the coupling coefficient can be taken as , obtaining then the decoupled state observed numerically. On the other hand, for lower capillary contrasts () the variations of occur at scales larger than the correlation length and can be taken as constant, obtaining then the regime described by .

Scaling class | ||||||||
---|---|---|---|---|---|---|---|---|

Superroughness | Non-local description | |||||||

Intrinsic A. | ||||||||

Intrinsic A. | Local | |||||||

Intrinsic A. | Description | |||||||

Decoupled state |

### v.2 Spontaneous Imbibition

In Fig. 2b there is plotted four interface profiles evaluated at the same time but different capillary contrast. The velocity profiles of the interface are shown in Fig. 6b. Since now there is not any imposed velocity, the velocity of the interface follows the Washburn’s law in each channel of different noise value.

When we increase the parameter a transition to a decoupled state (, ) also appears in the case of spontaneous imbibition. However, there is an important difference with the forced-flow case. As before, we calculate the local width in a small window size for different capillary contrasts (Fig. 8b). We observe that now, the local growth exponent changes suddenly to , indicating that the interface advances completely decoupled. The transition to the decoupled state can be now discontinuous in time for each value of the capillary contrast. It seems that there exists a length above which the fluctuations become decoupled. Therefore, for a given value of the capillary contrast, the clusters of interface with a size will become decoupled between each other at the time , when correlation length reaches the length . Above , the local description is not valid anymore.

## Vi Conclusions

By means of numerical integrations of a phase-field model we find that there are strong differences between the dynamics of fluctuations in the cases of low and high capillary contrast with columnar disorder. Although these differences also exist in a case of quenched noise composed of squares intrinsic-superrough (), they are more dramatic when the quenched noise is of columnar type. It turns out that the persistence of the noise forces the interface to adopt the same geometry, changing the nature of the interface motion.

In a low capillary contrast case, interfaces are superrough with a dynamics dominated by surface tension with exponents when the velocity is nearly constant or when the velocity varies with Washburn’s law. Furthermore, as interfaces are smooth and can be linearized around their mean value, a simple non local model for the interfacial evolution can be used to explain the observed dynamical scaling of fluctuations. On the other hand, when the capillary contrast is increased, interfaces are sharper and the correlation with disorder is more evident. The observed dynamical scaling corresponds then to an anomalous scaling description with a clear temporal shift at the power spectrum, and a subdiffusive behavior with dynamical exponents ranging from to , depending on the strength of capillary forces. A prominent point to remark is that this behavior can be explained by a local model made of coupled channels with a fluctuating force following Darcy’s law. One can interpret that in the high contrast case, the columnar disorder induces the existence of channels with more or less coupling, eliminating completely the non local character of imbibition in homogeneous geometries. Finally, the difference between forced-flow and spontaneous imbibition has also been elucidated in the high capillary contrast regime. While in the forced-flow case the interface gets completely decoupled above a critical capillary contrast, in spontaneous imbibition the same decoupled state can be suddenly achieved for a given value of capillary contrast.

## Vii Acknowledgments

We acknowledge financial support from DGI of the Ministerio de Educación y Ciencia (Spain), projects FIS2006-12253-C06-04 and FIS2006-12253-C06-05.

## References

- (1) M. Sahimi, Flow and transport in porous media and fractured rock, John Wiley and Sons, New York (1995).
- (2) M. Alava, M. Dubé, and M. Rost, Adv. Phys. 53, 83 (2004).
- (3) A.-L. Barabási and H.E. Stanley, Fractal Concepts in Surface Growth, Cambridge University Press, Cambridge (1995).
- (4) F. Family and T. Vicsek, J. Phys. A 18, L75 (1985).
- (5) J.J. Ramasco, J.M. López, and M.A. Rodríguez, Phys. Rev. Lett 84, 2199 (2000).
- (6) H. Leschhorn and L.-H. Tang, Phys. Rev. Lett. 70, 2973 (1993).
- (7) J.M. López and J. Schmittbuhl, Phys. Rev. E 57, 6405 (1998).
- (8) J. Maunuksela, M. Myllys, O.-P. Kähkönen, J. Timonen, N. Provatas, M.J. Alava, and T. Ala-Nissila, Phys. Rev. Lett. 79, 1515 (1997).
- (9) J. Soriano, J.J. Ramasco, M.A. Rodríguez, A. Hernández-Machado, and J. Ortín, Phys. Rev. Lett. 89, 026102 (2002).
- (10) J. Soriano, J. Ortín, and A. Hernández-Machado Phys. Rev. E 67, 056308 (2003).
- (11) J. Soriano, A. Mercier, R. Planet, A. Hernández-Machado, M.A. Rodríguez, and J. Ortín, Phys. Rev. Lett. 95, 104501 (2005).
- (12) S. V. Buldyrev, A.-L. Barabási, F. Caserta, S. Havlin, H. E. Stanley, and T. Vicsek, Phys. Rev. A, 45, R8313 (1992).
- (13) V.K. Horváth and H. E. Stanley, Phys. Rev. E, 52, 5166 (1995).
- (14) T. H. Kwon, A. E. Hopkins, and S. E. O’Donnell, Phys. Rev. E, 54, 685 (1996).
- (15) O. Zik, E. Moses, Z. Olami, and I. Webman, Euro. Phys. Lett., 38, 509 (1997).
- (16) A. S. Balankin, A. Bravo-Ortega, and D. Morales, Phylos. Mag. Lett., 80, 503 (2000).
- (17) M.A. Rubio, C.A. Edwards, A. Dougherty, and J.P. Gollub, Phys. Rev. Lett., 63, 1685 (1989).
- (18) V.K. Horváth, F. Family, and T. Vicsek, J. Phys. A, 24, L25 (1991).
- (19) S. He, G.L.M.K.S. Kahanda, and P.-Z. Wong, Phys. Rev. Lett., 69, 3731 (1992).
- (20) J. Soriano, J. Ortín, and A. Hernández-Machado Phys. Rev. E 66, 031603 (2002).
- (21) D. Geromichalos, F. Mugele, and S. Herminghaus, Phys. Rev. Lett., 89, 104503 (2002).
- (22) M. Dubé, M. Rost, K.R. Elder, M. Alava, S. Majaniemi, and T. Ala-Nissila, Phys. Rev. Lett. 83, 1628 (1999).
- (23) K.R. Elder, M. Grant, N. Provatas, and J.M. Kosterlitz, Phys. Rev. E 64, 021604 (2001).
- (24) A. Hernández–Machado, J. Soriano, A.M. Lacasta, M.A. Rodríguez, L. Ramírez–Piscina, and J. Ortín, Europhys. Lett. 55, 194 (2001).
- (25) A. Hernández-Machado, A.M. Lacasta, E. Mayoral, and E. C. Poiré, Phys. Rev. E, 68, 046310 (2003).
- (26) T. Laurila, C. Tong, I. Huopaniemi, S. Majaniemi, and T. Ala-Nissila, Eur. Phys. J. B 46, 553 (2005).
- (27) G. Barton Elements of Green’s functions and propagation: potentials, diffusion and waves, Oxford University Press, Oxford (1989).
- (28) E. Pauné and J. Casademunt, Phys. Rev. Lett. 90, 144504 (2003).
- (29) M. Pradas, J. M. López, and A. Hernández–Machado, Phys. Rev. E 76, 010102(R) (2007).
- (30) M. Pradas and A. Hernandez-Machado, Phys. Rev. E74 041608 (2006).
- (31) J.M. López and M.A. Rodríguez, Phys. Rev. E 52, 6442 (1995).