Return to Answer

This post is a work in progress, I will try to add the missing formulas and images later today.

In order to evaluate the resulting force we have to sum up all the tiny forces per surface element, the stresses, by means of integration (I have created an animated version of this figure that you can access here).

Our goal is now to find an analytical description for the drag. As we can imagine for highly turbulent (high Reynolds number) flow this is almost impossible to obtain as the flow will be chaotic, with no true steady state and instead depend on forming eddies. On the other side, for very low Reynolds numbers where the flow is smooth and follows the shape without separation we might be more lucky. In this flow regime we can use the symmetry of the flow and the corresponding forces it suffices to integrate the forces over tiny cylinders surfaces. This can be done by looking at a cylindrical segment whose area is given by a circle with radius $R*\sin(\theta)$ and a height equal to the differential arc length $R d\theta$ (here a YouTube video about it and see this post on Archimedes' Hat-Box Theorem).

$$ F_D = - \int\limits_{0}^{\pi} \tau_{r \theta} |_{r=R} \sin(\theta) dA + \int\limits_{0}^{\pi} \tau_{r r} |_{r=R} \cos(\theta) dA $$

$$ F_D = - 2 \pi R^2 \int\limits_{0}^{\pi} \tau_{r \theta} |_{r=R} \sin^2(\theta) d\theta + 2 \pi R^2 \int\limits_{0}^{\pi} \tau_{r r} |_{r=R} \sin(\theta)\cos(\theta) d\theta $$

(Approach, formulas and results will follow, they are quite complicated anyways)

The pressure in front is clearly the highest similar to the impulse you feel when sticking the hand out of a driving vehicle. On the other side, on the back, it is significantly lower.

With these distributions for pressure and velocity the drag force can be calculated to

$$F_D = 4 \pi \mu U R + 2 \pi \mu U R = 6 \pi \mu U R = 3 \pi \mu U D.$$

As can be seen this indeed takes the form that we predicted by smart dimensional analysis.

In order to evaluate the resulting force we have to sum up all the tiny forces per surface element, the stresses, by means of integration (I have created an animated version of this figure that you can access here).

Our goal is now to find an analytical description for the drag. As we can imagine for highly turbulent (high Reynolds number) flow this is almost impossible to obtain as the flow will be chaotic, with no true steady state and instead depend on forming eddies. On the other side, for very low Reynolds numbers where the flow is smooth and follows the shape without separation we might be more lucky. In this flow regime we can use the rotational symmetry of the flow ($\Phi$ in this case is an angle around the horizontal rotation axis) $$ u_{\Phi}=0, \quad \frac{\partial}{\partial \Phi}=0 $$ and the corresponding forces it suffices to integrate the forces over tiny cylinders surfaces. This can be done by looking at a cylindrical segment whose area is given by a circle with radius $R*\sin(\theta)$ and a height equal to the differential arc length $R d\theta$ (here a YouTube video about it and see this post on Archimedes' Hat-Box Theorem).

$$ F_D = - \int\limits_{0}^{\pi} \sigma_{r \theta} |_{r=R} \sin(\theta) dA + \int\limits_{0}^{\pi} \sigma_{r r} |_{r=R} \cos(\theta) dA $$

$$ F_D = - 2 \pi R^2 \int\limits_{0}^{\pi} \sigma_{r \theta} |_{r=R} \sin^2(\theta) d\theta + 2 \pi R^2 \int\limits_{0}^{\pi} \sigma_{r r} |_{r=R} \sin(\theta)\cos(\theta) d\theta $$

Vector identities and rotational symmetry

We can use the vector identity

$$\vec \nabla ^2 \vec u = \vec \nabla ( \underbrace{ \vec \nabla \cdot \vec u }_{= 0}) - \vec \nabla \times ( \vec \nabla \times \vec u)$$

to simplify the Stokes' momentum equation to

$$\vec \nabla p = - \mu \vec \nabla \times ( \vec \nabla \times \vec u).$$

Now we can take the cross-product of this equation and apply the vector identity

$$\vec \nabla \times \vec \nabla p = 0.$$

in order to eliminate the pressure and obtain the linear equation

$$\vec \nabla \times [ \vec \nabla \times ( \vec \nabla \times \vec u) ].$$

Further we will now switch to spherical coordinates as it is more convenient for a rotational symmetric problem

$$\vec u = \left(\begin{array}{c} u_r \\ u_\Theta \\ u_\Phi \end{array}\right)$$

The corresponding operators take the following form

$$\require{cancel} \vec \nabla \cdot \vec u =\frac{1}{r^{2}} \frac{\partial}{\partial r}\left(r^{2} u_{r}\right)+\frac{1}{r \sin \Theta} \frac{\partial}{\partial \Theta}\left(u_{\Theta} \sin \Theta\right)+\cancel{\frac{1}{r \sin \Theta} \frac{\partial u_{\Phi}}{\partial \Phi}},$$

$$\require{cancel} \vec \nabla \times \vec u =\left(\begin{array}{c} \cancel{\frac{1}{r \sin \Theta} \frac{\partial}{\partial \Theta}\left(u_{\Phi} \sin \Theta\right)}-\cancel{\frac{1}{r \sin \Theta} \frac{\partial u_{\Theta}}{\partial \Phi}} \\ \cancel{\frac{1}{r \sin \Theta} \frac{\partial u_{r}}{\partial \Phi}}-\cancel{\frac{1}{r} \frac{\partial}{\partial r}\left(r u_{\Phi}\right)} \\ \frac{1}{r} \frac{\partial}{\partial r}\left(r u_{\Theta}\right)-\frac{1}{r} \frac{\partial u_{r}}{\partial \Theta} \end{array}\right).$$

Stokes' stream function and solving the PDE

The derivation from here on can be found in a similar way in this NYU lecture notes, section 7.3. We can now introduce Stokes' stream function for rotational symmetric bodies

$$\frac{\partial \Psi}{\partial \Theta}:=u_{r} r^{2} \sin \Theta, \quad \frac{\partial \Psi}{\partial r}:=-u_{\Theta} r \sin \Theta.$$

This allows us to rewrite the Stokes' momentum equation as follows:

$$\vec \nabla \times \vec u =\left(\begin{array}{c} 0 \\ 0 \\ \frac{1}{r} \frac{\partial}{\partial r}\left(r u_{\Theta}\right)-\frac{1}{r} \frac{\partial u_{r}}{\partial \Theta} \end{array}\right)= \left(\begin{array}{c} 0 \\ 0 \\ -\frac{1}{r sin \Theta} \underbrace{ \left[ \frac{\partial^2 \Psi}{\partial r^2} + \frac{\sin \Theta}{r^2} \frac{\partial}{\partial\Theta} \left( \frac{1}{\sin \Theta} \frac{\partial \Psi}{\partial \Theta} \right) \right] }_{:= \mathcal{L}(\Psi)} \end{array}\right),$$

$$\vec \nabla \times(\vec \nabla \times \vec u)=\left(\begin{array}{c} - \frac{1}{r^2 \sin \Theta} \frac{\partial [\mathcal{L}(\Psi)]}{\partial \Theta} \\ \frac{1}{r \sin \Theta} \frac{\partial [ \mathcal{L}(\Psi)]}{\partial r} \\ 0 \end{array}\right),$$

$$\vec \nabla \times[\vec \nabla \times(\vec \nabla \times \vec u)]=\left(\begin{array}{c} 0 \\ 0 \\ \frac{1}{r \sin \theta} \underbrace{ \left\{\frac{\partial^{2}(\mathcal{L} \Psi)}{\partial r^{2}}+\frac{\sin \theta}{r^{2}} \frac{\partial}{\partial \Theta}\left[\frac{1}{\sin \theta} \frac{\partial(\mathcal{L} \Psi)}{\partial \Theta}\right]\right\} }_{= \mathcal{L}[\mathcal{L}(\Psi)]} \\ \end{array}\right).$$

This can also be written as

$$\vec \nabla \times[\vec \nabla \times(\vec \nabla \times \vec u)] =\frac{1}{r \sin \Theta} \mathcal{L}[\mathcal{L} (\Psi)] \vec e_{\Phi} =\frac{1}{r \sin \Theta} \mathcal{L}^{2} \Psi \vec e_{\Phi}$$

where the operator $\mathcal{L}$ is given as

$$\mathcal{L}=\frac{\partial^{2}}{\partial r^{2}}+\frac{\sin \theta}{r^{2}} \frac{\partial}{\partial \Theta}\left(\frac{1}{\sin \Theta} \frac{\partial}{\partial \Theta}\right).$$

Additionally the flow has to fulfill the boundary conditions. The fluid velocity at the wall must be zero (no-slip condition)

$$u_r = 0: \frac{\partial \Psi}{\partial \Theta} = 0, \quad u_\Theta = 0: \frac{\partial \Psi}{\partial r} = 0$$

and in the far-field the velocity is given by the unperturbed velocity $U_\infty$

$$u_{r} =\frac{1}{r^{2} \sin \Theta} \frac{\partial \Psi}{\partial \Theta}=U_\infty \cos \Theta,$$

$$u_{\Theta} =-\frac{1}{r \sin \Theta} \frac{\partial \Psi}{\partial r}=-U_\infty \sin \Theta.$$

Solving each of the equations for $\Psi$ by means of integration we find

$$\Psi=U_\infty \frac{r^{2}}{2} \sin ^{2} \Theta+C \quad \forall \Theta, r \rightarrow \infty$$

for the solution indefinitely far from the cylinder. Thus, the resulting partial differential equation could be solved by a product ansatz of the form

$$\Psi(r, \Theta)=f(r) \sin ^{2} \Theta.$$

We insert this ansatz into the partial differential equation $$\mathcal{L} (\Psi)=\underbrace{\left(f^{\prime \prime}-\frac{2}{r^{2}} f\right)}_{:= F(r)} \sin ^{2} \Theta,$$

$$\mathcal{L}^{2} (\Psi)=\mathcal{L}[\mathcal{L} (\Psi)]=\left(F^{\prime \prime}-\frac{2}{r^{2}} F\right) \sin ^{2} \Theta,$$

$$\mathcal{L}^{2} \Psi=0 \Leftrightarrow(\underbrace{F^{\prime \prime}-\frac{2}{r^{2}} F}_{= 0}) \underbrace{\sin ^{2} \Theta}_{=\neq 0}=0.$$

As the last term can't be assumed to vanish the differential equation that must be fulfilled is the simple Eulerian differential equation

$$F^{\prime \prime}-\frac{2}{r^{2}} F=0$$

that can be solved with the ansatz $F = C r^\lambda$ resulting in the algebraic equation

$$\cancel{r^2} \lambda (\lambda -1) \cancel{Cr^{\lambda-2}} - 2 \cancel{Cr^\lambda} = 0,$$

$$\lambda (\lambda - 1) - 2 = (\lambda - 2) (\lambda + 1) = 0,$$

which can be solved in a similar manner additionally considering the particular solution to

$$F = f^{\prime \prime}-\frac{2}{r^{2}} = A r^2 + \frac{B}{r}$$

which again results in $$f(r)=\frac{A}{10} r^{4}-\frac{B}{2} r+C r^{2}+\frac{D}{r}.$$

Thus we find for the stream function $\Psi$ and the radial $u_r$ and tangential $u_\Theta$ velocities

$$\Psi=\frac{1}{4} U_\infty R^{2}\left(\frac{R}{r}-3 \frac{r}{R}+2 \frac{r^{2}}{R^{2}}\right) \sin ^{2} \Theta,$$

$$u_{r}=U_\infty \left(1+\frac{1}{2} \frac{R^{3}}{r^{3}}-\frac{3}{2} \frac{R}{r}\right) \cos \Theta,$$

$$u_{\Theta}=U_\infty \left(-1+\frac{1}{4} \frac{R^{3}}{r^{3}}+\frac{3}{4} \frac{R}{r}\right) \sin \Theta.$$

Finally we can determine the pressure by integrating the radial Stokes' momentum equation

$$\frac{\partial p}{\partial r} =\mu \frac{1}{r^{2} \sin \Theta} \frac{\partial[\mathcal{L} (\Psi)]}{\partial \Theta},$$

to

$$p =-\frac{3}{2} \frac{\mu U_{\infty} R}{r^{2}} \cos \theta+D$$

Considering the pressure in the far field $r \to \infty$ we finally yield

$$p = p_{\infty}-\frac{3}{2} \frac{\mu U_{\infty} R}{r^{2}} \cos \Theta.$$

The pressure in front is clearly the highest similar to the impulse you feel when sticking the hand out of a driving vehicle. Behind the sphere it is anti-symmetrically lower.

Integrating the force

With these distributions for pressure and velocity finally the stresses and therefore the drag force can be evaluated by integration

$$\left. \sigma_{rr} \right|_{r=R} = - \left. p \right|_{r=R} + \left. \underbrace{ 2 \mu \frac{\partial u_r}{\partial r} }_{\tau_{rr}} \right|_{r=R} = - p_{\infty} + \frac{3}{2} \frac{\mu U_{\infty}}{R} \cos \Theta,$$

$$\left. \sigma_{r \Theta} \right|_{r=R} = \left. \underbrace{ \mu \left( \frac{1}{r} \frac{\partial u_r}{\partial \Theta} + \frac{\partial u_\Theta}{\partial r} - \frac{u_\Theta}{r} \right)}_{\tau_{r \Theta}} \right|_{r=R} = - \frac{3}{2} \frac{\mu U_\infty}{R} \sin \Theta,$$

$$F_D = 3 \pi \mu U_{\infty} R \int\limits_{\Theta = 0}^{\pi} \underbrace{(\sin^2 \Theta + \cos^2 \Theta)}_{=1} \sin \Theta d\Theta = 3 \pi \mu U_{\infty} R \left. cos \overline \Theta \right|_{\overline \Theta = 0}^{\pi} = 6 \pi \mu U R = 3 \pi \mu U D$$

As can be seen this indeed takes the form that we predicted by smart dimensional analysis. Evaluating each contribution independently we find

$$F_D = \underbrace{ 4 \pi \mu U R}_\text{viscous contribution} + \underbrace{2 \pi \mu U R}_\text{pressure contribution}.$$

$$ Re := \frac{ U L }{ \nu } = \frac{\rho \, U \, L}{\mu} = \overbrace{\rho A U^2}^\text{inertia} \cdot \underbrace{\frac{L}{\mu A U}}_\text{viscous force} \propto \frac{\text{inertial forces}}{\text{viscous damping forces}}.$$

Dimensional analysis

where $A$ is the projected area with units $[A] = m^2$. Thus, we will follow a more formal way of deriving the drag force that is valid for any flow regime of continuum flows, the integration of surface stresses.

In order to evaluate the resulting force we have to sum up all the tiny forces per surface element, the stresses, by means of integration (I have an animated figure prepared here).

Our goal is now to find an analytical description for the drag. As we can imagine for highly turbulent (high Reynolds number) flow this is almost impossible to obtain as the flow will be chaotic, with no true steady state and instead depend on forming eddies. On the other side, for very low Reynolds numbers where the flow is smooth and follows the shape without separation we might be more lucky. In this flow regime we can use the symmetry of the flow and the corresponding forces it suffices to integrate the forces over tiny cylinders surfaces. This can be done by looking at a cylindrical segment whose area is given by a circle with radius $R*\sin(\theta)$ and a height equal to the differential arc length $R d\theta$ (here a youtube video about it and see this post on Archimedes' Hat-Box Theorem).

$$ Re := \frac{ U L }{ \nu } = \frac{\rho \, U \, L}{\mu} = \overbrace{\rho A U^2}^\text{inertia} \cdot \underbrace{\frac{L}{\mu A U}}_\text{viscous forces} \propto \frac{\text{inertial forces}}{\text{viscous damping forces}}.$$

Educated guess and dimensional analysis

where $A$ is the projected area with units $[A] = m^2$. While we could try to conduct a set of controlled experiments to verify our claims, we will follow a more formal way of deriving the drag force that is valid for any flow regime of continuum flows, the integration of surface stresses.

In order to evaluate the resulting force we have to sum up all the tiny forces per surface element, the stresses, by means of integration (I have created an animated version of this figure that you can access here).

Our goal is now to find an analytical description for the drag. As we can imagine for highly turbulent (high Reynolds number) flow this is almost impossible to obtain as the flow will be chaotic, with no true steady state and instead depend on forming eddies. On the other side, for very low Reynolds numbers where the flow is smooth and follows the shape without separation we might be more lucky. In this flow regime we can use the symmetry of the flow and the corresponding forces it suffices to integrate the forces over tiny cylinders surfaces. This can be done by looking at a cylindrical segment whose area is given by a circle with radius $R*\sin(\theta)$ and a height equal to the differential arc length $R d\theta$ (here a YouTube video about it and see this post on Archimedes' Hat-Box Theorem).

Forces on structures: Integration of stresses

$$dA = \underbrace{2 R \sin(\theta) \pi}_\text{circle} \underbrace{R d\theta}_\text{height}$$

$$m \vec a = F_D$$

If this was a constant value this equation could be solved very easily by double integration. Generally the drag for though depends on the velocity $F_D = - k |\vec u| \vec u$ though and thus acceleration $\vec a = \frac{d \vec u}{dt}$ and velocity $u$ are coupled which makes them very complicated to solve analytically. For this reason one often turns to numerical solutions.

$$\sigma_{ij} = - p \delta_{ij} + 2 \mu S_{ij} - \frac{2}{3} \mu \sum\limits_{k \in \mathcal{D}} S_{kk} \delta_{ij}.$$

Stokes’ equations: The PDEs of creepying flows ($Re \ll 1$)

Now our goal is to reduce complexity by checking if we can neglect some terms in these equations as others are more important and will dominate the flow for very low Reynolds number anyways. For very low Reynolds number $Re \ll 1$ clearly the viscosity will dominate over the inertia. In order to determine precisely which terms can be neglected we convert the equations to their dimensionless form.

This simpler form allows us now to find closed descriptions for the pressure and the velocity. The pressure in front is clearly the highest similar to the impulse you feel when sticking the hand out of a driving vehicle. On the other side, on the back, it is significantly lower.

Dimensional analysis

The first thing we can attempt to do is trying to make an educated guess what the relationship between drag force and the different parameters might look like using dimensional analysis. The unit of force is Newton

$$ [F_D] = N = \frac{kg \, m}{s^2}.$$

A correct relationship must yield the same units. We might assume that parameters influencing the outcome will be the diameter of the sphere $[D] = m$, the velocity $[U] = \frac{m}{s}$, the viscosity $[\nu] = \frac{m^2}{s}$ and the density $[\rho] = \frac{kg}{m^3}$. Larger values of all of the variables states above should result in higher drag force. Thus it seems plausible that such a relationship might take the form

$$ F_D \propto \nu \rho U D = \mu U D.$$

Yet we can already see that this is only one potential correlation. For higher Reynolds numbers actually the correlation looks like

$$ F_D \propto \rho U^2 A $$

where $A$ is the projected area with units $[A] = m^2$. Thus, we will follow a more formal way of deriving the drag force that is valid for any flow regime of continuum flows, the integration of surface stresses.

Forces on structures: Integration of stresses

$$dA = \underbrace{2 R \sin(\theta) \pi}_\text{circle} \overbrace{R d\theta}^\text{height}$$

$$m \vec a = \vec F_D$$

If this force $\vec F_D$ was a constant value this equation could be solved very easily by double integration. Generally the drag for though depends on the velocity $\vec F_D = - k |\vec u| \vec u$ though and thus acceleration $\vec a = \frac{d \vec u}{dt}$ and velocity $u$ are coupled which makes them very complicated to solve analytically. For this reason one often turns to numerical solutions.

$$\sigma_{ij} = - p \delta_{ij} + 2 \mu S_{ij} - \frac{2}{3} \mu \sum\limits_{k \in \mathcal{D}} S_{kk} \delta_{ij}. \tag{4}\label{4}$$

Incompressible Navier-Stokes-equations

A first simplification is the assumption of an incompressible fluid which basically implies that density is a constant throughout the flow field. We furthermore silently assume that body forces can be neglected and the temperature does not change significantly throughout the flow field. As a consequence we only need the continuity and momentum equation which degenerate to:

$$\sum\limits_{j \in \mathcal{D}} \frac{\partial u_j}{\partial x_j} = 0,$$

$$\frac{\partial u_i}{\partial t} + \sum\limits_{j \in \mathcal{D}} u_j \frac{\partial u_i}{\partial x_j} = - \frac{1}{\rho} \frac{\partial p}{\partial x_i } + \frac{1}{\rho} \sum\limits_{j \in \mathcal{D}} \frac{\partial \tau_{ij}}{\partial x_j }.$$

and furthermore the last term in the stress tensor vanishes in equation \eqref{4} (for more details see).

Stokes’ equations: The PDEs of creepying flows ($Re \ll 1$)

Now our goal is to reduce complexity by checking if we can neglect some terms in these equations as others are more important and will dominate the flow for very low Reynolds number anyways. For very low Reynolds number $Re \ll 1$ clearly the viscosity will dominate over the inertia. In order to determine precisely which terms can be neglected we convert the equations to their dimensionless form by dividing through characteristic measures as follows:

$$x_i^*=\frac{x_i}{L}, \phantom{abc} u_i^*=\frac{u_i}{U}, \phantom{abc} \rho^*=\frac{\rho}{\rho_0}, \phantom{abc} t^*=\frac{t}{\frac{L}{U}},\phantom{abc} p^*=\frac{p}{\rho_0 \nu \frac{U}{L}} \phantom{ab}$$

This results in the two dimensionless equations

$$\sum\limits_{j \in \mathcal{D}} \frac{\partial u_j^*}{\partial x_j^* }=0$$

$$\underbrace{\overbrace{Re}^{\ll 1} \left( \frac{\partial u_i^*}{\partial t^*} + \sum\limits_{j \in \mathcal{D}} u_j^* \frac{\partial u_i^*}{\partial x_j^*} \right)}_{\ll 1} = - \frac{\partial p^*}{ \partial x_i^* } + \sum\limits_{j \in \mathcal{D}} \frac{\partial \tau_{ij}^*}{\partial x_j^* }$$

where the last one can be rewritten with equation \eqref{4} and neglecting the terms on the left-hand side due to their magnitude to

$$\frac{\partial p}{\partial x_i} = \mu \sum\limits_{j \in \mathcal{D}} \frac{\partial^2 u_i}{\partial x_j^2}.$$

This set of equations is often termed Stokes' equations and can be written symbolically as

$$\vec \nabla \cdot \vec u = 0,$$

$$\vec \nabla p = \mu \vec \nabla^2 \vec u = \mu \vec \Delta \vec u.$$

Now it is possible to find a closed descriptions for the pressure and the velocity by solving this much simpler system of partial differential equations.

(Approach, formulas and results will follow, they are quite complicated anyways)

The pressure in front is clearly the highest similar to the impulse you feel when sticking the hand out of a driving vehicle. On the other side, on the back, it is significantly lower.

With these distributions for pressure and velocity the drag force can be calculated to

$$F_D = 4 \pi \mu U R + 2 \pi \mu U R = 6 \pi \mu U R = 3 \pi \mu U D.$$

As can be seen this indeed takes the form that we predicted by smart dimensional analysis.