A rigorous bound on the scaling of dissipation with velocity amplitude in flow past a sphere

Abstract

A rigorous bound on the energy dissipation in the flow around a sphere is derived analytically. The bound is allowed to depend on the viscosity of the fluid, the radius of the sphere, the free-stream velocity and the amplitude of the difference between the actual flow field and the irrotational flow around the sphere. The bound becomes asymptotically independent of viscosity as the Reynolds number tends to infinity. Provided that the velocity fluctuations are of the same order as the free-stream velocity, the scaling of the bound is consistent with a Reynolds-number-independent drag coefficient at high Reynolds numbers.

1. Introduction

A fundamental problem of engineering fluid mechanics is the computation of the drag force on an obstacle at high Reynolds numbers. A paradigmatic example of this class of problems is the flow of fluid with kinematic viscosity $\nu$ and density $\rho$ around a sphere of radius $R$ with a flow velocity $U_\infty$ at large distances from the sphere. According to common reasoning, the momentum transport by turbulent eddies is so much more efficient than viscous friction that viscosity becomes an irrelevant quantity if only the Reynolds number is high enough. It follows from this assumption that the time-averaged drag force $F$ should become independent of $\nu$ and only depend on $\rho$, $U_\infty$ and $R$. In that situation, $F \propto \rho U_\infty ^2 R^2$ because this is the only combination of parameters with the dimension of a force, and the power dissipated in the fluid, $F U_\infty$, must scale as $\rho U_\infty ^3 R^2$.

The same dimensional argument pervades all empirical results on turbulence. Most importantly, a relation dating back to Kolmogorov (1941) states that energy dissipated per unit mass, $\epsilon$, in homogeneous and isotropic turbulence is determined by a mixing length $l$ and the amplitude of velocity fluctuations $U$ as $\epsilon \propto U^3/l$, again because this is the only plausible combination of parameters independent of $\nu$ with the right dimension. This relation will be called Kolmogorov scaling for short in this paper.

While the Kolmogorov scaling has proven very successful in practice and beyond homogeneous and isotropic turbulence, there is little theory based on the Navier–Stokes equation to support it. The only exact statements available come in the form of inequalities providing bounds on the energy dissipation. These bounds indeed have the expected form for turbulence driven by a volume force (Childress, Kerswell & Gilbert 2001; Doering & Foias 2002; Rollin, Dubief & Doering 2011; Tilgner 2017).

We also know techniques for deriving bounds in confined geometry. These use either the optimisation procedure demonstrated by Howard (1972) and Busse (1979), or a technique known as the background field method. This second option was applied to a variety of flows, including purely hydrodynamic settings such as Couette or channel flow (Doering & Constantin 1992; Constantin & Doering 1995), but also convection (Doering & Constantin 1996; Kerswell 2001; Siggers, Kerswell & Balmforth 2004; Plasting & Ierley 2005; Balmforth et al. 2006; Doering, Otto & Reznikoff 2006; Otto & Seis 2011; Whitehead & Doering 2012) and magnetohydrodynamics (Alboussière 2009). The bounds in these references are derived by analytical means with little or no numerical help, whereas much of the more recent work relied on the numerical solution of some optimisation problem. In the purely hydrodynamic examples, these results allow us, as a by-product, to verify whether the bounds on energy dissipation conform with the Kolmogorov scaling. No contradiction was found.

The situation is less clear for flows around obstacles. Kumar & Garaud (2020) recently presented the calculation of a bound on the energy dissipation in the flow around a thin plate of finite extent at zero angle of attack. Their result is again compatible with the Kolmogorov scaling, but their computation cannot be extended to the flow around bluff bodies. On the other hand, the oldest known application of the background field method goes back to work by Hopf (1941) in which he examines the flow caused by moving boundaries in a viscous fluid at rest at large distances from the boundaries. His paper contains a recipe for constructing a background field in nearly arbitrary geometries and he notes that if the geometry remains unchanged during the course of the motion, such as for a rotating cylinder, his method leads to bounds on the dissipation that vary algebraically with $\nu ^{-1}$, whereas in general, they depend exponentially on $\nu ^{-1}$.

Returning to a frame of reference in which the obstacle is at rest, Hopf's method allows us to find a bound on the drag force on a sphere which depends exponentially on the Reynolds number $U_\infty R/\nu$ and therefore does not tend to a limit independent of $\nu$. This raises the question of whether a state of boundary-driven turbulence exists in which the Kolmogorov scaling does not hold.

The bounds summarised above reflect two types of expressions with which the dissipation is computed empirically in turbulent flows. Both forms are independent of viscosity and answer two related, yet different questions. The two questions can be phrased as follows in the context of the flow around a sphere: suppose one is given the control parameters $\rho$, $R$ and $U_\infty$, but not $\nu$, can one deduce the dissipation in the flow assuming a high Reynolds number? The empirical answer to this first question is a dissipation proportional to $\rho U_\infty ^3 R^2$, but we have no bound compatible with this expression because of the viscosity dependence in Hopf's bound. On the other hand, we do have bounds in terms of control parameters independent of viscosity for the other flows mentioned previously. The second form parallels the Kolmogorov phenomenology of homogeneous and isotropic turbulence and answers another question: suppose one is given $\rho$, $R$ and the ‘velocity at large scales’ $U$, but not $\nu$, can one deduce the dissipation in the flow, assuming again a high Reynolds number? There is some slack in what exactly is meant by the velocity at large scales. Surely, $U_\infty$ itself is acceptable for $U$, in which case the first and second questions are identical. However, another common interpretation of $U$ typically found in turbulence models (Pope 2000) is that $U$ is computed from the square root of kinetic energy. Hopf's result leaves the possibility that a bound of this second form does not exist for the flow around bluff bodies. A flow whose dissipation violates a bound of the form $\rho U^3 R^2$ would contradict the universality of the Kolmogorov scaling.

The present paper deals with the second of these questions and computes a bound on the dissipation in the flow around a sphere that is allowed to depend on the energy of the flow. This is different from the problem of finding a bound that depends only on $U_\infty R/\nu$ because the relation between the kinetic energy and $U_\infty$ does not have to be determined. The calculation will employ the background field method. The background field will be chosen to be the irrotational flow past the obstacle, together with a boundary-layer correction to adapt the potential flow to the no-slip boundary conditions on the surface of the obstacle. It will be shown that the idea behind the Kolmogorov scaling remains intact, and that the bound found by Hopf can only be realised if the kinetic energy in the velocity field relative to the irrotational flow grows exponentially with the Reynolds number. The calculation is explicitly performed for a sphere, but the general procedure is applicable to bodies of different shape as well.

2. Calculation of the bound

A sphere of radius $R$ centred at the origin is surrounded by fluid of viscosity $\nu$ and density $\rho$ which moves in the $z$-direction with velocity $U_\infty$ at infinity. The velocity and pressure fields, $\boldsymbol {v}(\boldsymbol {r},t)$ and $p(\boldsymbol {r},t)$, are determined as a function of position $\boldsymbol {r}$ and time $t$ by the Navier–Stokes equation

(2.1a,b)\begin{equation} \partial_t \boldsymbol{v} + (\boldsymbol{v} \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v} ={-} \frac{1}{\rho} \boldsymbol{\nabla} p +\nu {\nabla}^2 \boldsymbol{v} , \quad \boldsymbol{\nabla} \boldsymbol{\cdot} \boldsymbol{v} = 0 \end{equation}

with the boundary conditions that $\boldsymbol {v} = 0$ on $r=|\boldsymbol {r}|=R$ and $\boldsymbol {v} \rightarrow U_\infty \hat {\boldsymbol {z}}$ for $r \rightarrow \infty$, where a hat denotes a unit vector. We choose a background field $\boldsymbol {v}_B(\boldsymbol {r})$ as the sum $\boldsymbol {v}_B = \boldsymbol {v}_P + \boldsymbol {v}_b$ of the potential flow around the obstacle, $\boldsymbol {v}_P(\boldsymbol {r})$, and a boundary-layer contribution $\boldsymbol {v}_b(\boldsymbol {r})$ such that $\boldsymbol {v}_P + \boldsymbol {v}_b = 0$ on $r=R$ and $\boldsymbol {v}_b = 0$ for $r \geqslant R+\delta$, where $\delta$ is a boundary-layer thickness whose size is available for later optimisation. These fields have the properties

(2.2)\begin{equation} \boldsymbol{\nabla} \boldsymbol{\cdot} \boldsymbol{v}_B = \boldsymbol{\nabla} \boldsymbol{\cdot} \boldsymbol{v}_P = \boldsymbol{\nabla} \boldsymbol{\cdot} \boldsymbol{v}_b =0 \end{equation}

(2.3)\begin{equation}\boldsymbol{\nabla} \times\boldsymbol{v}_P =0 \end{equation}

and fulfil the boundary conditions $\boldsymbol {v}_B=0$ and $\hat {\boldsymbol {r}} \boldsymbol {\cdot } \boldsymbol {v}_P=0$ on $r=R$, whereas at infinity, $\boldsymbol {v}_B=\boldsymbol {v}_P=U_\infty \hat {\boldsymbol {z}}$.

The total velocity field is written as $\boldsymbol {v}(\boldsymbol {r},t)=\boldsymbol {v}_B(\boldsymbol {r})+ \boldsymbol {u}(\boldsymbol {r},t)$, so that the fluctuation field $\boldsymbol {u}$ obeys $\boldsymbol {\nabla } \boldsymbol {\cdot } \boldsymbol {u} = 0$ and $\boldsymbol {u} = 0$ both on $r=R$ and at infinity. Inserting this decomposition into the Navier–Stokes equation and taking the dot product with $\boldsymbol {u}$ leads to

(2.4)\begin{align} &\frac{1}{2} \partial_t |\boldsymbol{u}|^2 + \frac{1}{2} \boldsymbol{\nabla} \boldsymbol{\cdot} (\boldsymbol{u} |\boldsymbol{u}|^2) +\boldsymbol{u} \boldsymbol{\cdot} (\boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_B + \frac{1}{2} \boldsymbol{\nabla} \boldsymbol{\cdot} (\boldsymbol{v}_B |\boldsymbol{u}|^2) +\boldsymbol{u} \boldsymbol{\cdot} (\boldsymbol{v}_B \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_B \nonumber\\ &\quad ={-}\frac{1}{\rho} \boldsymbol{\nabla} \boldsymbol{\cdot} (\boldsymbol{u} p) + \nu \boldsymbol{u} \boldsymbol{\cdot} {\nabla}^2 \boldsymbol{v}. \end{align}

If we integrate this equation over the entire fluid volume, all divergence terms disappear thanks to the boundary conditions, and

(2.5)\begin{align} \int \boldsymbol{u} \boldsymbol{\cdot} (\boldsymbol{v}_B \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_B \,\mathrm{d}V &= \int \boldsymbol{u} \boldsymbol{\cdot} \left[(\boldsymbol{\nabla} \times \boldsymbol{v}_B)\times \boldsymbol{v}_B + \tfrac{1}{2} \boldsymbol{\nabla} |\boldsymbol{v}_B|^2 \right] \,\mathrm{d}V \nonumber\\ &= \int \boldsymbol{u} \boldsymbol{\cdot} [(\boldsymbol{\nabla} \times \boldsymbol{v}_b)\times \boldsymbol{v}_B] \,\mathrm{d}V. \end{align}

In addition, applying a few times the vector identity $\boldsymbol {b} \boldsymbol {\cdot } (\boldsymbol {\nabla } \times \boldsymbol {a}) = \boldsymbol {\nabla } \boldsymbol {\cdot } (\boldsymbol {a} \times \boldsymbol {b}) + \boldsymbol {a} \boldsymbol {\cdot } (\boldsymbol {\nabla } \times \boldsymbol {b})$, one finds

(2.6)\begin{align} \int \boldsymbol{v} \boldsymbol{\cdot} \boldsymbol{\nabla}^2 \boldsymbol{v} \,\textrm{d}V & ={-} \int \boldsymbol{v} \boldsymbol{\cdot} \boldsymbol{\nabla} \times \boldsymbol{\nabla} \times \boldsymbol{v} \,\mathrm{d}V ={-} \int |\boldsymbol{\nabla} \times \boldsymbol{v}|^2 \,\mathrm{d}V\nonumber\\ & = \int (\boldsymbol{u} \boldsymbol{\cdot} {\nabla}^2 \boldsymbol{u} + \boldsymbol{v}_B \boldsymbol{\cdot} {\nabla}^2 \boldsymbol{v}_B + \boldsymbol{u} \boldsymbol{\cdot} {\nabla}^2 \boldsymbol{v}_B + \boldsymbol{v}_B \boldsymbol{\cdot} {\nabla}^2 \boldsymbol{u}) \,\mathrm{d}V\nonumber\\ & ={-} \int |\boldsymbol{\nabla} \times \boldsymbol{u}|^2 \,\textrm{d}V - \int |\boldsymbol{\nabla} \times \boldsymbol{v}_B|^2 \,\mathrm{d}V + 2 \int \boldsymbol{u} \boldsymbol{\cdot} {\nabla}^2 \boldsymbol{v}_B \,\mathrm{d}V \end{align}

so that

(2.7)\begin{align} \int \boldsymbol{u} \boldsymbol{\cdot} {\nabla}^2 \boldsymbol{v} \,\mathrm{d}V & ={-} \int |\boldsymbol{\nabla} \times \boldsymbol{u}|^2 \,\mathrm{d}V + \boldsymbol{u} \boldsymbol{\cdot} {\nabla}^2 \boldsymbol{v}_B \,\mathrm{d}V\nonumber\\ & ={-} \tfrac{1}{2} \int |\boldsymbol{\nabla} \times \boldsymbol{v}|^2 \,\mathrm{d}V - \tfrac{1}{2} \int |\boldsymbol{\nabla} \times \boldsymbol{u}|^2 \,\mathrm{d}V + \tfrac{1}{2} \int |\boldsymbol{\nabla} \times \boldsymbol{v}_B|^2 \,\mathrm{d}V. \end{align}

Let $V_\delta$ denote the volume $R \leqslant r \leqslant R+\delta$. Inserting $\boldsymbol {v}_B = \boldsymbol {v}_P + \boldsymbol {v}_b$, (2.5) and (2.7) into (2.4), and keeping track of integrands which are zero outside $V_\delta$ and which need to be integrated only over $V_\delta$, we obtain

(2.8)\begin{align} &\partial_t \int \frac{1}{2} |\boldsymbol{u}|^2 \,\mathrm{d}V + \int \boldsymbol{u} \boldsymbol{\cdot}(\boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_P \,\mathrm{d}V + \frac{\nu}{2} \int |\boldsymbol{\nabla} \times \boldsymbol{v}|^2 \,\mathrm{d}V - \frac{\nu}{2} \int_{V_\delta} |\boldsymbol{\nabla} \times \boldsymbol{v}_b|^2 \,\mathrm{d}V\nonumber\\ &\quad ={-} \frac{\nu}{2} \int |\boldsymbol{\nabla} \times \boldsymbol{u}|^2 \,\mathrm{d}V - \int_{V_\delta} \boldsymbol{u} \boldsymbol{\cdot}(\boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_b \,\mathrm{d}V - \int_{V_\delta} \boldsymbol{u} \boldsymbol{\cdot} [(\boldsymbol{\nabla} \times \boldsymbol{v}_b) \times \boldsymbol{v}_B] \,\mathrm{d}V\nonumber\\ &\quad \leqslant - \frac{\nu}{2} \int |\boldsymbol{\nabla} \times \boldsymbol{u}|^2\, \mathrm{d}V + \left|\int_{V_\delta} \boldsymbol{u} \boldsymbol{\cdot}(\boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_b \,\mathrm{d}V\right|\nonumber\\ & \qquad + \sqrt{ \int_{V_\delta} |\boldsymbol{u}|^2 \,\mathrm{d}V} \sqrt{ \int_{V_\delta} |(\boldsymbol{\nabla} \times \boldsymbol{v}_b) \times \boldsymbol{v}_B|^2 \,\mathrm{d}V}, \end{align}

where the Schwarz inequality was used in the last line.

We now have to specify expressions for $\boldsymbol {v}_P$ and $\boldsymbol {v}_b$ in order to evaluate the integrals in the last equation. The potential flow around a sphere is given by

(2.9)\begin{equation} \boldsymbol{v}_P = \hat{\boldsymbol{r}} U_\infty \cos \theta \left[ 1-\left( \frac{R}{r} \right)^3 \right] - \hat{\boldsymbol{\theta}} U_\infty \sin \theta \left[ 1 + \frac{1}{2}\left( \frac{R}{r} \right)^3 \right] \end{equation}

in spherical polar coordinates $(r, \theta , \varphi )$. Let us choose $\boldsymbol {v}_b$ as

(2.10)\begin{equation} \boldsymbol{v}_b = \left\{\begin{array}{@{}ll} \boldsymbol{\nabla} \times \boldsymbol{\nabla} \times \left[ -\dfrac{3}{2} U_\infty \cos \theta R(r-R) \left( 1 - \dfrac{r-R}{\delta} \right)^2\hat{\boldsymbol{r}} \right], & R \leqslant r \leqslant R+\delta \\ 0, & r \geqslant R+\delta \end{array} \right. . \end{equation}

With this choice, $\boldsymbol {v}_b$ is continuous everywhere, and all the integrations in the following can be carried out in closed form. However, we will finally be interested in the limit of high Reynolds numbers and the limit $\delta /R \ll 1$. The results of all algebraic manipulations will therefore be given only to the leading order in $\delta /R$ because this leads to considerably more transparent expressions. This procedure is sufficient for the ultimate goal of the calculation and alleviates the at any rate tedious algebra. For instance,

(2.11)\begin{equation} \int_{V_\delta} |\boldsymbol{\nabla} \times \boldsymbol{v}_b|^2 \,\mathrm{d}V = 24 {\rm \pi}U_\infty^2 \frac{R^2}{\delta} \left[1+{O}\left( \frac{\delta^2}{R^2} \right) \right] \end{equation}

and

(2.12)\begin{equation} \int_{V_\delta} |(\boldsymbol{\nabla} \times \boldsymbol{v}_b) \times \boldsymbol{v}_B|^2 \,\mathrm{d}V = \frac{3816}{175} {\rm \pi}U_\infty^4 \frac{R^2}{\delta} \left[1+{O}\left( \frac{\delta}{R} \right) \right]. \end{equation}

Let us call $\boldsymbol{\mathsf{e}}_b$ the deformation-rate tensor of the boundary-layer flow $\boldsymbol {v}_b$ so that ${e}_{b,ij} = (\partial _j v_{b,i} + \partial _i v_{b,j})/{2}$ and $\boldsymbol {u} \boldsymbol {\cdot }(\boldsymbol {u} \boldsymbol {\cdot } \boldsymbol {\nabla }) \boldsymbol {v}_b = \boldsymbol {u} \boldsymbol {\cdot } \boldsymbol{\mathsf{e}}_b \boldsymbol {u}$. If we write down $\boldsymbol{\mathsf{e}}_b$ in spherical coordinates and order the components within $\boldsymbol {u}$ in the sequence $(u_r, u_\theta , u_\varphi )$, we find

(2.13)\begin{equation} \boldsymbol{\mathsf{e}}_b = f_1(r) \boldsymbol{\mathsf{M}}_1 + f_2(r) \boldsymbol{\mathsf{M}}_2 \end{equation}

with matrices $\boldsymbol{\mathsf{M}}_1$ and $\boldsymbol{\mathsf{M}}_2$ defined by

(2.14)\begin{equation} \boldsymbol{\mathsf{M}}_1 = \begin{bmatrix} 0 & \sin \theta & 0 \\ \sin \theta & 0 & 0 \\ 0 & 0 & 0 \end{bmatrix} \end{equation}

and

(2.15)\begin{equation} \boldsymbol{\mathsf{M}}_2 = \begin{bmatrix} 2 \cos \theta & \sin \theta & 0 \\ \sin \theta & - \cos \theta & 0 \\ 0 & 0 & - \cos \theta \end{bmatrix} . \end{equation}

The functions $f_1(r)$ and $f_2(r)$ are of course zero for $r \geqslant R+\delta$ and, using the shorthand notation $s=r-R$, they are given within $V_\delta$ by

(2.16)\begin{equation} f_1(r)={-} \frac{3 R U_\infty}{2 \delta^{2} \left(R + s\right)^{3}} \left[R^{2} \left(2 \delta - 3 s\right) + R \left(4 \delta s - 6 s^{2}\right) + \delta^{2} s - 2 s^{3}\right] \end{equation}

and

(2.17)\begin{equation} f_2(r)= \frac{3 R U_\infty}{\delta^{2} \left(R + s\right)^{3}} \left[R \left(- \delta^{2} + 4 \delta s - 3 s^{2}\right) + \delta^{2} s - s^{3}\right] . \end{equation}

The eigenvalues $\lambda _1$ of $\boldsymbol{\mathsf{M}}_1$ obey $|\lambda _1| \leqslant 1$, whereas for the eigenvalues $\lambda _2$ of $\boldsymbol{\mathsf{M}}_2$, the inequality $|\lambda _2| \leqslant 2$ holds, so that

(2.18)\begin{equation} |\boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\mathsf{e}}_b \boldsymbol{u} | \leqslant |\boldsymbol{u}|^2 \left\{ |f_1(r)| + 2 |f_2(r)| \right\} \leqslant |\boldsymbol{u}|^2 3 \frac{U_\infty}{\delta} \left[1+\textit{O}\left( \frac{\delta}{R} \right) \right] \end{equation}

and

(2.19)\begin{equation} \left|\int_{V_\delta} \boldsymbol{u} \boldsymbol{\cdot}(\boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_b \,\mathrm{d}V\right| \leqslant 3 \frac{U_\infty}{\delta} \int_{V_\delta} |\boldsymbol{u}|^2 \,\mathrm{d}V \left[1+{O}\left( \frac{\delta}{R} \right) \right] . \end{equation}

We now want to find the maximum of the last two lines in the sequence (2.8). To this end, we first relate $\int _{V_\delta } |\boldsymbol {u}|^2 \,\mathrm {d}V$ to $\int |\boldsymbol {\nabla } \times \boldsymbol {u}|^2 \,\mathrm {d}V = \int |\boldsymbol {\nabla } \boldsymbol {u}|^2 \,\mathrm {d}V$ with the notation $|\boldsymbol {\nabla } \boldsymbol {u}|^2 = \sum _{i,j} |\partial _j u_i|^2$. A simple approach to achieve this goal is to first note that

(2.20)\begin{equation} |\boldsymbol{u}(\boldsymbol{r},t)|^2 = \left|\int_0^s \partial_r \boldsymbol{u}(\boldsymbol{r}',t) \,\mathrm{d}s' \right|^2, \end{equation}

where $\boldsymbol {r} = (R+s, \theta , \varphi )$ and $\boldsymbol {r}' = (R+s', \theta , \varphi )$. Inserting a factor 1 in the integral in (2.20) and applying Schwarz's inequality yields

(2.21)\begin{equation} |\boldsymbol{u}(\boldsymbol{r},t)|^2 \leqslant \left( \sqrt{s} \sqrt{\int_0^s |\partial_r \boldsymbol{u}(\boldsymbol{r}',t)|^2 \,\mathrm{d}s'} \right)^2 \leqslant s \int_0^s |\boldsymbol{\nabla} \boldsymbol{u}|^2 \,\mathrm{d}s' . \end{equation}

The integral of $|\boldsymbol {u}|^2$ over $V_\delta$ can therefore be estimated as follows:

(2.22)\begin{align} \int_{V_\delta} |\boldsymbol{u}|^2 \,\mathrm{d}V & = \int_0^\delta \,\mathrm{d}s \int_0^{2{\rm \pi}} \,\mathrm{d} \varphi \int_0^{\rm \pi} \,\mathrm{d} \theta \sin \theta (R+s)^2 |\boldsymbol{u}(\boldsymbol{r},t)|^2\nonumber\\ & \leqslant \int_0^\delta \,\mathrm{d}s s (R+s)^2 \int_0^{2{\rm \pi}} \,\mathrm{d} \varphi \int_0^{\rm \pi} \,\mathrm{d} \theta \sin \theta \int_0^s \,\mathrm{d}s' \frac{(R+s')^2}{(R+s')^2} |\boldsymbol{\nabla} \boldsymbol{u}|^2\nonumber\\ & \leqslant \int_0^\delta \,\mathrm{d}s s (R+s)^2 \frac{1}{R^2} \int_{V_s} |\boldsymbol{\nabla} \boldsymbol{u}|^2 \,\mathrm{d}V, \end{align}

where $V_s$ denotes the volume $R \leqslant r \leqslant R+s$. As the integrand $|\boldsymbol {\nabla } \boldsymbol {u}|^2$ is always positive, we obtain another valid inequality by extending the volume of integration from $V_s$ to the entire fluid volume:

(2.23)\begin{align} \int_{V_\delta} |\boldsymbol{u}|^2 \,\mathrm{d}V& \leqslant \left\{ \int |\boldsymbol{\nabla} \boldsymbol{u}|^2 \,\mathrm{d}V \right\} \left\{ \int_0^s \,\mathrm{d}s s (R+s)^2 \frac{1}{R^2} \right\}\nonumber\\ & \leqslant \left\{ \int |\boldsymbol{\nabla} \times \boldsymbol{u}|^2 \,\mathrm{d}V \right\} \left\{ \frac{1}{2} \delta^2 \left[1+{O}\left( \frac{\delta}{R} \right) \right] \right\} . \end{align}

The leading order in $\delta /R$ of the last right-hand side in (2.8) is therefore bounded from above by

(2.24)\begin{equation} \left( -\frac{\nu}{2} + \frac{3}{2} U_\infty \delta \right) \int |\boldsymbol{\nabla} \times \boldsymbol{u}|^2 \,\mathrm{d}V + \frac{6}{35} \sqrt{371 {\rm \pi}} R U_\infty^2 \sqrt{\delta} \sqrt{\int |\boldsymbol{\nabla} \times \boldsymbol{u}|^2 \,\mathrm{d}V} . \end{equation}

We now regard $\int |\boldsymbol {\nabla } \times \boldsymbol {u}|^2 \,\mathrm {d}V$ as the variable in expression (2.24) and consider what the maximum of this expression is if $\int |\boldsymbol {\nabla } \times \boldsymbol {u}|^2\,\mathrm {d}V$ is varied. In order for a maximum to exist, the condition

(2.25)\begin{equation} 3 U_\infty \delta \leqslant \nu \end{equation}

needs to be satisfied. If this is the case, the maximum of (2.24) is achieved for

(2.26)\begin{equation} \sqrt{\int |\boldsymbol{\nabla} \times \boldsymbol{u}|^2 \,\mathrm{d}V} = \frac{6}{35} \sqrt{371 {\rm \pi}} R U_\infty^2 \sqrt{\delta} \frac{1}{\nu-3 U_\infty \delta} . \end{equation}

Inserting the maximum of (2.24) together with (2.11) into (2.8) yields

(2.27)\begin{align} &\partial_t \int \frac{1}{2} |\boldsymbol{u}|^2 \,\mathrm{d}V + \int \boldsymbol{u} \boldsymbol{\cdot}(\boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_P \,\mathrm{d}V + \frac{\nu}{2} \int |\boldsymbol{\nabla} \times \boldsymbol{v}|^2 \,\mathrm{d}V\nonumber\\ &\quad \leqslant 12 {\rm \pi}R^2 U_\infty^3 \left( \frac{\nu}{U_\infty \delta} - \frac{159}{350} \dfrac{1}{3-\dfrac{\nu}{U_\infty \delta}} \right) \left[ 1 + {O}\left( \dfrac{\delta}{R} \right) \right] . \end{align}

The choice of $\delta$ that minimises the right-hand side of this inequality, yet respects condition (2.25), is

(2.28)\begin{equation} \frac{U_\infty \delta}{\nu} = \frac{5}{2991}(210-\sqrt{2226}) \end{equation}

so that (2.8) now turns into

(2.29)\begin{align} &\partial_t \int \frac{1}{2} |\boldsymbol{u}|^2 \,\mathrm{d}V + \frac{\nu}{2} \int |\boldsymbol{\nabla} \times \boldsymbol{v}|^2 \,\mathrm{d}V \nonumber\\ &\quad \leqslant \left|\int \boldsymbol{u} \boldsymbol{\cdot}(\boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_P \,\mathrm{d}V\right| + 12 {\rm \pi}\left(3+\frac{\sqrt{2226}}{35} \right) R^2 U_\infty^3 \left[ 1 + {O}\left( \frac{\delta}{R} \right) \right] . \end{align}

Finally, if we compute the deformation-rate tensor for the potential flow $\boldsymbol {v}_P$, we may conclude that

(2.30)\begin{align} \left|\int \boldsymbol{u} \boldsymbol{\cdot}(\boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\nabla}) \boldsymbol{v}_P \,\mathrm{d}V\right| & = \left|\int \frac{3}{2} \frac{U_\infty}{r} \left( \frac{R}{r} \right)^3 \boldsymbol{u} \boldsymbol{\cdot} \boldsymbol{\mathsf{M}}_2 \boldsymbol{u} \,\mathrm{d}V\right|\nonumber\\ & \leqslant 3 \frac{U_\infty}{R} \int \left( \frac{R}{r} \right)^4 |\boldsymbol{u}|^2 \,\mathrm{d}V . \end{align}

Let us use angle brackets $\langle ~ \rangle$ to denote the time average of any function $f(t)$:

(2.31)\begin{equation} \langle f(t) \rangle = \lim_{T \rightarrow \infty} \frac{1}{T} \int_0^T f(t) \,\mathrm{d}t. \end{equation}

Assuming a statistically stationary state in which the time average of a time derivative disappears, we deduce from (2.29) that

(2.32)\begin{align} \nu \left\langle \int |\boldsymbol{\nabla} \times \boldsymbol{v}|^2 \,\mathrm{d}V \right\rangle &\leqslant 6 \frac{U_\infty}{R} \int \left( \frac{R}{r} \right)^4 \langle |\boldsymbol{u}|^2 \rangle\, \mathrm{d}V \nonumber\\ &\quad + 24 {\rm \pi}\left(3+\frac{\sqrt{2226}}{35} \right) R^2 U_\infty^3 \left[ 1 + {O}\left( \frac{\delta}{R} \right) \right] . \end{align}

The viscosity does not appear explicitly on the right-hand side of this inequality, but there still is a viscosity dependence hidden in the definition of $\boldsymbol {u}$. To make this dependence apparent, we note that

(2.33)\begin{equation} |\boldsymbol{u}|^2 =|(\boldsymbol{v} - \boldsymbol{v}_P) - \boldsymbol{v}_b|^2 \leqslant (|\boldsymbol{v} - \boldsymbol{v}_P| + |\boldsymbol{v}_b|)^2 \leqslant 2 |\boldsymbol{v} - \boldsymbol{v}_P|^2 + 2 |\boldsymbol{v}_b|^2 \end{equation}

and

(2.34)\begin{equation} \int \frac{|\boldsymbol{u}|^2}{r^4} \,\mathrm{d}V \leqslant 2 \int \frac{|\boldsymbol{v} - \boldsymbol{v}_P|^2}{r^4} \,\mathrm{d}V + \frac{2}{R^4} \int |\boldsymbol{v}_b|^2 \,\mathrm{d}V \end{equation}

together with

(2.35)\begin{equation} \int |\boldsymbol{v}_b|^2 \,\mathrm{d}V = \frac{4{\rm \pi}}{5} R^2 U_\infty^2 \delta \left[ 1 + {O}\left( \frac{\delta^2}{R^2} \right) \right]. \end{equation}

If we insert these expressions into (2.32), we find that the term in $\int |\boldsymbol {v}_b|^2 \,\mathrm {d}V$ is smaller than the last term in (2.32) by a factor $\delta /R$ and is of the same order as previously neglected terms, so that this term has to be omitted for consistency as well. Let us simplify a numerical prefactor according to $327 < 24 {\rm \pi}(3+{\sqrt {2226}}/{35} ) < 328$ to finally arrive at

(2.36)\begin{equation} \nu \left\langle \int |\boldsymbol{\nabla} \times \boldsymbol{v}|^2 \,\mathrm{d}V \right\rangle \leqslant 12 \frac{U_\infty}{R} \int \left( \frac{R}{r} \right)^4 \langle |\boldsymbol{v} - \boldsymbol{v}_P|^2 \rangle \,\mathrm{d}V + 328 R^2 U_\infty^3 + {O}\left( \frac{\nu}{U_\infty R} \right) . \end{equation}

Relation (2.36) now has the characteristic features that Kolmogorov scaling demands. For any given velocity field, the right-hand side of (2.36) becomes asymptotically independent of viscosity as the Reynolds number tends to infinity, and it depends only on the velocity field itself, not on any of its derivatives. Deviations from the potential flow, such as the wake, contribute to the first term in the bound. This term weights the deviations such that the region close to the sphere contributes more than the far-away regions. The second term in the bound depends solely on $U_\infty$, which sets the amplitude of the irrotational flow. A purely irrotational flow does not cause any drag, a fact known as d'Alembert paradox. However, a purely irrotational flow is not compatible with the no-slip boundary conditions assumed here. Any non-zero $U_\infty$ thus implies a boundary layer or other rotational component, which explains the presence of the second term.

The Kolmogorov phenomenology originates from the study of homogeneous and isotropic turbulence. In this context, the dissipation per volume plays a central role and is proportional to $U^3/l$, where $U$ and $l$ are a characteristic velocity and length of the flow. The flow around a sphere is neither homogeneous nor isotropic, and its dissipation is finite and occurs in an infinite volume, so that the ratio of dissipation over flow volume is zero. Nonetheless, for the purpose of comparison with the Kolmogorov scaling, we may divide (2.36) by the reference volume $R^3$ to find that in the limit of $\nu /(U_\infty R)$ tending to zero,

(2.37)\begin{equation} \frac{\nu}{R^3} \left\langle \int |\boldsymbol{\nabla} \times \boldsymbol{v}|^2 \,\mathrm{d}V \right\rangle \leqslant 328 \frac{U^3}{R} \end{equation}

with

(2.38)\begin{equation} U = \left[ U_\infty^2 + \frac{3}{82}\frac{1}{R^3} \int \left( \frac{R}{r} \right)^4 \langle |\boldsymbol{v} - \boldsymbol{v}_P|^2 \rangle \,\mathrm{d}V \right]^{1/2} \end{equation}

as a possible choice for the characteristic velocity.

The empirical argument that the drag force should be proportional to $\rho U_\infty ^2 R$ requires, in addition, that $U$ is less than a bound proportional to $U_\infty$. The calculation presented here does not allow such a conclusion.

There are different ways in which the background field could be altered to produce relations similar to (2.36). The boundary-layer contribution $\boldsymbol {v}_b$ can obviously be optimised. However, there are also different possible choices for the irrotational part of the background field. For example, the potential flow around any obstacle enclosing the sphere is a viable option. Such a field could emulate the presence of a wake. If one aims for a short derivation, one might select the unperturbed velocity $U_\infty \hat {\boldsymbol {z}}$ as the irrotational part of the background. All the terms in (2.8) still allow for estimates of the required form if $\delta$ is fixed at a viscosity-independent value as for example $R$, so that $\boldsymbol {v}_b$ and $\boldsymbol {\nabla } \times \boldsymbol {v}_b$ have viscosity-independent supremum norms.

3. Conclusion

The relation (2.36) has the form expected from the Kolmogorov phenomenology. It bounds the first derivatives of the velocity field in the dissipation term in terms of the velocity field itself, and this bound asymptotically becomes independent of viscosity at high Reynolds numbers.

The prefactor in (2.36) is larger than what one can hope for from experimental data. Even if the velocity field outside the boundary layer was exactly irrotational, equation (2.36) would constrain the dissipation to $\rho \nu \langle \int |\boldsymbol {\nabla } \times \boldsymbol {v}|^2 \,\mathrm {d}V \rangle \leqslant 328 \rho R^2 U_\infty ^3$, whereas one finds experimentally that $\rho \nu \langle \int |\boldsymbol {\nabla } \times \boldsymbol {v}|^2 \,\mathrm {d}V \rangle \approx 0.6 \rho R^2 U_\infty ^3$ at $U_\infty R/\nu \approx 10^5$. This large discrepancy suggests that there are more optimal choices of the background field.

The result (2.36) is not enough to support the assumption mentioned in the introduction that the drag force should behave as $\rho U_\infty ^2 R^2$ at high Reynolds numbers. This scaling follows from (2.36) only under the additional assumption that a characteristic velocity such as (2.38) stays on the order of $U_\infty$ as $\nu$ tends to zero. We know from the work by Hopf that there is a bound on the dissipation that contains a viscosity dependence of the form $\exp (U_\infty R/\nu )$. The dissipation could become this large in the real flow either because the motion develops small length scales, or because the velocity becomes large. The inequality (2.36) excludes the former possibility. If indeed there exists a flow regime which exhausts Hopf's bound, it will have to do so by generating flows for which the integral $\int ({R}/{r} )^4 \langle |\boldsymbol {v} - \boldsymbol {v}_P|^2 \rangle \, \mathrm {d}V$ depends exponentially on $U_\infty R/\nu$. This could be realised by velocity amplitudes or the size of the wake growing exponentially with $U_\infty R/\nu$.