3.1.1. First convective correction: uniform flow

To see why a straightforward expansion of the form $\varTheta = \varTheta _0 + {\textit{Pe}}\,\varTheta _1 + \cdots$ , $\varTheta _0$ being the diffusion solution previously obtained, does not work, we examine the relative importance of the convection and diffusion terms in (2.1). For small ${\textit{Pe}}$ , convection only becomes important at distances from the particle much larger than its largest dimension. At such distances, the particle appears as a point source, so one may use $\varTheta _0 \approx {1}/{r}$ for obtaining estimates. Doing so, one finds ${\nabla} ^2 \varTheta \sim {1}/{r^3}$ and $\boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta \sim \hat {\boldsymbol{U}} \boldsymbol{\cdot }\boldsymbol{\nabla } ({1}/{r}) \sim {1}/{r^2}$ , where we have used $\boldsymbol{u} \approx \hat {\boldsymbol{U}}$ , $\hat {\boldsymbol{U}}$ being the unit vector along the ambient flow. The ratio of convection to diffusion terms is given by ${{\textit{Pe}}(\boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta )}/({{\nabla} ^2 \varTheta }) \sim ({{\textit{Pe}}/r^2})/({1/r^3}) = {\textit{Pe}}\,r$ , implying that convection becomes important when $r \sim O({\textit{Pe}}^{-1})$ or, in dimensional terms, at distances of $O(D/U)$ , this being the convective screening length. This contradicts the assumption implicit in the small- ${\textit{Pe}}$ expansion given previously, of convection being a perturbative influence everywhere in the domain. The implication is that the diffusion solution is a valid leading order approximation and convection a perturbative effect, only for $r \ll O({\textit{Pe}}^{-1})$ , corresponding to the so-called inner region. However, one has to account for convection at leading order for $r \gtrsim O({\textit{Pe}}^{-1})$ , the so-called outer region. The characteristic scale in the outer region is not $a$ (the particle or drop is a point source) but $D/U$ . Not accounting for this change in scale and using the above-mentioned regular expansion would, at $O({\textit{Pe}})$ , lead to $\varTheta _1$ satisfying

(3.13) \begin{equation} {\nabla} ^2 \varTheta _1 = \boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta _0. \end{equation}

Here, $\varTheta _1 = 0$ at the particle (drop) surface and must also decay in the far field. As noted already, ${\boldsymbol{u}} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta _0 \sim O(1/r^2)$ for $r \gg 1$ , and therefore, ${\nabla} ^2 \varTheta _1 \sim O(1/r^2)$ , implying that $\varTheta _1$ asymptotes to an $r$ -independent limiting form at large distances. The resulting conflict with the far-field decay requirement is Whitehead’s paradox.

The paradox is resolved in a matched asymptotic expansions approach, since the expansion $\varTheta = \varTheta _0 + Pe\,\varTheta _1 + \cdots$ is now restricted to the inner region and the far-field decay requirement replaced by the requirement of matching with the near-field limiting form of an outer-region expansion (described later), valid at large distances. The interpretation as an inner-region expansion does not, however, mean that the expansion proceeds in powers of ${\textit{Pe}}$ , as for a regular perturbation problem – the previous matching requirement leads to logarithmic corrections and the next term ( $\varTheta _2$ ) in the expansion is, in fact, $O({\textit{Pe}}^2\ln {\textit{Pe}})$ ; see § 3.1.3. Now, for all fore–aft symmetric particles, $\varTheta _0$ is an even function of $\boldsymbol{x}$ , satisfying $\varTheta _0(\boldsymbol{x}) = \varTheta _0(-\boldsymbol{x})$ . The solution of (3.13) may be written as $\varTheta _1= \varTheta _{1h} + \varTheta _{1p}$ , where the particular solution, $\varTheta _{1p}$ , is an odd function of $\boldsymbol{x}$ due to $\boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta _0$ being odd, which in turn is due to the evenness of both $\varTheta _0$ and the Stokes velocity field $\boldsymbol{u}$ associated with a translating fore–aft symmetric particle. Integration over the particle (drop) surface then leads to $\varTheta _{1p}$ not contributing to the transport rate, despite being the source of the above-mentioned paradox. The $O({\textit{Pe}})$ correction comes from the homogeneous solution, $\varTheta _{1h}$ , which, on account of the trivial surface boundary condition, is driven by the matching requirement with the outer-region expansion.

To analyse the outer region, we recognise the change in characteristic scale from $a$ to $D/U$ . Using $D/U$ to scale distances, we have the position vector in the outer region as $\tilde {\boldsymbol{x}} = {\textit{Pe}}\,{\boldsymbol{x}}$ , $\boldsymbol{x}$ being the position vector in (2.1) scaled using $a$ . Thus, $O(1)$ distances in the outer region correspond to asymptotically large distances of $O({\textit{Pe}}^{-1})$ in inner region units. Rewriting (2.1) in outer variables, the outer-region scalar field, $\tilde {\varTheta }(\tilde {x})$ , satisfies

(3.14) \begin{equation} \bigl[\hat {\boldsymbol{U}} + O({\textit{Pe}})\bigr] \boldsymbol{\cdot }\tilde {\boldsymbol{\nabla }} \tilde {\varTheta } - \tilde {\boldsymbol{\nabla} }^2 \tilde {\varTheta } ={\textit{Pe}}\left[\bar {Q}_0 + O({\textit{Pe}})\right]\delta (\tilde {\boldsymbol{x}}), \end{equation}

where the tilde’s denote outer-region quantities. Now, $\boldsymbol{u}$ in (2.1) comprises the ambient flow ( $\hat {\boldsymbol{U}}$ ) and the $O(1/r)$ Stokeslet field, the latter being $O({\textit{Pe}})$ at distances of $O({\textit{Pe}}^{-1})$ , as indicated. Similar to (3.5), the particle (drop) has been incorporated as a point source. Strictly speaking, the source strength is unknown in the presence of convection. However, as indicated on the right-hand side of (3.14), the convective correction is only $O({\textit{Pe}})$ and may be neglected at leading order, so that the source strength is the non-dimensional diffusive transport rate ( $\bar {Q}_0$ ), whose expressions for spheres and spheroids were already obtained.

Equation (3.14) suggests an outer-region expansion of the form $\tilde {\varTheta } = Pe\, \tilde {\varTheta }_0 + {\textit{Pe}}^2\,\tilde {\varTheta }_1 + \cdots$ , with neither the explicit form of the Stokeslet field nor the convective correction to the source strength being needed to determine $\tilde {\varTheta }_0$ , which satisfies

(3.15) \begin{equation} \hat {\boldsymbol{U}} \boldsymbol{\cdot }\tilde {\boldsymbol{\nabla }} \tilde {\varTheta }_0 - \tilde {\boldsymbol{\nabla} }^2 \tilde {\varTheta }_0 = \bar {Q}_0\delta (\boldsymbol{x}). \end{equation}

Calculation of the scalar transport rate proceeds more conveniently in Fourier space. Defining $\hat {\varTheta }_0(\boldsymbol{k}) = \int \tilde {\varTheta }_0 e^{-2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}} \,\text{d}\boldsymbol{\tilde x}$ and Fourier transforming (3.15), one obtains

(3.16) \begin{align} & 2\pi \mathrm{i} (\boldsymbol{k} \boldsymbol{\cdot }\hat {\boldsymbol{U}}) \hat {\varTheta }_0 + 4\pi ^2 k^2 \hat {\varTheta }_0 = \bar {Q}_0 \end{align}

(3.17) \begin{align} \Rightarrow & \hat {\varTheta }_0 = \frac {\bar {Q}_0}{2\pi \mathrm{i} k_U + 4\pi ^2 k^2}, \end{align}

with $k_U \equiv \boldsymbol{k} \boldsymbol{\cdot }\hat {\boldsymbol{U}}$ . The outer-region scalar field in physical space, at leading order, is

(3.18) \begin{align} \tilde {\varTheta }(\tilde {\boldsymbol{x}}) = {\textit{Pe}} \displaystyle \int e^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}}\frac {\bar {Q}_0}{2\pi \mathrm{i} k_U + 4\pi ^2 k^2} \,\text{d}\boldsymbol{k}, \end{align}

where $\boldsymbol{k}$ is scaled by $a^{-1}{\textit{Pe}}$ . The interest is in the limiting form of $\tilde {\varTheta }$ for $\tilde {\boldsymbol{x}} \rightarrow 0$ , which must match up to the far-field ( $\boldsymbol{x} \rightarrow \infty$ ) form of the inner-region expansion. Since the limiting form, at leading order, is just the far-field diffusion solution ${\bar {Q}_0}/({4\pi r})$ , or ${{\textit{Pe}} \bar {Q}_0}/({4\pi ^2k^2})$ in Fourier space, (3.18) may be rewritten in the form

(3.19) \begin{align} \tilde {\varTheta }(\tilde {\boldsymbol{x}}) &= \frac {\bar {Q}_0}{4\pi r} + {\textit{Pe}}\,\bar {Q}_0 \displaystyle \int e^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}}\left [ \frac {1}{2\pi \mathrm{i} k_U + 4\pi ^2 k^2} - \frac {1}{4\pi ^2 k^2} \right ]\text{d}\boldsymbol{k} \\[-12pt] \nonumber \end{align}

(3.20) \begin{align} &= \frac {\bar {Q}_0}{4\pi r} - {\textit{Pe}}\,\bar {Q}_0 \displaystyle \int e^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}}\frac {2\pi \mathrm{i} k_U}{(2\pi \mathrm{i} k_U + 4\pi ^2 k^2)4\pi ^2 k^2}\, \text{d}\boldsymbol{k}. \end{align}

One may now set $\tilde {\boldsymbol{x}} = 0$ in the ‘difference’ Fourier integral, so that

(3.21) \begin{align} \lim _{\tilde {\boldsymbol{x}} \rightarrow 0} \tilde {\varTheta }(\tilde {\boldsymbol{x}}) &= \frac {\bar {Q}_0}{4\pi r} - {\textit{Pe}} \frac {\mathrm{i}\bar {Q}_0}{2\pi } \displaystyle \int \frac {k_U}{(2\pi \mathrm{i} k_U + 4\pi ^2 k^2) k^2}\, \text{d}\boldsymbol{k} \\[-12pt] \nonumber \end{align}

(3.22) \begin{align} &= \frac {\bar {Q}_0}{4\pi r} - {\textit{Pe}}\frac {\bar {Q}_0}{8\pi }, \end{align}

on evaluating the difference integral. The leading term in (3.22) matches up to the far-field form of $\varTheta _0$ in the inner expansion; for a sphere, it reduces to $\varTheta _0$ exactly owing to the sphere-point-source equivalence. The $O({\textit{Pe}})$ term serves as a far-field boundary condition for $\varTheta _{1h}$ . Thus, for a prolate spheroid, for instance, one has to solve

(3.23) \begin{equation} {\nabla} ^2 \varTheta _{1h} = 0, \end{equation}

with $\varTheta _{1h} = 0$ at $\xi = \xi _0$ and $\varTheta _{1h} \rightarrow - {\bar {Q}_0(\xi _0)}/({8\pi })$ for $\xi \rightarrow \infty$ ; recall that $\bar {Q}_0 ={4\pi }/{(\xi _0\ln (({\xi _0 + 1})/({\xi _0-1}) ))}$ . The spherical symmetry of the far-field boundary condition implies that, to within an additive constant, $\varTheta _{1h}$ is identical to $\varTheta _0$ , being given by

(3.24) \begin{equation} \varTheta _{1h} = \frac {\bar {Q}_0(\xi _0)}{8\pi }\left [ \frac {\ln \left(\dfrac{\xi + 1}{\xi -1}\right)}{\ln \left(\dfrac{\xi _0 + 1}{\xi _0-1}\right)} -1 \right ]\!. \end{equation}

With the above-mentioned $O({\textit{Pe}})$ homogeneous solution, the inner-region scalar field is

(3.25) \begin{equation} \varTheta (\boldsymbol{x}) = \biggl [1+{\textit{Pe}} \frac {\bar {Q}_0(\xi _0)}{8\pi }\biggr ]\frac {\ln \left(\dfrac{\xi + 1}{\xi -1}\right)}{\ln \left(\dfrac{\xi _0 + 1}{\xi _0-1}\right)} + {\textit{Pe}}( - \frac {\bar {Q}_0}{8\pi } + \varTheta _{1p}), \end{equation}

where we have grouped the logarithmic term in $\varTheta _{1h}$ with $\varTheta _0$ . The coefficient of the combined term is proportional to the enhanced transport rate, given by $Q = \bar {Q}_0[1 + {\textit{Pe}}{\bar {Q}_0}/({8\pi })]$ . Dividing by $\bar {Q}_0$ leads to

(3.26) \begin{equation} {\textit{Nu}} = 1 + {\textit{Pe}} \frac {\bar {Q}_0}{8\pi } \end{equation}

for the small- ${\textit{Pe}}$ limiting form of the Nusselt number for uniform flow.

Since an ambient uniform flow enhances transport regardless of direction, ${\textit{Nu}}$ must be invariant to flow reversal ( $\hat {\boldsymbol{U}} \leftrightarrow - \hat {\boldsymbol{U}}$ ), as may be directly verified from implementing the change of variable ${\boldsymbol{k}} \leftrightarrow -{\boldsymbol{k}}$ in the integral in (3.21). Thus, ${\textit{Pe}}$ in (3.26) must be interpreted as $|{\textit{Pe}}|$ , leading to the ${\textit{Nu}}$ -curve not being differentiable at ${\textit{Pe}} = 0$ ! The usual argument of invariance to flow reversal would imply an $O({\textit{Pe}}^2)$ correction – clearly, this argument fails for a singular correction in an unbounded domain. Even seasoned physicists, at times, overlook this rather subtle point – in an otherwise insightful article on chemotactic reception of bacteria, Berg & Purcell (Reference Berg and Purcell1977) have a footnote arguing for the swimming-induced enhancement of transport (to a bacterium) being $O({\textit{Pe}}^2)$ ; Purcell (Reference Purcell1978) makes the same argument in the context of a linear flow. An analogous situation prevails in the momentum-transport scenario – the well-known Oseen correction to the Stokes drag is $3|{\textit{Re}}|/8$ , and not $3{\textit{Re}}/8$ (Benjamin Reference Benjamin1993). Note that ${\textit{Nu}}$ remains a differentiable function of ${\textit{Pe}}$ for any finite domain because, for ${\textit{Pe}}$ sufficiently small, the convective screening length ( $a{\textit{Pe}}^{-1})$ becomes larger than the domain size (the macro-scale $L_m$ ). The notion of an outer region no longer applies now and the correction comes from the inner-region, being $O({\textit{Pe}}^2)$ . In effect, the discontinuity in slope of the ${\textit{Nu}}$ -curve is smoothed out in a small but finite ${\textit{Pe}}$ -interval of $O(a/L_m)$ .

Determination of the convective correction to ${\textit{Nu}}$ in (3.26) only involved the limiting form of the outer-region scalar field. It is instructive to examine its complete structure since even weak convection accelerates the decay of the disturbance, while also profoundly altering its angular dependence at distances of $O(a{\textit{Pe}}^{-1})$ . This eventually leads to the disturbance being confined to a wake sufficiently far downstream, being exponentially small elsewhere. A signature of this fore–aft asymmetry may already be seen from neglecting the diffusive contribution (the term $4\pi ^2 k^2$ in the denominator in (3.18)), leading to $\tilde {\varTheta }(\tilde {\boldsymbol{x}}) \approx {\textit{Pe}}\bar {Q}_0\int ({e^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}}})/({2\pi \mathrm{i} k_U}) \,\text{d}\boldsymbol{k}$ . Using $\text{d}\boldsymbol{k} = \text{d}\boldsymbol{k}^\perp \text{d}k_U$ with $\boldsymbol{k}^\perp$ spanning the plane perpendicular to $\hat {\boldsymbol{U}}$ , $\tilde {\varTheta }(\tilde {\boldsymbol{x}}) \approx ({\bar {Q}_0}/{{\textit{Pe}}}) \delta (\boldsymbol{x}^\perp )\! \int \text{d}k_U {\sin (2\pi k_U\tilde {x}_U)}/({2\pi k_U}) = {\bar {Q}_0}/ ({2{\textit{Pe}}})\delta (\boldsymbol{x}^\perp )\text{sign}(x_U)$ , where $\text{sign}(x_U) = 1\,(-1)$ for $x_U = \boldsymbol{x} \boldsymbol{\cdot }\hat {\boldsymbol{U}} \gt 0\,(\lt 0)$ . Adding ${\bar {Q}_0}/ ({2{\textit{Pe}}})\delta ({\boldsymbol{x}}^\perp )$ to eliminate upstream influence, one obtains $\tilde {\varTheta }(\tilde {\boldsymbol{x}})=({\bar {Q}_0}/{\textit{Pe}})\delta (\boldsymbol{x}^\perp ){\mathcal H}(x_U)$ ( $\mathcal H$ being the Heaviside function), which represents an infinitely concentrated wake along the streamline downstream of the origin ( ${\boldsymbol{x}}^\perp = 0, x_U \gt 0$ ). The above addition may also be seen directly in Fourier space where, in the singular non-diffusive limit, $\hat {\varTheta }$ satisfies $2\pi \mathrm{i}k_U \hat {\varTheta } ={\textit{Pe}}\bar {Q}_0$ , whose general solution is of the form $\hat {\varTheta } = {\mathcal C}\delta (k_U) + Pe {\bar {Q}_0}/({2\pi \mathrm{i}k_U})$ . The constant multiplying the homogeneous solution, $\mathcal C$ , is determined from imposing the absence of an upstream influence, that is, from $\varTheta (\tilde {x}) = 0\,\forall \,x_U \lt 0$ (Nambiar, Garg & Subramanian Reference Nambiar, Garg and Subramanian2021).

The detailed structure of the outer-region scalar field, including the above-mentioned wake region, requires a full solution. Towards this end, it is convenient to include the unsteady term in (3.15) and start from the equation governing the impulse-response function:

(3.27) \begin{equation} \frac {\partial \tilde {\varTheta }^i_0}{\partial t} + \hat {\boldsymbol{U}} \boldsymbol{\cdot }\tilde {\boldsymbol{\nabla }} \tilde {\varTheta }^i_0 - \tilde {\boldsymbol{\nabla} }^2 \tilde {\varTheta }^i_0 = \bar {Q}_0\delta (\tilde {\boldsymbol{x}})\delta (t), \end{equation}

where the right-hand side also includes a delta function in time and corresponds to an impulsive forcing at the initial instant ( $t=0$ ). The solution of (3.15), corresponding to a steady forcing, is obtained via $\tilde {\varTheta }_0(\tilde {\boldsymbol{x}}) = \int _{0}^\infty \tilde {\varTheta }_0^i(\tilde {\boldsymbol{x}},s)\, \text{d}s$ . Equation (3.27) is identical to the Smoluchowski equation for the probability density of a Brownian particle starting from the origin at $t = 0$ , being subject to a uniform flow thereafter. To within a normalisation constant, $\tilde {\varTheta }^i_0$ may therefore be identified with this probability density. Without flow, the latter is a Gaussian with a variance that increases linearly in time. The ambient flow convects the centre of this expanding Gaussian envelope with unit speed along $\hat {\boldsymbol{U}}$ , so that

(3.28) \begin{equation} \tilde {\varTheta }^i_0(\tilde {\boldsymbol{x}},t) = \frac {\bar {Q}_0}{(4\pi t)^{{3}/{2}}} \exp \!\left [ - \frac {|\tilde {\boldsymbol{x}} - \hat {\boldsymbol{U}}t |^2}{4t} \right ]\!, \end{equation}

and therefore,

(3.29) \begin{align} \tilde {\varTheta }(\tilde {\boldsymbol{x}}) = \boldsymbol{Pe}\,\tilde {\varTheta }_0(\tilde {\boldsymbol{x}}) &= \boldsymbol{Pe}\,\bar {Q}_0\displaystyle \int _0^\infty \frac{\text{d}s}{(4\pi s)^{{3}/{2}}} \exp \!\left [ - \frac {|\tilde {\boldsymbol{x}} - \hat {\boldsymbol{U}}s |^2}{4s} \right ] \end{align}

(3.30) \begin{align} &= \frac {\bar {Q}_0}{4\pi r} \exp \left [ -\frac {\textit{Pe}}{2} ( r - x_U) \right ]\!, \\[10pt] \nonumber \end{align}

where we have reverted to inner-region coordinates in (3.30). This expression confirms the exponential smallness of the disturbance, on scales larger than $O(a{\textit{Pe}}^{-1})$ , almost everywhere in the domain, except where $r -x_U \ll x_U,\boldsymbol{x}^\perp$ , corresponding to the scalar wake. Since $r = (x_U^2 + {r^\perp }^2)^{{1}/{2}} \geqslant 0$ , the above-mentioned inequality is only satisfied downstream ( $x_U \gt 0$ ) when $r^\perp \ll x_U$ , implying $r -x_U \approx x_U [ (1+ {{r^\perp }^2}/{x_U^2})^{{1}/{2}} - 1] \approx {{r^\perp }^2}/({2x_U})$ . This leads to the approximate form $\tilde {\varTheta }_0^{w\textit{ake}}(\tilde {x}) = {\bar {Q}_0}/({4\pi x_U}) \exp [-\textit{Pe}({{r^\perp }^2}/({4x_U}))]$ within a paraboloidal wake growing in the downstream direction as $r^\perp \sim ({x_U}/{\textit{Pe}})^{{1}/{2}}$ , and within which the scalar continues to decay as the inverse of the distance (a behaviour that prevailed in all directions without flow); see figure 2. Note that $\lim \nolimits _{Pe \to \infty } ({\textit{Pe}}/({4\pi x_U})) \exp [-{\textit{Pe}}({{r^\perp }^2}/({4x_U})) ] = \delta ({\boldsymbol{x}}^\perp )$ with $x_U \gt 0$ , which recovers the generalised function expression for the infinitely localised wake given earlier. At large downstream distances, the entire transport occurs via the wake – $\int \text{d}S_{W\textit{ake}} \tilde {\varTheta }_0^{W\textit{ake}} = 2\pi \int _0^\infty r^\perp\, \text{d}r^\perp {\bar {Q}_0}/({4\pi x_U}) \exp [-{\textit{Pe}}({{r^\perp }^2}/({4x_U}))] = \bar {Q}_0$ ; the decay of $\tilde {\varTheta }_0^{W\textit{ake}}$ compensates for linear growth of the wake cross-section, the product of the two yielding $\bar {Q}_0$ .

Figure 2.

(a) Schematic of the paraboloidal wake behind a point source moving with speed $U$ ; it is the envelope of spheres centred at $(Ut,0,0)$ , with radii growing as $(Dt)^{{1}/{2}}$ , $\forall \,t \gt 0$ . (b) Convergence of the downstream scaled scalar profiles to the wake similarity form for ${\textit{Pe}}=0.1$ .

Figure 2. Long description

The image contains two graphs: a schematic diagram and a line graph. The schematic diagram on the left illustrates the paraboloidal wake behind a point source moving with speed U. It shows the envelope of spheres centered at different points, with radii growing as the square root of the product of the diffusion coefficient and time. The line graph on the right depicts the convergence of downstream scaled scalar profiles to the wake similarity form for various values of x sub U. The x-axis represents the scaled radial distance, while the y-axis represents the scaled scalar value. Different lines correspond to different values of x sub U, showing how the profiles converge to a similar form as x sub U increases. The schematic and the line graph together illustrate the behavior of scalar disturbances in the wake of a moving point source.

Using $t$ as the upper limit in (3.29) would reveal the manner in which the outer-region scalar field relaxes to its steady form, (3.30), after the ‘hot’ sphere is placed at the origin. The diffusive transport rate approaches its steady value in a time of $O(a^2/D)$ , with the unsteady contribution decaying as $O(t^{-{1}/{2}})$ for long times. The latter is the scalar analogue of the Basset–Boussinesq-type history contribution, which controls the approach to the steady Stokes drag in the momentum transport case. In contrast, the $O({\textit{Pe}})$ correction requires a steady state to be established in the outer region, which takes an asymptotically longer time of $O[(a{\textit{Pe}}^{-1})^2/D] \equiv O(a/U)$ – the outer region disturbance relaxes to (3.30) on a convective time scale, with the eventual approach towards steady state, for times much longer than $O(a/U)$ , being faster than $O(t^{-{1}/{2}})$ . This approach has been analysed in detail for momentum transport case (Mei, Lawrence & Adrian Reference Mei, Lawrence and Adrian1991; Lovalenti & Brady Reference Lovalenti and Brady1993), and more recently, for heat transport (Balachandar & Ha Reference Balachandar and Ha2021).

3.1.2. First convective correction: ambient linear flows

We now move to the more non-trivial case of an ambient linear flow, $\boldsymbol{u}^\infty = \boldsymbol{\varGamma } \boldsymbol{\cdot }\boldsymbol{x}$ , $\boldsymbol{\varGamma }$ being the transpose of the (non-dimensional) velocity gradient tensor. Use of the naive expansion $\varTheta = \varTheta _0 + Pe\,\varTheta _1 + \cdots$ again leads to $\varTheta _1$ satisfying (3.13). Now, the forcing function $\boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta _0 \approx (\boldsymbol{\varGamma } \boldsymbol{\cdot }\boldsymbol{x}) \boldsymbol{\cdot }\boldsymbol{\nabla } ({1}/{r}) \sim O(1/r)$ for $r \gg 1$ , implying $\varTheta _1 \sim O(r)$ , which leads to a greater (than uniform flow) conflict with the far-field decay requirement, and one therefore expects a stronger influence of convection. Indeed, comparing the convection and diffusion terms gives ${{\textit{Pe}}(\boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta )}/{{\nabla} ^2 \varTheta } \approx {{\textit{Pe}}(\boldsymbol{\varGamma } \boldsymbol{\cdot }\boldsymbol{x}) \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta _0}/{{\nabla} ^2 \varTheta _0} \sim ( {{\textit{Pe}}/r})/({1/r^3} )\sim {\textit{Pe}}\,r^2$ for $r \gg 1$ , so convective effects become important when $r \sim {\textit{Pe}}^{-({1}/{2})}$ , an asymptotically smaller distance compared with uniform flow (for the same ${\textit{Pe}}$ ), and that corresponds to a dimensional screening length of $(D/\dot {\gamma })^{{1}/{2}}$ . Defining the outer-region position vector as $\tilde {\boldsymbol{x}} = {\textit{Pe}}^{{1}/{2}} \boldsymbol{x}$ , one obtains

(3.31) \begin{equation} \left[\boldsymbol{\varGamma } \boldsymbol{\cdot }\tilde {\boldsymbol{x}} + O\big({\textit{Pe}}^{{3}/{2}}\big)\right] \boldsymbol{\cdot }\tilde {\boldsymbol{\nabla }} \tilde {\varTheta } - \tilde {\boldsymbol{\nabla} }^2 \tilde {\varTheta } = {\textit{Pe}}^{{1}/{2}} \left[\bar {Q}_0 + O\big({\textit{Pe}}^{{1}/{2}}\big)\right]\delta (\tilde {\boldsymbol{x}}), \end{equation}

in the outer region. The $O({\textit{Pe}}^{{3}/{2}})$ correction to the ambient flow, on the left-hand side, denotes the $O(1/r^2)$ stresslet velocity field induced by the neutrally buoyant particle (drop) at distances of $O({\textit{Pe}}^{-{1}/{2}})$ , while the $O({\textit{Pe}}^{{1}/{2}})$ correction to the source strength on the right-hand side anticipates the scaling of the correction to the diffusive transport rate. Unlike uniform flow, the fractional scaling makes the non-analyticity of the ${\textit{Nu}}$ -curve, at ${\textit{Pe}} = 0$ , obvious. Equation (3.31) suggests an outer-region expansion, $\tilde {\varTheta } = {\textit{Pe}}^{{1}/{2}}\tilde {\varTheta }_0 + Pe\,\tilde {\varTheta }_1 + \cdots$ , with

(3.32) \begin{equation} (\boldsymbol{\varGamma } \boldsymbol{\cdot }\tilde {\boldsymbol{x}}) \boldsymbol{\cdot }\tilde {\boldsymbol{\nabla }} \tilde {\varTheta }_0 - \tilde {\boldsymbol{\nabla} }^2 \tilde {\varTheta }_0 = \bar {Q}_0\delta (\tilde {\boldsymbol{x}}), \end{equation}

at leading order.

To analyse $\tilde {\varTheta }_0$ , we again consider the associated impulse-response field, $\tilde {\varTheta }^i_0$ , satisfying

(3.33) \begin{equation} \frac {\partial \tilde {\varTheta }_0^i}{\partial t} + (\boldsymbol{\varGamma } \boldsymbol{\cdot }\tilde {\boldsymbol{x}}) \boldsymbol{\cdot }\tilde {\boldsymbol{\nabla }} \tilde {\varTheta }_0^i - \tilde {\boldsymbol{\nabla} }^2 \tilde {\varTheta }_0^i = \bar {Q}_0 \delta (\tilde {\boldsymbol{x}})\delta (t). \end{equation}

The analogy with a Brownian particle in a linear flow again helps, since the Smoluchowski equation with a linear drift and constant diffusion coefficient, as in (3.33), leads to a Gaussian probability density, albeit an anisotropic one (Foister & Van de Wen Reference Foister and Van de Wen1980; Subramanian & Brady Reference Subramanian and Brady2004). One may write the formal solution of (3.33) as

(3.34) \begin{equation} \tilde {\varTheta }_0^i = \frac {\bar {Q}_0}{(2\pi )^{{3}/{2}}[\textit{det}\langle \tilde {\boldsymbol{x}}\tilde {\boldsymbol{x}}(t) \rangle ]^{{1}/{2}}} \exp \left [ -\frac {1}{2} \tilde {\boldsymbol{x}} \boldsymbol{\cdot }\langle \tilde {\boldsymbol{x}}\tilde {\boldsymbol{x}}\rangle ^{-1}(t) \boldsymbol{\cdot }\tilde {\boldsymbol{x}} \right ]\!, \end{equation}

where $\langle \tilde {\boldsymbol{x}} \tilde {\boldsymbol{x}} \rangle (t) = \int \tilde {\boldsymbol{x}}\tilde {\boldsymbol{x}}\,\tilde {\varTheta }_0^i\, \text{d}\tilde {\boldsymbol{x}}$ is the variance matrix. Inversion symmetry ( $\tilde {\boldsymbol{x}} \leftrightarrow -\tilde {\boldsymbol{x}})$ of the linear flow ensures that the Gaussian is always centred at the origin ( $\langle \tilde {\boldsymbol{x}}(t) \rangle = 0\,\forall \,t$ ). While $\langle \tilde {\boldsymbol{x}} \tilde {\boldsymbol{x}} \rangle = 2t{\boldsymbol{I}}$ in the quiescent (and uniform flow) cases, $\boldsymbol{I}$ being the identity matrix, it is anisotropic for a general linear flow, with the time dependence of the variance matrix elements depending sensitively on streamline geometry. This dependence reflects the differing nature of shear-enhanced Brownian motion: linear flows with open hyperbolic streamlines lead to an exponential-in-time spread of the probability density along certain directions, while those with closed streamlines lead to an oscillatory variation superposed on the linear-in-time (diffusive) growth already present in the absence of a shear (Foister & Van de Wen Reference Foister and Van de Wen1980). The marginal case of simple shear flow, with rectilinear streamlines, leads to variance elements growing algebraically with time. As will be seen, this range of behaviour has a direct impact on the structure of the scalar disturbance, including the existence of a wake and its geometry.

The above-mentioned aspects are best illustrated by considering a specific example – the one-parameter family of (canonical) planar linear flows. Before doing so, we introduce the parametrisation for a general incompressible linear flow that is convenient from a scalar transport viewpoint and will be used later for explicit calculations. We first note that $\boldsymbol{\varGamma }$ is traceless due to the incompressibility constraint, eliminating one of the original nine parameters. Further, $\boldsymbol{\varGamma }$ is the sum of the rate of strain ( $\boldsymbol{E} = ({1}/{2})({\boldsymbol \varGamma } + {\boldsymbol \varGamma }^\dagger )$ , symmetric) and vorticity ( $\boldsymbol{\varOmega } = ({1}/{2})({\boldsymbol \varGamma } - {\boldsymbol \varGamma }^\dagger )$ , antisymmetric) tensors, and choosing a coordinate frame aligned with the principal axes of $\boldsymbol{E}$ reduces the number of parameters to five. Finally, since the overall magnitude of the velocity gradient ( $\dot {\gamma }$ ) is already included in ${\textit{Pe}}$ and does not change streamline geometry, $\{ P \}$ in (2.2) and (2.3) comprises four flow-type parameters that define the geometry of an incompressible linear flow. A convenient choice for these parameters is $\epsilon =E_3/E_2$ measuring the non-axisymmetry of $\boldsymbol{E}$ , $E_i$ being the principal extensions ( $\sum _{i=1}^3 E_i = 0$ ); $\hat {\alpha }=|\boldsymbol{\omega }|/E_2$ measuring the relative magnitudes of vorticity and extension; and the polar ( $\theta _\omega$ ) and azimuthal ( $\phi _\omega$ ) angles that define the orientation of $\boldsymbol{\omega }$ in $\boldsymbol E$ -aligned coordinates (Krishnamurthy & Subramanian Reference Krishnamurthy and Subramanian2018a ; Narayanan & Subramanian Reference Narayanan and Subramanian2022). Thus, $\{ P \} \equiv \{ \epsilon ,\hat {\alpha },\theta _\omega ,\phi _\omega \}$ , with $\boldsymbol{\varGamma }$ given by

(3.35) \begin{align} \boldsymbol{\varGamma } = \begin{bmatrix} -(1+\epsilon ) & -\dfrac {\hat {\alpha } \cos \theta _{\omega }}{2} & -\dfrac {\hat {\alpha } \sin \theta _{\omega } \sin \phi _{\omega }}{2}\\ \dfrac {\hat {\alpha } \cos \theta _{\omega }}{2} & 1 & -\dfrac {\hat {\alpha } \sin \theta _{\omega } \cos \phi _{\omega }}{2}\\ \dfrac {\hat {\alpha } \sin \theta _{\omega } \sin \phi _{\omega }}{2} & \dfrac {\hat {\alpha } \sin \theta _{\omega } \cos \phi _{\omega }}{2} & \epsilon \end{bmatrix}. \end{align}

Accounting for axes relabelling and invariance to an overall sign change, all linear flow topologies may be shown to be covered for $\epsilon \!\in \![-2,0]$ , $\hat {\alpha }\!\geqslant \!0$ , $\theta _\omega ,\phi _\omega \!\in \![0,\!\pi /2]$ ; a schematic of an incompressible linear flow, with these ‘structural’ parameters, appears in figure 3(a). Figures 3(b) and 3(c) define two-parameter subsets of the 4-D space in figure 3(a). The former corresponds to non-axisymmetric extensional flows with aligned vorticity ( $\theta _\omega = 0$ ), so $\{ P\} \equiv \{\epsilon ,\hat {\alpha }\}$ ; the latter to axisymmetric extensional flows ( $\epsilon = -2, -0.5$ ) with inclined vorticity, so $\{ P\} \equiv \{\hat {\alpha }, \theta _\omega \}$ . Figure 3(d) depicts the one-parameter family ( $\hat {\alpha } = 0$ , $\{ P \} \equiv \{ \epsilon \}$ ) of extensional flows and figure 3(e) the one-parameter family ( $\epsilon = \theta _\omega =0, \{ P\} \equiv \{ \hat {\alpha } \}$ ) of planar linear flows. The parameters in figure 3(a) are an alternative to the four ratios that may be constructed from independent invariants associated with $\boldsymbol{E}$ and $\boldsymbol{\omega }$ ( $(\boldsymbol{E}:\boldsymbol{E})^{{1}/{2}}$ , $(\boldsymbol{\omega }\boldsymbol{\cdot }\boldsymbol{\omega })^{{1}/{2}}$ , $det(\boldsymbol{E})^{{1}/{3}}$ , $[(\boldsymbol{E}\boldsymbol{\cdot }\boldsymbol{E}):\boldsymbol{\omega }\boldsymbol{\omega }]^{({1}/{4})}$ and $(\boldsymbol{E}:\boldsymbol{\omega }\boldsymbol{\omega })^{{1}/{3}}$ ; see Poe & Acrivos (Reference Poe and Acrivos1975)).

Figure 3.

(a) Parameters defining an incompressible linear flow; (b) non-axisymmetric extensional flows with aligned vorticity; (c) axisymmetric extensional flows with inclined vorticity; (d) extensional flows; (e) (canonical) planar linear flows. Second row shows the streamline patterns for different planar flows; from left to right: solid-body rotation( $\alpha =-1$ , $\hat {\alpha } = \infty$ ), elliptic linear flow ( $\alpha =-0.5, \hat {\alpha } =6$ ), simple shear ( $\alpha = 0$ , $\hat {\alpha }=2$ ), hyperbolic linear flow ( $\alpha =0.5, \hat {\alpha }=2/3$ ) and planar extension ( $\alpha =1$ , $\hat {\alpha } = 0$ ).

Figure 3. Long description

The diagram illustrates various types of planar flows, including incompressible linear flow, non-axisymmetric extensional flows with aligned vorticity, axisymmetric extensional flows with inclined vorticity, extensional flows, and canonical planar linear flows. The second row shows the streamline patterns for different planar flows: solid-body rotation, elliptic linear flow, simple shear, hyperbolic linear flow, and planar extension.

While planar linear flows only constitute a one-dimensional (1-D) slice of the space in figure 3(a), they play an important role, being the local approximation of any incompressible 2-D velocity field. They are the simplest flows to illustrate the effect of a varying vorticity-to-extension ratio ( $\hat {\alpha }$ ) on different physical phenomena (scalar transport herein). They may be realised experimentally using a four-roll mill (Taylor Reference Taylor1934), variants of which have been used to study drop breakup and coalescence (Bentley & Leal Reference Bentley and Leal1986; Leal Reference Leal2001), and to generate the planar extensional flow that led to the ${\textit{Sh}}$ measurements in figure 1(a). In light of earlier efforts in the literature, we also use $\alpha ={(1-\hat {\alpha }/2)}/{(1+\hat {\alpha }/2)}$ to parametrise these flows. This choice and a $\pi /4$ -rotation of axes gives $\boldsymbol{\varGamma } = \Big[\begin{smallmatrix} 0 & 1 & 0\\ \alpha & 0 & 0\\ 0 & 0 & 0 \end{smallmatrix}\Big] $ ; $\alpha = 1$ , $0$ and $-1$ denote planar extension, simple shear and solid-body rotation, respectively, with the intervals $0 \lt \alpha \lt 1$ ( $\hat {\alpha } \lt 2$ ) and $-1 \lt \alpha \lt 0$ ( $\hat {\alpha } \gt 2$ ) corresponding respectively to hyperbolic and elliptic planar linear flows with open and closed streamlines; see second row of figure 3.

The absence of a third velocity component for planar linear flows implies that the scalar evolves diffusively along the third (vorticity) direction, this evolution being decoupled from that in the ( $x_1x_2$ -) plane of the flow. Equation (3.34) then reduces to the more specific form:

(3.36) \begin{align} \tilde {\varTheta }_0^i &= \frac {\bar {Q}_0}{(4\pi t)^{{1}/{2}}}\exp \left [- \frac {\tilde {x}_3^2}{4t} \right ] \frac {1}{2\pi \left(\left\langle \tilde {x}_1^2 \right\rangle \left\langle \tilde {x}_2^2 \right\rangle - \left\langle \tilde {x}_1\tilde {x}_2 \right\rangle ^2\right)^{{1}/{2}}}\nonumber\\&\quad \times\exp \left [ -\frac { \tilde {x}_1^2 \left\langle \tilde {x}_2^2 \right\rangle + \tilde {x}_2^2 \left\langle \tilde {x}_1^2 \right\rangle -2\tilde {x}_1\tilde {x}_2 \left\langle \tilde {x}_1\tilde {x}_2 \right\rangle }{2\left(\left\langle \tilde {x}_1^2 \right\rangle \left\langle \tilde {x}_2^2 \right\rangle - \left\langle \tilde {x}_1\tilde {x}_2 \right\rangle ^2\right)} \right ]\!, \end{align}

where the in-plane variances $\langle \tilde {x}_1^2 \rangle$ , $\langle \tilde {x}_2^2 \rangle$ and $\langle \tilde {x}_1\tilde {x}_2 \rangle$ are functions of both $\alpha$ and $t$ . The ordinary differential equations (ODEs) governing these variances may be obtained by operating both sides of (3.33) with $\int \tilde {\boldsymbol{x}}\tilde {\boldsymbol{x}}(.) \,\text{d}\tilde {\boldsymbol{x}}$ and using repeated integration by parts, while neglecting contributions from the surface at infinity, leading to

(3.37) \begin{equation} \frac{\text{d}}{\text{d}t}\left\langle \tilde {x}_i \tilde {x}_j \right\rangle = \varGamma _{il}\left\langle \tilde {x}_j \tilde {x}_l \right\rangle + \varGamma _{\textit{jl}} \left\langle \tilde {x}_i \tilde {x}_l \right\rangle + 2\delta _{\textit{ij}}, \end{equation}

in the general case. For planar linear flows, one only need consider the subset for $i,j = 1,2$ , solving which one obtains

(3.38) \begin{align} \langle \tilde {x}_1^2 \rangle (t) &= \frac {(1+\alpha )\left(\textit{e}^{2\alpha ^{{1}/{2}}t}-\textit{e}^{-2\alpha ^{{1}/{2}}t}\right)+4\alpha ^{{1}/{2}}(\alpha -1)t}{4\alpha ^{{3}/{2}}}, \end{align}

(3.39) \begin{align} \langle \tilde {x}_2^2 \rangle (t) &= \frac {(1+\alpha )\left(\textit{e}^{2\alpha ^{{1}/{2}}t}-\textit{e}^{-2\alpha ^{{1}/{2}}t}\right)-4\alpha ^{{1}/{2}}(\alpha -1)t}{4\alpha ^{{1}/{2}}}, \end{align}

(3.40) \begin{align} \langle \tilde {x}_1\tilde {x}_2 \rangle (t) &= \frac {(1+\alpha )}{4\alpha }\left[\textit{e}^{2\alpha ^{{1}/{2}}t}+\textit{e}^{-2\alpha ^{{1}/{2}}t}-2\right]\!. \end{align}

For hyperbolic flows, the variance along the extensional axis ( $\tilde {x}_1 =\alpha ^{{1}/{2}}\tilde {x}_2$ ) grows at an exponential rate for times longer than $O(\alpha ^{-{1}/{2}}\dot {\gamma }^{-1})$ , while that in the orthogonal direction ( $\tilde {x}_1 =-\alpha ^{{1}/{2}}\tilde {x}_2$ ) approaches a finite value due to a balance between the confining effect of flow and spreading action of diffusion. In contrast, for elliptic linear flows, one has complex exponentials and, thence, an oscillatory-cum-linear dependence on time; for solid-body rotation, diagonal variance components grow linearly in time as for the quiescent case. Evolution in simple shear is obtained as the limit of (3.38)–(3.40) for $\alpha \rightarrow 0$ , and yields $\langle \tilde {x}_1^2 \rangle = 2(t+{t^3}/{3})$ , $\langle \tilde {x}_2^2 \rangle = 2t$ and $\langle \tilde {x}_1\tilde {x}_2\rangle = t^2$ .

On using (3.36), the outer-region disturbance field satisfying (3.32) is given by

(3.41) \begin{align} \tilde {\varTheta }_0 &= \frac {\bar {Q}_0}{(2\pi )^{{3}/{2}}}\!\displaystyle \int _0^\infty \!\frac{\text{d}s}{\left[2s\left(\left\langle \tilde {x}_1^2\right\rangle (s) \left\langle \tilde {x}_2^2 \right\rangle (s) - \left\langle \tilde {x}_1\tilde {x}_2 \right\rangle ^2(s)\right)\right]^{{1}/{2}}}\exp \!\left [- \frac {\tilde {x}_3^2}{4s} \right ]\nonumber\\&\quad\times \exp \!\left [ -\frac { \tilde {x}_1^2 \left\langle \tilde {x}_2^2 \right\rangle (s) + \tilde {x}_2^2 \left\langle \tilde {x}_1^2 \right\rangle (s) -2\tilde {x}_1\tilde {x}_2 \left\langle \tilde {x}_1\tilde {x}_2 \right\rangle (s)}{2\left(\left\langle \tilde {x}_1^2 \right\rangle (s) \left\langle \tilde {x}_2^2 \right\rangle (s) - \left\langle \tilde {x}_1\tilde {x}_2 \right\rangle ^2(s)\right)} \right ]\!, \end{align}

with the variances as defined in (3.38)–(3.40). For obtaining the $\alpha$ -dependent transport enhancement, one proceeds along the lines of the previous uniform flow case. Thus, one adds and subtracts the diffusion solution ${\bar {Q}_0}/({4\pi r})$ in (3.41), with the subtracted solution written in the integral form $\bar {Q}_0\int _0^\infty ({\text{d}s}/{(4\pi s)^{{3}/{2}})}\textit{e}^{-{\tilde {r}^2}/({4s})}$ . Setting $\tilde {\boldsymbol{x}} = 0$ in the resulting difference integral leads to the following limiting form of the outer-region scalar field:

(3.42) \begin{align} \lim _{\tilde {x} \rightarrow 0} \tilde {\varTheta } = \frac {\bar {Q}_0}{4\pi r} - {\textit{Pe}}^{{1}/{2}} \frac {\bar {Q}_0}{4\pi ^{{3}/{2}}}\displaystyle \int _0^\infty \frac{\text{d}s}{s^{{1}/{2}}}\!\!\left [ \frac {1}{2s} - \frac {1}{\left(\left\langle \tilde {x}_1^2\right\rangle (s) \left\langle \tilde {x}_2^2 \right\rangle (s) - \left\langle \tilde {x}_1\tilde {x}_2 \right\rangle ^2(s)\right)^{{1}/{2}}}\right ]\!. \end{align}

The ambient shear enhances $\langle \tilde {x}_1^2\rangle \langle \tilde {x}_2^2 \rangle \!- \!\langle \tilde {x}_1\tilde {x}_2 \rangle ^2$ over its quiescent value ( $4s^2$ ), so the integral ( $\mathcal I$ , say) in (3.42) is positive. The term involving $\mathcal I$ serves as the far-field condition for the inner region, implying that the original expansion in powers of ${\textit{Pe}}$ was indeed naive. Unlike uniform flow, the inner-region expansion must proceed as $\varTheta = \varTheta _0 + {\textit{Pe}}^{{1}/{2}}\varTheta _{1h} + {\textit{Pe}}\,\varTheta _2 + \cdots$ . Again considering a prolate spheroid, ${\nabla} ^2 \varTheta _{1h} = 0$ with $\varTheta _{1h} = 0$ at $\xi = \xi _0$ and $\varTheta _{1h} \rightarrow -{\bar {Q}_0{\mathcal I}}/({4\pi ^{{3}/{2}}})$ for $\xi \rightarrow \infty$ . This gives $\varTheta _{1h} = ({\bar {Q}_0{\mathcal I}}/({4\pi ^{{3}/{2}}))} ({\ln [({\xi + 1})/({\xi -1})}]/{\ln [({\xi _0 + 1})/({\xi _0-1})]}-1)$ , whence the inner expansion takes the form

(3.43) \begin{equation} \varTheta (\boldsymbol{x}) = \left [1+ {\textit{Pe}}^{{1}/{2}}\frac {\mathcal I}{\pi ^{{1}/{2}}} \right ]\!\frac {\ln \left(\dfrac{\xi + 1}{\xi -1}\right)}{\ln \left(\dfrac{\xi _0 + 1}{\xi _0-1}\right)} - {\textit{Pe}}^{{1}/{2}} \frac {\bar {Q}_0{\mathcal I}}{4\pi ^{{3}/{2}}} + O({\textit{Pe}}), \end{equation}

with the enhanced transport identified as the coefficient of the logarithmic term. Thus,

(3.44) \begin{equation} {\textit{Nu}} = 1 + \frac {{\textit{Pe}}^{{1}/{2}}}{\pi ^{{1}/{2}}}\displaystyle \int _0^\infty \frac{\text{d}s}{s^{{1}/{2}}}\!\!\left [ \frac {1}{2s} - \frac {1}{\left (\left\langle \tilde {x}_1^2\right\rangle (s) \left\langle \tilde {x}_2^2 \right\rangle (s) - \langle \tilde {x}_1\tilde {x}_2 \rangle ^2(s)\right )^{{1}/{2}}}\right ] \end{equation}

for a particle (drop) in a planar linear flow.

The enhanced transport predicted by (3.44), and by (3.26) for uniform flow, arises from a depression of the scalar value (equalling $\lim \nolimits _{r \rightarrow \infty } {\textit{Pe}}\varTheta _{1h}$ and $\lim \nolimits _{r \rightarrow \infty } {\textit{Pe}}^{{1}/{2}}\varTheta _{1h}$ for uniform and linear flows, respectively) at the screening length. This depression may be interpreted as convection causing the disturbance to ‘vanish’ at the screening length itself, rather than at infinity. One might expect a linear flow calculation for a spheroid to be more difficult than the uniform flow one, since the flow causes an anisotropic particle to rotate; a spheroid in simple shear, for instance, rotates along Jeffery orbits for ${\textit{Re}} = 0$ (see figure 14(b) in § 4.2). However, the outer region, at leading order, only depends on the diffusive transport rate (the amplitude, $\bar {Q}_0$ , of the delta functions in (3.14) and (3.32)) which is evidently independent of the (changing) particle orientation. (The linear-flow calculation would indeed be more involved for a drop deforming in a time dependent manner. For small ${\textit{Ca}}$ , however, this corresponds to a very short initial transient.) Unsteadiness does come into play in the inner region, with the $O({\textit{Pe}})$ disturbance field being additionally forced by the unsteady term ( $({\partial \varTheta _0}/{\partial t})$ ). However, the $O({\textit{Pe}}^{{1}/{2}})$ outer-region contribution in (3.44) is evidently larger than any $O({\textit{Pe}})$ contribution arising from unsteady effects.

The $\alpha$ -dependent integral in (3.44), evaluated numerically, is plotted in figure 4, and shows the convective enhancement to be maximum for planar extension. For this flow, ${\textit{Nu}} = 1 + 0.43{\textit{Pe}}^{{1}/{2}}$ , which is one of the small- ${\textit{Pe}}$ correlations in figure 1(a) (the difference in coefficient, $0.36$ as opposed to $0.43$ , is due to the shear rate of Karp-Boss et al. (Reference Karp-Boss, Boss and Jumars1996) being defined as $({\boldsymbol{E}}\!:\!{\boldsymbol{E}})^{{1}/{2}} = \sqrt {2}\dot {\gamma }$ ). The enhancement decreases monotonically with increasing vorticity, the numerical coefficient equalling $0.257$ for simple shear ( $0.305$ on using $({\boldsymbol{E}}\!:\!{\boldsymbol{E}})^{{1}/{2}} = \dot {\gamma }/\sqrt {2}$ ), and eventually approaches zero as $(1+\alpha )^{{3}/{2}}$ for solid-body rotation. The kink at $\alpha = 0$ is indicative of the abrupt change in ambient flow topology.

Figure 4.

The small- ${\textit{Pe}}$ convective enhancement for a sphere in a planar linear flow, defined by $f(\alpha ) = {\textit{Pe}}^{-({1}/{2})}({\textit{Nu}}-1)$ , as a function of the flow-type parameter $\alpha$ . The bottom inset provides a magnified view of the kink at $\alpha = 0$ ; the top inset emphasises the ${\textit{Nu}}$ -scalings for $\alpha \rightarrow 0^\pm$ . Black dashed lines on the right denote the values for planar extension ( $\alpha = 1$ ) and simple shear ( $\alpha =0$ ).

Figure 4. Long description

The line graph illustrates the small-convective enhancement for a sphere in a planar linear flow, defined by a specific parameter, as a function of the flow-type parameter. The x-axis represents the flow-type parameter, while the y-axis represents the small-convective enhancement. The graph includes two insets: the bottom inset provides a magnified view of the kink at a specific point, and the top inset emphasizes the scaling behaviors for certain values. Black dashed lines on the right denote the values for planar extension and simple shear. The main graph shows a curve that increases and then levels off, with specific values highlighted for planar extension and simple shear. All values are approximated.

The outer-region origin of the convective correction implies its magnitude is sensitive to relatively distant bounding surfaces, potentially vitiating a comparison of theory with experiments (or computations). The predictions (3.26) and (3.44) are valid in the intervals $a/L_m \ll {\textit{Pe}} \ll 1$ and $(a/L_m)^2 \ll {\textit{Pe}} \ll 1$ , respectively; $L_m$ here being the effective size of the domain around the particle in experiments (or the periodic box dimensions in computations). As mentioned earlier, in the context of the analyticity of the ${\textit{Nu}}{-}{\textit{Pe}}$ relation, both predictions morph into a regular ${\textit{Pe}}^2$ -scaling for ${\textit{Pe}}$ smaller than the lower limits of the said intervals, due to the screening length exceeding $L_m$ . The latter seems to be the case in the experiments of Purcell (Reference Purcell1978) in figure 1(b), involving a sphere in planar extension, which show a quadratic (as opposed to a square-root) scaling for small ${\textit{Pe}}$ . The measurements were done in a four-roll mill, comprising symmetrically located rotating cylinders, with the effective domain around the central sphere (away from cylinder surfaces) estimated to be less than ten sphere radii.

The structure of the outer-region disturbance differs profoundly for positive and negative $\alpha$ , and for simple shear ( $\alpha = 0$ ). The disturbance field for the latter flow is

(3.45) \begin{align} \tilde {\varTheta }_0(\tilde {\boldsymbol{x}}) = \frac {\bar {Q}_0}{(4\pi )^{{3}/{2}}}\!\!\displaystyle \int _0^\infty \!\!\frac{\text{d}s}{s^{{3}/{2}}\left(1+ {s^2}/{12}\right)^{{1}/{2}}}\exp \!\left [\!- \frac {\tilde {x}_3^2}{4s} \right ]\! \exp \!\left [\! -\frac { \tilde {x}_1^2 + \tilde {x}_2^2\left(1+{s^2}/{3}\right) - \tilde {x}_1\tilde {x}_2 s}{4s(1+ {s^2}/{12})} \!\right ]\!. \end{align}

At distances much greater than $a{\textit{Pe}}^{-{1}/{2}}$ , corresponding to $\tilde {r} \gg 1$ , one must have $s \gg 1$ , so the integrand in (3.45) is not exponentially small, leading to the simpler form

(3.46) \begin{equation} \lim _{\tilde {r} \gg 1} \tilde {\varTheta }_0 = \frac {2\boldsymbol{\cdot }3^{{1}/{2}}\bar {Q}_0}{(4\pi )^{{3}/{2}}}\!\displaystyle \int _0^\infty \!\frac{\text{d}s}{s^{{5}/{2}}}\exp \!\left [- \frac {\tilde {x}_3^2}{4s} \right ]\! \exp \!\left [ -\left \{ \frac {3\tilde {x}_1^2}{s^3} + \frac {\tilde {x}_2^2}{s} - \frac {3\tilde {x}_1\tilde {x}_2}{s^2} \right \} \right ]\!. \end{equation}

Equation (3.46) may be examined in two distinguished limits. The first is where $\tilde {x}_1/\tilde {x}_2, \tilde {x}_1/\tilde {x}_3 \sim O(1)$ and corresponds to the region outside the wake. One now has $\tilde {x}_i \sim O(s^{{1}/{2}})$ and, therefore, $\tilde {x}_1^2/s^3, \tilde {x}_1\tilde {x}_2/s^2 \ll 1$ in the argument of the exponential. This leads to the simplified expression $\tilde {\varTheta }_0^{\textit{ outer}} \approx ({3^{{1}/{2}}}/({4\pi ^{{3}/{2}}}))\int _0^\infty ({\text{d}s}/{s^{{5}/{2}})}\exp [-(({\tilde {x}_3^2}/{4})+\tilde {x}_2^2) /{s}]$ , which gives $\tilde {\varTheta }_0^{\textit{ outer}} \approx ({3^{{1}/{2}}}/({8\pi }))(({\tilde {x}_3^2}/{4})+\tilde {x}_2^2)^{-{3}/{2}}$ . Of more interest is the wake region, where convection by simple shear and diffusion in the transverse directions are in balance. This requires $\tilde {x}_1^2/s^3 \sim \tilde {x}_2^2/s,\tilde {x}_3^2/s$ , implying $\tilde {x}_2,\tilde {x}_3 \sim |\tilde {x}_1|^{{1}/{3}}$ . The wake grows at a slower rate than in uniform flow – only as the cube root of distance along the flow direction. In the time $t_c \sim O[x_1/(\dot {\gamma }x_2)]$ that it moves through a distance of $O(x_1)$ along flow, the scalar diffuses through a distance $x_2,x_3 \sim O(D t_c)^{{1}/{2}}$ ; eliminating $t_c$ gives the dimensional rate of growth as $(Dx_1/\dot {\gamma })^{{1}/{3}}$ , consistent with the above-mentioned scaling. One may therefore define similarity variables $\eta = \tilde {x}_2/\tilde {x}_1^{{1}/{3}}$ , $\zeta = \tilde {x}_3/\tilde {x}_1^{{1}/{3}}$ , with the wake-disturbance field taking the form

(3.47) \begin{equation} \lim _{\tilde {r} \gg 1} \tilde {\varTheta }_0^{W\textit{ake}}(\eta ,\zeta ) = \frac {3^{{1}/{2}}}{(4\pi ^{{3}/{2}})|\tilde {x}_1|} \displaystyle \int _0^\infty \text{d}\hat {s}\frac {\textit{e}^{-\frac {3}{\hat {s}^3}}}{\hat {s}^{{5}/{2}}}\exp \!\left [-\left \{ \frac {\dfrac {\zeta ^2}{4}+\eta ^2}{\hat {s}} - \dfrac {3\eta }{\hat {s}^2} \right \} \right ]\!, \end{equation}

where $\hat {s}= s/\tilde {x}_1^{{2}/{3}}$ . The $|\tilde {x}_1|$ in (3.47) accounts for the wake along both positive and negative $x_1$ -axes, which arise due to the flow reversing across $x_2=0$ . The disturbance decays algebraically within and outside the wake, albeit at a faster rate outside ( $(\tilde {r}^\perp )^{-3}$ as opposed to $|\tilde {x}_1|^{-1}$ within the wake). Since the disturbance arises as the superposition of anisotropic Gaussian envelopes over an infinite time interval (see (3.45)), the algebraic decay can be understood by noting that the envelopes remain centred at the origin, allowing them to extend to infinity in all directions for sufficiently long times. The wake terminology for the region with the slower algebraic decay is nevertheless appropriate – the scalar is entirely transported through this region for sufficiently large $|x_1|$ ; $2\int \text{d}S_{wake} ({\boldsymbol{\varGamma }} \boldsymbol{\cdot }{\boldsymbol{x}})\,\tilde {\varTheta }_0 = 2 \int _{-\infty }^\infty \int _{-\infty }^\infty \text{d}\eta\, \text{d}\zeta \,\eta |\tilde {x}_1| \tilde {\varTheta }_0^{W\textit{ake}}= \bar {Q}_0$ , where the factor $2$ accounts for identical contributions from $\tilde {x}_1 \gt 0$ and ${\lt } 0$ . The structure of the outer field and convergence of the scaled large- $|\tilde {x}_1|$ scalar profiles to (3.47) are illustrated in figures 5(a) and 5(b).

Figure 5.

Outer-region disturbance fields in (a) simple shear and (c) planar extension; panels (b) and (d) illustrate the convergence of the appropriately scaled disturbance profiles, with increasing downstream distance, to the asymptotic self-similar forms given by (3.47) and (3.49).

Figure 5. Long description

The image contains four panels. Panels (a) and (c) show diagrams of disturbance fields in simple shear and planar extension, respectively. Panels (b) and (d) display line graphs illustrating the convergence of scaled disturbance profiles with increasing downstream distance to asymptotic self-similar forms. Each graph has an x-axis and y-axis with various labeled curves representing different downstream distances. The graphs in panels (b) and (d) show how the profiles converge to specific forms as the distance increases. All values are approximated.

For hyperbolic linear flows, the wake extends to infinity along the extensional axis ( $\tilde {x}_1 =\alpha ^{{1}/{2}}\tilde {x}_2$ ), while remaining finite in extent, orthogonal to it, in the plane of the flow; this being consistent with the behaviour of the variance components mentioned earlier. The structure of the far-field disturbance is most easily illustrated for planar extension, although the essential features remain the same for any $\alpha \gt 0$ at sufficiently large distances. For $\alpha = 1$ and $s \gg 1$ , (3.41) simplifies to

(3.48) \begin{equation} \lim _{\tilde {r} \gg 1} \tilde {\varTheta }_0 = \frac {\bar {Q}_0}{2^{{1}/{2}}(2\pi )^{{3}/{2}}} \exp \!\left [- \frac {\tilde {x}_2^2}{2} \right ]\!\displaystyle \int _0^\infty \!\frac{\text{d}s}{s^{{1}/{2}}\textit{e}^s}\exp \!\left [- \frac {\tilde {x}_3^2}{4s} \right ]\!\exp \!\left [ -\frac {\tilde {x}_1^2}{2\textit{e}^{2s}}\right ]\!, \end{equation}

where the exponential prefactor denotes the steady profile along $\tilde {x}_2$ ; unlike simple shear, the disturbance for large $\tilde {x}_2$ is exponentially small. Using $\tilde {x}_1/\tilde {x}_3 \sim O(1)$ , $\tilde {x}_3 \sim s^{{1}/{2}}$ leads to $\tilde {\varTheta }_0^{\textit{ outer}} \approx ({1}/({4\pi ^{{3}/{2}}}))\exp [- ({\tilde {x}_2^2}/{2}) ] \int _0^\infty({\text{d}s}/{s^{{1}/{2}}})\exp [- (s+ {\tilde {x}_3^2}/({4s})) ]= {1}/({4\pi })\exp [- ({\tilde {x}_2^2}/{2} +4|\tilde {x}_3|) ]$ . In the wake, one has $\tilde {x}_3 \sim s^{{1}/{2}}$ , $\tilde {x}_1^2 \sim \textit{e}^{2s}$ with $\tilde {x}_2 \sim O(1)$ , implying $\tilde {x}_3 \sim (\ln |\tilde {x}_1|)^{{1}/{2}}$ , a slow logarithmic spread in the vorticity direction. Using $s = \ln (\hat {s}|\tilde {x}_1|)$ , one obtains at leading order,

(3.49) \begin{align} \lim _{\tilde {r} \gg 1} \tilde {\varTheta }_0^{W\textit{ake}} = \frac {1}{4\pi ^{{3}/{2}}|\tilde {x}_1|(\ln |\tilde {x}_1|)^{{1}/{2}}} \exp \!\left [ - \frac {\tilde {x}_2^2}{2} \right ]\exp \!\left [ -\frac {\tilde {x}_3^2}{4\ln |\tilde {x}_1|}\right ] \displaystyle \int _0^\infty \frac{\text{d}\hat {s}}{\hat {s}^{2}}{e}^{- {1}/{(2\hat {s}^2)}}, \end{align}

with the $\hat {s}$ -integral evaluating to $({1}/{2})\varGamma ({1}/{4})$ . As for simple shear, the scalar is entirely transported through the wake: the $\tilde {x}_1^{-1}(\ln |\tilde {x}_1|)^{-{1}/{2}}$ -scaling of $\tilde {\varTheta }_0^{W\textit{ake}}$ compensating for the linear increase of the convecting flow and the logarithmic spread of the wake in the vorticity direction. The outer-region structure in planar extension is shown in figure 5(c) and the slow (logarithmic) convergence towards the self-similar form in figure 5(d).

For elliptic flows, as seen from (3.38)–(3.40), the diagonal variance elements grow linearly in time, with $\langle \tilde {x}_1^2 \rangle \approx {(1+|\alpha |)}/{|\alpha |}t$ , $\langle \tilde {x}_2^2 \rangle \approx (1+|\alpha |)t$ and $\langle \tilde {x}_1\tilde {x}_2 \rangle$ remaining bounded. Long-time evolution is therefore diffusive, albeit with an anisotropic diffusivity tensor, $\boldsymbol{D} = \Bigg[ \begin{smallmatrix} ({1+|\alpha |})/({2|\alpha |}) & 0 & 0 \\ 0 & ({1+|\alpha |})/{2} & 0 \\ 0 & 0& 1 \end{smallmatrix} \Bigg]$ , whose largest component is along the $x_1$ -direction (which, expectedly, diverges for $\alpha \rightarrow 0$ ). One has a diffusing ellipsoidal cloud in the far-field with $\lim \nolimits _{\tilde {r} \gg 1} \tilde {\varTheta }_0 = ( {\bar {Q}_0}/({4\pi })){(\boldsymbol{D}^{-1}:\tilde {\boldsymbol{x}}\tilde {\boldsymbol{x}})^{-{1}/{2}}}$ , there being no wake. For solid-body rotation, ${\boldsymbol{D}}^{-1}:\tilde {\boldsymbol{x}}\tilde {\boldsymbol{x}}= \tilde {r}^2$ , and the resulting limiting form, ${\bar {Q}_0}/({4\pi \tilde {r}})$ , is now exact and identical to pure diffusion.

To summarise, ${\textit{Nu}} = 1+ f({\{\mathcal P}\}){\textit{Pe}}^{{1}/{2}}$ for small ${\textit{Pe}}$ , for a particle (drop) in an ambient linear flow. For planar linear flows, $f(\{ {\mathcal P}\}) \equiv f(\alpha )$ and is defined by (3.44). Now, the multiplicative coupling of $\boldsymbol{u}$ and $\varTheta$ in (2.1) implies $f(\{ {\mathcal P} \})$ , for a given flow, cannot be calculated as a superposition of its values for simpler component flows that add up to the flow of interest; the implication being that, for a general linear flow, $f(\{ \mathcal P \})$ is not related in any simple way to the invariants of $\boldsymbol{E}$ and $\boldsymbol{\varOmega }$ . This fact is already evident from the radically different structures of the disturbance fields in planar extension, simple shear and solid-body rotation (simple shear results from superposing planar extension and solid-body rotation in equal proportions). Nevertheless, Batchelor (Reference Batchelor1979) found the first convective correction for an axisymmetric extension to have virtually the same numerical pre-factor as for planar extension, on using ${\textit{Pe}}$ based on the second invariant $({\boldsymbol{E}}\!:\!{\boldsymbol{E}})^{{1}/{2}}$ . Thus, ${\textit{Nu}} \approx 1 + 0.36 [({\boldsymbol{E}}:{\boldsymbol{E}})^{{1}/{2}}a^2/D]^{{1}/{2}}$ is expected to provide a good approximation for linear extensional flows (figure 3 d), with $f({\{\mathcal P}\}) \equiv f(\epsilon ) \approx 0.72(1+\epsilon + \epsilon ^2)$ for this family. For ambient linear flows with dominant vorticity ( $\hat {\alpha } \gt 2$ ), Batchelor (Reference Batchelor1979) proposed ${\textit{Nu}} \approx 1+0.4(a^2E_\omega /D)^{{1}/{2}}$ , provided $E_\omega = {\boldsymbol{E}}:\hat {\boldsymbol{\omega }}\hat {\boldsymbol{\omega }}$ , the projection of the extension rate along $\boldsymbol{\omega }$ ( $\hat {\boldsymbol{\omega }}=\boldsymbol{\omega }/|\boldsymbol{\omega }|$ being the associated unit vector), is larger than its magnitude in the orthogonal plane; the approximation fails for elliptic linear flows with $E_\omega = 0$ .

We end this subsection by noting that the small- ${\textit{Pe}}$ analysis may, in principle, be generalised to a nonlinear ambient flow, with $\boldsymbol{u}^\infty$ given by an $\boldsymbol{x}$ -polyad of an arbitrary order. The screening length is now $O(a{\textit{Pe}}^{-{1}/{n}})$ , with $n\geqslant 3$ for quadratic and higher order flows; ${\textit{Pe}}$ will itself be based on (the magnitude of) the analogue of $\boldsymbol{\varGamma }$ for such flows (for a quadratic flow, this would be a third-order tensor with units of $\text{length}^{-1}\,\text{time}^{-1}$ ) Boundaries, however, play an essential role in generating such nonlinear flows and must be included in a small- ${\textit{Pe}}$ transport rate calculation.

3.1.3. Higher order convective corrections

We now briefly analyse corrections beyond the leading order one discussed previously, which better highlight the singular nature of the small- ${\textit{Pe}}$ ${\textit{Nu}}$ -expansion, even for uniform flow, and also help increase its domain of validity. We start from Green’s second identity involving the scalar field of interest satisfying (2.1) and an auxiliary field ( $\varPhi$ ) satisfying Laplace’s equation with the same boundary conditions,

(3.50) \begin{equation} \displaystyle \int (\varPhi {\nabla} ^2 \varTheta - \varTheta {\nabla} ^2 \varPhi )\, \text{d}V = \displaystyle \int _{S_p} \{ \varTheta (\boldsymbol{n} \boldsymbol{\cdot }\boldsymbol{\nabla } \varPhi ) - \varPhi (\boldsymbol{n} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta ) \} \,\text{d}S. \end{equation}

The left-hand side integral is over the fluid domain (with ${\nabla} ^2 \varPhi = 0$ ) and $S_p$ on the right-hand side is the particle (drop) surface ( $\varTheta = \varPhi = 1$ on $S_p$ ). We neglect the contribution of the surface at infinity, based on the more rapid decay (of $\varTheta$ ) induced by convective effects in general; the analysis in earlier subsections shows this to be true except for ambient elliptic linear flows. Further, although the following arguments are based on the sphere solution, the generic structure of the ${\textit{Nu}}$ -expansion, for uniform flow in particular, will apply to other geometries too.

Using (2.1) and the surface boundary conditions, and identifying the first integral on the right-hand side of (3.50) with the diffusive transport rate, leads to

(3.51) \begin{align}{\textit{Nu}}&=1+ \frac {{\textit{Pe}}}{4\pi } \displaystyle \int \frac {1}{r}(\boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta ) \,\text{d}V, \end{align}

on using $\varPhi ={1}/{r}$ and $\bar {Q}_0 = 4\pi$ for a sphere. We now consider the scaling of the volume integral in (3.51) for ambient uniform and linear flows.

For uniform flow, on using the inner expansion to $O({\textit{Pe}})$ , the integral in (3.51) becomes ${\textit{Pe}}\int ({1}/{r})(\boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta _0)\,\text{d}V + {\textit{Pe}}^2\int ({1}/{r})(\boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta _1) \,\text{d}V$ . Writing $\boldsymbol{u} = \hat {\boldsymbol{U}} + \boldsymbol{u}'$ , $\boldsymbol{u}' \sim O(1/r)$ being the Stokes disturbance field at large distances, the $O({\textit{Pe}})$ integral takes the form ${\textit{Pe}}\int ({1}/{r})[\hat {\boldsymbol{U}} +O( {1}/{r})].O({1}/{r^2})r^2 \,\text{d}r \sim {\textit{Pe}}\int [O({1}/{r^3}) + O({1}/{r^4})]r^2\,\text{d}r$ – the ambient flow term appears to lead to a logarithmically divergent integral! Cutting off this divergence at $r \sim O({\textit{Pe}}^{-1})$ would imply a ${\textit{Nu}}$ contribution of $O({\textit{Pe}}\ln {\textit{Pe}})$ from the matching region ( $1 \ll r \ll {\textit{Pe}}^{-1}$ ), with logarithmically smaller $O({\textit{Pe}})$ contributions from inner and outer regions. However, $\varTheta _0$ being even in $\boldsymbol{x}$ , the matching and inner-region contributions vanish due to symmetry; the outer-region contribution survives, since $\tilde {\varTheta }_0$ only conforms to the basic requirement of axisymmetry about $\hat {\boldsymbol{U}}$ . This is consistent with the previous detailed considerations, where the leading correction to diffusive transport was an $O({\textit{Pe}})$ outer-region contribution – see (3.26). Next, considering the $O({\textit{Pe}}^2)$ integral, recall from (3.13) that the forcing term $\hat {\boldsymbol{U}} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta _0$ led to $\varTheta _{1p} \sim O(1)$ for $1 \ll r \ll {\textit{Pe}}^{-1}$ ; ${\boldsymbol{u}}' \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta _0$ leads to an additional $O(1/r)$ contribution in $\varTheta _{1p}$ . Thus, one may write the $O({\textit{Pe}}^2)$ integral in the form ${\textit{Pe}}^2\int O({1}/{r})[\hat {\boldsymbol{U}} +O({1}/{r})].[O({1}/{r})+O( {1}/{r^2})]r^2 \,\text{d}r$ – the dominant contribution is now linearly divergent. This points to an $O({\textit{Pe}})$ outer-region contribution, but one that has already been accounted for in solving for the leading order outer-region field. This may be understood by recognising that (3.18) for the outer-region scalar field captures the effect of the ambient flow alone to all orders in ${\textit{Pe}}$ ; a fact seen by expressing (3.18) in inner variables and expanding ${e}^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}} = {e}^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }{\textit{Pe}}{\boldsymbol{x}}}$ for small ${\textit{Pe}}$ . One only needed $\textit{e}^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}} \approx 1$ in (3.20) to capture the $O({\textit{Pe}})$ correction in (3.26). The next order term in $\tilde {\varTheta }$ requires the term linear in $\tilde {\boldsymbol{x}}$ in the previous expansion. Naïve expansion leads to divergent integrals, however; appropriate subtraction of divergent contributions leads to the following expression:

\begin{align*}&\lim _{\tilde {\boldsymbol{x}} \rightarrow 0} \int e^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}}\frac {1}{2\pi \mathrm{i} k_U + 4\pi ^2 k^2} \text{d}\boldsymbol{k} = \frac {1}{4\pi r} - {\textit{Pe}} \left[\int \frac {2\pi \mathrm{i} k_U}{\left(2\pi \mathrm{i} k_U + 4\pi ^2 k^2\right)4\pi ^2 k^2} \text{d}\boldsymbol{k}\right.\nonumber\\&\quad\left. - \int \big(e^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}} -1\big)\frac {2\pi \mathrm{i} k_U}{\big(4\pi ^2 k^2\big)^2}\text{d}\boldsymbol{k} \right] + {\textit{Pe}}\int \big(e^{2\pi \mathrm{i} \boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}}-1\big)2\pi \mathrm{i}k_U\nonumber\\&\quad\times \left[ \frac {1}{\big(4\pi ^2 k^2\big)^2} - \frac {1}{\big(4\pi ^2k^2\big)\left(2\pi \mathrm{i}k_U + 4\pi ^2 k^2\right)} \right]\!.\end{align*}

Table 2.

The small- ${\textit{Pe}}$ ${\textit{Nu}}$ -expansions for ambient uniform and linear flows.

Table 2. Long description

The equation presents the small expansions for ambient uniform and linear flows. It is divided into two main sections: one for ambient uniform flow and another for ambient linear flow. Each section is further divided into regimes based on the Peclet (Pe) and Reynolds (Re) numbers. For ambient uniform flow, the equation shows different contributions such as leading outer-region, second-order matching-region, and third outer-region contributions. For ambient linear flow, it includes leading outer-region contributions, outer-region contributions with source-strength corrections, and inner and outer-region contributions. The equation uses various mathematical terms like fractions, exponents, and functions to describe these contributions.

The first and second terms within brackets, at $O({\textit{Pe}})$ , are the $O(1)$ far-field forms of $\varTheta _{1h}$ and $\varTheta _{1p}\,= {\hat {\boldsymbol{U}} \boldsymbol{\cdot }{\boldsymbol{x}}}/({2r})$ . Using $e^{2\pi \mathrm{i}\boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}} - 1 = 2\pi \mathrm{i}\boldsymbol{k} \boldsymbol{\cdot }\tilde {\boldsymbol{x}}$ in the final integral leads to the $O(\tilde {\boldsymbol{x}})$ correction that must match up to a linearly diverging contribution, at $O({\textit{Pe}}^2)$ , in the inner expansion (the term $\varTheta _2$ ) – this term led to the linearly divergent integral in ${\textit{Nu}}$ .

The new contributions at the next order are thus logarithmically divergent integrals arising from coupling of the $O(1/r)$ Stokeslet field with the $O(1)$ term in $\varTheta _{1p}$ , and the ambient flow with the $O(1/r)$ term in $\varTheta _{1p}$ . Both integrands are even in $\boldsymbol{x}$ , so the divergence is real and there is an $O({\textit{Pe}}^2\ln {\textit{Pe}})$ ${\textit{Nu}}$ -contribution from the matching region. Using the far-field form, $\varTheta _{1p} \approx ({\hat {\boldsymbol{U}} \boldsymbol{\cdot }{\boldsymbol{x}}}/({2r}))(1- {3}/({2r}))$ , this contribution comes out to be $({1}/{2}){\textit{Pe}}^2 \ln {\textit{Pe}}$ ; smaller $O({\textit{Pe}}^2)$ contributions arise from the inner and outer regions, but require a more detailed analysis. The first entry in table 2 is the ${\textit{Nu}}$ expansion, to $O({\textit{Pe}}^3)$ , in the limit where inertia is sufficiently weak for the inertial screening length ( $a{\textit{Re}}^{-1}$ ) to be much larger than the convective screening length ( $a{\textit{Pe}}^{-1}$ ); it starts off with the aforementioned contributions, being accurate to within $10 \,\%$ for ${\textit{Pe}} \lesssim 0.3$ (Bell et al. Reference Bell, Byrne, Whiteley and Waters2013). The table also shows the ${\textit{Nu}}$ expansion in the complementary limit ${\textit{Pe}} \ll {\textit{Re}} \ll 1$ . The $O({\textit{Pe}})$ correction is unaffected since it arises from the ambient flow alone; the logarithmic divergence at the next order, described previously, is now cut off due to the more rapid decay of $\boldsymbol{u}'$ (rather than of $\varTheta$ ) at distances larger than $a{\textit{Re}}^{-1}$ , leading to a term of $O({\textit{Pe}}^2 \ln {\textit{Re}})$ .

For an ambient linear flow, using $\varTheta = \varTheta _0 + {\textit{Pe}}^{{1}/{2}}\varTheta _{1h} + {\textit{Pe}}\varTheta _2$ with $\boldsymbol{u} = \boldsymbol{\varGamma } \boldsymbol{\cdot }\boldsymbol{x} + \boldsymbol{u}'$ , $\boldsymbol{u}'$ being the $O(1/r^2)$ stresslet velocity field, one proceeds along lines similar to uniform flow. The integral in (3.51) involving the ambient flow is now linearly divergent and cutting off the divergence at $r \sim O({\textit{Pe}}^{-{1}/{2}})$ points to the $O({\textit{Pe}}^{{1}/{2}})$ outer-region contribution, given by (3.44), for planar linear flows. The first inner-region contribution arises only at $O({\textit{Pe}}^2)$ . While the absence of an $O({\textit{Pe}})$ correction again owes its origin to flow-reversal invariance, the invariance is not obvious since the original ( ${\boldsymbol{\varGamma }} \boldsymbol{\cdot }\boldsymbol{x}$ ) and time-reversed ( $-{\boldsymbol{\varGamma }} \boldsymbol{\cdot }\boldsymbol{x}$ ) linear flows are not in general the same; details of the invariance argument are given in § 4.1.1, and are, in fact, independent of ${\textit{Pe}}$ . Corrections at the intermediate orders, viz. $O({\textit{Pe}})$ and $({\textit{Pe}}^{{3}/{2}})$ , therefore, still arise from the outer region. Since the stresslet field in (3.31) is $O({\textit{Pe}}^{{3}/{2}})$ smaller and enters the outer-region expansion only at $O({\textit{Pe}}^2)$ , the above-mentioned corrections at lower orders arise simply from using the ${\textit{Pe}}$ -corrected estimates for the source term ( $\bar {Q}_0[1+f(\alpha ){\textit{Pe}}^{{1}/{2}}]$ and $\bar {Q}_0[1+f(\alpha ){\textit{Pe}}^{{1}/{2}}+f^2(\alpha ){\textit{Pe}}]$ ) on the right-hand side of (3.31), rather than $\bar {Q}_0$ ). The ${\textit{Nu}}$ expansions for an ambient linear flow, pertaining to both ${\textit{Re}} \ll {\textit{Pe}} \ll 1$ and ${\textit{Pe}} \ll {\textit{Re}} \ll 1$ , are also given in table 2.

The previous discussion summarises the nature of weak convective effects for a particle (drop) in a steady linear flow. It is natural to enquire about changes that arise when the ambient flow is unsteady, a scenario relevant to the examples in § 1: the role of bubble shape oscillations in enhancing transport rates and, thence, the performance of bubble column reactors, have been recognised (Martin et al. Reference Martin, Galan, Cerro and Montes2011); determining nutrient uptake rates by planktonic organisms requires analysis of transport by unsteady flow fields that are a natural consequence of ciliary kinematics (Magar & Pedley Reference Magar and Pedley2005; Gilpin, Prakash & Prakash Reference Gilpin, Prakash and Prakash2017). Unsteadiness may also arise as a result of adopting a particle-fixed reference frame for an anisotropic particle rotating in an otherwise steady shearing flow, a problem briefly considered in the next subsection. When the time scale characterising the unsteadiness is larger than $O(a^2/D)$ , the inner region may be regarded as quasi-steady and the ${\textit{Nu}}$ calculation proceeds along lines similar to the above-mentioned steady linear flows, but with $Q$ being time dependent. The Brownian motion analogy continues to hold. The outer-region disturbance field is still a Gaussian, but with components of the mean-squared displacement tensor now having to account for the ambient flow time dependence. For the case that the ambient linear flow varies periodically in time, matching the inner and outer limiting forms should, for instance, allow for the calculation of a time-averaged $Q$ , as a function of $\alpha$ and $Sl$ .

Figure 6.

(a) $\langle x_1^2 \rangle$ versus time for a Brownian particle (a diffusing scalar packet) in an oscillatory simple shear flow. The smallest $Sl$ plots transition from an initial linear growth to an intermediate cubic one, as in the steady case; a second transition to an eventual linear growth, across $t \sim O(Sl^{-1})$ , occurs due to oscillatory time dependence of the ambient flow. (b) $\langle x_1^2 \rangle$ versus time in an oscillatory hyperbolic linear flow, with $\alpha = 0.1$ , exhibiting a transition from an intermediate exponential growth to a linear one across $t \sim O(Sl^{-1})$ .

Figure 6. Long description

The image contains two line graphs labeled (a) and (b). Graph (a) shows the behavior of a Brownian particle in an oscillatory simple shear flow, with the y-axis representing a quantity and the x-axis representing the ratio of S to t. The graph illustrates transitions from initial linear growth to intermediate cubic growth and eventually to linear growth. Graph (b) depicts the behavior in an oscillatory hyperbolic linear flow, with a similar y-axis and x-axis. It shows a transition from intermediate exponential growth to linear growth. Both graphs highlight the complex dynamics of particle behavior under different flow conditions.

Armed with the knowledge that weak convection profoundly alters the structure of the scalar disturbance at distances of the order of the convective screening length, one may also enquire about the nature of the farfield transport in an unsteady linear flow; specifically, how are the schematics in figures 5(a) and 5(c) altered by ambient flow unsteadiness? Since the contrasting wake configurations in these figures arise from differing time dependencies of the mean squared displacement components, one needs to examine how these dependencies are altered by unsteadiness. We will only comment briefly on this aspect since it is not central to determining $Q$ . It is straightforward to determine the $\langle \tilde {x}_i\tilde {x}_j \rangle$ numerically, since the governing equations are still given by (3.37) with $\boldsymbol{\varGamma }$ now including a scalar time dependent amplitude. For a purely oscillatory linear flow, with $\boldsymbol{\varGamma } = {e}^{\mathrm{i}Sl t} \Big[ \begin{smallmatrix} 0 & 1 & 0\\ \alpha & 0 & 0\\ 0 & 0 & 0 \end{smallmatrix} \Big] $ , the variance elements grow at most linearly in time, for $Sl t \gg 1$ , regardless of $\alpha$ . Thus, while $\langle x_1^2 \rangle \sim \dot {\gamma }^2 D t^3$ for steady simple shear (in dimensional terms), the cubic growth in time only persists until $t \sim \omega ^{-1}$ for an oscillatory simple shear, with $\langle x_1^2 \rangle \sim \dot {\gamma }^2(D/\omega ^2)t$ for $t \gg \omega ^{-1}$ ; similar arguments apply to hyperbolic linear flows (see figures 6 a and 6 b). The eventual transition to a linear-in-time growth implies the existence of a larger secondary screening length beyond which the wake structures depicted in figure 5 are wiped out, the transport at larger distances being driven by an anisotropic effective diffusivity. This anisotropy is already imprinted in the above-mentinoed secondary screening length, which scales differently with $\omega$ along different directions – for simple shear, it is $O(\dot {\gamma }^2D/\omega ^3)^{{1}/{2}}$ and $O(D/\omega )^{{1}/{2}}$ along the flow and gradient directions; for a hyperbolic flow, one expects it to scale exponentially with $\omega ^{-1}$ along the extensional axis.

Of greater relevance is a time periodic ambient linear flow with a non-zero mean ( $\boldsymbol{\varGamma } = \boldsymbol{\varGamma }_\infty (\alpha ) + \boldsymbol{\varGamma }_o(\alpha ) \textit{e}^{\mathrm{i}Sl t}$ ) where, for long times, one intuitively expects the superposed oscillatory flow to manifest as an $Sl$ -dependent eddy diffusivity that adds to the bare Brownian one. While the idea of long-time diffusive behaviour, emerging from small-scale convection by an underlying time periodic or stochastic velocity field, has been dealt with in detail (Majda & Kramer Reference Majda and Kramer1999), using the method of multiple scales, the applicability of this paradigm to the case where the time periodic component is a linear flow (with no characteristic length scale) is not obvious.

4.1.1. Ambient linear flows with open streamlines; the Stokes limit

In § 3, the close relation between the small- ${\textit{Re}}$ and small- ${\textit{Pe}}$ scenarios led to an analogue of Whitehead’s paradox confounding a straightforward calculation of the first effects of convection on the transport rate. The matched asymptotic expansions approach in this regard was developed (Acrivos & Taylor Reference Acrivos and Taylor1962; Frankel & Acrivos Reference Frankel and Acrivos1968) along lines of the pioneering effort that resolved the paradox for small but finite ${\textit{Re}}$ (Proudman & Pearon Reference Proudman and Pearon1957). It should therefore come as no surprise that the singular nature of the large- ${\textit{Pe}}$ limit is similar to that of the large- ${\textit{Re}}$ limit for a streamlined body, where the momentum boundary layer remains attached even as ${\textit{Re}} \rightarrow \infty$ (a flat plate aligned with the free stream corresponding to the Blasius solution, or flow past wedge-shaped objects modelled by the Falkner–Skan family; see Lighthill (Reference Lighthill1950)). One faces the issue of a reduction of the order of the governing equation, (2.1), on a naive application of the infinite- ${\textit{Pe}}$ limit and, thence, the inability to satisfy all boundary conditions. The resolution is similar to that originally proposed by Prandtl in 1905 for the large- ${\textit{Re}}$ scenario – allowance for a thin boundary layer on the sphere surface where diffusion and convection proceed at comparable rates, and that approaches a jump discontinuity for ${\textit{Pe}} \rightarrow \infty$ .

We begin with the argument that points to the necessity of a boundary layer on the sphere surface, for ${\textit{Pe}} \gg 1$ , when one has an open-streamline topology; the conditions under which such a topology arises will be discussed in § 4.1.3. Naively setting ${\textit{Pe}} = \infty$ in (2.1), with the unsteady term neglected, leads to $\boldsymbol{u} \boldsymbol{\cdot }\boldsymbol{\nabla } \varTheta = 0$ , a first-order hyperbolic partial differential equation (PDE) with streamlines as the characteristic curves. The inference is that $\varTheta$ remains constant along any given streamline. Since open streamlines originate infinitely far upstream (where $\varTheta = 0$ ), one concludes $\varTheta = 0$ everywhere in the domain, leading to an incompatibility with the surface boundary condition. The error lies in assuming that $\boldsymbol{\nabla }, {\nabla} ^2$ in (2.1) remain bounded even as ${\textit{Pe}} \rightarrow \infty$ . This is not true – the gradients become steeper in a manner that the diffusion term in (2.1) remains comparable to the convective one within a thin boundary layer on the sphere surface. While the boundary layer thickness $\delta$ may be obtained formally by defining rescaled variables and implementing the correct leading order balance (Leal Reference Leal1992), it is also arrived at by noting that the convection time scale associated with the near-surface simple shear flow is $O(a/(\dot {\gamma }\delta ))$ . Here, we have  used that the near-surface shear rate is of the same order as the ambient one for small ${\textit{Re}}$ (when the sphere rotates, the near-surface shear is one transverse to the rotation-induced slip, as is the case for ambient vortical linear flows examined next). The time scale for diffusion across the boundary layer is $O(\delta ^2/D)$ . Equating the two time scales gives $\delta \sim a {\textit{Pe}}^{-{1}/{3}}$ – the scalar only diffuses a distance of $O(a {\textit{Pe}}^{-{1}/{3}})$ away from the sphere in the time of $O(\dot {\gamma }^{-1}{\textit{Pe}}^{{1}/{3}})$ in which near-surface shear sweeps it into the wake; while the inlet and wake locations are obvious for simple flows, they may in general be identified with invariant fixed points/limit cycles associated with the surface-streamline topology. The above-mentioned time scale of $O(\dot {\gamma }^{-1}{\textit{Pe}}^{{1}/{3}})$ is also that needed for the development of the boundary layer. If the sphere is suddenly immersed in the ambient flow at $t = 0$ , the transport for $t \ll \dot {\gamma }^{-1}{\textit{Pe}}^{{1}/{3}}$ is purely diffusive, across a boundary layer of thickness $O(Dt)^{{1}/{2}}$ , with the transport rate being $O(t^{-{1}/{2}})$ . Convective effects come into play for $t \sim O(\dot {\gamma }^{-1}{\textit{Pe}}^{{1}/{3}})$ , leading to a monotonic transition from the $t^{-{1}/{2}}$ -decrease to a plateau whose magnitude is determined by the structure of the steady boundary layer analysed later. Interestingly, the monotonic approach to steady state is only true when ${\textit{Re}}$ is small enough, or the geometry is such as to ensure the absence of attached recirculating eddies. Sharp corners, either due to an imposed geometry (Hasimoto & Sano Reference Hasimoto and Sano1980), or arising from a flow-transport coupling (Huang, Moore & Ristroph Reference Huang, Moore and Ristroph2015), lead to closed-streamline regions even in Stokes flows. In the presence of such eddies, the large- ${\textit{Pe}}$ dynamics of alignment of the iso-scalar contours with closed streamlines (within the eddies) leads to an oscillatory approach instead (Juncu Reference Juncu2010) and is discussed in the context of the drop-interior problem in § 5.2.

Figure 9.

(a) Scalar boundary-layer-wake configuration for a spherical particle in an ambient uniaxial extensional flow; the wake scalings are obtained in § 4.1.2. (b) Ensemble of theoretical predictions defining the large- ${\textit{Pe}}$ transport landscape.

Figure 9. Long description

The diagram consists of two parts. Part (a) illustrates the scalar boundary-layer-wake configuration for a spherical particle in an ambient uniaxial extensional flow. It shows the particle with labeled regions including the boundary layer and wake, along with arrows indicating flow directions and scalings. Part (b) presents an ensemble of theoretical predictions defining the large-transport landscape, with branches for spherical and non-spherical particles, and various flow conditions and equations.

For $t \gg \dot {\gamma }^{-1}{\textit{Pe}}^{{1}/{3}}$ , $\varTheta$ varies from $1$ to 0 across a steady $O(a{\textit{Pe}}^{-{1}/{3}})$ boundary layer attached to the sphere, only leaving it in an asymptotically narrow wake, as shown in figure 9(a). Here, $\varTheta$ is exponentially small outside the boundary layer and wake regions. Owing to the smaller gradients in the wake, the leading order transport occurs across the attached boundary layer. The flux, being proportional to the scalar gradient across the boundary layer, must scale inversely as its thickness and ${\textit{Nu}}$ must therefore scale as ${\textit{Pe}}^{{1}/{3}}$ . Thus, for ${\textit{Re}} \ll 1$ and with $\{ S \} \equiv \varnothing$ , (2.2) reduces to ${\textit{Nu}} = \mathcal{F}_{p\infty }(\{ P \}){\textit{Pe}}^{{1}/{3}}$ for ${\textit{Pe}} \gg 1$ , for a spherical particle in an open-streamline configuration. In the exceptional cases where boundary layer separation does not occur, one may infer the asymptotic ${\textit{Nu}}$ -scaling form even when ${\textit{Re}}$ is large. For ${\textit{Re}}, {\textit{Pe}} \gg 1$ , there is both a momentum and a scalar boundary layer on the surface of the body. Assuming laminar flow, the former has a thickness of $O(a{\textit{Re}}^{-{1}/{2}})$ with an associated shear rate $\dot {\gamma } \sim (U/a){\textit{Re}}^{{1}/{2}}$ , $a$ being an appropriate geometric scale; for the aligned flat plate, $a$ is the plate length. For $P\!r \gg 1$ , when the scalar boundary layer is much thinner than the momentum boundary layer, convection is driven by the near-surface simple shear, but with the ${\textit{Re}}$ -dependent $\dot {\gamma }$  of the momentum boundary layer; the arguments in the previous paragraph give ${\textit{Nu}} \propto {\textit{Re}}^{{1}/{6}}{\textit{Pe}}^{{1}/{3}} \equiv {\textit{Re}}^{{1}/{2}}P\!r^{{1}/{3}}$ . For $P\!r \ll 1$ , convection within the much thicker scalar boundary layer is largely driven by the uniform flow outside the momentum boundary layer, this being similar to the large- ${\textit{Pe}}$ transport from a drop analysed in § 5; one then obtains ${\textit{Nu}} \propto {\textit{Pe}}^{{1}/{2}}$ . At this point, it is pertinent to go back to the correlations in the second column of the table in figure 1(a) – the first and second entries in the third and second rows, respectively, correspond to the above-mentioned small- and large- ${\textit{Re}}$ (large- $P\!r$ ) correlations. The first entry in the third row linearly combines the small- ${\textit{Pe}}$ and large- ${\textit{Pe}}$ –large- ${\textit{Re}}$ ${\textit{Nu}}$ -asymptotes; for $A= 0.6$ , this empirical combination is termed the RanzMarshall correlation in the combustion literature.

Consideration of all possible linear flow topologies, to determine $\mathcal{F}_{p\infty }(\{ P \})$ in the above-mentioned ${\textit{Nu}}$ functional, might seem daunting in light of the 4-D parameter space in figure 3(a). However, the situation simplifies for a spherical particle, with one only needing to discriminate between two classes of ambient linear flows for ${\textit{Pe}} \rightarrow \infty$ – those with and without vorticity. This is seen from the velocity field induced by a torque-free sphere in a linear flow:

(4.1) \begin{align} \boldsymbol{u}(\boldsymbol{x}) &=\,\boldsymbol{\varOmega } \boldsymbol{\cdot }\boldsymbol{x} + \boldsymbol{E} \boldsymbol{\cdot }\boldsymbol{x} \left (1-\frac {1}{r^5}\right ) - \frac {5}{2}\boldsymbol{x}(\boldsymbol{E}:\boldsymbol{x}\boldsymbol{x})\left (\frac {1}{r^5}-\frac {1}{r^7}\right ). \end{align}

As a result, $\lim\nolimits _{r \rightarrow 1}\boldsymbol{u} = \boldsymbol{\varOmega }\boldsymbol{\cdot }\boldsymbol{n} + 5(\boldsymbol{E}\boldsymbol{\cdot }\boldsymbol{n})\boldsymbol{\cdot }(\boldsymbol{I}-\boldsymbol{n}\boldsymbol{n})(r-1) + ({15}/{2})\boldsymbol{n}(\boldsymbol{E}:\boldsymbol{n}\boldsymbol{n})(r-1)^2$ , with the first term corresponding to sphere rotation with an angular velocity equalling half the ambient vorticity. Using this limiting form in (2.1), and retaining only the dominant radial diffusion term on the right-hand side, leads to

(4.2) \begin{equation} \frac {15}{2}(\boldsymbol{E}:\boldsymbol{n}\boldsymbol{n})y^2\frac {\partial \varTheta }{\partial y} + (\boldsymbol{\varOmega }\boldsymbol{\cdot }\boldsymbol{n} + 5y\boldsymbol{E}\boldsymbol{\cdot }\boldsymbol{n})\boldsymbol{\cdot }\boldsymbol{\nabla }_s \varTheta = \frac {1}{\textit{Pe}}\frac {\partial ^2 \varTheta }{\partial y^2}, \end{equation}

governing the leading order scalar distribution within the boundary layer. Here, $y = r-1$ , so $y=0$ is the sphere surface and $\boldsymbol{\nabla }_s = (\boldsymbol{I}-\boldsymbol{n}\boldsymbol{n}) \boldsymbol{\cdot }\boldsymbol{\nabla }$ denotes the tangential component of the gradient operator. The near-surface streamline topology depends on the relative magnitudes of the tangential velocity contributions involving $\boldsymbol{E}$ and $\boldsymbol{\varOmega }$ on the left-hand side of (4.2). Since $y \sim \delta /a \sim O({\textit{Pe}}^{-{1}/{3}})$ and $|{\boldsymbol{\varOmega }}|/|{\boldsymbol{E}}| \equiv \hat {\alpha }$ , a comparison identifies two limiting scenarios: $\hat {\alpha } \ll {\textit{Pe}}^{-{1}/{3}}$ and $\hat {\alpha } \gg {\textit{Pe}}^{-{1}/{3}}$ .

For ${\textit{Pe}} \rightarrow \infty$ , the first scenario reduces to linear extensional flows (figure 3 d), so $\boldsymbol{\varGamma } = \boldsymbol{E} = \Bigg[\begin{smallmatrix} -(1+\epsilon ) & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & \epsilon \end{smallmatrix} \Bigg]$ from (3.35). Since a sphere does not rotate in these flows, the surface velocity is identically zero. The direction of the near-surface shear within the boundary layer is defined by surface stress-lines, curves tangential to the shear stress field on the sphere surface. The principal axes of $\boldsymbol{E}$ mark out three pairs of antipodal fixed points on the unit sphere – stable and unstable nodes corresponding to the largest and smallest principal rates of strain, and saddles corresponding to the intermediate one – which divide the surface into eight identical sectors contributing equally to the transport rate. As shown in figure 10(a–c), each sector is bounded by arcs connecting three fixed points (denoted by dashed lines). Barring the bounding arcs that start from, or end up at, the saddles, the stress-lines run from an unstable node to the neighbouring stable ones. For axisymmetric extension ( $\epsilon = -2$ or $-{1}/{2}$ ) alone, the great circle comprising stress-lines that connect the stable (for $\epsilon = -2$ ) or unstable (for $\epsilon = - {1}/{2}$ ) nodes with the saddles, degenerates to a fixed-point ring. While (4.2) has been analysed in Von Mises coordinates (Batchelor Reference Batchelor1979), it is better solved in a non-orthogonal coordinate system (Krishnamurthy & Subramanian Reference Krishnamurthy and Subramanian2018 a) that generalises to the more complicated drop surface-streamline topologies examined in § 5.1. Here, $C$ and $\tau$ are used to denote these latter coordinates, since a special instance of the above-mentioned system was first used to analyse the orientation distribution of a Brownian spheroid in simple shear, in which case the coordinates correspond to the Jeffery orbit constant ( $C$ ) and phase ( $\tau$ ) (Leal & Hinch Reference Leal and Hinch1971); see §§ 4.2 and 5.1.2. In spherical coordinates with the polar axis along $E_3$ , the tangential components of the near-surface shear flow, given in the paragraph following (4.1), are:

(4.3) \begin{align} u_\theta = -\frac {5}{4}y \sin 2\theta [(2+\epsilon )\cos 2\phi + 3\epsilon ];& u_\phi = 5y(2+\epsilon ) \sin \theta \sin \phi \cos \phi , \end{align}

Figure 10.

Surface stress-lines on a spherical-particle octant in a linear extensional flow: (a) $\epsilon =-1$ ; (b) $\epsilon =-1.5$ ; (c) $\epsilon =-1.95$ . Panel (b) clarifies the $C\!-\!\tau$ labels of the stress-lines; panel (c) shows the approach towards a fixed ring in the $x_1\!-\!x_2$ plane for $\epsilon \rightarrow -2$ ; panel (d) shows ${\textit{Nu}}/{\textit{Pe}}^{{1}/{3}}$ for a spherical particle, and ${\textit{Nu}}/\tilde {{\textit{Pe}}}^{{1}/{2}}$ for a spherical drop ( $\tilde {{\textit{Pe}}}={\textit{Pe}}/(1+\lambda )$ ), as a function of $\epsilon$ , as predicted from a boundary layer analysis; the drop calculation is in § 5.1.1; panel (e) compares the electrochemically measured transport rate, from a fixed sphere in simple shear, to the prediction (solid black line) for an ambient axisymmetric extension; ${\textit{Pe}}^* = ({\boldsymbol{E}}\!:\!{\boldsymbol{E}})^{{1}/{2}}(2a)^2/D$ . The inset in panel (e) shows the non-canonical stress-line topology in the experiments (contrast with that in panel a–c); reproduced with permission from Kutateladze, Nakoryakov & Iskakov (Reference Kutateladze, Nakoryakov and Iskakov1982).

Figure 10. Long description

The image contains multiple graphs and diagrams related to transport-rate correlations and surface stress-lines on a spherical particle in various flows. Panel (a) shows a table of transport-rate correlations from Perry’s Chemical Engineers Handbook. Panel (b) presents a plot of the Sherwood number versus Péclet number, demonstrating convectively enhanced scalar transport for a sphere in an ambient shearing flow. Equations (32) and (34) are small correlations, while equations (35) and (36) are large correlations. Filled symbols denote computations for a sphere in a uniaxial extensional flow, and open circles denote experiments for heat transport from a sphere in a planar extension. Panels (c), (d), and (e) show surface stress-lines on a spherical-particle octant in a linear extensional flow, the Sherwood number for a spherical particle and a spherical drop as a function of a parameter, and a comparison of electrochemically measured transport rate to the prediction for an ambient axisymmetric extension, respectively. The inset in panel (e) shows the non-canonical stress-line topology in the experiments.

where the coefficients of $y$ denote the angle-dependent surface shear stresses ( $\tau _{\theta y}$ and $\tau _{\phi y}$ ). Regarding $s = yt$ as a variable that measures time in units of $(y\dot {\gamma })^{-1}$ (so the $s$ -duration remains finite, even as the actual time taken to traverse any finite distance, along the sphere surface, diverges for $y\rightarrow 0$ ), leads to the equations

(4.4) \begin{align} \frac{\text{d}\theta }{\text{d}s} =-\frac {5}{4}\sin 2\theta [(2+\epsilon )\cos 2\phi + 3\epsilon ];\,\frac{\text{d}\phi }{\text{d}s}=5(2+\epsilon ) \sin \phi \cos \phi , \end{align}

governing the surface stress-lines. Integrating (4.4) yields the definitions of $C$ and $\tau$ :

(4.5) \begin{align} \tan \phi =\, \textit{e}^\tau , \tan \theta = C\big(1+\textit{e}^{2\tau }\big)^{{1}/{2}}{e}^{-( {(1+2\epsilon )\tau })/{(2+\epsilon )}} \end{align}

for linear extensional flows; $\tau = 5(2+\epsilon )s$ . The stress-line labels in a single octant ( $0 \leqslant \theta ,\phi \leqslant {\pi }/{2}$ ) are depicted in figure 10(b) for $\epsilon = -1.5$ . Here, $C = 0$ corresponds to the arc $\phi = {\pi }/{2}$ , while $C=\infty$ is the compound arc comprising $\phi = 0$ and $\theta = {\pi }/{2}$ linked at the saddle point $(\theta ,\phi )\equiv ({\pi }/{2},0)$ . The stress-lines in the octant interior correspond to $0 \lt C \lt \infty$ . Additionally, $\tau$ ranges from $-\infty$ (unstable node at $\theta = 0$ ) to $\infty$ (stable node at $\theta = \phi = {\pi }/{2}$ ) along each of these, its evolution becoming increasingly skewed for large $C$ – in this limit, $\tau$ remains large and negative until $\theta$ increases to ${\pi }/{2}$ , with $\phi$ still being small; it increases thereafter with increasing $\phi$ , $\theta$ remaining close to ${\pi }/{2}$ . Using (4.5), (4.2) takes the form

(4.6) \begin{equation} \frac {15}{2}\left ( \frac {C^2\left[\left(\textit{e}^{2\tau }-1\right)-\epsilon \right]+\textit{e}^{({2(1+2\epsilon )\tau }/{2+\epsilon })}\epsilon }{\left[\textit{e}^{({2(1+2\epsilon )\tau }/{2+\epsilon })}+C^2\left(1+\textit{e}^{2\tau }\right)\right]} \right )y^2 \frac{\partial \varTheta }{\partial y} + 5(2+\epsilon )y \frac{\partial \varTheta }{\partial \tau } = \frac {1}{{\textit{Pe}}} \frac {\partial ^2 \varTheta }{\partial y^2}, \end{equation}

The tangential convection only has a $\tau$ -component and that is $\tau$ -independent. Use of $C{-}\tau$ coordinates achieves, in some sense, the simplification one would normally expect for an axisymmetric ambient flow in spherical polar coordinates, with the polar axis aligned with the axis of symmetry.

Hereon, the analysis involving the derivation of a self-similar solution of (4.6) proceeds along standard lines (Chapter 9 of Leal Reference Leal1992). In terms of a re-scaled coordinate $Y = {\textit{Pe}}^{{1}/{3}}y$ , which takes $O(1)$ values within the boundary layer, the boundary conditions are $\varTheta = 1$ at $Y = 0$ and $\varTheta \rightarrow 0$ for $Y \rightarrow \infty$ , the latter denoting the exponential smallness of the disturbance outside the boundary layer. Next, defining a similarity variable, $\eta = Y/g(C,\tau )$ , with $g(C,\tau )$ denoting the boundary layer thickness field, one arrives at the following expression for the scalar field:

(4.7) \begin{equation} \varTheta (\eta ) = \frac {1}{\varGamma ({4}/{3})}\displaystyle \int _\eta ^\infty \textit{e}^{-q^3}\,\text{d}q. \end{equation}

Here, g is obtained from solving $({1}/{3})h_\tau ({\text{d}(g^3)}/{\text{d}\tau )} - h_r g^3 = 3$ , along with the constraint of finiteness at the inlet which, from figure 10(b), corresponds to $\tau \rightarrow -\infty$ for all finite $C$ ; $h_C(C,\tau ;\epsilon )$ and $h_\tau (\epsilon )$ in the above-mentioned ODE are the coefficients of $y^2$ and $y$ in the convective terms in (4.6). The boundary layer thickness comes out to be $g =[{9{\mathcal Q}(C,\tau ;\epsilon )}/({5(2+\epsilon ))}\int _{-\infty }^\tau Q(C,\tau ;\epsilon )^{-1}\text{d}\tau ']^{{1}/{3}}$ with $\mathcal Q(C,\tau _1;\epsilon ) = \textit{e}^{-{9(1+\epsilon )}/({2(2+\epsilon )})\tau }[\textit{e}^{({2(1+2\epsilon )}/({2+\epsilon )})\tau }+ (1+\textit{e}^{2\tau })C^2]^{{9}/{4}}$ . Now, ${\textit{Nu}}$ may be calculated as

(4.8) \begin{align} {\textit{Nu}} &= {\textit{Pe}}^{{1}/{3}}\frac {8}{4\pi \varGamma ({4}/{3})}\displaystyle \int _0^\infty \displaystyle \int _{-\infty }^\infty (h\,\text{d}C) (k\,\text{d}\tau )\sin \theta _{C\tau }\frac {1}{g} \\[-12pt]\nonumber \end{align}

(4.9) \begin{align} &= {\textit{Pe}}^{{1}/{3}}\frac {[15(2+\epsilon )]^{{1}/{3}}}{\pi \varGamma ({4}/{3})}\displaystyle \int _0^\infty C\left [\displaystyle \int _{-\infty }^{\infty } \frac{\text{d}\tau }{Q(C,\tau ;\epsilon )} \right ]^{{2}/{3}}\,\text{d}C. \end{align}

Consistent with earlier arguments that led to ${\textit{Nu}} \propto {\textit{Pe}}^{{1}/{3}}$ , ${\textit{Nu}}$ in (4.8) is proportional to the areal integral of the inverse boundary layer thickness. The areal element, $(h\,\text{d}C)(k\,\text{d}\tau )\sin \theta _{C\tau }$ , denotes an infinitesimal parallelogram on the sphere surface, with $h$ and $k$ the (diagonal) metric factors associated with $C$ and $\tau$ , and $\theta _{C\tau }\,(\neq {\pi }/{2}$ ) the angle between constant- $C$ and constant- $\tau$ lines. Although $\epsilon \in [-2,0]$ for all linear extensional flows, one only need evaluate (4.9) for $\epsilon \in (-2,1]$ . (Biaxial axisymmetric extension ( $\epsilon = -2$ ) is a singular limit as evident from the $C\text{-}$ relation in (4.5). One calculates ${\textit{Nu}}$ for this flow as a limiting process, rather than setting $\epsilon = - 2$ .) This is because the flows for $\epsilon \,\in \,[-1,0]$ are time-reversed versions of those in $[-2,-1]$ , and remarkably, ${\textit{Nu}}$ is invariant to ambient flow reversal, despite the reversal leading to a qualitative change in the structure of the disturbance field, especially for large ${\textit{Pe}}$ (Brenner Reference Brenner1967). This is readily proven starting from (3.50), where we now take $\varTheta$ and $\varPhi$ to be the scalar fields for the original ( $\boldsymbol{u}$ ) and time reversed ( $-{\boldsymbol{u}}$ ) flows, satisfying the same boundary conditions. The surface boundary condition implies that the right-hand side in (3.50) is the difference between the original and flow-reversed Nusselt numbers, leading to

(4.10) \begin{align} {\textit{Nu}}_\varTheta - {\textit{Nu}}_\varPhi &= \displaystyle \int \varPhi ({\nabla} ^2\varTheta - \varTheta {\nabla} ^2\varPhi ) \,\text{d}V, \end{align}

(4.11) \begin{align} &={\textit{Pe}} \displaystyle \int \boldsymbol{u} \boldsymbol{\cdot }[ \varPhi \boldsymbol{\nabla } \varTheta + \varTheta \boldsymbol{\nabla } \varPhi ] \,\text{d}V, \quad(\mbox{using}\,({2.1})) \end{align}

(4.12) \begin{align} &={\textit{Pe}}\displaystyle \int \boldsymbol{\nabla } \boldsymbol{\cdot }({\boldsymbol{u}}\varTheta \varPhi )\, \text{d}V,\quad \mbox{(Incompressibility)} \end{align}

(4.13) \begin{align} &={\textit{Pe}}\displaystyle \int _{S_p} ({\boldsymbol{u}} \boldsymbol{\cdot }{\boldsymbol{n}})\varTheta \varPhi\, \text{d}S =\,{\textit{Pe}}\displaystyle \int _{S_p} ({\boldsymbol{u}} \boldsymbol{\cdot }{\boldsymbol{n}}) \,\text{d}S = 0, \end{align}

provided there is no change in the particle or drop volume.

Evaluation of the integral in (4.9) leads to ${\textit{Nu}}/{\textit{Pe}}^{{1}/{3}}$ , which is plotted in figure 10(d) for $\epsilon \in (-2,-1]$ – the transport rate decreases slightly in going from axisymmetric to planar extension. Similar to the small- ${\textit{Pe}}$ scenario, this decrease can almost entirely be eliminated on using $(\boldsymbol{E}\!:\!\boldsymbol{E})^{{1}/{2}}$ as the shear rate; one then finds ${\textit{Nu}} \approx 0.9[({\boldsymbol{E}}\!:\!{\boldsymbol{E}})^{{1}/{2}}a^2/D]^{{1}/{3}}$ (Batchelor Reference Batchelor1979). Both the (heat-transport) experiments and computations in figure 1(b) are broadly consistent with the prediction. While the former are limited to ${\textit{Pe}} \lesssim 100$ , the large- ${\textit{Pe}}$ scaling has also been verified in later mass transport experiments (Kutateladze et al. Reference Kutateladze, Nakoryakov and Iskakov1982) which, for reasons mentioned in § 1, access larger ${\textit{Pe}}$ values (with ${\textit{Re}}$ small). The comparison of theory with these later experiments is shown in figure 10(e), where one observes a small but systematic discrepancy. While a discrepancy may arise from multiple sources (experimental errors, higher-order corrections to ${\textit{Nu}}$ , ${\textit{Re}}$ not being small enough), the one here likely arises from differences in the theoretical and experimental flow configurations. The experiments measure by electrochemical means (in the diffusion-controlled regime) the transport from a stationary sphere in a simple shear, realised in the narrow gap between co-rotating cylinders. As will be seen in § 4.1.3, the absence of sphere rotation is crucial; the fixed-sphere velocity field requires adding a term, $\boldsymbol{\varOmega } \boldsymbol{\cdot }\boldsymbol{x}/r^3$ to (4.1), reflective of an imposed torque, and the resulting surface stress-line topology is shown as an inset. The theoretical line is the above-mentioned ${\textit{Nu}}$ -expression, involving the second invariant of $\boldsymbol{E}$ , and can correspond to any of the stress-line topologies in figure 10(a–c). All of these topologies, however, differ from the experimental pattern with asymmetrical octants. The correct, but as yet unknown, prediction for the asymmetrical-octant topology ought to lie below the measurements in figure 10(e) – this would be consistent with the sign of higher-order corrections to be discussed in § 4.1.2.

We now examine the second scenario, $\hat {\alpha } \gg {\textit{Pe}}^{-{1}/{3}}$ , corresponding to virtually all linear flows with vorticity in the limit ${\textit{Pe}} \rightarrow \infty$ . For such flows, $\boldsymbol{\varOmega }$ in (4.1) is non-zero and the sphere rotates, with the surface streamlines being circles transverse to the rotation axis (along $\hat {\boldsymbol{\omega }}$ ). The $O({\textit{Pe}}^{-{1}/{3}})$ straining component in (4.1) leads to tightly spiralling near-surface streamlines (see first and third rows in figure 16); the spiralling direction, whether along or against the rotation axis, depends on the sign of $E_\omega$ . Denoting the polar and azimuthal angles, in a spherical coordinate system with axis along $\hat {\boldsymbol{\omega }}$ , as $(\theta ,\phi )$ , and using the expansion $\varTheta = \varTheta _0 + {\textit{Pe}}^{-{1}/{3}}\varTheta _1 + o({\textit{Pe}}^{-{1}/{3}})$ , (2.1) at leading order reduces to

(4.14) \begin{equation} \hat {\alpha } \frac {\partial \varTheta _0}{\partial \phi } = 0, \end{equation}

at steady state. Thus, $\varTheta _0$ is independent of $\phi$ , and at $O({\textit{Pe}}^{-{1}/{3}})$ , one therefore obtains

(4.15) \begin{equation} \frac {15}{2}(\boldsymbol{E}:\boldsymbol{n}\boldsymbol{n})Y^2\frac {\partial \varTheta _0}{\partial Y} + 5Y\big(\boldsymbol{E}: \boldsymbol{n}\hat {\boldsymbol{\theta }}\big) \frac {\partial \varTheta _0}{\partial \theta } + \hat {\alpha } \frac {\partial \varTheta _1}{\partial \phi } = \frac {\partial ^2 \varTheta _0}{\partial Y^2}, \end{equation}

where $\boldsymbol{n} = \sin \theta \cos \phi \hat {\boldsymbol{\omega }}^\perp + \sin \theta \cos \phi \hat {\boldsymbol{\omega }}^{\perp \prime } + \cos \theta \hat {\boldsymbol{\omega }}$ , $\hat {\boldsymbol{\theta }} = \cos \theta \cos \phi \hat {\boldsymbol{\omega }}^\perp + \cos \theta \sin \phi \hat {\boldsymbol{\omega }}^{\perp \prime }{}-\sin \theta \hat {\boldsymbol{\omega }}$ are the unit radial and (polar) tangential vectors in terms of the orthonormal unit vector triplet $[\hat {\boldsymbol{\omega }}^{\perp },\hat {\boldsymbol{\omega }}^{\perp \prime },\hat {\boldsymbol{\omega }}]$ . Integrating (4.15) over $\phi \,\in \,[0,2\pi ]$ , and using the single-valuedness of $\varTheta _1$ , yields the following $\phi$ -averaged equation:

(4.16) \begin{equation} \frac {15}{2}Y^2\left(\frac {1}{2\pi }\int _0^{2\pi }\boldsymbol{E}:\boldsymbol{n}\boldsymbol{n}\,\text{d}\phi \right)\frac {\partial \varTheta _0}{\partial Y} + 5Y\left(\frac {1}{2\pi } \int _0^{2\pi }\boldsymbol{E}: \boldsymbol{n}\hat {\boldsymbol{\theta }} \,\text{d}\phi \right) \frac {\partial \varTheta _0}{\partial \theta } = \frac {\partial ^2 \varTheta _0}{\partial Y^2}, \end{equation}

where $({1}/({2\pi }))\int _0^{2\pi }\boldsymbol{E}:\boldsymbol{n}\boldsymbol{n}\,\text{d}\phi =( {E_\omega }/{2})(3\cos ^2\theta - 1)$ , $({1}/({2\pi }))\int _0^{2\pi }\boldsymbol{E}:\boldsymbol{n}\hat {\boldsymbol{\theta }}\,\text{d}\phi = E_\omega \cos \theta \sin \theta$ , with $E_\omega =\epsilon -(1+2\epsilon )\sin ^2\theta _\omega + (2+\epsilon )\sin ^2\theta _\omega \sin ^2\phi _\omega$ . Thus, $\varTheta _0$ and the leading order transport rate are determined by a $\phi$ -averaged or ‘axisymmetrised’ extensional flow with the rate of strain tensor $E_\omega \hat {\boldsymbol{\omega }}\hat {\boldsymbol{\omega }} -{E_\omega }/{2}({\boldsymbol{I}} -\hat {\boldsymbol{\omega }}\hat {\boldsymbol{\omega }})$ . The $\phi$ -variations seeded by the original straining component are smeared out by a shear-enhanced diffusion mechanism arising from vorticity-induced rotation; the time scale for this axisymmetrisation being given in § 5.2. Herein, we only note the resulting non-monotonic approach to the eventual steady state, in contrast to the previous extensional flows.

The steps from hereon are again standard and simpler than the extensional flow case, owing to the axisymmetric $\phi$ -averaged stress-line topology. In terms of the similarity variable, $\eta = Y/g(\cos \theta )$ , $g(\cos \theta )$ denoting the axisymmetric boundary layer thickness field, $\varTheta _0(\eta )$ is still given by (4.7). One finds ${\textit{Nu}} = {1}/({2\varGamma ({4}/{3})}) \int _{-1}^1 \text{d}(\cos \theta )/{g}$ . The expression for $g$ depends on $\text{sign}(E_\omega )$ , being given by $[ ({3}/{25})(({2}/{\sin ^2\theta }) + ({{\mathcal B}(\cos ^2\theta ;{3}/{4},{1}/{2})}/(\sin ^3\theta |\cos \theta |^{{3}/{2}})))]^{{1}/{3}}$ for $E_\omega \gt 0$ (uniaxial extension) and $ [ ({3}/{25})( ({-2}/{\sin ^2\theta )} + ({{\mathcal B}({3}/{4},{1}/{2})- ({\mathcal B}(\cos ^2\theta ;{3}/{4},{1}/{2}))}/({\sin ^3\theta |\cos \theta |^{{3}/{2}}}))) ]^{{1}/{3}}$ for $E_\omega \lt 0$ (biaxial extension); ${\mathcal B}(a,b)$ and ${\mathcal B}(z;a,b)$ being the complete and incomplete Beta functions. The resulting ${\textit{Nu}}$ is the same for both flows, as expected from flow-reversal invariance, and is given by

(4.17) \begin{equation} {\textit{Nu}} = 0.97(E_\omega {\textit{Pe}})^{{1}/{3}} \end{equation}

for any non-zero $\hat {\alpha }$ , for ${\textit{Pe}} \rightarrow \infty$ . Interpreting sub-Kolmogorov turbulence as a fluctuating linear flow, the previous result (with $E_\omega$ replaced by $\langle E_\omega \rangle$ , $\langle .\rangle$ denoting an average over turbulent statistics), was used by Batchelor (Reference Batchelor1980) to evaluate the mass transport rate from a sub-Kolmogorov spherical particle in turbulence. While we comment on transport from particles (drops) in turbulent flows towards the end in § 6, it is worth pointing out that (4.17), when used for turbulence, relies on a non-trivial correlation between $\boldsymbol{E}$ and $\boldsymbol{\omega }$ at the Kolmogorov scale. There is indeed such a correlation – vortex stretching constrains $E_\omega$ to be positive. In fact, the preponderance of sub-Kolmogorov linear flows have a biaxial extensional character, with $\boldsymbol{\omega }$ aligned, on average, with the principal axis corresponding to the intermediate (positive) eigenvalue of $\boldsymbol{E}$ (Meneveau Reference Meneveau2011).

In summary, ${\textit{Nu}}$ for a spherical particle, in almost any ambient linear flow, does not directly depend on $\hat {\alpha }$ , at leading order, in the limit ${\textit{Pe}} \rightarrow \infty$ . Rather, this flow-type parameter only serves to demarcate the linear flows into two classes: (i) linear extensional flows where ${\textit{Nu}}$ is a function of $\epsilon$ alone; and (ii) vortical linear flows where ${\textit{Nu}}$ only depends on the combination of $\epsilon$ , $\theta _\omega$ and $\phi _\omega$ appearing in $E_\omega$ . The dependence on $E_\omega$ , rather than $E$ , may be interpreted as a rotation-induced suppression of transport. This is because, for (4.17) to be valid, one requires $E_\omega {\textit{Pe}} \gg 1$ , a condition that is restrictive when $E_\omega$ is small. As will be seen in § 4.1.3, rotation-induced suppression is complete for ambient planar linear flows with $E_{\omega } = 0$ .

4.1.2. Higher-order corrections

Returning to figure 1(b), the fit proposed for the computations therein (Karp-Boss et al. Reference Karp-Boss, Boss and Jumars1996) has the general form ${\textit{Nu}}=b {\textit{Pe}}^c+a$ , with $c$ equalling $1/3$ for large ${\textit{Pe}}$ , and the values of $b$ and $a$ depending on the ${\textit{Pe}}$ -interval examined; for example, $b=0.9$ (consistent with use of $(\boldsymbol{E:E})^{{1}/{2}}$ as the shear rate) and $a=0.56$ for $50\lt {\textit{Pe}}\lt 90$ , which is suggestive of an $O(1)$ correction to the leading order predictions obtained previously. To understand this, note that the analysis leading to (4.9) or (4.17) breaks down in the wake region, owing to a divergence of the boundary layer thickness ( $g$ ) at the angular location corresponding to the wake at ${\textit{Pe}} = \infty$ . Now, for any finite ${\textit{Pe}}$ , the wake has a finite angular spread. Further, the divergence leads to the scalar gradient becoming vanishingly small in the wake, as opposed to being small but finite. Both these point to the boundary-layer-based ${\textit{Nu}}$ predictions being an underestimate and to the next order correction being positive. Apart from the nature of the wake contribution, analysing corrections to (4.9) or (4.17) would help ascertain the range of validity of the boundary-layer-based prediction. For instance, knowing the first two terms in the large- ${\textit{Pe}}$ expansion would enable calculation of a (lower) threshold ${\textit{Pe}}$ at which they become comparable in magnitude.

It is easier to analyse corrections to (4.17) on account of the axisymmetry. For the uniaxial case, the expression for $g$ in the text before (4.17) diverges at the poles $\theta =0,\pi$ . Due to symmetry, one only need consider $\theta =0$ , in which case, $\lim _{\theta \rightarrow 0} g \approx (1-\cos \theta )^{-{1}/{2}}$ , with $g$ diverging as $\theta ^{-1}$ . The boundary layer analysis breaks down when the radial scale of $O({\textit{Pe}}^{-{1}/{3}}g) \equiv O({\textit{Pe}}^{-{1}/{3}}\theta ^{-1})$ , associated with the thickening boundary layer, becomes comparable to the shrinking angular scale $\theta$ . This happens when $\theta \sim O({\textit{Pe}}^{-{1}/{6}})$ , which gives the angular extent of the wake regions around $\theta = 0$ and $\pi$ . Their contribution to the transport rate is $\boldsymbol{\nabla }\varTheta _{wake}.S_{wake} \sim 1/[{\textit{Pe}}^{{1}/{3}}g(\theta )].\theta ^2$ with $\theta \sim O({\textit{Pe}}^{-{1}/{6}})$ , which comes out to be $O({\textit{Pe}}^{-{1}/{6}})$ . Thus, the ${\textit{Nu}}$ expansion for uniaxial extension proceeds as ${\textit{Nu}} = O({\textit{Pe}}^{{1}/{3}}) + O(1) + O({\textit{Pe}}^{-{1}/{6}})$ , with the $O(1)$ correction also originating within the boundary layer; the associated boundary-layer-wake configuration is shown in figure 9(a). For biaxial extension, the wake is in the neighbourhood of $\theta = {\pi }/{2}$ , collapsing to the equatorial fixed-point ring in figure 10(c) for ${\textit{Pe}} \rightarrow \infty$ (as opposed to a point for uniaxial extension). Now, $g$ diverges as $({\pi }/{2}-\theta )^{-{1}/{2}}$ for $\theta \rightarrow {\pi }/{2}$ , and arguments along the previous lines lead to a wake width of $O({\textit{Pe}}^{-{2}/{9}})$ , which makes an $O(1)$ contribution to the transport rate. The differing origin of the $O(1)$ corrections to ${\textit{Nu}}$ , for the uniaxial and biaxial cases, clearly point to the differing nature of the surface-flux distribution for large ${\textit{Pe}}$ . Yet, flow-reversal invariance implies that ${\textit{Nu}}$ remains unchanged! The $O(1)$ correction has only been calculated for a uniform flow and is positive (Acrivos & Goddard Reference Acrivos and Goddard1965). In terms of the local flow in the wake region, uniform flow is similar to uniaxial extension, and the associated ${\textit{Nu}}$ expansions must be similar, with the $O(1)$ correction for uniform flow also originating from the boundary layer region. As argued in the context of figure 10(e), the positive sign expected for the $O(1)$ correction implies that, in contrast to what is shown in the figure, the experiments likely approach the (correct) boundary-layer prediction from above.