Flow constrained by concentric spheres

Posted in math on Sunday, May 11 2014

The last few flow problems I toyed with used a simple momentum balance as the starting point, time to move on to other ways to solve flow problems such as the continuity equation and the equations of motion for the fluid (e.g. the Navier-Stokes equations ).

Today I'm going to examine the (laminar) flow between two concentric spheres (Transport Phenomena 3B.4). Suppose the fluid is incompressible, Newtonian, and the flow is creeping (i.e. $ \mathbf{ v \centerdot \nabla v} \approx 0 $)

concentric spheres

First off the choice of coordinate system is the normal spherical coordinates centered at the center of the inner sphere. The entrance and exit are assumed to be the same angular size, such that the flow enters at θ = ε and exits at θ = π - ε.

Suppose that the velocity is entirely in the θ direction: $$ \begin{array}{rl} v_{\theta} = v_{\theta}(r,\theta) & v_{\phi} = v_{r} = 0 \end{array} $$

From the continuity equation: $$ \begin{array}{rl} 0 =& \mathbf{ \nabla \centerdot v} \ 0 =& \frac{1}{r^2} \frac{\partial}{\partial r} \left( r^2 v_r \right) + \frac{1}{r \sin \theta} \frac{\partial}{\partial \theta} \left( v_{\theta} \sin \theta \right) + \frac{1}{r \sin \theta} \frac{\partial}{\partial \phi} v_{\phi} \ 0 =& \frac{1}{r \sin \theta} \frac{\partial}{\partial \theta} \left( v_{\theta} \sin \theta \right) \ 0 =& \frac{\partial}{\partial \theta} \left( v_{\theta} \sin \theta \right) \end{array} $$

This last bit implies that: $$ v_{\theta} \sin \theta = u(r) = \mbox{constant w.r.t. }\theta $$

Where u(r) is some function of r that has yet to be determined.

Since this is an incompressible Newtonian fluid we can skip along to using the Navier-Stokes equations in the form: $$ \frac{\partial}{\partial t} \rho \mathbf{v} = -\rho \left( \mathbf{v \centerdot \nabla v} \right) - \nabla p + \mu \nabla^2 \mathbf{v} + \rho \mathbf{g} $$

Assuming steady-state creeping flow, and combining gravity and pressure into the modified pressure..

$$ 0 = - \nabla \mathcal{P} + \mu \nabla^2 \mathbf{v} $$

in the θ direction

$$ 0 = -\frac{1}{r} \frac{\partial \mathcal{P}}{\partial \theta} + \mu \left[ \frac{1}{r^2} \frac{\partial}{\partial r} \left( r^2 \frac{\partial}{\partial r} v_{\theta}\right) + \frac{1}{r^2} \frac{\partial}{\partial \theta} \left( \frac{1}{\sin \theta} \frac{\partial}{\partial \theta} v_{\theta} \sin \theta \right) \right] $$

Recalling that: $$ \begin{array}{rl} \frac{\partial}{\partial \theta} v_{\theta} \sin \theta = 0 & v_{\theta} = \frac{u(r)}{\sin \theta} \end{array} $$

We end up with $$ 0 = -\frac{1}{r} \frac{d \mathcal{P}}{d \theta} + \frac{\mu}{r^2 \sin \theta} \frac{d}{dr} \left( r^2 \frac{d}{dr} u(r) \right) $$

Separating out variables: $$ \sin \theta \frac{d \mathcal{P}}{d \theta} = \frac{\mu}{r} \frac{d}{dr} \left( r^2 \frac{d}{dr} u(r) \right) = B = \mbox{constant w.r.t }r \mbox{ and } \theta $$

Let's solve for B by looking at the θ side of the equation $$ \sin \theta \frac{d \mathcal{P}}{d \theta} = B $$

$$ \Delta \mathcal{P} = \int \limits_{\epsilon}^{\pi - \epsilon} {B \over {\sin \theta}} d \theta = B \int \limits_{\epsilon}^{\pi - \epsilon} \csc \theta d \theta $$

Taking advantage of the result: $$ \int \csc \theta d \theta = \ln \tan \frac{\theta}{2} $$

As well as some basic trig manipulations...

$$ \Delta \mathcal{P} = 2 B \ln \left( \cot \left( \frac{\epsilon}{2} \right) \right) $$


$$ B = { {\Delta \mathcal{P}} \over {2 \ln \left( \cot \left( \frac{\epsilon}{2} \right) \right)} } $$

Now we can solve for the r side of the equation, and substitute in for B at the end.

$$ \frac{\mu}{r} \frac{d}{dr} \left( r^2 \frac{d}{dr} u(r) \right) = B $$

integrating once $$ \mu \left( r^2 \frac{d}{dr} u(r) \right) = \frac{B r^2}{2} + C_1 $$

integrating again $$ 2 \mu u(r) = B r - \frac{2 C_1}{r} + C_2 $$

If we take u(r) = 0 at r = kR and r = R (no-slip boundary condition), and do a bunch of algebra to get rid of the two constants of integration, we get:

$$ u(r) = \frac{B R}{2 \mu} \left( \left( \frac{r}{R} -1 \right) + k \left( \frac{R}{r} - 1 \right) \right) $$

Substituting in for B: $$ u(r) = \frac{\Delta \mathcal{P} R}{4 \mu \ln \left( \cot \left( \frac{\epsilon}{2}\right) \right)} \left( \left( \frac{r}{R} -1 \right) + k \left( \frac{R}{r} - 1 \right) \right) $$

We can determine the mass flow-rate around the spheres by integrating the velocity over a plane of constant θ $$ w = \int \limits_{0}^{2\pi} \int \limits_{kR}^{R} \rho v_{\theta} \sin \theta r dr d\phi = 2\pi \rho \int \limits_{kR}^{R} u(r) r dr $$

Subbing in for u(r) and integrating:

$$ w = { {\pi \rho R \Delta \mathcal{P}} \over {2 \mu \ln \left( \cot \left( \frac{\epsilon}{2}\right) \right)} } \left[ \frac{r^3}{3R} - \frac{1+k}{2} r^2 + kRr \right]_{kR}^{R} $$

After some manipulations, especially noting that k3 - 3k2 + 3k -1 = -(1-k)3, we end up at: $$ w = - { {\pi \rho R \Delta \mathcal{P} \left( 1-k\right)^3} \over {12 \mu \ln \left( \cot \left( \frac{\epsilon}{2}\right) \right)} } $$