CC BY
29
FRACTIONAL DIFFUSION EQUATION FOR AGING AND EQUILIBRATED RANDOM WALKS
V.Yu. Zaburdaev, I.M. Sokolov
We consider continuous time random walks and discuss situations pertinent to aging. These correspond to the case when the initial state of the system is known not at preparation (at t = 0) but at the later instant of time ti > 0 (intermediate-time initial condition). We derive the generalized aging diffusion equation for this case and express it through a single memory kernel. The results obtained are applied to the practically relevant case of the equilibrated random walks. We moreover discuss some subtleties in the setup of the aging subdiffusion problem and show that the behavior of the system depends on what was taken as the intermediate-time initial condition: whether it was coordinate of one particle given by measurement or the whole probability distribution. The two setups lead to different predictions for the evolution of a system. This fact stresses the necessity of a precise definition of aging statistical ensembles.
Keywords: Continuous time random walks, generalized diffusion equation, aging, statistical ensemble.
Introduction
Fractional diffusion equations (FDE) nowadays can be considered as an established and common mathematical tool with applications widely distributed in natural and social sciences, biology, and finance. Usually, one successfully uses them to explain anomalous scaling behavior of some quantities of interest. However, when finer details of anomalous transport processes come to a question, such as aging for example, an investigator, equipped with FDE but with the logic of classical diffusion, often fails to describe it correctly. An intrinsically asymptotic character of FDE is more crucial than in classical diffusion, since it involves nonlocal time/space operators. Therefore, some particular tasks require the knowledge of the dynamics underlying the FDE and not only their final form. In many cases, such dynamics can be modeled by continuous time random walks (CTRW) [1-4]. It allows to follow all derivation steps in detail and obtain useful analytical results. In the present work, we will also use CTRW as a basis model.
The problem of aging is one of the fundamental questions relevant in different fields of physics. The word «aging» applies to all situations when the result of a measurement depends explicitly on the time when this measurement was performed (i.e. what is the time elapsed between the instant when the system was prepared in its original state and the instant when the measurement started). Aging is an intrinsic feature of glasses
[5-7], colloidal systems [8,9], granular materials [10,11], diffusion processes in random environment [12, 13], etc. Measurement setups and physical reasons for such a time-inhomogeneous behavior can be quite different, e.g. the physical system can age due to some slow internal processes or due to external impacts.
One of the mechanisms of aging is described by CTRW with power-law waiting time distribution, where its nature is connected with the overall nonstationarity of the process. Different aspects of aging in continuous time random walks were addressed e.g. in [14-18]. The extensive discussion of properties of aging in trap models and their representation within the CTRW picture is given by Montus and Bouchaud [19]. We state however that the applications of CTRW are not exhausted by a mean-field description of diffusion in systems with traps. The approach can be applied in many other cases like diffusion on a spine of a comb and in some other loopless structures [20], as well as the one generated by chaotic maps (see Ref. [16] and references therein), where the diffusion anomalies are not caused by energetic disorder.
The problem discussed in the present article is related to the one of Ref. [16], where the author considered the aging in the subdiffusion process generated by a deterministic dynamical system. We generalize the results of the Ref. [16] and explicitly express the generalized diffusion equation describing aging CTRW by using a single memory kernel. Further, we consider the situation of random walks with truncated power-law waiting time distributions, which leads to equilibration. We show that the final form of the corresponding diffusion equation depends on whether this waiting time distribution possesses the second moment or not. In the latter case, the convergence of the probability density function (PDF) for the density of particles to a Gaussian is very slow.
We moreover stress that the results obtained in the Ref. [16] and generalized in the present work are pertinent to a very specific variant of the aging problem. To demonstrate this, we discuss a seemingly similar (but in reality very different) approach to aging. Although both approaches are intimately related, they lead to different predictions for the evolution of the system. We show that these approaches also correspond to different experimental situations and discuss their applicability and limitations.
The structure of the paper is as follows. In the Section 1 we formulate the aging problem corresponding to the setup of Ref. [16], introduce the generalized master equation (GME) approach, and derive the corresponding fractional diffusion equation for aging walks. We analyze the asymptotic form of the aging diffusion equation and express it through a single memory kernel. In the Section 2, we apply these results to the equilibrated random walks. We derive the corresponding transport equations with distributed order fractional time derivative and discuss particularly interesting behavior of the moments of the density distribution. In the Section 3, we discuss another approach to aging, the relevant experimental situation, and the corresponding results. The last section is reserved for conclusions.
1. Master equation for aging random walks
1.1. General considerations. Consider a system of particles created at t0 = 0. After the system was created, we let it evolve according to its internal dynamical laws during an interval of time ti > 0. The initial distribution of particles at preparation, n(x, 0) is unknown. At the time t1 we labeled some particles and thus created a known
density n(x,t1) of such marked particles (performed a measurement). For example, we could irradiate our system with light and create excited states of atoms or molecules of the interest inside the region of the order of the light beam radius or its penetration depth. Equivalently, we may consider a single particle that has started at time t = 0 from an unknown location and at the time t1 , was detected at the point x1 and labeled. The quantity of interest is n(x,t2), the density at a given instant of time t2 > t1, which we assume to be experimentally accessible. Then the following question can be put: If we only know the duration of the previous evolution, t1 (called aging time), and the position of the particle at t = t1, what can we say about the future evolution of the particle's position? How precise can we predict the position of this particle at the time t2 > t1, and how does this prediction depend on t1 and t2? In the case of many labeled particles, the same questions apply to the profile at t2 which can be obtained as a convolution of the concentration profile at t1 and the PDF of the single-particle displacements. We note that the density of points in a configuration space of some system can correspond not only to coordinates of some real particles but to whatever other coordinates characterizing e.g. temperature or magnetic field. In our explanation we, however, confine ourselves to the picture of particles.
The problem of aging in this setup can be considered as the «intermediate-time initial condition» problem: The system was created at t = 0 but the initial condition to the corresponding transport equation is posed at a later instant of time t = t1. For normal, Markovian diffusion this does not change the overall form of the transport equation. In non-Markovian cases, especially in the ones with long enough memory, it does, as we proceed to show. This aging problem is schematically illustrated in the Fig. 1. Assume that the density of labeled particles at t = t1 at point x = x1 is known. Starting from this point the particles diffuse, and at t2 = t1 +1 acquire some distribution (not shown). By n(x,t1,t) we denote our theoretical prediction for this distribution. In the case when the whole profile n(x, t = t1) is known, the concentration profile nA(x,t = t2) at t2 is given by a convolution of the former one and n(x, t1,t).
For simplicity we start with the CTRW model on a one-dimensional lattice; generalizations to the continuous case and to higher dimensions are quite evident. Consider a discrete set of sites marked by an index ¿.By m we denote the occupation probability of each site (density of particles). After a certain waiting time at the site i, a particle can jump to the two neighboring sites i — 1, i + 1 with equal probability 1/2. The waiting time distribution at a site is governed by the probability density ^(x). The properties of this function determine the regime of diffusion. If the mean waiting time, (t) = f^ x^(x)dx, is finite, the resulting transport process will be a normal diffusion. However, even in this case some aging effects can still emerge, with the only exception of the exponential
starts at t = 0 from an unknown initial distribution. At t = ti the particles with a profile n(x, t = ti) are labeled. At t2 = ti +1 for the diffusion from a single point we obtain the probability density n(x, ti, t) and for the full profile the density nA(x,t = t2). The dashed line and the grayscale shadow in the (x,t) plane are used to guide the eye for the diffusion from a single point
waiting time distribution [21]. If the mean waiting time is infinite, we are in the situation of anomalous subdiffusive behavior, where the aging effects are the strongest. It is known that in the subdiffusive regime asymptotic transport equations have a form of the fractional diffusion equation, where fractional differentiation appears in the temporal part. It is easy to see that in the case of the fractional time derivative, the solutions of the corresponding equation do not posses the semi-group property. This is a strong indication that the choice of the starting point or an intermediate initial condition may noticeably affect the following dynamics.
To consider this point in detail, we need to derive the transport equation which adequately describes aging. Our approach is based on the generalized master equation and is similar to the phenomenological derivation of the diffusion or Fokker-Planck equations using a combination of a continuity equation and the equation for currents [22]. The continuity assumption contains essentially two balance conditions guaranteeing the probability conservation: a local one (giving the balance between the probability gain and loss at one site) and the one for transitions between the two sites (representing particle conservation during the jumps). A balance equation at each site reads
ni(t)= j+(t) - j-(t), (1)
where j-(t) is the loss current, i.e. the probability for a particle to leave the site per unit time at time t, and j+(t) is the gain current at a site.
A particle arriving to the site i at time t comes either from the left or from the right. Probability conservation for transitions between sites then reads
1 - ^ 1
By combining (1) and (2) we get a continuity equation
j+(t) = 2 j--i(t) + 2 j~+i(t). (2)
ni(t) = 1 j-_1(t) + 2j-+i(t) - j-(t). (3)
According to the waiting time distribution, the loss current at time t is connected to the gain current at the site at all previous times: the particles which leave the site i at the time t (making a step from i to one of its neighbors) either were at the same site i from the very beginning, or came there at some later time 0 < t < t. A probability density to make a step at time t when arriving at t is given by the waiting time distribution ^(t — t). Then for the loss current j- (t) we can write:
j-(t) = V(1)(t, tiKo + / y(t — T)j+(T)dT. (4)
Jo
However, for the aging initial conditions the probability to make the first step after the measurement at t1 is different from ^(t). Here we postulated it to be independent of the spatial position i and denoted by ^(1)(t, t1). We also assumed that the system was created at the time t = 0, and a measurement was performed at t = t1. From here on we denote by t the time elapsed from t1 and denote by nii0 = ni(t1) the initial condition for the system's further evolution. By using Eqs. (2)-(3) we can rewrite (4) as:
j-(t)= ^(1)(t,t1 Ko +/ y(t — t) [fii(t)+ j-(t)] dx. (5)
o
Now we have to find the expression for the forward waiting time distribution of the first step t, t1) [23]. The forward waiting time t is counted starting from the observation point. Let us assume that the jump preceding t1 (numbered j — 1), took place at time tj-1 = z. The forward waiting time distribution V^t, t1) is
rt i
^(1)(T,ti)= f(z)ty(ti - z + T)dz, Jo
where f (z)dz is the probability to make a jump within the time interval between z and z + dz, so that f (z) is the time-dependent density of steps. This can be presented in the following form:
f (Z)=£ fn(z).
(6)
n=0
Here fn (t) is the probability density that it is exactly n-th jump that takes place at time t. This one is given by an n-fold convolution of the waiting time probability density ty(t) with itself. Under Laplace transform with respect to z Eq. (6) reads:
1
fs = 1+ ^s + ^ + ... =
1 - ^s
(7)
By indexes p, s we will denote the Laplace components corresponding to t and t1, and by k the Fourier component, respectively. The Laplace transform of V^t, t1) with respect to t1 is then:
^(t) =
^s -JO e-sT' ^(x')dx'
1 - ^s
The double Laplace transform of this function with respect to both times t and t1 has the following form [16,23]:
,(1) = ^ —
s
(8)
^ (p - S)(1 - ^s)' With this information at hand, we take the Laplace transform of (5) with respect to t to obtain:
j-p = ^ (ti)ni,0 + ^p
p nip - Ui,0 + j
i,p
From the above equation we find the connection between loss current and the occupation probability:
p ^p . ^(h) — ^p
ji,p = 1 -' ni,p +
-ni,0.
(9)
Fp 1 — ^p
Now we insert (9) into the Laplace transform of (3) and convert it back to the time domain:
ni(t) = Î 0(t - T) 0
1 n-i(T) + 1 ni+i(T) - ni(T)
dT +
+ ®(1)(t,ti)
11
2 ni-1,0 + 2 ni+1,0 - ni,0
where O(t) and <S>(1)(t,t1) are defined through their Laplace transforms:
p 1 - V p V 7 1 - ^p
(10)
(11)
e
1.2. The second memory term. We note now that Eq. (10) can be rewritten in the form
^ = dît jf M (t -x) x
d rt /0
1 Ui-i(x) + 1 Ui+i(x) — Ui{x)
(12)
dx +
+ $(l)(Mi)
1 Ui-l(tl) + 1 Ui+l(tl) - Ui(tl)
where the memory kernel M(t) is the inverse Laplace transform of the function Mp = (1 — which is connected with the density of steps f (t) (having the Laplace representation fp = 1/(1 — ^p)) via M(t) = f (t) — S(t). Thus, M(t) is the density of all steps excluding the first one. Evidently, the asymptotic properties of the functions M and f are the same.
Now we concentrate on the temporal asymptotic behavior of the source term. By using Eqs. (8) and (11), for the double Laplace transform of (t, t1) we can write:
p,s
1
p — s 1
p — s
1
1
%p
1 — %s 1 — %p J 1 1 J
1 — %s 1 — %p_
s(1 — %p) 1
s(1 — %p) ' s
(13)
1
+ - •
In order to understand the structure of the above expression we use one remarkable property of the Laplace transform. For any function h, such that its Laplace transform exists, the following property for double Laplace transform holds [21]:
hh [h(t + tl)]p,s = ^^
s — p
(14)
where s and p are Laplace space variables corresponding to t and t1 respectively. By using this fact we can rewrite the expression for $(1)(t,ti) (13) in time domain:
$(l)(t,tl)= f (t + tl) — f (t)+ à(t),
(15)
where f (t) is the density of steps given by the inverse Laplace transform of 1/(1 — This one, as we have already seen, is connected with the first memory kernel M(t) via M(t) = f (t) — S(t), so that it can be rewritten in the form (t, t1) = M(t + t1) + +S(t + t1) — M(t). However, for whatever aged system (t1 > 0) the delta function vanishes, so that the final result
$(l)(t,tl) = M (t + tl ) — M (t)
(16)
follows. This is a new expression which connects the source term with the memory kernel of the generalized master equation for aging continuous time random walks (12), which thus gets the form
dtUi(t) = dt lM(t—T)
+ [M (t + tl) — M (t)]
1 Ui-l(x) + 1 Ui+l(x) — Ui (t)
dx +
1 Ui-l(tl) + 1 Ui+l (tl) — Ui(tl)
X
In the more general continuous case, when the jumps' lengths are no more discrete but rather distributed according to some probability density g(x), Eq. (17) turns to an integral equation of the form:
d d f*
dtn(x,t) = dt J0 M(t - x) J g(y)[n(x - y, x) - n(x, x)\ dydx + (18)
g(y) [n(x - y,ti) - n(x,ti)] dy.
The discrete case, Eq. (17), corresponds to the choice g(x) = [S(x - a) + S(x + a)]/2 with a being the lattice spacing. The generalized master equations, Eq. (17) and Eq. (18), contain the standard random walk part represented by the first term and an additional term representing the memory on the initial conditions. This memory term is expressed through the same CTRW memory kernel as encountered in the first term on the right hand side and vanishes when either t\ = 0, or ^(x) is an exponential. Before considering the asymptotic form of the above Eqs. (17), (18) we give their solution in Fourier-Laplace space:
+ ti) - M(t)]p (gk - 1) + i} n(ti)k
nk'p = p (1 - (gk - 1)MP) . (19)
Indexes k and p denote the Fourier and Laplace transforms of the corresponding functions.
1.3. Asymptotic form of the aging equation. Now we would like to obtain the asymptotic form of the above transport equation corresponding to large values of x and t (small k and p respectively). Assuming ni(t) to change slowly enough as the function of the site number i one can change to the continuous description introducing the density n(x,t) with x = ia. The difference operators in Eq. (17) can then be considered as a discrete approximation to a Laplacian, so that the corresponding master equations takes the form of generalized diffusion equation with the additional memory term
d , a2 d /■*„„., w ,
d^n(x,t) = — dt M(t - x)An(x, x)dx +
a2
+ [M(t + ti) - M(t)]2An(x,ti). (20)
The same form follows from the more general continuous form, Eq. (18), provided the jump length distribution g(x) has a finite second moment ((x2) < to). The corresponding proof follows either the standard Kramers-Moyal procedure or can easily be obtained as the small-k expansion in the Fourier space. Note that the continuous limit is obtained not by tending the lattice spacing a to zero, but using the smallness of the ratio between typical observation scales and the microscopic scale given by a: x/a » 1. One usually takes the first leading term in the expansion with respect to this small parameter, finds the solution, and then estimates the remaining correction terms. A more detailed discussion of the transition to the continuous limit in the anomalous diffusion problems can be found in Ref. [24].
For exponential waiting time distribution ty(t) = t-1 exp(-t/t0) the memory kernel M (t) = 1/t0 and the equation takes form of the ordinary diffusion equation
d a2
dtn(x,t) = 2t0 An(x,t) (21)
lacking the second memory term. The combination K = a2/2t0 is the usual diffusion coefficient of the process. Another common choice of the waiting time distribution (determined by both, real practical situations and mathematical convenience) is a power law function, e.g. ^(t) = Y/to(1 + t/to)1+Y, with 0 < y < 1, where to again gives us the characteristic temporal scale of waiting times. The expansion of its Laplace transform for a small p is: = 1 — r(1 — Y)t0PY + O(p), and the corresponding memory kernel M(t) is given by
1 tY_1
M (t) ~ —---,-(22)
U r(Y)t0 r(1 — y)' ^ j
in which we recognize an expression proportional to the integral kernel of the fractional derivative. After substituting this expansion into (20) and some algebra we obtain:
d
Trn(x,t) = Ky oDl-YAn(x,t) + (23)
dt
+ K 1
r(1 — y)
1
(t + ti)1-Y t1-Y
An(x, t1),
where the generalized diffusion coefficient KY stands for KY = a2/[2r(Y)t0]. The Eq. (23) generalizes the aging subdiffusion equation derived in Ref. [16] (cf. Eq. (18) there). It is now applicable to any initial density of particles n(x,t1) whereas the setup of Ref. [16] assumes a localized initial condition as a delta peak at zero.
Let us briefly discuss the properties of the above equation. It is easy to check that the total number of particles is conserved by setting k = 0 in Eq. (19) since nk=0(t) = n(x,t)dx. Just from the normalization condition it follows that
(gk — 1)lk=0 = 0 and we immediately obtain:
n(t1)k=0
np,k=0 = -,
p
which after the inverse Laplace transform demonstrates the conservation of the total number of particles. The behavior of the mean square displacement, M2(t) = (x2(t, t1)), can also be computed with a help of the answer in Fourier-Laplace space (19) by the following formula: M2(t) = — (d2/dk2)nk(t)|k=0. By taking the initial condition as a delta function at zero n(x, t1) = S(x) and for large t we obtain (cf. [16,18]):
a2
<x2<"1» = t0r(1 + Y)r(1 — y) !(' + ^ — ^ (24)
By setting t1 = 0 we recover the standard answer for the subdiffusion. For t1 > 0 we can see that the diffusion gets slower for larger values of t1 . This is an intrinsic feature of the nonstationary subdiffusion process since the frequency of jumps decays with time. After a longer preceding evolution t1 it takes longer time to make the first jump after the observation started. Only at times t ^ t1 the mean square displacement approaches the standard behavior tY. For the whole class of processes described by Eq. (20) the nonnegativity of the probability density for the power-law kernels considered in guaranteed due to the fact that the process discussed is subordinated to the Brownian motion. In the case of the localized initial condition can be expressed in the integral form through the Levy stable distribution [16].
1
2. Equilibrated random walks
Let us now consider the system which relaxes to a true equilibrium. Such a system would correspond to CTRW with the waiting time distribution possessing the first moment, (x) = fj x^(x)dx, which is however so large that the intermediate power-law asymptotics of the function is still seen. The examples are a theta-truncated power-law ^(x) = A/(1 + x)1+a6(T — t) or an exponentially truncated one ^(x) = A/(1 + x)1+a exp(—t/T) (A is the normalization constant), or an exponentially truncated one-sided Levy-distribution with the Laplace transform = exp[—A(p + 1/X)a + A/Xa]. Here A is the parameter with the dimension [T]a and of the absolute value of unity (the typical step time is set to one), and X = (Aa/x)1/'1-"', where (x) is the mean waiting time. The equilibrated case is interesting because it is experimentally relevant: systems are typically created much earlier than measurements are performed and have enough time to equilibrate. The time lag between the labeling particles at t = t1 and final measurement at t2 is however small enough to probe nonequilibrium dynamics. The general result for the aging problems can be easily applied to this concrete example. The equilibrated case corresponds to the limit t1 — o of aging CTRW. Analogously to the aging problem, we illustrate the equilibrated problem setup in the Fig. 2.
The system prepared at t = 0 evolves during the time t1 , with t1 - o (in
practice it means that the system was created a long time before the measurement took place and develops a uniform distribution of particles). At the time t = t1 we select some profile of the particles n(x,t = t1) and follow only these particles until the time
t2 = t1 + t.
We start from the transport equation (20) by taking its Fourier-Laplace transform:
Fig. 2. Problem setup for the equilibrated random walks. At t = ti we label some of the particles from their uniform distribution and create the profile n(x,t = ti). Later we follow only the diffusion of those marked particles and obtain at t2 = ti + t the density profile nA(x, t = t2)
pnk,p — nk (t1) = —
a_k2nk _ a_k20(D(t1) 1 — 2 knk'p 2 k (t1)'
(25)
where we returned to our earlier notation for the second memory kernel In the above equation, the waiting time distribution has a finite first moment, (x), and therefore = 1 — (x) p + o(p). In addition, the limit t1 — o should be applied. In this case the expression for the second memory kernel $1(t) is simpler than before. Its Laplace transform is given by (see Eqs. (8) and (11)):
$(1) p
p (x) 1 — ^p
+ 1
(26)
We note that in the time domain we now have $(1)(t) = 1/ (x) — f (t) + S(t) and that f (t) still has a physical meaning of the (time-dependent) density of steps. For a pure power-law
1
1
CTRW this function monotonously decays to zero with time, while for equilibrated walks it decays to 1/ (x). The transition between the both regimes takes place at at time tc ~ (x). However, this does not mean that the second memory function vanishes at times which exceed tc: the difference 1/ (x) - M(t) still can behave as a pure power-law in some cases. It is important that the second memory function is governed by the subleading term in the expansion of 1/(1 - ^p), whose behavior depends on whether the function ^(x) possesses the second moment (we shall call this situation «a sharp cut-off») or does not.
In the case of a sharp cut-off and for a small p, = 1 - (x) p + M2p2 + o(p2), where M2 is the second moment of waiting times. Now we have:
11
1 - ^p {t) P + M2p2 + o(p2)
= 71P[1 + (M2/ {T))p + o(p)] {T) p
IT- + T2 + o(l). (27)
{t)p {t)2
This means that for small p
^ = 1 - M + o(T), (28)
{t)
so that the function $(1)(t) is integrable, with the integral being equal to 1 — M2/ {t)2. This in turn means that <^1(t) decays faster than t-1 and its exact decay form depends on how many moments does ^(t) actually possess. We note that the existence of the second moment is necessary for such a behavior. The overall fractional diffusion equation for such a process in the long-time limit tends to a normal diffusion equation of the type
dn(xt = k An(x, t) + Cb(t)Año, (29)
where K stands for the diffusion coefficient K = a2/ {t) and the nonnegative constant C is given by C = (1 — M2/ {T)2)a2, i.e. to the ordinary diffusion equation with some correction to the initial condition. This correction vanishes only for the Markovian random walk process with exponential waiting time distribution for which M2 = {t)2.
In the case when the second moment does not exist, the situation can be vastly different. In order to present this type of behavior let us consider a special example of the waiting time distribution for which we can convert Eq. (25) in the Fourier-Laplace domain into a distributed-order fractional diffusion equation with a source in space and time. Let us take in a form
^p = exp
p {t)
1 + A(p {t))
1-a
(30)
This function is completely monotonic, which, according to the Bernstein's theorem, means that it is a Laplace transform of a probability density function. To prove this it is enough to note that has a form e-hp with the function hp > 0, hp\p=0 = 0 and possesses a monotonous derivative (see [25]). For p 0, behaves as ~ 1 - p (x), i.e. the corresponding PDF ^(x) has the mean value of (x). The parameter A is a free parameter of the order of unity which governs the precise form of the crossover.
Provided that p is small enough, we can take œ 1-p (t) /[1+A(p (t))1-"] which corresponds to 1/(1 — ^p) = 1/(p (t)) + A(p (r))-a. By substituting this expression into (25) we obtain:
pnk,p - no,k — -K [1+ A(p (t))1 a] k2nk,p + + K (t) A(p (T))-ak2no,k.
(31)
In normal space-time domain this equation corresponds to
dn(x, t) df~
=K
1 + A (t)
1—a ni—a
0Dt
An(x, t) —
- K ~r7~—~ t—i+aArio-r(a)
(32)
Note the distributed order derivative (1 + Ax1_a 0Dj-a) on the r.h.s. [25], which is typical for systems showing crossover behavior, and the extremely slowly dissolving of the initial condition!
The difference between the cases of the equilibrated random walks corresponding to the waiting time distributions with or without second moment (mirrored in the difference of the form of the corresponding generalized Fokker-Planck equation) is intimately connected with the question whether the mean forward waiting time (the mean waiting time of the first step after the beginning of observations) exists. Indeed, looking at the limit of the Laplace transform in the first argument of the forward waiting time distribution for the equilibrated case (t1 — to) we see that ^P1^ — (1 - )/ (x) p. Therefore, the
first moment of ^(1)(t), (ti) = L°° ^(1)(t, oc)tdt = - dp ^P0 is equal to M2/2 (x)2
p p=0
provided M2, the second moment of the waiting time distribution, exists, and diverges otherwise.
One of the interesting peculiarities of the equilibrated random walks is found in the behavior of its moments. It is known that the second moment of equilibrated CTRW behaves diffusively at all times [26,27]. However, the higher moments of the density profile are much more exotic when the second moment of ^(x) does not exist. From Eq. (25) we obtain:
nk,p —
1 — ^P
a2k
p2 (t) 1 — % + %a2k2/2
1 _ 1 — %p a2k2/2 p p2 (t) 1 — %p + a2k2
1
p
1
1 a2k2
p (t) 2
+(—1)
2n
1
+ ... +
2n 2n
a2nk
p (t) (1 — %p)n—i
+
(33)
Now we easily find the even moments of the density by the following formula:
M2n = (—1)ra(d2ra/dk2ra)nk|k=n. We thus have:
Mo p _ 1/p
Mo(t) = 1
M2 p _ 2Kp
—2
M2(t) _ 2Kt
M _ 24K2(t)
M4 'p _ P2(—p)
M4(t) _24K2 (t) /4 (t)
(34)
M2np _ Pnff^Cl1 M2n(t) _ (2n)!Kn (x)n-1 /2n(t).
All odd moments are zero due to the the symmetry of the case considered here. The functions f2n(t) are the inverse Laplace transforms of f2n,p = 1/p2(1 — ^p)n-1. For the waiting time distribution with infinite second moment considered above, its asymptotic Laplace transform is as follows:
/2n P _
1 + A(p (t)) pn+1 (T)
1-a
n- 1
(35)
In the time domain it corresponds to the following behavior:
An
/2n(t) ^
n— 1 (T)-a(n-1)
r[(n - 1)a + 2] 1 tn
n!t
n—1
t(n— 1)a+1 for t < (t)
for t > (t)
(36)
This means, for example, that the distribution of the particles' positions at an intermediate time is considerably platykurtotic, and only slowly tends to a Gaussian in course of time with a typical transition time (x). We note here that such a slow convergence to a Gaussian is typical also for other random walks models with truncated power-laws, as exemplified by truncated Levy flights [28,29].
3. Green's function approach and the problem of aging ensembles
In the previous sections, we considered the results for a specific setup of an aging problem corresponding especially to a displacement of a single particle labeled at t = t1 in course of its further temporal evolution or of the ensemble of such labeled particles. Let us consider another setup for aging problem. Assume that the full profile of the density of particles (or the full probability density for a single particle) at t = t1 is known (for example, we started with a thermodynamically large ensemble of particles and were able to measure their density at t = t1 with high enough precision) or is known in principle as it is assumed when calculating, say, the correlation functions [30]. Does it give us additional information which might improve the prediction of the particles' positions? What would be the following evolution of this profile and how will it depend on the aging time?
The GME approach used in Sec. 1. relies on a special definition of the aging ensemble. The system was created at t = 0 and evolves according to the CTRW dynamics
up to the time t1, where the first measurement of the particle's position (or the the particles' positions) takes place. No explicit notion about the initial distribution of the particles' positions at t = 0 is assumed or used. In any case, we follow up only the particles tagged at t = t1. The situation considered in the present section corresponds to a different ensemble.
Imagine now that we have two experiments (probably with different initial conditions), which are such that the exact measured density in one of them is the same as the density of labeled particles in the another one. Will the further development of the densities be the same? If not, what is the difference? As we proceed to show, the two situations lead to different results for the higher moments of the distribution. It is worth mentioning that the distribution of the first exit times after the measuring event is the only quantity affected by aging and it determines the future evolution of the system for both definitions [16]. It depends on the initial conditions through the correlation between the particle's position and the forward waiting time, as discussed in Ref. [31].
The difference between the two setups, the one of Fig. 1 and the one of Fig. 3, on the qualitative level, can be understood already now. In the given point x1 at the time t = t1, there are particles which have arrived different times ago and, therefore, will make their next steps also at different times [32, 33]. In the first aging setup, we virtually ignore the dependence of these exit times on x1 (the relative position in the full density profile, if the latter is known). Therefore, we assume that in all points the distribution of exit times depends only on t1 . In the second aging problem, we take into consideration the relative position of the point x1. It is clear that particles which occur to be found at the wing of the density profile typically have made much more steps than those found close to the center of the profile. Moreover, the particles which had to make a lot of jumps probably have arrived at their actual position very recently, i.e. start their last waiting period not long ago before t1. This implies that the distribution of exit times depends on x1 which affects the diffusion of particles out of this point. Therefore corresponding profiles in the Fig. 1 and 3, n(x,t1,t) and n(x,x1,t1,t) are different. By repeating the same argumentation for each point in the profile n(x,t = t1) we may conclude that nA and nG will be different as well. We will come back to this question later in the text.
Let us now return to the case when the PDF of particles' positions at t = t1 is known exactly and note that the temporal evolution of the PDF from true initial conditions to the state at time t is described by the linear evolution operator S(t). If the time evolution operator has an inverse defined at least on a set of relevant PDFs n(x, t), one can propagate this PDF back in time to t = 0, thus getting n(x, 0) = nn(x) and then forward in time up
the full density profile, n(x,t = ti) is measured. By using the Green's function it is possible to propagate it back and reproduce the initial distribution at t = 0, n(x,t = 0). Then it is forwarded until the time t2 = ti + t to obtain the prediction for the future evolution nG(x,t = t2). In contrast to the first aging problem, diffusion from a single point xi will depend on both, the starting point xi and the aging time ti : n(x, xi,ti,t). The grayscale shadows and dashed line in the (x,t) plane show the diffusion from the initial distribution and a single point.
to t2 > ti to obtain n(x, t2)*:
n(x, t2) = S(t2)S-l(ti)n(x, ti). (37)
This defines the linear time evolution operator S(t2; ti) = S(t2)S-i(ti) which gives us the PDF at the time t2 provided the PDF at the time ti is known.
Let us first find the explicit form of the S(t2; ti) operator and then discuss the difference between this aging problem, and the first one, considered above. The time evolution operator, S(t2; ti), follows from the explicit form of the solution of the transport equation via the Green's functions method: it is the integral operator containing the Green's functions of the standard CTRW equation, G(x,x',t). Let us recall the form of the transport equation for CTRW without aging, i.e. all particles were introduced at t = 0 without history (the first part of Eq. (10) or Eq. (18)):
r t
n(x,t)= <$(t — t) g(y)[n(x — y, t) — n(x, t)] dydx. (38)
Jo J-<x
The Green's function is defined as the solution of the above equation with n(x, t = 0) = = S(x — x'). Then, for arbitrary initial distribution n0(x) we can write:
G(x,x',t)no(x')dx'. (39)
For the homogeneous situation considered here, the Green's function depends only on the difference of its spatial variables: G(x,x',t) = G(x — x',t). Therefore, the solution can be found as a convolution of the initial condition with the Green's function. By denoting the convolution operation with *, we can write S(t) = G(x,t)*, i.e.
n(x, t) = S(t)n0(x) = G(x, t) * n0(x). (40)
The PDFs of the particles' distributions at time ti and t2 = ti +1 are thus given by
n(x,ti +1) = G(x,ti + t) * n0(x),
n(x,ti) = G(x,ti) * n0(x). (41)
Under the Fourier transforms, convolutions are changed into simple products of Fourier components:
nk(ti + t) = Gk (ti + t2)no,k,
nk (ti) = Gk (ti)no,k ■ (42)
*This is essentially the idea used in Ref. [29]. However, the technical implementation of this idea is faulty. This can be seen when considering Eq. (18) of Ref. [29]. We note that taking t2 = ti = t in this equation has to give W(X2,t; Xi,t) = 8(x2 — xi): the particle cannot be at two different points at the same time. On the other hand we note that due to the completeness of the system of eigenfunctions of the corresponding Hermitian operators we have Y1 n ^n(Xi)^n(X2) = S(x2 — xi). Due to the uniqueness of the eigenfunction expansion this expression and Eq. (18) with t2 = ti = t are compatible only in the case
when Ea(z)Ea(—z) = 1, which can be proven wrong for whatever a = 1. A trivial example is given by
2
Ei/2(z) = ez [1 + erf(z)]. Thus the assumption that backward propagator is the solution of the fractional backwards Kolmogorov equation seems to be not a quite transparent choice in all cases except for the trivial Markovian one where Ei(z) = exp(z).
From these two we formally find:
nk(tl + t) = Gk(t) + ^ • nk(il), Gk (t1)
and now write the time evolution operator S>(t2,t1) = S(t2)S as an integral
convolution operator T(x, t2, t1)* with the integral kernel being the inverse Fourier transform of Gk (ti + t)/Gk (ti).
For the CTRW model its Green's function can be easily found analytically [3], so that also the corresponding operator S(t2, t1) can be immediately calculated. The Green's function for a CTRW equation (the solution of (38)) in the Fourier-Laplace space is
Gk ,p =
1 - %
P (1 - %p9k ) '
(44)
We can now compare the resulting density distributions obtained in the two aging problems, namely Eq. (18) and (43). To do this we analyze the difference of the two densities Sn(x, t) = nG — nA, with nG being the PDF predicted by the Green's function approach (43), and nA being the one given by the by the probabilistic approach (18) of Section 1.2. In particular we concentrate on the behavior of moments of the density difference defined as (xn)bn = f^ xnbndx.
As an example we consider the following situation. Imagine that the system evolved starting from a sharp (delta-function) initial condition n(x,t = 0) = = S(x). In this case the particles' distribution at t = ti is given by the inverse Fourier transform of nk(t1) = Gk(t1). This distribution is now used as the intermediate-time initial condition at t = t1 in both approaches, in the one of Sec. 1 and in the one of the present Section. Then the results of these two approaches are compared. Eq. (43) reduces in this case to a
simple formula n£ (t +11) = Gk (t +11).
Fig. 4. The forth moment of 8n(x, t) = n — n , as a function of time, t, for different aging times ti and Y = 1/2. Its deviation from zero shows the difference between the prediction given by Eq. (18) as compared to the Green's function result (43).
Omitting the details of rather tedious calculations, we can show that the second
moment of the density difference,
(x
bn-
is equal to zero. However, already the next
even moment,
bn
deviates from zero. The final expression for the forth moment of
the density difference reads:
(X4) bn
= 24
sin ny f dx
nyr(y)2 / (t - t)1-^ 0
(45)
x —
t^1-y 1
T
t1 + t
[t1 - (t1 - T1)Y] dT1.
In Fig. 4 we plot (x4)bn given by (45) as a function of time for three different aging
times, for the case y = 1/2 (y is the exponent in the power law tail of the waiting
4
x
t
x
time distribution). These findings show that the approach based on the backward-forward propagation delivers a result which is different from the one given by the approach of Sec. 1. It is also necessary to stress, that the aged propagator n(x,t,ti) given by (18) is not the Green's function of aged anomalous diffusion, since it does not reproduce the PDF of the particles exactly. The Gaussian situation showing no aging is the only one when the both approaches are equivalent (and give the same result): in this case Gk(t) a exp(—Dk2t) so that Gk(ti + t)/Gk(ti) = Gk(t).
Let us discuss some implications of the two approaches. The Green's function approach utilizes more information about the intermediate stage of the system (the full knowledge of the whole distribution) and corresponds to a different nonequilibrium statistical ensemble than the one of Sec. 1. On the other hand, the approach of Sec. 1 assumes no knowledge about the initial condition, i.e. essentially starts from the assumption that at t = 0 the distribution of particles in space was uniform (which is the most reasonable assumption about the state of the statistical system provided no additional information is given and corresponds exactly to the Jaynes' information approach to statistical mechanics [34]). It is clear that the homogeneous distribution does not evolve in course of the time and stays uniform up to the time t = ti. At t = ti we choose some profile n(x, t = ti) from this distribution (see the Fig. 2) and follow its evolution.
As we have mentioned above, the physical reason for being able to predict exactly the dynamics of particles after a measurement, is the precise knowledge of their exit times. The Green's function approach uses this information only indirectly, just propagating the profile back and then forward and automatically provides the corresponding distributions. This approach is not literally applicable in a whatever case when the intermediate time density distribution is not a full profile resulting from the previous evolution of some initial density. There exist, however methods which work with microscopic distributions directly and can be applied to such cases as well. The first one is based on the generalized transport equation with microscopic details taken into account [21], another one considers two-point probability distributions [31,35]. Both of these more general approaches (not restricted to the aging problems) give answers corresponding to the Green's function approach when starting from a concentrated initial condition (and therefore for a whatever exactly known one). Moreover, these approaches allow for the full solution of the aging problem in the setup of the Sec. 1, where the exit times after the measurement tend to be independent of coordinate. The results for the equilibrated walks provided by methods of Refs. [21,31] and by that of the Sec. 1 coincide since the information about the initial state of the system is forgotten.
It is important to state that the two ensembles corresponding to the two approaches discussed in this work are not the only two alternatives, but two of a quite broad spectrum of possibilities. They provide however the two limiting situations: the full knowledge and no knowledge about this distribution. There are various situations possible, when the distribution of the particles' positions at the intermediate time is known to some extent or with some uncertainty, in which case the results for aging will depend on what exactly was measured at t = ti and what the precision of this measurement was.
Conclusions
Summarizing our findings, we have derived a general equation for the PDF of the particle's position in aging CTRW, corresponding to the «intermediate-time initial condition» setup, in which the coordinate of the particle is measured at some time ti after
the preparation of the system considered to take place at t = 0. We were able to express the second memory function, describing the dissolving of the initial condition, through the memory kernel of the CTRW equation. We discussed the ensuing forms of the transport equation for the case of CTRW with power-law waiting time distributions lacking the first moment, as well for the equilibrating situations, where the first moment exists. We have shown that the exact form of the corresponding equation for equilibrating walks depends on whether the second moment of the waiting time distribution exists or not. In the first case, the normal diffusion equation appears, in the second case, the evolution is described by the diffusion equation with distributed-order temporal derivative and with a very slowly dissolving second memory function following the power law. The asymptotic behavior of the higher moments of the density profile for this case was also calculated and indicated a slow convergence to the Gaussian behavior.
Moreover, as we tried to show, the problem of aging is a very delicate task even in the framework of exactly solvable model of CTRW. To illustrate this we have considered and solved another possible aging setup corresponding to the full knowledge of the PDF of the particles' positions at time ti after the preparation of the system at t = 0. The two setups give different predictions for the evolution of the particles density and correspond to different experimental realizations. This shows that there is still a room for further developments in the area. One of those is the question of the prediction of the evolution based on an incomplete intermediate-time initial condition. This work is currently in progress.
We would like to thank Profs. J. Klafter and E. Barkai for useful discussions. IMS thankfully acknowledges the financial support by DFG within the SFB555 research project.
References
1. Montroll E.W. and Shlesinger M.F. Nonequilibrium Phenomena II: From Stochastic to Hydrodynamics. In «Studies in Statistical Mechanics». / Eds. J. Leibowitz and E.W. Montroll. North-Holland, Amsterdam, 1984. Vol. 11. 1.
2. HausJ.W. and Kehr K.W. Phys. Rep. 1987. Vol. 150, № 5-6. 263.
3. MetzlerR. and Klafter J. Phys. Rep. 2000. Vol. 339, № 1. 1.
4. BouchaudJ.P. and Georges A. Phys. Rep. 1990. Vol. 195, № 4-5. 127.
5. Struick L.C.E. Physical Aging in Amorphous Polymers and Other Materials. Elsevier, Houston, 1978.
6. Bouchaud J.P., Cugliandolo L.F., Mezard M., and Kurchan J. In «Spin Glasses and Random Fields». Edited by A.P. Young. Singapore: World Scientific, 1997.
7. Dupuis V., Bert F., Bouchaud J.P., Hammann J., Ladieu F., Parker D., and Vincent E. PRAMANA-J. Phys. 2006. Vol. 64, № 6. 1109.
8. Bellour M., Knaebel A., Harden J.L., Lequeux F. and Munch J.P. Phys. Rev. E. 2003. Vol. 67, № 3. 031405.
9. Abou B., Bonn D, and Meunier J.Phys. Rev. E. 2001. Vol. 64, № 2. 021510.
10. Ovarlez G., Clement E. Phys. Rev. E. 2003. Vol. 68, № 3. 031302.
11. Josserand C., Tkachenko A.V., Mueth D.M., and Jaeger H.M. Phys. Rev. Lett. 2000. Vol. 85, № 17. P. 3632.
12. Havlin S., Kiefer J.E., and Weiss G.H. Phys.Rev. A. 1987. Vol. 36, № 3. 1403.
13. P. Le Doussal, Monthus C. and Fisher D.S. Phys. Rev. E. 1999. Vol. 59, № 5. 4795.
14. Aquino G., Bologna M., Grigolini P., and West B.J. Phys. Rev. E. 2004. Vol. 70, № 3. 036105.
15. Allegrini P., Aquino G., Grigolini P., Palatella L., and Rosa A. Phys. Rev. E. 2003. Vol. 68, № 5. 056123.
16. Barkai E. Phys. Rev. Lett. 2003. Vol. 90, № 10. 104101.
17. Barkai E. and Cheng Y.C., J. Chem Phys. 2003. Vol. 118, № 14. 6167.
18. Sokolov I.M., Blumen A. and Klafter J.Europhys. Lett. 2001. Vol. 56. 175; Sokolov I.M., Blumen A. and Klafter J.Physica A. 2001. Vol. 302. 268.
19. Monthus C. and Bouchaud J.P. J. Phys. A: Math. and Gen. 1996. Vol. 29, № 14. 3847.
20. Havlin S. and Ben-Avraham D. Adv. Phys. 2002. Vol. 51, № 1. 187.
21. Zaburdaev V.Yu. and Chukbar K.V. JETP Lett. 2003. Vol. 77. 551.
22. Chechkin A.V., Gorenflo R. and Sokolov I.M. J. Phys. A: Math. and Gen. 2005. Vol. 38. L679-L684.
23. Gordeche C. and Luck J.M. J. Stat. Phys. 2001. Vol. 104. 489.
24. Scalas E, Gorenflo R., and Mainardi F. Phys. Rev. E. 2004. Vol. 69. 011107.
25. Sokolov I.M., Chechkin A.V., and Klafter J.Acta Physica Polonica B. 2004. Vol. 35. 1323.
26. Tunaley J.K.E. Phys. Rev. Lett. 1974. Vol. 33. 1037.
27. Tunaley J.K.E. J. Stat. Phys. 1976. Vol. 14. 461.
28. Mantegna R.N. and Stanley H.E. Phys. Rev. Lett. 1994. Vol. 73. 2946.
29. Sokolov I.M., Chechkin A.V and Klafter J. Physica A. 2004. Vol. 336. 245.
30. Barsegov V. and Mukamel S. J. Phys. Chem. A. 2004. Vol. 108. 15.
31. Barkai E. and Sokolov I.M. J. Stat. Mech. 2007. P08001.
32. Chukbar K.V. JETP. 1995. Vol. 81. 1025.
33. Zaburdaev V.Yu. and Chukbar K.V. JETP. 2002. Vol. 94. 252.
34. Jaynes E.T. Phys. Rev. 1957. Vol. 106. 620.
35. BauleA. and Friedrich R. Europhys. Lett. 2007. Vol. 77. 10002.
Поступила в редакцию 13.07.2009
ДРОБНОЕ УРАВНЕНИЕ ДИФФУЗИИ ДЛЯ СТАРЕЮЩИХ И РАВНОВЕСНЫХ СЛУЧАЙНЫХ БЛУЖДАНИЙ
В.Ю. Забурдаев, И.М. Соколов
В настоящей работе рассматриваются случайные блуждания с непрерывным временем в ситуациях, описывающих старение процесса. Такие ситуации встречаются в том случае, если начальные условия известны не в момент приготовления системы (Ь = 0), а в более поздний промежуточный момент времени при Ь > 0. Для этого случая нами выведено обобщенное уравнение диффузии, содержащее одно ядро памяти. Полученное уравнение использовано для описания важного специального
случая блужданий в состоянии равновесия. Кроме того, обсуждены особенности различных постановок задач для проблем субдиффузии со старением и показано, что поведение системы зависит от того, как именно ставятся начальные условия в промежуточный момент времени: задана ли координата одной частицы (определенная измерением) или полная форма распределения вероятностей положений частиц. Такие две постановки задач ведут к различным предсказаниям о дальнейшей эволюции системы. Полученный результат подчеркивает важность правильного определения статистического ансамбля для стареющих систем.
Ключевые слова: Случайные блуждания с непрерывным временем, обобщенное уравнение диффузии, старение, статистический ансамбль.
Забурдаев Василий Юрьевич - родился (1979) в Москве, в 2002 году окончил Институт естественных наук и экологии (в настоящее время - факультет нано-, био-, информационных и когнитивных технологий МФТИ). По окончании института работал в Российском научном центре «Курчатовский Институт». В 2003 году защитил диссертацию на соискание ученой степени кандидата физико-математических наук в области теории аномальных процессов переноса в плазме. После защиты диссертации 3 года работал в институте Макса Планка в Гёттингене (Германия). Занимался исследованием динамики мокрых гранулярных сред. С 2007 года работает в институте теоретической физики Технического университета Берлина. Является автором 20 научных статей по проблемам транспортных процессов в плазме, нелинейной физики, аномальной диффузии и стохастического транспорта. Стипендиат Президента РФ и Фонда Сороса. Лауреат премии имени И.В. Курчатова. Hardenbergstrasse 36 D-10623 Berlin Germany Institut für Theoretische Physik, Technische Universität Berlin E-mail: vasily.zaburdaev@tu-berlin.de
Соколов Игорь Михайлович - родился (1958) в Москве, окончил физический ф-т МГУ (1981). Работал в отделе теоретической физики ФИАН, в настоящее время является профессором Гумбольдтского университета в Берлине (Германия). Защитил диссертацию на соискание ученой степени кандидата физико-математических наук (1984, МГУ), получил докторскую степень во Фрайбургском университете (1995). Научные интересы в области теории неупорядоченных систем, нелинейной динамики, физики полимеров, кинетики. Автор более 200 научных статей, монографии «Статистическая термодинамика и стохастическая теория неравновесных систем» (в соавторстве с В. Эбелингом), редактор монографии «Аномальный транспорт» (совместно с Р. Клагесом и Г. Радонсом) Член редакционной коллегии «Среднеевропейского физического журнала».
Newtonstrasse 15 D-12161 Berlin Germany Institut für Physik Humboldt-Universitat zu Berlin E-mail: igor.sokolov@physik.hu-berlin.de
