Next, we examine how to solve the heat equation in a circular region:
\(\{ (x,y) | x^2+y^2 \leq R^2 \}\) as shown in the following figure.
Figure 7.8.1. A circular region where a point is written in polar coordinates. Instead of solving the equation in cartesian coordinates, we look to write the heat equation in \emph{polar coordinates}. A point in polar coordinates is labelled
\((r,\theta)\) where
\begin{align*}
r \amp = \sqrt{x^2+y^2} \amp \theta \amp = \tan^{-1} \frac{y}{x}
\end{align*}
or written as
\(x\) and
\(y\) in terms of
\(r\) and
\(\theta\text{,}\)
\begin{align*}
x \amp = r \cos \theta \amp y \amp = r \sin \theta .
\end{align*}
To convert the heat equation to polar coordinates, we need to write the right hand side of (
(7.6.1) ) or
\begin{equation*}
\frac{\partial^2 u}{\partial {x}^2} + \frac{\partial^2 u}{\partial {y}^2}
\end{equation*}
in terms of
\(r\) and
\(\theta\text{.}\) This is basically an exercise in using the chain rule with multiple independent variables. We start by finding the first partial derivatives of
\(u\) with respect to
\(x\) and
\(y\text{.}\)
\begin{align*}
u_x = \frac{\partial u}{\partial x} \amp = \frac{\partial u}{\partial r} \frac{\partial r}{\partial x} + \frac{\partial u}{\partial \theta} \frac{\partial \theta}{\partial x} = u_r r_x + u_{\theta} \theta_x \\
u_y = \frac{\partial u}{\partial y} \amp = \frac{\partial u}{\partial r} \frac{\partial r}{\partial y}+ \frac{\partial u}{\partial \theta} \frac{\partial \theta}{\partial y} = u_r r_y + u_{\theta} \theta_y
\end{align*}
and differentiating again, we get:
\begin{align*}
u_{xx} \amp = (u_r r_x)_x + (u_{\theta} \theta_x)_x \qquad\text{using the product rule}\\
\amp = (u_r)_x r_x + u_r r_{xx} + (u_\theta)_x \theta_x + u_{\theta} \theta_{xx} \\
\amp = (u_{rr} r_x + u_{r\theta} \theta_x) r_x + u_r r_{xx} + (u_{\theta r} r_x + u_{\theta\theta} \theta_x) \theta_x + u_{\theta} \theta_{xx}
\end{align*}
\begin{align*}
u_{yy} \amp =(u_rr_y)_y +(u_{\theta}\theta_y) _y \\
\amp = (u_r)_y r_y + u_r r_{yy} + (u_{\theta})_y \theta_y + u_{\theta} \theta_{yy} \\
\amp = (u_{rr} r_y + u_{r\theta} \theta_y )r_y + u_r r_{yy} + (u_{r\theta} r_y + u_{\theta\theta} \theta_y) \theta_y + u_{\theta} \theta_{yy}
\end{align*}
To complete this, we need to find
\(r_x, \theta_x, r_{xx}, \theta_{xx}\)
\begin{align*}
r_x \amp = \frac{x}{\sqrt{x^2+y^2}} = \frac{x}{r} \\
r_{xx} \amp = \frac{\partial }{\partial x} \frac{x}{r} = \frac{r -x r_x}{r^2} = \frac{r-x(x/r)}{r^2} = \frac{r^2-x^2}{r^3} = \frac{y^2}{r^3} \\
\theta_{x} \amp = \frac{1}{1+(y/x)^2} \frac{\partial }{\partial x} \biggl(\frac{y}{x}\biggr) = \frac{x^2}{x^2+y^2} \biggl( -\frac{y}{x^2} \biggr) = -\frac{y}{r^2} \\
\theta_{xx} \amp = \frac{\partial }{\partial x} \biggl( -\frac{y}{r^2} \biggr) = -y \frac{-2}{r^3} \frac{\partial r}{\partial x} =
\frac{2y}{r^3} \frac{x}{r} = \frac{2xy}{r^4}
\end{align*}
as well as
\(r_y, \theta_y, r_{yy}\) and
\(\theta_{yy}\text{.}\)
\begin{align*}
r_y \amp = \frac{y}{\sqrt{x^2+y^2}} = \frac{y}{r} \\
r_{yy} \amp = \frac{\partial }{\partial y} \biggl( \frac{y}{r} \biggr) = \frac{r - y\frac{\partial r}{\partial y} }{r^2} =
\frac{r-y(y/r)}{r^2} = \frac{r^2-y^2}{r^3} = \frac{x^2}{r^3} \\
\theta_y \amp = \frac{1}{1+(y/x)^2} \frac{\partial }{\partial y} \biggl(\frac{y}{x}\biggr) = \frac{x^2}{x^2+y^2} \biggl( \frac{1}{x} \biggr) = \frac{x}{r^2} \\
\theta_{yy} \amp = \frac{\partial }{\partial y}\biggl( \frac{x}{r^2} \biggr) = x \biggl( \frac{-2}{r^3}\biggr) \frac{\partial r}{\partial y} = \frac{-2x}{r^3} \cdot \frac{y}{r} = -\frac{2xy}{r^4}
\end{align*}
So now we write
\(u_{xx}+u_{yy}\)
\begin{align*}
u_{xx} \amp = \biggl(u_{rr} \frac{x}{r} + u_{r\theta} \biggl(-\frac{y}{r^2} \biggr) \biggr)\frac{x}{r} + u_r \frac{y^2}{r^3} + \biggl(u_{r\theta}\frac{x}{r} + u_{\theta\theta}\biggl(-\frac{y}{r^2} \biggr) \biggr)\biggl(-\frac{y}{r^2} \biggr) \\
\amp \qquad \qquad + u_{\theta}\frac{2xy}{r^4} \\
u_{yy} \amp = \biggl(u_{rr} \frac{y}{r} + u_{r\theta} \frac{x}{r^2} \biggr)\frac{y}{r} + u_r \frac{x^2}{r^3} + \biggl(u_{r \theta}\frac{y}{r} + u_{\theta\theta} \frac{x}{r^2} \biggr) \frac{x}{r^2} + u_{\theta}\frac{-2xy}{r^4} \\
u_{xx}+u_{yy} \amp = u_{rr} \frac{x^2+y^2}{r^2} + u_{r\theta} \biggl(-\frac{xy}{r^3} + \frac{xy}{r^3} \biggr) + u_r \biggl( \frac{y^2+x^2}{r^3}\biggr) \\
\amp \qquad \qquad + u_{r \theta} \biggl(-\frac{xy}{r^3} + \frac{xy}{r^3} \biggr) + u_{\theta \theta}\biggl(\frac{y^2}{r^4} + \frac{x^2}{r^4} \biggr) + u_{\theta} \biggl( \frac{2xy}{r^4} - \frac{2xy}{r^4} \biggr) \\
\amp = u_{rr} + \frac{1}{r} u_r + \frac{1}{r^2} u_{\theta\theta}
\end{align*}
Thus the heat equation in polar coordinates is
\begin{equation}
\frac{1}{\kappa} \frac{\partial u}{\partial t} = \frac{\partial^2 u}{\partial {r}^2} + \frac{1}{r} \frac{\partial u}{\partial r}+ \frac{1}{r^2} \frac{\partial^2 u}{\partial {\theta}^2} \label{eq:heat:eqn:polar}\tag{7.8.1}
\end{equation}
Before solving this equation in general, we will examine a simpler equation that is related.
Subsection 7.8.1 Rotationally Symmetric Solutions
The general case is a bit hard to deal with, so we first start with a
rotationally symmetric solution , which means that there is no
\(\theta\) dependence. Thus the term
\(\frac{\partial^2 u}{\partial {\theta}^2} =0\) and the heat equation becomes:
\begin{equation*}
\frac{1}{\kappa} \frac{\partial u}{\partial t} = \frac{\partial^2 u}{\partial {r}^2} + \frac{1}{r} \frac{\partial u}{\partial r}
\end{equation*}
To solve this using separation of variables, let
\(u(r,t)=T(t) R(r)\) and substitution into the PDE:
\begin{align*}
\frac{1}{\kappa} T'R \amp = TR'' + \frac{1}{r} R' T \amp\amp \text{divide by $RT$} \\
\frac{1}{\kappa} \frac{T'}{T} \amp = \frac{r R'' + R'}{r R}
\end{align*}
As before, the only way that the left hand side can equal the right hand side is if each side only depends on a constant or
\begin{align*}
r R'' + R' \amp = -\lambda r R \amp T' + \lambda \kappa T \amp = 0
\end{align*}
The boundary conditions become
\(R(L)=0\) and
\(R(0)\) is finite, so we get the Sturm-Liouville problem:
\begin{align*}
r R'' + R' + \lambda r R \amp = 0, \amp R(L)\amp = 0, \qquad R(0) \text{ is finite}
\end{align*}
This is a difficult problem to solve in that letting
\(R(r) = e^{pr}\) does not work. There are techniques to solve this, but let’s use Maple to get:
\begin{equation*}
R(r) = c_1 J_0(\sqrt{\lambda} r) + c_2 Y_0(\sqrt{\lambda} r)
\end{equation*}
where
\(J_0\) and
\(Y_0\) are bessel functions.
Subsection 7.8.2 Solving the Sturm-Liouville Problem
If we apply the first ``boundary’’ condition that
\(R(0)\) must be finite, we take that
\(c_2=0\) since
\(Y_0(x)\) is not finite at
\(x=0\text{.}\) The next boundary condition is
\(R(L)=0\) or
\begin{equation*}
R(L) = c_1 J_0(\sqrt{\lambda} L) = 0
\end{equation*}
and let
\(\sigma_n\) be the
\(n\) th root of
\(J_0(x)\text{.}\) so
\begin{align*}
\sqrt{\lambda} L \amp = \sigma_n \\
\lambda \amp = \frac{\sigma_n^2}{L^2}
\end{align*}
are the eigenvalues of the problem with eigenfunctions:
\begin{equation*}
R_n(r) = J_0\biggl(\frac{\sigma_n x}{L}\biggr)
\end{equation*}
\begin{equation*}
T' + \kappa \lambda T = 0
\end{equation*}
\begin{equation*}
T_n = e^{-\kappa \sigma_n t/L}
\end{equation*}
\begin{align*}
u \amp = \sum_{n=1}^{\infty} u_n(r,t) \\
\amp = \sum_{n=1}^{\infty} c_n J_0\biggl(\frac{\sigma_n x}{L} \biggr) e^{-\kappa \sigma_n t/L}
\end{align*}
Lastly, using the initial condition, we get:
\begin{equation*}
u(r,0) = f(r) =\sum_{n=1}^{\infty} c_n J_0\biggl(\frac{\sigma_n r}{L} \biggr)
\end{equation*}
which is a Fourier-type series with
\begin{equation*}
c_n = \frac{\int_0^L r f(r) J_0(\sigma_n r/L) \, dr}{\int_0^{L} r J_0^2(\sigma_n r/L) \, dr}
\end{equation*}
Example 7.8.2 .
Find the solution using the initial condition
\(f(r)=1-r\) and let
\(L=1\text{.}\)
Solution .
Again, we need only to find the Fourier Coefficients. The first three are:
\begin{align*}
c_1 \amp = \frac{\int_0^1 r(1-r) J_0(\sigma_1 r) \, dr}{\int_0^1 r J_0^2(\sigma_1 r) \, dr} = 0.7845194227 \\
c_2 \amp = \frac{\int_0^1 r(1-r) J_0(\sigma_2 r) \, dr}{\int_0^1 r J_0^2(\sigma_2 r) \, dr} = 0.06868885655\\
c_3 \amp = \frac{\int_0^1 r(1-r) J_0(\sigma_3 r) \, dr}{\int_0^1 r J_0^2(\sigma_3 r) \, dr} = 0.05311413893
\end{align*}
\begin{equation*}
u(r,t) = \sum_{n=1}^{\infty}c_n J_0(\sigma_n r) e^{-\kappa \sigma_n t}
\end{equation*}
Subsection 7.8.4 Sturm-Liouville Problems from the heat equation
The first Sturm-Liouville problem is:
\begin{equation*}
\Theta'' + \lambda \Theta = 0
\end{equation*}
and
\(\Theta(0) = \Theta(2\pi)\) as well as
\(\Theta'(0) = \Theta'(2\pi)\text{.}\) If
\(\lambda < 0\text{,}\) then
\(\Theta(\theta) = c_1 e^{-\sqrt{-\lambda} \theta} + c_2 e^{\sqrt{-\lambda} \theta}\text{.}\) No solution of this exists that satsifies the boundary conditions. If
\(\lambda = 0\text{,}\) then
\(\Theta(\theta) = c_1 + c_2 \theta\text{.}\) To satisfy the boundary conditions
\(\Theta(\theta) = c_1\text{.}\) If
\(\lambda > 0\text{,}\) then
\(\Theta(\theta) = c_1 \cos \sqrt{\lambda} \theta + c_2 \sin \sqrt{\lambda} \theta\text{.}\)
\begin{equation*}
\Theta (0) = c_1 = \Theta(2\pi) = c_1 \cos 2\pi\sqrt{\lambda} + c_2 \sin 2\pi \sqrt{\lambda}
\end{equation*}
which is satisfied when
\(\lambda = n^2\text{,}\) for
\(n=1,2,3,\ldots\text{.}\) The derivative of
\(\Theta\) is
\begin{equation*}
\Theta(\theta) = -n c_1 \sin n \theta + n c_2 \cos n\theta
\end{equation*}
The second boundary condition:
\begin{equation*}
\Theta'(0) = n c_2 \cos (0) = \Theta'(2\pi) = -n c_1(0) + n c_2 \cos 2\pi n
\end{equation*}
\begin{align*}
\Theta(\theta) \amp = 1 \amp
\Theta(\theta) \amp= \sin n \theta \amp \amp \text{and} \amp \Theta(\theta) \amp = \cos n \theta
\end{align*}
each satisfy the boundary condition. The next differential equation is
\begin{equation*}
r^2 R''+ r R' + (r^2 \nu-n^2) R = 0
\end{equation*}
\begin{equation*}
R(r) = c_1 J_n(\sqrt{\nu} r) + c_2 Y_n(\sqrt{\nu} r)
\end{equation*}
and the boundary conditions are
\(R(0)\) is finite and
\(R(L)=0\text{.}\) The condition at
\(r=0\) sets
\(c_2=0\) and the other condition:
\begin{equation*}
0 = c_1 J_n(\sqrt{\nu}L)
\end{equation*}
\begin{equation*}
\sqrt{\nu} L = \sigma_{n,i}
\end{equation*}
where
\(\sigma_{n,i}\) is the
\(i\) th root of
\(J_n(x)\text{.}\) Thus the eigenvalue is
\begin{equation*}
\nu = \frac{\sigma_{n,i}^2}{L^2}
\end{equation*}
\begin{equation*}
J_n\biggl(\frac{\sigma_{n,i} r}{L} \biggr)
\end{equation*}
\begin{equation*}
\frac{1}{\kappa} \frac{T'}{T} = \nu = \frac{\sigma_{n,i}^2}{L^2}
\end{equation*}
of which the solution is:
\begin{equation*}
T = e^{-\kappa (\sigma_{n,i}^2/L^2) t}
\end{equation*}
\begin{align*}
u_{0,i} \amp = J_0\biggl(\frac{\sigma_{0,i} r}{L} \biggr) e^{-\kappa (\sigma_{0,i}^2/L^2) t} \\
u_{n,i} \amp = R_{n,i}(r) \Theta_n T_{n,i} \\
\amp = J_n\biggl(\frac{\sigma_{n,i} r}{L} \biggr) (a_n \cos n \theta + b_n \sin \theta) e^{-\kappa (\sigma_{n,i}^2/L^2) t}
\end{align*}
and using the principle of superposition the full solution is:
\begin{align*}
u(r,\theta,t) \amp = \sum_{i=1}^{\infty} a_0 J_0\biggl(\frac{\sigma_{0,i} r}{L} \biggr) e^{-\kappa (\sigma_{0,i}^2/L^2) t} + \\
\amp \qquad \sum_{i=1}^{\infty} \sum_{n=1}^{\infty} J_n\biggl(\frac{\sigma_{n,i} r}{L} \biggr) ( a_n \cos n \theta + b_n \sin \theta) e^{-\kappa (\sigma_{n,i}^2/L^2) t}
\end{align*}
Finding the coefficients. In this case we use the initial condition that
\begin{equation*}
u(r,\theta,0) = f(r,\theta)
\end{equation*}
and substituting into the solution:
\begin{equation*}
f(r,\theta) = \sum_{i=1}^{\infty} a_0 J_0\biggl(\frac{\sigma_{0,i} r}{L} \biggr) +
\sum_{i=1}^{\infty} \sum_{n=1}^{\infty} J_n\biggl(\frac{\sigma_{n,i} r}{L} \biggr) ( a_n \cos n \theta + b_n \sin \theta)
\end{equation*}
\begin{align*}
a_{0,i} \amp = \frac{\int_0^L \int_0^{2\pi} r f(r,\theta) J_0(\sigma_{0,i} r/L) \,d\theta \, dr}
{2\pi \int_0^L r J_0(\sigma_{0,i} r/L)^2 \, dr}\\
a_{n,i} \amp = \frac{\int_0^L \int_0^{2\pi} r f(r,\theta) J_n(\sigma_{n,i} r/L)\cos n\theta \,d\theta \, dr}
{\pi \int_0^L r J_n(\sigma_{0,i} r/L)^2 \, dr} \\
b_{n,i} \amp = \frac{\int_0^L \int_0^{2\pi} r f(r,\theta) J_n(\sigma_{n,i} r/L)\sin n\theta \,d\theta \, dr}
{\pi \int_0^L r J_n(\sigma_{0,i} r/L)^2 \, dr}
\end{align*}
Example 7.8.3 .
Solve the equation above when
\(f(r,\theta) = J_1(\sigma_{1,1} r) \sin \theta\text{.}\) Use
\(\kappa=0.04\) and
\(L=1\text{.}\)
Solution .
Again, we just need to compute the coefficients above. Use can either use Maple or make a symmetry argument to see that
\begin{align*}
a_{0,i} \amp = 0 \\
a_{n,i} \amp = 0 \\
b_{1,1} \amp = 1 \\
a_{n,i} \amp = 0 \qquad n \neq 1, i \neq 1.
\end{align*}
\begin{equation*}
u(r, \theta,t) = J_1(\sigma_{1,1} r) \sin \theta e^{-0.04 \sigma{1,1}^2 t}
\end{equation*}
