2.1 Porous medium, liquid and tracers

The PMMA grains composing the porous medium have been specifically chosen for their smooth surface, good transparency and monodispersity. The solid volume fraction is close to the random packing fraction of frictional spheres $\simeq 0.58$ and the porosity (void fraction) is therefore $\simeq 0.42$. The fluid employed, Triton X-100, is Newtonian and has a large viscosity, $\unicode[STIX]{x1D702}=0.270~\text{Pa}~\text{s}$. The typical average fluid velocity is $\langle U\rangle \sim 100~\unicode[STIX]{x03BC}\text{m}~\text{s}^{-1}$. This yields a Reynolds number $\unicode[STIX]{x1D70C}\langle U\rangle d/\unicode[STIX]{x1D702}\sim 10^{-3}$, ensuring that inertial effects are negligible. At room temperature, the fluid refractive index ($\simeq 1.492$) is close to that of the PMMA grains ($\simeq 1.491$), which makes the porous medium nearly transparent. To optimize the transparency, we finely match both indices of refraction down to the fourth digit (Souzy et al. Reference Souzy, Lhuissier, Villermaux and Metzger2017) by adjusting the temperature of the water bath surrounding the porous medium to $T=28.1\pm 0.5\,^{\circ }\text{C}$. The flow visualizations and the dispersion experiments have been performed by seeding the fluid with small fluorescent particles (PMMA B-particles from MF-Rhodamine) as shown in figure 1(b). The particles have a diameter $d_{t}=3.23~\unicode[STIX]{x03BC}\text{m}$, which is much smaller than the pore size ($d_{t}\simeq d/620$). This allows an unprecedented resolution for passively sampling the fluid flow – as a comparison $d_{t}/d$ was around $1/63$, $1/32$ and $1/35$ for the measurements by Moroni & Cushman (Reference Moroni and Cushman2001), Datta et al. (Reference Datta, Chiang, Ramakrishnan and Weitz2013) and Carrel et al. (Reference Carrel, Morales, Dentz, Derlon, Morgenroth and Holzner2018), respectively. The particle diffusivity is negligible (the Péclet number based on their size, $\unicode[STIX]{x1D702}\langle U\rangle dd_{t}/k_{B}T$, where $k_{B}$ is the Boltzmann constant, is larger than $10^{10}$).

2.2 Imaging

The solid structure and fluid velocity inside the porous medium are imaged in a two-dimensional (2-D) plane formed by a planar vertical laser sheet, which is displaced across the porous medium to reconstruct the 3-D fields (see figure 1a). Each section of the porous medium illuminated by the sheet is imaged with a high-resolution camera (ORCA-Flash4.0, $2048\times 1048~\text{pixels}$, 16 bit) mounted with a macroscopic lens (Sigma APO-Macro-180mm-F3.5-DG). This camera has been chosen for its linear response and its high signal-to-noise ratio. The laser sheet is formed by reflecting a laser beam (2 W, 532 nm) on a rotating mirror (${\sim}10^{4}~\text{r.p.m.}$), which yields a highly uniform sheet. It is focused down to a thickness of ${\approx}60~\unicode[STIX]{x03BC}\text{m}$ using plano-convex lenses. Finally, scattered light is discarded with the help of a high-pass filter (590 nm) located before the lens; this ensures that only the light emitted by the fluorescent tracers is collected by the camera.

3.1 Velocity measurements

Figure 1(b) shows a zoom on a typical image where the grain sections are highlighted in grey. All the grains have the same size. The apparent size differences only arise from their different distances to the laser sheet plane. The fluid velocity is measured by tracking the fluorescent tracers (bright dots in the image) using the open-source particle tracking velocimetry software TracTrac developed by Heyman (Reference Heyman2019). The volume fraction of tracers is sufficiently small (${\sim}10^{-5}$) not to affect the flow. The algorithm tracks each individual tracer and provides its displacement between successive frames. Within this stationary flow, the velocity field is obtained by accumulating data over 1000 successive images (separated by 0.1 s) such that, in each pixel within the fluid, approximately 10 particles were detected and tracked. The strength of the present TracTrac algorithm is to provide the fluid velocity within each pixel (here the pixel size is $d/228$), thus largely outperforming any particle image velocimetry techniques in terms of spatial resolution. Figure 2(a) shows a map of the raw 2-D velocity amplitude $U_{2\text{-}\text{D}}/\langle U_{2\text{-}\text{D}}\rangle$, with $U_{2\text{-}\text{D}}\equiv \sqrt{u^{2}+v^{2}}$, obtained with the algorithm. The raw data are post-processed to successively (i) mask regions inside the grains, (ii) fill the few pixels where no tracer has been detected with the surrounding pixels’ median value, and (iii) smooth the velocity field using a Gaussian filter of width 5 pixels. The post-processed field is shown in figure 2(b). As shown in figure 2(c), the post-processing has a negligible effect on the velocity distribution.

Figure 3.

(a) Experimental 3-D velocity field and microstructure. The colour scale encodes the fluid velocity magnitude, $U=\sqrt{u^{2}+v^{2}+w^{2}}$. (b,c) Distributions of the velocity components (b) $u$ along the mean flow direction ($\boldsymbol{e}_{\boldsymbol{x}}$) and (c) $v,w$ along the transverse directions ($\boldsymbol{e}_{\boldsymbol{y}},\boldsymbol{e}_{\boldsymbol{z}}$). (d) Distribution of the velocity magnitude.

The 3-D velocity field is obtained by repeating the above 2-D velocity field measurement over two perpendicular sets of 50 vertical sections spaced by $d/10$, as depicted in figure 1(a). The sampled volume of porous medium is $9d\times 5d\times 5d$. The two sets of 2-D fields are linearly interpolated to provide the three velocity components ($u,v,w$) in the three directions ($\boldsymbol{e}_{\boldsymbol{x}},\boldsymbol{e}_{\boldsymbol{y}},\boldsymbol{e}_{\boldsymbol{z}}$) over the whole sampled domain ($2048\times 1048\times 1048$  voxels). Finally, a weighted divergence correction scheme (Wang et al. Reference Wang, Gao, Wei, Li and Wang2017) is used to ensure that the interpolated field is divergence-free. Again, both the interpolation and the divergence-free algorithm have been checked to have a negligible effect on the velocity distributions.

The reconstructed 3-D velocity field is shown together with the porous medium microstructure in figure 3(a). It provides the three velocity components with a higher spatial resolution ($d/228$) than the most recent numerical simulations performed at the pore scale (Maier et al. Reference Maier, Kroll, Kutsovsky, Davis and Bernard1998; Kang et al. Reference Kang, Anna, Nunes, Bijeljic, Blunt and Juanes2014; Jin et al. Reference Jin, Langston, Pavlovskaya, Hall and Rigby2016).

3.2 Velocity distributions

The distribution of the velocity components ($u,v,w$) and magnitude $U=\sqrt{u^{2}+v^{2}+w^{2}}$ are presented in figure 3(b–d). The velocity distributions are similar to those reported from previous experiments (Datta et al. Reference Datta, Chiang, Ramakrishnan and Weitz2013; Holzner et al. Reference Holzner, Morales, Willmann and Dentz2015; Carrel et al. Reference Carrel, Morales, Dentz, Derlon, Morgenroth and Holzner2018) and numerical simulations (Dentz et al. Reference Dentz, Icardi and Hidalgo2018). The distribution of the streamwise component ($u$) is skewed towards positive values and has an exponential tail. Small negative velocities, indicating the presence of back-flow, are present, although with a very weak probability. The distributions of the transverse velocities ($v,w$) are symmetric and have exponential shapes. Last, the distribution of the velocity magnitude is also exponentially tailed but it is found to follow, at low velocity, a power law, $\text{p.d.f.}~(U)\propto U^{\unicode[STIX]{x1D6FC}}$, with a positive exponent, $\unicode[STIX]{x1D6FC}\simeq 0.89$. As discussed in the introduction, the characteristics of the velocity distribution at small velocities control the longitudinal dispersion behaviour (Berkowitz & Scher Reference Berkowitz and Scher1995; Bijeljic, Mostaghimi & Blunt Reference Bijeljic, Mostaghimi and Blunt2011). In order to clarify these crucial characteristics, we have refined the velocity measurements using a dynamic-range technique.

3.3 Refined measurements at low velocities

The tracking algorithm accurately detects the particle displacements from 0.1 to 5 pixels only. Below 0.1 pixel, the displacement signal compares with the noise level, and above 5 pixels, the algorithm is likely to swap neighbouring particles. This experimental limitation actually biases the velocity field, especially at low velocity. To circumvent this problem, a dynamic-range measurement is used. It consists in performing three measurements of the velocity field with a constant frame rate (as described in § 3.1) but with a 10-fold increase of the flow rate between each measurement. Note that this dynamic range measurement was performed only over a single vertical plane since the statistics of the fluid velocity over that plane is found to be identical and thus representative of that obtained over the whole sampled volume. For each measurement, only the velocities corresponding to particle displacements in the confidence range of 0.5–5 pixels are retained.

Figure 4.

Dynamic-range velocity measurements. (a) Map showing the regions of the three ranges of measurement: (i) yellow, fastest regions of the flow measured at low flow rate; (ii) light blue, intermediate regions; (iii) dark blue, slowest regions measured at large flow rate. (b) Comparison between the velocity distributions obtained with the dynamic-range technique and with a single measurement (same data as in figure 2c).

Figure 5.

(a) Distribution of the streamline orientation $\cos \unicode[STIX]{x1D703}$ and average velocity magnitude for a given orientation $\langle U\rangle _{\cos \unicode[STIX]{x1D703}}$. Here $A$ and $B$ are normalizing constants ($A=a/\!\sinh a$, $B=2(a+b)\sinh a\sinh b/[ab\sinh (a+b)]$). Inset: definition of the streamline orientation $\unicode[STIX]{x1D703}$. (b) Relative distribution of the velocity magnitude for a given orientation. The black line represents (3.2). Inset: joint distribution of the magnitude and orientation.

Figure 4(a) shows the regions corresponding to the three velocity ranges of measurement, with the fastest regions of the flow (measured with the slowest flow rate) in yellow, the intermediate regions in light blue and the slowest regions (largest flow rate) in dark blue. The unbiased distribution of the 2-D velocity magnitude obtained with this dynamic-range technique is found to be flat at low velocity ($\text{p.d.f.}(U_{2\text{-}\text{D}})\propto {U_{2\text{-}\text{D}}}^{0}$); see figure 4(b). This refined measurement procedure was not extended to the 3-D velocity field because of time constraints. However, since both $\text{p.d.f.}(u)$ and $\text{p.d.f.}(v)$ are flat at low velocity (see appendix A), the fact that $\text{p.d.f.}(U_{2\text{-}\text{D}}\equiv \sqrt{u^{2}+v^{2}})$ is also flat indicates that $u$ and $v$ are strongly correlated for small velocities (small values of $v$ are found where $u$ is small). Moreover, by symmetry, $w$ must be similarly correlated to $u$. Therefore, since the three velocity components are correlated and since their distribution is flat at low velocity, the 3-D velocity magnitude distribution, $\text{p.d.f.}(U_{3\text{-}D}\equiv \sqrt{u^{2}+v^{2}+w^{2}})$, must be flat as well at low velocity, in agreement with the numerical simulations performed recently by AlAdwani (Reference AlAdwani2017) and Dentz et al. (Reference Dentz, Icardi and Hidalgo2018).

Before discussing the implications of the shape of the velocity distribution on the dispersion process (see next section), we provide below an empirical correlation between the magnitude of the velocity and its orientation.

3.4 Velocity magnitude and orientation

Figure 5 presents the statistics of the streamline local orientations $\text{p.d.f.}(\cos \unicode[STIX]{x1D703})$ and the mean velocity magnitude for a given orientation $\langle U\rangle _{\cos \unicode[STIX]{x1D703}}$ obtained from the 3-D experimental field over the velocity range for which measurements are accurate ($U\geqslant 0.3\langle U\rangle$). The streamline orientation is found to be well fitted by $\text{p.d.f.}(\cos \unicode[STIX]{x1D703})\propto \text{e}^{a\cos \unicode[STIX]{x1D703}}$, with a strong dependence on $\cos \unicode[STIX]{x1D703}$ embedded in the value of $a\simeq 4$; see inset of figure 5(a). This trend simply reflects that streamlines are strongly aligned with the longitudinal direction ($\boldsymbol{e}_{\boldsymbol{x}}$). Besides, we find that the magnitude of the velocity is also strongly correlated with its orientation. The average velocity magnitude for a given orientation $\langle U\rangle _{\cos \unicode[STIX]{x1D703}}$ increases when streamlines are more aligned with the flow, following the exponential trend $\langle U\rangle _{\cos \unicode[STIX]{x1D703}}\propto \text{e}^{b\cos \unicode[STIX]{x1D703}}$, with $b\simeq 1.5$.

These two experimental observations suggest a compact, 3-D characterization of the one-point statistics of the velocity field. Assuming for simplicity that the relative distribution of the velocity for each orientation is independent of the orientation (which is presumably a reasonable ansatz given the isotropy of the medium structure) imposes that the joint distribution of $U$ and $\unicode[STIX]{x1D703}$ follows

(3.1)$$\begin{eqnarray}\text{p.d.f.}(U,\cos \unicode[STIX]{x1D703})=\frac{1}{\langle U\rangle _{\cos \unicode[STIX]{x1D703}}}\,\text{p.d.f.}(\cos \unicode[STIX]{x1D703})\,\text{p.d.f.}(U/\langle U\rangle _{\cos \unicode[STIX]{x1D703}}),\end{eqnarray}$$

where, as shown in figure 5(b), the relative velocity distribution is found experimentally to be close to

(3.2)$$\begin{eqnarray}\text{p.d.f.}(U/\langle U\rangle _{\cos \unicode[STIX]{x1D703}})=\text{e}^{-U/\langle U\rangle _{\cos \unicode[STIX]{x1D703}}}.\end{eqnarray}$$

Equation (3.1) summarizes all the one-point information about the 3-D velocity field. It captures the main dependences observed experimentally, including an accurate description of the statistics at small velocity, since it embeds $\text{p.d.f.}(U)\propto U^{0}$ for $U\ll \langle U\rangle _{\cos \unicode[STIX]{x1D703}}$.

4.1 Dispersion measurements

Figure 6.

Dispersion measurements: successive experimental images of a blob of small non-diffusive tracers flowing through the index-matched 3-D bead-pack (see also supplementary movie 1, available at https://doi.org/10.1017/jfm.2020.113). Only those particles contained in the $z=0$ plane are visible (the field of view is $20d$ long in the streamwise direction and the tracers have a diameter $d_{t}=d/620$).

Figure 7.

Dispersion measurements: trajectories of non-diffusive point particles advected numerically in the 3-D experimental velocity field (see also supplementary movies 2 and 3 for an immersive exploration of the 3-D velocity field (Pore Aventura 2019)). The six thick lines highlight six specific trajectories.

The dispersion process induced by this flow field is characterized in two ways: (i) directly, by monitoring the evolution of a blob of small non-diffusive particles (see figure 6), and (ii) indirectly, by numerically advecting non-diffusive point particles in the experimental 3-D velocity field (see figure 7).

Figure 8.

Dispersion results. (a) Overlay of 13 blobs of non-diffusive tracers injected at 13 different locations in the porous medium at $t/\unicode[STIX]{x1D70F}=0$ (red), 3 (green) and 10.5 (blue), with $\unicode[STIX]{x1D70F}=d/\langle U\rangle$. (b) Trajectories of 1250 non-diffusive point particles advected numerically in the experimental velocity field (grey) and their positions at $t/\unicode[STIX]{x1D70F}=0$ (red), 43 (green) and 136 (blue). (c) Variances of the position ${\unicode[STIX]{x1D70E}_{i}}^{2}\equiv \langle i^{2}\rangle -\langle i\rangle ^{2}$, with $i=x$, $y$ or $z$, obtained with both methods. The error bars represent the typical standard deviation between experimental runs. (d) Same data over a longer time.

The first method consists in injecting with a thin needle a small blob of fluid, typically of the size of the grains, seeded with non-diffusive fluorescent tracers (see § 2) in the bulk of the porous medium. Once the injection is completed, the flow is turned on and the evolution of the blob is captured within one plane illuminated by the laser sheet (see figure 6). The dispersion process is quantified through image processing by measuring the spatial variance $\unicode[STIX]{x1D70E}^{2}$ of the tracer locations in both the longitudinal and transverse directions. The measurement can be achieved until $t/\unicode[STIX]{x1D70F}=10$, since all tracers have to remain in the field of view. The experiment is repeated 13 times by injecting the blob at different locations within the porous medium.

The second method consists in using the 3-D experimental fluid velocity field to advect numerically non-diffusive point particles. Some examples of trajectories are provided in figure 7. A total of 1250 tracers are positioned randomly within the fluid phase and advected using an explicit Runge–Kutta scheme. In order to allow a long-time advection ($t/\unicode[STIX]{x1D70F}>100$), each point particle that reaches a boundary of the fluid domain is re-injected upwards at a location where the fluid velocity is the closest (in both magnitude and orientation) to its current velocity. This re-injection procedure was found to have no effect on the dispersion statistics. As shown in figure 3(d), the velocity distribution of these Lagrangian point particles is identical to that of the Eulerian velocity field, as expected in the present case of a volume-conserving flow (Dentz et al. Reference Dentz, Kang, Comolli, Le Borgne and Lester2016).

The dispersion results are gathered in figure 8. Figure 8(a) shows overlay images of the 13 blob experiments at $t/\unicode[STIX]{x1D70F}=0$, 3 and 10.5, where $\unicode[STIX]{x1D70F}=d/\langle U\rangle$. The typical dispersion envelope explored by these non-diffusive particles is indicated with a dotted line. Figure 8(b) shows the trajectories of the 1250 numerical point particles (in grey) as well as their positions at $t/\unicode[STIX]{x1D70F}=0$, 43 and 136 (in colour). The dispersion process is quantified by measuring the variance of the tracers’ position, in both the longitudinal and transverse directions, i.e. ${\unicode[STIX]{x1D70E}_{i}}^{2}\equiv \langle i^{2}\rangle -\langle i\rangle ^{2}$, where $i$ stands for $x$, $y$ or $z$. The variances obtained with the two methods are compared in figure 8(c) until $t/\unicode[STIX]{x1D70F}=10$. After a short ballistic regime, where $\unicode[STIX]{x1D70E}^{2}\propto t^{2}$ for $t/\unicode[STIX]{x1D70F}\lesssim 1$, the variance becomes seemingly anomalous for $1\lesssim t/\unicode[STIX]{x1D70F}\leqslant 10$. It is super-diffusive in the streamwise direction ($\unicode[STIX]{x1D70E}_{x}^{2}\propto t^{1.4}$) and sub-diffusive transversely (${\unicode[STIX]{x1D70E}_{y}}^{2}\propto t^{0.7}$). Note that similar exponents have also been reported for the same range of advection times by Kang et al. (Reference Kang, Anna, Nunes, Bijeljic, Blunt and Juanes2014), who simulated the fluid flow in the 3-D pore space of a real rock. Moreover, although the 3-D velocity field is slightly biased at small velocities (as discussed in § 3.3), the streamwise and transverse dispersions obtained from the numerical advection are very similar to those measured directly in the set-up, which validates the numerical advection procedure.

To characterize the dispersion over longer times, we rely on numerical advection measurements. These measurements, presented in figure 8(d), show that the anomalous trends are only transient. For $t/\unicode[STIX]{x1D70F}\gtrsim 10$, both the longitudinal and transverse dispersions eventually relax to a Fickian regime with variances following

(4.1a,b)$$\begin{eqnarray}\unicode[STIX]{x1D70E}_{x}^{2}=2D_{\Vert }t\quad \text{and}\quad {\unicode[STIX]{x1D70E}_{y}}^{2}={\unicode[STIX]{x1D70E}_{z}}^{2}=2D_{\bot }t,\end{eqnarray}$$

respectively, with dispersion coefficients equal to

(4.2a,b)$$\begin{eqnarray}D_{\Vert }\simeq (0.45\pm 0.05)d^{2}/\unicode[STIX]{x1D70F}\quad \text{and}\quad D_{\bot }\simeq (0.027\pm 0.002)d^{2}/\unicode[STIX]{x1D70F}.\end{eqnarray}$$

4.2 Link between dispersion and velocity distributions

The above measurements require some interpretation. The ballistic dispersion regime lasts typically until $t\sim \unicode[STIX]{x1D70F}$, which indicates that the velocity of the tracers remains correlated over a displacement of order $d$. Beyond that time, the dispersion becomes anomalous.

In the longitudinal direction, velocities are essentially positive and the (super-diffusive) anomaly is controlled by the occurrence of the smallest velocities, which significantly delay some of the tracers. Following Berkowitz et al. (Reference Berkowitz, Cortis, Dentz and Scher2006), we model the displacement of each tracer as a succession of random longitudinal steps with a fixed length $d$ and a random velocity $u$ that is conserved over the step but uncorrelated with the previous velocities. This CTRW is fully constrained by the velocity distribution. It is computed numerically for $10^{4}$ uncorrelated tracers (the first step length is chosen randomly between 0 and $d$, to account for the random initial position of the tracers, and the tracers’ displacement is sampled with a substep scale at desired times). For an exponential velocity distribution, which has the same flat shape ($p(u)\propto u^{0}$) at low velocity and the same mean value $\langle u\rangle /\langle U\rangle =0.79$ as the experimental distribution, the CTRW yields $\unicode[STIX]{x1D70E}_{x}^{2}/d^{2}\sim \log (t/\unicode[STIX]{x1D70F})t/\unicode[STIX]{x1D70F}$. As shown in figure 9, this trend is in good agreement with the experimental data over the range $1\lesssim t/\unicode[STIX]{x1D70F}\lesssim 10$ and agrees with the prediction that the dispersion is super-diffusive when the low-velocity exponent of the distribution ($p(u)\propto u^{\unicode[STIX]{x1D6FC}}$) satisfies $\unicode[STIX]{x1D6FC}\leqslant 0$ (Berkowitz et al. Reference Berkowitz, Cortis, Dentz and Scher2006).

To understand the final relaxation towards a Fickian regime, the finite minimal velocity of the tracers has to be considered. Indeed, in a real flow, advected tracers always have a finite minimal velocity. This finite minimal velocity can arise from the finite size $d_{t}$ of the tracers, which prevents them from sampling the smallest fluid velocities close to the walls. In the present case, this minimal velocity can be estimated as $U_{min}/\langle U\rangle \approx 10d_{t}/d\approx 0.02$, since the typical pore size is equal to $d/10$; see figure 3(a). The finite minimal velocity can also result from the finite thermal diffusivity of the tracers. Such tracers can escape from a region of low velocity more quickly through the effect of molecular diffusion than through the effect of advection alone (Saffman Reference Saffman1959). In such a case, the minimal velocity is $U_{min}\sim (10\langle U\rangle /d)\times \unicode[STIX]{x1D6FF}$, where $\unicode[STIX]{x1D6FF}\sim \sqrt{Dd/U_{min}}$ is the transverse diffusive displacement and $D$ is the tracers’ (thermal) diffusion coefficient. This yields $U_{min}/\langle U\rangle \approx (10^{2}D/\langle U\rangle d)^{1/3}$. This minimal velocity has a crucial impact on the longitudinal dispersion, precisely because the dispersion anomaly is controlled by the low velocities. To illustrate this point, we simulate a CTRW based on the same exponential distribution but discarding all velocities below $U_{min}$. As shown in figure 9, when $U_{min}/\langle U\rangle$ is set to 0.01, 0.03, 0.1 or 0.2, the variance relaxes to a Fickian regime beyond $t_{F}/\unicode[STIX]{x1D70F}\simeq 100$, 33, 10 and 5, respectively, instead of following $\log (t/\unicode[STIX]{x1D70F})t/\unicode[STIX]{x1D70F}$ as previously. The Fickian regime is reached once the velocities of all the tracers become decorrelated. The CTRW dispersion matching the experiments is obtained with $U_{min}/\langle U\rangle =0.1$. This is consistent with the fact that the 3-D velocity distribution, as a result of the limited experimental resolution, departs from a flat distribution below $U/\langle U\rangle \simeq 0.1$ (see figure 3d).

Figure 9.

Streamwise and transverse variances, $\unicode[STIX]{x1D70E}_{x}^{2}/d^{2}$ and ${\unicode[STIX]{x1D70E}_{z}}^{2}/d^{2}$, obtained from the numerical advection in the 3-D experimental velocity field compared to the CTRW simulations. CTRW based on a flat velocity distribution ($\text{p.d.f.}(u)\propto u^{0}$ at small $u$) yields a persistent anomalous dispersion $\unicode[STIX]{x1D70E}_{x}^{2}/d^{2}\sim \log (t/\unicode[STIX]{x1D70F})t/\unicode[STIX]{x1D70F}$ (green line), which matches the experiments until $t/\unicode[STIX]{x1D70F}\sim 10$. For a flat distribution truncated below $U_{min}/\langle U\rangle =0.01$, 0.03, 0.1 and 0.2, the dispersion eventually becomes Fickian ($\unicode[STIX]{x1D70E}_{x}^{2}\propto t$) after $t_{F}/\unicode[STIX]{x1D70F}\approx 100$, 33, 10 and 5, respectively (grey lines). The experimental transverse dispersion is well captured by a CTRW with anticorrelated velocities (see text). The black dashed lines indicate the slope 1.

By contrast, the sub-diffusive anomaly in the transverse direction is not controlled by the smallest velocities but is a consequence of the meandering trajectories inside the porous medium, which makes the transverse velocities alternately positive and negative with similar magnitudes. This can be shown by simulating a unidimensional CTRW with a large probability $p$ that the tracers change direction between steps. Using an exponential velocity distribution with a mean value $\langle v\rangle /\langle U\rangle =0.34$ identical to the experimental distribution (figure 3c) and setting $p=0.86$ results in a sub-diffusive dispersion in good agreement with the experiments. As shown in figure 9, the transient sub-diffusive trend is recovered over the typical duration $\unicode[STIX]{x1D70F}/(1-p)$, over which velocities remain anticorrelated.

With the above information at hand, we can explain the values of the dispersion coefficients ($D_{\Vert }$ and $D_{\bot }$) obtained in the (long-time) Fickian regime. For tracers with a constant velocity over the porous scale $d$, the duration of each step is $\unicode[STIX]{x1D70F}^{\prime }=d/U$. The probability of $\unicode[STIX]{x1D70F}^{\prime }$ is fully constrained by the tracer velocity distribution $p_{U}(U)$ according to $p_{\unicode[STIX]{x1D70F}^{\prime }}(\unicode[STIX]{x1D70F}^{\prime })=(d^{2}/\langle U\rangle \unicode[STIX]{x1D70F}^{\prime 3})p_{U}(d/\unicode[STIX]{x1D70F}^{\prime })$, which is defined over $0\leqslant \unicode[STIX]{x1D70F}^{\prime }\leqslant d/U_{min}$ (see details in appendix B). Now, the probability that a tracer has moved over a large curvilinear distance $s$, i.e. that $s/d$ steps have been covered at a given time $t$, is given by the $s/d$-fold convolute of $p_{\unicode[STIX]{x1D70F}^{\prime }}$, i.e. ${p_{\unicode[STIX]{x1D70F}^{\prime }}}^{\ast \,s/d}(t)$. Using the central limit theorem and $\langle x\rangle /s=\langle u\rangle /\langle U\rangle$, the streamwise spatial variance $\unicode[STIX]{x1D70E}_{x}^{2}$ can be expressed in terms of the variance of the step durations, ${\unicode[STIX]{x1D70E}_{\unicode[STIX]{x1D70F}^{\prime }}}^{2}\equiv \langle \unicode[STIX]{x1D70F}^{\prime 2}\rangle -\langle \unicode[STIX]{x1D70F}^{\prime }\rangle ^{2}$, as $\unicode[STIX]{x1D70E}_{x}^{2}/d^{2}\approx \langle u\rangle ^{2}/\langle U\rangle ^{2}{\unicode[STIX]{x1D70E}_{\unicode[STIX]{x1D70F}^{\prime }}}^{2}t/\unicode[STIX]{x1D70F}^{3}$. For a flat velocity distribution ($p_{U}(U)\propto U^{0}$), the step duration variance is prescribed by the minimal velocity as ${\unicode[STIX]{x1D70E}_{\unicode[STIX]{x1D70F}^{\prime }}}^{2}/\unicode[STIX]{x1D70F}^{2}\approx \ln (\langle U\rangle /U_{min})$. Altogether, one obtains

(4.3)$$\begin{eqnarray}D_{\Vert }=\frac{\unicode[STIX]{x1D70E}_{x}^{2}}{2t}\approx \frac{1}{2}\ln \left(\frac{\langle U\rangle }{U_{min}}\right)\frac{\langle u\rangle ^{2}}{\langle U\rangle ^{2}}\frac{d^{2}}{\unicode[STIX]{x1D70F}},\end{eqnarray}$$

which makes precise how the dispersion coefficient depends on the minimal velocity $U_{min}$. With the experimental values $\langle u\rangle \simeq 0.79\langle U\rangle$ and $U_{min}\simeq 0.1\langle U\rangle$, equation (4.3) yields $D_{\Vert }\simeq 0.7d^{2}/\unicode[STIX]{x1D70F}$, which is within the same order of magnitude as the measured value $0.45d^{2}/\unicode[STIX]{x1D70F}$ reported in (4.2).

In the transverse direction, the same arguments together with the probability, $1-p$, that a step is not anticorrelated with the previous ones, give similarly

(4.4)$$\begin{eqnarray}D_{\bot }=\frac{{\unicode[STIX]{x1D70E}_{y}}^{2}}{2t}=\frac{{\unicode[STIX]{x1D70E}_{z}}^{2}}{2t}\approx \frac{1}{2}\ln \left(\frac{\langle U\rangle }{U_{min}}\right)(1-p)\frac{\langle v\rangle ^{2}}{\langle U\rangle ^{2}}\frac{d^{2}}{\unicode[STIX]{x1D70F}}.\end{eqnarray}$$

For the experimental values $\langle v\rangle =\langle w\rangle \simeq 0.34\langle U\rangle$ and $U_{min}\simeq 0.1\langle U\rangle$, this yields $D_{\bot }\simeq 0.02d^{2}/\unicode[STIX]{x1D70F}$, which is also found to be close to the measurement $0.027d^{2}/\unicode[STIX]{x1D70F}$ reported in (4.2).