# Peeter Joot's (OLD) Blog.

• ## Archives

 Adam C Scott on avoiding gdb signal noise… Ken on Scotiabank iTrade RESP …… Alan Ball on Oops. Fixing a drill hole in P… Peeter Joot's B… on Stokes theorem in Geometric… Exploring Stokes The… on Stokes theorem in Geometric…

• 294,114

## Spin down of coffee in a bottomless cup.

Posted by peeterjoot on April 25, 2012

[Click here for a PDF of this post with nicer formatting (especially if my latex to wordpress script has left FORMULA DOES NOT PARSE errors.)]

# Motivation.

Here’s a variation of a problem outlined in section 2 of [1], which looked at the time evolution of fluid with initial rotational motion, after the (cylindrical) rotation driver stops, later describing this as the spin down of a cup of tea. I’ll work the problem in more detail than in the text, and also make two refinements.

• I drink coffee and not tea.
• I stir my coffee in the interior of the cup and not on the outer edge.

Because of the second point I’ll model my stir stick as a rotating cylinder in the cup and not by somebody spinning the cup itself to stir the tea. This only changes the solution for the steady state part of the problem.

# Guts

We’ll work in cylindrical coordinates following the conventions of figure (1).

Figure 1: Fluid flow in nested cylinders.

We’ll assume a solution that with velocity azimuthal in direction, and both pressure and velocity that are only radially dependent.

\begin{aligned}\mathbf{u} = u(r) \hat{\boldsymbol{\phi}}.\end{aligned} \hspace{\stretch{1}}(2.1)

\begin{aligned}p = p(r)\end{aligned} \hspace{\stretch{1}}(2.2)

Let’s first verify that this meets the non-compressible condition that eliminates the $\mu \boldsymbol{\nabla} (\boldsymbol{\nabla} \cdot \mathbf{u})$ term from Navier-Stokes

\begin{aligned}\boldsymbol{\nabla} \cdot \mathbf{u}&=\left(\hat{\mathbf{r}} \partial_r + \frac{\hat{\boldsymbol{\phi}}}{r} \partial_\phi + \hat{\mathbf{z}} \partial_z\right) \cdot \left(u \hat{\boldsymbol{\phi}}\right) \\ &=\hat{\boldsymbol{\phi}} \cdot\left(\hat{\mathbf{r}} \partial_r u + \frac{\hat{\boldsymbol{\phi}}}{r} \partial_\phi u + \hat{\mathbf{z}} \partial_z u\right) +u\left(\hat{\mathbf{r}} \cdot \partial_r \hat{\boldsymbol{\phi}} + \frac{\hat{\boldsymbol{\phi}}}{r} \cdot \partial_\phi \hat{\boldsymbol{\phi}} + \hat{\mathbf{z}} \cdot \partial_z \hat{\boldsymbol{\phi}}\right) \\ &=\hat{\boldsymbol{\phi}} \cdot \hat{\mathbf{r}} \partial_r u +u\frac{\hat{\boldsymbol{\phi}}}{r} \cdot \left(-\hat{\mathbf{r}}\right) \\ &= 0.\end{aligned}

Good. Now let’s express each of the terms of Navier-Stokes in cylindrical form. Our time dependence is

\begin{aligned}\rho \partial_t u(r, t) \hat{\boldsymbol{\phi}}=\rho \hat{\boldsymbol{\phi}} \partial_t u.\end{aligned} \hspace{\stretch{1}}(2.3)

Our inertial term is

\begin{aligned}\begin{aligned}\rho (\mathbf{u} \cdot \boldsymbol{\nabla}) \mathbf{u}&=\frac{\rho u}{r} \partial_\phi (u \hat{\boldsymbol{\phi}}) \\ &=\frac{\rho u^2}{r} (-\hat{\mathbf{r}}).\end{aligned}\end{aligned} \hspace{\stretch{1}}(2.4)

Our pressure term is

\begin{aligned}-\boldsymbol{\nabla} p=-\hat{\mathbf{r}} \partial_r p,\end{aligned} \hspace{\stretch{1}}(2.5)

and our Laplacian term is

\begin{aligned}\begin{aligned}\mu \boldsymbol{\nabla}^2 \mathbf{u}&=\mu \left( \frac{1}{{r}} \partial_r ( r \partial_r) + \frac{1}{{r^2}} \partial_{\phi\phi} + \partial_{z z}\right)u(r) \hat{\boldsymbol{\phi}} \\ &=\mu \left( \frac{\hat{\boldsymbol{\phi}}}{r} \partial_r ( r \partial_r u) + \frac{-\hat{\mathbf{r}} u}{r^2} \right).\end{aligned}\end{aligned} \hspace{\stretch{1}}(2.6)

Putting things together, we find that Navier-Stokes takes the form

\begin{aligned}\rho \hat{\boldsymbol{\phi}} \partial_t u+\frac{\rho u^2}{r} (-\hat{\mathbf{r}})=-\hat{\mathbf{r}} \partial_r p+\mu \left( \frac{\hat{\boldsymbol{\phi}}}{r} \partial_r ( r \partial_r u) + \frac{-\hat{\boldsymbol{\phi}} u}{r^2} \right),\end{aligned} \hspace{\stretch{1}}(2.7)

which nicely splits into an separate equations for the $\hat{\boldsymbol{\phi}}$ and $\hat{\mathbf{r}}$ directions respectively

\begin{aligned}\frac{1}{{\nu}} \partial_t u=\frac{1}{r} \partial_r ( r \partial_r u) - \frac{u}{r^2}\end{aligned} \hspace{\stretch{1}}(2.8a)

\begin{aligned}\frac{\rho u^2}{r}=\partial_r p.\end{aligned} \hspace{\stretch{1}}(2.8b)

Before $t = 0$ we seek the steady state, the solution of

\begin{aligned}r \partial_r ( r \partial_r u) - u = 0.\end{aligned} \hspace{\stretch{1}}(2.9)

We’ve seen that

\begin{aligned}u(r) = A r + \frac{B}{r}\end{aligned} \hspace{\stretch{1}}(2.10)

is the general solution, and can now fit this to the boundary value constraints. For the interior portion of the cup we have

\begin{aligned}{\left.{{A r + \frac{B}{r}}}\right\vert}_{{r = 0}} = 0\end{aligned} \hspace{\stretch{1}}(2.11)

so $B = 0$ is required. For the interface of the “stir-stick” (moving fast enough that we can consider it having a cylindrical effect) at $r = R_1$ we have

\begin{aligned}A R_1 = \Omega R_1,\end{aligned} \hspace{\stretch{1}}(2.12)

so the interior portion of our steady state coffee velocity is just

\begin{aligned}\mathbf{u} = \Omega r \hat{\boldsymbol{\phi}}.\end{aligned} \hspace{\stretch{1}}(2.13)

Between the cup edge and the stir-stick we have to solve

\begin{aligned}A R_1 + \frac{B}{R_1} &= \Omega R_1 \\ A R_2 + \frac{B}{R_2} &= 0,\end{aligned} \hspace{\stretch{1}}(2.14)

or

\begin{aligned}A R_1^2 + B &= \Omega R_1^2 \\ A R_2^2 + B &= 0.\end{aligned} \hspace{\stretch{1}}(2.16)

Subtracting we find

\begin{aligned}A = -\frac{\Omega R_1^2}{R_2^2 - R_1^2}\end{aligned} \hspace{\stretch{1}}(2.18a)

\begin{aligned}B = \frac{\Omega R_1^2 R_2^2}{R_2^2 - R_1^2},\end{aligned} \hspace{\stretch{1}}(2.18b)

so our steady state coffee flow is

\begin{aligned}\mathbf{u} =\left\{\begin{array}{l l}\Omega r \hat{\boldsymbol{\phi}}& \quad \mbox{latex r \in [0, R_1]} \\ \frac{\Omega R_1^2}{R_2^2 – R_1^2} \left( \frac{R_2^2}{r} -r \right)\hat{\boldsymbol{\phi}}& \quad \mbox{$r \in [R_1, R_2]$} \\ \end{array}\right.\end{aligned} \hspace{\stretch{1}}(2.19)

## Time evolution.

We can use a separation of variables technique with $u(r, t) = R(r) T(t)$ to find

\begin{aligned}\frac{1}{{\nu}} \frac{T'}{T} = \frac{1}{{R}} \left( \frac{1}{r} \partial_r ( r \partial_r R) - \frac{R}{r^2}\right)= -\lambda^2,\end{aligned} \hspace{\stretch{1}}(2.20)

which gives us

\begin{aligned}T \propto e^{-\lambda^2 \nu t},\end{aligned} \hspace{\stretch{1}}(2.21)

and $R$ specified by

\begin{aligned}0 = r^2 \frac{d^2 R}{dr^2} + r \frac{d R}{dr} + R \left( r^2 \lambda^2 - 1 \right).\end{aligned} \hspace{\stretch{1}}(2.22)

Checking [2] (9.1.1) we see that this can be put into the standard form of the Bessel equation if we eliminate the $\lambda$ term. We can do that writing $z = r \lambda$, $\mathcal{R}(z) = R(z/\lambda)$ and noting that $r d/dr = z d/dz$ and $r^2 d^2/dr^2 = z^2 d^2/dz^2$, which gives us

\begin{aligned}0 = z^2 \frac{d^2 \mathcal{R}}{dr^2} + z \frac{d \mathcal{R}}{dr} + \mathcal{R} \left( z^2 - 1 \right).\end{aligned} \hspace{\stretch{1}}(2.23)

The solutions are

\begin{aligned}\mathcal{R}(z) = J_{\pm 1}(z), Y_{\pm 1}(z).\end{aligned} \hspace{\stretch{1}}(2.24)

From (9.1.5) of the handbook we see that the plus and minus variations are linearly dependent since $J_{-1}(z) = -J_1(z)$ and $Y_{-1}(z) = -Y_1(z)$, and from (9.1.8) that $Y_1(z)$ is infinite at the origin, so our general solution has to be of the form

\begin{aligned}\mathbf{u}(r, t) = \hat{\boldsymbol{\phi}} \sum_\lambda c_\lambda e^{-\lambda^2 \nu t} J_{1}(r \lambda).\end{aligned} \hspace{\stretch{1}}(2.25)

In the text, I see that the transformation $\lambda \rightarrow \lambda/a$ (where $a$ was the radius of the cup) is made so that the Bessel function parameter was dimensionless. We can do that too but write

\begin{aligned}\mathbf{u}(r, t) = \hat{\boldsymbol{\phi}} \sum_\lambda c_\lambda e^{-\frac{\lambda^2}{R_2^2} \nu t} J_{1}\left(\lambda \frac{r}{R_2}\right).\end{aligned} \hspace{\stretch{1}}(2.26)

Our boundary value constraint is that we require this to match 2.19 at $t = 0$. Let’s write $R_2 = R$, $R_1 = a R$, $z = r/R$, so that we are working in the unit circle with $z \in [0, 1]$. Our boundary problem can now be expressed as

\begin{aligned}\frac{1}{{\Omega R}} \sum_\lambda c_\lambda J_{1}\left(\lambda z\right)=\left\{\begin{array}{l l}z & \quad \mbox{latex z \in [0, a]} \\ \frac{1}{\frac{R^2}{a^2} – 1} \left( \frac{1}{{z}} – z\right)& \quad \mbox{$z \in [a, 1]$} \\ \end{array}\right.\end{aligned} \hspace{\stretch{1}}(2.27)

Let’s pull the $\Omega R$ factor into $c_\lambda$ and state the problem to be solved as

\begin{aligned}\mathbf{u}(r, t) = \Omega R \hat{\boldsymbol{\phi}} \sum_{i=1}^n c_i e^{-\frac{\lambda_i^2}{R^2} \nu t} J_{1}\left(\lambda_i \frac{r}{R}\right)\end{aligned} \hspace{\stretch{1}}(2.28a)

\begin{aligned}\sum_{i = 1}^n c_i J_{1}\left(\lambda_i z\right) = \phi(z)\end{aligned} \hspace{\stretch{1}}(2.28b)

\begin{aligned}\phi(z) = \left\{\begin{array}{l l}z & \quad \mbox{latex z \in [0, a]} \\ \frac{a^2}{1 – a^2} \left( \frac{1}{{z}} – z\right)& \quad \mbox{$z \in [a, 1]$} \\ \end{array}\right..\end{aligned} \hspace{\stretch{1}}(2.28c)

Looking at section 2.7 of [3] it appears the solutions for $c_i$ can be obtained from

\begin{aligned}c_i = \frac{\int_0^1 r\phi(z) J_1(\lambda_i z) dz}{\int_0^1 r J_1^2(\lambda_i z) dz},\end{aligned} \hspace{\stretch{1}}(2.29)

where $\lambda_i$ are the zeros of $J_1$.

To get a feel for these, a plot of the first few of these fitting functions is shown in figure (2).

Figure 2: First four zero crossing Bessel functions.

Using Mathematica in bottomlessCoffee.cdf, these coefficients were calculated for $a = 0.6$. The $n = 1, 3, 5$ approximations to the fitting function are plotted with a comparison to the steady state velocity profile in figure (3).

Figure 3: Bessel function fitting for the steady state velocity profile for n = 1, 3, 5.

As indicated in the text, the spin down is way too slow to match reality (this can be seen visually in the worksheet by animating it).

# References

[1] D.J. Acheson. Elementary fluid dynamics. Oxford University Press, USA, 1990.

[2] M. Abramowitz and I.A. Stegun. {\em Handbook of mathematical functions with formulas, graphs, and mathematical tables}, volume 55. Dover publications, 1964.

[3] H. Sagan. Boundary and eigenvalue problems in mathematical physics. Dover Pubns, 1989.