2006 International Conference in Honor of Jacqueline Fleckinger. Electronic Journal of Differential Equations, Conference 16, 2007, pp. ?? ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu ftp ejde.math.txstate.edu (login: ftp)

FRACTIONAL FLUX AND NON-NORMAL DIFFUSION

ALI ABDENNADHER, MARIE-CHRISTINE NEEL

Dedicated to Jacqueline Fleckinger on the occasion of an international conference in her honor

Abstract. Fractional diffusion equations are widely used for mass spreading in heterogeneous media. The correspondence between fractional equations and random walks based upon stable Lévy laws, keeps in analogy with that between heat equation and Brownian motion. Several definitions of fractional derivatives yield operators, which coincide on a wide domain and can be used in fractional partial differential equations. Then, the various definitions are useful in different purposes: they may be very close to some physics, or to numerical schemes, or be based upon important mathematical properties. Here we present a definition, which enables us to describe the flux of particles, performing a random walk. We show that it is a left inverse to fractional integrals. Hence it coincides with Riemann-Liouville and Marchaud’s derivatives when applied to functions, belonging to suitable domains.

1. Introduction

Fick’s law is a basic tool for the transport of dissolved matter. When it holds, the concentration of solute evolves according to heat equation. Nevertheless, there exist heterogeneous media where experimental evidence indicates that dispersion does not obey Fick’s law and heat equation: data from electronics [17] [18] and passive tracer experiments [3][4] show heavy tails, apparently connected with statistics giving some importance to extreme events, in similarity with densities of stable probability Lévy laws. Gaussian statistics are a limiting case of the latters, and distribute successive jump lengths of Brownian motions, which serve as small scale models for mass transport. On the macroscopic level, the concentration of a cloud of particles performing Brownian motions, satisfies Fick’s law and heat equation.

Data with heavy tails correspond to Lévy flights, which also are random walks, with jump lengths distributed according to stable Lévy laws. The correspondence between small and large scales is obtained by means of time and length references which we let tend to zero. If, moreover, they satisfy some scaling relation, the concentration of walkers tends to a limit, evolving according to a variant of heat

2000 Mathematics Subject Classification. 26A33, 45K05, 60E07, 60J60, 60J75, 83C40. Key words and phrases. Diffusion equations; fractional calculus; random walks; stable probability distributions. \begin{matrix}c \\ ◯ \end{matrix}2007 Texas State University - San Marcos. Published May 15, 2007.

equation, with Laplacean being replaced by derivatives of fractional order. They are non-local operators.

For unbounded domains, these results were derived from Generalized Master Equation via Fourier’s transform [5] [16]. They were extended to semi-infinite and bounded domains after some adaptations, under hypotheses which we want to escape [10]. Here we present a novel definition of Riemann-Liouville or Marchaud’s derivatives. It seems to us that it makes it possible to interpret fluxes directly for random walks observed from the macroscopic point of view, without passing through Generalized Master Equation, Fourier’s transform and space-fractional heat equation. Later, the gain in simplicity will allow us to tackle problems mixing, for instance, boundary conditions and skewness.

2. A new expression for Riemann-Liouville’s and Marchaud’s fractional derivatives

Riemann-Liouville and Marchaud’s fractional derivatives interpolate between integer orders of differentiation, and generalize many aspects of the notion. Real and complex orders of differentiation can be defined in this context. Here we aim at introducing a new way of presenting fractional calculus, in connection with particle counting and random walks. In this purpose, we only need to consider real valued orders of differentiation. And for the moment, in a sake of simplicity, we focus on one-dimensional problems, future works will be devoted to more general dimensions.

After having recalled essentials of the most widely used basic tools of fractional calculus, we show that a novel definition yields more or less similar objects. Then, we outline the connection with fluxes of particles.

2.1. Riemann-Liouville Fractional integrals and derivatives. Let α be positive (the definition applies to complex numbers with positive real part): the leftsided fractional integral of order α, computed over [a, x] is

I_{a,+}^{α}φ(x) = \frac{1}{Γ(α)} ∫_{a}^{x}(x − y)^{α−1}φ(y)dy  (2.1)

according to [15]. Here we will mainly consider the case a = −∞, with the simplified notation I_{+}^{α}φ(x) = I_{−∞,+}^{α}φ(x) of [14]. Right-sided integrals are

I_{b,−}^{α}φ(x) = \frac{1}{Γ(α)} ∫_{x}^{b}(y − x)^{α−1}φ(y)dy  (2.2)

with I_{−}^{α}φ(x) = I_{+∞,−}^{α}φ(x).

The corresponding left-sided Riemann-Liouville derivative of order α is

D_{+}^{α}φ(x) = (\frac{d}{dx})^{[α]+1}I_{+}^{1−\{α\}} = (\frac{d}{dx})^{[α]+1} \frac{1}{Γ([α] + 1 − α)} ∫_{−∞}^{x}(x − y)^{−\{α\}}φ(y)dy,

(2.3)

where [.] denotes integer part, while {.} is defined by α = [α]+{α}. The right-sided Riemann-Liouville derivative is

D_{−}^{α}φ(x) = (−\frac{d}{dx})^{[α]+1}I_{−}^{1−\{α\}} = (−\frac{d}{dx})^{[α]+1} \frac{1}{Γ([α] + 1 − α)} ∫_{x}^{+∞}(y−x)^{−\{α\}}φ(y)dy.

(2.4)

When α is a positive integer, D_{−}^{α} and D_{+}^{α} are usual right and left-sided derivatives or order α. A natural question is of whether fractional derivatives defined by (2.3) and (2.4) share with derivatives of integer order the property of being left inverses to the corresponding integrals. In fact, when φ is in L_{loc}^{1}(ℝ), if moreover the integrals I_{±}^{[α]+1} is absolutely convergent, we have (D_{±}^{α}I_{±}^{α}φ)(x) = φ(x) a.e., due to [14, Lemma 4.7]. When the above hypotheses are satisfied, D_{±}^{α} can be thought of as being a left inverse to I_{±}^{α}, which fails to hold if I_{±}^{1}φ and I_{±}^{[α]+1}φ do not belong to L^{1}(ℝ), even for φ in L^{p}(ℝ) with 1 < p < 1/α. Marchaud’s definition seems to give a more appropriate left inverse to fractional integrals.

2.2. Marchaud’s Fractional derivatives. Marchaud’s definition combines generalized finite differences and fractional integrals. A rather general definition of finite differences was given by [14] in the form

(Δ_{t}^{n}f)(x) = \frac{1}{d_{n}} |\begin{matrix}f(x − λ_{0}t) & 1 & λ_{0} & … & λ_{0}^{n−1} \\ f(x − λ_{1}t) & 1 & λ_{1} & … & λ_{1}^{n−1} \\ . & . & . & … & . \\ f(x − λ_{n}t) & 1 & λ_{n} & … & λ_{n}^{n−1} \end{matrix}| ,

with

d_{n} = |\begin{matrix}1 & λ_{1} & … & λ_{1}^{n−1} \\ . & . & … & . \\ 1 & λ_{n} & … & λ_{n}^{n−1} \end{matrix}| .

For λ_{0} = 0, λ_{1} = 1, … λ_{n} = n, and with T_{t} denoting the translation of amplitude t, the definition of Δ_{t}^{n} just becomes (Id−T_{t})^{n}f(x) = Σ_{j=0}^{n}(_{j}^{n})(−1)^{j}f(x−jt), which is f(x − t) − f(x) for n = 1. With these notations, Marchaud’s derivative D_{±}^{α} of function f is the limit, when ε → 0+ of

D_{±,ε}^{α}f(x) = \frac{1}{∫_{0}^{+∞} t^{−α−1}(1 − e^{−t})^{n}dt} ∫_{ε}^{+∞} t^{−α−1}Δ_{±t}^{n}f(x)dt,  (2.5)

with n > α. For 0 < α < 1, we have n = 1 and (5) becomes

D_{±,ε}^{α}f(x) = \frac{−1}{Γ(−α)} ∫_{ε}^{+∞} t^{−α−1}[f(x) − f(x ∓ t)]dt.

In this definition, the relative freedom let to Δ_{t}^{n}, is useful when α is a complex number. In view of our objective, we will focus on real valued orders for integrals and derivatives, hence λ_{0} = 0, λ_{1} = 1, … λ_{n} = n with Δ_{t}^{n}f(x) = (Id − T_{t})^{n}f(x), used by [15] will be enough. For α = 1, we have to put n = 2 in (5) if we want to use this expression, but we also can consider that D_{±}^{α} is the usual left or right-sided derivative of order α when α is a non-negative integer.

We thus have a left inverse for I_{±}^{α} in a wider domain, which in some sense is optimal, since it provides a characterization of I_{±}^{α}L^{p} for 1 < p < 1/α. Indeed, for 0 < α < 1, [15, Theorem 6.2] states the following: If the L^{p}(ℝ) limit of D_{±,ε}^{α}f exists when ε → 0+, or if sup_{ε>0}∥D_{±,ε}^{α}f∥_{L^{p}(ℝ)} is finite, if moreover, f belongs to L^{p}(ℝ) with 1 ≤ r < ∞, then f belongs to I_{±}^{α}L^{p}(ℝ) and there exists φ s.t. f(x) = I_{±}^{α}φ(x) almost everywhere. in ℝ. The theorem was stated in L^{p}(ℝ), but the proof adapts without any modification to L^{p}] − ∞, a] for D_{+}^{α} and to L^{p}[a, +∞[ for D_{−}^{α}. Derivatives D_{±}^{α} and D_{±}^{α} coincide for functions of the form I_{±}^{α}φ with φ in L_{loc}^{1} such that I_{±}^{[α]+1} converges absolutely [14].

Other expressions yield the left inverse of I_{±}^{α}. Among them, the Grünwald-Letnikov fractional derivative [15] or order α of f is the limit, when mesh h tends to zero, of h^{−α} times the series Σ_{k=0}^{∞}(−1)^{k}(_{k}^{α})f(x − kh). It provides useful approximations to Riemann-Liouville or Marchaud’s derivatives, connected with finite differences numerical schemes. Here we present a further expression for the left inverse of I_{±}^{α}, not very different from Grünwald-Letnikov operator, since it contains an integrals in place of the above evocated series. Then, we will discuss the physical meaning.

2.3. A new expression for the inverse of I_{±}^{α}. Here we consider 0 < α ≤ 1. With some modifications, the following adapts to all positive values of α.

Notation Let F be a function, while l is positive. Set

W_{l,±}^{α,F} f(x) = l^{−α−1} ∫_{0}^{+∞} f(x ∓ t)F (t/l)dt.  (2.6)

The limit of W_{l,±}^{α,F} f when l tends to zero, if it ever exists, will be denoted by W_{±}^{α,F} f.

• We will say that F satisfies Hypothesis (H1) if F belongs to L^{1}(ℝ^{+}) with ∫_{0}^{∞} F (t)dt = 0.

• We will say that F satisfies Hypothesis (H2) if, in a neighborhood [A, +∞[ of +∞, there exists a function F_{1} such that ∫_{A}^{+∞} y^{α}|F_{1}(y)|dy < ∞ and F (x) = F_{1}(x) + λx^{−α−1}, for 0 < α < 1 but F (x) = F_{1}(x) + λx^{−2−ε} with ε > 0 for α = 1.

We will see that (H1) and (H2) imply that W_{±}^{α,F} is a left inverse to I_{±}^{α}. In this purpose, let us consider W_{±}^{α,F} ∘ I_{±}^{α}.

Let φ belong to L^{p}[a, +∞[. We have

W_{l,−}^{α,F} ∘ I_{−}^{α}φ(x) = \frac{l^{−α−1}}{Γ(α)} ∫_{0}^{+∞} F (t/l) ∫_{x+t}^{+∞} φ(y)(y − x − t)^{α−1}dydt.

Setting t = lT yields

W_{l,−}^{α,F} ∘ I_{−}^{α}φ(x) = \frac{l^{−α}}{Γ(α)} ∫_{0}^{+∞} F (T) ∫_{x+lT}^{+∞} φ(y)(y − x − lT)^{α−1}dy dT

= \frac{1}{Γ(α)} ∫_{0}^{+∞} F (T) ∫_{T}^{+∞} φ(x + lθ)(θ − T)^{α−1}dθ dT,

with y = x + lθ. Then, Fubini’s theorem yields

W_{l,−}^{α,F} ∘ I_{−}^{α}φ(x) = \frac{1}{Γ(α)} ∫_{0}^{+∞} φ(x + lθ) ∫_{0}^{θ} F (T)(θ − T)^{α−1} dT dθ,  (2.7)

as soon as I_{+}^{α}(HF)(θ) = ∫_{0}^{θ} F (T)(θ − T)^{α−1} dT is integrable in ℝ^{+}. Let us use this point, which will be stated in Lemma 2.2 below. Let H denote Heaviside’s function: on the right-hand side of (2.7) we have ∫_{ℝ} φ(x + lθ)(I_{+}^{α}(HF))(θ)dθ which, by [15, Theorem 1.3] is an approximation of ∫_{0}^{+∞} I_{+}^{α}(HF)(θ)dθ times Identity in L^{p}.

For φ in L^{p}] − ∞, a], instead of (2.7) we have

W_{l,+}^{α,F} ∘ I_{+}^{α}φ(x) = \frac{1}{Γ(α)} ∫_{0}^{+∞} φ(x − lθ) ∫_{0}^{θ} F (T)(θ − T)^{α−1} dT dθ.  (2.8)

Hence the following Theorem holds.

Theorem 2.1. Suppose F satisfies hypotheses (H1) and (H2), with 0 < α ≤ 1.

(i) For φ in L^{p}[a, +∞[, W_{l,−}^{α} ∘I_{−}^{α}φ tends in L^{p}[a, +∞[ to ∫_{0}^{+∞} I_{+}^{α}HF (t)dt×φ when l tends to zero, and pointwise everywhere φ is right continuous.

(ii) For φ in L^{p}]−∞, a], W_{l,+}^{α} ∘I_{+}^{α}φ tends in L^{p}]−∞, a] to ∫_{0}^{+∞} I_{+}^{α}HF (t)dt×φ when l tends to zero, and pointwise everywhere φ is left continuous.

It remains to prove the following lemma.

Lemma 2.2. If F satisfies (H1) and (H2), with 0 < α ≤ 1, then ∫_{0}^{θ} F (T)(θ − T)^{α−1} dT is integrable in ℝ^{+}.

Proof. If F is as F_{1} in hypothesis (H2)’s statement, [14, Lemma 4.12] shows that

Γ(α)I_{−}^{α}(HF)(θ) = ∫_{0}^{θ} F (T)(θ − T)^{α−1} dT

is in L^{1}. It is enough to prove the lemma for F = −\frac{1}{α}χ_{[0,1]} + x^{−α−1}χ_{[1,+∞[} if α is less than 1, for F = −\frac{1}{1+ε}χ_{[0,1]} + x^{−2−ε}χ_{[1,+∞[} if α is equal to 1, since modifying F_{1} will immediately lead to the general case. For α = 1, the result is obvious, for α less than 1, we have

∫_{0}^{x}(x − y)^{α−1}χ_{[0,1]}(y)dy = \frac{x^{α} − (x − 1)^{α}}{α}

for x > 1, and

∫_{0}^{x}(x − y)^{α−1}y^{−α−1}χ_{[1,+∞[}(y)dy = x^{−1}(G(1) − G(1/x) + \frac{x^{α} − 1}{α})

when x is large enough, with G being defined by G(X) = ∫_{0}^{X}[(1−z)^{α−1}−1]z^{−α−1}dz. From this we deduce

\begin{matrix}∫_{0}^{x} F (t)(x − t)^{α−1}dt \\ = α^{−1}(x^{α−1} − α^{−1}x^{α}(1 − (1 − 1/x)^{α})) + x^{−1}(G(1) − α^{−1}) − x^{−1}G(1/x). \end{matrix} (2.9)

Function \frac{(1−t)^{α−1}−1}{t} is continuous and integrable in [0, 1[. In the neighborhood of

0, \frac{(1−t)^{α−1}−1}{t}t^{−α} is equivalent to (1 − α)t^{−α}, hence G(1/x) is equivalent to x^{α−1} when x is large. Hence x^{−1}G(1/x) is integrable in a neighborhood of +∞. It also is the case for α^{−1}[x^{α−1} − α^{−1}x^{α}(1 − (1 − 1/x)^{α}). We now will check that G(1) − α^{−1} is zero. To see this, set g(p, q) = ∫_{0}^{1}((1 − t)^{q−1} − 1)t^{p−1}dt. For complex valued p and q satisfying Re(p) > 0 and Re(q) > 0, ∫_{0}^{1}(1 − t)^{q−1}t^{p−1}dt is a beta function [1] and we have

g(p, q) = \frac{Γ(p)Γ(q)}{Γ(p + q)} − \frac{1}{p}.  (2.10)

Let us fix q = α, and vary the complex number p: t^{p} is a function of p, whose derivative t^{p}Ln(t) is dominated by the L^{1}]0, 1[ function t^{p}|Ln(t)| for Re(p) ≥ p_{0} > −1, so that, by dominated convergence, g(p, α) is derivable with respect to p. Hence it is analytic for Re(p) ≥ p_{0} > −1. Since \frac{Γ(q)}{Γ(p+q)} is also analytic in the neighborhood of 0 while Γ(p) has a simple pole with residuum 1, the right-hand side of (2.10) is holomorphic for Re(p) ≥ p_{0} > −1. Hence relation (2.10) holds for p = −α, and Lemma 2.2 is proved. □

Therefore Theorem 2.1 holds. It states that the operator W_{l,−}^{α,F}, which is defined on I_{−}^{α}L^{p}[a; +∞[, has a limit in L^{p}[a; +∞[ when l tends to zero. Up to multiplication by a function of F and α, the limit is a left inverse to I_{−}^{α}, hence it coincides with D_{−}^{α}. Similarly, W_{l,+}^{α,F}, defined in I_{+}^{α}L^{p}] − ∞, a], tends to D_{+}^{α}, times a function of F and α. Theorem 2.1 adapts to higher values of α, provided Hypothesis 1 is made stronger.

For values of α between 0 and 1, will see that Theorem 2.1 allows us to represent the flux of particles within the frame work of a wide class of Random Walks.

3. Particles flux for Lévy Flights, in the macroscopic limit

Brownian Motion is a particular case of Lévy Flights. The latters are Continuous Time Random Walks: a large number of particles perform a succession of independent jumps, whose lengths X_{i} are identically distributed. To be more precise, with l being a length scale, the density of X_{i}/l is the normalized αstable Lévy density L_{α,θ} of exponent α between 1 and 2 and with skewness parameter θ (see Appendix A). Waiting times T_{i} between successive jumps are such that the independent random variables T_{i}/τ have density ψ, whose average is 1. Here, for definiteness, we set ψ(t) = e^{−t}. Looking at the cloud of particles from the macroscopic point of view means that we let length and time scales l and τ tend to zero. Then, if the scaling relation l^{α}/τ = K holds [5] [16], the probability of finding a particle in a given interval tends to a limit, which has a density satisfying a space-fractional diffusion equation such as (3.7). This implies that the flux of particles satisfies a fractional generalization of Fick’s law [13]. All these results are based upon Generalsized Master Equation and Fourier’s analysis.

In fact, we will see that Theorem 2.1 connects more directly particle flux and fractional derivatives.

3.1. Computing the flux for Lévy flights with length scale l and mean waiting time τ satisfying l^{α}/τ = K. For a given particle, the location after n^{th} jump is Σ_{i=0}^{n}X_{i}, and it happens at time Σ_{i=0}^{n}T_{i}. Let us denote by μ(., t) the measure giving the probability μ(I, t) that the particle be in interval I at time t. With this notation, the balance of particles crossing abscissa x during [t, t + dt[ is the difference of two expressions. The first one is the probability ∫_{−∞}^{x} F_{α,θ}^{d}(\frac{x−y}{l})dμ(y, t)\frac{ψ(0)}{τ}dt of crossing x to the right, with F_{α,θ}^{d}(y/l) being the probability ∫_{y}^{+∞} \frac{1}{l}L_{α,θ}(z/l)dz = ∫_{y/l}^{+∞} L_{α,θ}(z)dz for a jump to have an amplitude

of more than y. The second is the probability ∫_{x}^{+∞} F_{α,θ}^{g}(\frac{x−y}{l})dμ(y, t)\frac{ψ(0)}{τ}dt of

crossing x to the left, with F_{α,θ}^{g}(−y/l) being the probability ∫_{−∞}^{−y/l} \frac{1}{l}L_{α,θ}(z/l)dz for a jump to have an amplitude of more than y, but to the left. The flux is the probability rate, hence the following difference:

Kl^{−α}[ ∫_{−∞}^{x} F_{α,θ}^{d}(\frac{x − y}{l})dμ(y, t) − ∫_{x}^{+∞} F_{α,θ}^{g}(\frac{x − y}{l})dμ(y, t)].

When μ(., t) has density C(., t), the flux Q_{l}^{α,θ}C(., t)(x) is given by

Q_{l}^{α,θ}C(., t)(x) = Kl^{−α}[ ∫_{0}^{+∞} C(x − y, t)F_{α,θ}^{d}(\frac{y}{l})dy − ∫_{0}^{+∞} C(x + y)F_{α,θ}^{g}(\frac{−y}{l})dy].

(3.1)

Both integrals are similar to (2.6), except that F_{α,θ}^{d} and F_{α,θ}^{g}(−.) satisfy (H2) with α − 1 instead of α, according to Appendix A, but of course not (H1). Appendix B shows that ∫_{0}^{+∞} F_{α,θ}^{d}(y)dy = ∫_{0}^{+∞} F_{α,θ}^{g}(−y)dy = I_{α,θ}. Hence, with f_{α,θ} being a compactly supported function of class L^{1} s.t. ∫_{0}^{+∞} f_{α,θ}(y)dy = I_{α,θ}, setting F̃_{α,θ}^{d} = F_{α,θ}^{d} − f_{α,θ} and F̃_{α,θ}^{g}(−y) = F_{α,θ}^{g}(−y) − f_{α,θ}(y) yields functions, satisfying (H1) and (H2) with α − 1 instead of α. Then, we have Q_{l}^{α,θ}C(., t)(x) = Q_{+,l}^{α,θ}C(., t)(x) − Q_{−,l}^{α,θ}C(., t)(x), with

\begin{aligned}Q_{+,l}^{α,θ}f(x) & = Kl^{−α} ∫_{0}^{+∞} F_{α,θ}^{d}(y/l)(f(x − y) − f(x))dy \\ & = K[(W_{l,+}^{α−1,F̃_{α,θ}^{d}}f)(x) + l^{−α} ∫_{0}^{+∞} f_{α,θ}(y/l)(f(x − y) − f(x))dy] \end{aligned} (3.2)

at the left of x, and

\begin{aligned}Q_{−,l}^{α,θ}f(x) & = Kl^{−α} ∫_{0}^{+∞} F_{α,θ}^{g}(−y/l)(f(x + y) − f(x))dy \\ & = K[(W_{l,−}^{α−1,F̃_{α,θ}^{g}(−.)}f))(x) + l^{−α} ∫_{0}^{+∞} f_{α,θ}(y/l)(f(x + y) − f(x))dy] \end{aligned}

(3.3)

at the right. Since F̃_{α,θ}^{d} and F̃_{α,θ}^{g}(−.) satisfy (H1) and (H2) with α − 1 instead of α,

hence (W_{l,+}^{α−1,F̃_{α,θ}^{d}}f)(x) tends to ∫_{0}^{+∞} I_{+}^{α−1}(HF̃_{α,θ}^{d})(y)dyD_{+}^{α−1}(f)(x) in L^{p}] − ∞, a] when f belongs to I_{+}^{α−1}L^{p}] − ∞, a] and

(W_{l,−}^{α−1,F̃_{α,θ}^{g}(−.)}f)(x) tends to ∫_{0}^{+∞} I_{+}^{α−1}(HF̃_{α,θ}^{g}(−.))(y)dyD_{−}^{α−1}(f)(x)

in L^{p}[a, +∞[ when f belongs to I_{−}^{α−1}L^{p}[a, +∞[. We will see that appropriately choosing f_{α,θ} allows us to see on the right-hand sides of (3.2) and (3.3) expressions which are “local fractional derivatives”, in the sense of Kolwankar and Gangal.

3.2. Kolwankar and Gangal’s local fractional derivatives. The notion of “a local fractional derivative” was introduced [9] in view of building a tool, designed for the study of continuous but nowhere differentiable functions frequently occurring in Nature and economics. Those fractional derivatives share some properties with previously defined ones, such as chain rule or generalized Leibniz rule [2]. They are very useful for to compute fractal dimensions of graphs. In fact, they vanish for smooth enough functions, and hence can become “invisible”.

For q between 0 and 1, the right-sided Kolwankar and Gangal’s [9] fractional derivative of order q of function f, computed at x, will be denoted by

D_{+}^{KG,q}f(x) = lim_{h→0+} \frac{d}{dh}I_{x,+}^{1−q}(f(.) − f(x))(x + h).

Let us suppose that f is continuous is [x, x + ε], with positive ε. When the limit exists, it is equal to the limit, when h tends to 0+, of h^{−1}I_{x,+}^{1−q}(f(.) − f(x))(x + h), due to l’Hôpital’s rule and to lim_{h→0}(I_{x,+}^{1−q}(f(.) − f(x))(x + h)) = 0 . Moreover, we have h^{−1}I_{x,+}^{1−q}(f(.) − f(x))(x + h) = \frac{h^{−q}}{Γ(1−q)} ∫_{0}^{1}(1 − t)^{−q}(f(x + th) − f(x))dt.

At the left, we have

D_{−}^{KG,q}f(x) = lim_{h→0+} \frac{d}{dh}I_{x,−}^{1−q}(f(x) − f(.))(x − h),

also equal to the limit, when h tends to 0+, of h^{−1}I_{x,−}^{1−q}(f(x) − f(.))(x − h). If, for positive and finite b − a and with q < q + ε < 1, function f belongs to Hölder space H^{q+ε}[a, b], then sup_{x,y∈[a,b]}(\frac{|f(y)−f(x)|}{|y−x|^{q+ε}}) is finite, hence Kolwankar and Gangal’s local derivatives of order q are zero in [a, b].

For α < 2, an appropriate choice of f_{α,θ} yields that the second expressions on the right-hand sides of (3.2) and (3.3) tend to right and left-sided local fractional derivatives of order α − 1. The choice is f_{α,θ}(t) = I_{α,θ}(2 − α)(1 − t)^{1−α}χ_{[0,1]}; then, we have l^{−α} ∫_{0}^{+∞} f_{α,θ}(y/l)(f(x+y)−f(x))dy = I_{α,θ}l^{1−α}(2−α) ∫_{0}^{1}(1−t)^{1−α}(f(x+ lt) − f(x))dt. Consequently, when f has a local derivative of order α − 1 at point x+, l^{−α} ∫_{0}^{+∞} f_{α,θ}(y/l)(f(x + y) − f(x))dy has a limit when l tends to zero, and the limit is I_{α,θ}Γ(3 − α) times the right-sided local derivative of order 1 − α. The same holds at the left of x: l^{−α} ∫_{0}^{+∞} f_{α,θ}(y/l)(f(x − y) − f(x))dy tends to I_{α,θ}Γ(3 − α) times the left-sided local derivative of order 1 − α. When f is differentiable at x, the limit is zero.

For α = 2, we have the usual derivative instead of local Kolwankar and Gangal’s derivatives, the above choice of f_{α,θ} is no longer relevant, but the end of next subsection will show that the method yields Fick’s law simply and directly.

We now are ready for looking at the limit of operator flux Q_{l}^{α,θ} when l tends to zero.

3.3. The flux, in the l → 0 limit. We know from (2.1) that the flux of particles, passing trough x at time t, is Q_{l}^{α,θ}C(., t)(x) = Q_{+,l}^{α,θ}C(., t)(x)−Q_{−,l}^{α,θ}C(., t)(x). Moreover, according to (3.2) and (3.3), Q_{±,l}^{α,θ}f(x) splits into two terms. For 0 < α < 2,

with moreover D_{±}^{α−1}f ∈ L^{p}, the first one, KW_{l,±}^{α−1,F̃_{α,θ}^{d}}f, tends to Kλ_{±}D_{±}^{α−1}f in L^{p}, according to Theorem 2.1, with

λ_{+} = ∫_{0}^{+∞} I_{+}^{α−1}HF̃_{α,θ}^{d}(y)dy,  (3.4)
λ_{−} = ∫_{0}^{+∞} I_{+}^{α−1}HF̃_{α,θ}^{g}(−y)dy.  (3.5)

The second term Kl^{−α} ∫_{0}^{+∞} f_{α,θ}(y/l)(f(x∓y)−f(x))dy in (3.2) and (3.3) tends to KI_{α,θ}D_{±}^{KG,α−1}f at points where the local fractional derivative exists. Parameter I_{α,θ} can be computed numerically, with the help of integral expressions [20] for stable Lévy distributions. It also can be deduced from Appendix B. Oppositely, direct computation of λ_{±} from (3.4) and (3.5) is not easy, but can be avoided and

replaced by directly checking W_{±}^{α−1,F̃_{αθ}^{d,g}}f and D_{±}^{α−1}f, with particular functions f in I_{+}^{α−1}] − ∞, a] or I_{−}^{α−1}[a, +∞[ such that the local derivative is zero.

Let us show that considering f = χ_{[1,2[} yields λ_{−}. Indeed, for x ∈]1, 2[, the local derivative exists and is equal to zero, while we have

Q_{−,l}^{α,θ}χ_{[1,2[}(x) = l^{−α} ∫_{2−x}^{+∞} F_{α,θ}^{g}(−y/l)dy = C_{α,−θ} \frac{(2 − x)^{1−α}}{α − 1} + O(l)

for 1 < α < 2, with C_{α,−θ} being defined by (3.7). But

D_{−}^{α−1}χ_{[1,2[}(x) = \frac{−1}{1 − α} ∫_{2−x}^{+∞} y^{−α}dy = −\frac{(2 − x)^{1−α}}{Γ(2 − α)} .

For x ≥ 2, we have D_{−}^{α−1}χ_{[1,2[}(x) = 0, hence D_{−}^{α−1}χ_{[1,2[} belongs to L^{p}[1, +∞[ for 1 ≤ p < \frac{1}{α−1}, so that χ_{[1,2[} belongs to I_{−}^{α−1}L^{p}[1, +∞[. Hence, according to [15, Theorem 6.2], we have λ_{−} = −\frac{Γ(2−α)}{α−1}C_{α,−θ} = \begin{matrix}sin & \frac{π}{2} (α+θ) \\ sin & πα \end{matrix}, and similarly λ_{+} =

\begin{matrix}sin & \frac{π}{2} (α−θ) \\ sin & πα \end{matrix}. Hence for sufficiently well-behaved functions (in L^{p}(ℝ) ∩ H^{α−1+ε}(ℝ)) the limit of the operator, giving the flux is

f ↦ K(\begin{matrix}sin & \frac{π}{2} (α − θ) \\ sin & πα \end{matrix}D_{+}^{α−1}f − \begin{matrix}sin & \frac{π}{2} (α + θ) \\ sin & πα \end{matrix}D_{−}^{α−1}f),  (3.6)

which is a fractional variant of Ficks law. In fact, for α = 2, the method has to be slightly adapted but yields Fick’s law itself.

In subsection 2.2, we pointed out that case α = 2 has to be considered separately. To do this, take A(l) a function, tending to +∞ when l tends to zero, with lA(l) tending to zero: for instance, we can choose A(l) = l^{−1/2}. Parameter θ is equal to zero, L_{α}^{θ} is even and superscripts d and g in F_{2,0}^{g,d} are of no use: instead we put F_{2,0}. We have

Q_{−,l}^{2,0}f(x) = Kl^{−2} ∫_{0}^{+∞} F_{2,0}(y/l)(f(x + y) − f(x))dy

= Kl^{−1} ∫_{0}^{+∞} F_{2,0}(y)(f(x + ly) − f(x))dy

= K ∫_{0}^{A(l)} F_{2,0}(y)\frac{f(x + ly) − f(x)}{ly}dy

+ l^{−1}K ∫_{A(l)}^{+∞} F_{2,0}(y)(f(x + ly) − f(x))dy.

If f is differentiable at point x,

∫_{0}^{A(l)} F_{2,0}(y)\frac{f(x + ly) − f(x)}{ly}dy

tends to the usual derivative f′(x), times ∫_{0}^{+∞} F_{2,0}(y)ydy, itself equal to 1/2 due to F_{2,0}(x) = ∫_{x}^{+∞} \frac{1}{2\sqrt{π}} e^{−y^{2}/4}dy. And l^{−1}| ∫_{A(l)}^{+∞} F_{2,0}(y)(f(x + ly) − f(x))dy| is less

than l^{−2}F_{2,0}(A(l)) ∫_{lA(l)}^{+∞} |f(x + y) − f(x)|dy, which tends to 0 when f is fixed in L^{1}. Similar results are obtained on the left side of x, hence for α = 2, in the limit “l tends to zero” operator flux tends to K \frac{d}{dx}f(x) which satisfies Fick’s law.

The fractional version (3.6) implies a space-fractional variant of heat equation.

3.4. Space-fractional heat equation. For functions f of the form I_{±}^{α−1}φ with φ in L^{p}(ℝ) ∩ L^{1}(ℝ), if, moreover, f belongs to H^{α−1+ε}, (3.6) is of the form

\frac{K}{Γ(2 − α)} \frac{d}{dx}[\begin{matrix}sin & \frac{π}{2} (α − θ) \\ sin & πα \end{matrix} ∫_{−∞}^{x}(x − y)^{1−α}f(y)dy

+ \begin{matrix}sin & \frac{π}{2} (α + θ) \\ sin & πα \end{matrix} ∫_{x}^{+∞}(y − x)^{1−α}f(y)dy],

which yields operator flux Q^{∗} for particles performing Lévy Flights in the diffusive limit (l and τ tend to zero with l^{α}/τ = K), when f is the concentration at time t. Casting the expression for Q^{∗}C^{∗} into mass conservation law ∂_{t}C^{∗} = −∂_{x}Q^{∗}C^{∗} yields for the evolution of the concentration C^{∗}(x, t)

∂_{t}C^{∗}(x, t) = K∇_{x,α}^{θ}C^{∗}(x, t)  (3.7)

with ∇_{x,α}^{θ} being the Riesz-Feller fractional derivative of order α and skewness parameter θ [8], defined by

(∇_{x,α}^{θ}f)(x) = \frac{−1}{Γ(2 − α)} \frac{d^{2}}{dx^{2}} [\begin{matrix}sin & \frac{π}{2} (α − θ) \\ sin & πα \end{matrix} ∫_{−∞}^{x}(x − y)^{1−α}f(y)dy

+ \begin{matrix}sin & \frac{π}{2} (α + θ) \\ sin & πα \end{matrix} ∫_{x}^{+∞}(y − x)^{1−α}f(y)dy],

(3.8)

in agreement with [16].

Two points were essential in the reasoning, leading to (3.6) and (3.8). The first one is the asymptotic behavior, for x → +∞, of the cumulated probabilities P {X > x} and P {X < −x} for a jump to have an amplitude X larger than x, and directed to the right or to the left. Both probabilities have to satisfy (H2) with α − 1 instead of α. An other property of stable laws was used, in view of (H1): it is the fact that the integrals, over ℝ^{+}, of P {X > x} and P {X < −x}, are equal. This allowed us to subtract this integral times f(x) from both sides of the difference, giving the flux, without any net change. In fact, any Continuous Time Random Walk made of successive independent jumps, identically distributed according to a random variable lX satisfying both conditions, has a flux whose limit is (3.6) when l tends to zero, provided mean waiting time τ exists with also l^{α}/τ = K.

Conclusion

Among many objects, interpolating between derivatives of integer orders, several tools termed fractional derivatives, were designed for various purposes. Some of them are connected with the idea that integration and derivation are inverses of each other.

Within this frame work, there are several ways for to define fractional derivatives, which are more or less similar to each other. They are more or less interesting, according to the sets of functions, which we want them to operate on. Among them, Grünwald-Letnikov derivatives led to performing numerical schemes. Theorem 2.1 indicates a novel definition of fractional derivatives, not so far from Grünwald-Letnikov’s: an integral replaces a series. It seems to be appropriate for to represent fluxes of particles performing Continuous Time Random Walks satisfying some hypotheses. Among them, Lévy flights play an important role, since stable laws are ubiquitous in Nature. We developped this point for random walks in a free one-dimensional space, also using the local derivatives invented by Kolwankar and Gangal for fractal graphs.

Combining those objects also applies to situations with boundary conditions, for instance in a half space {x ∈ ℝ/x > 0} limited by a wall at x = 0. There are several possibilities for the interaction between wall and particles. For instance, we can imagine that they do not exchange any energy: particles bouncing on the wall continue the distance, they had to fly if there were no wall, but they stay on the same side. Then, when writing down the balance of particles crossing abscissa x, we have to take into account that, among random walkers flying to the left (in the direction of the wall) some of them bounce and come back to the right of x: they have to be excluded from balance W_{l,−}. If a particle located in y ∈]0, x[ has to jump a length larger than |y − 2x| to the left, it arrives at the right of the wall and has to be taken into account in W_{l,+}. By doing this, we obtain that the flux through x is the sum of two terms. The first one is the flux corresponding to a concentration profile equal to the even extension of the actual one, to a free space without any wall. The second term is proportional to the left-sided Riemann-Liouville or Marchaud’s derivative of order α − 1, if α denotes the stability exponent of the jump length distribution of the Lévy flight. It contains a factor, which becomes zero when the distribution is symmetric, in agreement with results [10], previously obtained by an other method.

Appendix A: Densities of alpha stable Lévy laws

Stable laws are a generalization of Gaussian statistics. In many occasions, and here also, the word “stable” refers to some property, invariant under a definite set of transformations, as in the following definition. Definition: Let X be a random variable, distributed according to the probability law F. Random variable X and law F are said to be stable if [11] for every (a_{1}, a_{2}) ∈ ℝ^{+2} and (b_{1}, b_{2}) ∈ ℝ^{2}, there exist a ∈ ℝ^{+} and b ∈ ℝ such that F (a_{1}x + b_{1}) ∗ F (a_{2}x + b_{2}) (the law of the random variable a_{1}X_{1} + b_{1} + a_{2}X_{2} + b_{2}, with X_{1} and X_{2} being independent and distributed according to F) be equal to F (aX + b).

When F is as in the above definition, for any sequence of independent random variables X_{i} identically distributed according to F, there exists a sequence c_{n} of positive numbers such that \frac{X_{1}+...X_{n}}{c_{n}} be distributed according to F itself for any positive integer n [6] [7]. Moreover, c_{n} is a power of n, and the inverse α of the exponent belongs to ]0, 2] and serves as a label for the law: it is called the stability exponent of the law, which is said to be α stable. For α = 2 we have normal law, which is symmetric. For α ∈]0, 2[, stable laws may be symmetric or skewed. Stable laws play an important role in Nature because they are attractors, which are defined below. Definition: Let F be the probability of a sequence of independent random variables X_{n}. The probability law G is an attractor for F if there exists sequences A_{n} and B_{n}, with B_{n} > 0, such that the law of \frac{X_{1}+···+X_{n}}{B_{n}} − A_{n} tends to G when n tends to ∞ [6].

Loosely speaking, α stable laws are attractors for probability laws whose density behaves asymptotically as x^{−α−1} if α belongs to ]0, 2[, normal law (with α = 2)is an attractor for probability laws whose asymptotics is x^{−α′−1} with α′ ≥ 2 [6] [7].

Except for some values (e.g. α = 1 or 2), the density of a stable law cannot be given in closed form. But, up to translations and dilatations, the Fourier transform is e^{−|k|^{α}e^{isign(k)πθ/2}}. The corresponding density L_{α}^{θ} satisfies L_{α}^{θ}(−x) = L_{α}^{−θ}(x). Up to dilatations and translations, two labels determine stable densities: the stability exponent α, and the skewness parameter θ, which belongs to [α − 2, 2 − α].

In neighborhoods of ∞, except for α = 2, L_{α}^{θ} behaves as a negative power of the variable [19] [12]. For 1 < α < 2, α − 2 < θ ≤ 2 − α and x > A > 0, we have

L_{α}^{θ}(x) = \frac{1}{πx}Σ_{n=1}^{+∞}(−x^{−α})^{n} \frac{Γ(1 + nα)}{n!} sin \frac{nπ}{2} (θ − α).  (3.9)

We will denote by C_{α}^{θ} = \frac{−1}{π}Γ(1 + α) sin \frac{π}{2} (θ − α) the coefficient of the leading term in expansion (3.9).

Appendix B: Integrals of cumulated alpha stable Lévy laws

Due to symmetry, the integrals ∫_{0}^{+∞} F_{α,θ}^{d}(y)dy and ∫_{0}^{+∞} F_{α,θ}^{g}(−y)dy are equal for θ = 0. In fact, and this point is important for us, this equality holds for all admissible values of θ. Let us prove the claim.

First, notice that F_{α,θ}^{g}(−x) = ∫_{−∞}^{−x} L_{α}^{θ}(y)dy = ∫_{x}^{+∞} L_{α}^{−θ}(−y)dy = F_{α,−θ}^{d}(x). Then, we will uses Mellin’s transform, defined by Mω(z) = ∫_{0}^{+∞} t^{z−1}ω(t)dt for function ω. With z = 1 we see that ∫_{0}^{+∞} F_{α,θ}^{d}(y)dy = MF_{α,θ}^{d}(1), while we have F_{α,θ}^{d}(x) = I_{−}^{1}L_{α}^{−θ}(x), hence ∫_{0}^{+∞} F_{α,θ}^{d}(y)dy = (MI_{−}^{1}L_{α}^{−θ})(1).

For z > 1 and sufficiently good-behaved functions in neighborhoods of ∞, such as L_{α}^{θ}, we have

(MI_{−}^{1}ω)(z) = \frac{Γ(z)}{Γ(z + 1)}(Mω)(z + 1),

for z < α according to [14] page 44. From this, due to F_{α,θ}^{d}(x) = ∫_{x}^{+∞} L_{α}^{θ}(y)dy = I_{−}^{1}L_{α}^{θ}(x), we deduce

(MF_{α,θ}^{d})(z) = \frac{Γ(z)}{Γ(z + 1)}(ML_{α}^{θ})(z + 1).

The Mellin transform ML_{α}^{θ} is given in [19]:

(ML_{α}^{θ})(z) = \frac{1}{α} \frac{Γ(z)Γ((1 − z)α^{−1})}{Γ((1 − z)\frac{α−θ}{2α})Γ(1 − (1 − z)\frac{α−θ}{2α})},

which is of the form

(ML_{α}^{θ})(z) = \frac{1}{πα}Γ(z)Γ(\frac{1 − z}{α}) sin ((1 − z)π \frac{α − θ}{2α})  (3.10)

due to complements formula for Gamma functions [1]. In fact, [19] proved (20) for 0 < Re(z) < 1. Nevertheless, ML_{α}^{θ}(z), as a function of z, is holomorphic for 0 < Re(z) < α + 1, due to the behavior of L_{α}^{θ}(x) for large real values of x. On the right-hand side of (3.10), Γ(z)Γ((1 − z)α^{−1}) is holomorphic also except at poles of Γ((1−z)α^{−1}), which means that we have to exclude 1 from {z ∈ ℂ/0 < Re(z) < α+ 1}. Then, analytic continuation extends (3.10) to {z ∈ ℂ/0 < Re(z) < α+1}−{1}.

¿From this we deduce

∫_{0}^{+∞} F_{α,θ}^{d}(y)dy = (MF_{α,θ}^{d})(2)

= \frac{Γ(2)Γ(−1/α)}{απ} sin π \frac{θ − α}{2α}

= −\frac{Γ(−1/α)}{απ} cos π \frac{θ}{2α}.

We see that ∫_{0}^{+∞} F_{α,θ}^{d}(y)dy is an even function of θ, hence the claimed result.

References

[1] M. Abramowitz and I. Stegun; Handbook of mathematical functions, Dover, New York 1965.

[2] A. Babakhani and Varsha Daftardar-Gejji; On calculus of local fractional derivatives, J. Math. Anal. Appl. 270, 66-79 (2002).

[3] D. A. Benson, S. W. Wheatcraft and M. M. Meerschaert; The fractional-order equation of Lévy motion, Water Resour. Res. 36, 1413-1424 (2000).

[4] D. A. Benson, R. Schumer, M. M. Meerschaert and S. W. Wheatcraft, Fractional dispersion, Lévy motion, and the MADE tracer tests, Transp. Por. Med. 42, 211-240 (2001).

[5] A. Compte; Stochastic foundations of fractional dynamics, Phys. Rev. E, 53, 4191 (1996).

[6] W. Feller; An Introduction to Probability Theory and its Applications, vol. II, Page 581, Wiley series in probability and mathematical statistics, J. Wiley and sons, New York, Chichester, Brisbane, Toronto (1970).

[7] B.V. Gnedenko and A.N. Kolmogorov; Limit distributions for sums of independent variables, Addison Wesley (1968).

[8] R. Gorenflo and F. Mainardi; Random Walks models for space-fractional diffusion processes, Fractional Calculus and Applied Analysis, 1 (2), 167, (1999).

[9] K. M. Kolwankar and A. D. Gangal; Fractional differentiation of nowhere differentiable functions and dimension, CHAOS 6, 505-513 (1996).

[10] N. Krepysheva, L. Di Pietro and M.C. Néel; Fractional diffusion and reflective boundary condition, Physica A(368), 355-361 (2006).

[11] P. Lévy; Théorie de l’addition des variables aléatoires, Gauthier-Villars, Paris (1937).

[12] F. Mainardi, Y. Luchko, G. Pagnini; The fundamental solution of the space-time fractional diffusion equation, Fractional Calculus and Applied Analysis, vol. 4,(2) pp 153-192 (2001).

[13] P. Paradisi, R. Cesari, F. Mainardi, and F. Tampieri; The fractional Fick’s law for non-local transport processes, Physica A 293(1-2), 130 (2001).

[14] B. Rubin, Fractional integrals and potentials, Harlow: Longman (1996).

[15] S. G. Samko, A. A. Kilbas and O. I. Marichev; Fractional integrals and derivatives: theory and applications (Gordon and Breach, New York (1993)).

[16] E. Scalas, R. Gorenflo, F. Mainardi; Uncoupled continuous-time random walks: solutions and limiting behaviour of the master equation, Phys. Rev E, 692, 011107 (2004).

[17] H. Scher and M. Lax; Stochastic transport in a disorered solid. I. Theory, Phys. Rev B, 710, 4491-4502 (1973).

[18] H. Scher and M. Lax; Stochastic transport in a disorered solid. II. Theory, Phys. Rev B, 710, 44502-4519 (1973).

[19] W. R. Schneider; Stable distributions: Fox function representation and generalization, in Stochastic process in classical and quantum systems , Ed. S. Albeveiro, G. Casatti, D. Merlini, Springer, Berlin, pp 497-511, (1986).

[20] R. Weron; On the Chambers-Mallows-Stuck method for simulating skewed stable random variables, Stat. and Prob. Lett. 28 pp 165-171, (1996).

Ali Abdennadher Department of mathematics, Institut National des Sciences Appliquées et de Technologie, Centre Urbain Nord, BP 676 Cedex 1080 Charguia Tunis, Tunisie

E-mail address: ali abdennadher@yahoo.fr

Marie-Christine Neel UMRA 1114 Climat Sol Environnement, University of Avignon, 33, rue Pasteur, 84000 Avignon, France

E-mail address: marie-christine.neel@univ-avignon.fr