Effect of Ergodic and Non-Ergodic Fluctuations on a Charge Diffusing in a Stochastic Magnetic Field

Abstract

In this paper, we study the basic problem of a charged particle in a stochastic magnetic field. We consider dichotomous fluctuations of the magnetic field where the sojourn time in one of the two states are distributed according to a given waiting-time distribution either with Poisson or non-Poisson statistics, including as well the case of distributions with diverging mean time between changes of the field, corresponding to an ergodicity breaking condition. We provide analytical and numerical results for all cases evaluating the average and the second moment of the position and velocity of the particle. We show that the field fluctuations induce diffusion of the charge with either normal or anomalous properties, depending on the statistics of the fluctuations, with distinct regimes from those observed, e.g., in standard Continuous-Time Random Walk models.

1. Introduction

Diffusive processes occur in many physical, chemical and engineering applications. When the diffusion process is taking place, the quantities related to the spreading species take random values. Since Einstein and Smoluchowski’s work on Brownian motion [1,2], diffusive phenomena have been a fundamental subject of intense research. Both derivations (Einstein and Smoluchowski’s) lead to the well-known diffusion relationship, in the one-dimensional case, $\langle x^2\rangle =2Dt$, with D the diffusion coefficient. Several relevant physical and biological phenomena have been discovered in the last few decades, showing an anomalous relationship between mean-squared displacement and time, $\langle x^2\rangle \propto t^{\alpha }$. For example, diffusion through porous media or within a crowded cellular environment, making anomalous diffusion a relevant subject of research work [3,4,5,6] and Refs. [7,8] for a review.

This paper presents a detailed study of a particle moving in a fluctuating magnetic field that is directly connected to plasma physics. In particular, it is relevant in many technological applications such as, for example, plasma confinement [9]. We focus on the diffusion of the particles caused by the magnetic field fluctuations, which can destroy the plasma confinement. The fluctuations of physical quantities are an almost inevitable occurrence and generate diffusion processes of the related quantity [10,11,12]. It is important, therefore, to have an adequate theoretical framework to model their effect. We aim to fill the gap relative to the case of diffusion in a fluctuating magnetic field, which, to our knowledge, has not yet been explored and modeled so far in the case of non-ordinary statistics.

The paper is organized as follows: In Section 2 we introduce and formally define the problem for generic fluctuations. In Section 3, we introduce the case of dichotomous fluctuations, and we afford a complete analytical and numerical treatment of Poisson and non-Poisson statistics. In Section 4, we look at the average squared displacement and, starting from the Poisson case, we move to non-Poisson, power-law distributed fluctuations exploring, therefore also the non-ergodic regime when the average time of the fluctuations diverges (see Ref. [13] for an extended discussion, and Ref. [14] and references therein for a review and a historical perspective). The diffusion properties in this regime are characterized by the analytical and numerical derivation of the average squared displacement. Section 5 and Section 6 summarize results and draw final conclusions.

2. Stochastic Equation for the Magnetic Force

Let us consider the following classical equation [15,16]

$$\begin{array}{c}\hfill m\frac{d\mathbf{v}}{dt}=q\left(\mathbf{v}\times \mathbf{B}+\mathbf{E}\right)\end{array}$$

(1)

where m is the mass of the charge q, $\mathbf{v}$ is the velocity, $\mathbf{B}$ the magnetic field, and $\mathbf{E}$ the electric field. We focus on the case where the particle is traveling in a region with a uniform magnetic field randomly fluctuating, i.e., $\mathbf{B}={\mathbf{B}}_0+{\mathbf{B}}_1(t)$ where ${\mathbf{B}}_1(t)=B_1\mathbf{\xi }(t)$ with $\mathbf{\xi }(t)$ the stochastic fluctuation of the magnetic field, where $\mathbf{\xi }(t)=\xi (t){\mathbf{n}}_{\xi }$. We shall consider a dichotomous fluctuation, with values of $\pm 1$ where the sojourn time in one of the two states are distributed according to the distribution $\psi (t)$. We will consider the case when $\psi (t)$ is an exponential, $\psi (t)=\gamma exp[-\gamma t]$, or Poisson case, and the case when $\psi (t)$ is a power law, $\psi (t)\propto t^{-\alpha -1}$ with $0<\alpha <2$, or non-Poisson case.

The main reason for choosing dichotomous fluctuations rests on the fact that taking the finite values fluctuations may better represent a physical system. Additionally, we stress that in the appropriate limit, the most well-known noises in the literature, such as the gaussian white noise and the white shot noise are recovered within this framework [17,18]. Furthermore, we consider the case where no external electric field is applied. Taking the magnetic field in the z direction, $\mathbf{B}\equiv B(t)\mathbf{k}=[B_0+B_1\xi (t)]\mathbf{k}$, and using cartesian coordinates, we have

$$\begin{array}{ccc}\hfill & & \frac{d{v}_x}{dt}=\omega (t)v_y+\frac{q}{m}{E}_x\hfill \end{array}$$

(2)

$$\begin{array}{ccc}\hfill & & \frac{d{v}_y}{dt}=-\omega (t)v_x+\frac{q}{m}{E}_y\hfill \end{array}$$

(3)

with $\omega (t)=qB(t)/m$ the time-dependent Larmor frequency and $E_x$, $E_y$ the induced electrical field. As a further simplification, we consider the time scale of the fluctuation much larger than the time associated with the unperturbed Larmor frequency ${\omega }_0=q{B}_0/m$. It is worthy to note that Equations (2) and (3) can be reduced to a second order stochastic differential equation for the function $z=x+iy$. The resulting equation has a strong similarity with the mechanical system studied in Refs. [19,20] where the authors study a linear damped oscillator with a noise perturbing both the oscillator mass and the friction. A detailed study of this equation is out of the scope of this paper, and it is left to an upcoming publication.

Neglecting the induced electrical field, $\mathbf{E}\approx 0$, (see also [21]), we have

$$\begin{array}{ccc}\hfill & & \frac{d{v}_x}{dt}=\omega (t)v_y\hfill \end{array}$$

(4)

$$\begin{array}{ccc}\hfill & & \frac{d{v}_y}{dt}=-\omega (t)v_x\hfill \end{array}$$

(5)

For the sake of completeness, we end this section showing the equation for the density probability $P(v_x.v_y,t)=P(\mathbf{v},t)$ associated with the stochastic Equations (4) and (5). To obtain a closed equation for $P(\mathbf{v},t)$ we will assume that $\xi$ is a Poisson process although, in the next sections, the analysis of the relevant quantiles, $\langle v_x(t)\rangle$, $\langle v_y(t)\rangle$, $\langle x{(t)}^2+y{(t)}^2\rangle$, will include also non-Poisson processes. Using the Liouville approach we write the following continuity equation

$$\begin{array}{c}\hfill \frac{\partial \rho }{\partial t}=-\nabla ·\left[\rho \mathbf{v}\times \mathit{\omega }\right]=-\nabla ·\left[\rho \mathbf{v}\times {\mathit{\omega }}_0\right]-\nabla ·\left[\xi \rho \mathbf{v}\times {\mathit{\omega }}_1\right]\end{array}$$

(6)

where for the sake of compactness we dropped the function arguments and we introduced the symbols

$$\nabla \equiv \frac{\partial }{\partial v_x}\mathbf{i}+\frac{\partial }{\partial v_y}\mathbf{j},\mathit{\omega }=\left[\frac{q{B}_0}{m}+\frac{q{B}_1}{m}\xi (t)\right]\mathbf{k}\equiv {\mathit{\omega }}_0+{\mathit{\omega }}_1\xi (t).$$

The stochastic density $\rho$ is related to $P(v_x.v_y,t)$ via the Van Kampen’s lemma [22] $\langle \rho \rangle =P(\mathbf{v},t)$ where the average is performed on the $\xi$ realizations. Additionally, we will use the Shapiro–Loginov formulae of differentiation [23]

$$\frac{\partial }{\partial t}\langle \xi (t)\rho (t)\rangle =-\gamma \langle \xi (t)\rho (t)\rangle +\langle \xi (t)\frac{\partial }{\partial t}\rho (t)\rangle ,$$

(7)

that holds true for processes with a nth correlation function fulfilling the condition

$$\frac{\partial }{\partial t}\langle \xi (t)\xi (t_1)\cdots \xi (t_n)\rangle =-\gamma \langle \xi (t)\xi (t_1)\cdots \xi (t_n)\rangle .$$

In particular, this applies to Poisson, Gaussian and Markov jump processes, with a correlation function given by $\langle \xi (t_1)\xi (t_2)\rangle \sim exp[-\gamma |t_2-t_1|]$. Taking the average of Equation (6), defining $P_1(\mathbf{v},t)\equiv \langle \xi \rho \rangle$ and using Equation (7), we may write the system

$$\begin{array}{ccc}\hfill & & \frac{\partial P(\mathbf{v},t)}{\partial t}=-\mathbf{v}\times {\mathit{\omega }}_0·\nabla P-\mathbf{v}\times {\mathit{\omega }}_1·\nabla P_1,\hfill \end{array}$$

(8)

$$\begin{array}{ccc}\hfill & & \frac{\partial P_1(\mathbf{v},t)}{\partial t}=-\gamma P_1\times {\mathit{\omega }}_0·\nabla P_1+\mathbf{v}\times {\mathit{\omega }}_1·\nabla P\hfill \end{array}$$

(9)

Taking the time derivative of Equation (8), and after some algebra we obtain the following equation for the probability density

$$\begin{array}{ccc}\hfill \frac{{\partial }^{2}P(\mathbf{v},t)}{\partial t^2}& =& -\gamma \frac{\partial P}{\partial t}-\gamma \mathbf{v}\times {\mathit{\omega }}_0·\nabla P+\mathbf{v}\times {\mathit{\omega }}_0·\nabla \left[\mathbf{v}\times {\mathit{\omega }}_0·\nabla P\right]+\hfill \\ \hfill & & -\mathbf{v}\times {\mathit{\omega }}_1·\nabla \left[\mathbf{v}\times {\mathit{\omega }}_1·\nabla P\right].\hfill \end{array}$$

(10)

3. Dichotomous Processes

As stated in Section 2, in this section we will consider a magnetic field with a fluctuating component which is assumed to be dichotomous. Dichotomous fluctuations have the nice property of assuming finite values but, despite their relative simplicity, they can be shown to allow one to recover both gaussian white noise and white shot noise [17] within an appropriate limit procedure. Formally the exact solution of Equations (4) and (5) is

$$\begin{array}{ccc}\hfill & & v_x(t)=v_{0}sin\left[{\int }_0^t\omega (u)du\right]=v_0\mathrm{Im}\left[exp\left[i{\omega }_{0}t\right]exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]\right],\hfill \end{array}$$

(11)

$$\begin{array}{ccc}\hfill & & v_y(t)=v_{0}cos\left[{\int }_0^t\omega (u)du\right]=v_0\mathrm{Re}\left[exp\left[i{\omega }_{0}t\right]exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]\right],\hfill \end{array}$$

(12)

where $v_0$ is the initial velocity along the y axis. For the average we have

$$\begin{array}{ccc}\hfill & & \langle v_x(t)\rangle =v_0\mathrm{Im}\left[exp\left[i{\omega }_{0}t\right]〈exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]〉\right],\hfill \end{array}$$

(13)

$$\begin{array}{ccc}\hfill & & \langle v_y(t)\rangle =v_0\mathrm{Re}\left[exp\left[i{\omega }_{0}t\right]〈exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]〉\right].\hfill \end{array}$$

(14)

As we infer from Equations (13) and (14), we need to evaluate the average of the exponential of the noise integral. For this purpose, we consider the stochastic equation

$$\frac{dU}{dt}=\xi (t)$$

(15)

with $\xi (t)$ a dichotomous fluctuation where the sojourn time in one of the two states are distributed according to the distribution $\psi (t)$. We assume that the event, occurring at each random time $t_i$, changes the $\xi (t)$ sign. When these events occur with a constant rate $\gamma$ this corresponds to a Poisson process, and is characterized by an exponential distribution. As stated in the Introduction, we will consider here both the case of exponential (Poisson) distribution with

$$\psi (t)=\gamma exp[-\gamma t]$$

(16)

and the case of non-Poisson process with power-law distribution, characterized by the following asymptotic behavior

$$\psi (t)\propto {\left(\frac{t}{T}\right)}^{-\alpha -1},t\gg T$$

(17)

with $0<\alpha <2$, which corresponds to a regime with a diverging second moment ($1<\alpha <2$) and diverging first and second moment ($0<\alpha <1$). The latter case, characterized by the absence of a finite time scale, corresponds to a condition of ergodicity breaking. The formal solution of Equation (15) is

$$U(t)={\int }_0^t\xi (u)du.$$

(18)

Considering the formal solution for the velocity of the charge, Equations (13) and (14), we need to evaluate the quantity $〈exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]〉$. For this purpose, we use the exact formula in the Laplace representation [24,25] (see Ref. [26] for detailed calculations)

$$\begin{array}{ccc}\hfill & & \mathcal{L}\left[〈exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]〉\right]=\frac{1}{2}\left[\frac{\left(1+\widehat{\psi }\left(s-i{\omega }_1\right)\right)\widehat{\mathsf{\Psi }}\left(s+i{\omega }_1\right)}{1-\widehat{\psi }\left(s-i{\omega }_1\right)\widehat{\psi }\left(s+i{\omega }_1\right)}+\right.\hfill \\ \hfill & & \left.\frac{\left(1+\widehat{\psi }\left(s+i{\omega }_1\right)\right)\widehat{\mathsf{\Psi }}\left(s-i{\omega }_1\right)}{1-\widehat{\psi }\left(s-i{\omega }_1\right)\widehat{\psi }\left(s+i{\omega }_1\right)}\right]=\mathrm{Re}\left[\frac{\left(1+\widehat{\psi }\left(s+i{\omega }_1\right)\right)\widehat{\mathsf{\Psi }}\left(s-i{\omega }_1\right)}{1-\widehat{\psi }\left(s-i{\omega }_1\right)\widehat{\psi }\left(s+i{\omega }_1\right)}\right]\hfill \end{array}$$

(19)

where $\widehat{\psi }(s)$ is the Laplace transform of $\psi (t)$ and $\mathsf{\Psi }\left(t\right)$ (and consequently $\widehat{\mathsf{\Psi }}(s)$ is its Laplace transform) is the probability that no switch occurs for a generic interval of time t, i.e.,

$$\mathsf{\Psi }\left(t\right)=1-{\int }_0^t\psi \left(t_1\right)d{t}_1={\int }_t^{\infty }\psi \left(t_1\right)d{t}_1$$

(20)

The Poisson case does not present difficulties and, in the time representation, gives the expression

$$\begin{array}{ccc}\hfill & & 〈exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]〉=\hfill \\ \hfill & & exp\left[-\frac{\gamma }{2}t\right]\left[\frac{\gamma sinh\left(t\sqrt{\frac{{\gamma }^2}{4}-{\omega }_1^2}\right)}{2\sqrt{\frac{{\gamma }^2}{4}-{\omega }_1^2}}+cosh\left(t\sqrt{\frac{{\gamma }^2}{4}-{\omega }_1^2}\right)\right].\hfill \end{array}$$

(21)

We are now in position to write a closed expression for Equations (13) and (14). The average of the velocity components is

$$\begin{array}{ccc}\hfill & & \langle v_x(t)\rangle =v_{0}sin\left[{\omega }_{0}t\right]\hfill \\ \hfill & & exp\left[-\frac{\gamma }{2}t\right]\left[\frac{\gamma sinh\left(t\sqrt{\frac{{\gamma }^2}{4}-{\omega }_1^2}\right)}{2\sqrt{\frac{{\gamma }^2}{4}-{\omega }_1^2}}+cosh\left(t\sqrt{\frac{{\gamma }^2}{4}-{\omega }_1^2}\right)\right]\hfill \\ \hfill & & \langle v_y(t)\rangle =v_{0}cos\left[{\omega }_{0}t\right]\hfill \end{array}$$

(22)

$$\begin{array}{ccc}\hfill & & exp\left[-\frac{\gamma }{2}t\right]\left[\frac{\gamma sinh\left(t\sqrt{\frac{{\gamma }^2}{4}-{\omega }_1^2}\right)}{2\sqrt{\frac{{\gamma }^2}{4}-{\omega }_1^2}}+cosh\left(t\sqrt{\frac{{\gamma }^2}{4}-{\omega }_1^2}\right)\right]\hfill \end{array}$$

(23)

We now study Equation (19) for the non-Poisson case with a power-law distribution of the type of Equation (16) with $0<\alpha <2$ where the non-ergodicity of the process plays an important role. Some difficulties arise in inverting the Laplace transform, in particular in the region defined by $1<\alpha <2$. This is partially due to the fact that while for $0<\alpha <1$ the calculation in the Laplace transform can be carried out using as waiting-time distribution the derivative of the Mittag–Leffler function (see for example Ref. [8] and references therein), for the region $0<\alpha <2$ the Laplace transform is usually a complicated function, and the inversion of the final result is not an easy task. Traditionally the inversion of a Laplace transform for a large value of time t is performed using the Tauberian theorem, i.e., taking the development for the Laplace parameter $s\to 0$. If the function to invert is a hard-to-handle function, it is not always clear where to stop the development (see [27] for detailed examples). To overcome this difficulty, we may use as waiting-time distribution [28]

$$\psi (t)=\frac{sin\left(\frac{\pi \alpha }{2}\right)cos\left(\frac{t}{T}\right)+cos\left(\frac{\pi \alpha }{2}\right)sin\left(\frac{t}{T}\right)}{Tcos\left(\frac{\pi \alpha }{2}\right)}-\frac{sin\left(\frac{\pi \alpha }{2}\right){cos}_{\alpha }\left(\frac{t}{T}\right)+cos\left(\frac{\pi \alpha }{2}\right){sin}_{\alpha }\left(\frac{t}{T}\right)}{Tcos\left(\frac{\pi \alpha }{2}\right)}.$$

(24)

where T is a time-scale parameter and, by definition [29,30],

$$\begin{array}{c}\hfill {cos}_{\alpha }t\equiv \frac{E_{\alpha }^{it}+E_{\alpha }^{-it}}{2},{sin}_{\alpha }t\equiv \frac{E_{\alpha }^{it}-E_{\alpha }^{-it}}{2i}\end{array}$$

(25)

and

$$E_{\alpha }^t\equiv D_t^{\alpha }exp[t]=\sum _{n=0}^{\infty }\frac{t^{n-\alpha }}{\mathsf{\Gamma }(n+1-\alpha )},$$

(26)

where $D_t^{\alpha }$ is the Riemann-Liouville fractional derivative. The functions ${cos}_{\alpha }t$ and ${sin}_{\alpha }t$ compensate the oscillatory behavior of the ordinary trigonometric functions, and what remains is a positive power law, i.e., $t^{-\alpha -1}$. A rigorous proof is given in Ref. [28]. Despite its complicated structure in time representation, its Laplace transform is

$$\widehat{\psi }(s)=\frac{1+sTtan\left(\frac{\pi \alpha }{2}\right)-sec\left(\frac{\pi \alpha }{2}\right){(sT)}^{\alpha }}{{(sT)}^2+1},0<\alpha <2,\alpha \ne 1.$$

(27)

For $\alpha =1$ we must take the limit and we obtain

$$\widehat{\psi }(s)=\underset{\alpha \to 1}{lim}\frac{1+sTtan\left(\frac{\pi \alpha }{2}\right)-sec\left(\frac{\pi \alpha }{2}\right){(sT)}^{\alpha }}{{(sT)}^2+1}=\frac{1+\frac{2sT}{\pi }log(sT)}{{(sT)}^2+1}.$$

(28)

The proposed distribution has a simple structure based on power law and its validity is in the non-Poisson ranges $0<\alpha <2$. Using the property $\mathcal{L}\left[exp[\pm i{\omega }_{1}t]f(t)\right]=\widehat{f}(s\mp i{\omega }_1)$, defining the new Laplace variables, $v=s-i{\omega }_1$, we can reduce the inverse Laplace problem of Equation (19) the following inversion Laplace transform

$$\begin{array}{c}\hfill \mathrm{Re}\left[\frac{\left(1+\widehat{\psi }\left(v+2i{\omega }_1\right)\right)\widehat{\mathsf{\Psi }}\left(v\right)}{1-\widehat{\psi }\left(v+2i{\omega }_1\right)\widehat{\psi }\left(v\right)}\right].\end{array}$$

(29)

Since we are interested in the asymptotic limit, taking the limit for $v\to 0$ we may invert Equation (29). We find for the asymptotic expression

$$\begin{array}{c}\hfill 〈exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]〉\approx A\frac{exp\left[i{\omega }_{1}t+\varphi \right]}{t^{\alpha }}\end{array}$$

(30)

where A and $\varphi$ are constant depending on ${\omega }_1,T,\alpha$. In the region $0<\alpha <1$ there is no dependence on the parameter T and we obtain

$$\begin{array}{ccc}\hfill & & 〈exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]〉\approx 2\mathrm{Re}\left[exp(i{\omega }_{1}t)E_{\alpha }\left[-{(2i{\omega }_{1}t)}^{\alpha }\right]\right]\hfill \\ \hfill & & \approx \frac{\sqrt{\pi }{J}_{\alpha -\frac{1}{2}}({\omega }_{1}t)}{\mathsf{\Gamma }(1-\alpha ){(2{\omega }_{1}t)}^{\alpha -\frac{1}{2}}}\hfill \end{array}$$

(31)

where $E_{\alpha }\left[-{(2i{\omega }_{1}t)}^{\alpha }\right]$ is the Mittag–Leffler function defined as

$$E_{\alpha }\left[z\right]=\sum _{n=0}^{\infty }\frac{z^n}{\mathsf{\Gamma }(n\alpha +1)},$$

and $J_{\nu }(z)$ is the Bessel function of the first kind. From an asymptotic point of view, all the expressions contained in Equations (30) and (31) have the same accuracy. The advantage of the expression written as in the last line of Equation (31) is that in the case $\alpha =1/2$, it provides an exact expression, i.e.

$$\begin{array}{ccc}\hfill & & \langle v_x(t)\rangle =v_0\mathrm{Im}\left[exp\left[i{\omega }_{0}t\right]〈exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]〉\right]=\hfill \\ \hfill & & v_{0}sin\left[{\omega }_{0}t\right]J_0\left[{\omega }_{1}t\right],\hfill \end{array}$$

(32)

$$\begin{array}{ccc}\hfill & & \langle v_y(t)\rangle =v_0\mathrm{Re}\left[exp\left[i{\omega }_{0}t\right]〈exp\left[i{\omega }_1{\int }_0^t\xi (u)du\right]〉\right]=\hfill \\ \hfill & & v_{0}cos\left[{\omega }_{0}t\right]J_0\left[{\omega }_{1}t\right].\hfill \end{array}$$

(33)

This can be directly checked using the distribution $P(x,t)$ derived by Lamperti [31,32], and which describes the distribution associated with Equation (15) for $0<\alpha <1$. Integrating Equations (32) and (33) the result with respect to time, we obtain

$$\begin{array}{ccc}\hfill & & \langle x(t)\rangle =v_0{\int }_0^{t}sin\left[{\omega }_{0}t\right]J_0\left[{\omega }_{1}t\right]dt+x_0,\hfill \end{array}$$

(34)

$$\begin{array}{ccc}\hfill & & \langle y(t)\rangle =v_0{\int }_0^{t}cos\left[{\omega }_{0}t\right]J_0\left[{\omega }_{1}t\right]dt+y_0\hfill \end{array}$$

(35)

Asymptotically we have

$$\begin{array}{ccc}\hfill & & x(\infty )=\left\{\begin{array}{cc}\frac{v_0}{\sqrt{{\omega }_0^2-{\omega }_1^2}}+x_0,\hfill & \mathrm{if}{\omega }_1<{\omega }_0,\hfill \\ x_0,\hfill & \mathrm{if}{\omega }_0<{\omega }_1.\hfill \end{array}\right.\hfill \end{array}$$

(36)

$$\begin{array}{ccc}\hfill & & y(\infty )=\left\{\begin{array}{cc}\frac{v_0}{\sqrt{{\omega }_1^2-{\omega }_0^2}},\hfill & \mathrm{if}{\omega }_0<{\omega }_1,\hfill \\ y_0,\hfill & \mathrm{if}{\omega }_1<{\omega }_0.\hfill \end{array}\right.\hfill \end{array}$$

(37)

The resonant case ${\omega }_1={\omega }_0$ generates diverging average positions

$$\begin{array}{ccc}\hfill & & \langle x(t)\rangle =v_{0}t[sin({\omega }_{0}t)J_0({\omega }_{0}t)-cos({\omega }_{0}t)J_1({\omega }_{0}t)]+x_0,\hfill \end{array}$$

(38)

$$\begin{array}{ccc}\hfill & & \langle y(t)\rangle =v_{0}t[sin({\omega }_{0}t)J_1({\omega }_{0}t)+cos({\omega }_{0}t)J_0({\omega }_{0}t)]+y_0.\hfill \end{array}$$

(39)

Figure 1 and Figure 2 show the comparison between analytical results and numerical simulations for the quantities $\langle v_x(t)\rangle$ and $\langle v_y(t)\rangle$. Figure 3, Figure 4, Figure 5 and Figure 6 show single realizations of the stochastic trajectories, and Figure 7 and Figure 8 show the comparison between analytical results and numerical simulations for the quantities $\langle x(t)\rangle$ and $\langle y(t)\rangle$. Finally, Figure 9 shows the percent error as a function of the number of the realizations. The error decreases starting from $10\%$ (blu line, 10k realizations) to $<1\%$ (yellow line, 200k realizations)

4. Diffusion

In this section, we will evaluate $\langle r^2\rangle =\langle x^2+y^2\rangle$ in the Poissonian and non-Poisson case. Using Equations (11) and (12) we have

$$\begin{array}{ccc}\hfill & & x^2(t)={\int }_0^t{\int }_0^t{v}_x(t_1)v_x(t_2)d{t}_{1}d{t}_2=v_0^2{\int }_0^t{\int }_0^{t}sin{\varphi }_{t_1}sin{\varphi }_{t_2}d{t}_{1}d{t}_2\hfill \end{array}$$

(40)

$$\begin{array}{ccc}\hfill & & y^2(t)={\int }_0^t{\int }_0^t{v}_y(t_1)v_y(t_2)d{t}_{1}d{t}_2=v_0^2{\int }_0^t{\int }_0^{t}cos{\varphi }_{t_1}cos{\varphi }_{t_2}d{t}_{1}d{t}_2\hfill \end{array}$$

(41)

where we set

$${\varphi }_t={\int }_0^t\omega (u)du={\omega }_{0}t+{\omega }_1{\int }_0^t\xi (u)du.$$

Consequently

$$x^2+y^2=r^2=v_0^2{\int }_0^t{\int }_0^{t}cos({\varphi }_{t_1}-{\varphi }_{t_2})d{t}_{1}d{t}_2.$$

(42)

Our goal is to evaluate the quantity $\langle r^2\rangle$. For the sake of compactness let us define the complex quantity

$$r_2(t)=v_0^2{\int }_0^t{\int }_0^{t}exp\left[i({\varphi }_{t_1}-{\varphi }_{t_2})\right]d{t}_{1}d{t}_2,$$

(43)

and take its time derivative

$$\begin{array}{c}\hfill \frac{\partial }{\partial t}{r}_2(t)=2v_0^2{\int }_0^{t}exp\left[i{\omega }_0(t-t_1)\right]〈exp\left[i{\omega }_1{\int }_{t_1}^t\xi (u)du\right]〉d{t}_1.\end{array}$$

(44)

For a Poisson process, with exponential waiting times distribution and correlation, the distribution of the first observed jump/event is the same as that of any other following event [24], this means that when averaging over the fluctuations, shifting the time origin, will not affect the result. For a generic non-Poissonian process, but with a finite time scale, this remains a good approximation in the long-time limit, so that we may re-write Equation (44) as

$$\begin{array}{c}\hfill \frac{\partial }{\partial t}{r}_2(t)\approx 2v_0^2{\int }_0^{t}exp\left[i{\omega }_0(t-t_1)\right]〈exp\left[i{\omega }_1{\int }_0^{t-t_1}\xi (u)du\right]〉d{t}_1.\end{array}$$

(45)

Formally, for $t\to \infty$, Equation (45) is the Laplace transform of the averaged function where $s\to i{\omega }_0$. Using the result of Equation (19), we may write

$$\begin{array}{ccc}\hfill & & \frac{\partial }{\partial t}{r}_2(t)\approx v_0^2\left[\frac{\left(1+\widehat{\psi }\left(i{\omega }_0-i{\omega }_1\right)\right)\widehat{\mathsf{\Psi }}\left(i{\omega }_0+i{\omega }_1\right)}{1-\widehat{\psi }\left(i{\omega }_0-i{\omega }_1\right)\widehat{\psi }\left(i{\omega }_0+i{\omega }_1\right)}+\right.\hfill \\ \hfill & & \left.\frac{\left(1+\widehat{\psi }\left(i{\omega }_0+i{\omega }_1\right)\right)\widehat{\mathsf{\Psi }}\left(i{\omega }_0-i{\omega }_1\right)}{1-\widehat{\psi }\left(i{\omega }_0-i{\omega }_1\right)\widehat{\psi }\left(i{\omega }_0+i{\omega }_1\right)}\right],\hfill \end{array}$$

(46)

namely a constant. We deduce that for the Poisson and the non-Poisson case, but with $1<\alpha <2$, we have ordinary diffusion, i.e.,

$$\begin{array}{c}\hfill \langle r^2\rangle =\mathrm{Re}[r_2(t)]\propto t\end{array}$$

(47)

Please note that in the Poisson case the approximation (45) is actually an exact expression and $〈exp\left[i{\omega }_1{\int }_0^{t-t_1}\xi (u)du\right]〉$ is given by Equation (21). The integration of (45) and the subsequent integration to obtain $\langle r^2\rangle$ does not present difficulties being the integral functions exponential functions. Neglecting the transient due to the exponentials with negative real part, the asymptotic expression for $\langle r^2\rangle$ reads as

$$\langle r^2\rangle \approx \frac{2\gamma v_0^2{\omega }_1^2}{{\gamma }^2{\omega }_0^2+{\left({\omega }_0^2-{\omega }_1^2\right)}^2}t,t\to \infty .$$

(48)

The diffusion is faster at resonance, when ${\omega }_0={\omega }_1$. The numerical check for the quantity $\langle r^2\rangle$, obtained integrating the uniform circular motion between switches of the magnetic field value, is shown in Figure 10 (Poisson case) and Figure 11 ($1<\alpha <2$). The agreement with the analytical result, Equation (47), is remarkable.

When $0<\alpha <1$ there is not a finite time scale, and we cannot use the approximation given by Equation (45). It is important to notice that while in the Poissonian case the distribution of the first observed event/jump coincides with the distribution $\psi (t)$ of any other event, in the non-Poissonian case it will be different from $\psi (t)$ and it will be function also of the time $t_a<0$, the time at which the system is prepared, in other words, the distribution of the first observed event will have a two-times dependence $\psi (t,t_a)$. The distribution of the first jump must be considered, and we are forced to use the full formula of Ref. [24], i.e.,

$$\begin{array}{ccc}\hfill & & \mathcal{L}\left[{\int }_0^{t}exp\left[i{\omega }_0(t-t_1)\right]〈exp\left[i{\omega }_1{\int }_{t_1}^t\xi (u)du\right]〉d{t}_1\right]=\hfill \\ \hfill & & \frac{1}{2}\frac{\left[{\widehat{f}}_{-}(s)+{\widehat{f}}_{+}(s)\widehat{\psi }\left(s-i{\delta }_{-}\right)\right]\widehat{\mathsf{\Psi }}\left(s-i{\delta }_{+}\right)}{1-\widehat{\psi }\left(s-i{\delta }_{+}\right)\widehat{\psi }\left(s-i{\delta }_{-}\right)}+\frac{1}{2}{\widehat{F}}_{+}\left(s\right)+\hfill \\ \hfill & & \frac{1}{2}\frac{\left[{\widehat{f}}_{+}(s)+{\widehat{f}}_{-}(s)\widehat{\psi }\left(s-i{\delta }_{+}\right)\right]\widehat{\mathsf{\Psi }}\left(s-i{\delta }_{-}\right)}{1-\widehat{\psi }\left(s-i{\delta }_{+}\right)\widehat{\psi }\left(s-i{\delta }_{-}\right)}+\frac{1}{2}{\widehat{F}}_{-}\left(s\right)\hfill \end{array}$$

(49)

where ${\delta }_{\pm }={\omega }_0\pm {\omega }_1$, and

$${\widehat{f}}_{\pm }(s)=\mathcal{L}\left[{\int }_0^{t}f(t-t_1,t_1)exp[i{\delta }_{\pm }(t-t_1)]d{t}_1\right]=\frac{\widehat{\psi }\left(s\right)-\widehat{\psi }\left(s-i{\delta }_{\pm }\right)}{-i{\delta }_{\pm }\left[1-\widehat{\psi }\left(s-i{\delta }_{+}\right)\right]}$$

(50)

being $f(\tau ,t_1)$ the conditional probability density that, fixed at time $t_1$, the first next switching event of the variable $\xi (t)$ occurs at time $t_1+\tau$. It is important to notice that differently from the Poisson case, this distribution is different from the distribution $\psi (t)$ of any other event. It coincides with $\psi$ only for the Poisson case. Analogously

$${\widehat{F}}_{\pm }(s)=\mathcal{L}\left[{\int }_0^{t}F(t-t_1,t_1)exp[i{\delta }_{\pm }(t-t_1)]d{t}_1\right]=\frac{1/s-{\widehat{f}}_{\pm }(s)}{s-i{\delta }_{\pm }}$$

(51)

where

$$F(\tau ,t_1)=1-{\int }_0^{\tau }f({\tau }^{\prime },t_1)d{\tau }^{\prime }$$

is the conditional probability that, fixed $t_1$, no switch occurs between $t=t_1$ and $t=t_1+\tau$. Using distribution (27) or in alternative the Mittag–Leffler distribution, we may write for the inverse Laplace transform of ${\widehat{f}}_{\pm }(s)$ the following

$$\begin{array}{ccc}\hfill & & f_{\pm }(t)={\mathcal{L}}^{-1}\left[{\widehat{f}}_{\pm }(s)\right]=-\frac{1}{{\delta }_{\pm }}\frac{d}{dt}{E}_{\alpha }\left[-{\left(\frac{t}{T}\right)}^{\alpha }\right]-\frac{1}{{\delta }_{\pm }T}\frac{1}{\mathsf{\Gamma }(\alpha )}{\left(\frac{t}{T}\right)}^{\alpha -1}-\hfill \\ \hfill & & \frac{1}{{\delta }_{\pm }T\mathsf{\Gamma }(\alpha )}{\int }_0^{t}exp\left[i{\delta }_{\pm }(t-t_1)\right]{\left(\frac{t-t_1}{T}\right)}^{\alpha -1}\frac{d}{d{t}_1}{E}_{\alpha }\left[-{\left(\frac{t_1}{T}\right)}^{\alpha }\right]d{t}_1\hfill \end{array}$$

(52)

From Equation (52), one can also understand why the contribution of the first jump/ event becomes crucial for a non-Poissonian process for $0<\alpha <1$. In this regime, differently from the Poisson regime, not only is such distribution different from that of any following event [24], but it also becomes dominant asymptotically. Analyzing Equation (51) from an asymptotic point of view, namely $sT\to 0$ and $\delta T\ll 1$ we deduce that contribution to the diffusion is given by $F_{\pm }$.

$$\begin{array}{c}\hfill F_{\pm }(t)=\frac{i}{{\delta }_{\pm }}-\frac{iexp[i{\delta }_{\pm }t]}{{\delta }_{\pm }}+{\int }_0^{t}exp[i{\delta }_{\pm }(t-t_1)]f(t_1)d{t}_1\end{array}$$

(53)

As $f_{\pm }(t)$ is a function decaying with time we may approximate $F_{\pm }(t)$ for $t\to \infty$ (see Ref. [27] for more details)

$$\begin{array}{c}\hfill F_{\pm }(t)\approx \frac{i}{{\delta }_{\pm }}-\frac{iexp[i{\delta }_{\pm }t]}{{\delta }_{\pm }}+exp[i{\delta }_{\pm }t]{\widehat{f}}_{\pm }(i{\delta }_{\pm }).\end{array}$$

(54)

The contribution of $F_{\pm }(t)$ to $\langle r^2(t)\rangle$ generates constant and oscillating terms while the dominant term is a power law. Indeed, we have

$$\begin{array}{c}\hfill \langle r^2(t)\rangle \sim 2v_0^{2}Re[A]\frac{1}{\mathsf{\Gamma }[\alpha ]\mathsf{\Gamma }[\alpha -1]}{\left(\frac{t}{T}\right)}^{\alpha },\mathrm{for}t\to \infty .\end{array}$$

(55)

The analytical result is obtained expanding for small s the Laplace transform of the expression for the derivative of $\langle r^2(t)\rangle$, neglecting the terms ${\widehat{F}}_{\pm }(s)$ for the reason stated above. The following expression gives the coefficient A

$$A=\frac{1}{1-{\widehat{\psi }}_{+}{\widehat{\psi }}_{-}}\left[\frac{{\widehat{\mathsf{\Psi }}}_{+}(1-{\widehat{\psi }}_{+}){\widehat{\psi }}_{-}}{i{\delta }_{+}}+\frac{{\widehat{\mathsf{\Psi }}}_{+}(1-{\widehat{\psi }}_{-})}{i{\delta }_{-}}\right]-\frac{1-{\widehat{\psi }}_{+}}{{\delta }_{+}^2}+\left(+\leftrightarrow -\right)$$

(56)

with ${\widehat{\mathsf{\Psi }}}_{\pm }$ and ${\widehat{\psi }}_{\pm }$ the Laplace transform of the respective functions $\mathsf{\Psi }(t)$ and $\psi (t)$ evaluated in $s=-i{\delta }_{\pm }$ and the last term within round brackets is obtained from previous terms by exchanging + and − subscripts. The numerical check is shown in Figure 11.

5. Results

To summarize, we have introduced a complete framework to describe the motion of a charged particle in a fluctuating magnetic field. We considered both ergodic and non-ergodic fluctuations of the magnetic field. We find that in the case of ergodic fluctuations the diffusion is asymptotically normal, while for non-ergodic fluctuations, we find anomalous diffusion properties. The diffusion is characterized by the mean-squared displacement, which we derived analytically and confirmed numerically for both regimes as illustrated in Figure 10 and Figure 11. In the case of non-Poisson fluctuations, we provide analytical formulae in Laplace transform from which we extract the asymptotic time behavior for the mean-squared displacement.

6. Discussion and Conclusions

The problem of the motion of a charged particle in fluctuating magnetic field has not been investigated when fluctuations of the field have non-ordinary statistical properties despite its physical interest and possible applications e.g., to plasma models. This paper starts filling this gap by considering dichotomic fluctuations with both Poisson and non-Poisson properties. For the Poisson case, we find the driven equation for the probability density. On the contrary, the equation for probability density function in the non-ordinary statistics case is still an open problem. We developed a theoretical framework for the first and second moment of the particle position and afforded both an analytical and numerical description. We neglected in this framework the effect of induced electric field due to the variations of the magnetic field, which we leave out for further investigation in an upcoming publication. Interestingly, we find that diffusion is normal either when fluctuations are Poisson or non-Poisson but with a finite mean time, differently from the standard case of continuous-time random walk (CTRW) which shows anomalous diffusion behavior for power-law distribution with finite mean time (either in the velocity or jump model) [33]. When the fluctuations time scale diverges, i.e., for non-ergodic fluctuations, an anomalous diffusion regime emerges, again differently from standard CTRW where a ballistic regime applies for power-law distributions with diverging mean time. This difference is related to the fact that the length of the CTRW jumps is proportional to the elapsed time while in the case studied in the paper, the particle is forced to move along a circular trajectory. Thus, the jumps cannot exceed the radius of the circumference, causing the relevant differences between the two diffusion processes.