Ultraslow scaled Brownian motion

We define and study in detail \emph{utraslow scaled Brownian motion (USBM)\/} characterised by a time dependent diffusion coefficient of the form $D(t)\simeq 1/t$. For unconfined motion the mean squared displacement (MSD) of USBM exhibits an ultraslow, logarithmic growth as function of time, in contrast to the conventional scaled Brownian motion. In an harmonic potential the MSD of USBM does not saturate but asymptotically decays inverse-proportionally to time, reflecting the highly non-stationary character of the process. We show that the process is weakly non-ergodic in the sense that the time averaged MSD does not converge to the regular MSD even at long times, and for unconfined motion combines a linear lag time dependence with a logarithmic term. The weakly non-ergodic behaviour is quantified in terms of the ergodicity breaking parameter. The USBM process is also shown to be ageing: observables of the system depend on the time gap between initiation of the test particle and start of the measurement of its motion. Our analytical results are shown to agree excellently with extensive computer simulations.


Introduction
In the wake of the development of modern particle tracking techniques strong deviations of the time dependence of the mean squared displacement (MSD) from the linear law x 2 (t) t derived by Einstein [1] and Smoluchowski [2] have been observed in a variety of complex fluidic environments [3,4,5,6,7]. Typically, anomalous diffusion of the power-law form is observed, where, depending on the value of the anomalous diffusion exponent α, we distinguish subdiffision with 0 < α < 1 and superdiffusion with α > 1 [8,9]. Accordingly, subdiffusion was observed in the cytoplasm of living cells [10,11], in artificially crowded liquids [12,13], and in structured or functionalised environments [14]. Also superdiffusive motion was found in living cells [15,16].
Recently, interest in ultraslow diffusion processes with the logarithmic form x 2 (t) log γ (t) (2) of the MSD with different values for the exponent γ has been revived [6]. Ultraslow diffusion may be generated by periodically iterated maps [17] and observed for random walks on bundled structures [18]. A prototype model for ultraslow diffusion is provided by Sinai diffusion in quenched landscapes with random force field, for which γ = 4 [19,20,21,22]. In the context of Sinai diffusion ultraslow continuous time random walks with super heavy-tailed waiting times with γ > 0 [22,23,24,25] were discussed. Ultraslow scaling of the MSD of the form (2) were obtained in aperiodic environments (variable γ) [26] and vacancy induced motion (γ = 1) [27]. Moreover, it occurs in heterogeneous diffusion processes with exponentially varying diffusivity (γ = 2) [28], or interacting many-body systems in low dimensional disordered environments with γ = 1/2 [29], the dynamics of the latter being governed by an ultraslow, ageing counting processes [30]. The logarithmic time dependence (2) with γ = 1 of the MSD is also observed for the self diffusion of particles in free cooling granular gases with constant, subunity restitution coefficient in the homogeneous cooling state [31]. Granular gases are rarefied granular systems, in which particles move along ballistic trajectories between instantaneous collisions [31]. They are common in Space, for instance, in protoplanetary discs, interstellar clouds and planetary rings [32]. At terrestrial conditions granular gases may be obtained by placing granular matter into containers with vibrating [33] or rotating [34] walls. If no net external forces (gravitation, etc.) are acting on the granular system, the motion of granular particles gradually slows down due to dissipative collisions between them [31]. This microgravity condition can be achieved, inter alia, with parabolic airplane flights or satellites [35,36,37] or by the use of diamagnetic levitation [38]. We note that in very dense two-dimensional lattice gas systems, ultraslow diffusion emerges, as well [39]. Figure 1 shows the crossover from the ballistic to the ultraslow form (2) of the MSD of a granular gas with constant restitution coefficient ε = 0.8 in the homogeneous cooling state. Haff's law demonstrates that the kinetic temperature of such a free granular gas with constant restitution coefficient decays inverse-proportionally with time, T (t) 1/t [40]. For the effective self diffusion of the gas particles-mediated by particle-particle collisions-this property translates into the time dependent diffusion coefficient D(t) 1/t [41,42,43]. We note that a diffusivity of the form D(t) = D 0 + D 1 /t with a component decaying inverse-proportionally with time was used in the modelling of the motion of molecules in porous environments [44] as well as of water diffusion in brain tissue measured by magnetic resonance imaging [45].
Here we study in detail the process of ultraslow scaled Brownian motion (USBM) with time dependent diffusion coefficient D(t) 1/t. Starting from the Langevin equation for USBM and a summary of the simulations procedure we present analytical and numerical results for the MSD and the time averaged MSD for the cases of Figure 1. Time dependence of the ensemble averaged MSD x 2 (t) obtained from event driven molecular dynamics simulations of three-dimensional force-free granular gases [43]. At short times the particles follow ballistic trajectories, while for longer times the ensemble averaged MSD has a logarithmic time dependence. The inset focuses on the logarithmic long time behaviour. unconfined (Section 2) and confined (Section 3) motion. We analyse in detail the disparity between the ensemble and time averaged MSD and quantify the statistical scatter of the amplitude of the time averaged MSD of individual realisations of the USBM process. Moreover we study the ageing properties of USBM, that is, the explicit dependence of the physical observables on the time difference between the initiation of the system and the start of the observation. In Section 4 we present our Conclusions. In the Appendix we present details of the calculation of higher order moments and the ergodicity breaking parameter.

Overdamped Langevin equation for ultraslow scaled Brownian motion
Anomalous diffusion processes with power-law form (1) of the MSD are often modelled in terms of scaled Brownian motion (SBM) characterised by an explicitly time dependent diffusivity of the power-law form D(t) t α−1 with 0 < α < 2, see, for instance, references [46,47,48,49,50,51] as well as the study by Saxton [52] and further references Figure 2. Schematic of the motion of a Brownian particle in a bath with decreasing temperature T (t) t 2α−2 for 0 ≤ α < 1. The diffusion coefficient of the Brownian particle decays with time as D(t) t α−1 . USBM corresponds to the case α = 0, while standard SBM is strictly limited to 0 < α < 2 [54].
therein. In SBM this form of D(t) is combined with the regular Langevin equation [53] in which ζ(t) represents white Gaussian noise with the normalised covariance and zero mean ζ(t) = 0. While for a system connected to a thermal reservoir a description in terms of a time dependent temperature underlying SBM is unphysical [54], time dependent diffusion coefficients appear naturally in systems that are open or dissipate energy into other degrees of freedom such as the granular gases discussed above, see the schematic in figure 2. In fact, granular gases with a viscoelastic, relative particle speed-dependent restitution coefficient correspond to SBM with α = 1/6 [31,43]. Diffusion in media with explicitly time dependent temperature can, for instance, also be observed in snow melt dynamics [55,56]. A diffusion equation with a time dependent diffusivity proportional to t 2 was originally introduced by Batchelor [57] to describe the anomalous Richardson relative diffusion [58] in turbulent atmospheric systems. SBM with diffusivity D(t) t α−1 was studied extensively during the last few years [59,60,61,54,62]. In particular, the weakly non-ergodic disparity between ensemble and time averages in SBM as well as its ageing behaviour were analysed [60,61,54,62], see also below. Processes with both time and position dependent diffusion coefficients were also reported [63]. SBM is a Markovian process with stationary increments ζ(t), however, it is rendered non-stationary by the time dependence of the coefficient D(t). SBM is therefore fundamentally different [6,54] from seemingly similar processes such as fractional Brownian motion or fractional Langevin equation motion [64].
Following the motivation from our studies of granular gases with constant restitution coefficient [43] we here consider USBM with the time dependent diffusion coefficient The time scale τ 0 defines the characteristic time beyond which the long time scaling D(t) ∼ D 0 τ 0 /t sets in. We here introduce τ 0 to avoid a divergence of D(t) at t = 0. The case (5) is explicitly excluded in the allowed range for the scaling exponent α in SBM and, as we will see, constitutes a new class of stochastic processes. In the following we solve the overdamped Langevin equation (3) with the time dependent diffusion coefficient (5) analytically and perform extensive computer simulations of the corresponding finitedifference analogue of the Langevin equation. In this procedure, at each time step the increment of the particle position takes on the value where W i+1 − W i is the increment of the standard Wiener process and D(i) is the value of the time dependent diffusivity (5) at the time instant i. We simulated N = 10 3 independent particles (runs) with the parameters τ 0 = 1 and D 0 = 1/2 in all graphs presented below.

Ensemble and time averaged mean squared displacements
From direct integration of the Langevin equation (3) with the time dependent diffusivity (5) we find the ultraslow, logarithmic growth of the ensemble averaged MSD. USBM therefore reproduces the asymptotic behaviour of the MSD for granular gases in the homogeneous cooling state and with constant restitution coefficient [43], as shown in figure 1.
In addition to the ensemble averaged MSD x 2 (t) of the particle motion, it is often useful to compute the time averaged MSD Here, the lag time ∆ defines the width of the averaging window slid over the time series x(t) of the particle position of overall length t (the measurement time). Time averages of the form (8) are often used in experiments and large scale simulations studies based on single particle tracking approaches, in which typically few but long trajectories are available [10,11,65]. The careful analysis of the time averaged MSD (8) provides additional important information on the studied process as compared to the ensemble averaged MSD x 2 (t) , see, for instance, the analyses in references [11,65]. Often one takes the additional average over N individual particle traces δ 2 i (∆), For ergodic processes ‡ such as Brownian motion, fractional Brownian motion, and fractional Langevin equation motion the time averaged MSD converges to the ensemble averaged MSD in the limit of sufficiently long times, lim t→∞ δ 2 (∆) = x 2 (∆) [6]. This property is due to the stationarity of the increments of these processes [66]. The ergodic behaviour lim t→∞ δ 2 (∆) = x 2 (∆) of these processes holds for unconfined motion when the system is in fact out-of-equilibrium, an advantage of the particular definition (8).
Moreover, ergodic systems fulfil the equivalence even at finite t [6]. Systems in which we observe the disparity δ 2 (∆) = x 2 (∆) and therefore also lim t→∞ δ 2 (∆) = x 2 (∆) are called weakly non-ergodic [4,5,6,7,67,68]. § To calculate the time averaged MSD (9) for USBM we do not need to consider the mixed position autocorrelations in the definition of the time averaged MSD, as the expression in the angular brackets simplify as follows, This is due to the property for stochastic processes whose increments are independent random variables. We thus find the exact form for the time averaged MSD of USBM, ‡ We consider processes ergodic in the Boltzmann-Khinchin sense when the long time average of a physical observable converges to the associated time average. § Note that also transiently non-ergodic behaviour may become relevant as it may mask intrinsic relaxation times when time averages are measured [12,69].
In contrast, this is not valid in the case of granular gases, where particles move ballistically in between instantaneous collisions [43], or for processes driven by long-range correlated increments such as fractional Brownian motion or fractional Langevin equation motion [6,64,70].
where we introduced the auxiliary function The time averaged MSD (13) thus crosses over from the limiting behaviour at short lag times τ 0 ∆ t combining a linear with a logarithmic ∆ dependence, to the purely logarithmic law We see that as the lag time ∆ approaches the measurement time t, the time average MSD approaches the MSD (7), δ 2 (t) → x 2 (t) . The results of our simulations of the USBM process for both ensemble and time averaged MSDs agree very well with the above analytical results, as demonstrated in figure 3. In that plot the thin grey curves depict the simulations results for the time averaged MSD for individual trajectories. The amplitude spread between different trajectories is fairly small for ∆ t and increases when the lag time ∆ approaches the trace length t due to worsening statistics.

Stochasticity of the time averaged mean squared displacement and ergodicity breaking parameter
Even ergodic processes such as Brownian motion exhibit a certain degree of stochasticity of time averaged observables for shorter measurement times. The amplitude fluctuations at a given lag time ∆ of the time averaged MSD as compared to the trajectory average (9) is quantified in terms of the ergodicity breaking parameter [6,70,71,72] where in the second equality we introduced the relative deviation [71] ξ = δ 2 (∆) The necessary condition for ergodicity of a stochastic process is that the ergodicity breaking parameter vanishes in the limit of infinitely long trajectories. Brownian motion provides the basal level for the approach to ergodicity according to [70] Fractional Brownian motion and fractional Langevin equation motion are ergodic [64,70]. Weakly non-ergodic processes, which are characterised by the disparity δ 2 (∆) = x 2 (∆) [4,5,6,7,68,71] include continuous time random walks with scalefree distributions of waiting times [4,5,6,68,71] and heterogeneous diffusion processes [74,82]. In the limit of long traces, the value of their ergodicity breaking parameter remains finite, which is indicative of the intrinsic randomness of time averages of these processes. In contrast, the ergodicity breaking parameter for SBM vanishes in the limit of long trajectories [61]. The ergodicity breaking parameter for USBM is derived in the Appendix. The final expression in the relevant limit τ 0 ∆ t reads where the constant C = π 2 /6 − 1 0.645. Thus, the time averaged MSD for USBM becomes increasingly reproducible as the length of the time traces is extended, albeit the approach to zero is logarithmically slow. We demonstrate the functional form of the ergodicity breaking parameter as function of the lag time ∆ for two different measurement times and the approach of EB to its asymptotic behaviour (20) in figure  4.  The ergodicity breaking parameter quantifies the statistical spread of the time averaged MSD. An important indicator for different types of stochastic processes is also the complete distribution φ(ξ) [6,68,71,73]. As shown in figure 5 this distribution has an asymmetric bell-shaped curve approximately centred around the ergodic value ξ = 1. The tail at larger ξ values appears somewhat longer compared to the tail at shorter ξ. ¶ For longer lag times at fixed overall length t of the time series the width of the distribution φ(ξ) grows. This is consistent with the fact that at larger value of ∆/t the time averages become more random. In figure 5 we also show a fit to the function which appears to capture the functional behaviour reasonably well. We note that the shape of φ(ξ) appears narrower compared to the one of heterogeneous diffusion processes with power-law space dependence of the diffusivity [74] which was fitted by a three-parameter Gamma distribution [74,75]. In comparison, the distribution φ(ξ) for standard SBM is quite narrow, although it widens as the exponent α approaches zero and particularly as the lag time ∆ grows [54].

Ageing ultraslow scaled Brownian motion
For processes with stationary increments such as Brownian motion or fractional Brownian motion, if we initiate the system at t = 0 but start recording it only at some later time t a , the physical observables will not explicitly depend on the ageing time t a . + However, for several anomalous processes pronounced ageing effects are found. These include continuous time random walk processes with scale free distributions of waiting times [77,78], correlated continuous time random walks [79], non-linear maps generating subdiffusion [80], systems with annealed and quenched disorder [81], heterogeneous diffusion processes [82], or standard SBM [62]. In contrast to subdiffusive continuous time random walk processes, in which ageing emerges due to the divergence of a characteristic waiting time [77], in ultraslow SBM the non-stationarity of the system stems from the explicit time dependence of the diffusion coefficient. When the recording of the particle position starts at a finite time t a , this ageing time explicitly appears in the particle's MSD. For the aged MSD [77,6] in analogy to equation (7) we find that In the limit of strong ageing, t a t, this expression yields the linear scaling  of the MSD with time t, the ageing time t a rescaling the effective particle diffusivity. The transition between this ageing-dominated linear scaling for the MSD and the anomalous logarithmic time dependence in the weak ageing limit t t a is clearly seen in figure 6. For the aged time averaged MSD [6,77] we obtain the result where the auxiliary function (t) was defined in equation (14). In the limit τ 0 ∆ t and ∆ t a the aged time averaged MSD factories into a term containing all the information on the ageing and measurement times t a and t, and another capturing the physically relevant dependence on the lag time ∆ and the measurement time t, This factorisation is analogous to that of heterogeneous diffusion processes [82], scalefree subdiffusive continuous time random walks [77], and standard SBM [62]. However, in contrast to these processes the aged time averaged MSD for short lag times does not factorise into the product of the non-aged time averaged MSD (16) and a factor containing the ageing time. For strong ageing t a t we obtain the linear scaling In this limit, that is, the system becomes apparently ergodic and we observe the equality δ 2 a (∆, t a ) = x 2 a (∆, t a ) , as can be seen from comparison with equations (23) and (26). Figure 7 shows the convergence of the time averaged MSD to the limiting behaviour (26). Such a behaviour was previously observed for aged subdiffusive SBM [62], heterogeneous diffusion processes [82], and continuous time random walk processes [77]. In the case of USBM this phenomena has a clear physical explanation: at the beginning of the experiment the diffusion coefficient D(t) significantly decreases during the measurement time t τ 0 from D(0) = D 0 to D(t) ∼ D 0 τ 0 /t, and the system is strongly nonstationary. In contrast, after a long ageing period t a t the diffusion coefficient remains practically unchanged during the measurement time, D(t a + t) D(t a ) = D 0 τ 0 /t a . Figure 7 explicitly shows how the amplitude of the time averaged MSD is reduced due to ageing in the system. How do the fluctuations of individual time averaged MSD traces change in the presence of ageing? The derivation of the ergodicity breaking parameter for the aged process is provided in the Appendix. The final result in the limit ∆ t and ∆ t a assumes the form . In the strong ageing limit t a t the ergodicity breaking parameter EB a is independent of the ageing time t a , and it asymptotically converges to the result (19) of Brownian diffusion. Our theoretical results agree well with the simulations, as witnessed by figure 8. For weak ageing t a ∆, t the result (13) of the non-aged USBM process is recovered.

Confined ultraslow scaled Brownian motion
The motion of particles in external confinement is an important physical concept for applications of stochastic processes, and it is also relevant from an experimental point of view. Namely, the motion of particles in cells may repeatedly hit the cell wall, or the tracer particles may experience a restoring force in particle tracing experiments by help of optical tweezers. Here we consider the generic case of confinement in an harmonic potential. USBM in the presence of such a linear restoring force is governed by the overdamped Langevin equation with additional Hookean force term −kx,

Ensemble and time averaged mean squared displacements
The ensemble averaged MSD follows directly from this stochastic equation, and we obtain Here we defined the auxiliary function where in the second line we used the definition of the exponential integral The asymptotic behaviour of the MSD for long times t 1/k has the time dependence Reflecting the temporal decay of the temperature encoded in the time dependent diffusion coefficient (5) we observe the 1/t scaling of the MSD in confinement. This underlines the highly non-stationary and athermal character of this process [60,54,62]. The time averaged MSD for confined USBM is obtained from the relation The covariance of the position for ultraslow SBM in confinement can no longer be simplified according to equation (12) but has the time dependence Introducing relations (29) and (34) into equation (33) we obtain For long times and strong external confinement, {t, t a , ∆} {1/k, τ 0 } this expression simplifies to The time averaged MSD has a pronounced plateau for ∆ t, t a ,  (29), while the blue line denotes equation (35). The red line shows the asymptotic behaviour (36), and the horizontal dashed line the leading term (37). The symbols correspond to the simulations of equation (28).
that is, in this regime the time averaged MSD is independent of the lag time, compare the discussion in references [54,62]. Simulations based on the Langevin equation with the Hookean forcing are in excellent agreement with these analytical results, as shown in figure 9.

Ageing ultraslow scaled Brownian motion in confinement
3.2.1. Ensemble averaged mean squared displacement. For confined ageing USBM, in which we measure the MSD starting from the ageing time t a until time t, the result for the MSD becomes where E (x) is defined in equation (30). Expression (38) reduces to equation (29) for vanishing ageing, t a = 0. However, even in the presence of weak ageing, t a 1/k, at long times t 1/k the behaviour of the MSD reads Figure 10. Ensemble averaged MSD x 2 a (t, t a ) for confined ageing USBM at different ageing times: t a = 0 (no ageing, black line), t a = 0.1 (weak ageing, red line), and t a = 10 5 (strong ageing, blue line). Note that for better visibility the curve for t a = 10 5 was multiplied by a factor of 10 3 .
contrasting the behaviour in equation (32). The ensemble averaged MSD for ageing USBM at different ageing times is depicted in figure 10. At short times t < 1/k the weakly aged MSD follows the non-aged behaviour. Eventually it attains the plateau given by the first term in equation (39), instead of decaying towards zero as in the non-aged case. In the analysis of experimental data times the exact moment of the system's initiation may often not be known, for instance, when measuring biological cells. The apparent plateau revealed here for confined ageing USBM dynamics may thus erroneously be mistaken as a signature of a stationary process.
Expanding the exponential integral in equation (38), in the strong ageing limit t a {τ 0 , 1/k} we find For t 1/k we recover the unconfined result (23). In the opposite limit t 1/k the behaviour of equation (40) crosses over to In this case we recover a transition between two plateaus, as it was observed for subdiffusive SBM [62]. Namely, for short measurement times t t a we find from Figure 11. Ensemble and time averaged MSDs x 2 a (t, t a ) and δ 2 a (∆) for confined ageing USBM. The symbols depict simulations of equation (28). The blue line corresponds to the theoretical result (44), and the red line shows the asymptotic (45). The horizontal dashed line shows the leading term (46).
while at long measurement times t t a this turns to This behaviour, which appears unique for USBM and SBM, is depicted in figure 10).

Time averaged mean squared displacement.
The time averaged MSD for ageing confined USBM is derived analogously to the non-aged case, yielding Figure 12. Ergodicity breaking parameter EB as function of t/∆ in the non-aged (t a = 0) and aged (t a = 3000) cases.
In the limit of strong confinement 1/k {t a , t, ∆} this expression can be significantly simplified to obtain For ∆ t, t a we again find an apparent plateau, In the case of strong ageing t a t we find Comparison to equation (42) shows that the time averaged MSD becomes equal to the ensemble MSD in this strong ageing regime, and ergodicity is apparently restored as in the unconfined case. The behaviour of the ensemble and time averaged MSDs for confined ageing USBM are depicted in figure 11. The ergodicity breaking parameter EB for confined USBM is depicted in figure 12 for both absence and presence of ageing. It is a decreasing function of the ratio t/∆ for large t/∆, while at small values of t/∆ it remains practically unchanged.

Conclusions
We proposed and studied ultraslow scaled Brownian motion, a new anomalous stochastic process with a time dependent diffusion coefficient of the form D(t) 1/t. Formally USBM corresponds to the lower bound α = 0 of scaled Brownian motion with diffusivity D(t) t α−1 (0 < α) [59,60,54,61,62], yet its dynamical behaviour is significantly different. We showed that USBM yields a logarithmic time dependence of the MSD rather than the power-law scaling of SBM. USBM's time averaged MSD was shown to acquire a combination of power-law and logarithmic lag time dependence. USBM is weakly non-ergodic and ageing. The ergodicity breaking parameter quantifying the random character of time averages of the MSD has a weak logarithmic dependence on the ratio ∆/t of lag time ∆ and length t of the recorded trajectories, tending to zero in the limit of infinitely long traces and/or short lag times. In the case of strong ageing the system tends to usual Brownian motion and the behaviour of the system becomes apparently ergodic. Under external confinement the behaviour of the USBM dynamics exhibits an apparent plateau for the time averaged MSD, while the ensemble averaged MSD decays proportionally to 1/t at longer times, reflecting the highly non-stationary character of USBM. Ageing produces an apparent plateau for the ensemble averaged MSD and a crossover between two plateaus for the time averaged MSD. USBM adds to the rich variety of ultraslow processes with logarithmic growth of the ensemble averaged MSD yet displays several unique features in comparison to other ultraslow processes.
Potential applications of USBM are foremost in the description of random particle motion in intrinsically non-equilibrium system such as free cooling granular gases or systems coupled to explicitly time dependent thermal reservoirs. On a more general level we hope that the discussion of ultraslow processes will lead to a rethinking of claims in diffusion studies that certain particles appear immobile. Namely, one often observes a population splitting into a (growing) fraction of immobile particles and another fraction of particles performing anomalous diffusion of the form (1) [83]. Ageing continuous time random walks [77] or heterogeneous diffusion processes [28,82] give rise to such a behaviour. However, given the tools provided here on ultraslow diffusion it might be worthwhile checking whether the observe "immobile" particles may in fact perform logarithmically slow diffusion.
In the case of ageing EB a (∆) = lim Expanding the integrand for t 1 > t a 1, we get Evaluating the integral for t a ∆, t ∆, we obtain , (A. 12) and the ergodicity breaking parameter EB a (∆) is then given by equation (27).