Heat Equation in Cylindrical Coordinates

We now solve the full heat equation in cylindrical coordinates without assuming radial symmetry.

$$ \dot{T} = \alpha \nabla^2 T $$

where \( T = T(r,\theta,z,t). \) The domain is \( 0 \leq r \leq R, \quad 0 \leq \theta \leq \pi, \quad 0 \leq z \leq H, \quad t \geq 0. \)

The boundary conditions are

$$ T(R,\theta,z,t) = 0 $$ $$ \partial_{\theta} \, T(r,0,z,t) = 0, \quad \partial_{\theta} \, T(r,\pi,z,t) = 0 $$ $$ \partial_z \, T(r,\theta,0,t) = 0, \quad \partial_z \, T(r,\theta,H,t) = 0 $$

with initial condition

$$ T(r,\theta,z,0) = T_i $$

\( T(r,\theta,z,t) = f(r) g(\theta) h(z) p(t) \)

$$ \frac{p'(t)}{p(t)} = \alpha \frac{\nabla^2 \big(f(r) g(\theta) h(z)\big)}{f(r) g(\theta) h(z)} = \alpha \lambda \quad \Rightarrow \quad p(t) = e^{\alpha \lambda t} $$

Since the left-hand side depends only on \( t \), while the right-hand side depends only on spatial variables, both sides must equal a constant.

In cylindrical coordinates

$$ \nabla^2 = \frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial}{\partial r}\right) + \frac{1}{r^2}\frac{\partial^2}{\partial \theta^2} + \frac{\partial^2}{\partial z^2} $$ $$ \frac{1}{f(r)}\frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial}{\partial r}f(r)\right) + \frac{1}{r^2}\frac{g''(\theta)}{g(\theta)} + \frac{h''(z)}{h(z)} = \lambda $$ $$ g'(0) = g'(\pi) = 0 $$ $$ g''(\theta) = A g(\theta), \quad g_n(\theta) = \cos(n\theta), \quad A_n = -n^2, \quad n = 0,1,2,\dots $$

These are the angular modes compatible with the insulating boundary conditions.

$$ h'(0) = h'(H) = 0 $$ $$ h''(z) = B h(z), \quad h_m(z) = \cos\left(\frac{m\pi z}{H}\right), \quad B_m = -\left(\frac{m\pi}{H}\right)^2, \quad m = 0,1,2,\dots $$

These are the axial standing-wave modes.

$$ \frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial}{\partial r}f(r)\right) + \frac{A_n}{r^2}f(r) + B_m f(r) = \lambda f(r) $$ $$ \frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial}{\partial r}f(r)\right) - \frac{n^2}{r^2}f(r) = (\lambda - B_m)f(r) $$

Let

$$ \lambda - B_m = \lambda + \left(\frac{m\pi}{H}\right)^2 = - k^2 $$

The axial eigenvalue therefore shifts the effective radial eigenvalue.

$$ \frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial}{\partial r}f(r)\right) - \frac{n^2}{r^2}f(r) = - k^2 f(r) $$

Rescale \( r \) by \( k \):

$$ x = kr \quad \frac{d}{dx} = \frac{1}{k} \frac{d}{dr} $$

Then the radial equation becomes

$$ \frac{k^2}{x} \frac{d}{dx} \left( x \frac{d}{dx} f(r) \right) - \frac{n^2 k^2}{x^2} f(r) = - k^2 f(r) $$ $$ \frac{d^2}{dx^2}f(r) + \frac{1}{x}\frac{d}{dx}f(r) + \left(1 - \frac{n^2}{x^2}\right)f(r) = 0 $$

This is the Bessel equation. Solutions finite at \( r = 0 \) are

$$ f(r) = J_n(kr) $$

The cylindrical geometry therefore leads naturally to Bessel eigenfunctions.

$$ f(R) = J_n(kR) = 0 $$

\( n = 0 \), \( \, J_0(0) = 1 \) so \( k = 0 \) is not a solution of \( J_0(kR) = 0\)

\( n > 0 \), \( \, J_n(0) = 0 \) so the boundary condition is solved but the corresponding solution vanishes everywhere and cannot be an eigenvector

Hence we must solve \( J_n(kR) = 0 \) for \( k > 0 \), the solutions are labeled by a positive integer \( \ell \), \( k_{n,\ell} \)

The equation above then has eigenvectors \( f_{n,\ell}(r) = J_n(k_{n,\ell} r) \) and eigenvalues \( - k_{n,\ell}^2 \)

$$ \lambda_{n,m,\ell} = - k_{n,\ell}^2 - \left(\frac{m\pi}{H}\right)^2 $$

The general solution is then

$$ T(r,\theta,z,t) = \sum_{n=0}^{\infty} \sum_{m=0}^{\infty} \sum_{\ell=1}^{\infty} c_{n,m,\ell} f_{n,\ell}(r) g_n(\theta) h_m(z) e^{\alpha \lambda_{n,m,\ell} t} $$

The full solution is therefore built from radial, angular, and axial eigenmodes evolving independently in time.

The operators

$$ \partial_{\theta}^2, \quad \partial_z^2, \quad \frac{1}{r} \partial_r (r \partial_r) - \frac{n^2}{r^2} $$

are symmetric under the corresponding inner products. Therefore the following are orthogonal

$$ \langle f_{n,\ell}, f_{n,\tilde{\ell}} \rangle = \delta_{\ell\tilde{\ell}}\|f_{n,\ell}\|^2, \quad \langle u,v \rangle = \int_0^R u(r)v(r)r\,dr $$ $$ \langle g_n, g_{\tilde{n}} \rangle = \delta_{n\tilde{n}}\|g_n\|^2, \quad \langle u,v \rangle = \int_0^\pi u(\theta)v(\theta)\,d\theta $$ $$ \langle h_m, h_{\tilde{m}} \rangle = \delta_{m\tilde{m}}\|h_m\|^2, \quad \langle u,v \rangle = \int_0^H u(z)v(z)\,dz $$

With the initial condition

$$ T(r,\theta,z,0) = \sum_{n=0}^{\infty} \sum_{m=0}^{\infty} \sum_{\ell=1}^{\infty} c_{n,m,\ell} f_{n,\ell}(r) g_n(\theta) h_m(z) = T_i $$

The coefficients are

$$ c_{n,m,\ell} = \frac{\langle f_{n,\ell},1\rangle}{\langle f_{n,\ell},f_{n,\ell}\rangle} \, \frac{\langle g_n,1\rangle}{\langle g_n,g_n\rangle} \, \frac{\langle h_m,1\rangle}{\langle h_m,h_m\rangle} \, T_i $$ $$ \langle g_n,1\rangle = \int_0^\pi \cos(n\theta)d\theta = \begin{cases} \pi & n=0 \\ 0 & n>0 \end{cases} $$ $$ \langle g_n,g_n\rangle = \int_0^\pi \cos^2(n\theta)d\theta = \begin{cases} \, \pi & n=0 \\[0.4em] \dfrac{\pi}{2} & n>0 \end{cases} $$ $$ \langle h_m,1\rangle = \int_0^H \cos\left(\frac{m\pi z}{H}\right)dz = \begin{cases} H & m=0 \\ \, 0 & m>0 \end{cases} $$ $$ \langle h_m,h_m\rangle = \int_0^H \cos^2\left(\frac{m\pi z}{H}\right)dz = \begin{cases} \, H & m=0 \\[0.4em] \dfrac{H}{2} & m>0 \end{cases} $$

\( c_{n,m,\ell} \) is only nonzero for \( n=m=0 \). The constant initial condition therefore excites only the radially symmetric modes.

$$ T(r,\theta,z,t) = \sum_{\ell=1}^{\infty} c_{0,0,\ell} f_{0,\ell}(r) e^{- \alpha k_{0,\ell}^2 t} $$ $$ T(r,\theta,z,t) = T_i \sum_{\ell=1}^{\infty} \frac{\int_0^R J_0(k_{0,\ell} r)\,r\,dr}{\int_0^R J_0^2(k_{0,\ell} r)\,r\,dr} J_0(k_{0,\ell} r) e^{- \alpha k_{0,\ell}^2 t} $$

where \( J_0 \) is a Bessel function and \( k_{0,\ell} \) is the \( \ell \)-th positive solution of \( J_0(k_{0,\ell} R) = 0 \).