\input{aefm-macro}

\def\f{\frac}
\def\erf{\rm erf}
\def\p{\partial}
\def\be{\begin{equation}}
\def\ee{\end{equation}}
\def\a{\alpha}
\def\n{\nabla}
\def\z{\zeta}
\def\la{\langle}
\def\ra{\rangle}
\def\lb{\left[}
\def\rb{\right]}
\def\lcb{\left\{}
\def\rcb{\right\}}
\def\lp{\left(}
\def\rp{\right)}
\def\rhu{\rightharpoonup}
\def\o{\over}
\def\ep{\epsilon}
\def\gl{\stackrel{>}{\scriptstyle <}}
\def\lg{\stackrel{<}{\scriptstyle >}}
\def\lap{\nabla^2}
\def\no{\noindent}
\def\h{\hat{}}
\def\hs{\hspace{0.3 cm}}
\def\definition{\stackrel{\rm def}{=}}
\def\erf{{\rm erf}}
\def\erfc{{\rm erfc}}
\def\pvi{{\int\mskip -33 mu - \quad}}
\def\Real{\mbox{Re}\,}
\def\Imag{\mbox{Im}\,}
\def\pvint{{\int\!\!\!\!\!\!-}}
\def\2int{\int\!\!\!\int}
\def\3int{\int\!\!\!\int\!\!\!\int}
 \def\sgn{\mbox{sgn}}
\def\bfq{{\bf q}}
\def\bfn{{\bf n}}
\def\bfx{{\bf x}}
\def\o{\over}
\def\cd{{\cal D}}
2-5Stokes.tex 

%%%%%%%%%%%%%%%%%%%%%%%%%%
  \setcounter{chapter}{2}
\setcounter{section}{4}
 
 
 
\section{Stokes flow past a sphere}
\setcounter{equation}{0}

[Refs] \\
Lamb: {\em Hydrodynamics}\\
Acheson : {\em Elementary Fluid Dynamics}, p. 223 ff\\

 One of the fundamental results in low   Reynolds hydrodynamics is the  Stokes solution  for steady flow  past    a small sphere. The apllicatiuon range widely form the determination of electron charges to the physics of aerosols. 
 
 The continuity equation reads
\be \nabla\cdot\vec q = 0 \label{eqno(1)}\ee
With inertia neglected, the approximate momentum  equation is 
\be 0=-\f{\nabla p}{ \rho} + \nu\nabla^2\vec
q\label{eqno(2)}\ee
Physically, the presssure gradient drives the flow by overcoming viscous resistence, but does affect the fluid inertia significantly. 

Refering to Figure \ref{fig:sphericalcoord} for the spherical coordinate system $(r,\theta,\phi)$. Let the 
ambient velocity  be upward  and along the polar $(z)$ axis: $(u,v,w) = (0,0,W)$. 
Axial symmetry demands
\[ {\p\o \p \phi}=0, \quad \mbox{and}\quad \vec q = (q_r(r,\theta),
q_\theta(r,\theta), 0)\] 
 Using a known formula for the divergence in spherical polar coordinates, Eq. (\ref{eqno(1)}) becomes
\be 
\f{1}{ r^2} \f{\p}{\p r }(r^2 q_r) +\f{1}{r\sin\theta}{\p\o\p \theta}\left(q_\theta\sin
\theta\right)=0\label{eqno(7)}\ee
An equivalent and  physically more revealing  way is to write
\be 
  \f{\p}{\p r }(r^2 q_r \sin \theta) +\f{\p}{\p \theta}\left(r q_\theta\sin
\theta\right)=0\label{eqno(7)}\ee

As in the case of rectangular coordinates, we define  the stream function $\psi$ to satisify the continuity equation (\ref{eqno(7)}) identically
\be 
q_r = {1\o r^2\sin \theta}{\p\psi\o\p\theta},\quad
q_\theta=-{1\o r\sin\theta}{\p\psi\o\p r} \label{eqno(8)}\ee
   
 
At infinity, the uniform velocity $W$ along $z$ axis can be decomposed into radial and polar components 
\be q_r = W\cos\theta ={1\o r^2\sin\theta}{\p\psi\o \p \theta}, \quad
  q_\theta= -W\sin\theta 
  =-{1\o r\sin\theta} {\p\psi\o\p r}, \quad r \sim \infty\ee
The corresponding stream function at infinity  follows by integration 
\be \psi = {W\o 2}r^2\sin^2\theta, \quad r \sim \infty\label{eqno(9)}\ee

 Using the vector identity  
\be \nabla \times (\nabla\times \vec q) =
\nabla(\nabla\cdot\vec q) -\nabla^2 \vec
q\label{eqno(5)}\ee
and (\ref{eqno(1)}), we get
\be 
\nabla^2\vec q = -\nabla \times(\nabla \times\vec q) = -\nabla
\times \vec  
\zeta\label{eqno(4)}\ee
Taking the curl  of (\ref{eqno(2)}) and using (\ref{eqno(4)}) we get
\be 
\nabla \times (\nabla \times \vec \zeta)=0\label{eqno(5)}\ee
%%%%%%%%%
\begin{figure}
\label{fig:sphericalcoord}
\begin{center}
\includegraphics[scale=1]{sphcoor.eps}\end{center}
 \caption{The spherical coordinates} \end{figure}
%%%%%%%%%%%



After some  straightforward algebra  given in the Appendix, we can show   that
\be 
\vec q =\nabla \times \left({\psi \vec e_\phi \o r\sin\theta}
\right)\label{eqno(10)}\ee
and \be \vec \zeta = \nabla \times \vec q =\nabla\times \nabla \times \left({\psi \vec e_\phi \o r\sin\theta}
\right)=
 -\f{\vec e_\phi }{r\sin \theta}\lp  \f{\p^2\psi}{\p
r^2}+\f{\sin\theta}{ r^2}\f{\p  }{\p\theta}\left({1\o\sin\theta}
\f{\p \psi  }{ \p\theta}\rp\rp\label{eqno(11)}\ee

Now from (\ref{eqno(5)})
\[
\nabla \times \nabla \times(\nabla\times \vec q )=
\nabla\times \nabla\times\lb\nabla\times\lp\nabla\times\f{\psi \vec
e_\phi}{r\sin\theta} \rp\rb   = 0\]
hence, the momentum equation (\ref{eqno(5)}) becomes a scalar equation for $\psi$.
\be \lp  \f{\p^2 }{\p
r^2}+\f{\sin\theta}{ r^2}\f{\p  }{\p\theta}\left({1\o\sin\theta}
\f{\p    }{ \p\theta}\rp\rp^2\psi=0 
\label{eqno(13)}\ee
The boundary conditions on the sphere are 
\be 
q_r = 0\quad q_\theta= 0\quad \hbox{on}\quad r=a\label{eqno(14)}\ee
The boundary conditions at $\infty$ is
\be 
\psi\to {W\o 2} r^2\sin^{2}\theta\label{eqno(15)}\ee

Let us try a solution of the form: 
\be \psi(r, \theta)=f (r) \sin^2 \theta\label{eqno(16)}\ee
then $f$ is governed by the equi-dimensional differential equation: 
\be 
\left[\f{d^2}{ dr^2}-\f{2}{ r^2}\right]^2 f = 0\ee
whose solutions are of the form $f(r)\propto r^n$, 
It is easy to verify that $n=-1,1,2,4$ so that
\[ f(r) = \f{A}{ r} + Br+ Cr^2 + Dr^4  \]
or \[\psi = \sin^2\theta\left[\f{A}{ r} + Br + Cr^2 +
Dr^4\right]\]

To satisfy (\ref{eqno(15)}) we set $ D=0, C={W/ 2}$.
To satisfy (\ref{eqno(14)})  we use (\ref{eqno(8)}) to get 
\[
q_r = 0=\f{W}{ 2} + \f{A}{ a^3} + \f{B}{ a} = 0, \quad 
 q_\theta=0=W-\f{A}{ a^3} + \f{B}{ a}=0 \]
Hence
\[
A=\f{1}{ 4}Wa^3, \qquad B = -\f{3}{ 4} Wa\]
Finally the stream function  is
\be \psi=\f{W}{2} \left[ { r^2}  + 
  \f{a^3}{ 2r}  -  \f{3ar }{ 2 }\right]\sin^2 \theta  \label{eqno(17)}\ee
Inside the  parentheses, the  first term corresponds to the uniform
flow, and the second term to the doublet; together they represent an
inviscid flow past a sphere. The third term is called the Stokeslet,
representing the viscous correction.

 The 
velocity components   in the fluid are: (cf. (\ref{eqno(8)}) :
\begin{eqnarray} 
 q_r  &=& W\cos\theta\left[1 +  \f{a^3}{ 2r^3} -
  \f{3a}{ 2r}\right]\\
 q_\theta &=& -W\sin\theta \left[1 -   \f{a^3}{ 4 r^3} -
\f {3a}{ 4 r}\right]  \label{eqno(18)}\end{eqnarray}

\subsection{Physical Deductions}
\begin{enumerate}
\item {\em Streamlines}:   With respect to the  the equator along $\theta
= {\pi/ 2}$,  $\cos \theta$  and $q_r$  are  odd  while    
$\sin\theta$ and  $q_\theta$ are even.  Hence the 
streamlines (velocity vectors) are symmetric fore and aft.
\item  {\em Vorticity}: 
\[ \vec \zeta=  \zeta_\phi \vec
e_\phi=
\left(\f{1}{ r} \f{\p (r q_\theta
)}{\p r}- \f{1}{ r} \f{\p q_r}{ \p\theta}\right) \vec
e_\phi
  = -\f{3}{ 2} Wa \f{\sin\theta}{ r^2 }\vec e_\phi \]

 \item  {\em Pressure }:  From the $r$-component of
momentum equation
\[ \f{\p p}{ \p r} = \f{\mu  Wa}{ r^3} \cos\theta  (=-\mu \nabla
\times(\nabla
\times \vec q )) \]
Integrating  with respect to $r$ from $r$ to $\infty$, we get 
\be
p = p_\infty - \f{3}{ 2}\f{\mu Wa}{ r^3}\cos \theta\label{eqno(20)}\ee

 \item  {\em Stresses and strains}:
\[
\f{1}{ 2}e_{rr} = \f{\p q_r}{ \p r} = W\cos\theta\left( \f{3a}{ 2r^2}- \f{3a^3}{
2r^4}\right)\]
 On the sphere,  $r = a$,  $e_{rr}
= 0$ hence $ \tau_{rr} =0$ and 
\be 
\sigma_{rr} = - p + \tau_{rr} =-p_\infty + {3\o
2}{\mu W\o a} \cos\theta\label{eqno(21)}\ee
On the other hand 
\[ e_{r\theta}= r\f{\p }{\p r} \left(\f{ q_\theta}{ 
r}\right) +\f{1}{ r}\f{\p q_r}{ \p \theta} = -\f{3}{ 2}
\f{Wa^3}{ r^4}\sin\theta \]
Hence at $r = a$:
\be
\sigma_{r\theta}  =\tau_{r\theta}= \mu e_{r\theta} =
-\f{3}{ 2}\f{\mu W}{ a} \sin \theta\label{eqno(22)}\ee
The resultant stress on the sphere is parallel to the $z$ axis. 
 \[
\Sigma_z = \sigma_{rr} \cos\theta
-\tau_{r\theta}\sin\theta = -
p_\infty \cos\theta + \f{3}{ 2} \f{\mu W}{ a}\]
 

The constant part exerts a   net drag in $z$
direction  
\be 
D= \int^{2\pi}_o a d \phi\int^{\pi}_o d\theta \sin
\theta \Sigma_z   =
=  \f{3}{ 2} \f{\mu W}{ a} 4\pi a^2 = 6\pi\mu
Wa\label{eqno(23)}\ee
This is the celebrated Stokes formula. 

   A drag coefficient can be defined as
\be 
C_D= \f{D}{ \f{1}{ 2}\rho W^2   \pi a^2} = \f{6 \pi\mu Wa}{ {1\o 2}
\rho W^2 \pi a^2}=
\f{24}{\f{\rho W(2a)}{ \mu}} =  \f{24}{ Re_d}\label{eqno(24)}\ee

\item  {\em Fall velocity } of a particle  through a
fluid. Equating the drag and the   buoyant weight of the eparticle
\[
6\pi\mu W_o a = \f{4\pi}{ 3}a^3(\rho_s-\rho_f)g\] 
hence 
\[ W_o = \f{2}{ 9} g \left(\f{a^2}{\nu}\f{\Delta\rho}{
\rho_f}\right) = 217.8 \left(\f{a^2}{ \nu}\f{\Delta\rho}{
\rho_f}\right) \]
in cgs units. For  a  sand grain in water,
\[
{\Delta\rho\o\rho_f}= \f{2.5-1}{ 1} = 1.5, \quad 
\nu = 10^{-2} {\rm cm}^2/s\]
\be 
W_o = 32,670\  a^2 {\rm cm}/s\ee
To have some quantitative ideas, let us consider two sand of two sizes  : \begin{eqnarray*}&&  a = 10^{-2} {\rm cm} = 10^{-4} {\rm m}:  ~~~~  W_o=   3.27
{\rm cm}/s;\\
&&  a = 10^{-3} {\rm cm} = 10^{-5} = 10 \mu m, ~~~~ W_o= 0.0327 {\rm
cm}/s = 117 {\rm cm/hr}\end{eqnarray*}

  For a  water droplet in air,
\[\f{\Delta \rho}{ \rho_f}=\f{1}{ 10^{-3}} = 10^{3},\quad
 \nu = 0.15\  {\rm cm}^2/{\rm sec}\]
then 
\be W_o ={ (217.8) 10^3\o 0.15} a^2\ee 
in cgs units. If $a= 10^{-3} \mbox{ cm} = 10 \mu\mbox{m}$, then $W_o=  1.452 \ 
{\rm cm/sec}$. 


 

\end{enumerate}
\section*{ Details of derivation}
\underline{Details  of (\ref{eqno(10)}).}
\[
\nabla \times\left(\f{\psi}{ r\sin\theta}\vec
e_\phi\right)= \f{1}{ r^2\sin\theta}\left| 
\begin{array}{ccc}
 \vec  e_r & \vec  e_\theta & r \sin\theta\vec  e_\phi \\
						\f{ \p}{ \p r} &\f{\p}{ \p \theta} &\f{\p }{ \p\phi}\\
							0 &  0 & \psi\end{array} \right|\]
\[= \vec e_r\left(\f{1}{ r^2\sin\theta} \f{\p\psi}{\p
\theta}\right)-\vec e_\theta \left(\f{1}{
r\sin\theta} \f{\p\psi}{\p r}\right)\]

\noindent \underline{ Details of (\ref{eqno(11)}).}
\begin{eqnarray*}\lefteqn{
\nabla \times \nabla \times\f{\psi\vec e_\phi }{ r\sin\theta}=\nabla
\times\vec q} \\
&& =\f{1}{ r^2\sin\theta}\left| \begin{array}{ccc}
\vec e_r & r\vec e_\theta &
r\sin\theta\vec e_\phi\\
\f{\p}{\p r} & \f{\p}{\p\theta} & \f{\p}{\p \phi}\\
\f{1}{ r^2\sin\theta} \f{\p\psi}{\p\theta} &
\f{-1}{\sin\theta}\f{\p\psi}{\p r} & 0 \end{array}\right|\\&&
=\f {\vec e_\theta}{
r\sin\theta}\left[\f{\p^2\phi}{\p r^2} + \f{\sin\theta}{ r^2}
\f{\p}{\p\theta} \left( \f{1}{\sin\theta}
\f{\p\psi}{\p\theta}\right)\right] \end{eqnarray*}
 
 \end{document}