2.1 Irreversibility

Following on from the work of Bragg et al. (Reference Bragg, Ireland and Collins2016), we shall consider statistically stationary, homogeneous, isotropic turbulence, in which the inertial particle statistics are also stationary. From a practical perspective, this is perhaps not the most relevant system in which to study the irreversibility of particle-pair dispersion. For example, in many practical dispersion problems, particles will often be injected into the turbulent flow with velocities such that the as the particles begin to disperse, they will undergo an initial transitory phase before their single-time statistics become stationary. However, from a fundamental perspective, a statistically stationary system is the most natural starting point for examining irreversibility, since in such a system, trivial causes of irreversibility (e.g. global energy decay/growth) are removed, and the mechanisms generating irreversibility must arise solely from the intrinsic dynamics of the system. In future work we will consider non-stationary effects, in order to assess more carefully the implications of irreversible dispersion for practical problems.

For the system of interest, the statistics of the pair separation only depend upon the evolution of their separation vector $\boldsymbol{r}^{p}$ and not their centre of mass. The FIT and BIT PDFs describing the particle-pair dispersion process may then be defined as

(2.1) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}^{{\mathcal{F}}}(\boldsymbol{r},t|\unicode[STIX]{x1D743},t^{\prime })\equiv \langle \unicode[STIX]{x1D6FF}(\boldsymbol{r}^{p}(t)-\boldsymbol{r})\rangle _{\boldsymbol{ r}^{p}(t^{\prime })=\unicode[STIX]{x1D743}},\quad t^{\prime }\in [0,t], & \displaystyle\end{eqnarray}$$

(2.2) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}^{{\mathcal{B}}}(\boldsymbol{r},t^{\prime }|\unicode[STIX]{x1D743},t)\equiv \langle \unicode[STIX]{x1D6FF}(\boldsymbol{r}^{p}(t^{\prime })-\boldsymbol{r})\rangle _{\boldsymbol{ r}^{p}(t)=\unicode[STIX]{x1D743}},\quad t^{\prime }\in [0,t], & \displaystyle\end{eqnarray}$$

where $\boldsymbol{r}$ denotes the time-independent configuration-space coordinate, and $\langle \cdot \rangle _{\boldsymbol{r}^{p}(t^{\prime })=\unicode[STIX]{x1D743}}$ denotes an ensemble average conditioned upon $\boldsymbol{r}^{p}(t^{\prime })=\unicode[STIX]{x1D743}$ , and similarly for the BIT case. Due to isotropy, we need only consider the PDFs for $\Vert \boldsymbol{r}^{p}\Vert$ , that is

(2.3) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}^{{\mathcal{F}}}(r,t|\unicode[STIX]{x1D709},t^{\prime })\equiv \langle \unicode[STIX]{x1D6FF}(\Vert \boldsymbol{r}^{p}(t)\Vert -r)\rangle _{\Vert \boldsymbol{r}^{p}(t^{\prime })\Vert =\unicode[STIX]{x1D709}},\quad t^{\prime }\in [0,t], & \displaystyle\end{eqnarray}$$

(2.4) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,t^{\prime }|\unicode[STIX]{x1D709},t)\equiv \langle \unicode[STIX]{x1D6FF}(\Vert \boldsymbol{r}^{p}(t^{\prime })\Vert -r)\rangle _{\Vert \boldsymbol{r}^{p}(t)\Vert =\unicode[STIX]{x1D709}},\quad t^{\prime }\in [0,t], & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D709}\equiv \Vert \unicode[STIX]{x1D743}\Vert$ . One way to analyse the behaviour of these PDFs is to consider their evolution equations $\unicode[STIX]{x2202}_{t}{\mathcal{P}}^{{\mathcal{F}}}$ and $\unicode[STIX]{x2202}_{t^{\prime }}{\mathcal{P}}^{{\mathcal{B}}}$ . We find it more insightful, however, to consider integral formulations for the PDFs that may be constructed by substituting into (2.3) and (2.4) the integral representations of the solutions for $\boldsymbol{r}^{p}(t)$ and $\boldsymbol{r}^{p}(t^{\prime })$ , respectively.

From the kinematic equation $\dot{\boldsymbol{r}}^{p}(t)\equiv \boldsymbol{w}^{p}(t)$ we obtain

(2.5) $$\begin{eqnarray}\displaystyle d_{t}\Vert \boldsymbol{r}^{p}(t)\Vert ^{2}\equiv 2\Vert \boldsymbol{r}^{p}(t)\Vert w_{\Vert }^{p}(t), & & \displaystyle\end{eqnarray}$$

where $w_{\Vert }^{p}(t)\equiv \Vert \boldsymbol{r}^{p}(t)\Vert ^{-1}\boldsymbol{r}^{p}(t)\boldsymbol{\cdot }\boldsymbol{w}^{p}(t)$ , leading to the result

(2.6) $$\begin{eqnarray}\displaystyle \Vert \boldsymbol{r}^{p}(t)\Vert \equiv \left[\Vert \boldsymbol{r}^{p}(t^{\prime })\Vert ^{2}+2\int _{t^{\prime }}^{t}\Vert \boldsymbol{r}^{p}(t^{\prime \prime })\Vert w_{\Vert }^{p}(t^{\prime \prime })\,\text{d}t^{\prime \prime }\right]^{1/2}. & & \displaystyle\end{eqnarray}$$

Using this result in (2.3) and (2.4) we obtain

(2.7) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}^{{\mathcal{F}}}(r,t|\unicode[STIX]{x1D709},t^{\prime })=\left\langle \unicode[STIX]{x1D6FF}\left(\left[\unicode[STIX]{x1D709}^{2}+2\int _{t^{\prime }}^{t}\Vert \boldsymbol{r}^{p}(t^{\prime \prime })\Vert w_{\Vert }^{p}(t^{\prime \prime })\,\text{d}t^{\prime \prime }\right]^{1/2}-r\right)\right\rangle _{\Vert \boldsymbol{r}^{p}(t^{\prime })\Vert =\unicode[STIX]{x1D709}}, & \displaystyle\end{eqnarray}$$

(2.8) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,t^{\prime }|\unicode[STIX]{x1D709},t)=\left\langle \unicode[STIX]{x1D6FF}\left(\left[\unicode[STIX]{x1D709}^{2}-2\int _{t^{\prime }}^{t}\Vert \boldsymbol{r}^{p}(t^{\prime \prime })\Vert w_{\Vert }^{p}(t^{\prime \prime })\,\text{d}t^{\prime \prime }\right]^{1/2}-r\right)\right\rangle _{\Vert \boldsymbol{r}^{p}(t)\Vert =\unicode[STIX]{x1D709}}. & \displaystyle\end{eqnarray}$$

Since we are interested in the statistically stationary state of the system, for which the dispersion only depends upon the time difference $s\equiv t-t^{\prime }$ , we may set the conditioning time to zero, and note the time ordering $t^{\prime }\leqslant t$ , to obtain

(2.9) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)=\left\langle \unicode[STIX]{x1D6FF}\left(\left[\unicode[STIX]{x1D709}^{2}+2\int _{0}^{s}\Vert \boldsymbol{r}^{p}(s^{\prime })\Vert w_{\Vert }^{p}(s^{\prime })\,\text{d}s^{\prime }\right]^{1/2}-r\right)\right\rangle _{\unicode[STIX]{x1D709}}, & \displaystyle\end{eqnarray}$$

(2.10) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)=\left\langle \unicode[STIX]{x1D6FF}\left(\left[\unicode[STIX]{x1D709}^{2}-2\int _{-s}^{0}\Vert \boldsymbol{r}^{p}(s^{\prime })\Vert w_{\Vert }^{p}(s^{\prime })\,\text{d}s^{\prime }\right]^{1/2}-r\right)\right\rangle _{\unicode[STIX]{x1D709}}, & \displaystyle\end{eqnarray}$$

where $s\geqslant 0$ , and $\langle \cdot \rangle _{\unicode[STIX]{x1D709}}\equiv \langle \cdot \rangle _{\Vert \boldsymbol{r}^{p}(0)\Vert =\unicode[STIX]{x1D709}}$ . Note that, to be precise, whereas in the FIT case, $s\in [0,\infty ]$ , in the BIT case $s\in [0,T_{0}]$ , where $-s=-T_{0}$ corresponds to the initial time of the dynamical system, since in the original variables $t^{\prime }\in [0,t]$ not $t^{\prime }\in [-\infty ,t]$ . However, since we are considering the stationary state, we could effectively set $T_{0}\rightarrow \infty$ .

In Bragg et al. (Reference Bragg, Ireland and Collins2016) we presented arguments for the physical mechanisms that generate irreversible pair dispersion in turbulence. Here we will present the arguments in alternative form, especially to bring out in more detail some aspects of the problem that were not considered thoroughly in Bragg et al. (Reference Bragg, Ireland and Collins2016).

There are essentially two ways that (2.9) and (2.10) could differ. First, they will differ trivially if the statistics of $w_{\Vert }^{p}$ (conditioned on $\Vert \boldsymbol{r}^{p}(0)\Vert =\unicode[STIX]{x1D709}$ ) differ for $s$ and $-s$ . In one sense this effect would appear to be absent for the system we are considering because of statistical stationarity. However, the conditioning $\langle \cdot \rangle _{\unicode[STIX]{x1D709}}$ actually means that (2.9) and (2.10) depend upon multi-time statistics of $w_{\Vert }^{p}$ , and as such they could differ for $s$ and $-s$ . The second effect is non-trivial: equation (2.9) reveals that only particle pairs whose motion is dominated by $w_{\Vert }^{p}>0$ (i.e. such that $\int _{0}^{s}\Vert \boldsymbol{r}^{p}(s^{\prime })\Vert w_{\Vert }^{p}(s^{\prime })\,\text{d}s^{\prime }>0$ ) will go to larger separations $r>\unicode[STIX]{x1D709}$ , whereas only particle pairs whose motion is dominated by $w_{\Vert }^{p}<0$ (i.e. such that $\int _{0}^{s}\Vert \boldsymbol{r}^{p}(s^{\prime })\Vert w_{\Vert }^{p}(s^{\prime })\,\text{d}s^{\prime }<0$ ) will go to smaller separations $r<\unicode[STIX]{x1D709}$ . Conversely, the BIT result in (2.10) reveals that only particle pairs whose motion was dominated by $w_{\Vert }^{p}<0$ were at larger separations $r>\unicode[STIX]{x1D709}$ in the past, whereas only particle-pairs whose motion was dominated by $w_{\Vert }^{p}>0$ were at smaller separations $r<\unicode[STIX]{x1D709}$ in the past. It follows then that if the PDFs of $w_{\Vert }^{p}$ are themselves asymmetric, the separation PDF will be irreversible, (this is actually a necessary but not sufficient criteria for dispersion irreversibility. The other condition required is that the correlation time scales of $w_{\Vert }^{p}$ must be finite, but this condition is always satisfied in real turbulent flows. See Bragg et al. (Reference Bragg, Ireland and Collins2016) for further details) i.e. ${\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)\neq {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)$ .

These observations are purely kinematic; the dynamical origin of any irreversibility depends upon the dynamics governing $\boldsymbol{w}^{p}$ . Following on from the work in Bragg et al. (Reference Bragg, Ireland and Collins2016), we shall consider heavy, point particles whose motion is governed by a Stokes drag force. In this case

(2.11) $$\begin{eqnarray}\displaystyle \dot{\boldsymbol{w}}^{p}(t)=\frac{1}{\unicode[STIX]{x1D70F}_{p}}(\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{x}^{p}(t),\boldsymbol{r}^{p}(t),t)-\boldsymbol{w}^{p}(t)), & & \displaystyle\end{eqnarray}$$

where $\boldsymbol{x}^{p}(t)$ and $\boldsymbol{x}^{p}(t)+\boldsymbol{r}^{p}(t)$ are the locations of the two particles, $\boldsymbol{w}^{p}(t)$ is their relative velocity, and $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{x}^{p}(t),\boldsymbol{r}^{p}(t),t)$ is the difference in the fluid velocity evaluated at the two particle positions. To avoid any confusion, we note that in (2.11), the time label ‘ $t$ ’ is generic and does not necessarily coincide with the ‘ $t$ ’ appearing in the expressions (2.3) and (2.4). Furthermore, since we are considering homogeneous turbulence, we will usually neglect the $\boldsymbol{x}^{p}(t)$ argument when discussing statistical properties of the dispersion.

In the absence of turbulence (i.e. setting $\unicode[STIX]{x0394}\boldsymbol{u}=\mathbf{0}$ ), (2.11) would reduce to $\dot{\boldsymbol{w}}^{p}(t)=-\boldsymbol{w}^{p}(t)/\unicode[STIX]{x1D70F}_{p}$ , and the dispersion would be trivially irreversible simply because of dissipation, giving $d_{t}\Vert \boldsymbol{w}^{p}(t)\Vert ^{2}<0$ (where $d_{t}\equiv d/dt$ ). However, in the presence of turbulence, $d_{t}\Vert \boldsymbol{w}^{p}(t)\Vert ^{2}$ is not strictly negative since the fluid does work on the particles. In this case, understanding the irreversibility is much more complex since the behaviour of $\boldsymbol{w}^{p}$ then depends upon how the inertial particle pairs interact with the field $\unicode[STIX]{x0394}\boldsymbol{u}$ , and the properties of $\unicode[STIX]{x0394}\boldsymbol{u}$ itself. More specifically, since we are interested in statistical irreversibility, we must understand how $w_{\Vert }^{p}$ in (2.9) and (2.10) depends upon the statistical properties of $\unicode[STIX]{x0394}\boldsymbol{u}$ . A crucial point in this regard is that in (2.9) and (2.10), the behaviour of $w_{\Vert }^{p}(s^{\prime })$ is not only controlled by its dynamical equation, but also by the fact that the trajectories (along which $w_{\Vert }^{p}(s^{\prime })$ is measured) that contribute to the ensembles in (2.9) and (2.10) are conditioned, i.e. the contributing trajectories have to satisfy $\Vert \boldsymbol{r}^{p}(s^{\prime }=0)\Vert =\unicode[STIX]{x1D709}$ . This trajectory conditioning is the reason why FIT and BIT separating pairs experience values of $w_{\Vert }^{p}(s^{\prime })$ with different signs, as discussed earlier.

The forward solution to (2.11), satisfying the condition $\Vert \boldsymbol{r}^{p}(s^{\prime }=0)\Vert =\unicode[STIX]{x1D709}$ , may be written as

(2.12) $$\begin{eqnarray}\displaystyle \boldsymbol{w}^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)=\text{e}^{-s^{\prime }/\unicode[STIX]{x1D70F}_{p}}\boldsymbol{w}^{p}(0|\unicode[STIX]{x1D709},0)+\frac{1}{\unicode[STIX]{x1D70F}_{p}}\int _{0}^{s^{\prime }}\text{e}^{-(s^{\prime }-s^{\prime \prime })/\unicode[STIX]{x1D70F}_{p}}\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(s^{\prime \prime }|\unicode[STIX]{x1D709},0),s^{\prime \prime })\,\text{d}s^{\prime \prime }. & & \displaystyle\end{eqnarray}$$

The backward solution for $\boldsymbol{w}^{p}$ can be represented by (2.12) with $s^{\prime }\rightarrow -s^{\prime }$ and $s^{\prime }\in [0,T_{0}]$ , in which case the backward trajectory is represented as a terminal value problem. In this time-reversed solution, the term $\text{e}^{s^{\prime }/\unicode[STIX]{x1D70F}_{p}}\boldsymbol{w}^{p}(0|\unicode[STIX]{x1D709},0)$ does not cause the solution $\boldsymbol{w}^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)$ to blow up. This is because $\boldsymbol{w}^{p}(0|\unicode[STIX]{x1D709},0)$ is not an arbitrary end condition for the solution (unlike the forward solution where $\boldsymbol{w}^{p}(0|\unicode[STIX]{x1D709},0)$ can be chosen arbitrarily), but implicitly contains the information on the trajectory history of the particle, and this regularizes the solution to the time-reversed equation. This can be seen by constructing $\boldsymbol{w}^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)$ as an initial value solution, with initial time $-s^{\prime }=-T_{0}$

(2.13) $$\begin{eqnarray}\boldsymbol{w}^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)=\text{e}^{(s^{\prime }-T_{0})/\unicode[STIX]{x1D70F}_{p}}\boldsymbol{w}^{p}(-T_{0}|\unicode[STIX]{x1D709},0)+\frac{1}{\unicode[STIX]{x1D70F}_{p}}\int _{-T_{0}}^{-s^{\prime }}\text{e}^{(s^{\prime }+s^{\prime \prime })/\unicode[STIX]{x1D70F}_{p}}\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(s^{\prime \prime }|\unicode[STIX]{x1D709},0),s^{\prime \prime })\,\text{d}s^{\prime \prime }.\end{eqnarray}$$

Setting $s^{\prime }=0$ in this solution shows how $\boldsymbol{w}^{p}(0|\unicode[STIX]{x1D709},0)$ depends upon the trajectory history of the particle, and by inserting this expression for $\boldsymbol{w}^{p}(0|\unicode[STIX]{x1D709},0)$ into the time-reversed form of (2.12), one indeed finds that the resulting solution does not blow up. We will find both the ‘terminal value’ and ‘initial value’ representations of the backward solution $\boldsymbol{w}^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)$ useful in our discussion.

If $\boldsymbol{w}^{p}(0|\unicode[STIX]{x1D709},0)\neq \mathbf{0}$ , then in the regime $s^{\prime }/\unicode[STIX]{x1D70F}_{p}\ll 1$ we have

(2.14) $$\begin{eqnarray}\displaystyle \Vert \text{e}^{-s^{\prime }/\unicode[STIX]{x1D70F}_{p}}\boldsymbol{w}^{p}(0|\unicode[STIX]{x1D709},0)\Vert \gg \left\Vert \frac{1}{\unicode[STIX]{x1D70F}_{p}}\int _{0}^{s^{\prime }}\text{e}^{-(s^{\prime }-s^{\prime \prime })/\unicode[STIX]{x1D70F}_{p}}\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(s^{\prime \prime }|\unicode[STIX]{x1D709},0),s^{\prime \prime })\,\text{d}s^{\prime \prime }\right\Vert , & & \displaystyle\end{eqnarray}$$

which shows that in this regime, the evolution of $\boldsymbol{w}^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)$ is determined by dissipation. Inserting the leading-order behaviour $\boldsymbol{w}^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)\approx \text{e}^{-s^{\prime }/\unicode[STIX]{x1D70F}_{p}}\boldsymbol{w}^{p}(0|\unicode[STIX]{x1D709},0)$ into (2.9) and (2.10) leads to ${\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)\neq {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)$ . In this regime the dispersion is irreversible simply because of dissipation, and particles will separate faster BIT than FIT simply because the particle kinetic energy is decaying in time. Formally, this dissipation effect will dominate the irreversibility in the short-time regime $s/\unicode[STIX]{x1D70F}_{p}\ll 1$ . However, its effect can persist up to $s/\unicode[STIX]{x1D70F}_{p}\lesssim O(1)$ in certain regimes of the flow, a point we shall discuss further in § 2.2. In Bragg et al. (Reference Bragg, Ireland and Collins2016) we did not properly consider this short-time, dissipation-induced irreversibility mechanism, and as we shall now show, the irreversibility mechanisms identified in Bragg et al. (Reference Bragg, Ireland and Collins2016) actually apply when $s/\unicode[STIX]{x1D70F}_{p}\not \ll 1$ .

After the initial transient stage, i.e. for $s^{\prime }/\unicode[STIX]{x1D70F}_{p}>1$ , the forward and backward solutions for $\boldsymbol{w}^{p}$ become (using (2.13) for the backward case)

(2.15) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{w}^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)\approx \frac{1}{\unicode[STIX]{x1D70F}_{p}}\int _{0}^{s^{\prime }}\text{e}^{-(s^{\prime }-s^{\prime \prime })/\unicode[STIX]{x1D70F}_{p}}\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(s^{\prime \prime }|\unicode[STIX]{x1D709},0),s^{\prime \prime })\,\text{d}s^{\prime \prime }, & \displaystyle\end{eqnarray}$$

(2.16) $$\begin{eqnarray}\displaystyle & \displaystyle \boldsymbol{w}^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)\approx \frac{1}{\unicode[STIX]{x1D70F}_{p}}\int _{-T_{0}}^{-s^{\prime }}\text{e}^{(s^{\prime }+s^{\prime \prime })/\unicode[STIX]{x1D70F}_{p}}\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(s^{\prime \prime }|\unicode[STIX]{x1D709},0),s^{\prime \prime })\,\text{d}s^{\prime \prime }. & \displaystyle\end{eqnarray}$$

As discussed earlier, if the PDFs of $\boldsymbol{w}^{p}$ are asymmetric, then ${\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)\neq {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)$ . In Bragg et al. (Reference Bragg, Ireland and Collins2016), by considering the behaviour of $\boldsymbol{w}^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)$ and $\boldsymbol{w}^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)$ we argued that there are two distinct physical mechanisms that generate such asymmetry in the PDFs of $\boldsymbol{w}^{p}$ . We refer the reader to that paper for detailed explanations; here we summarize. Let us first define $St_{r}(s^{\prime })\equiv \unicode[STIX]{x1D70F}_{p}/\unicode[STIX]{x1D70F}_{r}(s^{\prime })$ , where $\unicode[STIX]{x1D70F}_{r}(s^{\prime })$ is the eddy turnover time based upon the particle separation at time $s^{\prime }$ . If we project (2.15) and (2.16) onto $\boldsymbol{r}^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)$ and $\boldsymbol{r}^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)$ , respectively, then in the regime $St_{r}(s^{\prime })\ll 1$ we obtain

(2.17) $$\begin{eqnarray}\displaystyle & \displaystyle w_{\Vert }^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)\approx \unicode[STIX]{x0394}u_{\Vert }(\boldsymbol{r}^{p}(s^{\prime }|\unicode[STIX]{x1D709},0),s^{\prime }), & \displaystyle\end{eqnarray}$$

(2.18) $$\begin{eqnarray}\displaystyle & \displaystyle w_{\Vert }^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)\approx \unicode[STIX]{x0394}u_{\Vert }(\boldsymbol{r}^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0),-s^{\prime }). & \displaystyle\end{eqnarray}$$

As noted earlier, because of the trajectory conditioning, particles that are separating FIT at time $s^{\prime }$ have $w_{\Vert }^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)>0$ , whereas particles that are separating BIT at time $-s^{\prime }$ have $w_{\Vert }^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)<0$ . However, due to the dynamical fluxes in the turbulent velocity field, the PDF of $\unicode[STIX]{x0394}u_{\Vert }$ is negatively skewed in three dimensions. Consequently, we expect (statistically) that $|w_{\Vert }^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)|>|w_{\Vert }^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)|$ , i.e. the energy flux in three-dimensional turbulence causes particle pairs with $St_{r}\ll 1$ to be pushed together more strongly than apart. In Bragg et al. (Reference Bragg, Ireland and Collins2016) we referred to this as the local irreversibility mechanism (LIM), since it arises from the local (in a temporal sense) properties of $\unicode[STIX]{x0394}u_{\Vert }$ experienced by the particle pair. For fluid particles, the LIM is the only irreversibility mechanism since they do not experience the short-time dissipation source of irreversibility (because they do not experience a drag force) described earlier.

When $St_{r}(s^{\prime })\not \ll 1$ , the path-history contribution in (2.15) and (2.16) is important, and the inertial particle dynamics is temporally non-local. Since particles that are separating FIT at time $s^{\prime }$ have $w_{\Vert }^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)>0$ , then their separation would have been smaller in the past. On the other hand, since particles that are separating BIT at time $-s^{\prime }$ have $w_{\Vert }^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)<0$ , then their separation would have been larger in the past. Since $\unicode[STIX]{x0394}\boldsymbol{u}$ on average increases with increasing separation (true instantaneously in the dissipation range) then the path integral in (2.16) will be (statistically speaking) larger than that in (2.15). Consequently, we expect (statistically) that $|w_{\Vert }^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)|>|w_{\Vert }^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)|$ , just as for the LIM case. In Bragg et al. (Reference Bragg, Ireland and Collins2016) this was referred to as the non-local irreversibility mechanism (NLIM), since it arises from the non-local dynamics of the inertial particles, and operates when $St_{r}(s^{\prime })\not \ll 1$ provided that the particle separation lies in the range where the statistics of $\unicode[STIX]{x0394}\boldsymbol{u}$ depend upon separation (i.e. at sub-integral scales). As emphasized in Bragg et al. (Reference Bragg, Ireland and Collins2016), the NLIM does not arise from skewness in the underlying field $\unicode[STIX]{x0394}\boldsymbol{u}$ , and would operate even in a kinematic field where $\unicode[STIX]{x0394}\boldsymbol{u}$ has a symmetric PDF.

Whereas in the regime $s\ll \unicode[STIX]{x1D70F}_{p}$ , dissipation plays an explicit role in the irreversibility, its effect is implicit outside of this regime. Dissipation is implicit in the NLIM since the memory of the inertial particles (on which the NLIM depends) comes from their dissipative dynamics. However, dissipation is not the actual cause of irreversibility outside the regime $s\ll \unicode[STIX]{x1D70F}_{p}$ . This is demonstrated by the following argument: if the statistics of the field $\unicode[STIX]{x0394}\boldsymbol{u}$ were independent of separation, then the system described by (2.11) would be mathematically identical to the case of single inertial particles moving in a homogeneous turbulent flow field, in which the particle velocity PDFs are necessarily symmetric. For such a system, the statistics of (2.15) and (2.16) would be identical, and through (2.9) and (2.10) this would imply reversible dispersion. This is consistent with the known fact that the statistics of single-particle velocities in stationary, homogeneous turbulence are reversible (Falkovich et al. Reference Falkovich, Xu, Pumir, Bodenschatz, Biferale, Boffetta, Lanotte and Toschi2012). Thus, dissipation alone cannot break the symmetry of the PDF of $\boldsymbol{w}^{p}$ and therefore it is not, in and of itself, the cause of irreversibility outside the regime $s\ll \unicode[STIX]{x1D70F}_{p}$ .

Another important point is that the asymmetry of the PDF of $\boldsymbol{w}^{p}$ also plays a role in the regime $s\ll \unicode[STIX]{x1D70F}_{p}$ , where the explicit source of irreversibility comes from the particle dissipation. In the regime $s\ll \unicode[STIX]{x1D70F}_{p}$ we have the FIT behaviour $w_{\Vert }^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)\approx \text{e}^{-s^{\prime }/\unicode[STIX]{x1D70F}_{p}}w_{\Vert }^{p}(0|\unicode[STIX]{x1D709},0)$ , and BIT we have $w_{\Vert }^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)\approx \text{e}^{s^{\prime }/\unicode[STIX]{x1D70F}_{p}}w_{\Vert }^{p}(0|\unicode[STIX]{x1D709},0)$ . However, since FIT separating particles must have $w_{\Vert }^{p}(s^{\prime }|\unicode[STIX]{x1D709},0)>0$ , and BIT have $w_{\Vert }^{p}(-s^{\prime }|\unicode[STIX]{x1D709},0)<0$ , then even in the short-time limit, the irreversibility will be affected by the asymmetry of the PDF of $w_{\Vert }^{p}(0|\unicode[STIX]{x1D709},0)$ .

In addition to the role of dissipation in the short-time limit, another issue that was not thoroughly considered in Bragg et al. (Reference Bragg, Ireland and Collins2016) is the role of ‘preferential sampling’ on the irreversibility. Preferential sampling relates to the fact that due to the way inertial particles interact with the topology of the fluid velocity field, they do not uniformly sample the field $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r},t)$ , showing a tendency to avoid rotation dominated regions of the flow (Maxey Reference Maxey1987; Ireland, Bragg & Collins Reference Ireland, Bragg and Collins2016). This effect is implicit in both LIM and NLIM, because in either case, the fluid velocity difference driving the particle relative motion is $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t),t)$ , the statistics of which (conditioned on $\boldsymbol{r}^{p}(t)=\boldsymbol{r}$ ) will deviate from those of $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r},t)$ due to preferential sampling. Both the LIM and the NLIM would still generate irreversible dispersion even in the absence of preferential sampling, however, their strength will be quantitatively affected by preferential sampling. For example, in the local case, the strength of the dispersion irreversibility will be affected since the asymmetry in the PDF of $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t),t)$ may not be the same as the asymmetry in the PDF of $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r},t)$ . Similarly, in the non-local case, the strength of the dispersion irreversibility will be affected by the fact that the values of $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t),t)$ along the path history of the particle pair will be biased due to preferential sampling (similar to how the non-local clustering mechanism in the dissipation range is affected by the preferential sampling effect, see Bragg, Ireland & Collins Reference Bragg, Ireland and Collins2015b ). The only assumption we are making is that the skewness of the PDF of $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t),t)$ (more precisely, its longitudinal component) remains negative in three dimensions $\forall St$ , and that the moments of $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t),t)$ increase with increasing separation $\forall St$ . In other words, we assume that preferential sampling causes the statistics of $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t),t)$ to differ from those of $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r},t)$ quantitatively, but not qualitatively. These assumptions seem very reasonable, and will be validated using DNS data in § 3.

In Bragg et al. (Reference Bragg, Ireland and Collins2016), we used the arguments for the irreversibility mechanisms to predict that the mean-square particle separation should be faster BIT than FIT $\forall St_{r}$ , which was confirmed by DNS results. When considering the full PDFs ${\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)$ and ${\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)$ , the effect is more subtle. In particular, we must consider the behaviour for $r<\unicode[STIX]{x1D709}$ and $r>\unicode[STIX]{x1D709}$ separately.

As discussed earlier, for $r>\unicode[STIX]{x1D709}$ , ${\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)$ is governed by particles whose motion is dominated by $w_{\Vert }^{p}>0$ , whereas for $r>\unicode[STIX]{x1D709}$ , ${\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)$ is governed by particles whose motion was dominated by $w_{\Vert }^{p}<0$ . Since the PDFs for $w_{\Vert }^{p}$ are negatively skewed then we would expect that for $r>\unicode[STIX]{x1D709}$ , ${\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)>{\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)$ . However, for $r<\unicode[STIX]{x1D709}$ the relationship is inverted: for $r<\unicode[STIX]{x1D709}$ , ${\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)$ is governed by particles whose motion is dominated by $w_{\Vert }^{p}<0$ , whereas for $r<\unicode[STIX]{x1D709}$ , ${\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)$ is governed by particles whose motion was dominated by $w_{\Vert }^{p}>0$ . Since the PDFs for $w_{\Vert }^{p}$ are negatively skewed then we would expect that for $r<\unicode[STIX]{x1D709}$ , ${\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)>{\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)$ .

With respect to this prediction, three remarks are in order. First, $r>\unicode[STIX]{x1D709}$ and $r<\unicode[STIX]{x1D709}$ in the previous argument should not be taken too literally, since we would expect a transition region between the two predicted behaviours. Second, the prediction that for $r>\unicode[STIX]{x1D709}$ , ${\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)>{\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)$ is consistent with the findings in Bragg et al. (Reference Bragg, Ireland and Collins2016) for the mean-square separations, namely that the BIT mean-square separation should be greater than the FIT counterpart. The connection is that the mean-square separation is itself dominated by the behaviour of particles at $r>\unicode[STIX]{x1D709}$ , i.e. $\langle \Vert \boldsymbol{r}^{p}(\pm s)\Vert ^{2}\rangle _{\unicode[STIX]{x1D709}}\geqslant \unicode[STIX]{x1D709}^{2}$ . Third, for $\unicode[STIX]{x1D709}$ sufficiently small and/or $s$ sufficiently large, the prediction for $r<\unicode[STIX]{x1D709}$ may not hold due to the effect of particle pairs whose separation $\Vert \boldsymbol{r}^{p}\Vert$ has passed through a minimum on the interval $s^{\prime }\in [0,s]$ .

So far we have considered how, for a given $St$ , ${\mathcal{P}}^{{\mathcal{F}}}$ and ${\mathcal{P}}^{{\mathcal{B}}}$ might differ. Another important question is how ${\mathcal{P}}^{{\mathcal{F}}}$ and ${\mathcal{P}}^{{\mathcal{B}}}$ are each affected by changes in $St$ . To understand this, we consider the relative velocities of inertial particle governed by (2.11) at a given separation $\boldsymbol{r}$ (ignoring the initial condition and the $\boldsymbol{x}$ coordinate)

(2.19) $$\begin{eqnarray}\displaystyle \boldsymbol{w}^{p}(t|\boldsymbol{r})=\frac{1}{\unicode[STIX]{x1D70F}_{p}}\int _{0}^{t}\text{e}^{-(t-s)/\unicode[STIX]{x1D70F}_{p}}\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(s|\boldsymbol{r},t),s)\,\text{d}s, & & \displaystyle\end{eqnarray}$$

where $\boldsymbol{w}^{p}(t|\boldsymbol{r})$ denotes the value of $\boldsymbol{w}^{p}$ at time $t$ conditioned upon $\boldsymbol{r}^{p}(t)=\boldsymbol{r}$ . The solution (2.19) may be split up into its local and non-local contributions. (By local and non-local we are referring to the way $\boldsymbol{w}^{p}$ depends dynamically upon the field $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r},t)$ . There is another kinematic kind of non-locality: $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t|\boldsymbol{r},t),t)$ is kinematically non-local since $\boldsymbol{r}^{p}(t|\boldsymbol{r},t)$ itself depends upon a time integral of $\boldsymbol{w}^{p}$ because $\dot{\boldsymbol{r}}^{p}(t)\equiv \boldsymbol{w}^{p}(t)$ . In this kinematic sense, even fluid particle dynamics would be non-local (their dynamics is also non-local due to the non-local pressure in the Navier–Stokes equation, however this is a spatial non-locality, not temporal, which is the kind we are discussing here). In all of our discussion, by local and non-local we only refer to the temporal dependence of $\boldsymbol{w}^{p}$ on $\unicode[STIX]{x0394}\boldsymbol{u}$ , which is only non-local for inertial, not fluid particles.) The local contribution is simply the part that comes from the replacement $\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(s|\boldsymbol{r},t),s)\rightarrow \unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t|\boldsymbol{r},t),t)$ in (2.19), allowing us to write

(2.20) $$\begin{eqnarray}\boldsymbol{w}^{p}(t|\boldsymbol{r})=\underbrace{\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t|\boldsymbol{r},t),t)}_{\mathit{Local}}+\underbrace{\frac{1}{\unicode[STIX]{x1D70F}_{p}}\int _{0}^{t}\text{e}^{-(t-s)/\unicode[STIX]{x1D70F}_{p}}(\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(s|\boldsymbol{r},t),s)-\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(t|\boldsymbol{r},t),t))\,\text{d}s}_{\mathit{Non}\text{-}\mathit{Local}}.\end{eqnarray}$$

The local contribution is the same, regardless of whether the pair separation was larger or smaller in the past, but not so for the non-local contribution. Similar to the arguments for the non-local irreversibility mechanism, the non-local contribution in (2.20) will be larger for particles whose separation satisfies $\Vert \boldsymbol{r}^{p}(s|\boldsymbol{r},t)\Vert >r$ than for those satisfying $\Vert \boldsymbol{r}^{p}(s|\boldsymbol{r},t)\Vert <r$ for $s<t$ . The local contribution in (2.20) will therefore (statistically) be more significant for separating particles than approaching particles, and in this sense the former is less influenced by the effects of particle inertia than the latter (noting that only the local part survives in the limit of weak inertia). Coupled with the previous irreversibility arguments, this implies that for $r>\unicode[STIX]{x1D709}$ , ${\mathcal{P}}^{{\mathcal{B}}}$ should be more strongly affected by particle inertia than ${\mathcal{P}}^{{\mathcal{F}}}$ , but the opposite for $r<\unicode[STIX]{x1D709}$ .

2.2 Dependence of the PDFs on $s,\unicode[STIX]{x1D709},r$

We now turn to consider how ${\mathcal{P}}^{{\mathcal{F}}}$ and ${\mathcal{P}}^{{\mathcal{B}}}$ depend on $s,\unicode[STIX]{x1D709},r$ . An important quantity in this respect is the time scale $\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D709}}^{p}$ , defined as the time scale corresponding to the autocovariance $\langle w_{\Vert }^{p}(0)w_{\Vert }^{p}(s^{\prime })\rangle _{\unicode[STIX]{x1D709}}$ .

The first regime to consider is the ‘short-time regime’. Although formally this regime is $s/\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D709}}^{p}\ll 1$ , theoretical results based on this asymptotic limit are often accurate up to $s/\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D709}}^{p}\leqslant O(1)$ . From a physical perspective, the significance of this regime is that it corresponds to the regime where the dispersion behaviour is dominated by the particle relative velocities at the conditioning time $s=0$ .

We begin by considering the regime $\unicode[STIX]{x1D709}\ll \unicode[STIX]{x1D702}$ and $St\geqslant O(1)$ , for which a question of significant interest pertains to understanding how the presence of caustics in the inertial particle relative velocities in the dissipation range (Wilkinson & Mehlig Reference Wilkinson and Mehlig2005; Wilkinson, Mehlig & Bezuglyy Reference Wilkinson, Mehlig and Bezuglyy2006; Gustavsson & Mehlig Reference Gustavsson and Mehlig2011) influences the behaviour of ${\mathcal{P}}^{{\mathcal{B}}}$ . As discussed earlier, in the short-time regime $s\ll \unicode[STIX]{x1D70F}_{p}$ we have $\boldsymbol{w}^{p}(-s^{\prime })\approx \text{e}^{s^{\prime }/\unicode[STIX]{x1D70F}_{p}}\boldsymbol{w}^{p}(0)$ , and from this we obtain

(2.21) $$\begin{eqnarray}\displaystyle \boldsymbol{r}^{p}(-s)\approx \boldsymbol{r}^{p}(0)+{\mathcal{G}}(-s)\boldsymbol{w}^{p}(0), & & \displaystyle\end{eqnarray}$$

where ${\mathcal{G}}(-s)\equiv \unicode[STIX]{x1D70F}_{p}(1-\text{e}^{s/\unicode[STIX]{x1D70F}_{p}})$ . This approximation is valid provided that

(2.22) $$\begin{eqnarray}\displaystyle \Vert \text{e}^{s/\unicode[STIX]{x1D70F}_{p}}\boldsymbol{w}^{p}(0)\Vert \gg \left\Vert \frac{1}{\unicode[STIX]{x1D70F}_{p}}\int _{-s}^{0}\text{e}^{(s^{\prime }+s^{\prime \prime })/\unicode[STIX]{x1D70F}_{p}}\unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(s^{\prime }),s^{\prime })\,\text{d}s^{\prime }\right\Vert . & & \displaystyle\end{eqnarray}$$

Although this regime is typically only realized for $s\ll \unicode[STIX]{x1D70F}_{p}$ , in the caustic regions it may be valid for much longer times, say for $s\leqslant O(\unicode[STIX]{x1D70F}_{p})$ , because in the caustic regions $\Vert \boldsymbol{w}^{p}\Vert \gg \Vert \unicode[STIX]{x0394}\boldsymbol{u}\Vert$ . As such, (2.21) may be valid beyond the ‘short-time’ regime $s\ll \unicode[STIX]{x1D70F}_{p}$ , especially when $\boldsymbol{w}^{p}(0)$ corresponds to the large values described by the heavy tails of its PDF, tails that are much heavier than those for the PDF of $\unicode[STIX]{x0394}\boldsymbol{u}$ (see Bec et al. (Reference Bec, Biferale, Cencini, Lanotte and Toschi2010a ), Ireland et al. (Reference Ireland, Bragg and Collins2016)).

Using (2.21) we obtain

(2.23) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)\approx \langle \unicode[STIX]{x1D6FF}(\Vert \unicode[STIX]{x1D743}+{\mathcal{G}}(-s)\boldsymbol{w}^{p}(0)\Vert -r)\rangle _{\unicode[STIX]{x1D709}}. & & \displaystyle\end{eqnarray}$$

Although only formally valid for the regime $St\gg 1$ , DNS results show that for particle pairs in the dissipation range with $St\geqslant O(1)$ , the caustics in their relative velocities lead to the statistics of the longitudinal and perpendicular projections of $\boldsymbol{w}^{p}(0)$ being approximately equal (Ireland et al. Reference Ireland, Bragg and Collins2016), and using this we obtain

(2.24) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)\approx \langle \unicode[STIX]{x1D6FF}([\unicode[STIX]{x1D709}^{2}+2\unicode[STIX]{x1D709}{\mathcal{G}}(-s)w_{\Vert }^{p}(0)+3({\mathcal{G}}(-s)w_{\Vert }^{p}(0))^{2}]^{1/2}-r)\rangle _{\unicode[STIX]{x1D709}}. & & \displaystyle\end{eqnarray}$$

Since the random variable in the argument of $\unicode[STIX]{x1D6FF}(\cdot )$ in (2.24) is $w_{\Vert }^{p}(0)$ , let us define ${\mathcal{F}}(w_{\Vert }^{p}(0))\equiv [\unicode[STIX]{x1D709}^{2}+2{\mathcal{G}}(-s)\unicode[STIX]{x1D709}w_{\Vert }^{p}(0)+3({\mathcal{G}}(-s)w_{\Vert }^{p}(0))^{2}]^{1/2}-r$ . Then, using standard results on the properties of $\unicode[STIX]{x1D6FF}(\cdot )$ we write

(2.25) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)=\langle \unicode[STIX]{x1D6FF}({\mathcal{F}}(w_{\Vert }^{p}(0)))\rangle _{\unicode[STIX]{x1D709}}={\mathcal{C}}\mathop{\sum }_{{\mathcal{Z}}_{n}}|{\mathcal{F}}^{\prime }({\mathcal{Z}}_{n})|^{-1}\langle \unicode[STIX]{x1D6FF}(w_{\Vert }^{p}(0)-{\mathcal{Z}}_{n})\rangle _{\unicode[STIX]{x1D709}}, & & \displaystyle\end{eqnarray}$$

where ${\mathcal{C}}$ is a normalization constant to ensure $\int _{0}^{\infty }{\mathcal{P}}^{{\mathcal{B}}}\,\text{d}r=1$ , and ${\mathcal{Z}}_{n}$ is the $n\text{th}$ root of ${\mathcal{F}}$ , i.e. ${\mathcal{F}}({\mathcal{Z}}_{n})\equiv 0$ , of which there are two

(2.26) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{Z}}_{1}=\frac{-\unicode[STIX]{x1D709}+\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}{3{\mathcal{G}}(-s)}, & \displaystyle\end{eqnarray}$$

(2.27) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{Z}}_{2}=\frac{-\unicode[STIX]{x1D709}-\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}{3{\mathcal{G}}(-s)}. & \displaystyle\end{eqnarray}$$

The following argument shows that only ${\mathcal{Z}}_{1}$ is viable: in the limit $s\rightarrow 0$ , $|2{\mathcal{G}}(-s)\unicode[STIX]{x1D709}w_{\Vert }^{p}(0)|\gg |3({\mathcal{G}}(-s)w_{\Vert }^{p}(0))^{2}|$ , and solutions with $r>\unicode[STIX]{x1D709}$ must therefore correspond to $w_{\Vert }^{p}(0)<0$ , and solutions with $r<\unicode[STIX]{x1D709}$ must correspond to $w_{\Vert }^{p}(0)>0$ (recall that ${\mathcal{G}}(-s)\leqslant 0$ ). The root ${\mathcal{Z}}_{2}$ is not consistent with this; it implies that $r>\unicode[STIX]{x1D709}$ corresponds to $w_{\Vert }^{p}(0)>0$ . We also note that both roots imply that solutions for ${\mathcal{P}}^{{\mathcal{B}}}$ are only real for $r\geqslant \sqrt{2/3}\unicode[STIX]{x1D709}$ . This is because we assumed in (2.24) that the statistics of the longitudinal and perpendicular projections of $\boldsymbol{w}^{p}(0)$ are approximately equal in the caustic regions, and this gave rise to the factor $3$ . Let us now suppose that in a given realization they are not equal, and so we would replace $3$ in (2.24) with $\unicode[STIX]{x1D70E}$ , where $\unicode[STIX]{x1D70E}\equiv \Vert \boldsymbol{w}^{p}(0)\Vert ^{2}/[w_{\Vert }^{p}(0)]^{2}$ , such that $\unicode[STIX]{x1D70E}\geqslant 1$ . If we then consider the case $[\unicode[STIX]{x1D709}^{2}+2\unicode[STIX]{x1D709}{\mathcal{G}}(-s)w_{\Vert }^{p}(0)+\unicode[STIX]{x1D70E}({\mathcal{G}}(-s)w_{\Vert }^{p}(0))^{2}]^{1/2}=r$ , we obtain the solution for $w_{\Vert }^{p}(0)$

(2.28) $$\begin{eqnarray}\displaystyle w_{\Vert }^{p}(0)=\frac{-\unicode[STIX]{x1D709}\pm \sqrt{\unicode[STIX]{x1D709}^{2}(1-\unicode[STIX]{x1D70E})+\unicode[STIX]{x1D70E}r^{2}}}{\unicode[STIX]{x1D70E}{\mathcal{G}}(-s)}. & & \displaystyle\end{eqnarray}$$

This can be interpreted as the value(s) of $w_{\Vert }^{p}(0)$ that generate a particle-pair separation equal to $r$ at time $-s$ , conditioned on the separation being $\unicode[STIX]{x1D709}$ at time $s=0$ . The solution for $w_{\Vert }^{p}(0)$ is real only if $\sqrt{\unicode[STIX]{x1D709}^{2}(1-\unicode[STIX]{x1D70E})+\unicode[STIX]{x1D70E}r^{2}}\geqslant 0$ , i.e. $r\geqslant \unicode[STIX]{x1D709}\sqrt{(\unicode[STIX]{x1D70E}-1)/\unicode[STIX]{x1D70E}}$ . This shows that unless $\unicode[STIX]{x1D70E}=1$ , the minimum separation that a particle pair can reach is ${>}0$ , i.e. the particle trajectories cannot cross. But $\unicode[STIX]{x1D70E}=1$ corresponds to the case where the perpendicular component of the particle relative velocity at $s=0$ is zero. This makes sense; only if the perpendicular component of $\boldsymbol{w}^{p}(0)$ is identically zero can the two-particle trajectories intersect. If the perpendicular component of $\boldsymbol{w}^{p}(0)$ is finite, i.e. $\unicode[STIX]{x1D70E}>1$ , then the particle pair will reach a minimum separation, before they begin to move apart. The implication of this is that to correctly predict ${\mathcal{P}}^{{\mathcal{B}}}$ over the range $r\in [0,\unicode[STIX]{x1D709}]$ , we would need to account for the range of possible values of $\unicode[STIX]{x1D70E}$ (we will return to this point momentarily when we consider the case of fluid particles). We leave the incorporation of this additional complexity to future work, and focus here on the behaviour of separating particle pairs, i.e. ${\mathcal{P}}^{{\mathcal{B}}}$ for $r\geqslant \unicode[STIX]{x1D709}$ .

Following the previous arguments, using only the root ${\mathcal{Z}}_{1}$ , the result in (2.25) leads to

(2.29) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)={\mathcal{C}}|{\mathcal{G}}(-s)|^{-1}\frac{r}{\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}\langle \unicode[STIX]{x1D6FF}(w_{\Vert }^{p}(0)-{\mathcal{Z}}_{1})\rangle _{\unicode[STIX]{x1D709}}. & & \displaystyle\end{eqnarray}$$

With this approach, we have translated the problem of predicting ${\mathcal{P}}^{{\mathcal{B}}}$ to the problem of predicting $\langle \unicode[STIX]{x1D6FF}(w_{\Vert }^{p}(0)-{\mathcal{Z}}_{1})\rangle _{\unicode[STIX]{x1D709}}$ . However, the latter is simply the PDF of $w_{\Vert }^{p}(0)$ , concerning which a number of asymptotic predictions have been derived. The theoretical work by Gustavsson & Mehlig (Reference Gustavsson and Mehlig2011, Reference Gustavsson and Mehlig2014) predicts the asymptotic result

(2.30) $$\begin{eqnarray}\displaystyle \langle \unicode[STIX]{x1D6FF}(w_{\Vert }^{p}(0)-w_{\Vert })\rangle _{\unicode[STIX]{x1D709}}\propto (\unicode[STIX]{x1D709}/\unicode[STIX]{x1D70F}_{p})^{3-d_{2}}|w_{\Vert }|^{d_{2}-4}, & & \displaystyle\end{eqnarray}$$

where $d_{2}(St)$ is the correlation dimension and $w_{\Vert }$ is the phase-space variable corresponding to $w_{\Vert }^{p}(0)$ . The result (2.30) applies when caustics are present, except in the limits $|w_{\Vert }|\rightarrow 0$ and $|w_{\Vert }|\rightarrow \infty$ , the latter because the algebraic tail of the PDF is cut off at large $|w_{\Vert }|$ due to the finitude of the large-scale fluid velocity fluctuations. The result in (2.30) has recently been confirmed using DNS (Perrin & Jonker Reference Perrin and Jonker2015), for $St\geqslant O(1)$ , and separations in the dissipation range.

We may use (2.30) in (2.29) and obtain

(2.31) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)\propto \left(\frac{\unicode[STIX]{x1D709}|{\mathcal{G}}(-s)|}{\unicode[STIX]{x1D70F}_{p}}\right)^{3-d_{2}}\frac{r}{\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}|\unicode[STIX]{x1D709}-\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}|^{d_{2}-4}. & & \displaystyle\end{eqnarray}$$

Equation (2.31) predicts that ${\mathcal{P}}^{{\mathcal{B}}}$ grows with $s$ as $|{\mathcal{G}}(-s)|^{3-d_{2}}\equiv |\unicode[STIX]{x1D70F}_{p}(1-\text{e}^{s/\unicode[STIX]{x1D70F}_{p}})|^{3-d_{2}}$ , which is a very non-trivial behaviour since $d_{2}(St)$ takes on non-integer values. Indeed, in the limit $s/\unicode[STIX]{x1D70F}_{p}\rightarrow 0$ , the leading-order behaviour is $|{\mathcal{G}}(-s)|^{3-d_{2}}\sim s^{3-d_{2}}$ . An important point to note, however, is that the extent of $r$ -space over which (2.31) will be accurate will itself be a function of time, since (2.30), upon which (2.31) depends, does not describe the entire distribution of $w_{\Vert }^{p}(0)$ . The implication of this is that the proportionality constant in (2.31) is then itself implicitly a function of time. Therefore (2.31) can only be used to predict the $r$ -dependence, not the $s$ dependence of ${\mathcal{P}}^{{\mathcal{B}}}$ . However, as a test, we simulated (2.24) for the case where $w_{\Vert }^{p}(0)$ is governed entirely by (2.30), and the results matched (2.31) exactly, confirming that in this regime, the time dependence of ${\mathcal{P}}^{{\mathcal{B}}}$ is indeed given by $|{\mathcal{G}}(-s)|^{3-d_{2}}$ .

It is straightforward to show that the same analysis for the FIT case leads to the result

(2.32) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)\propto \left(\frac{\unicode[STIX]{x1D709}{\mathcal{G}}(s)}{\unicode[STIX]{x1D70F}_{p}}\right)^{3-d_{2}}\frac{r}{\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}|\unicode[STIX]{x1D709}-\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}|^{d_{2}-4}. & & \displaystyle\end{eqnarray}$$

The only difference between (2.31) and (2.32) is that the former contains $|{\mathcal{G}}(-s)|^{3-d_{2}}$ while the latter contains ${\mathcal{G}}^{3-d_{2}}(s)$ . These functions satisfy $|{\mathcal{G}}(-s)|^{3-d_{2}}/{\mathcal{G}}^{3-d_{2}}(s)\geqslant 1$ , and their difference captures the short-time irreversibility produced by the dissipation arising from the drag force on the particles. Consistent with our earlier arguments, these results predict that particles should separate faster BIT than FIT. However, the difference between (2.31) and (2.32) does not fully account for irreversibility in the short-time regime since we have used (2.30) which does not capture the asymmetry in the distribution $\langle \unicode[STIX]{x1D6FF}(w_{\Vert }^{p}(0)-w_{\Vert })\rangle _{\unicode[STIX]{x1D709}}$ . Unfortunately, there is at present no known way to predict the asymmetry of $\langle \unicode[STIX]{x1D6FF}(w_{\Vert }^{p}(0)-w_{\Vert })\rangle _{\unicode[STIX]{x1D709}}$ .

The results (2.31) and (2.32) predict that ${\mathcal{P}}^{{\mathcal{B}}}$ and ${\mathcal{P}}^{{\mathcal{F}}}$ will change with $Re_{\unicode[STIX]{x1D706}}$ through $d_{2}$ . Since it is known that for $St\gtrsim 1$ , $d_{2}$ decreases with increasing $Re_{\unicode[STIX]{x1D706}}$ (see Ireland et al. (Reference Ireland, Bragg and Collins2016), where $c_{1}\equiv 3-d_{2}$ ), then the tails of ${\mathcal{P}}^{{\mathcal{B}}}$ and ${\mathcal{P}}^{{\mathcal{F}}}$ will decay faster with increasing $Re_{\unicode[STIX]{x1D706}}$ . This is the opposite of what we would expect for fluid particles, as we shall soon discuss. It is important to note, however, that the functional forms in (2.31) and (2.32) are in another sense universal. That is, these functional forms apply even in the absence of any intermittency in the underlying fluid velocity field. This follows from the universality of (2.30) itself, which holds in both simple Gaussian flow fields (Gustavsson & Mehlig Reference Gustavsson and Mehlig2011) and also in Navier–Stokes turbulence (Perrin & Jonker Reference Perrin and Jonker2015). These results therefore suggest that extreme events (here characterized by large separations over short time scales) in the pair separation of inertial particles with $St\geqslant O(1)$ occur, not as a consequence of the intermittency in the small-scale turbulence itself, but as a consequence of their strongly non-local dynamics.

Next, we consider $St=0$ and arbitrary $\unicode[STIX]{x1D709}$ . In this case, following Batchelor (Reference Batchelor1952b ), we may use a short-time expansion in $s$ and obtain

(2.33) $$\begin{eqnarray}\displaystyle \boldsymbol{r}^{p}(-s)\approx \boldsymbol{r}^{p}(0)-s\unicode[STIX]{x0394}\boldsymbol{u}^{p}(0)+O(s^{2}), & & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x0394}\boldsymbol{u}^{p}(0)\equiv \unicode[STIX]{x0394}\boldsymbol{u}(\boldsymbol{r}^{p}(0),0)$ , leading to

(2.34) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)\approx \langle \unicode[STIX]{x1D6FF}(\Vert \unicode[STIX]{x1D743}-s\unicode[STIX]{x0394}\boldsymbol{u}^{p}(0)\Vert -r)\rangle _{\unicode[STIX]{x1D709}}. & & \displaystyle\end{eqnarray}$$

A significant complication compared to the $St\geqslant O(1)$ case is that the parallel and perpendicular components of $\unicode[STIX]{x0394}\boldsymbol{u}^{p}(0)$ are not statistically equivalent in general, or even approximately so. Therefore, if we write

(2.35) $$\begin{eqnarray}\displaystyle \Vert \unicode[STIX]{x0394}\boldsymbol{u}^{p}(0)\Vert ^{2}=\unicode[STIX]{x1D707}^{p}[\unicode[STIX]{x0394}u_{\Vert }^{p}(0)]^{2}, & & \displaystyle\end{eqnarray}$$

then the coefficient $\unicode[STIX]{x1D707}^{p}\geqslant 1$ will be random, varying significantly from one realization of the ensemble to another. In principle, we could capture the effect of $\unicode[STIX]{x1D707}^{p}$ on ${\mathcal{P}}^{{\mathcal{B}}}$ in the short-time limit by using (2.35) in (2.34) to obtain

(2.36) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)\approx \langle \unicode[STIX]{x1D6FF}([\unicode[STIX]{x1D709}^{2}-2s\unicode[STIX]{x1D709}\unicode[STIX]{x0394}u_{\Vert }^{p}(0)+\unicode[STIX]{x1D707}^{p}(s\unicode[STIX]{x0394}u_{\Vert }^{p}(0))^{2}]^{1/2}-r)\rangle _{\unicode[STIX]{x1D709}}, & & \displaystyle\end{eqnarray}$$

and then apply to this the decomposition

(2.37) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)\equiv \int _{\unicode[STIX]{x1D707}}{\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0,\unicode[STIX]{x1D707})\unicode[STIX]{x1D711}(\unicode[STIX]{x1D707})\,\text{d}\unicode[STIX]{x1D707}, & & \displaystyle\end{eqnarray}$$

where

(2.38) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0,\unicode[STIX]{x1D707})=\langle \unicode[STIX]{x1D6FF}([\unicode[STIX]{x1D709}^{2}-2s\unicode[STIX]{x1D709}\unicode[STIX]{x0394}u_{\Vert }^{p}(0)+\unicode[STIX]{x1D707}^{p}(s\unicode[STIX]{x0394}u_{\Vert }^{p}(0))^{2}]^{1/2}-r)\rangle _{\unicode[STIX]{x1D709},\unicode[STIX]{x1D707}^{p}=\unicode[STIX]{x1D707}} & \displaystyle\end{eqnarray}$$

(2.39) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D711}(\unicode[STIX]{x1D707})\equiv \langle \unicode[STIX]{x1D6FF}(\unicode[STIX]{x1D707}^{p}-\unicode[STIX]{x1D707})\rangle , & \displaystyle\end{eqnarray}$$

$\unicode[STIX]{x1D707}$ being the sample-space variable corresponding to $\unicode[STIX]{x1D707}^{p}$ . In general we cannot evaluate (2.37) since we do not know $\unicode[STIX]{x1D711}(\unicode[STIX]{x1D707})$ . However, from (2.33) we find that to leading order in $s$ , $\Vert \boldsymbol{r}(-s)\Vert ^{2}-\unicode[STIX]{x1D709}^{2}\propto s\unicode[STIX]{x0394}u_{\Vert }^{p}(0)$ , and we might expect that $\unicode[STIX]{x0394}u_{\Vert }^{p}(0)$ is typically larger in regions where $\unicode[STIX]{x1D707}=O(1)$ than in regions where $\unicode[STIX]{x1D707}\gg 1$ . For example in regions of strong, quasi-solid body rotation where $\unicode[STIX]{x1D707}\gg 1$ , $\unicode[STIX]{x0394}u_{\Vert }^{p}(0)$ would be very small compared with its typical values in regions dominated by extensional straining where $\unicode[STIX]{x1D707}=O(1)$ . This then implies that ${\mathcal{P}}^{{\mathcal{B}}}$ should be dominated by the contributions with $\unicode[STIX]{x1D707}=O(1)$ in the limit $s\rightarrow 0$ . As an approximation we assume that the statistics of the parallel and perpendicular components of $\unicode[STIX]{x0394}\boldsymbol{u}^{p}(0)$ are equal, corresponding to assuming $\unicode[STIX]{x1D711}(\unicode[STIX]{x1D707})=\unicode[STIX]{x1D6FF}(3-\unicode[STIX]{x1D707})$ in (2.37). Then, following the same procedure as was used to derive the result for the $St\geqslant O(1)$ case, we obtain

(2.40) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)={\mathcal{C}}\frac{r}{s\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}\left\langle \unicode[STIX]{x1D6FF}(\unicode[STIX]{x0394}u_{\Vert }^{p}(0)-{\mathcal{Z}})\right\rangle _{\unicode[STIX]{x1D709}}, & & \displaystyle\end{eqnarray}$$

where

(2.41) $$\begin{eqnarray}\displaystyle {\mathcal{Z}}=\frac{\unicode[STIX]{x1D709}-\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}{3s}. & & \displaystyle\end{eqnarray}$$

Since we are considering a statistically stationary system where the fluid particles are fully mixed (globally), then $\langle \unicode[STIX]{x1D6FF}(\unicode[STIX]{x0394}u_{\Vert }^{p}(0)-{\mathcal{Z}})\rangle _{\unicode[STIX]{x1D709}}$ is identical to the PDF of the fluid velocity differences at scale $\unicode[STIX]{x1D709}$ , i.e. $\langle \unicode[STIX]{x1D6FF}(\unicode[STIX]{x0394}u_{\Vert }^{p}(0)-{\mathcal{Z}})\rangle _{\unicode[STIX]{x1D709}}=\langle \unicode[STIX]{x1D6FF}(\unicode[STIX]{x0394}u_{\Vert }(\unicode[STIX]{x1D709},0)-{\mathcal{Z}})\rangle$ , and this can be described using the multifractal formalism. In order to capture the asymmetry of this PDF, we use the multifractal model from Chevillard et al. (Reference Chevillard, Castaing, Arneodo, Lvque, Pinton and Roux2012)

(2.42) $$\begin{eqnarray}\left\langle \unicode[STIX]{x1D6FF}\left(\unicode[STIX]{x0394}u_{\Vert }(\unicode[STIX]{x1D709},0)-{\mathcal{Z}}\right)\right\rangle _{\unicode[STIX]{x1D709}}=\int _{h_{min}}^{h_{max}}\frac{1}{\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}}{\mathcal{P}}_{h,\unicode[STIX]{x1D709}}{\mathcal{P}}_{\unicode[STIX]{x1D6FF}}\left(\frac{{\mathcal{Z}}}{\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}}\right)\,\text{d}h,\end{eqnarray}$$

where

(2.43) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}(h)\equiv (\unicode[STIX]{x1D709}/L)^{h}(1+(\unicode[STIX]{x1D702}_{h}/\unicode[STIX]{x1D709})^{2})^{(h-1)/2}, & \displaystyle\end{eqnarray}$$

(2.44) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}_{h,\unicode[STIX]{x1D709}}(h)\equiv A^{-1}(\unicode[STIX]{x1D709}/L)^{1-{\mathcal{D}}_{h}}(1+(\unicode[STIX]{x1D702}_{h}/\unicode[STIX]{x1D709})^{2})^{(1-{\mathcal{D}}_{h})/2}, & \displaystyle\end{eqnarray}$$

(2.45) $$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}_{\unicode[STIX]{x1D6FF}}(y)\equiv \frac{1}{\sqrt{2\unicode[STIX]{x03C0}}}(\text{e}^{-y^{2}/2}-\unicode[STIX]{x1D706}_{\unicode[STIX]{x1D709}}y(y^{2}-1)\text{e}^{-y^{2}/2a^{2}}), & \displaystyle\end{eqnarray}$$

$A$ is a normalization constant, $h_{min}=0$ , $h_{max}=1$ , ${\mathcal{D}}_{h}=1-(h-c_{1})^{2}/2c_{2}$ , $a=\sqrt{9/10}$ , $c_{2}=1/40$ , $c_{1}=(1/3)+(3c_{2}/2)$ , $\unicode[STIX]{x1D70E}=\sqrt{2}u^{\prime }$ , $\unicode[STIX]{x1D702}_{h}=L(Re/R^{\ast })^{-1/(h+1)}$ , $L$ is the integral length scale, $Re=\unicode[STIX]{x1D70E}L/\unicode[STIX]{x1D708}$ and $R^{\ast }=52$ . The parameter $\unicode[STIX]{x1D706}_{\unicode[STIX]{x1D709}}$ allows the model to capture the asymmetry of $\langle \unicode[STIX]{x1D6FF}(\unicode[STIX]{x0394}u_{\Vert }(\unicode[STIX]{x1D709},0)-{\mathcal{Z}})\rangle$ , and is given by $\unicode[STIX]{x1D706}_{\unicode[STIX]{x1D709}}=(2\langle \unicode[STIX]{x1D716}\rangle \unicode[STIX]{x1D709}/15-\unicode[STIX]{x1D708}\unicode[STIX]{x1D735}_{\unicode[STIX]{x1D709}}\langle [\unicode[STIX]{x0394}u_{\Vert }(\unicode[STIX]{x1D709},0)]^{2}\rangle )/\langle \unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}^{3}\rangle$ , where $\langle \unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}^{3}\rangle \equiv \int _{h_{min}}^{h_{max}}\unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}^{3}(h){\mathcal{P}}_{h,\unicode[STIX]{x1D709}}(h)\,\text{d}h$ .

In the multifractal theory, $h$ plays the role of the scaling exponent, describing how $\unicode[STIX]{x0394}u_{\Vert }(\unicode[STIX]{x1D709},0)$ scales with $\unicode[STIX]{x1D709}$ on a given fractal subset of the system. The exponent $h$ has the PDF ${\mathcal{P}}_{h,\unicode[STIX]{x1D709}}(h)$ , and the associated singularity spectrum is ${\mathcal{D}}_{h}$ . The result in (2.42) captures the transition from an essentially Gaussian PDF for $\unicode[STIX]{x0394}u_{\Vert }(\unicode[STIX]{x1D709},0)$ at $\unicode[STIX]{x1D709}=L$ , to a highly non-Gaussian PDF at $\unicode[STIX]{x1D709}=\unicode[STIX]{x1D702}$ , exhibiting heavy, stretched exponential-type decaying tails. It is important to note that (2.42) is based upon the ‘log-normal’ assumption (Kolmogorov Reference Kolmogorov1962), which is known to lead to unphysical descriptions of the most extreme fluctuations of the fluid velocity field (see Frisch Reference Frisch1995). As a result, using (2.42) in (2.40) could lead to a failure to correctly describe the dispersion of particle pairs that separate extremely fast. Nevertheless, the use of the multifractal model of Chevillard et al. (Reference Chevillard, Castaing, Arneodo, Lvque, Pinton and Roux2012) is just a choice; if future multifractal models for $\langle \unicode[STIX]{x1D6FF}(\unicode[STIX]{x0394}u_{\Vert }(\unicode[STIX]{x1D709},0)-{\mathcal{Z}})\rangle$ emerge that capture both the asymmetry of the PDF and go beyond the log-normal assumption, they could be used in (2.40) instead. We refer the reader to Benzi et al. (Reference Benzi, Biferale, Paladin, Vulpiani and Vergassola1991), Frisch (Reference Frisch1995), Boffetta, Mazzino & Vulpiani (Reference Boffetta, Mazzino and Vulpiani2008) for detailed discussions of the multifractal formalism in turbulence.

Using (2.42) in (2.40) we obtain

(2.46) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)=\int _{h_{min}}^{h_{max}}\frac{{\mathcal{C}}r{\mathcal{P}}_{h,\unicode[STIX]{x1D709}}}{\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}s\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}{\mathcal{P}}_{\unicode[STIX]{x1D6FF}}\left(\frac{\unicode[STIX]{x1D709}-\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}{3\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}s}\right)\,\text{d}h. & & \displaystyle\end{eqnarray}$$

As expected, the result in (2.46) shows that for $r=O(\unicode[STIX]{x1D709})$ , the shape of ${\mathcal{P}}^{{\mathcal{B}}}$ is affected by the initial condition, but for $r\gg \unicode[STIX]{x1D709}$ , the shape of ${\mathcal{P}}^{{\mathcal{B}}}$ reflects the shape of the underlying PDF for $\unicode[STIX]{x0394}u_{\Vert }(\unicode[STIX]{x1D709},t)$ . It therefore also describes the variation in the behaviour of ${\mathcal{P}}^{{\mathcal{B}}}$ in the short-time regime as $\unicode[STIX]{x1D709}$ is varied, transitioning from an approximately Gaussian PDF for $\unicode[STIX]{x1D709}=O(L)$ to a highly non-Gaussian PDF for $\unicode[STIX]{x1D709}\leqslant O(\unicode[STIX]{x1D702})$ .

The corresponding prediction for ${\mathcal{P}}^{{\mathcal{F}}}$ is

(2.47) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)=\int _{h_{min}}^{h_{max}}\frac{{\mathcal{C}}r{\mathcal{P}}_{h,\unicode[STIX]{x1D709}}}{\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}s\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}{\mathcal{P}}_{\unicode[STIX]{x1D6FF}}\left(\frac{-\unicode[STIX]{x1D709}+\sqrt{3r^{2}-2\unicode[STIX]{x1D709}^{2}}}{3\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6FD}_{\unicode[STIX]{x1D709}}s}\right)\,\text{d}h. & & \displaystyle\end{eqnarray}$$

Since ${\mathcal{P}}_{\unicode[STIX]{x1D6FF}}$ is asymmetric, it follows that (2.46) and (2.47) predict ${\mathcal{P}}^{{\mathcal{B}}}\neq {\mathcal{P}}^{{\mathcal{F}}}$ , consistent with our earlier arguments that the irreversibility of fluid particle-pair dispersion arises because the PDF of $\unicode[STIX]{x0394}u_{\Vert }$ is asymmetric. More specifically, ${\mathcal{P}}_{\unicode[STIX]{x1D6FF}}$ gives rise to a negatively skewed PDF for $\unicode[STIX]{x0394}u_{\Vert }$ , such that (2.46) and (2.47) predict the fluid particles separate faster BIT than FIT, in agreement with the arguments of § 2.1.

The multifractal prediction for the PDF of $\unicode[STIX]{x0394}u_{\Vert }$ captures the enhanced intermittency as $Re_{\unicode[STIX]{x1D706}}$ is increased. Through (2.46) this leads to the prediction that the tails of ${\mathcal{P}}^{{\mathcal{B}}}$ decay more slowly as $Re_{\unicode[STIX]{x1D706}}$ is increased, as expected. What is interesting is that this is in contrast to the behaviour of $St\gtrsim 1$ particles, for which, according to (2.31), the tails of ${\mathcal{P}}^{{\mathcal{B}}}$ decay faster as $Re_{\unicode[STIX]{x1D706}}$ is increased.

We now consider the regime $s/\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D709}}^{p}\gg 1$ (but finite), in which most of the particles will be found in the inertial range of the turbulence if $Re_{\unicode[STIX]{x1D706}}\rightarrow \infty$ . For fluid particles in this regime, a famous prediction of Richardson is that ${\mathcal{P}}^{{\mathcal{F}}}$ should be independent of $\unicode[STIX]{x1D709}$ with the distribution (see Salazar & Collins (Reference Salazar and Collins2009))

(2.48) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)\propto \frac{r^{2}}{\langle \Vert \boldsymbol{r}^{p}(s)\Vert ^{2}\rangle _{\unicode[STIX]{x1D709}}^{3/2}}\exp \left[-\frac{{\mathcal{A}}^{{\mathcal{F}}}r^{2/3}}{\langle \Vert \boldsymbol{r}^{p}(s)\Vert ^{2}\rangle _{\unicode[STIX]{x1D709}}^{1/3}}\right], & & \displaystyle\end{eqnarray}$$

where $\langle \Vert \boldsymbol{r}^{p}(s)\Vert ^{2}\rangle _{\unicode[STIX]{x1D709}}\propto \langle \unicode[STIX]{x1D716}\rangle s^{3}$ , and ${\mathcal{A}}^{{\mathcal{F}}}>0$ is a constant. The most recent tests of the accuracy of (2.48) have demonstrated that it describes ${\mathcal{P}}^{{\mathcal{F}}}$ quite well for fluid particles at long times when the particles are in the inertial range (Bitane, Homann & Bec Reference Bitane, Homann and Bec2013; Biferale et al. Reference Biferale, Lanotte, Scatamacchia and Toschi2014), except for very small and very large $r$ , which correspond to extreme events in the pair-separation process (Scatamacchia, Biferle & Toschi Reference Scatamacchia, Biferle and Toschi2012; Biferale et al. Reference Biferale, Lanotte, Scatamacchia and Toschi2014). Under the phenomenology employed by Richardson we also have

(2.49) $$\begin{eqnarray}\displaystyle {\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)\propto \frac{r^{2}}{\langle \Vert \boldsymbol{r}^{p}(-s)\Vert ^{2}\rangle _{\unicode[STIX]{x1D709}}^{3/2}}\exp \left[-\frac{{\mathcal{A}}^{{\mathcal{B}}}r^{2/3}}{\langle \Vert \boldsymbol{r}^{p}(-s)\Vert ^{2}\rangle _{\unicode[STIX]{x1D709}}^{1/3}}\right], & & \displaystyle\end{eqnarray}$$

where $\langle \Vert \boldsymbol{r}^{p}(-s)\Vert ^{2}\rangle _{\unicode[STIX]{x1D709}}\propto \langle \unicode[STIX]{x1D716}\rangle s^{3}$ , and ${\mathcal{A}}^{{\mathcal{B}}}>0$ is a constant.

The form of Richardson’s PDF follows from the assumption that the pair separation is governed by a diffusion process, with a diffusion coefficient

(2.50) $$\begin{eqnarray}\displaystyle {\mathcal{D}}(r)\propto r^{4/3}, & & \displaystyle\end{eqnarray}$$

which can also be arrived at using Kolmogorov’s 1941 (K41) scaling. Let us now consider how the shape of the Richardson PDF might be modified for inertial particles through the effects of the inertia on the scaling of ${\mathcal{D}}(r)$ . In the regime $St_{r}\ll 1$ , it is straightforward to show that ${\mathcal{D}}(r)\propto r^{4/3}+O[St_{r}]$ (Bragg, Ireland & Collins Reference Bragg, Ireland and Collins2015a ). For $St_{r}\geqslant O(1)$ , it is more difficult to determine the scaling of ${\mathcal{D}}$ since the inertial particle relative velocity statistics in the inertial range do not in general scale with $r$ as simple power laws. Results from PDF theories for inertial particles predict that (Zaichik & Alipchenkov Reference Zaichik and Alipchenkov2009; Bragg & Collins Reference Bragg and Collins2014a )

(2.51) $$\begin{eqnarray}\displaystyle {\mathcal{D}}(r)=\unicode[STIX]{x1D70F}_{p}\unicode[STIX]{x1D706}+\unicode[STIX]{x1D70F}_{p}\langle w_{\Vert }^{p}(t)w_{\Vert }^{p}(t)\rangle _{r}, & & \displaystyle\end{eqnarray}$$

where

(2.52) $$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D706}(r)=\frac{1}{St_{r}(1+St_{r})}\langle \unicode[STIX]{x0394}u_{\Vert }(r,t)\unicode[STIX]{x0394}u_{\Vert }(r,t)\rangle , & & \displaystyle\end{eqnarray}$$

is a diffusion coefficient that describes the effect of $\unicode[STIX]{x0394}u_{\Vert }(r,t)$ on the transport of particles in $r$ -space (Bragg & Collins Reference Bragg and Collins2014a ). When $St_{r}\ll 1$ , $\unicode[STIX]{x1D706}\propto r^{4/3}$ , but when $St_{r}\gg 1$ $\unicode[STIX]{x1D706}\propto r^{2}$ .

In the inertial range, the effect of the non-local inertial particle dynamics on $\langle w_{\Vert }^{p}w_{\Vert }^{p}\rangle _{r}$ is quite weak when $St_{r}\leqslant O(1)$ compared with the effect of inertial filtering (Ireland et al. Reference Ireland, Bragg and Collins2016). We may then approximate $\langle w_{\Vert }^{p}w_{\Vert }^{p}\rangle _{r}$ based on its local behaviour (Zaichik & Alipchenkov Reference Zaichik and Alipchenkov2009; Bragg & Collins Reference Bragg and Collins2014b )

(2.53) $$\begin{eqnarray}\displaystyle \langle w_{\Vert }^{p}(t)w_{\Vert }^{p}(t)\rangle _{r}\approx \frac{1}{1+St_{r}}\langle \unicode[STIX]{x0394}u_{\Vert }(r,t)\unicode[STIX]{x0394}u_{\Vert }(r,t)\rangle . & & \displaystyle\end{eqnarray}$$

Provided then that $St_{r}\leqslant O(1)$ in the inertial range, substituting (2.53) in (2.51) gives

(2.54) $$\begin{eqnarray}\displaystyle {\mathcal{D}}(r)\approx \unicode[STIX]{x1D70F}_{r}\langle \unicode[STIX]{x0394}u_{\Vert }(r,t)\unicode[STIX]{x0394}u_{\Vert }(r,t)\rangle , & & \displaystyle\end{eqnarray}$$

which becomes ${\mathcal{D}}(r)\propto r^{4/3}$ under K41 scaling. These considerations therefore suggest that in the regime for which Richardson’s PDF shape holds true for fluid particles, the same PDF shape should also apply even up to $St_{r}=O(1)$ . Of course we have neglected to consider here the fact that Richardson’s Markovian assumption will in principle be even less applicable for finite $St_{r}$ than it is for fluid particles. We will not consider this here, but it is an issue that we shall address in a future publication.

Finally, in the limit $s/\unicode[STIX]{x1D70F}_{\unicode[STIX]{x1D709}}^{p}\rightarrow \infty$ , we expect the state of the system to be independent of its state at $s=0$ , giving the time-independent distributions (In (2.56) we have assumed that the single time probability measures of the system are time invariant on the time interval $(-\infty ,+\infty )$ . This is not necessarily true since any real dynamical system will pass through an initial transient regime. However, in reality we expect the asymptotic scaling (2.56) to be realized for $-s\ll -\unicode[STIX]{x1D70F}_{I}$ , where $\unicode[STIX]{x1D70F}_{I}$ is the integral time scale for the system. Therefore, if we denote the initial time of the system by $-T_{0}$ (e.g. the time at which the experiment began), then provided $-T_{0}\lll -\unicode[STIX]{x1D70F}_{I}$ , the result in (2.56) makes sense for $-T_{0}\ll -s\ll -\unicode[STIX]{x1D70F}_{I}$ .)

(2.55) $$\begin{eqnarray}\displaystyle & \displaystyle \lim _{s\rightarrow \infty }{\mathcal{P}}^{{\mathcal{F}}}(r,s|\unicode[STIX]{x1D709},0)\propto r^{2}g^{\infty }(r), & \displaystyle\end{eqnarray}$$

(2.56) $$\begin{eqnarray}\displaystyle & \displaystyle \lim _{s\rightarrow \infty }{\mathcal{P}}^{{\mathcal{B}}}(r,-s|\unicode[STIX]{x1D709},0)\propto r^{2}g^{\infty }(r), & \displaystyle\end{eqnarray}$$

where $g^{\infty }(r)$ is the stationary form of the radial distribution function (RDF). Since $g^{\infty }(r)$ is strongly dependent upon $St$ , the asymptotic state of the particle pairs will vary with $St$ . Furthermore, for $St>0$ , the asymptotic forms of ${\mathcal{P}}^{{\mathcal{F}}}$ and ${\mathcal{P}}^{{\mathcal{B}}}$ will in principle depend upon $Re_{\unicode[STIX]{x1D706}}$ through $g^{\infty }(r)$ , e.g. see Ireland et al. (Reference Ireland, Bragg and Collins2016).