Heat Equation in Cylindrical Coordinates

Consider the heat equation in cylindrical coordinates for the temperature distribution inside a half-cylinder of radius \( R \) and height \( H \), with insulating boundary conditions along the rectangular face and zero temperature on the curved and half-disk faces.

$$ \dot{T} = \alpha \nabla^2 T $$ $$ T(R,\theta,z,t) = 0 $$ $$ T(r,\theta,0,t) = T(r,\theta,H,t) = 0 $$ $$ \partial_{\theta} \, T(r,0,z,t) = \partial_{\theta} \, T(r,\pi,z,t) = 0 $$ $$ T(r,\theta,z,0) = T_i $$

With insulating boundary conditions on the face bisection of the cylinder, this problem can be extended symmetrically to an entire cylinder. Therefore, we expect the solutions to be independent of \( \theta \), and we can solve for \( T(r,z,t) \).

$$ T(r,z,t) = f(r) h(z) p(t) $$ $$ \frac{p'(t)}{p(t)} = \alpha \frac{\nabla^2 f(r) h(z)}{f(r) h(z)} = \alpha\lambda $$

Since the left-hand side depends only on time and the right-hand side only on spatial variables, both sides must equal a constant.

$$ p(t) = e^{\alpha\lambda t} $$ $$ \frac{1}{r}\frac{\partial}{\partial r} \left( r\frac{\partial f(r)}{\partial r} \right) + \frac{h''(z)}{h(z)}f(r) = \lambda f(r) $$ $$ h(0) = h(H) = 0 $$ $$ h''(z) = B h(z) $$ $$ h_m(z) = \sin\left(\frac{m\pi z}{H}\right), \quad B_m = -\left(\frac{m\pi}{H}\right)^2, \quad m = 1,2,3,\dots $$

These are the standing-wave modes compatible with the vanishing boundary conditions

$$ \frac{1}{r}\frac{\partial}{\partial r} \left( r\frac{\partial f(r)}{\partial r} \right) = (\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 f(r)}{\partial r} \right) = - k^2 f(r) $$

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

$$ x = kr $$ $$ \frac{1}{x}\frac{\partial}{\partial x} \left( x\frac{\partial}{\partial x}f(r(x)) \right) = -f(r(x)) $$ $$ \frac{\partial^2}{\partial x^2}f(r(x)) + \frac{1}{x}\frac{\partial}{\partial x}f(r(x)) + f(r(x)) = 0 $$

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

$$ f(r) = J_0(kr) $$ $$ f(R) = J_0(kR) = 0 $$

Since \( J_0(0) = 1 \), \( k = 0 \) is not a solution.

We must solve

$$ J_0(kR) = 0, \quad k > 0 $$

Solutions are labelled by a positive integer \( \ell \), \( k_\ell \)

The radial eigenvectors are

$$ f_\ell(r) = J_0(k_\ell r) $$

and the eigenvalues are

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

The general solution is then

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

That is,

$$ T(r,z,t) = \sum_{m=1}^{\infty} \sum_{\ell=1}^{\infty} c_{m,\ell} J_0(k_\ell r) \sin\left(\frac{m\pi z}{H}\right) e^{\alpha\lambda_{m,\ell}t} $$

The solution is therefore built from products of radial and axial eigenmodes, each decaying exponentially in time.

As shown before, with the indicated inner products, the operators

$$ \partial_z^2, \quad r^{-1} \partial_r (r \partial_r) $$

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

$$ \langle f_\ell,f_{\tilde{\ell}}\rangle = \delta_{\ell\tilde{\ell}}\|f_\ell\|^2, \quad \langle u,v\rangle = \int_0^R u(r)v(r) \, r\,dr $$ $$ \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 $$

Using these orthogonality relations, the coefficients \( c_{m,\ell} \) can be determined from the initial condition.

$$ T(r,z,0) = \sum_{m=1}^{\infty} \sum_{\ell=1}^{\infty} c_{m,\ell} f_\ell(r) h_m(z) = T_i $$ $$ c_{m,\ell} = T_i \frac{\langle f_\ell,1\rangle}{\langle f_\ell,f_\ell\rangle} \frac{\langle h_m,1\rangle}{\langle h_m,h_m\rangle} $$ $$ \langle h_m,1\rangle = \int_0^H \sin\left(\frac{m\pi z}{H}\right)\,dz $$ $$ = \left[ -\frac{H}{m\pi} \cos\left(\frac{m\pi z}{H}\right) \right]_0^H $$ $$ = \frac{H}{m\pi}\left(1-\cos(m\pi)\right) $$ $$ = \begin{cases} \dfrac{2H}{m\pi} & m \text{ odd} \\[0.4em] \: \: \, 0 & m \text{ even} \end{cases} $$ $$ \langle h_m,h_m\rangle = \int_0^H \sin^2\left(\frac{m\pi z}{H}\right)\,dz = \frac{H}{2} $$

Therefore,

$$ c_{m,\ell} = T_i \frac{4}{m\pi} \frac{\langle f_\ell,1\rangle}{\langle f_\ell,f_\ell\rangle}, \quad m \text{ odd} $$

Thus the constant initial condition excites only the odd axial modes.

Hence, we have

$$ T(r,z,t) = \frac{4T_i}{\pi} \sum_{\substack{m=1 \\ m\text{ odd}}}^{\infty} \sum_{\ell=1}^{\infty} \frac{1}{m} \frac{ \int_0^R J_0(k_\ell r)\,r\,dr }{ \int_0^R J_0^2(k_\ell r)\,r\,dr } J_0(k_\ell r) \sin\left(\frac{m\pi z}{H}\right) e^{\alpha\lambda_{m,\ell}t} $$

where \( J_0 \) is a Bessel function and \( k_\ell \) is the \( \ell \)-th positive solution of \( J_0(k_\ell R) = 0 \, \) and \( \, \displaystyle \lambda_{m,\ell} = - k_\ell^2 - \left(\frac{m\pi}{H}\right)^2 \)

Higher radial and axial modes decay more rapidly, so the temperature distribution gradually approaches the equilibrium state \( T = 0 \).

See the full solution with angular dependence.