Generalized fractional Poisson process and related stochastic dynamics

TL;DR

This paper surveys the generalized fractional Poisson process (GFPP), a flexible renewal process involving Prabhakar functions, and develops related stochastic models including CTRWs and fractional diffusion equations, extending classical processes.

Contribution

It introduces the GFPP with two parameters, derives explicit formulas, and develops new stochastic dynamics models including fractional Kolmogorov equations and diffusion limits.

Findings

GFPP generalizes classical Poisson and Erlang processes.

Derived explicit formulas for GFPP and related stochastic models.

Established fractional diffusion equations from GFPP-based random walks.

Abstract

We survey the 'generalized fractional Poisson process' (GFPP). The GFPP is a renewal process generalizing Laskin's fractional Poisson counting process and was first introduced by Cahoy and Polito. The GFPP contains two index parameters with admissible ranges $0<\beta\leq 1$, $\alpha >0$ and a parameter characterizing the time scale. The GFPP involves Prabhakar generalized Mittag-Leffler functions and contains for special choices of the parameters the Laskin fractional Poisson process, the Erlang process and the standard Poisson process. We demonstrate this by means of explicit formulas. We develop the Montroll-Weiss continuous-time random walk (CTRW) for the GFPP on undirected networks which has Prabhakar distributed waiting times between the jumps of the walker. For this walk, we derive a generalized fractional Kolmogorov-Feller equation which involves Prabhakar generalized fractional operators governing the stochastic motions on the network. We analyze in $d$ dimensions the 'well-scaled' diffusion limit and obtain a fractional diffusion equation which is of the same type as for a walk with Mittag-Leffler distributed waiting times. The GFPP has the potential to capture various aspects in the dynamics of certain complex systems.

Equations

$$prob(\Delta t\le t)=\Psi (t)={\int }_0^t\chi (\tau )\mathrm{d}\tau ,t\ge 0,\underset{t\to \infty }{\lim }\Psi (t)=1-0,$$

$$prob(\Delta t>t)={\Phi }^{(0}(t)=1-\Psi (t)={\int }_t^{\infty }\chi (\tau )\mathrm{d}\tau .$$

$${\chi }^{(2)}(t)$$

$$={\int }_0^{\infty }{\int }_0^{\infty }\chi ({\tau }_1)\chi ({\tau }_2)\delta \left(t-({\tau }_1+{\tau }_2)\right)\mathrm{d}{\tau }_1\mathrm{d}{\tau }_2,$$

$${\chi }^{(n)}(t)={\int }_0^{\infty }\dots {\int }_0^{\infty }\chi ({\tau }_1)\dots \chi ({\tau }_n)\delta \left(t-\sum _{j=1}^n{\tau }_j\right)\mathrm{d}{\tau }_1\dots \mathrm{d}{\tau }_n.$$

$${\Phi }^{(n)}(t)={\int }_0^t\left(1-\Psi (t-\tau )\right){\chi }^{(n)}(\tau )\mathrm{d}\tau ={\int }_0^t{\Phi }^{(0)}(t-\tau ){\chi }^{(n)}(\tau )\mathrm{d}\tau ,$$

$$\sum _{n=0}^{\infty }{\Phi }^{(n)}(t)=1,$$

$$\begin{array}[]{l}\displaystyle\chi^{(n)}(t)=\int_{0}^{t}\chi^{(n-1)}(\tau)\chi(t-\tau){\rm d}\tau,\hskip 28.45274pt\chi^{(0)}(t)=\delta(t),\\ \\ \displaystyle\Phi^{(n)}(t)=\int_{0}^{t}\Phi^{(n-1)}(\tau)\chi(t-\tau){\rm d}\tau,\hskip 28.45274pt\Phi^{(0)}(t)=1-\Psi(t),\end{array}$$

$${\tilde{\chi }}^{(n)}(s)$$

$$={\int }_0^{\infty }\dots {\int }_0^{\infty }\chi ({\tau }_1)\dots \chi ({\tau }_n)\mathrm{d}{\tau }_1\dots \mathrm{d}{\tau }_n{e}^{-s{\sum }_{j=1}^n{\tau }_j}$$

$$={\left({\int }_0^{\infty }\chi (\tau )e^{-s\tau }\mathrm{d}\tau \right)}^n=(\tilde{\chi }(s){)}^n,$$

$${\tilde{\Phi }}^{(n)}(s)=\frac{1-\tilde{\chi }(s)}{s}(\tilde{\chi }(s){)}^n=\frac{1}{s}\left({\tilde{\chi }}^{(n)}(s)-{\tilde{\chi }}^{(n+1)}(s)\right),n\in {\mathbb{N}}_0,$$

$${\Phi }^{(n)}(t)$$

$$={\mathcal{L}}^{-1}\left\{\frac{1}{s}\left({\tilde{\chi }}^{(n)}(s)-{\tilde{\chi }}^{(n+1)}(s)\right)\right\},$$

$$G(t,v)=\sum _{n=0}^{\infty }{\Phi }^{(n)}(t)v^n,|v|\le 1,$$

$$\begin{array}[]{l}\displaystyle{\tilde{G}}(s,v)=\sum_{n=0}^{\infty}{\tilde{\Phi}}^{(n)}(s)v^{n}=\frac{1-{\tilde{\chi}}(s)}{s}\sum_{n=0}^{\infty}({\tilde{\chi}}(s))^{n}v^{n}\\ \displaystyle\hskip 36.98866pt=\frac{1-{\tilde{\chi}}(s)}{s}\frac{1}{1-v{\tilde{\chi}}(s)},\hskip 28.45274pt|v|\leq 1.\end{array}$$

$$\displaystyle{\bar{n}}(t)=\sum_{n=0}^{\infty}n\Phi^{(n)}(t)=\frac{d}{dv}G(t,v){\Big{|}}_{v=1}=\mathcal{L}^{-1}\left\{\frac{{\tilde{\chi}}(s)}{s(1-{\tilde{\chi}}(s))}\right\},$$

$$\left(\widetilde{\stackrel{ˉ}{n}}(s)+\frac{1}{s}\right)\widetilde{\chi (s)}=\frac{\chi (s)}{s(1-\tilde{\chi }(s))}=\widetilde{\stackrel{ˉ}{n}}(s),$$

$$\bar{n}(t)={\int }_0^t(1+\bar{n}(\tau ))\chi (t-\tau )\mathrm{d}\tau =1-{\Phi }^{(0)}(t)+{\int }_0^t\bar{n}(\tau )\chi (t-\tau )\mathrm{d}\tau ,$$

$$\bar{n}(t)=\sum _{n=0}^{\infty }(n+1){\Phi }^{(n+1)}(t)=\sum _{n=0}^{\infty }n{\Phi }^{(n+1)}+\sum _{n=0}^{\infty }{\Phi }^{(n+1)},$$

$${\left(\frac{d}{dt}+\sigma \right)}^{\beta }{g}_{\beta ,\sigma }(t)=\delta (t),\beta \ge 0\sigma >0.$$

$$g_{\beta ,\sigma }(t)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }\frac{e^{i\omega t}}{(\sigma +i\omega {)}^{\beta }}\mathrm{d}\omega ,$$

$$\frac{1}{(\sigma +i\omega {)}^{\beta }}={\int }_{-\infty }^{\infty }\Theta (\tau )e^{-(\sigma +i\omega )\tau }\frac{{\tau }^{\beta -1}}{(\beta -1)!}\mathrm{d}\tau .$$

$$g_{\beta ,\sigma }(t)={\left(\frac{d}{dt}+\sigma \right)}^{-\beta }\delta (t)=e^{-\sigma t}\Theta (t)\frac{t^{\beta -1}}{(\beta -1)!}=e^{-\sigma t}{g}_{\beta }(t),$$

$$g_{\beta }(t)$$

$$=\Theta (t)\frac{t^{\beta -1}}{\Gamma (\beta )}={\mathcal{L}}^{-1}\{s^{-\beta }\}=\frac{1}{2\pi i}{\int }_{\sigma -i\infty }^{\sigma +i\infty }{s}^{-\beta }{e}^{st}\mathrm{d}s$$

$$\underset{ϵ\to 0+}{\lim }$$

$$G_{\beta ,\xi }(t)$$

$$=\xi t^{\beta -1}\sum _{n=0}^{\infty }\frac{(-1{)}^n{\xi }^n{t}^{n\beta }}{\Gamma (\beta n+\beta )}=\xi t^{\beta -1}{E}_{\beta ,\beta }(-\xi t^{\beta }),t\in {\mathbb{R}}_{+}$$

$$=\frac{d}{dt}[1-E_{\beta }(-\xi t^{\beta })],0<\beta \le 1,$$

$$E_{\beta }(z)=\sum _{n=0}^{\infty }\frac{z^n}{\Gamma (\beta n+1)},\Re \{\beta \}>0,z\in \mathbb{C},$$

$$E_{\beta ,\gamma }(z)=\sum _{n=0}^{\infty }\frac{z^n}{\Gamma (\beta n+\gamma )},\Re \{\beta \}>0,\Re \{\gamma \}>0,z\in \mathbb{C},$$

$$G_{\beta ,\xi }(t)$$

$$=e^{\sigma t}\xi {\left\{\xi +{\left(\frac{d}{dt}+\sigma \right)}^{\beta }\right\}}^{-1}\delta (t)$$

Full text

MSC 2010: Primary 60K05, 33E12, 26A33; Secondary 60J60, 65R10, 60K40

Key Words and Phrases: Prabhakar fractional calculus; renewal process; Mittag-Leffler type functions; generalized fractional Poisson process; fractional ordinary and partial differential equations; generalized fractional Kolmogorov-Feller equation; continuous-time random walk; Non-Markovian random walks

1 Introduction

A characteristic feature in many complex systems are non-exponential time patterns with asymptotic power-law behavior as occurring for instance in anomalous transport and diffusion phenomena [19, 27, 36, 53, 58, 65] (and many others). On the other hand, it has been found that random walk approaches have the capacity to capture various dynamic features of certain complex systems. Markovian random walk models have been established to investigate a variety of aspects of stochastic processes in continuous spaces and networks as various as human mobility patterns [62] and random motions with long-range jumps on undirected networks with the emergence of Lévy flights [40, 41, 63, 64].

However, Markovian approaches are not able to capture (long-time) memory features. For the description of the stochastic behavior of such systems non-Markovian models of time fractional diffusion and fractional Fokker Planck equations have been developed exhibiting non-Markovian long-time memory effects [22, 35, 36, 65].

Continuous-time random walks (CTRWs) are stochastic processes where the walker makes jumps at random times which follow a renewal process. The CTRW approach was introduced by Montroll and Weiss in their 1965 seminal paper [44]. The notion of *Continuous Time Random Walk - CTRW * was coined by Scher and Lax [56] who launched in this way the field of ‘anomalous transport’ [58]. The CTRW mathematically is equivalent to a compound renewal process or a random walk subordinated to a renewal process [8, 19]. The CTRW meanwhile has been applied to a large variety of problems including charge transport in disordered systems [57], the stochastic processes in protein folding [5], the dynamics of chemical reactions [61], glass phenomenology [43], transport phenomena in low dimensional chaotic systems [27, 60] and generally in anomalous transport models [1, 26].

On the other hand, a huge literature exists [36] with various generalizations of the CTRW approach including variants leading to a fractional Fokker Planck equation [3] and the so-called ‘Aging Continuous Time Random Walk’ represents a variant where the stochastic properties of the walk depend on the ‘age’ of the walk [2]. The CTRW approach also was generalized to the case that waiting times of successive jumps are correlated by memory kernels [7].

The underlying renewal process in classical CTRWs is the Poisson process with exponential waiting time distributions and the Markov property [10, 36]. However, as mentioned it has turned out that many ‘complex systems’ exhibit time patterns with heavy-tailed power-law behavior rather than exponential waiting time characteristics [29, 36].

A non-Markovian fractional generalization of the Markovian Poisson process was introduced and analyzed by Laskin [29] who called this process ‘fractional Poisson process’. He demonstrated the utmost importance of the fractional Poisson distribution in several applications [30]. Various properties of the fractional Poisson process were meanwhile analyzed by several authors [4, 18, 24, 32, 34, 52] and consult also [47, 49] where space-time variants of the fractional Poisson process are introduced. The Laskin fractional Poisson process is a renewal process with Mittag-Leffler distributed waiting times exhibiting asymptotic heavy-tailed power law arrival densities with long-time memory effects as observed in a wide class of complex systems [36, 53, 65].

The heavy-tailed power-law asymptotic features of the Mittag-Leffler waiting time law in CTRWs and asymptotic universality of Mittag-Leffler pattern for large observation times among other properties were thoroughly analyzed in the references [4, 17, 19, 21].

The subject of the present survey paper is a generalization of Laskin’s fractional Poisson counting process and the analysis of the resulting stochastic dynamics on graphs. This generalization first was proposed by Cahoy and Polito in a seminal paper [6]. We call this process the ‘generalized fractional Poisson process’ (GFPP). Both the GFPP like the fractional Poisson process exhibit non-Markovian and long-memory effects.

The paper is organized as follows. In the first part (Section 2) we recall general concepts of renewal theory employed in the Continuous Time Random Walk (CTRW) of Montroll and Weiss [44]. Based on this approach we recall in Section 3 the Laskin fractional Poisson process and distribution [29, 53]. Section 4.1 is devoted to the GFPP counting process. We derive its waiting time distribution with related quantities and analyze their asymptotic behaviors. In Section 4.2 the probability distribution for $n$ arrivals within a time interval (‘state-probabilities’) is deduced for the GFPP. We refer this distribution to as the ‘generalized fractional Poisson distribution’ (GFPD). We derive in Section 4.3 the expected number of arrivals in a GFPP and analyze asymptotic power-law features. The GFPP contains two index parameters with admissible ranges $0<\beta\leq 1$ and $\alpha>0$ and a parameter characterizing the time scale. The Laskin fractional Poisson process and the standard Poisson process are contained as special cases in the GFPP.

As an application we consider a CTRW on undirected networks (Section 5). In this part we subordinate a random walk on the network to the GFPP: We derive (Section 5.2) the ‘generalized fractional Kolmogorov-Feller equation’ that governs the resulting stochastic motion. We show that for $\alpha=1$ the Laskin fractional Kolmogorov-Feller equations are recovered as well as for $\alpha=1$, $\beta=1$ the classical Kolmogorov-Feller equations of the standard Poisson process. Finally in Section 6 we analyze stochastic motions subordinated to a GFPP in the infinite $d$-dimensional integer lattice and derive in a ‘well-scaled diffusion limit’ the governing diffusion equation which is of the same fractional type as for the fractional Poisson process. For $\alpha=1$, $\beta=1$ this equation reduces to the standard diffusion equation (Fick’s second law) of normal diffusion.

The entire demonstration is accompanied by supplementary materials/Appendices where detailed mathematical derivations are performed and properties are discussed.

2 Renewal processes and waiting time distributions

2.1 Preliminaries

In order to analyze the GFPP counting process let us recall some basic properties of renewal processes and of the Montroll-Weiss continuous-time random walk (CTRW) approach. In the present paper, we deal with causal generalized functions and distributions [13]. A function or distribution is referred to as causal if it has the representation $\Theta(t)f(t)$ where $\Theta(t)$ is the Heaviside step function $\Theta(t)=1$ for $t\geq 0$ and $\Theta(t)=0$ for $t<0$. Throughout the paper we utilize Laplace transforms of causal functions and distributions which we define by Eq. (A.2) in Appendix A where we briefly outline the properties relevant in our analysis. We employ notation $\mathcal{L}\{f(t)\}=:{\tilde{f}}(s)$ for the Laplace transform of a causal function $\Theta(t)f(t)$.

Let us now evoke the basic features of renewal processes and waiting time distributions as a basis of the Montroll-Weiss continuous time random walk (CTRW) approach [44, 56] and consult also [18, 19, 20, 32]. We assume that a walker makes jumps at random times $0<t_{1}<t_{2}<\ldots t_{n}<\ldots\infty$ on a set of accessible sites. The walk can be conceived as a stream of jump events that happen at ‘arrival times’ $t_{k}$. We start the observation of the walk at time $t=0$ where the walker is sitting on his initial position until making the first jump at $t=t_{1}$ (arrival of the first jump event) and so forth. The waiting times $\Delta t_{k}=t_{k}-t_{k-1}>0$ between two successive (jump-) events are assumed to be drawn from the same probability density distribution $\chi(t)$. The waiting times $\Delta t_{k}$ then are independent and identically distributed (IID) random variables. $\chi(\tau){\rm d}\tau$ indicates the probability that the first jump event (first arrival) happens exactly at time $\tau$. We refer $\chi(t)$ to as waiting time probability density function (waiting time PDF) or short ‘jump density’. The probability that the walker has jumped at least once within the time interval $[0,t]$ is given by the cumulative probability distribution

[TABLE]

with the obvious initial condition $\Psi(t=0)=0$. (2.1) is also called ‘failure probability’. The probability that no jump event happened within $[0,t]$, i.e. that the waiting time is greater than $t$ is given by

[TABLE]

The cumulative probability distribution (2.2) that the walker during $[0,t]$ still is waiting on its initial site also is called ‘survival probability’. We notice that densities such as the jump density $\chi(t)$ have physical dimension $[sec]^{-1}$, whereas the cumulative distributions (2.1), (2.2) are dimensionless probabilities. The jump density $\chi(\tau)$ is non-negative and normalized (See (2.1)), i.e. for $t\rightarrow\infty$ the walker (almost) surely has made (at least) one jump. This means that survival probability (2.2) tends to zero as time approaches infinity $\lim_{t\rightarrow\infty}\Phi^{(0}(t)\rightarrow 0$.

What is the PDF that the walker makes its second jump exactly at time $t$? Since the waiting times $\Delta t_{1}$ and $\Delta t_{2}$ are IID the probability that the walker makes its first jump exactly at $\Delta t_{1}$ and the second jump at the time $\Delta t_{1}+\Delta t_{2}$ is ${\rm d}\chi^{(2)}=\chi(\Delta t_{1})\chi(\Delta t_{2}){\rm dt}_{1}{\rm d}t_{2}$. Integrating this quantity over the line $t=\Delta t_{1}+\Delta t_{2}$ (and dividing by a time increment) yields the PDF (probability per unit time) that the second step is made exactly at time $t$, namely

[TABLE]

where $\delta(u)$ indicates Dirac’s $\delta$-distribution. In this relation we used that $\chi(t)=\chi(t)\Theta(t)$ is as a causal function. The PDF $\chi^{(n)}(t)$ that the walker makes its $n$th jump exactly at time $t$ then is

[TABLE]

We refer (2.4) to as $n$-jump probability density. The physical dimension of the $n$-step PDF $\chi^{(n)}(t)$ (for all $n$) is $[sec]^{-1}$ whereas the cumulative probability $\Psi(t)$ of (2.1) is dimensionless. It is important to note that the PDF $\chi^{(n)}(t)$ to make the $n$th jump exactly at time $t$ is not the same as the dimensionless probability $\Phi^{(n)}(t)$ that the walker within the time interval $[0,t]$ has made $n$ jumps. In the latter case the walker may have made the $n$th step at a time $\tau<t$ and waited for $t-\tau$ (with probability $\Phi^{(0)}(t-\tau)=1-\Psi(t-\tau)$). The probability for $n$ arrivals within the time interval $[0,t]$ is hence given by

[TABLE]

for $n=0,1,2,\dots$. The probabilities (2.5) also are referred to as ‘state-probabilities’. This convolution covers also $n=0$ when we account for $\chi^{(0)}(t)=\delta(t)$. It follows that $\Phi^{(n)}(t)$ ($n=0,1,2,\ldots$) is a dimensionless probability distribution and normalized on the ‘state-space’ $n\in\mathbb{N}_{0}$, namely

[TABLE]

i.e., at any time $t\geq 0$ the walker for sure either waits, i.e. has performed $n=0$ jumps, or has performed $n=1,2,\ldots,\infty$ jumps. The normalization (2.6) is proved subsequently. It is straightforward to see that both (2.4) and (2.5) fulfill the recursion

[TABLE]

where $n=1$ recovers $\chi_{1}(t)=\chi(t)=\int_{0}^{t}\chi(t-\tau)\delta(\tau){\rm d}\tau$. In all relations we account for causality $\chi(t)=\Theta(t)\chi(t)$. The $n$-jump PDF $\chi^{(n)}(t)$ of (2.4) has the Laplace transform

[TABLE]

reflecting the product rule for convolutions where ${\tilde{\chi}}(s)$ indicates the Laplace transform of the (single-) jump density $\chi(t)=\chi^{(1)}(t)$. The property ${\tilde{\chi}}^{(n)}(s)|_{s=0}=\int_{0}^{\infty}\chi^{(n)}(t){\rm d}t=1$ in (2.8) reflects the normalization condition, i.e. for $t\rightarrow\infty$ the walker has made $n=0,1,2,\ldots,\infty$ jumps almost surely. We then get for the Laplace transform of (2.5)

[TABLE]

where for $n=0$ we have for the Laplace transform of the survival probability ${\tilde{\Phi}}^{(0)}(s)=\frac{1-{\tilde{\chi}}(s)}{s}$ when we account for ${\chi}_{0}(t)=\delta(t)$ with ${\tilde{\chi}}^{(0)}(s)=\mathcal{L}\left\{\delta(t)\right\}=1$. It is important to notice that $|{\tilde{\chi}}(s)|\leq{\chi(s=0)}=1$ (equality only for $s=0$ reflecting normalization). The following interpretation of Eq. (2.9) appears mention worthy. Taking into account the initial condition of $\Phi^{(n)}(0)=\delta_{n0}$, i.e. zero for $n\geq 1$ then it follows that

[TABLE]

where $\int_{0}^{t}\chi^{(m)}(\tau){\rm d}\tau$ is the probability for at least $m$ arrivals within $[0,t]$. It follows that the difference in (2.10) gives the probability $\Phi^{(n)}(t)$ for exactly $n$ arrivals within $[0,t]$.

It is convenient now to introduce a (dimensionless) generating function $G(t,v)$ which allows us to recover many quantities of interest in elegant manner. We define this generating function as

[TABLE]

where the value $G(t,1)=1$ reflects normalization condition (2.6). In the Laplace domain we get

[TABLE]

We observe that ${\tilde{G}}(s,v=1)=s^{-1}$ ($s\neq 0$) and thus proves $G(t,1)=\sum_{n=0}^{\infty}\Phi^{(n)}(t)=1$, i.e. normalization relation (2.6). It follows from $|{\tilde{\chi}}(s)|\leq{\tilde{\chi}}(s=0)=1$ that the series (2.12) is converging for $|v|\leq 1$.

A quantity of great interest in the following analysis is the expected number of arrivals ${\bar{n}}(t)$ (expected number of jumps) that happen in the time interval $[0,t]$. This quantity is obtained by

[TABLE]

${\bar{n}}(t)$ is a dimensionless function of $t$. In view of relation (2.13) the following observation appears mention worthy. Taking into account the identity

[TABLE]

where in the time domain this relation reads

[TABLE]

and is referred to as renewal equation [8, 18, 19]. The renewal equation is equivalently obtained by

[TABLE]

where we identify $\sum_{n=0}^{\infty}n\Phi^{(n+1)}=\int_{0}^{t}{\bar{n}}(\tau)\chi(t-\tau){\rm d}\tau$ and $\sum_{n=0}^{\infty}\Phi^{(n+1)}=1-\Phi^{0}(t)=\Psi(t)$ as the two contributions of Eq. (2.15).

3 Fractional Poisson distribution and process

In this section, we recall the fractional Poisson process introduced by Laskin who coined this name [29]. The fractional Poisson process also was derived in different manners by several authors [4, 18, 32, 52]. This process is a natural generalization of the standard Poisson process. The subsequent Section 4 is devoted to a generalization of Laskin’s fractional Poisson process which we call the ‘generalized fractional Poisson process’ (GFPP).

In order to capture effects of causality we employ the machinery of causal Green’s functions (causal propagators) and their relations to Laplace transforms. Let us recall some basic properties of this machinery. We introduce the causal propagator $g_{\beta}(t)=\Theta(t)g_{\beta}(t)$ by

[TABLE]

By using $\delta(t)=(2\pi)^{-1}\int_{-\infty}^{\infty}e^{i\omega t}{\rm d}\omega$ for the Dirac $\delta$-distribution, we obtain the Fourier transform of this propagator as ${\bar{g}}_{\beta,\sigma}(\omega)=\left(\sigma+i\omega\right)^{-\beta}$. This propagator is then defined by the Fourier integral

[TABLE]

where we observe that

[TABLE]

We employ throughout this paper as equivalent notations $\zeta!=\Gamma(\zeta+1)$. Plugging relation (3.3) into (3.2) yields

[TABLE]

where $\sigma>0$. With these equations we can write

[TABLE]

with $s=\sigma+i\omega$ and $\sigma=\Re(s)>0$ where the least line shows the equivalence to Laplace inversion.

We identify the causal propagator $g_{\beta}(t)=\lim_{\sigma\rightarrow 0+}g_{\beta,\sigma}(t)=\Theta(t)\frac{t^{\beta-1}}{\Gamma(\beta)}$ with the kernel of the Riemann-Liouville fractional integral of order $\beta$. This kernel includes the limit $\lim_{\beta\rightarrow 0+}g_{\beta}(t)=\delta(t)$. Hence $g_{\beta}(t)$ is the causal Green’s function that solves $(\frac{d}{dt})^{\beta}g_{\beta}(t)=\delta(t)$. Applying $(\frac{d}{dt})^{-\beta}$ on a causal function $f(t)=\Theta(t)f(t)$ hence gives the convolution

[TABLE]

which we identify with the Riemann-Liouville fractional integral of order $\beta$ [17, 18, 46]. In Appendix B we show that the causal propagator $(\frac{d}{dt})^{\beta}\delta(t)$ defines the Riemann-Liouville fractional derivative kernel. These relations shed light on the connection of fractional operators with causal propagators.

With these general remarks, let us now introduce the causal propagator (the causal waiting time PDF) that defines Laskin’s fractional Poisson process (See Appendix B for some properties):

[TABLE]

which is the Mittag-Leffler density containing the Mittag-Leffler function of order $\beta$ [17, 19, 22, 29]

[TABLE]

and the generalized (two-parameter) Mittag-Leffler function

[TABLE]

where $E_{\beta,1}(z)=E_{\beta}(z)$ recovers the Mittag-Leffler function (3.8) of order $\beta$. Since the argument $-\xi t^{\beta}$ of the Mittag-Leffler density (3.7) is dimensionless, $\xi$ is a positive characteristic parameter of physical dimension $[sec]^{-\beta}$. Let us now rederive the Mittag-Leffler density in the following way establishing the connection of causal propagators with Laplace transforms

[TABLE]

The convergence of the series requires $\sigma=\Re\{s\}>\xi^{1/\beta}$, where $s=\sigma+i\omega$ can be identified with the Laplace variable. We observe then

[TABLE]

where the causal propagators $G_{\beta,\sigma,\xi}(t)$, and $G_{\beta,\xi}(t)$ fulfill

[TABLE]

One can then introduce the fractional Poisson process as the renewal process with the Mittag-Leffler jump density (3.7) which reduces for $\beta=1$ to the exponential density $\chi_{1}(t)=\xi e^{-\xi t}$ of the standard Poisson process.

The survival probability in the fractional Poisson process is of the form of a Mittag-Leffler function

[TABLE]

where for $\beta=1$ the standard Poisson process with survival probability $E_{1}(-\xi t)=e^{-\xi t}$ is recovered. The fractional Poisson distribution is defined by the probabilities $\Phi_{\beta}^{(n)}(t)$ for $n$ arrivals within $[0,t]$ with jump density (3.7) and is obtained as [29] (See Appendix B for a derivation):

[TABLE]

where $0<\beta\leq 1$. For $\beta=1$ the fractional Poisson distribution recovers the standard Poisson distribution $\Phi^{(n)}_{1}(t)=\frac{(\xi t)^{n}}{n!}e^{-\xi t}$.

4 Generalization of the fractional Poisson process

4.1 Jump density and survival probability distribution

Here we survey the ‘Generalized Fractional Poisson Process’ (GFPP) first introduced by Cahoy and Polito [6]. The GFPP is a generalization of the Laskin fractional Poisson process. The jump density (waiting time PDF) of the GFPP has the Laplace transform [6, 37]

[TABLE]

The characteristic dimensional constant $\xi$ in (4.1) has physical dimension $[sec]^{-\beta}$. Per construction {\tilde{\chi}}_{\beta,\alpha}(s){\big{|}}_{s=0}=\int_{0}^{\infty}\chi_{\beta,\alpha}(t){\rm d}t=1 reflects normalization of the jump density. The GFPP recovers for $\alpha=1$ with $0<\beta\leq 1$ the fractional Poisson process, for $\beta=1$ with $\alpha>0$ an Erlang type process where we allow any $\alpha>0$, i.e. also non-integer $\alpha$ [18], and for $\alpha=1$, $\beta=1$ the standard Poisson process. For simplicity we call subsequently the process with $\beta=1$, $\alpha>0$ (for any $\alpha\in\mathbb{R}_{+}$) ‘Erlang process’.

We observe in (4.1) that the jump density in the fractional range $0<\beta<1$ and for all $\alpha>0$ has diverging mean (diverging expected time of first arrival), namely -\frac{d}{ds}{\tilde{\chi}}_{\beta,\alpha}(s){\big{|}}_{s=0}=\int_{0}^{\infty}t\chi_{\beta,\alpha}(t)\rightarrow+\infty reflecting occurrence of very long waiting times. In order to determine the time domain representation of (4.1) it is convenient to introduce the Pochhammer symbol which is defined as [33]

[TABLE]

For $c\neq 0$ (4.2) can be represented as $(c)_{m}=\frac{\Gamma(c+m)}{\Gamma(c)}$. For $c=0$ this expression has be regularized to define its value at $c=0$ with $(0)_{m}=\delta_{m0}=\lim_{c\rightarrow 0+}\frac{\Gamma(c+m)}{\Gamma(c)}$ which occurs as an important limiting case for $\lim_{\alpha\rightarrow 0+}\chi_{\beta,\alpha}(t)=\mathcal{L}^{-1}(1)=\delta(t)$. With this regularization definition (4.2) can be extended to all $c\in\mathbb{C}$. We can then evaluate the jump density defined by (4.1) by taking into account ($0<\beta\leq 1$, $\alpha>0$)

[TABLE]

In the same way as above (See Eq. (3.10)) we get the causal time domain representation of the jump density of the GFPP in terms of causal propagators in form [6, 37]

[TABLE]

For $\alpha=1$ (4.4) reduces to the Mittag-Leffler density (3.7), then for $\beta=1$ with $\alpha>0$ to the so called Erlang-density $\chi_{1,\alpha}(t)=\frac{\xi^{\alpha}t^{\alpha-1}}{\Gamma(\alpha)}e^{-\xi t}$ and for $\alpha=1$, $\beta=1$ to the density $\chi_{1,1}(t)=\xi e^{-\xi t}$ of the standard Poisson process. In expression (4.4) appears a generalization of the Mittag-Leffler function, the so called Prabhakar function which was introduced by Prabhakar [51] and has the representation

[TABLE]

The Prabhakar function converges for all $z\in{\mathbb{C}}$. The Prabhakar function is analyzed by several authors [11, 33, 59] and has recently attracted much interest. We call (4.4) ‘Prabhakar density’ or also ‘GFPP (jump-) density’. The GFPP (Prabhakar generalization of fractional Poisson counting process) and the related Prabhakar generalization of fractional calculus is a newly emerging rapidly developing field. For a comprehensive review of properties and applications we refer to the recent review article by Giusti et al. [15] and for an outline of Prabhakar theory and related generalized fractional calculus consult [12, 14, 50]. For an outline of Mittag-Leffler functions of several types and related topics we refer to the references [23, 25, 42, 48, 54]. Physical applications of the Prabhakar theory can be found in the references [9, 16, 38, 55] and an analysis of complete monotonicity of the Prabhakar function with a non-Debye relaxation model is developed in [31]. Discrete-time variants of Prabhakar renewal processes are developed in a recent paper [39]. We observe that for $t$ small the Prabhakar density (4.4) behaves as

[TABLE]

The power-law asymptotic behavior $\sim t^{\alpha\beta-1}$ occurring for $t$ small is different from that of the fractional Poisson process which is weakly singular $\sim t^{\beta-1}$ and is reproduced for $\alpha=1$. The Prabhakar density is weakly singular only for $\alpha\beta<1$ and non-singular else. The Prabhakar density (4.4) is plotted in Fig. 1 for different values $0<\beta<1$ and $\alpha=0.5$, $\alpha=1$, respectively. $\alpha=1$ in Fig. 1(b) recovers the Mittag-Leffler density (3.7). We see in this plot the behavior $\lim_{t\rightarrow 0}\chi_{\beta,\alpha}(t)\sim t^{\alpha\beta-1}\rightarrow+\infty$ occurring for $\alpha\beta<1$. The jump frequency in this regime is extremely high for short observation times where the jump density takes $\delta(t)$-shape in the limit $\alpha\beta\rightarrow 0+$. We can directly integrate (4.4) to obtain the cumulative distribution defined in Eq. (2.1), namely

[TABLE]

The probability distribution (4.7) (failure probability) of at least one arrival within $[0,t]$ in a GFPP is plotted for different values of $0<\beta<1$ in Fig. 2 for $\alpha=0.5$ and $\alpha=1$ (Mittag-Leffler case), respectively. For small $\xi^{\frac{1}{\beta}}t$ this distribution is dominated by the order $n=0$, namely

[TABLE]

where we also see in this relation that initial condition $\Psi_{\beta,\alpha}(0)=0$ is always fulfilled. The $t^{\alpha\beta}$ power-law behavior can be observed in the Fig. 2. The smaller $\alpha\beta\rightarrow 0$ the more $\Psi_{\beta,\alpha}(t)$ approaches Heaviside step function shape. The GFPP survival probability yields

[TABLE]

with initial condition $\Phi_{\beta,\alpha}^{(0)}(0)=1$ and $\lim_{t\rightarrow\infty}\Phi_{\beta,\alpha}^{(0)}(t)=0$. For $\alpha=1$ (4.9) becomes

[TABLE]

reproducing the Mittag-Leffler survival probability (3.13) in the fractional Poisson process. On the other hand let us consider $\beta=1$ and $\alpha>0$ which gives

[TABLE]

i.e. for $\beta=1$ the jump density is light-tailed with exponentially evanescent behavior for $t\rightarrow\infty$, but for $t\rightarrow 0$ we get power law behavior $\chi_{1,\alpha}(t)\sim t^{\alpha-1}$. The density (4.11) is an Erlang type density (extending the traditional Erlang density to any $\alpha\in\mathbb{R}_{+}$). For $\alpha=1$, $\beta=1$ (4.11) recovers the density $\chi_{1,1}(t)=\xi e^{-\xi t}$ of the standard Poisson process. Let us now consider the behavior for large times $t\rightarrow\infty$. Expanding (4.1) for $|s|$ small yields

[TABLE]

Thus we get fat-tailed behavior, namely

[TABLE]

where the exponent is independent of $\alpha$ and with $-\Gamma(-\beta)=\frac{\Gamma(1-\beta)}{\beta}>0$ for $0<\beta<1$.

4.2 Generalized fractional Poisson distribution

In this section we deduce the generalized counterpart $\Phi^{(n)}_{\beta,\alpha}(t)$ to the fractional Poisson distribution (3.14). The probability distribution of $n$ arrivals within $[0,t]$ is with (2.9) and (4.1) defined by

[TABLE]

where ${\tilde{\chi}}_{\beta,n\alpha}(s)=\frac{\xi^{n\alpha}}{(\xi+s^{\beta})^{n\alpha}}$ and with the cumulative probabilities $\Psi_{\beta,n\alpha}(t)$ (Eq. (4.7) with $\alpha\rightarrow n\alpha$). We obtain

[TABLE]

With this result we obtain from (4.14) the probability of $n$ arrivals within $[0,t]$ (GFPP state probabilities)

[TABLE]

($n=0,1,2,\ldots$) which is a function of the ‘dimensionless time’ $\xi^{\frac{1}{\beta}}t$ and contains the two index parameters $0<\beta\leq 1$ and $\alpha>0$. Expression (4.16) for the GFPP state probabilities also was obtained in [6]. We call the distribution $\Phi^{(n)}_{\beta,\alpha}(t)$ (4.16) generalized fractional Poisson distribution (GFPD). For $n=0$ (4.16) yields $\Phi^{(0)}_{\beta,\alpha}(t)=1-\Psi_{\beta,\alpha}(t)$ of Eq. (4.9). We can write the GFPD compactly as

[TABLE]

with

[TABLE]

for $m=1,2,\ldots$, $n=0,1,2,\ldots$. The GFPD is a dimensionless probability distribution where $\xi^{-\frac{1}{\beta}}$ having physical dimension $[sec]$ defines a characteristic time scale in the GFPP. It follows that for small dimensionless times the GFPD behaves as ($n=0,1,2,\ldots$)

[TABLE]

covering order $m=0$ of representations (4.16) and (4.17). It follows that the GFPP state probabilities fulfill the initial conditions

[TABLE]

i.e. at $t=0$ per construction of the CTRW the walker is on his departure site. Further of interest is the asymptotic behavior for large $t\xi^{\frac{1}{\beta}}$. To see this behavior we expand the Laplace transform (4.14) up to the lowest non-vanishing order in $\frac{s^{\beta}}{\xi}$ to arrive at

[TABLE]

where $0<\beta<1$, $\alpha>0$, $\forall n=0,1,\ldots,\infty$. This power law holds universally for all $\alpha>0$ with $0<\beta<1$ and contains the Laskin case $\alpha=1$ and is also independent of $n$. The universal power law scaling in (4.21) indeed can be attributed to non-Markovianity and long-time memory [39]. We will come back to this important issue in Sect. 4.3.

Now let us consider the important limit $\alpha=1$ for the GFPD in more details. Then the functions (4.18) become time independent coefficients, namely

[TABLE]

and for the oder $m=0$ we have \frac{1}{\Gamma(n\beta+1)}=\left(\frac{(n+m)!}{n!m!}\frac{1}{\Gamma(\beta(m+n)+1)}\right){\big{|}}_{m=0} thus we get for the GFPD (4.17) for $\alpha=1$ the expression

[TABLE]

which we identify with Laskin’s fractional Poisson distribution (3.14) [29].

4.3 Expected number of arrivals in a GFPP

Here we analyze the asymptotic behavior of the average number of arrivals ${\bar{n}}(t)$ within the interval of observation $[0,t]$. To this end we consider the generating function (2.12)

[TABLE]

where $|v|\leq 1$. The expected number ${\bar{n}}_{\beta,\alpha}(t)$ of steps the walker makes in the interval $[0,t]$ is {\bar{n}}_{\beta,\alpha}(t)=\frac{d}{dv}\mathcal{G}_{\beta,\alpha}(v,t){\big{|}}_{v=1}=\sum_{n=0}^{\infty}n\Phi^{(n)}_{\beta,\alpha}(t) and has the Laplace transform

[TABLE]

with ${\tilde{\bar{n}}}_{\beta,\alpha}(s)\approx\frac{\xi}{\alpha}s^{-\beta-1}$ for $|s|\rightarrow 0$. We obtain hence the asymptotic behavior for $t\xi^{\frac{1}{\beta}}\rightarrow\infty$ large as

[TABLE]

This result includes the Erlang case $\beta=1$ where for all $\alpha>0$ the average number of arrivals ${\bar{n}}_{1,\alpha}(t)\approx\frac{\xi}{\alpha}t$ increases linearly for $t\xi$ large. Now the following observation appears mention worthy. In view of the asymptotic power law behavior of the GFPD (4.21) which is independent of $n$ and the power law (4.26) we can define

[TABLE]

which gives for $\xi^{\frac{1}{\beta}}t\rightarrow\infty$ a universal dimensionless positive constant $\mathcal{C}_{\beta}<1$ which is independent of $\alpha$. In this relation we utilized the Euler reflection formula $\Gamma(1-\beta)=\frac{\pi}{\Gamma(\beta)\sin(\pi\beta)}$. On the other hand we obtain the asymptotic behavior for $t\xi^{\frac{1}{\beta}}\rightarrow 0$ when we expand ${\tilde{\bar{n}}}_{\beta,\alpha}(s)$ for $|s|\rightarrow\infty$ to arrive at

[TABLE]

Both asymptotic expressions (4.26) and (4.28) coincide for $\alpha=1$ recovering the (exact) expression (B.11) for the fractional Poisson process.

5 Montroll-Weiss CTRW and GFPP on undirected networks

5.1 CTRW on undirected networks

In this section, we analyze random motions on undirected networks (graphs) subordinated to a GFPP. In such a walk the jumps from one node to another take place at GFPP arrival times. To develop this model we employ the CTRW approach by Montroll and Weiss [44] and consult also [56]. For general outlines and historical developments of the CTRW we refer to the references [28, 58]. We first recall the Montroll-Weiss CTRW for undirected networks. Then in Section 5.2 we derive the generalization of the fractional Kolmogorov-Feller equation for the walk subordinated to a GFPP. Finally, in Section 6 we explore this motion on infinite $d$-dimensional integer lattices and analyze the ‘well-scaled diffusion limit’.

We consider an undirected connected network with $N$ nodes which we denote with $p=1,\ldots,N$. For our convenience, we employ Dirac’s $|bra\rangle\langle ket|$ notation. In an undirected network the positive-semidefinite $N\times N$ Laplacian matrix ${\mathbf{L}}=(L_{pq})$ characterizes the network topology and is defined by [40, 45]

[TABLE]

where ${\mathbf{A}}=(A_{pq})$ denotes the adjacency matrix having elements one if a pair of nodes is connected and zero otherwise. Further we do not allow self-connections thus $A_{pp}=0$. In an undirected network adjacency and Laplacian matrices are symmetric. The diagonal elements $L_{pp}=K_{p}$ of the Laplacian matrix are referred to as the degrees of the nodes $p$ counting the number of connected neighbor nodes with a node $p$ with $K_{p}=\sum_{q=1}^{N}A_{pq}$. The one-step transition matrix ${\mathbf{W}}=(W_{pq})$ relating the network topology with the probabilities for one jump in the random walk then is defined by [40, 45]

[TABLE]

which generally is a non-symmetric matrix for networks with variable degrees $K_{i}\neq K_{j}$ $(i\neq j)$. The one-step transition matrix $W_{pq}$ defines the conditional probability that the walker which is on node $p$ reaches in one step node $q$ [40, 45]. We consider now a CTRW where the walker performs random steps (‘jumps’) governed by the one-step transition matrix $(W_{pq})$ at random arrival times $0<t_{1}<t_{2}<\ldots t_{n}<\ldots\infty$ in a renewal process (See Section 2). We start the observation at $t=0$ where walker sits on its departure node.

We introduce the $N\times N$ transition matrix ${\mathbf{P}}(t)=(P_{ij}(t))$ indicating the (dimensionless) probability to find the walker at time $t$ on node $j$ under the condition that the walker was sitting at node $i$ at time $t=0$ (start of the observation). The transition matrix fulfills the (row-) normalization condition $\sum_{j=1}^{N}P_{ij}(t)=1$ and stochasticity implies that $0\leq P_{ij}(t)\leq 1$. We restrict us on connected undirected networks and allow variable degrees. In such networks the transition matrix is non-symmetric $P_{ij}(t)\neq P_{ji}(t)$ (See [40] for a detailed discussion). Assuming an initial condition ${\mathbf{P}}(0)$, then the probability to find the walker on node $j$ at time $t$ is given by the transition matrix $P_{ij}(t)$ which has the representation [8]

[TABLE]

where $\Phi^{(n)}(t)$ are the probabilities that the walker within the time interval $[0,t]$ has made $n$ steps (Eq. (2.5)) and ${\mathbf{W}}$ indicates the one-step transition matrix (5.2). The convergence of the series (5.3) can be easily proved by using that ${\mathbf{W}}$ has uniquely eigenvalues $|\lambda_{m}|\leq 1$ [40] together with the normalization of the state probabilities to see that $|\sum_{n=0}^{\infty}\Phi^{(n)}(t)(\lambda_{m})^{n}|\leq 1$. We observe the relation

[TABLE]

with elements

[TABLE]

This relation can be seen as the ‘CTRW on networks’- counterpart of the renewal equation (2.15) (For an outline of related aspects in renewal theory consult [8, 18, 19]). We assume now the initial condition that the walker at $t=0$ sits on a known node $P_{ij}(0)=\delta_{ij}$. Then the transition matrix (5.3) has the Laplace transform

[TABLE]

It is useful to account for the canonical representation of the one-step transition matrix [40]

[TABLE]

where we have a unique eigenvalue $\lambda_{1}=1$ and $|\lambda_{m}|<1$ ($m=2,\ldots,N$) (Consult [40] for an analysis of the spectral structure). In Eq. (5.6) we have introduced the right- and left eigenvectors $|\Phi_{j}\rangle$ and $\langle{\bar{\Phi}}_{j}|$ of ${\mathbf{W}}$, respectively. The first term $|\Phi_{1}\rangle\langle{\bar{\Phi}}_{1}|$ corresponding to $\lambda_{1}=1$ is the stationary distribution. Hence the canonical representation of (5.5) has the form

[TABLE]

with the eigenvalues

[TABLE]

We mention that for $\lambda_{1}=1$ we have ${\tilde{P}}(1,s)=s^{-1}$ thus in the time domain the stationary amplitude yields $P(1,t)=\mathcal{L}^{-1}(s^{-1})=\Theta(t)$. Equation (5.8) is the celebrated Montroll-Weiss formula [44, 56, 58] and is of utmost importance in a wide range of fields.

5.2 Generalized fractional Kolmogorov-Feller equation

In this section, we specify the renewal process in the CTRW to a GFPP and derive a generalization of the Laskin’s fractional Kolmogorov-Feller equation [29] describing the time evolution of the transition matrix. The GFPD probabilities $\Phi^{(n)}_{\beta,\alpha}(t)$ have with Eq. (4.1) the Laplace transforms

[TABLE]

For $n=0$ we have the Laplace transform of the survival probability

[TABLE]

These Laplace transforms fulfill the following relations

[TABLE]

and for $n=0$ we have

[TABLE]

Now the goal is to derive the generalized counterpart of the fractional Kolmogorov-Feller equation which is obtained with the causal time domain representations of relations (5.11) and (5.12). To this end we introduce the convolution kernel

[TABLE]

In the causal time domain Eq. (5.11) takes the form

[TABLE]

where ${}_{0}\!\mathcal{D}_{t}^{\beta,\alpha}\cdot\Phi(t)=:\int_{0}^{t}\Phi(\tau)\mathcal{D}^{\beta,\alpha}(t-\tau){\rm d}\tau$ has to be read as a convolution. The time domain representation of this kernel is obtained explicitly as (See Appendix C for a derivation, especially Eqs. (C.11)-(C.14))

[TABLE]

and

[TABLE]

where $E_{a,b,c}(z)$ is a Prabhakar Mittag-Leffler type function defined in Eq. (5.23). We emphasize that the $\Theta(t)$-function must be taken ‘under the time derivative’ of integer order $\lceil\alpha\beta\rceil$ where the ceiling function $\lceil\gamma\rceil$ indicates the smallest integer greater or equal to $\gamma$. Above distributional representations for $\mathcal{D}^{\beta,\alpha}(t)$ define then convolution operators in the sense of relation (C.1) and related properties of Laplace transforms are outlined in Appendix A (See property (A.3) together with decomposition (C.3)). For $n=0$ from Eq. (5.12) we get

[TABLE]

The kernel

[TABLE]

is evaluated explicitly in Appendix D (Eq. (D.2)) where $K^{(0)}_{\beta,\alpha}(t)$ is obtained as

[TABLE]

and we have $\mathcal{D}^{\beta,\alpha}(t)=\frac{d}{dt}K^{(0)}_{\beta,\alpha}(t)$. In all these relations the expressions holding for $\alpha\beta\in\mathbb{N}$ are recovered by the expressions for $\alpha\beta\notin\mathbb{N}$ in the limit $\alpha\beta\rightarrow\lceil\alpha\beta\rceil-0$. These kernels contain the Prabhakar type function

[TABLE]

where $E_{u,v}^{w}(z)$ is defined in (4.5). The transition probability matrix (5.3) of the GFPP walk then writes

[TABLE]

Then we obtain with (5.14), (5.17) the equation

[TABLE]

which takes in the causal time domain the form ($0<\beta\leq 1$, $\alpha>0$)

[TABLE]

We call equation (5.26) ‘generalized fractional Kolmogorov-Feller equation’. This equation can be rewritten as

[TABLE]

to give

[TABLE]

This representation of the generalized fractional Kolmogorov-Feller equation emphasizes the connection with Prabhakar fractional calculus [9, 14, 15, 50]. Taking into account the explicit representations (5.15) and (5.22) this equation writes for $\alpha\beta\notin{\mathbb{N}}$

[TABLE]

Laskin fractional Poisson process limit $\alpha=1$:

Let us consider (5.26) with (5.29) for the Laskin fractional case $\alpha=1$ and $0<\beta<1$. Then we have (See (5.22))

[TABLE]

thus Eq. (5.28) takes the form

[TABLE]

The left-hand side of this equation is the Caputo-fractional derivative of order $\beta$. This equation can also be written as (corresponding to (5.29))

[TABLE]

where ${}_{0}\!D_{t}^{\beta}$ is the Riemann-Liouville fractional derivative of order $\beta$. (5.31) and (5.32) are equivalent representations of the fractional Kolmogorov-Feller equation derived by Laskin and others [4, 18, 19, 24, 29, 32, 36, 53, 65]. By accounting for $\lim_{\beta\rightarrow 1-0}\frac{t^{-\beta}}{\Gamma(1-\beta)}=\delta(t)$ in (5.31)) one can see that the fractional K-F equation for $\alpha=1$ and $\beta=1$ reduces to the classical K-F equation for the walk subordinated to the standard Poisson process. In this way we have shown that the generalized fractional K-F equation contains for $\alpha=1$ the well-known classical counterparts.

6 Generalized fractional diffusion in the infinite $d$-dimensional integer lattice

In this section, we analyze the effects of the GFPP on diffusive features in the infinite $d$-dimensional integer lattice ${\mathbb{Z}}^{d}$. Especially we derive here the ‘well-scaled’ diffusion limit. We denote with $\vec{q}=(q_{1},..,q_{d})$ and $q_{j}=0,\pm 1,\pm 2,\ldots\in{\mathbb{Z}_{0}}$ the lattice points (nodes of the network). The canonical representation of the one-step transition matrix (5.6) becomes then a Fourier integral over the $d$-dimensional Brillouin zone $k_{j}\in[-\pi,\pi]$ and has the eigenvalues [40, 41]

[TABLE]

where $\mu(\vec{k})$ indicates the eigenvalues of the Laplacian matrix (5.1) of the lattice. We assumed that any node is connected only with its $2d$ closest neighbor nodes. The Laplace transform (5.7) of the transition matrix then writes

[TABLE]

We assume here the initial condition that the walk starts in the origin $\vec{q}=\vec{0}$. We subordinate a normal random walk on the infinite $d$-dimensional integer lattice to a GFPP. In order to analyze diffusive motions we account for the behavior of the eigenvalues (6.1) for $k\rightarrow 0$ (where $k=|\vec{k}|$), namely

[TABLE]

Then we get for the Montroll-Weiss equation (5.8) for $k\rightarrow 0$ small

[TABLE]

This relation contains the Laplace transform ${\tilde{\bar{n}}}_{\beta,\alpha}(s)$ of the expected number of arrivals (4.25). Before we further exploit this relation let us first derive the diffusion equation which is emerging in a ‘well-scaled diffusion limit’. To this end we rewrite the exact Montroll-Weiss equation (6.4) as follows

[TABLE]

We notice that the matrix representation of this equation gives a generalized fractional K-F equation of the type (5.26). Eq. (6.5) takes for small $k$ the representation

[TABLE]

where we utilized that the Laplacian eigenvalues $\mu(\vec{k})\approx k^{2}$ for $k$ small (See (6.1)). Eq. (6.6) holds asymptotically for small (dimensionless wave numbers) $k={\bar{k}}h$ where we introduce a new wave vector ${\vec{\bar{k}}}=({\bar{k}}_{1},..,{\bar{k}}_{d})$ ($k_{j}\in\mathbb{R}$ having physical dimension $[cm]^{-1}$). Now we consider the ‘well-scaled diffusion limit’ of (6.6). Generally there are different ways (not necessarily well-scaled procedures) to define the diffusion limit leading to different types of diffusion equations [19, 37, 38, 40]. Let us introduce the Laplace operator $\Delta=\sum_{j=1}^{d}\frac{\partial^{2}}{\partial x_{j}^{2}}$ of the $\mathbb{R}^{d}$ with respect to the new rescaled quasi-continuous spatial coordinates $\vec{x}=h\vec{p}\in\mathbb{R}^{d}$ (having units $[cm]$). Then we introduce $P_{(\beta,\alpha)}(\vec{q},t)=h^{d}\mathcal{P}_{(\beta,\alpha)}(\vec{x},t)$ where $\mathcal{P}_{(\beta,\alpha)}(\vec{x},t)$ indicates the spatial transition probability kernel having physical units $[cm]^{-d}$. For any finite ${\bar{k}}$ we can choose $h>0$ sufficiently small such that $\mu(\vec{k})\approx k^{2}=({\bar{k}}h)^{2}\rightarrow 0$ thus the left-hand side of relation (6.6) tends to zero for $h\rightarrow 0$. To maintain the right-hand side in relation (6.6) also small requires that $\xi(h)\rightarrow\infty$ for $h\rightarrow 0$ thus we can expand in equation (6.6) $(1+\frac{1}{\xi}s^{\beta})^{\alpha}\approx 1+\frac{\alpha}{\xi}s^{\beta}$ leading to relation

[TABLE]

In order the left-hand side of (6.7) remains finite for $h\rightarrow 0$ requires the scaling $\xi(h)\sim\xi_{0}h^{-2}\rightarrow\infty$ (where $\xi_{0}>0$ is a new constant independent of $h$).

The constant $\frac{\xi_{0}}{2d\alpha}$ has physical dimension $[cm^{2}sec^{-\beta}]$ and can be interpreted as a generalized fractional diffusion constant (recovering for $\beta=1$ the units $[cm^{2}sec^{-1}]$ of normal diffusion). We observe that the index $\alpha$ enters Eq. (6.7) only as a scaling parameter. The well-scaled diffusion limiting equation (6.7) then leads to the fractional diffusion equation

[TABLE]

where $\mathcal{P}_{\beta,\alpha}(\vec{x},t)=\delta^{(d)}(\vec{x})$ is the initial condition in the form of a Dirac’s $\delta$-function in $\mathbb{R}^{d}$ and $\vec{x}=h\vec{q}\in\mathbb{R}^{d}$ are the rescaled quasi-continuous coordinates. ${}_{0}\!D_{t}^{\beta}\cdot P_{\beta,\alpha}(q,t)$ indicates the Riemann-Liouville fractional derivative of order $\beta$ (See Appendix C). The right-hand side of (6.8) can be read as Caputo-fractional derivative of order $\beta$ of $\mathcal{P}_{\beta,\alpha}(\vec{x},t)$ where the term $\,\frac{t^{-\beta}}{\Gamma(1-\beta)}$ reflects the non-Markovian feature with long-time memory of the walk (For a discussion consult [37, 38, 39]).

We notice that the well-scaled diffusive limit of generalized fractional diffusion leads to the same type of purely fractional diffusion equation (6.8) as for the fractional Poisson process ($\alpha=1$ with $0<\beta<1$) having Mittag-Leffler waiting times. (6.8) hence reflects ‘asymptotic Mittag-Leffler universality’ [21]. Indeed Eq. (6.8) coincides with the fractional diffusion equation widely used in fractional dynamics [35, 65] (and many others).

Finally we consider the standard Poisson limit $\alpha=1$ and $\beta=1$. Then we get with (6.7) Fick’s second law of normal diffusion

[TABLE]

where $\frac{\xi_{0}}{2d}$ indicates the diffusion constant of standard (normal-) diffusion. This equation is also obtained by the limit $\beta\rightarrow 1-0$ in Eq (6.8) taking into account $\lim_{\beta\rightarrow 1-0}\frac{t^{-\beta}}{\Gamma(1-\beta)}=\delta(t)$. This behavior reflects the Markovian memoryless nature of underlying the pure Poisson process.

Finally let us consider the mean squared displacement in generalized fractional diffusion. Then we obtain the Laplace-Fourier transform of the mean squared displacement ${\tilde{\sigma}^{2}}(s)$ by taking into account (6.4) for $\xi(h)\sim\xi_{0}h^{-2}$ large

[TABLE]

where $\Delta_{{\vec{\bar{k}}}}$ stands for the Laplace operator with respect to the new wave vector $\vec{{\bar{k}}}\in\mathbb{R}^{d}$. For the GFPP the asymptotic expression for the expected number of arrivals ${\bar{n}}_{\beta,\alpha}(t)$ was determined in (4.26). We hence obtain for the mean square displacement

[TABLE]

We have for $0<\beta<1$ and for all $\alpha>0$ including fractional Poisson $\alpha=1$ sublinear $t^{\beta}$ power-law behavior corresponding to fat-tailed jump density (4.13). The sublinear power law $\sigma_{\beta,\alpha}^{2}(t)\sim t^{\beta}$ ($0<\beta<1$) corresponds also to the universal limit for $t\xi^{\frac{1}{\beta}}\rightarrow\infty$ of fat-tailed waiting time PDFs where extremely long waiting times occur. Such behavior is well known in the literature for anomalous subdiffusion ($0<\beta<1$) [36].

In contrast for $\beta=1$ for all $\alpha>0$ (Erlang regime which includes the standard Poisson $\alpha=1$, $\beta=1$) we obtain normal diffusive behavior with linear increase $\sigma^{2}_{\beta,\alpha}(t)=\frac{\xi_{0}}{\alpha}t$ of the mean squared displacement. This behavior can be attributed to normal diffusion [36]. The emergence of universal scaling laws for $t\xi^{\frac{1}{\beta}}\rightarrow\infty$ indeed reflects the asymptotic Mittag-Leffler universality for time fractional dynamics [21].

7 Conclusions

We have analyzed a generalization of the Laskin fractional Poisson process, the generalized fractional Poisson process (GFPP) and derived the generalized fractional Poisson distribution (GFPD) (Eq.(4.16)). The GFPP is (unless in the standard Poisson case $\alpha=1$ and $\beta=1$) non-Markovian and introduces long-time memory effects. For $\alpha=1$, $0<\beta<1$ the fractional Poisson process, for $\beta=1$, $\alpha>0$ the Erlang process, and for $\alpha=1$, $\beta=1$ the (memoryless) standard Poisson process are recovered. We showed that for $\alpha=1$ the GFPD state probabilities reduce to Laskin’s fractional Poisson distribution (3.14) and for $\alpha=1$ with $\beta=1$ to the standard Poisson distribution. The GFPP contains two index parameters $\alpha>0$ and $0<\beta\leq 1$ and a characteristic time scale $\xi^{-\frac{1}{\beta}}$ controlling the dynamic behavior. For long dimensionless observation times the asymptotic universal scaling of the fractional Poisson process emerges in the range $0<\beta<1$, $\alpha>0$ (asymptotic Mittag-Leffler universality).

We developed the Montroll-Weiss CTRW with GFPP (Prabhakar distributed) waiting times between the jumps and analyzed the resulting stochastic motions in undirected networks. We derived the ‘generalized fractional Kolmogorov-Feller equation’ (Eqs. (5.26)-(5.28)) which governs the time-evolution of the transition matrix in undirected networks.

As an application we analyzed for the $d$-dimensional infinite integer lattice the ‘well-scaled diffusion limit’ and the resulting stochastic motion. We obtained the same type of fractional diffusion equation as for CTRWs with Mittag-Leffler distributed waiting times in a fractional Poisson process. All equations turn in the limits $\alpha=1$, $0<\beta<1$, and $\beta=1$, $\alpha=1$ into their classical counterparts of Laskin’s fractional Poisson process and standard Poisson process, respectively.

The GFPP has a huge potential of further applications in anomalous diffusion and transport phenomena including turbulence, non-Markovian stochastic motions, and in the dynamics of complex systems.

Appendix A Laplace transforms of causal functions

In this appendix, we discuss some properties related with causal distributions and generalized functions [13] and their Laplace transforms. First we introduce the Heaviside step function

[TABLE]

where we emphasize that $\Theta(0)=1$ and we have $\delta(t)=\frac{d}{dt}\Theta(t)$ for Dirac’s $\delta$-distribution. We consider uniquely causal functions and distributions $f(t)=\Theta(t)f(t)$ which may take non-zero values only for non-negative times $t\geq 0$. First of all we introduce the Laplace transform ${\tilde{f}}(s)$ of a causal function $f(t)=\Theta(t)f(t)$ by

[TABLE]

with $\sigma=\Re\{s\}>\sigma_{0}$. We choose the lower integration limit infinitesimally negative $0-$ thus $\Theta(0-)=\lim_{\epsilon\rightarrow 0+}\Theta(-\epsilon)=0$. We denote by $\mathcal{L}^{-1}\left\{{\tilde{f}}(s)\right\}=f(t)$ the inverse Laplace transform. As an important issue in our analysis let us compare the two Laplace transforms (i) $\mathcal{L}\{\frac{d^{m}}{dt^{m}}\left[\Theta(t)f(t)\right]\}$ and (ii) $\mathcal{L}\{\Theta(t)\frac{d^{m}}{dt^{m}}f(t)\}$ where $m\in\mathbb{N}_{0}$ is integer or zero. One obtains straightforwardly (i) by $m$ partial integrations

[TABLE]

On the other hand, we have for (ii)

[TABLE]

In the analysis of this paper we often suppress writing the Heaviside $\Theta(t)$-function in causal functions in cases where no derivatives are involved. Further we mention that we utilize both notations $\Gamma(\beta)=(\beta-1)!$ for the $\Gamma$-function.

Appendix B Some properties of the fractional Poisson process

Now we are interested in the state probabilities (2.5), i.e. the probabilities for $n$ arrivals within interval $[0,t]$ for the fractional Poisson process with Mittag-Leffler density (3.7). This probability is given by

[TABLE]

where the Laplace transform is with ${\tilde{\chi}}_{\beta}^{(n)}(s)=\chi_{\beta}^{n}(s)=\frac{\xi^{n}}{(\xi+s^{\beta})^{n}}$,

[TABLE]

For evaluation it is convenient to employ the generating function of relation (2.11)

[TABLE]

and its Laplace transform

[TABLE]

where (B.2) can be expressed by \frac{1}{n!}\frac{d^{n}}{dv^{n}}{\tilde{G}}_{\beta}(s,v){\big{|}}_{v=0}. Again from ${\tilde{G}}_{\beta}(s,1)=s^{-1}$ follows that the $\Phi^{(n)}_{\beta}(t)$ are a normalized probability distribution on $n=\{0,1,2,\ldots\infty\}$, i.e. G_{\beta}(t,v){\Big{|}}_{v=1}=\sum_{n=0}^{\infty}\Phi^{(n)}_{\beta}(t)=1. It follows from (3.13) that Laplace inversion of (B.4) yields again a Mittag-Leffler function where we have to replace in $\xi\rightarrow\xi(1-v)$ in (3.13). In this way we obtain for the generating function (B.3) the Mittag-Leffler function

[TABLE]

This expression was also obtained by Laskin [29]. The probabilities $\Phi^{(n)}_{\beta}(t)$ are then obtained by

[TABLE]

This relation yields the fractional Poisson distribution, i.e. the probability for $n$ arrivals within $[0,t]$ in a fractional Poisson process ($0<\beta\leq 1$) [29]

[TABLE]

where $n=0$ gives the Mittag-Leffler survival probability $\Phi^{(0)}(t)=E_{\beta}(-\xi t^{\beta})$. For $\beta=1$ relation (B.7) recovers the standard Poisson distribution

[TABLE]

with the generating function (B.5)

[TABLE]

The normalization of the fractional Poisson distribution is easily verified by

[TABLE]

The expected number of arrivals (See Eq. (2.13)) in a fractional Poisson process within time interval $[0,t]$ is obtained as

[TABLE]

Appendix C Evaluation of the kernel $\mathcal{L}^{-1}\left\{(s^{\beta}+\xi)^{\alpha}\right\}$

Let us make a brief general remark on notations. Convolution kernels that we write in the distributional representation $k(t)=\frac{d^{m}}{dt^{m}}\left[\Theta(t)f(t)\right]$, $m\in\mathbb{N}_{0}$ (having Laplace transform ${\tilde{k}}(s)=s^{m}{\tilde{f}}(s)$), where the Heaviside step function is included into the differentiation define convolutions as follows

[TABLE]

On the other hand if $k(t)$ appears as function and is not integrated we can read it in on $\mathbb{R}_{+}$ as $k(t)=\frac{d^{m}}{dt^{m}}\left[\Theta(t)f(t)\right]=\frac{d^{m}}{dt^{m}}f(t)$.

Now let us evaluate the time domain representation of the generalized fractional operator of Eq. (5.13) defined by

[TABLE]

The following simple representation will turn out to be highly useful to evaluate kernel (C.2), namely

[TABLE]

where in the evaluation the important property (A.3) will be employed. We have to consider two possible cases: (i) $\alpha\beta$ is not integer $\alpha\beta\notin\mathbb{N}$ and (ii) $\alpha\beta$ is integer $\alpha\beta\in\mathbb{N}$.

Let us first evaluate case (i) $\alpha\beta\notin\mathbb{N}$. To this end it is convenient to introduce the ceiling function $\lceil(..)\rceil$ where $m=\lceil\lambda\rceil$ indicates the smallest integer larger or equal to $\lambda$. In this way by putting $m=\lceil\alpha\beta\rceil$ in relation (C.3) we see that in the part $s^{\alpha\beta-\lceil\alpha\beta\rceil}(1+\xi s^{-\beta})^{\alpha}$ of the decomposition occur terms $s^{-\gamma}$ with negative powers $-\gamma<0$ leading to inverse Laplace transforms $\mathcal{L}^{-1}\{s^{-\gamma}\}=\Theta(t)\frac{t^{\gamma-1}}{\Gamma(\gamma)}$ giving (up to coefficients) Riemann-Liouville fractional integral kernels where $\gamma-1\geq\lceil\alpha\beta\rceil-\alpha\beta-1>-1$, i.e. the terms are either weakly singular or non-singular and hence in both cases integrable on $\mathbb{R}_{+}$. With these observations let us now derive the causal time domain representation of the kernel (C.2). We can write

[TABLE]

taking the form of a Fourier integral ($m=\lceil\alpha\beta\rceil$)

[TABLE]

Then for $\sigma>\xi^{\frac{1}{\beta}}$ we can expand (C.5) with respect to $(\sigma+i\omega)^{-\beta}\xi$ and integrate each term to arrive at

[TABLE]

where we notice that $-1<{\lceil\alpha\beta\rceil}-\beta\alpha-1<0$. This relation contains the Prabhakar type function (5.23) which converges in the entire $z$-plane. Let us consider order $n=0$ in (C.6) which can be evaluated in the same way as above ($\gamma=\alpha\beta$) where $\gamma>0$ with $\gamma\notin\mathbb{N}$, namely

[TABLE]

We identify this expression with the kernel of the Riemann-Liouville fractional derivative of order $\gamma$ (e.g. [42] and many others) which defines a convolution ${}_{0}\!D^{\gamma}_{t}\cdot f(t)$ in the sense of (C.1), namely

[TABLE]

One observes that integer-order derivatives are also covered

[TABLE]

where (C.7) then takes $\lim_{\gamma\rightarrow\lceil\gamma\rceil-0}\frac{d^{\lceil\gamma\rceil}}{dt^{\lceil\gamma\rceil}}\delta(t)$ ($\lceil\gamma\rceil=\gamma$).

Then to complete our demonstration let us finally evaluate also the part

[TABLE]

We notice that in this series the order $n=0$ yields Dirac’s $\delta$-distribution $\delta(t)$ and $E_{\alpha,\beta,1}(\xi t^{\beta})$ is a Prabhakar type function defined in Eq. (5.23). From above considerations (See Eq. (C.3)) it follows that we can represent the kernel (C.6) by the fractional Riemann-Liouville derivative of order $\alpha\beta$ as

[TABLE]

where

[TABLE]

We notice that ${\mathcal{B}}_{\beta,\alpha}(t)$ is a (weakly singular) ‘Prabhakar kernel’ [We refer to the recent paper of Giusti [14] for an introduction of ‘Prabhakar generalized fractional calculus’] where ${\mathcal{B}}_{\beta,\alpha}(t)=e^{-\alpha}_{\beta,\lceil\alpha\beta\rceil-\alpha\beta}(-\xi,t)$ (in his notation).

Now let us evaluate case (ii) where $\lceil\alpha\beta\rceil=\alpha\beta\in\mathbb{N}$ is integer. We obtain then

[TABLE]

Since for $\alpha\beta\in\mathbb{N}$ we have $\lceil\alpha\beta\rceil=\alpha\beta$ thus this expression coincides with (C.10). Then we get

[TABLE]

This result is consistent with (C.11) and is recovered from this expression in the limiting case $\alpha\beta\rightarrow\lceil\alpha\beta\rceil-0$ when we account for the limit of the fractional Riemann-Liouville derivative (C.9) reducing then to an integer order derivative.

Appendix D Evaluation of

$\mathcal{L}^{-1}\left\{\frac{(s^{\beta}+\xi)^{\alpha}-\xi^{\alpha}}{s}\right\}$

Here we evaluate the kernel

[TABLE]

that occurs in (5.17) which is now straightforward accounting for the results in previous Appendix C. Hence we arrive at

[TABLE]

For $\alpha\beta<1$ we have $\lceil\alpha\beta\rceil-1=0$ thus in that case $K^{0}_{\beta,\alpha}(t)={\mathcal{B}}_{\beta,\alpha}(t)$. We observe in these relations that $\mathcal{D}^{\beta,\alpha}(t)=\frac{d}{dt}K^{0}_{\beta,\alpha}(t)$.