Heat Equation in spherical coordinates

Consider the problem of a sphere of material that starts at a non-uniform temperature, $T = r^{2}$ and is covered with insulation on the outer surface so that no heat gets out. We take the coordinate $r$ to measure position radially from the centre of the sphere with the outer surface given at $ r = 1.$ We take t as time and the variable $ T (r, t)$ as the temperature. The equation governing the heat flow, is the heat equation $$\frac{\delta T}{\delta t}=\frac{1}{r^2}\frac{\delta}{\delta r}(r^2\frac{\delta T}{\delta r}), 0\leq r\leq 1, t>0$$

We insist that T remains finite as $ r → 0$ and the outer surface boundary condition is $$\frac{\delta T}{\delta t}(1,t)=0, t>0$$ with the initial condition $$T(r,0) = r^2, 0 < r < 1.$$

(a) Find the solution T (t, r) using separation of variables. Note: you can leave the answer in a form where eigenvalues are given by the roots of an equation. By exploiting the orthogonality of the eigenfunctions you should give the integral formulas necessary to compute the coefficients in the solutions. You do not need to evaluate the integrals.

(b) Find numerically, to five decimal places accuracy, the value of the smallest of the eigenvalues in (a), corresponding to the first non-constant term in the solution.

This is what I have obtained so far:

$$\frac{\delta T}{\delta t}=\frac{1}{r^2}\frac{\delta}{\delta r}(r^2\frac{\delta T}{\delta r})$$ $=\frac{\delta^2 T}{\delta r^2}+ \frac{2}{r} \frac{\delta T}{\delta r}$

Let $T=R(r)Q(t)$

$\frac{\delta T}{\delta t}= R(r)Q'(t)$

$\frac{\delta T}{\delta r}= R'(r)Q(t)$

$\frac{\delta^2 T}{\delta r^2}= R''(r)Q(t)$

Then putting back into the equation we get $$\frac{Q'(t)}{Q(t)}=\frac{R''(r)}{R(r)} + \frac{2}{r}\frac{R'(r)}{R(r)}= \mu$$ which gives the two equations: $$Q'(t)-\mu Q(t)=0$$ and $$rR''(r)+2R'(r)-\mu rR(r)=0$$

Then I have done

Let $o(r)=rR(r)$

$o'(r)=R(r)+rR'(r)$

$o''(r)=2R'(r)+rR''(r)$

Putting this back into the equation we get $o''(r)-\mu o(r)=0$

Now how do I proceed? using the boundary conditions how do I implement this into the equation.

Solution 1:

The equation $$ R''(r)+2rR'(r) - \mu rR(r) = 0 $$ can be written as $$ (rR(r))''-\mu(rR(r))=0. $$ Therefore, if you substitute $S(r) = rR(r)$, you get $S''-\mu S=0$, which is a familiar equation.

You have the wrong condition for insulated. Insulated means that the normal derivative of the heat distribution at the boundary is $0$. This is because heat flows according to the temperature gradient; it flows from hot to cold. So that means you want $$ R'(1) = 0. $$ In terms of $S(r)=(rR(r))$, that gives $$ \left.\frac{d}{dr}\left(\frac{1}{r}S(r)\right)\right|_{r=1}=0, \\ -S(1)+S'(1)=0. $$ The final equations for $Q(t)$ and $S(r)$ are $$ Q(t) = e^{\mu t},\\ S''(r)-\mu S(r) = 0,\\ S(0)=0,\;\;\;\; S(1)-S'(1)=0. $$ Note: There is a $\mu > 0$ that satisfies the equation and the condition at $1$: $\mu = 1$ with solution $S(r)=e^{r}$. However, this solution is eliminated because $S(0) \ne 0$. For $\mu < 0$, let $\rho = -\mu$, and the $S$ solutions which satisfy $S(0)=0$ are $$ S_{\rho}(r) = \frac{\sin(\sqrt{\rho}\,r)}{\sqrt{\rho}}. $$ The above is the correct solution rather than just $\sin(\sqrt{\rho}r)$; this is not obvious, but the general theory proves that the correct way is to choose the values of $S_{\rho}(0)$ and $S_{\rho}'(0)$ to be constants, such as the case above where $S_{\rho}(0)=0$ and $S_{\rho}'(0)=1$. Then, it turns out that the limit case where $\rho \rightarrow 0$ will automatically work without a special case. So, the correct $S_{0}$ is $$ \lim_{\rho\rightarrow 0}S_{\rho}(r)=r. $$ That trick eliminates special cases, and guarantees that you get the correct solutions. The resulting functions will have an everywhere convergent power series in $\rho$.

The values of $\rho \ge 0$ for which $S_{\rho}$ matches the right endpoint condition must satisfy $$ \frac{\sin(\sqrt{\rho}\,1)}{\sqrt{\rho}}-\cos(\sqrt{\rho}\,1)=0,\\ \sqrt{\rho} = \tan(\sqrt{\rho}). $$ Notice that $\rho = 0$ is a solution, and, in fact, $S_{0}(r)=r$ is a solution because $S_{0}'(1)=1$ and $S_{0}(1)=1$. The other solutions $\rho > 0$ are those satisfying the above transcendental equation $\sqrt{\rho}=\tan(\sqrt{\rho})$.

Remember that the final solutions are $R(r)=\frac{1}{r}S(r)$, which gives you $$ R_{0}(r)=1,\;\; R_{n}(r) =\frac{\sin(\lambda_{n}r)}{r},\;\; n=1,2,3,\cdots $$ and the separated solution $$ T(r,t) = C_{0}+\sum_{n=1}^{\infty} C_{n}e^{-\lambda_{n}^{2}t}\frac{\sin(\lambda_{n}r)}{r}. $$ where $\lambda_1 < \lambda_2 < \lambda_3 < \cdots$ are the positive solutions of $\tan \lambda = \lambda$. The functions $\sin(\lambda_{n}r)/r$ are orthogonal on $[0,1]$ with respect to the weight function $r^{2}$, which can be seen by putting the original equation for $R$ into selfadjoint form: $$ -(r^{2}R(r))' = \mu r^{2}R(r). $$ The correct weight function is the multiplier on the right: $r^{2}$. Using this orthogonality with respect to $r^{2}dr$, the constants $C_{n}$ are found from $$ r^{2} = F(r,0) = C_{0} + \sum_{n=1}^{\infty}C_{n}\frac{\sin(\lambda_n r)}{r} $$ through the orthogonality conditions: $$ \int_{0}^{1}r^{2}\cdot 1 \cdot r^{2}dr = C_{0}\int_{0}^{1}1\cdot 1\cdot r^{2}dr \\ \implies C_{0} = \frac{3}{4},\\ \int_{0}^{1}r^{2}\frac{\sin(\lambda_{n}r)}{r}r^{2}dr = C_{n}\int_{0}^{1}\frac{\sin(\lambda_{n}r)}{r}\frac{\sin(\lambda_{n}r)}{r}r^{2}dr. \\ \implies C_{n} = \frac{\int_{0}^{1}r^{3}\sin(\lambda_{n}r)dr} {\int_{0}^{1}\sin^{2}(\lambda_{n}r)dr} $$ You can integrate and use the fact that $\tan(\lambda_{n})=\lambda_{n}$ in order to express everything in terms of the roots $\lambda_{n}$.

Reference: Free Google E-Book: Joseph Fourier - Analaytical Theory of Heat.
Go to page 270 and look at his graph. You can back up to the theory above that, and more graphs after.

The details I presented here are basically straight out of Fourier's original work. Take a look at the link above, which is the English translation of Fourier's original manuscript. You can download it for free. Fourier solves this problem completely, and even describes how to numerically find the eigenvalues. This is for a sphere, which he used to analyze a thermometer.

Solution 2:

From a basic ODE class you should know that the solution to $y'' = \lambda y$ is either ( $\lambda <0$) $$y = A \sin (\sqrt {-\lambda} t) + B \cos ( \sqrt {-\lambda} t ) $$ ( $\lambda =0$ ) $$y = A t + B $$ ( $\lambda >0$) $$ y = A \exp (\sqrt \lambda t) + B \exp ( -\sqrt \lambda t )$$ Use the boundary conditions to see what case you fall into, this will give you your eigenvalue equation.

Hint: What does the condition $$ \frac{ \partial T }{\partial t } (1,t) = 0 \quad t>0 $$ imply about $\mu$? Then use Fourier series with the initial condition to deduce the solution.