Disorder-Induced Anomalous Mobility Enhancement in Confined Geometries

Transport in disordered and amorphous materials has attracted vast attention for many decades [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]. The study of systems’ response to external forces, particularly with an aim to optimize transport, constitutes an imperative focus of research [11]. The external force can be a result of an electric field pulling on an electron through a conductor [12, 13, 14] or a pressure gradient pushing on a molecule diffusing in a channel [15, 16, 17, 18]. The classical depiction of such dynamics is Drude’s model of current flow in a metal [19]. It describes, through the application of kinetic theory, the diffusion of electrons by repeated encounters with immobile hard scatterers (such as ions or impurities). When an external electric field is applied, the charge carriers experience a net drift velocity related to the mean free path between scattering events. The resulting picture is that the response to the external force, i.e., carrier mobility, is an intrinsic property of the medium. Therefore, for transport in restricted geometry, like a channel, the expectation is that the mobility will be independent of channel width (or cross-section) or sometimes will increase with channel width due to the availability of new pathways. In this work, we explore the mobility properties of transport inside a channel with the presence of strong and quenched disorder. Specifically, we aim to demonstrate that quenched and strong disorder can redefine our understanding of the dependence of mobility on geometry.

Experiments in amorphous materials [20, 21, 22, 23, 24] have shown that a packet of charge carries do not propagate in a Gaussian manner and instead exhibit a dispersion of carrier transit times. Scher and Montroll [25, 5] termed this phenomenon anomalous transport and suggested that carriers are affected by deep traps or local areas of arrest. When the duration $\tau$ of such events follows a power-law distribution, i.e. $\sim\tau^{-1-\alpha}$, and $0<\alpha<1$, the transport becomes subdiffusive [3]. Meaning that the mean squared displacement (MSD) is not proportional to time $t$ but rather grows sublinearly, i.e. MSD$\sim t^{\alpha}$, as observed in amorphous materials [20, 21, 24, 22, 23, 26, 27, 28, 25], biological cells [29, 6, 30, 31, 32, 8, 33], granular materials [34, 35, 10], non-Newtonian fluids [36] and other systems [37, 38, 39]. These power-law distributed waiting times ($\sim\tau^{-(1+\alpha)}$), as detected in various systems [40, 36, 41, 42, 9], can appear naturally due to the exponential distribution of the depths of energetic wells that give rise to the regions of local arrest. The strong disorder ($0<\alpha<1$) results in a diverging mean trapping time that disrupts regular diffusive properties and leads to aging, weak ergodicity breaking, and non-self averaging [43, 44, 45, 46, 47, 48]. Most theoretical studies address the annealed version of the strong disorder. Namely, the waiting times in the trapping regions are uncorrelated, and each time the particle returns to the same arrest region, it is localized for a different time. Such framework was termed the continuous time random walk (CTRW) [25, 3, 6, 49], a very popular model of anomalous transport. The quenched version with a strong disorder, termed the quenched trap model (QTM), treats the trapping times during revisits of the arrest region as correlated. For the QTM regular techniques and Stokes-Einstein–Smoluchowski theory do not apply due to strong correlations and memory effects [3, 50, 51]. Scaling arguments and renormalization group approach [3, 52, 53] among other works [54, 55, 56, 57, 58] suggest that for dimensions $d>2$, QTM behaves qualitatively as CTRW in the subdiffusion regime. But big differences can be witnessed as we show.

In this work we explore the effects of a strong quenched disorder on particle mobility under the geometrical constraint of a channel with width $w$. By utilizing the recently developed double-subordination technique [59, 60, 61, 58, 62], we obtain an analytical expression for the mobility and its dependence on the external driving force $f$, temperature $T$ and width $w$. For low temperatures, we find that the mobility is a decreasing function of $w$. Namely, the response to an external drive weakens as the channel cross-section grows. Such a counterintuitive enhancement with decreasing $w$ appears only when the disorder is quenched and strong. When one of these requirements is omitted, the mobility is independent of the channel width.

The quenched trap model. The physical picture behind QTM is a thermally activated particle jumping between energetic traps. When a particle is in a trap, the average escape time $\tau$ is provided by the Arrhenius law $\tau\propto\exp\left(E_{\mathbf{r}}/T\right)$, where $E_{\mathbf{r}}>0$ is the depth of trap at position $\mathbf{r}$ and $T$ is the temperature. When the distribution of the energetic traps $E_{\mathbf{r}}$ is exponential, $\phi(E_{\mathbf{r}})=\exp(-E_{\mathbf{r}}/T_{g})/T_{g}$, the distribution of the average escape time is

$A=|\Gamma(-T/T_{g})|T/T_{g}$ and $\Gamma(\dots)$ is the Gamma function. For $T<T_{g}$, the slow power-law decay of $\psi(\tau)$ leads to a diverging mean escape, when averaged over disorder. In the following, we set $\alpha=T/T_{g}$ and focus on the glassy regime $0<\alpha<1$, where QTM exhibits aging and non-self-averaging behavior [50, 51, 53]. The average escape times $\tau_{\mathbf{r}}$ serve in QTM as the waiting times. Each time the particle visits position $\mathbf{r}$, it spends there exactly the same time $\tau_{\mathbf{r}}$ hence the disorder is quenched. In [47], no difference was found between setting quenched waiting times or setting quenched average waiting times. The quenched variables $\{\tau_{\mathbf{r}}\}$ are positive, independent, identically distributed (i.i.d) random variables with probability density function (PDF) provided by Eq. (1). We consider the spatial process between different positions $\mathbf{r}$ as a random hop process on a two-dimensional square lattice with lattice spacing $a$ taken to be $1$ (a.u). At time $t=0$, the particle starts at $\mathbf{r}=0$ and stays at this position for the period $\tau_{\mathbf{0}}$ before jumping to some random site ${\mathbf{r}}^{\prime}$ where it waits for the period $\tau_{\mathbf{r}^{\prime}}$ and then the random jump + waiting period procedure continues. The probability of transition (jump) from $\mathbf{r}$ to $\mathbf{r^{\prime}}$ is provided by $p(\mathbf{r}^{\prime};\mathbf{r})$. We assume that the lattice is translationally invariant in space, i.e., $p(\mathbf{r}^{\prime};\mathbf{r})$ is a function of $\mathbf{r}^{\prime}-\mathbf{r}$: $p(\mathbf{r}^{\prime}-\mathbf{r})$. The disorder averaged positional PDF of finding the particle at position $\mathbf{r}$ at time $t$, $\langle P(\mathbf{r},t)\rangle$ ($\langle\cdots\rangle$ represents the averaging over disorder) is found by utilizing the double subordination technique [59, 60] that we briefly describe below.

The diffusion front. The effect of correlations imposed by quenched disorder is appreciated when the measurement time $t$ is written in terms of the local waiting times $\tau_{\mathbf{r}}$. Namely, $t=\sum_{\mathbf{r}}n_{\mathbf{r}}\tau_{\mathbf{r}}$, where $n_{\mathbf{r}}$ is the number of visits to $\mathbf{r}$ up to time $t$. Although all the different $\tau_{\mathbf{r}}$ are i.i.d, the $\{n_{\mathbf{r}}\}$ are correlated, like in our case of nearest-neighbor hopping on a lattice where $n_{\mathbf{r}}$ is very similar to the nearest-neighbour $n_{\mathbf{r^{\prime}}}$. By fixing the values of $\{n_{\mathbf{r}}\}$ and averaging over $\{\tau_{\mathbf{r}}\}$ (disorder averaging) the Laplace pair of $t$, i.e. $\langle e^{-ut}\rangle$, is $\sim e^{-AS_{\alpha}u^{\alpha}}$ where

Since the Laplace pair of the one-sided Lévy distribution $l_{\alpha,A,1}(\eta)$ is $\sim e^{-Au^{\alpha}}$ we obtain that $t\sim S_{\alpha}^{1/\alpha}\eta$, where $\eta$ is a random variable distributed according to $l_{\alpha,A,1}(\eta)$. This connection between $t$, $\eta$ and fixed $S_{\alpha}$ allows to obtain the PDF of $S_{\alpha}$ for fixed $t$, $\mathcal{N}_{t}\left(\mathcal{S}_{\alpha}\right)$, by changing variables from $\eta=t/{S_{\alpha}}^{1/\alpha}$ to $S_{\alpha}=(t/\eta)^{\alpha}$. Therefore, $\mathcal{N}_{t}\left(\mathcal{S}_{\alpha}\right)\sim\frac{t}{\alpha}\left(\mathcal{S}_{\alpha}\right)^{-1/\alpha-1}l_{\alpha,A,1}\left[\frac{t}{\left(\mathcal{S}_{\alpha}\right)^{1/\alpha}}\right]$. The explicit form of $\mathcal{N}_{t}\left(\mathcal{S}_{\alpha}\right)$ allows performing the first, what is commonly called, subordination [63] and express the disorder averaged $\langle P({\mathbf{r}},t)\rangle$ by using $S_{\alpha}$ as the local time of the process. Namely, for the conditional probabilty $P_{S_{\alpha}}({\mathbf{r}})$ of finding the particle at position $\mathbf{r}$ for a given $S_{\alpha}$ (i.e., at “time” $S_{\alpha}$), the law of total probability yields

The second subordination is applied to $P_{S_{\alpha}}(\mathbf{r})$ We use the number of jumps, $N$, to represent $P_{S_{\alpha}}(\mathbf{r})$ in terms of $W_{N}(r)$, the probability to find the particle at ${\mathbf{r}}$ after $N$ jumps, and $\mathcal{G}_{S_{\alpha},\mathbf{r}}(N)$ the probability of different values of $N$ for a prescribed $S_{\alpha}$ and $\mathbf{r}$. The law of total probability yields $P_{S_{\alpha}}(\mathbf{r})=\sum_{N=0}^{\infty}W_{N}(\mathbf{r})\mathcal{G}_{S_{\alpha},\mathbf{r}}(N)$ and then according to Eq. (3)

When $t\to\infty$, the probability of small $S_{\alpha}$ is negligible and for large $S_{\alpha}$, when the probability of eventual return to the origin $Q_{0}$ is $<1$, it was shown that [58], $\mathcal{G}_{S_{\alpha}}(N,\mathbf{r})\to\delta\left(S_{\alpha}-\Lambda N\right)$ where

$Li_{a}(b)=\sum_{j=1}^{\infty}b^{j}/j^{a}$ is the Polylogarithm function [64]. $Q_{0}$ is computed when the spatial process is treated as a function of $N$. Therefore, for QTM where the spatial process is defined solely by the jump probabilities $p(\mathbf{r}^{\prime}-\mathbf{r})$, Eq. (4) yields

Equation (6), first obtained in [58], presents the disorder averaged propagator of QTM in terms of the spatial process on a lattice as a function of $N$, and a transformation from $N$ to $t$. Two points are in place: (I) The distribution $W_{N}(r)$ is defined by the jump probabilities $p(\mathbf{r}^{\prime}-\mathbf{r})$ and found by the standard techniques for a random walk (RW) on a lattice [65]. (II) For $\Lambda=1$ Eq. (6) displays the propagator for the annealed version of the disorder (CTRW) [3]. Therefore, $\Lambda$, quantifies the difference between quenched and annealed disorder as a function of $Q_{0}$ (Eq. (5)) and depends on the geometry and the external force. Below, we utilize Eq. (6) and compute the response to external constant force $F$ acting on a single particle in a $2$-dimensional channel of width $w$. For this purpose, we first compute the average position along the longitudinal axis of the channel $\hat{x}$.

The average displacement. The motion occurs on top of a lattice in a two-dimensional channel and is unrestricted in the $\hat{x}$ direction. The width of the channel, in the $\hat{y}$ direction, is $w$, which is also the number of sites across $\hat{y}$ (the lattice spacing is $a=1$). Due to the translational invariance of the spatial process and transition probabilities $p(\mathbf{r}-\mathbf{r^{\prime}})$, we use periodic boundary conditions for $\hat{y}$. The case of reflecting boundary conditions will be addressed below. The strength of the force $f$, applied only along $\hat{x}$, is characterized by the dimensionless parameter $F=af/T$, where $k_{B}$ is set to $1$. The transition probabilities $p(\mathbf{r}-\mathbf{r^{\prime}})$ allow transitions only to the nearest neighbors on the square lattice. Namely, $p_{\rightarrow}$ ($p_{\leftarrow}$) is the probability for a single jump to the right (left) along $\hat{x}$ and $p_{\uparrow}$ ($p_{\downarrow}$) is the probability for a single jump up (down) along $\hat{y}$. The detailed balance condition dictates that $p_{\rightarrow}/p_{\leftarrow}=e^{F}$ and $p\uparrow/p_{\downarrow}=1$. Therefore, due to the normalization condition $p_{\uparrow}+p_{\downarrow}+p_{\leftarrow}+p_{\rightarrow}=1$, we obtain that $p_{\rightarrow}=Be^{F/2}$, $p_{\leftarrow}=Be^{-F/2}$ and $p_{\uparrow}=p_{\downarrow}=B=1/[2\cosh(F/2)+2]$. We are interested in the mean position $\langle x(t)\rangle=\sum_{\mathbf{r}}x\langle P(\mathbf{r},t)\rangle$. After one single jump the average displacement along $\hat{x}$ is $p_{\rightarrow}-p_{\leftarrow}=\tanh(F/4)$, therefore after $N$ jumps the average displacement is $\sum_{\mathbf{r}}xW_{N}(\mathbf{r})=N\tanh(F/4)$. Then according to Eq. (6) $\langle x(t)\rangle=\sum_{N}N\tanh(F/4)\frac{t/\Lambda^{1/\alpha}}{\alpha N^{\frac{1}{\alpha}-1}}l_{\alpha,A,1}\left(\frac{t/\Lambda^{1/\alpha}}{N^{1/\alpha}}\right)$. We take the limit $t\to\infty$, replace the summation by integration [3] and use the relation $\int_{0}^{\infty}y^{q}l_{\alpha,A,1}(y)dy=A^{q/\alpha}\Gamma(1-q/\alpha)/\Gamma(1-q)$ for $q/\alpha<1$ [66] and obtain the average displacement in $\hat{x}$

Eq. (7) shows that the average displacement is anomalous in time $\sim t^{\alpha}$. Such departure from the Einstein relation that predicts a linear dependence on time is a direct consequence of diverging mean waiting times, and for the annealed disorder was termed as Generalized Einstein relation [67]. The return probability $Q_{0}$ (Eq. (5)) depends on geometry, jump probabilities, and $F$. To finalize the calculation of $\langle x(t)\rangle$ we find the explicit form of $Q_{0}$.

The return probability $Q_{0}$ is computed in terms of $f_{N}(\mathbf{0})$, the first return probability to the origin after $N$ steps. Namely, $Q_{0}=\sum_{N=0}^{\infty}f_{N}(0)$. The probability $f_{N}(\mathbf{0})$ determines the probability $W_{N}(\mathbf{0})$ since according to the renewal equation [65], $W_{N}(\mathbf{0})=\delta_{N,0}+\sum_{i=1}^{N}f_{N}(\mathbf{0})W_{N-i}(\mathbf{0})$, which yields for the generating function of the first return probability, $\tilde{f}_{z}(\mathbf{0})=\sum_{N=0}^{\infty}f_{N}(\mathbf{0})z^{N}$, the result $\tilde{f}_{z}(\mathbf{0})=1-1/\tilde{W}_{z}(\mathbf{0})$, where $\tilde{W}_{z}(\mathbf{0})=\sum_{N=0}^{\infty}W_{N}(\mathbf{0})z^{N}$. By noting that $Q_{0}=\tilde{f}_{1}(\mathbf{0})$ the connection between $Q_{0}$ and $W_{N}(\mathbf{0})$ is finally established [65], namely $Q_{0}=1-1/\tilde{W}_{z=1}(\mathbf{0})$. Since $W_{N}(\mathbf{r})$ is a convolution of $N$ random variables, i.e., steps, its Fourier transform is the $N$th power of the single-step characteristic function $\lambda(\mathbf{k})=e^{ik_{x}}p_{\rightarrow}+e^{-ik_{x}}p\leftarrow+e^{ik_{y}}p_{\uparrow}+e^{-ik_{y}}p_{\downarrow}$ and therefore $\tilde{W}_{z=1}(\mathbf{0})=\frac{1}{4\pi^{2}}\int_{-\pi}^{\pi}\int_{-\pi}^{\pi}\frac{1}{1-\lambda(\mathbf{k})}dk_{x}dk_{y}$. In the supplemental material (SM) we show how this integral is computed and eventually obtain the explicit form of the return probability $Q_{0}$

When $w\to\infty$, the result in Eq. (8) converges to the known result [5] for an unbounded $2$-dimensional square lattice $\lim_{w\to\infty}Q_{0}=1-1\Big{/}\left[\frac{2}{\pi}\mathbf{K}\left(\frac{4}{(1+\cosh(\frac{F}{2}))^{2}}\right)\right]$, where $\bm{K}(k)=\int_{0}^{\pi/2}d\gamma\left(1-k\sin^{2}\gamma\right)^{-1/2}$ is the complete elliptic integral of the first kind [64].

Equations (7 - 8) provide $\langle x(t)\rangle$ as a function of time, arbitrarily large external force $F$ and the width of the channel $w$. For small forces Eq. (7) yields (see SM)

while $4/F\gg w/\sqrt{2}$ and $w$ is an integer $\geq 1$. Our main result is immediately apparent from Eq. (9): $\langle x(t)\rangle$ unexpectedly decays with growing channel width $w$! The dependence is $w^{\alpha-1}$, meaning the motion is faster for narrow channels than for wide channels. In Fig. 1 an excellent agreement between this analytical result and numerical simulation is displayed. In addition, we observe that the dependence on $F$ is complex and non-linear as indicated through Eq.(8) and seen in Fig. 1. For small $F$ the dependence is simplified to $\sim F^{\alpha}$. Namely, the strong disorder and it’s quenched nature that impose prolonged correlations make the usual assumption of linear response inapplicable in the quasi $2$d, as was found previously for $1$d QTM  [3, 68, 58, 59]. We note that for the case of strong annealed disorder, (CTRW), the dependence on $F$ (when $F\to 0$) is linear [59] (see the dashed line in Fig. 1).

To emphasize the mobility enhancement due to the channel width constraint, we calculate (see SM) the ratio of the $\langle x(t)\rangle$ for a given $w$, and the average displacement for unrestricted $2$d motion, $\langle x_{\infty}(t)\rangle$, i.e., $w\to\infty$. For $F\to 0$ we obtain

implying that imposing a geometrical constraint enhances the transport, and the stronger the constraint (narrower channel), the larger the enhancement! Note that the logarithmic term in $F$ enters Eq. (10) due to the critical properties of $Q_{0}$ in unrestricted $2$d (see [59] and SM for details). In Fig.2 we present the excellent agreement between the analytical results for $\langle x(t)\rangle\Big{/}\langle x_{\infty}(t)\rangle$ and numerical simulations when explored as a function of $F$, $w$ and $\alpha$. The enhancement associated with geometrical restriction repeats itself in all the presented cases and is preserved for smaller times (see SM). Fig. 2(a) shows that the effect disappears ($\langle x(t)\rangle\Big{/}\langle x_{\infty}(t)\rangle=1$) when the disorder is not strong ($\alpha>1$), or if the disorder is not quenched.

In our derivation, we assumed periodic boundary conditions in the channel. In Fig. 2(a), we show (numerically) that for reflecting boundary conditions, the obtained effects of non-linear dependence on $F$ and transport enhancement due to geometrical restriction are preserved. Additional details are provided in SM, and we intend to address this issue in future work. Equation 9 summarizes the unconventional effect of strong and quenched disorder. The regular expectation for $\langle x(t)\rangle$ is $\langle x(t)\rangle=\mu Ft$ where $\mu$ is the mobility. Strong disorder modifies the temporal dependence and the regular Einstein relation. Quenchness breaks linear response and introduces the non-linear dependence on $F$, and here we have shown that the properties of the mobility $\mu$ are counterintuitive. First of all, the mobility $\mu$ for QTM is anomalous since it can’t be defined as $\langle x(t)\rangle/tF$, but rather as $\mu=\langle x(t)\rangle/t^{\alpha}F^{\alpha}$. From Eq. (9) $\mu=w^{\alpha-1}\big{/}4^{\alpha}A\Gamma^{2}[\alpha+1]$ while $\alpha=T/Tg<1$. The mobility is enhanced as the channel width $w$ decreases. While the expectation is that additional pathways, which start to appear with relaxed geometrical constraints, will lead to a speed-up of the transport [69], we observe an opposite behavior. In the presence of strong and quenched disorder, stricter geometrical constraints improve mobility.

Mathematical reasons for such counterintuitive enhancement are rooted in the properties of the local time $S_{\alpha}$, transformation constraint $\Lambda$, and geometrical dependence of $Q_{0}$. The intuition behind the found effect is based on what is known as the “big jump principle” [70, 71, 72]. When scale-free waiting times (Eq. (1)) are in play, it is not the accumulation of many events but rather a single maximal event that governs the overall behavior. Naturally, when such a single event is excluded, for example, by replacing the site with maximal waiting time by significantly shorter $\tau$, it will lead to faster transport [73]. Our results suggest that narrowing the channel width decreases the number of possible sites the particle will sample during transport and effectively modifies this single dominant arrest time. The quenched nature of the disorder is a crucial ingredient for this to work. A further in-depth analysis and experimental research of this phenomenon is warranted. We expect that such enhancement will be useful to optimize transport in a channel media relevant for applications in nanotechnology and nanomedicine [11] and transport in porous media [73].