Description of diffusive and propagative behavior on fractals

Grup de Fı ́sica Estadı ́stica, Departmento de Fı ́sica, Universitat Auto `noma de Barcelona, E-08193 Bellaterrra, Spain Departamento de Medicina, Universitat Internacional de Catalunya, c./Gomera s/n, 08190-Sant Cugat del Valle ́s, Barc lona, Spain Departamento de Fı ́sica, Universitat de Girona, Campus de Montilivi, 17071 Girona, Catalonia, Spain ~Received 30 September 2003; revised manuscript received 15 December 2003; published 31 March 2004 !

͗r 2 ͘ϳt 2/d w ͑1͒ exhibited by this kind of spatially correlated structures, as well as many other characteristics of these processes, has been developed by scaling arguments and proved either analytically or numerically before ͑see Ref. ͓1͔ and references therein͒.Specifically, lots of efforts have been done in order to achieve an expression for the probability distribution function ͑PDF͒ P(r,t) for diffusion on fractals ͓2-5͔, i.e., the function giving the probability that a particle or individual stays at time t at a distance r from the origin.
As noted by some authors, the geometry of fractals, on its own, makes local and global diffusion behave in a different way ͓5-7͔.Nonetheless, for most applications we just need to know how the particles move at global scales.For this aim, the homogenization properties of fractals become essential.According to them, the averaged global diffusion behaves as a Brownian motion, although it does not locally ͑this involves important effects, as the ''restoration of isotropy'' analyzed in Ref. ͓7͔͒.It allows us to define an averaged PDF in the asymptotic regime properly ͓8͔.The usual form taken for this averaged PDF is, according to scaling arguments and computer simulations ͓1,9,10͔, where d f is the Hausdorff dimension of the fractal, d w is the dimension of the random walk therein, and c is a constant.
Recently, a great advance was done by Mosco ͓3͔, who considered an alternative approach by introducing an intrinsic metric for fractals, defined by the relation sϵr d w /2 , where s and r are the intrinsic and Euclidean distances, respectively.The importance of this formalism is that the anomalous dynamics of fractals is absorbed therein, so it gives an ideal framework to consider dynamic ideas that in nonfractal systems are developed in the usual Euclidean space.From his arguments, Mosco obtained an expression for the PDF which reads where d min is the scaling exponent of the minimum path possible within the fractal (l min ) between two points of the structure, namely l min ϳr d min .We think that Eq. ͑3͒ is a generalization of the habitual form ͑2͒, which is recovered in the case d min ϭ1 or by averaging in the appropriate regime, as shown previously both theoretically Ref. ͓11͔ and by simulations ͓12͔.Here, we begin from the general ideas presented in ͓3͔ and try to complete them ͑namely, we look for a generalized diffusion equation͒ in order to reach a physical description of diffusion and propagation processes on fractals, as well as for potential application purposes which we describe below.

I. DIFFUSION EQUATION
One of the first attempts to seek an equation able to describe diffusion on fractals was performed by O'Shaughnessy and Procaccia ͓2͔ before expressions ͑2͒ and ͑3͒ had been predicted.Their expression was deduced from the conservation and transport equations 1 with D(r) the diffusion coefficient ͑or conductivity͒ in fractals and J(r,t) the radial flux of particles.Although their arguments seem theoretically correct, they were not able to predict the form ͑2͒ for the PDF.In consequence, later generalizations have been proposed by different authors ͓9,10͔.
In these papers, the general idea is to replace the temporal derivative in Eq. ͑4͒ by a fractional derivative, so this leads to a solution which may reproduce the exponent d w /(d w Ϫ1) in Eq. ͑2͒.However, we note that this approach poses two problems.
͑i͒ The correct parameters are only achieved by means of a normalization, which is different from that used by O'Shaugnessy and Procaccia, and its physical meaning has not been argued anywhere.
͑ii͒ When fractional derivatives are considered into the diffusion equation, a new term on r appears in the solution for the PDF.Until now, this extra term has not been justified by scaling arguments or somehow.
In general, we believe that the physical justification for fractional derivatives is unclear for fractal geometries.We pretend here to show that there is no need for fractional calculus in order to find an equation giving the solution in Eq. ͑3͒, so we overcome both drawbacks ͑i͒ and ͑ii͒ by means of the following method.
First, the procedure from O'Shaugnessy and Procaccia can be directly adapted to the intrinsic metric of the fractal just by writing in their equation s and d S instead of r and where d S is the mass scaling exponent in the intrinsic metric ͑see Table I͒.The main reason that leads us to propose this approach is that undesirable dynamic consequences from fractal nature are absorbed into the intrinsic metric ͓3͔.We argue that these consequences would be the reason why previous methods did not reach the expected results.The approach in Ref. ͓2͔ can account for the fractional geometric properties of fractals, but not totally for their fractional dynamic properties, as the parameter d w did not appear there.
In consequence, the new metric gives an alternative which is more suitable in order to develop the idea of a diffusion process.
We pretend now to find the explicit dependence of D(s).This parameter is defined from Eq. ͑5͒, which is equivalent to Fick's law, but it does not have exactly the same properties as the classical diffusion coefficient, so we will rather refer to it as the conductivity.Likewise, we can define the resistivity of the media as ϭ1/D.This resistivity should be proportional to the number of steps that the particles must move to go through the media; then, arguments from random walks on fractals should be useful in order to find .
Walks on fractals are characterized by the existence of two behavior scales, as claimed in some mathematics works ͓5,13͔ and proved by simulations in different cases ͓12͔.According to this idea of two separate scales, the media is divided into small blocks of size ; within the small blocks the motion is classical, so 2 ϳt, but at scales larger than the effect of heterogeneities appears and motion depends on the fractal parameters ͓14͔.It means that the self-invariance properties of the fractal are not valid at short distances ͑as the ideal concept of self-similarity at all scales does not hold for fractals in practice͒.
From these arguments, it has been shown ͓3͔ that for large scales s, the number of blocks n crossed by walkers grows as nϳ(s 2 /t) ␥ , where ␥ϭd min S /(2Ϫd min S ) and d min

S
is the equivalent to d min in the intrinsic metric ͑Table I͒.In the asymptotic regime, the number of steps done by the walker is proportional to the number n of blocks crossed, so we predict that ϳn.Nonetheless, this result must be renormalized in order to recover constant resistivity for the appropriate ͑homoge-neous͒ case.According to this, we will divide n by the number of blocks n hom corresponding to motion on homogeneous media, as in that case D is the classical diffusion coefficient and is then a constant.Finally we find and hence DϭD 0 (s 2 /t) 1Ϫ␥ for the conductivity.Therefore, we obtain an unexpected temporal dependence on D, a surprising result that has been proposed before, in fact, from experience for specific systems that involve heterogeneous media ͓15,16͔.
It is clear from the arguments above that the temporal dependence on D comes from the existence of two separate scales and the relation between them: random walks ͑diffu-sion͒ take place through the inner scale, but our equations concern asymptotic scales.This is why this dual behavior must be taken into account.Nonetheless, because the idea of a temporal dependence for the conductivity may seem counterintuitive, we have performed numerical simulations for random walks on two-dimensional ͑2D͒ percolation clusters in order to prove Eq. ͑5͒.We chose a circle at a certain radius from the origin of the walk ͑to keep r, and so s, fixed͒ and measured the flux of particles and the spatial derivative of the density at that radius, as a function of time.Introducing the results into Eq.͑5͒ we confirmed that the conductivity D decays with t as a power law ͑Fig.1͒.The exponent of the power law is expected to be for percolation clusters ␥Ϫ1 ϭd min /(d w Ϫd min )Ϫ1ϭϪ0.35Ϯ0.01͓17͔, while the fitted line in Fig. 1 gives Dϳt Ϫ0.33Ϯ0.04 , so the agreement found is good.The little discrepancy and the error bar given above are due to the fact that the spatial derivative is difficult to estimate from discretized simulations exactly.
It is important to note that our simulations involve some limitations.In Fig. 1 the predicted scaling behavior is only found over barely two decades, while one would expect that this universal behavior would extend over many time scales.First of all, it is obvious that we must wait a transitory time till the particles reach the distance r where the observations are done.On the other hand, we must ensure that the asymptotic regime is achieved, i.e., that rt Ϫ1/d w ӷ1 ͓8͔; as we keep r fixed, it means that there is also a maximum threshold in time where our results stop holding ͓18͔.We could certainly choose a higher distance r to obtain a higher upper threshold, but then the transitory time would increase too.All this, in addition to restrictions concerning the running time of the simulations, explain why the power-law behavior predicted only holds in our simulations for a certain time interval.Maybe some experimental work could be useful in order to confirm our predictions for D.
Now, after we have achieved an explicit form for D(s), we can introduce it into Eq.͑6͒ and find the point-source solution for the PDF, by means of the normalization ͐d S s d s Ϫ1 P(s,t)dsϭ1 ͓19͔.This expression is exactly the same as that given by Mosco just by scaling arguments ͓3͔.In view of these results, we propose Eq. ͑6͒ as the general equation governing the asymptotic diffusive behavior in the natural ͑intrinsic͒ metric of the fractal.
This PDF expression has a particular interest since it describes the dynamics of the fractal only in terms of the geometric parameters d S and d min S .The fact that the dynamics of the fractal is determined just by its geometry is an idea not much emphasized before, although it is essential as it agrees with some basic characteristics of the fractal.We stress that the relevance of the geometry on itself only holds for the intrinsic metric, in agreement with the arguments by Mosco; this is the main reason that led him to define the new metric s ͓3͔.
Finally, the corresponding diffusion equation in the Euclidean space, in function of the distance rϵs 2/d w , can be written as where uϭd w d min /(d w Ϫd min ).The exact solution of Eq. ͑9͒ is, as expected, Eq. ͑3͒.It means that we have obtained, a partial differential equation that reproduces the main results of fractal diffusion ͑note that classical diffusion is recovered for the appropriate values d f ϭ1, d w ϭ2, d min ϭ1) and respects the arguments on the intrinsic metric presented by Mosco.

II. PROPAGATION PROPERTIES
Equation ͑9͒ describes diffusion on fractal media, i.e., how the probability density varies temporally according to the spatial behavior of individuals.But we are allowed to include into Eq.͑4͒ a reaction term that accounts for individuals appearing and disappearing at any occupied point of the media.Thus, we can easily generalize our discussion to reaction-diffusion processes on fractal media just by adding a new term, which yields ͑here we choose a logistic term for production, since it is a well-known case and very common in biological applications ͓20-23͔͒.
Reaction-diffusion systems have many different applications ͓24͔ due to their propagation properties ͑their solutions have in general the form of traveling fronts ͓25͔͒.We wonder whether these characteristics can also be found and analyzed for our heterogeneous systems.In order to find the asymptotic speed for the traveling wave fronts in Eq. ͑10͒ we will make use of Hamilton-Jacobi dynamics ͓25͔.The starting point is the hyperbolic scaling procedure t→t/, r →r/ and the representation of the rescaled PDF P (r,t) ϭ P(r/,t/) in WKB form 1. Time dependence of the parameter D(s) in Eq. ͑5͒ for a percolation cluster on a square lattice ͑empty circles͒, which leads to Dϳt Ϫ0.33Ϯ0.04 , and comparison with the simulations for homogeneous diffusion ͑full circles͒, which shows no time-dependence, as expected.
where the action functional G has to be found.It follows from Eq. ͑11͒ that, as long as the function G(r,t) ϭlim →0 G (r,t) is positive, the rescaled field P (r,t)→0 as →0.The boundary of the set where G(r,t)Ͼ0 can be regarded as a reaction front.From these ideas, the Hamilton-Jacobi method yields for the speed of the traveling front ͓26͔ It is interesting to note that the speed of the front does not depend on the geometric dimension d f and the front travels for large times with constant speed if d min ϭ1 and is decelerated if d min Ͼ1; thus, the only parameter determining whether the front is accelerated or not is the fractal dimension of the minimum path d min ͓26͔.All of these effects have been checked by performing numerical simulations of Eq. ͑10͒ ͑see Fig. 2͒ for values appropriate to two structures which are known to agree with the behavior of Eq. ͑2͒ ͓1͔: the Sierpinski gasket (d w ϭ2.32 and d f ϭ1.58) and percolation clusters on 2D near criticality ͓17͔.In both cases, the figure shows that the speed and the acceleration found are in good agreement with Eq. ͑12͒.
We want also to emphasize the fact that, when the usual conditions on fractals d w Ͼ2 and d min у1 are considered, the speed given by Eq. ͑12͒ is always lower than the well-known Fisher's expression for homogeneous media ͓23͔, as other authors have intuitively pointed out before ͓1͔.For the classical geometry ͑namely, d w ϭ2 and d min ϭ1), Eq. ͑13͒ is recovered from Eq. ͑12͒, as it should.

III. DISCUSSION
Our approach has allowed us to reach a diffusion equation in accordance with Eqs.͑1͒ and ͑2͒, which are results widely accepted for fractal dynamics.Previous attempts to do this, in spite of achieving some agreement with simulations ͓4͔, have been rarely supported by physical arguments.When fractional derivatives are introduced into these diffusion equations ͑see, for example, Ref. ͓10͔͒, the system is explicitly forced to show the anomalous behavior and the PDF form in Eq. ͑2͒.On the contrary, the intrinsic metric leads to the expected results more naturally.In fact, the intrinsic metric only transforms the spatial coordinate (sϵr d w /2 ).Then, there is no need to alter the temporal terms in our equation ͑in contrast to what happened in previous approaches͒, since the dynamic properties arise now directly from the geometry of the fractal itself.We consider that this approach not only agrees with the features of fractals, but is far more intuitive.
Besides the theoretical interest of our Eq.͑9͒ on the field of fractal dynamics, we have shown the potential interest of this kind of expressions for application to propagative processes.Specifically, classical reaction-diffusion equations have been widely used to the study of biological expansions and invasions ͓20-22,24͔ since the pioneering work by Fisher ͓23͔.In spite of the good results obtained from them, one of the main objections argued by biologists against such methods is that they assume homogeneous media.Systems characterized by a fractal spatial dimension as those considered here have been pointed out as a possible solution ͓27͔, and specific biological problems where spatial heterogeneities play an essential role have already been raised, such as epidemics on networks ͓28͔ and the coexistence of biological species ͓29͔.For these applications and many others, the expansion rate of the front found here becomes useful to determine the asymptotic dynamical behavior.Thus, our approach may offer a suitable generalization to these kinds of studies, by improving the classical expression by Fisher, as shown by Eq. ͑12͒.Likewise, a more general study about the wave front speeds derived from fractal diffusion equations is performed in ͓26͔.When d min ϭ1 ͑solid line and squares͒ we obtain that v is a constant, while for higher values (d min ϭ1.1; dotted line and circles͒ the front is clearly decelerated, as predicted by Eq. ͑12͒.All the variables plotted are adimensional.For simplicity, 4aD 0 ϭ1 was considered.

FIG. 2 .
FIG. 2. Comparison between the theoretical expression for wavefront speed ͑12͒ ͑lines͒ and that from numerical simulations of Eq. ͑10͒ ͑points͒.͑a͒ Dependence on the parameter u for two wellknown structures at a fixed time tϭ250: Sierpinski gasket (d w ϭ2.32, d f ϭ1.58; solid line and squares͒ and percolation on 2D near criticality (d w ϭ2.88, d f ϭ1.90; dotted line and circles͒.͑b͒ When d min ϭ1 ͑solid line and squares͒ we obtain that v is a constant, while for higher values (d min ϭ1.1; dotted line and circles͒ the front is clearly decelerated, as predicted by Eq. ͑12͒.All the variables plotted are adimensional.For simplicity, 4aD 0 ϭ1 was considered.

TABLE I .
Equivalence between the Euclidean and the new intrinsic metrics and their characteristic parameters ͑from Ref. ͓3͔͒.