Floquet Theory of Classical Relaxation in Time-Dependent Field

Abstract

The anomalous transport of particles in the presence of a time-dependent field is considered in the framework of a comb model. This turbulent-like dynamics consists of inhomogeneous time-dependent advection along the x-backbone and Brownian motion along the y-side branches. This geometrically constrained transport leads to anomalous diffusion along the backbone, which is described by a fractional diffusion equation with time-dependent coefficients. The time periodic process leads to localization of the transport and a particular form of relaxation. The analytical approach is considered in the framework of the Floquet theory, which is developed for the fractional diffusion equation with periodic in time coefficients. This physical situation is considered in detail and analytical expressions for both the probability density function and the mean squared displacement are obtained. The new analytical approach is developed in the framework of the fractional Floquet theory that makes it possible to investigate a new class of anomalous diffusion in the presence of time periodic fields.

1. Introduction

Diffusion in the presence of time-dependent fields is a long lasting issue, starting from turbulent diffusion [1,2], and it attracts much attention in contemporary studies; see, e.g., [3,4,5,6]. Among these studies, fractional relaxation [7] is one of the central issues. A typical example is anomalous diffusion in MRI processes [8,9].

In the paper, we suggest another example of anomalous transport of particles in the presence of a time periodic field and multiplicative noise. For the description of the phenomenon, we consider a two-dimensional Langevin equation in the Matheron de Marsily form [10],

(1)${\dot{x}}_1=v\left(x_2\right)cos\left(\omega t\right)x_1,{\dot{x}}_2=\xi \left(t\right),v\left(x_2\right)=v\delta \left(x_2\right).$

Here $\xi \left(t\right)$ is a random Gaussian $\delta$-correlated process (white noise) with the correlation function $⟨\xi \left(t\right)\xi \left(t^{\prime }\right)⟩=2D\delta (t-t^{\prime })$, where D is a diffusion coefficient and $⟨\cdots ⟩$ means functional averaging over the Gaussian process, and $\delta (t-t^{\prime })$ is the Dirac delta function. This specific form of Equation (1) corresponds exactly to the two-dimensional comb geometry, where the transport along the $x_1$ direction is possible at $x_2=0$ only, while Brownian motion takes place along the $x_2$ side branches (or fingers). When the oscillating frequency $\omega =0$, the comb model has been discussed in Refs. [11,12] as a geometrical realization of turbulent diffusion like geometric Brownian motion [13]; see also an extended discussion in Ref. [14]. In that case, an exponential in time spread $exp\left[{\int }_0^{t}v\left(x_2\left(\tau \right)\right)d\tau \right]$ takes place along the $x_1$ axis, while $x_2\left(t\right)=B\left(t\right)={\int }_0^t\xi \left(\tau \right)d\tau$ is the Brownian functional, which describes Brownian motion along the $x_2$ axis.

This topologically constrained dynamics of a particle in the presence of time-dependent fields is a new issue of anomalous diffusion inside the inhomogeneous non-stationary environment. So far, such a task was not treatable analytically. Here we present an exact analytical method, which makes it possible to treat the system in the framework of the Floquet theory. The Floquet theory describes linear differential equations with time periodic coefficients [15] and it is widely used in driven quantum systems, including the interaction of radiation with matter, quantum nonlinear resonances, quantum chaos, and so on [16,17,18,19]. The Floquet theorem states that the solution for this class of systems has the form $f\left(t\right)=e^{-iϵt}u\left(t\right)$, where $u(t+T)=u\left(t\right)$ is periodic in time with the period $T=2\pi /\omega$, while $ϵ$ is a quasienergy spectrum. The Floquet theory for the fractional time Schrödinger equation has been developed in Ref. [20]. Here we generalize this approach for fractional diffusion equations by considering the realistic system in the form of the topologically constrained Langevin Equation (1).

The paper is organized as follows. In Section 2, the time-dependent comb model is introduced. The corresponding fractional diffusion equation is obtained in Laplace space. Section 3 is devoted to the application of the fractional Floquet theorem to the time fractional diffusion equation. Analytical expressions for the backbone, marginal, and comb probability density functions are obtained in Section 4. A process of relaxation and the mean squared displacement of the backbone transport are discussed in Section 5. Discussion and conclusion are presented in Section 6. Some details of calculations are presented in Appendix A.

2. Time-Dependent Comb Model

In the present case, $\omega \ne 0$, the geometrical realization of the transport, is analogous to turbulent diffusion, discussed in Ref. [12]; however, the transport itself is more sophisticated and will be studied in this section. Introducing a probability density function (PDF)

$P(x,y,t)=⟨\delta (x_1\left(t\right)-x)\delta (x_2\left(t\right)-y)⟩,$

we consider $x_1$ as the x axis and $x_2$ as the y axis. Then following either the Furutsu–Novikov formula [21,22] or straightforward calculations of the functional integral [23], we obtain the Fokker–Planck equation in the form of a two-dimensional comb model,

(2)${\partial }_{t}P(x,y,t)=-v\delta \left(y\right)cos\left(\omega t\right){\partial }_{x}xP(x,y,t)+D{\partial }_y^{2}P(x,y,t).$

The specific property of the comb model is that the random spread along the x axis (backbone) is possible only at $y=0$ in the form of randomly inhomogeneous advection, while Brownian motion in the y axis is a homogeneous process. This Fokker–Planck equation conserves the probability flow, $\int P(x,y,t)dxdy=1$, which follows from the natural (zero) boundary conditions at infinity for the PDF and its first derivatives with respect to x and y. Symmetrical initial conditions $P_0(x,y)$ will be specified in sequel.

Following the standard prescription according to Refs. [24,25,26], used for $\omega =0$, the solution to Equation (2) is considered in the form of the inverse Laplace transform,

(3)$P(x,y,t)={\mathcal{L}}^{-1}\left[\tilde{P}(x,y,s)\right]=\frac{1}{2\pi i}{\int }_{c-i\infty }^{c+i\infty }exp\left[-\left|y\right|\sqrt{s/D}\right]\tilde{f}(x,s)e^{st}ds,$

where $\tilde{P}(x,y,s)=\mathcal{L}\left[P(x,y,t)\right]\equiv {\int }_0^{\infty }P(x,y,t)e^{-st}dt$ is the Laplace image of the PDF, while the Laplace image $\tilde{f}(x,s)$ determines the backbone PDF at $y=0$, namely $f(x,t)=P(x,0,t)$, which immediately follows from Equation (3) when the number of particles on the backbone is not conserved. Another characteristic of the backbone transport with the conservation of the number of diffusive particles is the marginal PDF $P_1(x,t)$, which is

(4)$P_1(x,t)={\int }_{-\infty }^{\infty }P(x,y,t)dy=\frac{2\sqrt{D}}{2\pi i}{\int }_{-i\infty }^{i\infty }{s}^{-\frac{1}{2}}\tilde{f}(x,s)e^{st}ds,$

and it relates to the backbone PDF in Laplace space as follows:

(5)${\tilde{P}}_1(x,s)=\mathcal{L}\left[P_1(x,t)\right]=2\sqrt{D/s}\tilde{f}(x,s).$

In the section, we consider the anomalous transport in terms of the backbone PDF $f(x,t)$.

Performing the Laplace transform of the comb Equation (2), we obtain

(6)$s\tilde{P}(x,y,s)-P_0(x,y)=-v\delta \left(y\right){\partial }_{x}x\mathcal{L}\left\{cos\left(\omega t\right){\mathcal{L}}^{-1}\left[\tilde{P}(x,y,s)\right]\right\}+D{\partial }_y^2\tilde{P}(x,y,s),$

where $P_0(x,y)=f_0\left(x\right)\delta \left(y\right)$. Taking into account the second derivative of the Laplace image in Equation (3) with respect to y, we obtain

(7)$D{\partial }_y^2\tilde{P}(x,y,s)=s\tilde{P}(x,y,s)-2\sqrt{Ds}\delta \left(y\right)\tilde{P}(x,y,s).$

Then taking into account Equations (3) and (7), we obtain Equation (6) for the Laplace image of the backbone PDF $\tilde{f}(x,s)$ as follows:

(8)$2\sqrt{sD}\tilde{f}(x,s)-f_0\left(x\right)=-v{\partial }_{x}x\mathcal{L}\left[cos\left(\omega t\right)f(x,t)\right],$

where we also take into account that

$\delta \left(y\right){\mathcal{L}}^{-1}\left[\tilde{P}(x,y,s)\right]={\mathcal{L}}^{-1}\left[\tilde{P}(x,0,s)\right]\delta \left(y\right)=f(x,t)\delta \left(y\right).$

3. Fractional Floquet Theorem

Equation (8) is the fractional diffusion equation in Laplace space. It can be easily verified by the Laplace inversion that it corresponds to the time fractional diffusion equation, where ${\mathcal{L}}^{-1}\left[s^{\frac{1}{2}}\tilde{f}\left(s\right)\right]\left(t\right)$ is the Riemann–Liouville fractional derivative of the order of $\frac{1}{2}$, ref. [27] (in sequel; see also Equation (16)). Therefore, the solution to Equation (8) can be presented in the form of the fractional Floquet theorem [20], which is

(9)$f(x,t)=\sum _{n=-\infty }^{\infty }{C}_n\left(x\right){\mathrm{\Lambda }}_n\left(t\right)e^{in\omega t},$

where ${\mathrm{\Lambda }}_n\left(t\right)$ is a function which will be determined, and the initial condition is $f_0\left(x\right)=\sum {\mathrm{\Lambda }}_n\left(0\right)C_n\left(x\right)$. Therefore, its Laplace image reads

(10)$\tilde{f}(x,s)=\sum _{n=-\infty }^{\infty }{C}_n\left(x\right)\tilde{\mathrm{\Lambda }}(s-in\omega ).$

The Laplace transform of the r.h.s. of Equation (8) yields the chain of transformations

(11)$\begin{array}{c}\mathcal{L}\left[cos\left(\omega t\right)f(x,t)\right]=\frac{1}{2}{\displaystyle \sum _n}{C}_n\left(x\right)\mathcal{L}\left[{\mathrm{\Lambda }}_n\left(t\right)\left(e^{i(n+1)\omega t}+e^{i(n-1)\omega t}\right)\right]\hfill \\ \hspace{1em}\hspace{1em}=\frac{1}{2}{\displaystyle \sum _n}{C}_n\left(x\right)\left[\tilde{\mathrm{\Lambda }}\left(s-i(n+1)\omega \right)+\tilde{\mathrm{\Lambda }}\left(s-i(n-1)\omega \right)\right]\\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\frac{1}{2}{\displaystyle \sum _n}\tilde{\mathrm{\Lambda }}(s-in\omega )\left[C_{n-1}\left(x\right)+C_{n+1}\left(x\right)\right],\end{array}$

where the shift of the indexes is performed in the last line.

In the next step, we consider the term $\frac{1}{2}\left[C_{n-1}\left(x\right)+C_{n+1}\left(x\right)\right]$, which is discussed in Appendix A. According to Equation (A4), it corresponds to the equation

(12)$-\frac{v}{2}{\partial }_{x}x\left[C_{n-1}\left(x\right)+C_{n+1}\left(x\right)\right]=(i\omega D/v)C_n\left(x\right).$

This expression yields Equation (8) for the Laplace image in the closed form, which is

(13)$\sum _n\left[\tilde{\mathrm{\Lambda }}(s-in\omega )\left(s^{\frac{1}{2}}-in\overline{\omega }\right)-{\left(4D\right)}^{-\frac{1}{2}}\mathrm{\Lambda }\left(0\right)\right]C_n\left(x\right)=0,\overline{\omega }\equiv \omega \frac{D^{\frac{1}{2}}}{2v}.$

In general case, the expression in the square brackets does not equal zero; however, in this specific case, it is reasonable to suppose that it does. Then, setting the brackets equal to zero and performing the Laplace inversion, we obtain for $\nu =\frac{1}{2}$

(14)$\begin{array}{c}2\sqrt{D}{\mathrm{\Lambda }}_n\left(t\right)e^{in\omega t}=\mathrm{\Lambda }\left(0\right)\frac{1}{2\pi i}{\int }_{c-i\infty }^{c+i\infty }\frac{e^{st}ds}{s^{\nu }-in\overline{\omega }}=\mathrm{\Lambda }\left(0\right)t^{\nu -1}{E}_{\nu ,\nu }\left(in\overline{\omega }{t}^{\nu }\right)\hfill \\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\mathrm{\Lambda }\left(0\right)t^{\nu -1}{\displaystyle \sum _{k=0}^{\infty }}\frac{{\left(in\overline{\omega }{t}^{\nu }\right)}^k}{\mathrm{\Gamma }(k\nu +\nu )},\nu =\frac{1}{2},\end{array}$

where $E_{\alpha ,\beta }\left(z\right)$ is the two-parameter Mittag-Leffler function (V. 3 [28]). Therefore, the backbone PDF is

(15)$f(x,t)=\sum _n{C}_n\left(x\right){\mathrm{\Lambda }}_n\left(t\right)e^{in\omega t}=\mathrm{\Lambda }\left(0\right){\left(4D\right)}^{-\frac{1}{2}}{t}^{\nu -1}\sum _n{C}_n\left(x\right)E_{\nu ,\nu }\left(in\overline{\omega }{t}^{\nu }\right).$

In the limit $\nu =1$, the PDF reduces to the standard result, as expected according to the Floquet theorem, $f(x,t)=\sum C_n{e}^{in\omega t}$; see Appendix A.

3.1. Verification of the Solution

Since the solution (15) is obtained by assumption on the square brackets in Equation (13), it should be verified by substitution in the equation in the task. To that end, we first perform the Laplace inversion of Equation (8) and then substitute the solution in the l.h.s. of the obtained equation. Without restriction of generality, we set (for a while) $\mathrm{\Lambda }\left(0\right)=1$. Then we have the following chain of transformations for $t>0$:

(16)$\begin{array}{c}\frac{1}{2\pi i}{\int }_{c-i\infty }^{c+i\infty }{s}^{\nu }\tilde{f}(x,s)e^{st}ds-{\left(4D\right)}^{\frac{1}{2}}{f}_0\left(x\right)\delta \left(t\right)=D_t{\int }_{c-i\infty }^{c+i\infty }{s}^{\nu -1}\tilde{f}(x,s)e^{st}ds\hfill \\ \hspace{1em}\hspace{1em}={\displaystyle \sum _n}{C}_n\left(x\right)D_t{\int }_{c-i\infty }^{c+i\infty }\frac{s^{\nu -1}{e}^{st}}{s^{\nu }-in\overline{\omega }}ds={\displaystyle \sum _n}{C}_n\left(x\right)D_t{E}_{\nu ,1}\left(in\overline{\omega }{t}^{\nu }\right)\\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}={\displaystyle \sum _n}{C}_n\left(x\right)D_t{\displaystyle \sum _k}\frac{{\left(in\overline{\omega }{t}^{\nu }\right)}^k}{\mathrm{\Gamma }(k\nu +1)}=t^{\nu -1}{\displaystyle \sum _n}in\overline{\omega }{C}_n\left(x\right)E_{\nu ,\nu }\left(in\overline{\omega }{t}^{\nu }\right),\end{array}$

where $D_t\equiv \frac{d}{dt}$ and $\overline{\omega }=\omega \frac{D^{\frac{1}{2}}}{2v}$.

Let us consider the r.h.s. of Equation (8). Performing the identity transformation, we have

(17)$\begin{array}{c}{\mathcal{L}}^{-1}\left\{\left(-\frac{v}{2D^{\frac{1}{2}}}\right){\partial }_{x}x\mathcal{L}\left[cos\left(\omega t\right)f(x,t)\right]\right\}\hfill \\ \hspace{1em}=\left(-\frac{v}{2D^{\frac{1}{2}}}\right){\partial }_{x}x{\displaystyle \sum _n}{C}_n\left(x\right){\mathcal{L}}^{-1}\left\{\mathcal{L}\left[t^{\nu -1}{E}_{\nu ,\nu }\left(in\overline{\omega }{t}^{\nu }\right)\right]\right\}\\ \hspace{1em}=\left(-\frac{v}{2D^{\frac{1}{2}}}\right){\partial }_{x}x{\displaystyle \sum _n}{C}_n\left(x\right)\frac{1}{2}{\mathcal{L}}^{-1}\left[\tilde{\mathrm{\Lambda }}\left(s-i(n+1)\omega \right)+\tilde{\mathrm{\Lambda }}\left(s-i(n-1)\omega \right)\right]\\ \hspace{1em}={\displaystyle \sum _n}{\mathcal{L}}^{-1}\tilde{\mathrm{\Lambda }}(s-in\omega )\left(-\frac{v}{2D^{\frac{1}{2}}}\right){\partial }_{x}x\frac{1}{2}\left[C_{n-1}\left(x\right)+C_{n+1}\left(x\right)\right]\\ \hfill \hspace{1em}={\displaystyle \sum _n}in\overline{\omega }{C}_n\left(x\right)t^{\nu -1}{E}_{\nu ,\nu }\left(in\overline{\omega }{t}^{\nu }\right),\end{array}$

which coincides exactly with the result in Equation (16). Therefore, Equation (15) is the solution of Equation (8) at the condition (12), which determines the coefficients $C_n\left(x\right)$ for $x\in {\mathbb{R}}_{+}$. Then the next step of the analysis is the solution of that equation.

3.2. Space Dependence of the PDF

Let us consider Equation (12), which results from the standard Floquet theorem (see Appendix A), and it reads

(18)$-\overline{v}{\partial }_{x}x\left[C_{n+1}\left(x\right)+C_{n-1}\left(x\right)\right]=2in\omega C_n\left(x\right),\overline{v}=\frac{v^2}{D}.$

Performing the variable change $x=x_0{e}^y$, and for $x\in {\mathbb{R}}_{+}$, $y=ln(x/x_0)$, where $y\in \mathbb{R}$ (not to be confused with the finger’s coordinate y), we have ${\partial }_x=\frac{1}{x}{\partial }_y$ and ${\partial }_{x}x=1+{\partial }_y$, and $C_n\left(x\right)\to C_n\left(y\right)$. Then performing the Fourier transform ${\mathcal{C}}_n\left(\xi \right)=\mathcal{F}\left[C_n\left(y\right)\right]$, we obtain Equation (18) as follows:

(19)$\frac{\overline{v}(-\xi +i)}{\omega }\left[{\mathcal{C}}_{n+1}\left(\xi \right)+{\mathcal{C}}_{n-1}\left(\xi \right)\right]=2n{\mathcal{C}}_n\left(\xi \right),$

which is the recurrence relation of the Bessel functions $J_n\left(z\right)$, where $z=\frac{\overline{v}(\xi +i)}{\omega }$ [29]. Performing the Fourier inversion, we have

(20)$C_n\left(y\right)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{J}_n\left(\frac{\overline{v}(-\xi +i)}{\omega }\right)e^{iy\xi }d\xi .$

3.2.1. Initial Condition

The result can be related to the initial condition as well. To this end, we use the integral representation of the Bessel function [29],

(21)$J_n\left(\frac{\overline{v}(-\xi +i)}{\omega }\right)\equiv J_n\left(z\right)=\frac{1}{2\pi }{\int }_{-\pi }^{\pi }{e}^{izsin\eta -in\eta }d\eta .$

The solution (9) for $t=0$ reads $f(x,0)\equiv f_0\left(x\right)={\sum }_n\mathrm{\Lambda }\left(0\right)C_n\left(x\right)$, where $\mathrm{\Lambda }\left(0\right)={\mathrm{\Lambda }}_n(t=0)$ is independent of n. Then we have from Equations (20) and (21)

(22)$\begin{array}{c}{f}_0\left(x\right)/\mathrm{\Lambda }\left(0\right)={\displaystyle \sum _n}{C}_n\left(x\right)\to {\displaystyle \sum _n}{C}_n\left(y\right)\hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\frac{1}{{\left(2\pi \right)}^2}{\int }_{-\infty }^{\infty }{e}^{iy\xi }d\xi {\int }_{-\pi }^{\pi }exp\left[\frac{i\overline{v}(-\xi +i)}{\omega }sin\eta \right]{\displaystyle \sum _n}{e}^{in\eta }d\eta \\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{e}^{iy\xi }d\xi =\delta \left(y\right)=\delta \left(ln(x/x_0)\right)=x_0\delta (x-x_0).\end{array}$

Now taking $\mathrm{\Lambda }\left(0\right)=\frac{1}{x_0}$, we obtain $f_0\left(x\right)=\delta (x-x_0)$. Note also that

$\sum _{n=-\infty }^{\infty }{e}^{in\eta }=\sum _m\delta (\eta -2\pi m)=\delta \left(\eta \right),\eta \in [-\pi ,\pi ].$

3.2.2. Fourier Inversion

To obtain the coefficients $C_n\left(y\right)$ we perform the Fourier inversion in Equation (20) using the integral representation (21) of the Bessel function. Then we have

(23)$\begin{array}{c}{C}_n\left(x\right)\to C_n\left(y\right)=\frac{1}{{\left(2\pi \right)}^2}{\int }_{-\infty }^{\infty }{\int }_{-\pi }^{\pi }{e}^{i\frac{\overline{v}(-\xi +i)}{\omega }sin\eta -in\eta }d\eta e^{iy\xi }d\xi \hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\frac{1}{{\left(2\pi \right)}^2}{\int }_{-\pi }^{\pi }d\eta e^{-in\eta }{e}^{-\frac{\overline{v}}{\omega }sin\eta }{\int }_{-\infty }^{\infty }{e}^{-i\xi \frac{\overline{v}}{\omega }sin\eta }{e}^{iy\xi }d\xi \\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\frac{1}{\left(2\pi \right)}{\int }_{-\pi }^{\pi }d\eta e^{-in\eta }{e}^{-\frac{\overline{v}}{\omega }sin\eta }\delta \left(-\frac{\overline{v}}{\omega }sin\eta +y\right).\end{array}$

Taking into account the property of the Dirac $\delta$ function

(24a)$\begin{array}{cc}\hfill & \delta \left(\frac{\overline{v}}{\omega }sin\eta -y\right)=\frac{\omega }{\overline{v}}{\left[cos\left({\eta }_0\right)\right]}^{-1}\delta (\eta -{\eta }_0),\hfill \end{array}$

(24b)$\begin{array}{cc}\hfill & {\eta }_0=arcsin(\omega y/\overline{v})=arcsin\left[(\omega /\overline{v})ln(x/x_0)\right],\hfill \end{array}$

where $\left|\omega y/\overline{v}\right|\le 1$, we obtain

(25)$C_n\left(x\right)=\frac{x_0}{x}{\left[{\left(\frac{\overline{v}}{\omega }\right)}^2-{ln}^2\left(\frac{x}{x_0}\right)\right]}^{-\frac{1}{2}}exp\left[inarcsin\left(\frac{\omega }{\overline{v}}ln(x/x_0)\right)\right].$

Here we also used the fact that $y=ln(x/x_0)$. Then the restriction condition is $\left|ln(x/x_0)\right|\le \overline{v}/\omega =v^2/D\omega$. Eventually, we obtain the dimensionless parameter $v^2/D\omega$, which defines the diffusive area on the backbone: $x_0{e}^{-v^2/D\omega }\le x\le x_0{e}^{v^2/D\omega }$. This also leads to the complete correspondence of dimensions of the parameters of the system.

4. Probability Density Functions

All ingredients for the construction of the backbone and marginal PDFs are obtained. First we consider the backbone PDF.

4.1. Backbone PDF

Substituting Equation (25) into expression (15) for the backbone PDF and taking into account that $\mathrm{\Lambda }\left(0\right)=1/x_0$, we have

(26)$f(x,t)=\frac{{\left(4D\right)}^{-\frac{1}{2}}{t}^{\nu -1}}{x{\left[{\left(\frac{\overline{v}}{\omega }\right)}^2-{ln}^2\left(\frac{x}{x_0}\right)\right]}^{\frac{1}{2}}}\sum _{n}exp\left[inarcsin\left(\frac{\omega }{\overline{v}}ln(x/x_0)\right)\right]E_{\nu ,\nu }\left(in\overline{\omega }{t}^{\nu }\right).$

It is also convenient to use the concise form of Equation (26), which is

(27)$f(x,t)={\left(4Dt\right)}^{-\frac{1}{2}}{\eta }_0^{\prime }\left(x\right)\sum _n{e}^{in{\eta }_0\left(x\right)}{E}_{\nu ,\nu }\left(in\overline{\omega }{t}^{\frac{1}{2}}\right),$

where we use $\nu =1/2$ and ${\eta }_0^{\prime }\left(x\right)\equiv \frac{d{\eta }_0\left(x\right)}{dx}$, and ${\eta }_0\left(x\right)$ is defined in Equation (24b).

It is reasonable to use this expression to describe the relaxation process on the backbone of the comb. To this end, we present the Mittag-Leffler functions as follows [28]:

(28)$E_{\alpha ,\beta }\left(z\right)=\frac{1}{\mathrm{\Gamma }\left(\beta \right)}+z{E}_{\alpha ,\beta +\alpha }\left(z\right),$

which yields $E_{\nu ,\nu }\left(z\right)=\frac{1}{\mathrm{\Gamma }\left(\nu \right)}+z{E}_{\nu ,2\nu }\left(z\right)$. If z limits to zero faster than $n^{1/\nu }$, then the second term can be neglected in the initial time limit. For the first term, we obtain

(29)$f(x,t\to 0)\approx \frac{{\left(4Dt\right)}^{-\frac{1}{2}}}{\mathrm{\Gamma }(1/2)}\delta (x-x_0),$

where details of the calculation are presented in Appendix A.2.

This expression describes the relaxation process on the backbone of the comb. However, this information is not complete, since the backbone PDF (26) describes relaxation accomplished also by loosing the number of particles on the backbone (leaking into the side branches). This immediately follows from the integration

$\overline{f}\left(t\right)={\int }_0^{\infty }f(x,t)dx=\frac{{\left(4Dt\right)}^{-\frac{1}{2}}}{\mathrm{\Gamma }(1/2)}.$

Therefore, for the complete description of the backbone relaxation, one should use the normalization condition [24,30], $f(x,t)\to f(x,t)/\overline{f}\left(t\right)$, which also yields the correct limit for the initial condition $f_0\left(x\right)=\delta (x-x_0)$. However, a more reasonable consideration is in the terms of the marginal PDF according to Equation (4), or the comb PDF (3).

4.2. Marginal PDF

Performing the Laplace transform of the backbone PDF (26) and substituting it into Equation (4), we obtain

(30)$\begin{array}{c}{P}_1(x,t)=\frac{1}{x{\left[{\left(\frac{\overline{v}}{\omega }\right)}^2-{ln}^2\left(\frac{x}{x_0}\right)\right]}^{\frac{1}{2}}}{\displaystyle \sum _n}exp\left[inarcsin\left(\frac{\omega }{\overline{v}}ln(x/x_0)\right)\right]\frac{1}{2\pi i}{\int }_{-i\infty }^{i\infty }\frac{s^{-\frac{1}{2}}{e}^{st}ds}{s^{\frac{1}{2}}-in\overline{\omega }}\hfill \\ \hfill \hspace{1em}=\frac{1}{x{\left[{\left(\frac{\overline{v}}{\omega }\right)}^2-{ln}^2\left(\frac{x}{x_0}\right)\right]}^{\frac{1}{2}}}{\displaystyle \sum _n}exp\left[inarcsin\left(\frac{\omega }{\overline{v}}ln(x/x_0)\right)\right]E_{\frac{1}{2}}\left(in\overline{\omega }{t}^{\frac{1}{2}}\right).\end{array}$

Now the number of particles on the backbone is conserved. Integrating Equation (30) with respect to x, we obtain

(31)$\begin{array}{c}{\int }_0^{\infty }\frac{dx/x}{{\left[{\left(\frac{\overline{v}}{\omega }\right)}^2-{ln}^2\left(\frac{x}{x_0}\right)\right]}^{\frac{1}{2}}}{\displaystyle \sum _n}exp\left[inarcsin\left(\frac{\omega }{\overline{v}}ln(x/x_0)\right)\right]E_{\frac{1}{2}}\left(in\overline{\omega }{t}^{\frac{1}{2}}\right)\hfill \\ \hspace{1em}\hspace{1em}={\int }_{-\infty }^{\infty }\frac{d\left(ln(x/x_0)\right)}{{\left[{\left(\frac{\overline{v}}{\omega }\right)}^2-{ln}^2\left(\frac{x}{x_0}\right)\right]}^{\frac{1}{2}}}{\displaystyle \sum _n}exp\left[inarcsin\left(\frac{\omega }{\overline{v}}ln(x/x_0)\right)\right]E_{\frac{1}{2}}\left(in\overline{\omega }{t}^{\frac{1}{2}}\right)\\ \hfill \hspace{1em}\hspace{1em}={\int }_{-\pi }^{\pi }d{\eta }_0{\displaystyle \sum _n}{e}^{in{\eta }_0}{E}_{\frac{1}{2}}\left(in\overline{\omega }{t}^{\frac{1}{2}}\right)={\displaystyle \sum _n}{\delta }_{n,0}{E}_{\frac{1}{2}}\left(in\overline{\omega }{t}^{\frac{1}{2}}\right)=1.\end{array}$

Note that the upper limit of the integration taken at infinity is formal, since we know that the region of the integration (diffusion area) is restricted for both x and y, and the exact limits are known for ${\eta }_0\left(x\right)\in [-\pi ,\pi ]$ to which we have eventually arrived along the change of the variables.

From Equation (31), the marginal PDF can be written in the compact form as follows:

(32)$P_1(x,t)=\frac{d{\eta }_0\left(x\right)}{dx}\sum _n{e}^{in{\eta }_0}{E}_{\frac{1}{2}}\left(in\overline{\omega }{t}^{\frac{1}{2}}\right)$

(cf. Equation (27)). Then, using the definition (A11a), it is convenient to present the Mittag-Leffler function for the real argument, that is,

(33)$\begin{array}{c}{P}_1(x,t)=\frac{d{\eta }_0\left(x\right)}{dx}{\displaystyle \sum _n}{e}^{in{\eta }_0}\left[{\displaystyle \sum _{k=0}^{\infty }}\frac{{(-n^2{\overline{\omega }}^{2}t)}^k}{\mathrm{\Gamma }\left(2k\frac{1}{2}+1\right)}+in\overline{\omega }{t}^{\frac{1}{2}}{\displaystyle \sum _{k=0}^{\infty }}\frac{{(-n^2{\overline{\omega }}^{2}t)}^k}{\mathrm{\Gamma }\left(2k\frac{1}{2}+\frac{3}{2}\right)}\right]\hfill \\ \hspace{1em}=\frac{d{\eta }_0\left(x\right)}{dx}\left[1+{\displaystyle \sum _n}{}^{\prime }{e}^{in{\eta }_0}{E}_{1,1}\left(-n^2{\overline{\omega }}^{2}t\right)+i\overline{\omega }{t}^{\frac{1}{2}}{\displaystyle \sum _n}{}^{\prime }{e}^{in{\eta }_0}n{E}_{1,\frac{3}{2}}\left(-n^2{\overline{\omega }}^{2}t\right)\right]\\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\frac{d{\eta }_0\left(x\right)}{dx}\left[1+2{\displaystyle \sum _{n=1}^{\infty }}cos\left(n{\eta }_0\right)e^{-n^2{\overline{\omega }}^{2}t}+i\overline{\omega }{t}^{\frac{1}{2}}{\displaystyle \sum _n}{}^{\prime }{e}^{in{\eta }_0}n{E}_{1,\frac{3}{2}}\left(-n^2{\overline{\omega }}^{2}t\right)\right].\end{array}$

This is the first central result of the paper, since the main properties of the fractional transport are described by the marginal PDF $P_1(x,t)$ and this situation will be discussed in Section 5. Here we admit the normalization condition of the asymptotic large time behavior. The asymptotic behavior of the second Mittag-Leffler function in Equation (33) can be presented as follows (see Equation (A11b)):

(34)$E_{1,\frac{3}{2}}(-n^2{\overline{\omega }}^{2}t)\approx \frac{{\left(n^2{\overline{\omega }}^{2}t\right)}^{-1}}{\mathrm{\Gamma }(1/2)},n^2{\overline{\omega }}^{2}t\to \infty .$

Note that $\mathrm{\Gamma }(1/2)=\sqrt{\pi }$. Taking into account the power law decay with time in Equation (34), we obtain Equation (32) in the large time asymptotic limit as follows (cf. Equation (A12)):

(35)$P_1(x,t)\approx \frac{d{\eta }_0\left(x\right)}{dx}\left[1+2\sum _{n=1}^{\infty }cos\left(n{\eta }_0\right)e^{-n^2{\overline{\omega }}^{2}t}-\frac{2}{\overline{\omega }{\left(\pi t\right)}^{\frac{1}{2}}}\sum _{n=1}^{\infty }\frac{sin\left[n{\eta }_0\left(x\right)\right]}{n}\right].$

Since this solution is the large time approximation, to keep the conservation of the probability, a normalization constant should be introduced, which is ${\overline{P}}_1^{st}\left(x\right)={\int }_{-\pi }^{\pi }d{\eta }_0=2\pi$, while integration of the second and third terms yields zero.

4.3. Two-Dimensional PDF

The complete description of the transport on the two-dimensional comb is according to the comb PDF (3). Performing the Laplace transform of the backbone PDF (26) and substituting it into Equation (3), we obtain

(36)$P(x,y,t)={\mathcal{L}}^{-1}\left[\tilde{P}(x,y,s)\right]={\left(4D\right)}^{-\frac{1}{2}}\frac{d{\eta }_0\left(x\right)}{dx}\sum _n{e}^{in{\eta }_0}\frac{1}{2\pi i}{\int }_{c-i\infty }^{c+i\infty }\frac{e^{-\left|y\right|\sqrt{s/D}}{e}^{st}}{s^{\frac{1}{2}}-in\overline{\omega }}ds,$

The Laplace inversion corresponds to the table integral [31], which reads

(37)$\frac{1}{2\pi i}{\int }_{c-i\infty }^{c+i\infty }\frac{e^{-\left|y\right|\sqrt{s/D}}{e}^{st}}{s^{\frac{1}{2}}-in\overline{\omega }}ds=\frac{e^{-\frac{y^2}{4Dt}}}{\sqrt{\pi t}}+in\overline{\omega }{e}^{-in\frac{\overline{\omega }\left|y\right|}{\sqrt{D}}-n^2{\overline{\omega }}^{2}t}\mathrm{erfc}\left(\frac{\left|y\right|}{\sqrt{4Dt}}-in\overline{\omega }{t}^{\frac{1}{2}}\right),$

where $\mathrm{erfc}\left(z\right)$ is the (complementary) error function [29]. Thus the comb PDF consists of two terms $P(x,y,t)=P^{\left(1\right)}(x,y,t)+P^{\left(2\right)}(x,y,t)$. The first one corresponds to Brownian motion inside the fingers at the fixed initial backbone’s position,

(38)$P^{\left(1\right)}(x,y,t)={\left(4D\right)}^{-\frac{1}{2}}\frac{d{\eta }_0\left(x\right)}{dx}\sum _n{e}^{in{\eta }_0}\frac{e^{-\frac{y^2}{4Dt}}}{\sqrt{\pi t}}=\frac{e^{-\frac{y^2}{4Dt}}}{\sqrt{4\pi Dt}}\delta (x-x_0).$

The second term describes the interplay of Brownian motion along the fingers and anomalous transport along the backbone,

(39)$P^{\left(2\right)}(x,y,t)=\frac{i\overline{\omega }}{2\sqrt{D}}\frac{d{\eta }_0\left(x\right)}{dx}\sum _{n}n{e}^{in{\eta }_0\left(x\right)}{e}^{-in\frac{\overline{\omega }\left|y\right|}{\sqrt{D}}-n^2{\overline{\omega }}^{2}t}\mathrm{erfc}\left(\frac{\left|y\right|}{\sqrt{4Dt}}-in\overline{\omega }{t}^{\frac{1}{2}}\right).$

Integration of the comb PDF with respect to x yields the first term in the form of Brownian motion along the fingers, while the second term vanishes according to the integral

${\int }_0^{\infty }{\eta }_0^{\prime }\left(x\right)sin\left({\eta }_0\right)dx={\int }_{-\pi }^{\pi }sin\left({\eta }_0\right)d{\eta }_0=0.$

Therefore, integration of the comb PDF with respect to x yields the finger’s marginal PDF $P_2(y,t)=\frac{e^{-\frac{y^2}{4Dt}}}{\sqrt{\pi t}}$ in the form of Brownian motion.

Let us estimate the second term by means of the asymptotic expansion of the error function [29]

$e^{z^2}\mathrm{erfc}\left(z\right)\approx \frac{1}{\sqrt{\pi }z}.$

Then we have

(40)$P^{\left(2\right)}(x,y,t)=\frac{i\overline{\omega }}{\sqrt{4D\pi }}\frac{d{\eta }_0\left(x\right)}{dx}\sum _{n}n{e}^{in{\eta }_0\left(x\right)}\frac{e^{-\frac{y^2}{4Dt}}}{\left(\frac{\left|y\right|}{\sqrt{4Dt}}-in\overline{\omega }{t}^{\frac{1}{2}}\right)}.$

To perform summation over n, the denominator in Equation (40) is presented in the integral form,

(41)$\frac{1}{\left(\frac{\left|y\right|}{\sqrt{4Dt}}-in\overline{\omega }{t}^{\frac{1}{2}}\right)}={\int }_0^{\infty }{e}^{-u\left(\frac{\left|y\right|}{\sqrt{4Dt}}-in\overline{\omega }{t}^{\frac{1}{2}}\right)}du.$

This yields

(42)$\begin{array}{c}{\displaystyle \sum _n}n{e}^{in{\eta }_0\left(x\right)}{e}^{iun\overline{\omega }{t}^{\frac{1}{2}}}=\frac{1}{i}\frac{d}{d{\eta }_0\left(x\right)}{\displaystyle \sum _{k=-\infty }^{\infty }}\delta \left(u\overline{\omega }{t}^{\frac{1}{2}}+{\eta }_0\left(x\right)-2\pi k\right)\hfill \\ \hfill ={\left(i{\overline{\omega }}^{2}t\right)}^{-1}\frac{d}{du}{\displaystyle \sum _{k=-\infty }^{\infty }}\delta \left(u+\frac{{\eta }_0\left(x\right)-2\pi k}{\overline{\omega }{t}^{\frac{1}{2}}}\right).\end{array}$

Note that this integration yields $u=\frac{2\pi k-{\eta }_0\left(x\right)}{\overline{\omega }{t}^{\frac{1}{2}}}\ge 0$, which means that k is a non-negative integer. Integrating with respect to u, we obtain

(43)$\begin{array}{c}{P}^{\left(2\right)}(x,y,t)=\frac{e^{-\frac{y^2}{4Dt}}{\left(\overline{\omega }t\right)}^{-1}}{\sqrt{4D\pi }}\frac{d{\eta }_0\left(x\right)}{dx}{\int }_0^{\infty }{e}^{-u\frac{\left|y\right|}{\sqrt{4Dt}}}\frac{d}{du}{\displaystyle \sum _{k=0}^{\infty }}\delta \left(u+\frac{{\eta }_0\left(x\right)-2\pi k}{\overline{\omega }{t}^{\frac{1}{2}}}\right)\hfill \\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}=\frac{d{\eta }_0\left(x\right)}{dx}\frac{\left|y\right|}{4D\overline{\omega }\sqrt{\pi t^3}}exp\left[-\frac{y^2}{4Dt}+{\eta }_0\left(x\right)\frac{\left|y\right|}{\sqrt{4D}\overline{\omega }t}\right]\frac{1}{1-e^{-2\pi \frac{\left|y\right|}{\sqrt{4D}\overline{\omega }t}}}.\end{array}$

It is interesting to admit that integration of the asymptotic expression (40) with respect to x yields $P^{\left(2\right)}(x,y,t)=0$ due to the oscillations in time. However, summation over n smooths out these oscillations. Moreover, the summation leads to entanglement of the x and y directions of anomalous diffusion, and the integration with respect to x does not lead to Brownian motion on the fingers. Therefore, for the approximate expression, summation and integration do not commute.

5. Mean Squared Displacement

In this section, we estimate the mean squared displacement (MSD) for the case of fractional diffusion along the backbone described by the marginal PDF $P_1(x,t)$.

Taking into account Equations (32) and (33), we consider the marginal PDF in the following form:

(44)$P_1(x,t)=\frac{d{\eta }_0\left(x\right)}{dx}\left[1+2\sum _{n=1}^{\infty }cos\left(n{\eta }_0\right)e^{-n^2{\overline{\omega }}^{2}t}+i\overline{\omega }{t}^{\frac{1}{2}}\sum _n{}^{\prime }{e}^{in{\eta }_0}n{E}_{1,\frac{3}{2}}\left(-n^2{\overline{\omega }}^{2}t\right)\right],$

where we used

(45)$\begin{array}{c}{E}_{\frac{1}{2}}\left(in\overline{\omega }{t}^{\frac{1}{2}}\right)={\displaystyle \sum _{k=0}^{\infty }}\frac{{(-n^2{\overline{\omega }}^{2}t)}^k}{\mathrm{\Gamma }\left(k+1\right)}+in\overline{\omega }{t}^{\frac{1}{2}}{\displaystyle \sum _{k=0}^{\infty }}\frac{{(-n^2{\overline{\omega }}^{2}t)}^k}{\mathrm{\Gamma }\left(k+\frac{3}{2}\right)}\hfill \\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=e^{-n^2{\overline{\omega }}^{2}t}+in\overline{\omega }{t}^{\frac{1}{2}}{E}_{1,\frac{3}{2}}\left(-n^2{\overline{\omega }}^{2}t\right).\end{array}$

For the summations, we consider the asymptotic expression (35)

(46)$P_1(x,t)\approx \frac{d{\eta }_0\left(x\right)}{dx}\left[1+2\sum _{n=1}^{\infty }cos\left(n{\eta }_0\right)e^{-n^2{\overline{\omega }}^{2}t}-\frac{2}{\overline{\omega }{\left(\pi t\right)}^{\frac{1}{2}}}\sum _{n=1}^{\infty }\frac{sin\left[n{\eta }_0\left(x\right)\right]}{n}\right].$

Summations in Equation (46) are the table series [32]. We start from the last term, which depends on the limits of ${\eta }_0\left(x\right)$; see Section 5.4.2 in Ref. [32]. Note that the limits of the integration follow from the definition of the Bessel function in Equation (23). However, the physical limits of the angle ${\eta }_0\left(x\right)$ are determined by the transport/diffusive area on the backbone, which is $x_0{e}^{-v^2/D\omega }\le x\le x_0{e}^{v^2/D\omega }$. Substituting the limits in Equation (24b), we obtain the physical limits, which are according to the expression $-1\le sin\left[{\eta }_0\left(x\right)\right]\le 1$. Therefore, taking ${\eta }_0\left(x\right)\in [\frac{\pi }{2},\frac{3\pi }{2}]$, we obtain that the summation in Equation (46) yields $\frac{\pi -{\eta }_0\left(x\right)}{2}$, where ${\eta }_0\left(x\right)=\pi$ for $x=x_0$.

Eventually, we obtain that the time in the second Mittag-Leffler function is separated from the space coordinate that results in the time decay $\sim t^{-\frac{1}{2}}$ for any averaged values. Therefore, this term can be neglected for the MSD.

Now we estimate the marginal PDF in Equation (33) according to the first term in Equation (45), which is

(47)$\begin{array}{c}{P}_1(x,t)=\frac{d{\eta }_0\left(x\right)}{dx}\left[1+{\displaystyle \sum _n}{}^{\prime }{e}^{in{\eta }_0}{E}_{1,1}\left(-n^2{\overline{\omega }}^{2}t\right)\right]\hfill \\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\frac{d{\eta }_0\left(x\right)}{dx}\left[{\displaystyle \sum _{n=0}^{\infty }}2e^{-n^2{\overline{\omega }}^{2}t}cos\left[n{\eta }_0\left(x\right)\right]-1\right].\end{array}$

The summation yields ( see Section 5.4.11.2 in Ref. [32])

(48)$P_1(x,t)=\frac{d{\eta }_0\left(x\right)}{dx}{\vartheta }_3\left(\frac{{\eta }_0}{2},\kappa \right),$

where ${\vartheta }_3\left(\frac{{\eta }_0}{2},\kappa \right)$ is the elliptic theta function [29,33] and $e^{-{\omega }^{2}t}\equiv e^{-\pi \kappa }$ and $\kappa ={\omega }^{2}t/\pi$, respectively, and ${\eta }_0\equiv {\eta }_0\left(x\right)$.

The specific property of the theta function is that it describes diffusion in the finite area ${\eta }_0\in [-\pi ,\pi ]$ according to the diffusion equation, which in our case reads [33]

(49)$\frac{\partial {\vartheta }_3\left(\frac{{\eta }_0}{2},\kappa \right)}{\partial \kappa }=\frac{1}{\pi }\frac{{\partial }^2{\vartheta }_3\left(\frac{{\eta }_0}{2}\kappa \right)}{\partial {\eta }_0^2}$

with the reflecting boundary conditions at ${\eta }_0=\pm \pi$ and the initial condition ${\vartheta }_3\left(\frac{{\eta }_0}{2},\kappa =0\right)=\delta \left({\eta }_0\right)$.

According to the method of images (see, e.g., Ref. [34]), a formal solution of the equation reads

(50)$\begin{array}{c}{\vartheta }_3\left(\frac{{\eta }_0}{2},\kappa \right)={\left(4\kappa \right)}^{-\frac{1}{2}}{\displaystyle \sum _{n=-\infty }^{\infty }}exp\left[-\frac{\pi {({\eta }_0+2n\pi )}^2}{4\kappa }\right]\hfill \\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}=\frac{e^{-\frac{\pi {\eta }_0^2}{4\kappa }}}{2\sqrt{\kappa }}\left\{1+2{\displaystyle \sum _{n=1}^{\infty }}{e}^{-n^2{\pi }^3/\kappa }cosh\left(\frac{n{\pi }^2{\eta }_0}{\kappa }\right)\right\}.\end{array}$

It follows from Equation (50) that the leading term is the normal distribution, which is the fundamental Gaussian solution to diffusion Equation (50) for ${\eta }_0\in \mathbb{R}$, which corresponds to the absence of the boundaries [34].

Therefore, the MSD of the angle ${\eta }_0$ is $⟨{\eta }_0^2\left(\kappa \right)⟩=\kappa /\pi ={\omega }^{2}t/{\pi }^2$. However, the angle MSD is restricted by the diffusive area, and therefore it behaves like a sawtooth function of time. Then, taking the sine function from both sides of the MSD and taking into account Equation (24b), we obtain

$sin\left(\sqrt{⟨{\eta }_0^2⟩}\right)\approx \frac{\omega }{2\overline{v}}ln\left(⟨x^2⟩/x_0^2\right)=sin\left(\omega t^{\frac{1}{2}}/\pi \right).$

Then, from Equation (24b), we also have $x^2=x_0^2{e}^{2ln(x/x_0)}=x_0^{2}exp\left[2\frac{\overline{v}}{\omega }sin\left({\eta }_0\right)\right]$, and the average yields the approximate expression for the MSD as follows:

(51)$⟨x^2\left(t\right)⟩\sim x_0^{2}exp\left[2\frac{\overline{v}}{\omega }sin\left(\omega t^{\frac{1}{2}}/\pi \right)\right].$

6. Conclusions

The anomalous transport of particles in the presence of a time-dependent field is considered in the framework of the two-dimensional $(x,y)$ diffusion equation (2), where the space dependence of the transport coefficients corresponds to the comb geometry. In this case, the transport along the x direction, which is the backbone, is possible only at $y=0$, while Brownian motion along the y direction, or fingers, takes place for any $x\in \mathbb{R}$. The diffusion equation is known as a comb model. There are various transport scenarios along the backbone [25,26].

Here a new scenario is suggested in the form of inhomogeneous advection with the periodic in time velocity, $v\left(t\right){\partial }_{x}x$, where $v\left(t\right)=vcos\left(\omega t\right)$. This dependence of inhomogeneous advection on time leads to an essential difference in the analysis and the results from the time-independent case. The main issue of the consideration is the marginal PDF $P_1(x,t)$, which describes the classical relaxation transport for the dilatation operator at the condition of the conservation of the number of particles (or the probability). However, for $\omega \ne 0$, the straightforward equation for the marginal PDF cannot be obtained. Therefore, the backbone PDF $f(x,t)$ is studied first and then using the relation (5), the marginal PDF is obtained.

One of the central results of the analysis is Equation (48), which is a concise asymptotic form of the marginal PDF $P_1(x,t)$, described by the $({\eta }_0,\kappa )$ variables, where ${\eta }_0={\eta }_0\left(x\right)\in [-\pi ,\pi ]$ is an effective coordinate, while $\kappa ={\omega }^{2}t/\pi \in {\mathbb{R}}_{+}$ is an effective time. The correspondence between the real $(x,t)$ variables and the effective $({\eta }_0,\kappa )$ variables follows immediately from the definition of averaged values. Indeed, for an arbitrary function $g\left(x\right)$, we have

$⟨g\left(x\right)⟩\left(t\right)={\int }_0^{\infty }g\left(x\right)P_1(x,t)dx={\int }_{-\pi }^{\pi }G\left({\eta }_0\right){\vartheta }_3({\eta }_0,\kappa )d{\eta }_0=⟨G\left({\eta }_0\right)⟩\left(\kappa \right),$

where $G\left({\eta }_0\right)=g\left(x_0{e}^{sin\left({\eta }_0\right)}\right)$, and the main contribution to the elliptic theta function ${\vartheta }_3({\eta }_0,\kappa )$ is according to the normal distribution in Equation (50) in the effective $({\eta }_0,\kappa )$ space. To some extent, this very specific normal distribution is an effective relaxation, described in the framework of the $({\eta }_0,\kappa )$ variables.

Another important point of the analysis is the comb parameters and correspondence of their dimensions. It should be pointed out that in the standard comb model, the transport is controlled by two parameters only; these are the backbone velocity v and the finger’s diffusion coefficient D (in the present case). Performing the scaling of combinations of these parameters, it is possible to describe the comb model in the framework of dimensionless variables and parameters. In the present case, which is performed in the dimension framework, the time-dependent comb model is controlled by an additional parameter, which is the frequency $\omega$. This leads to the appearance of additional control parameters $\overline{v}=\frac{v^2}{D}$ and $\overline{\omega }=\omega \frac{D^{\frac{1}{2}}}{2v}$, which control the correct dimension of the obtained fractional transport equation for the comb PDF and the Floquet theory.

In conclusion, it should be pointed out that this geometrically constrained transport leads to anomalous diffusion along the backbone. When inhomogeneous convection is taken in the form of a dilatation operator (see Appendix A.3), then the fractional diffusion equation, obtained for the backbone PDF, becomes the quantum fractional Schrödinger equation just by a simple multiplication by the imaginary unit and the Planck constant, $i\hslash$. In this case, the marginal PDF is also the Green function for the fractional quantum transport. This quantum case, however, deserves separate consideration, as well.

Data Availability Statement

No new data were created or analyzed in this study. Data sharing is not applicable to this article.

Conflicts of Interest

The author declares no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

References

1. Hentschel, H.G.E.; Procaccia, I. Relative diffusion in turbulent media: The fractal dimension of clouds. Phys. Rev. A; 1984; 29, pp. 1461-1470. [DOI: https://dx.doi.org/10.1103/PhysRevA.29.1461]

2. Csanady, G. Turbulent Diffusion in the Environment; Reidel: Dortdrecht, The Netherlands, 1973.

3. Khoury, M.; Lacasta, A.M.; Sancho, J.M.; Romero, A.H.; Lindenberg, K. Charged particle transport in antidot lattices in the presence of magnetic and electric fields: Langevin approach. Phys. Rev. B; 2008; 78, 155433. [DOI: https://dx.doi.org/10.1103/PhysRevB.78.155433]

4. Shkilev, V.P. Subdiffusion in a Time-Dependent Force Field. J. Exp. Theor. Phys.; 2012; 114, pp. 830-835. [DOI: https://dx.doi.org/10.1134/S1063776112030089]

5. Burioni, R.; Cassi, D.; Giusiano, G.; Regina, S. Anomalous diffusion and Hall effect on comb lattices. Phys. Rev. E; 2003; 67, 016116. [DOI: https://dx.doi.org/10.1103/PhysRevE.67.016116]

6. Wu, J.; Berland, K.M. Propagators and Time-Dependent Diffusion Coefficients for Anomalous Diffusion. Biophys. J.; 2008; 95, pp. 2049-2052. [DOI: https://dx.doi.org/10.1529/biophysj.107.121608]

7. Mainardi, F. Fractional relaxation-oscillation and fractional diffusion-wave phenomena. Chaos Solitons Fractals; 1996; 7, pp. 1461-1477. [DOI: https://dx.doi.org/10.1016/0960-0779(95)00125-5]

8. Zhou, X.; Gao, Q.; Abdullah, O.; Magin, R. Studies of anomalous diffusion in the human brain using fractional order calculus. Magn. Reson. Med.; 2010; 63, pp. 562-569. [DOI: https://dx.doi.org/10.1002/mrm.22285] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/20187164]

9. Magin, R.L.; Abdullah, O.; Baleanu, D.; Zhou, X.J. Anomalous diffusion expressed through fractional order differential operators in the Bloch–Torrey equation. J. Magn. Reson.; 2008; 190, pp. 255-270. [DOI: https://dx.doi.org/10.1016/j.jmr.2007.11.007]

10. Matheron, G.; De Marsily, G. Is transport in porous media always diffusive? A counterexample. Water Resour. Res.; 1980; 16, pp. 901-917. [DOI: https://dx.doi.org/10.1029/WR016i005p00901]

11. Iomin, A. Quantum dynamics and relaxation in comb turbulent diffusion. Chaos Solitons Fractals; 2020; 139, 110305. [DOI: https://dx.doi.org/10.1016/j.chaos.2020.110305]

12. Sandev, T.; Iomin, A.; Kocarev, L. Hitting times in turbulent diffusion due to multiplicative noise. Phys. Rev. E; 2020; 102, 042109. [DOI: https://dx.doi.org/10.1103/PhysRevE.102.042109]

13. Black, F.; Scholes, M. The Pricing of Options and Corporate Liabilities. J. Political Econ.; 1973; 81, pp. 637-654. [DOI: https://dx.doi.org/10.1086/260062]

14. Armitage, J.V. Chapter 10-The Riemann hypothesis and the Hamiltonian of a quantum mechanical system. Number Theory and Dynamical Systems; London Mathematical Society Lecture Note Series Dodson, M.M.; Vickers, J.A.G. Cambridge University Press: Cambridge, UK, 1989; Volume 134, pp. 153-172.

15. Teschl, G. Ordinary Differential Equations and Dynamical Systems; Graduate Studies in Mathematics AMS: Providence, RI, USA, 2012; Volume 140.

16. Grifoni, M.; Hänggi, P. Driven quantum tunneling. Phys. Rep.; 1998; 304, pp. 229-354. [DOI: https://dx.doi.org/10.1016/S0370-1573(98)00022-2]

17. Santoro, G.E. Introduction to Floquet; Lecture Notes SISSA: Trieste, Italy, Available online: https://www.ggi.infn.it/sft/SFT_2019/LectureNotes/Santoro.pdf (accessed on 18 May 2023).

18. Casati, G.; Molinari, L. “Quantum chaos” with time periodic Hamiltonians. Prog. Theor. Phys. Suppl.; 1989; 98, pp. 287-322. [DOI: https://dx.doi.org/10.1143/PTPS.98.287]

19. Satija, I.I. The Butterfly in the Quantum World: The Story of the Most Fascinating Quantum Fractal; Morgan & Claypool Publishers: San Rafael, CA, USA, 2016.

20. Iomin, A. Fractional Floquet theory. Chaos Solitons Fractals; 2023; 168, 113196. [DOI: https://dx.doi.org/10.1016/j.chaos.2023.113196]

21. Furutsu, K. On the theory of radio wave propagation over inhomogeneous earth. J. Res. Nat. Bur. Stand; 1963; 67, pp. 39-62.

22. Novikov, E.A. Functionals and the Random-Force Method in Turbulence Theory. Sov. Phys. JETP; 1965; 20, pp. 1290-1294.

23. Zinn-Justin, J. Quantum Field Theory and Critical Phenomena; Clarendon Press: Oxford, UK, 1990.

24. Arkhincheev, V.E.; Baskin, E.M. Anomalous diffusion and drift in the comb model of percolation clusters. Sov. Phys. JETP; 1991; 73, pp. 161-165.

25. Iomin, A.; Méndez, V.; Horsthemke, W. Fractional Dynamics in Comb-like Structures; World Scientific: Singapore, 2018.

26. Sandev, T.; Iomin, A. Special Functions of Fractional Calculus; World Scientific: Singapore, 2022.

27. Mainardi, F. Fractional Calculus and Waves in Linear Viscoelasticity; Imperial College Press: Singapore, 2010.

28. Bateman, H.; Erdélyi, A. Higher Transcendental Functions [Volumes I–III]; McGraw-Hill: New York, NY, USA, 1953.

29. Abramovitz, M.; Stegun, I.A. Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables; Dover Publications: New York, NY, USA, 1972.

30. Baskin, E.; Iomin, A. Superdiffusion on a Comb Structure. Phys. Rev. Lett.; 2004; 93, 120603. [DOI: https://dx.doi.org/10.1103/PhysRevLett.93.120603]

31. Bateman, H.; Erdélyi, A. Tables of Integral Transforms [Volumes I & II]; McGraw-Hill: New York, NY, USA, 1954.

32. Prudnikov, A.; Brychkov, Y.A.; Marichev, O.I. Integrals and Series; Nauka: Moscow, Russia, 1981.

33. Jahnke, E.; Emde, F.; Lösch, F. Tables of Higher Functions; McGraw-Hill: New York, NY, USA, 1960.

34. Balakrishnan, V. Elements of Nonequilibrium Statistical Mechanics; Springer: Cham, Switzerland, 2021.

35. Howland, J.S. Stationary scattering theory for time-dependent Hamiltonians. Math. Ann.; 1974; 207, pp. 315-335. [DOI: https://dx.doi.org/10.1007/BF01351346]

36. Howland, J.S. Two problems with time-dependent Hamiltonians. Mathematical Methods and Applications of Scattering Theory; Lecture Notes in Physics DeSanto, J.A.; Sáenz, A.W.; Zachary, W.W. Springer: Berlin, Germany, 1980; Volume 130, pp. 163-168.

37. Berman, G.; Vishik, M. Long time evolution of quantum averages near stationary points. Phys. Lett. A; 2003; 319, pp. 352-359. [DOI: https://dx.doi.org/10.1016/j.physleta.2003.10.038]

38. Berry, M.V.; Keating, J.P. H = xp and the Riemann zeros. Supersymmetry and Trace Formulae; Springer: Berlin/Heidelberg, Germany, 1999; pp. 355-367.

39. Sierra, G. H=xp with interaction and the Riemann zeros. Nucl. Phys. B; 2007; 776, pp. 327-364. [DOI: https://dx.doi.org/10.1016/j.nuclphysb.2007.03.049]

40. Landau, L.D.; Lifshitz, E.M. Quantum Mechanics: Non-Relativistic Theory; Course of Theoretical Physics Pergamon: Oxford, UK, 1977; Volume 3.