Find the normal modes of a rectangular drum with sides of length L_{x} and L_{y}
solution
The geometry of the problem is
Using Cartesian coordinates. Wave displacement is u\equiv u\left ( x,y,t\right ) (out of page).\begin{align*} \frac{\partial ^{2}u\left ( x,y,t\right ) }{\partial t^{2}} & =c^{2}\left ( \frac{\partial ^{2}u}{\partial x^{2}}+\frac{\partial ^{2}u}{\partial y^{2}}\right ) \\ 0 & <x<L_{x}\\ 0 & <y<L_{y} \end{align*}
Boundary conditions on x\begin{align*} u\left ( 0,y,t\right ) & =0\\ u\left ( L_{x},y,t\right ) & =0 \end{align*}
And boundary conditions on y\begin{align*} u\left ( x,0,t\right ) & =0\\ u\left ( x,L_{y},t\right ) & =0 \end{align*}
Solution
Let u=X\left ( x\right ) Y\left ( y\right ) T\left ( t\right ) . Substituting into the PDE gives\begin{align*} \frac{1}{c^{2}}T^{\prime \prime }XY & =X^{\prime \prime }YT+Y^{\prime \prime }XT\\ \frac{1}{c^{2}}\frac{T^{\prime \prime }}{T} & =\frac{X^{\prime \prime }}{X}+\frac{Y^{\prime \prime }}{Y} \end{align*}
Hence, using \lambda as first separation constant we obtain\begin{align*} \frac{1}{c^{2}}\frac{T^{\prime \prime }}{T} & =-\lambda \\ \frac{X^{\prime \prime }}{X}+\frac{Y^{\prime \prime }}{Y} & =-\lambda \end{align*}
The time ODE becomes T^{\prime \prime }+c^{2}\lambda T=0 And the space ODE becomes \frac{X^{\prime \prime }}{X}+\frac{Y^{\prime \prime }}{Y}=-\lambda Separating the space ODE again \frac{X^{\prime \prime }}{X}=-\lambda -\frac{Y^{\prime \prime }}{Y}=-\mu Where \mu is the new separation variable. This gives two new separate ODE’s\begin{align*} \frac{X^{\prime \prime }}{X} & =-\mu \\ -\lambda -\frac{Y^{\prime \prime }}{Y} & =-\mu \end{align*}
Or\begin{align*} X^{\prime \prime }+\mu X & =0\\ Y^{\prime \prime }+Y\left ( \lambda -\mu \right ) & =0 \end{align*}
Solving for X ODE first, and knowing that \mu >0 from nature of boundary conditions, we obtain X\left ( x\right ) =A\cos \left ( \sqrt{\mu }x\right ) +B\sin \left ( \sqrt{\mu }x\right ) Applying B.C. at x=0 0=A Hence X\left ( x\right ) =B\sin \left ( \sqrt{\mu }x\right ) . Applying B.C. at x=L_{x} 0=B\sin \left ( \sqrt{\mu }L_{x}\right ) Hence \begin{align} \sqrt{\mu }L_{x} & =n\pi \nonumber \\ \mu _{n} & =\left ( \frac{n\pi }{L_{x}}\right ) ^{2}\qquad n=1,2,3,\cdots \tag{1} \end{align}
Therefore the X_{n}\left ( x\right ) solution is\begin{equation} X_{n}\left ( x\right ) =B_{n}\sin \left ( \frac{n\pi }{L_{x}}x\right ) \qquad n=1,2,3,\cdots \tag{2} \end{equation} Solving the Y\left ( y\right ) ODE using the same eigenvalues found above Y^{\prime \prime }+Y\left ( \lambda -\left ( \frac{n\pi }{L_{x}}\right ) ^{2}\right ) =0 The solution is Y\left ( y\right ) =C\cos \left ( \sqrt{\lambda -\left ( \frac{n\pi }{L_{x}}\right ) ^{2}}y\right ) +D\sin \left ( \sqrt{\lambda -\left ( \frac{n\pi }{L_{x}}\right ) ^{2}}y\right ) Applying first B.C. Y\left ( 0\right ) =0 gives 0=C Hence Y\left ( y\right ) =D\sin \left ( \sqrt{\lambda -\left ( \frac{n\pi }{L_{x}}\right ) ^{2}}y\right ) Applying second B.C. Y\left ( L_{y}\right ) =0 0=D\sin \left ( \sqrt{\lambda -\left ( \frac{n\pi }{L_{x}}\right ) ^{2}}L_{y}\right ) Hence\begin{align*} \sqrt{\lambda -\left ( \frac{n\pi }{L_{x}}\right ) ^{2}}L_{y} & =m\pi \qquad m=1,2,3,\cdots \\ \lambda _{nm}-\left ( \frac{n\pi }{L_{x}}\right ) ^{2} & =\left ( \frac{m\pi }{L_{y}}\right ) ^{2}\\ \lambda _{nm} & =\left ( \frac{m\pi }{L_{y}}\right ) ^{2}+\left ( \frac{n\pi }{L_{x}}\right ) ^{2}\qquad n=1,2,3,\cdots ,m=1,2,3,\cdots \end{align*}
Hence the Y_{nm} solution is Y_{nm}=D_{nm}\sin \left ( \frac{m\pi }{L_{y}}y\right ) \qquad n=1,2,3,\cdots ,m=1,2,3,\cdots We notice that X_{n}\left ( x\right ) solution depends on n only, while Y_{nm}\left ( y\right ) solution depends on n and m. Now that we found \lambda we can we solve the time T\left ( t\right ) ode\begin{align*} T_{nm}^{\prime \prime }+c^{2}\lambda _{nm}T_{nm} & =0\\ T_{nm}\left ( t\right ) & =E_{nm}\cos \left ( c\sqrt{\lambda _{nm}}t\right ) +F_{nm}\sin \left ( c\sqrt{\lambda _{nm}}t\right ) \end{align*}
Combining all solution , and merging all constants into two, we find\begin{align*} u_{nm}\left ( x,y,t\right ) & =X_{n}\left ( x\right ) Y_{nm}\left ( y\right ) T_{nm}\left ( t\right ) \\ & =\left ( B_{n}X_{n}\right ) \left ( D_{nm}\sin \left ( \frac{m\pi }{L_{y}}y\right ) \right ) \left ( E_{nm}\cos \left ( c\sqrt{\lambda _{nm}}t\right ) +F_{nm}\sin \left ( c\sqrt{\lambda _{nm}}t\right ) \right ) \\ & =B_{n}X_{n}\sin \left ( \frac{m\pi }{L_{y}}y\right ) \left ( E_{nm}^{\prime }\cos \left ( c\sqrt{\lambda _{nm}}t\right ) +F_{nm}^{\prime }\sin \left ( c\sqrt{\lambda _{nm}}t\right ) \right ) \\ & =X_{n}\sin \left ( \frac{m\pi }{L_{y}}y\right ) \left ( E_{nm}^{\prime \prime }\cos \left ( c\sqrt{\lambda _{nm}}t\right ) +F_{nm}^{\prime \prime }\sin \left ( c\sqrt{\lambda _{nm}}t\right ) \right ) \end{align*}
Where E_{nm}^{\prime \prime },F_{nm}^{\prime \prime } are the new constants after merging them with the other constants. Renaming E_{nm}^{\prime \prime }=A_{nm},F_{nm}^{\prime \prime }=B_{nm} the above solution can be written as\begin{align} u\left ( x,y,t\right ) & =\sum _{n=1}^{\infty }\sum _{m=1}^{\infty }X_{n}\left ( x\right ) Y_{mn}\left ( y\right ) T_{mn}\left ( t\right ) \nonumber \\ & =\sum _{n=1}^{\infty }\sum _{m=1}^{\infty }A_{nm}\sin \left ( \frac{n\pi }{L_{x}}x\right ) \sin \left ( \frac{m\pi }{L_{y}}y\right ) \cos \left ( c\sqrt{\lambda _{nm}}t\right ) \nonumber \\ & +\sum _{n=1}^{\infty }\sum _{m=1}^{\infty }B_{nm}\sin \left ( \frac{n\pi }{L_{x}}x\right ) \sin \left ( \frac{m\pi }{L_{y}}y\right ) \sin \left ( c\sqrt{\lambda _{nm}}t\right ) \tag{3} \end{align}
To solve this completely, we apply initial conditions to find A_{nm},B_{nm}. But the problem is just asking for the normal modes. These are given by X_{n}\left ( x\right ) Y_{mn}\left ( y\right ) . Therefore for n=1, we have the modes \sin \left ( \frac{\pi }{L_{x}}x\right ) \sin \left ( \frac{\pi }{L_{y}}y\right ) ,\sin \left ( \frac{\pi }{L_{x}}x\right ) \sin \left ( \frac{2\pi }{L_{y}}y\right ) ,\sin \left ( \frac{\pi }{L_{x}}x\right ) \sin \left ( \frac{3\pi }{L_{y}}y\right ) ,\cdots and for n=2 we have \sin \left ( \frac{2\pi }{L_{x}}x\right ) \sin \left ( \frac{\pi }{L_{y}}y\right ) ,\sin \left ( \frac{2\pi }{L_{x}}x\right ) \sin \left ( \frac{2\pi }{L_{y}}y\right ) ,\sin \left ( \frac{2\pi }{L_{x}}x\right ) \sin \left ( \frac{3\pi }{L_{y}}y\right ) ,\cdots and so on.
n | m=1 | 2 | 3 | 4 |
1 | \sin \left ( \frac{\pi }{L_{x}}x\right ) \sin \left ( \frac{\pi }{L_{y}}y\right ) | \sin \left ( \frac{\pi }{L_{x}}x\right ) \sin \left ( \frac{2\pi }{L_{y}}y\right ) | \sin \left ( \frac{\pi }{L_{x}}x\right ) \sin \left ( \frac{3\pi }{L_{y}}y\right ) | \cdots |
2 | \sin \left ( \frac{2\pi }{L_{x}}x\right ) \sin \left ( \frac{\pi }{L_{y}}y\right ) | \sin \left ( \frac{2\pi }{L_{x}}x\right ) \sin \left ( \frac{2\pi }{L_{y}}y\right ) | \sin \left ( \frac{2\pi }{L_{x}}x\right ) \sin \left ( \frac{3\pi }{L_{y}}y\right ) | \cdots |
3 | \sin \left ( \frac{3\pi }{L_{x}}x\right ) \sin \left ( \frac{\pi }{L_{y}}y\right ) | \sin \left ( \frac{3\pi }{L_{x}}x\right ) \sin \left ( \frac{2\pi }{L_{y}}y\right ) | \sin \left ( \frac{3\pi }{L_{x}}x\right ) \sin \left ( \frac{3\pi }{L_{y}}y\right ) | \cdots |
\vdots | \vdots | \vdots | \vdots | \vdots |
To draw these modes, let us assume that L_{x}=1,L_{y}=1. This gives
n | m=1 | 2 | 3 | 4 |
1 | \sin \left ( \pi x\right ) \sin \left ( \pi y\right ) | \sin \left ( \pi x\right ) \sin \left ( 2\pi y\right ) | \sin \left ( \pi x\right ) \sin \left ( 3\pi y\right ) | \cdots |
2 | \sin \left ( 2\pi x\right ) \sin \left ( \pi y\right ) | \sin \left ( 2\pi x\right ) \sin \left ( 2\pi y\right ) | \sin \left ( 2\pi x\right ) \sin \left ( 3\pi y\right ) | \cdots |
3 | \sin \left ( 3\pi x\right ) \sin \left ( \pi y\right ) | \sin \left ( 3\pi x\right ) \sin \left ( 2\pi y\right ) | \sin \left ( 3\pi x\right ) \sin \left ( 3\pi y\right ) | \cdots |
\vdots | \vdots | \vdots | \vdots | \vdots |
The following is a plot of the above modes for illustrations with the code used to generate these plots.
The following is 3D view of the above modes.
Find the normal modes of an acoustic waves in a hollow sphere of radius R. The wave equation is \nabla ^{2}\psi \left ( r,\theta ,\phi ,t\right ) =\frac{1}{c^{2}}\psi _{tt} With boundary conditions \psi _{r}=0 at r=0 and at r=r_{0}. (I used r_{0} in place of R because wanted to use R\left ( r\right ) for separation of variables).
What is the lowest frequency?
solution
Let \psi \left ( r,\theta ,\phi ,t\right ) =u\left ( r,\theta ,\phi \right ) e^{-i\omega t} Substituting this back in the original PDE gives \nabla ^{2}u\left ( r,\theta ,\phi \right ) +\frac{\omega ^{2}}{c^{2}}u\left ( r,\theta ,\phi \right ) =0 Let k=\frac{\omega }{c} (wave number) and the above becomes\begin{equation} \nabla ^{2}u+k^{2}u=0 \tag{1} \end{equation} The above is called the Helmholtz PDE. In spherical coordinates it becomes \overset{\text{Radial part}}{\overbrace{u_{rr}+\frac{2}{r}u_{r}}}+\overset{\text{Angular part}}{\overbrace{\frac{1}{r^{2}}\left ( \frac{\cos \theta }{\sin \theta }u_{\theta }+u_{\theta \theta }\right ) +\frac{1}{r^{2}\sin ^{2}\theta }u_{\phi \phi }}}+k^{2}u=0 Let u\left ( r,\theta ,\phi \right ) =R\left ( r\right ) \Theta \left ( \theta \right ) \Phi \left ( \phi \right ) and the above becomes R^{\prime \prime }T\Theta \Phi +\frac{2}{r}R^{\prime }T\Theta \Phi +\frac{1}{r^{2}}\left ( \frac{\cos \theta }{\sin \theta }\Theta ^{\prime }RT\Phi +\Theta ^{\prime \prime }RT\Phi \right ) +\frac{1}{r^{2}\sin ^{2}\theta }\Phi ^{\prime \prime }R\Theta T+k^{2}R\Theta T=0 Dividing by R\Theta \Phi \neq 0 gives\begin{align*} \frac{R^{\prime \prime }}{R}+\frac{2}{r}\frac{R^{\prime }}{R}+\frac{1}{r^{2}}\left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) +\frac{1}{r^{2}\sin ^{2}\theta }\frac{\Phi ^{\prime \prime }}{\Phi }+k^{2} & =0\\ r^{2}\sin ^{2}\theta \frac{R^{\prime \prime }}{R}+r^{2}\sin ^{2}\theta \frac{2}{r}\frac{R^{\prime }}{R}+\sin ^{2}\theta \left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) +k^{2}r^{2}\sin ^{2}\theta & =-\frac{\Phi ^{\prime \prime }}{\Phi } \end{align*}
The left side depends only on r,\theta and the right side depends only on \phi . Let the second separation constant be m^{2} and the above becomes\begin{equation} r^{2}\sin ^{2}\theta \frac{R^{\prime \prime }}{R}+r^{2}\sin ^{2}\theta \frac{2}{r}\frac{R^{\prime }}{R}+\sin ^{2}\theta \left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) +k^{2}r^{2}\sin ^{2}\theta =-\frac{\Phi ^{\prime \prime }}{\Phi }=m^{2} \tag{2} \end{equation} Which gives the first angular ODE as\begin{equation} \Phi ^{\prime \prime }+m^{2}\Phi =0 \tag{2A} \end{equation} We now go back to (2) to obtain the rest of the solutions. We now have\begin{align*} r^{2}\sin ^{2}\theta \frac{R^{\prime \prime }}{R}+r^{2}\sin ^{2}\theta \frac{2}{r}\frac{R^{\prime }}{R}+\sin ^{2}\theta \left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) +k^{2}r^{2}\sin ^{2}\theta & =m^{2}\\ k^{2}r^{2}+r^{2}\left ( \frac{R^{\prime \prime }}{R}+\frac{2}{r}\frac{R^{\prime }}{R}\right ) +\left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) & =\frac{m^{2}}{\sin ^{2}\theta }\\ k^{2}r^{2}+r^{2}\left ( \frac{R^{\prime \prime }}{R}+\frac{2}{r}\frac{R^{\prime }}{R}\right ) & =-\left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) +\frac{m^{2}}{\sin ^{2}\theta } \end{align*}
The left side depends on r and the right side depends on \theta only. Let the separation constant be l\left ( l+1\right ) where l is integer which results in \begin{equation} k^{2}r^{2}+r^{2}\left ( \frac{R^{\prime \prime }}{R}+\frac{2}{r}\frac{R^{\prime }}{R}\right ) =-\left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) +\frac{m^{2}}{\sin ^{2}\theta }=l\left ( l+1\right ) \tag{3} \end{equation} Therefore the next angular ODE is\begin{align} -\left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) +\frac{m^{2}}{\sin ^{2}\theta } & =l\left ( l+1\right ) \nonumber \\ -\left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) +\frac{m^{2}}{\sin ^{2}\theta }-l\left ( l+1\right ) & =0\nonumber \\ \left ( \frac{\cos \theta }{\sin \theta }\frac{\Theta ^{\prime }}{\Theta }+\frac{\Theta ^{\prime \prime }}{\Theta }\right ) -\frac{m^{2}}{\sin ^{2}\theta }+l\left ( l+1\right ) & =0\nonumber \\ \Theta ^{\prime \prime }+\frac{\cos \theta }{\sin \theta }\Theta ^{\prime }+\left ( l\left ( l+1\right ) -\frac{m^{2}}{\sin ^{2}\theta }\right ) \Theta & =0 \tag{4} \end{align}
Let z=\cos \theta , then \frac{d\Theta }{d\theta }=\frac{d\Theta }{dz}\frac{dz}{d\theta }=-\frac{d\Theta }{dz}\sin \theta and \begin{align*} \frac{d^{2}\Theta }{d\theta ^{2}} & =\frac{d}{d\theta }\left ( -\frac{d\Theta }{dz}\sin \theta \right ) \\ & =-\frac{d^{2}\Theta }{dz^{2}}\frac{dz}{d\theta }\sin \theta -\frac{d\Theta }{dz}\cos \theta \\ & =\frac{d^{2}\Theta }{dz^{2}}\sin ^{2}\theta -\frac{d\Theta }{dz}\cos \theta \end{align*}
But \sin ^{2}\theta =1-\cos ^{2}\theta =1-z^{2} and the above becomes \frac{d^{2}\Theta }{d\theta ^{2}}=\frac{d^{2}\Theta }{dz^{2}}\left ( 1-z^{2}\right ) -\frac{d\Theta }{dz}z Using these in (4) gives\begin{align} \frac{d^{2}\Theta }{dz^{2}}\left ( 1-z^{2}\right ) -\frac{d\Theta }{dz}z+\frac{z}{\sin \theta }\left ( -\frac{d\Theta }{dz}\sin \theta \right ) +\left ( l\left ( l+1\right ) -\frac{m^{2}}{1-z^{2}}\right ) \Theta \left ( z\right ) & =0\nonumber \\ \left ( 1-z^{2}\right ) \Theta ^{\prime \prime }-2z\Theta ^{\prime }+\left ( l\left ( l+1\right ) -\frac{m^{2}}{1-z^{2}}\right ) \Theta \left ( z\right ) & =0 \tag{3A} \end{align}
And finally, we obtain the final ODE, which is the radial ODE from (3)\begin{align} k^{2}r^{2}+r^{2}\left ( \frac{R^{\prime \prime }}{R}+\frac{2}{r}\frac{R^{\prime }}{R}\right ) & =l\left ( l+1\right ) \nonumber \\ k^{2}r^{2}R+r^{2}\left ( R^{\prime \prime }+\frac{2}{r}R^{\prime }\right ) -l\left ( l+1\right ) R & =0\nonumber \\ r^{2}R^{\prime \prime }+2rR^{\prime }+\left ( k^{2}r^{2}-l\left ( l+1\right ) \right ) R & =0\nonumber \\ R^{\prime \prime }+\frac{2}{r}R^{\prime }+\left ( k^{2}-\frac{l\left ( l+1\right ) }{r^{2}}\right ) R & =0 \tag{4A} \end{align}
In summary we have obtained the following 4 ODE’s to solve (1A,2A,3A,4A)\begin{align} \Phi ^{\prime \prime }+m^{2}\Phi & =0\tag{2A}\\ \left ( 1-z^{2}\right ) \Theta ^{\prime \prime }-2z\Theta ^{\prime }+\left ( l\left ( l+1\right ) -\frac{m^{2}}{1-z^{2}}\right ) \Theta \left ( z\right ) & =0\tag{3A}\\ R^{\prime \prime }+\frac{2}{r}R^{\prime }+\left ( k^{2}-\frac{l\left ( l+1\right ) }{r^{2}}\right ) R & =0 \tag{4A} \end{align}
Solution to (2A) requires m to be integer due to periodicity requirements of solution. The solution is \Phi \left ( \phi \right ) =e^{\pm im\phi }. Equation (3A) is the associated Legendre ODE. Since we are taking l as integer then the solution is known to be \Theta \left ( z\right ) =P_{l}^{m}\left ( z\right ) +Q_{l}^{m}\left ( z\right ) where P_{l}^{m}\left ( z\right ) is called the associated Legendre polynomial and Q_{l}^{m} is the Legendre function of the second kind. Finally (4A) can be converted to Bessel ODE as shown in class notes using the transformation R\left ( r\right ) =\frac{u\left ( r\right ) }{\sqrt{r}} which results in u^{\prime \prime }+\frac{1}{r}u^{\prime }+\left ( k^{2}-\frac{\left ( l+\frac{1}{2}\right ) ^{2}}{r^{2}}\right ) u=0 Which has solution J_{l+\frac{1}{2}}\left ( kr\right ) . The second solution J_{-\left ( l+\frac{1}{2}\right ) }\left ( kr\right ) is rejected since it is not finite at zero and hence makes the solution blow up at center of sphere. Therefore solution to (4A) is \begin{align*} R\left ( r\right ) & =C\sqrt{\frac{\pi }{2kr}}J_{l+\frac{1}{2}}\left ( kr\right ) \\ & =Cj_{l}\left ( kr\right ) \end{align*}
Where C is arbitrary constant. Putting all the above together, then the final solution is \psi \left ( r,\theta ,\phi ,t\right ) =\left \{ \begin{array} [c]{c}e^{-i\omega t}\end{array} \right . \left \{ \begin{array} [c]{c}e^{im\phi }\\ e^{-im\phi }\end{array} \right . \left \{ \begin{array} [c]{c}P_{l}^{m}\left ( \cos \theta \right ) \\ Q_{l}^{m}\left ( \cos \theta \right ) \end{array} \right . \left \{ \begin{array} [c]{c}j_{l}\left ( kr\right ) \end{array} \right . Where j_{l}\left ( kr\right ) are the spherical Bessel functions. Now we need to satisfy the boundary conditions. Since only j_{l}\left ( kr\right ) depends on r, then \psi _{r}=0 at r=0 and at r=r_{0} are equivalent to looking at R^{\prime }\left ( r\right ) =0 at r=0 and r=r_{0}. Therefore we need to find the smallest l,k which satisfy both conditions. This will give the lowest frequency.
I found from DLMF that the series expansion of j_{l}\left ( kr\right ) is\begin{equation} j_{l}\left ( kr\right ) =\frac{\left ( kr\right ) ^{l}}{\left ( 2l+1\right ) !!}\left ( 1-\frac{\left ( kr\right ) ^{2}}{2\left ( 2l+3\right ) }+\frac{\left ( kr\right ) ^{4}}{8\left ( 2l+5\right ) \left ( 2l+3\right ) }+\cdots \right ) \tag{5} \end{equation} Hence for r\rightarrow 0, we can approximate the above as the following by ignoring all higher order terms \lim _{r\rightarrow 0}j_{l}\left ( kr\right ) =\frac{\left ( kr\right ) ^{l}}{\left ( 2l+1\right ) !!} Which means for small r, the derivative is \frac{d}{dr}j_{l}\left ( kr\right ) =\frac{l\left ( kr\right ) ^{l-1}}{\left ( 2l+1\right ) !!} At r=0 then setting \left [ \frac{d}{dr}j_{l}\left ( kr\right ) \right ] _{r\rightarrow 0}=0 is satisfied for all l. Now taking derivative of (5) gives
\frac{d}{dr}j_{l}\left ( kr\right ) =\frac{l\left ( kr\right ) ^{l-1}}{\left ( 2l+1\right ) !!}\left ( 1-\frac{\left ( kr\right ) ^{2}}{2\left ( 2l+3\right ) }+\frac{\left ( kr\right ) ^{4}}{8\left ( 2l+5\right ) \left ( 2l+3\right ) }+\cdots \right ) +\frac{\left ( kr\right ) ^{l}}{\left ( 2l+1\right ) !!}\left ( 1-\frac{2\left ( kr\right ) }{2\left ( 2l+3\right ) }+\frac{4\left ( kr\right ) ^{3}}{8\left ( 2l+5\right ) \left ( 2l+3\right ) }+\cdots \right )
At r=r_{0} the above becomes
\left [ \frac{d}{dr}j_{l}\left ( kr\right ) \right ] _{r\rightarrow r_{0}}=\frac{l\left ( kr_{0}\right ) ^{l-1}}{\left ( 2l+1\right ) !!}\left ( 1-\frac{\left ( kr_{0}\right ) ^{2}}{2\left ( 2l+3\right ) }+\frac{\left ( kr_{0}\right ) ^{4}}{8\left ( 2l+5\right ) \left ( 2l+3\right ) }+\cdots \right ) +\frac{\left ( kr_{0}\right ) ^{l}}{\left ( 2l+1\right ) !!}\left ( 1-\frac{2\left ( kr_{0}\right ) }{2\left ( 2l+3\right ) }+\frac{4\left ( kr_{0}\right ) ^{3}}{8\left ( 2l+5\right ) \left ( 2l+3\right ) }+\cdots \right )
Now we ask, for which values of l is the above zero? If we let l\rightarrow \infty then we obtain\begin{align*} \left [ \frac{d}{dr}j_{l}\left ( kr\right ) \right ] _{\substack{r\rightarrow r_{0}\\l\rightarrow \infty }} & =\lim _{l\rightarrow \infty }\frac{l\left ( kr_{0}\right ) ^{l-1}}{\left ( 2l+1\right ) !!}+\frac{\left ( kr_{0}\right ) ^{l}}{\left ( 2l+1\right ) !!}\\ & =0 \end{align*}
Therefore, to satisfy both \left [ \frac{d}{dr}j_{l}\left ( kr\right ) \right ] _{r\rightarrow 0}=0 and \left [ \frac{d}{dr}j_{l}\left ( kr\right ) \right ] _{r\rightarrow r_{0}}=0 we need l\rightarrow \infty . In other words, a very large integer. The larger l is, the lower the radial frequency. In addition, increasing k while keeping l fixed will increase the frequency. And decreasing k while keeping l fixed decreases the frequency. And for fixed k, increasing l decreases the frequency.
A sphere of radius R is at temperature u=0. At time t=0 it is immersed in a heat bath of temperature u_{0}. What is the temperature distribution u\left ( r,t\right ) as function of time?
solution
Note: I Used u\left ( r,t\right ) instead of T\left ( r,t\right ) as the dependent variable to allow using T\left ( t\right ) for separation of variables without confusing it with the original T\left ( r,t\right ) .
The PDE specification is, solve for u\left ( r,t\right ) u_{t}=k\nabla ^{2}u\qquad t>0,0<r<R With initial conditions u\left ( r,0\right ) =0 And boundary conditions\begin{align*} u\left ( R,t\right ) & =u_{0}\\ \left \vert u\left ( 0,t\right ) \right \vert & <\infty \end{align*}
Where the second B.C. above means the temperature u is bounded at origin (center of sphere). In spherical coordinates, the PDE becomes (There are no dependency on \theta ,\phi due to symmetry), and only radial dependency.\begin{equation} \frac{1}{k}u_{t}=\frac{1}{r}\left ( ru\right ) _{rr}\tag{1} \end{equation} To simplify the solution, let U\left ( r,t\right ) =ru\left ( r,t\right ) And we obtain a new PDE\begin{equation} \frac{1}{k}U_{t}=U_{rr}\tag{2} \end{equation} And the boundary conditions u\left ( R,t\right ) =u_{0} becomes U\left ( R,t\right ) =Ru_{0} and the initial conditions becomes U\left ( r,0\right ) =0. So we will solve (2) and not (1). But since the boundary conditions are not homogenous, we can not use separation of variables. We introduce a reference function w\left ( r\right ) which need to satisfy the nonhomogeneous boundary conditions only. Let w\left ( r\right ) =Br. When r=R then Ru_{0}=BR or B=u_{0} When r=0 then w=0 which is bounded. Hence w\left ( r\right ) =u_{0}r Therefore, the solution now can be written as\begin{equation} U\left ( r,t\right ) =v\left ( r,t\right ) +u_{0}r\tag{3} \end{equation} Where v\left ( r,t\right ) now satisfies the PDE but with homogenous B.C. Substituting (3) into (2) gives\begin{align} v_{t} & =k\frac{\partial ^{2}}{\partial r^{2}}\left ( v\left ( r,t\right ) +u_{0}r\right ) \nonumber \\ v_{t} & =kv_{rr}\left ( r,t\right ) \tag{4} \end{align}
We need to solve the above but with homogenous boundary conditions\begin{align*} v\left ( R,t\right ) & =0\\ \left \vert v\left ( 0,t\right ) \right \vert & <\infty \end{align*}
This is standard PDE, who can be solved by separation of variables. let v=F\left ( r\right ) T\left ( t\right ) , hence (4) becomes\begin{align*} T^{\prime }F & =kF^{\prime \prime }T\\ k\frac{T^{\prime }}{T} & =\frac{F^{\prime \prime }}{F}=-\lambda ^{2} \end{align*}
Which gives F^{\prime \prime }+\lambda ^{2}F=0 Due to boundary conditions only \lambda >0 is eigenvalues. Hence solution is F\left ( r\right ) =A\cos \left ( \lambda r\right ) +B\sin \left ( \lambda r\right ) At r=0, since bounded, say 0, then we can take A=0, leaving the solution F\left ( r\right ) =B\sin \left ( \lambda r\right ) At r=R 0=B\sin \left ( \lambda R\right ) For nontrivial solution\begin{align*} \lambda R & =n\pi \qquad n=1,2,3,\cdots \\ \lambda _{n} & =\frac{n\pi }{R} \end{align*}
Hence eigenfunctions are F_{n}\left ( r\right ) =\sin \left ( \frac{n\pi }{R}r\right ) \qquad n=1,2,3,\cdots The time ODE is therefore T^{\prime }+\lambda ^{2}kT=0 with solution T_{n}\left ( t\right ) =A_{n}e^{-\left ( \frac{n\pi }{R}\right ) ^{2}kt}. Hence the solution to (4) is v\left ( r,t\right ) =\sum _{n=1}^{\infty }A_{n}e^{-\left ( \frac{n\pi }{R}\right ) ^{2}kt}\sin \left ( \frac{n\pi }{R}r\right ) Therefore from (3) U\left ( r,t\right ) =\left ( \sum _{n=1}^{\infty }A_{n}e^{-\left ( \frac{n\pi }{R}\right ) ^{2}kt}\sin \left ( \frac{n\pi }{R}r\right ) \right ) +u_{0}r But U\left ( r,t\right ) =ru\left ( r,t\right ) , hence\begin{equation} u\left ( r,t\right ) =\left ( \frac{1}{r}\sum _{n=1}^{\infty }A_{n}e^{-\left ( \frac{n\pi }{R}\right ) ^{2}kt}\sin \left ( \frac{n\pi }{R}r\right ) \right ) +u_{0}\tag{5} \end{equation} Now we find A_{n} from initial conditions. At t=0\begin{align*} 0 & =u_{0}+\frac{1}{r}\sum _{n=1}^{\infty }A_{n}\sin \left ( \frac{n\pi }{R}r\right ) \\ -ru_{0} & =\sum _{n=1}^{\infty }A_{n}\sin \left ( \frac{n\pi }{R}r\right ) \end{align*}
Therefore A_{n} are the Fourier series coefficients of -ru_{0}\begin{align*} \frac{R}{2}A_{n} & =-\int _{0}^{R}ru_{0}\sin \left ( \frac{n\pi }{R}r\right ) dr\\ A_{n} & =-\frac{2u_{0}}{R}\int _{0}^{R}r\sin \left ( \frac{n\pi }{R}r\right ) dr\\ & =-\frac{2u_{0}}{R}\left ( -1\right ) ^{n+1}\frac{R^{2}}{n\pi }\\ & =\left ( -1\right ) ^{n}\frac{2R}{n\pi }u_{0} \end{align*}
Hence the solution (5) becomes\begin{align} u\left ( r,t\right ) & =u_{0}+u_{0}\frac{2R}{r\pi }\sum _{n=1}^{\infty }\left ( -1\right ) ^{n}\frac{1}{n}e^{-k\left ( \frac{n\pi }{R}\right ) ^{2}t}\sin \left ( \frac{n\pi }{R}r\right ) \nonumber \\ & =u_{0}\left ( 1+\frac{2R}{r\pi }\sum _{n=1}^{\infty }\left ( -1\right ) ^{n}\frac{1}{n}e^{-k\left ( \frac{n\pi }{R}\right ) ^{2}t}\sin \left ( \frac{n\pi }{R}r\right ) \right ) \tag{7} \end{align}
Verification of solution
Verification that (7) satisfies the PDE u_{t}=k\nabla ^{2}u. Taking time derivative of (7) gives\begin{equation} u_{t}=-u_{0}\frac{2R}{r\pi }k\sum _{n=1}^{\infty }\left ( -1\right ) ^{n}\frac{1}{n}\left ( \frac{n\pi }{R}\right ) ^{2}e^{-k\left ( \frac{n\pi }{R}\right ) ^{2}t}\sin \left ( \frac{n\pi }{R}r\right ) \tag{8} \end{equation} And taking space derivatives of (7) gives\begin{align*} u_{x} & =u_{0}\frac{2R}{r\pi }\sum _{n=1}^{\infty }\left ( -1\right ) ^{n}\frac{1}{n}e^{-k\left ( \frac{n\pi }{R}\right ) ^{2}t}\frac{n\pi }{R}\cos \left ( \frac{n\pi }{R}r\right ) \\ u_{xx} & =-u_{0}\frac{2R}{r\pi }\sum _{n=1}^{\infty }\left ( -1\right ) ^{n}\frac{1}{n}e^{-k\left ( \frac{n\pi }{R}\right ) ^{2}t}\left ( \frac{n\pi }{R}\right ) ^{2}\sin \left ( \frac{n\pi }{R}r\right ) \end{align*}
Hence ku_{xx} becomes\begin{equation} ku_{xx}=-u_{0}\frac{2R}{r\pi }k\sum _{n=1}^{\infty }\left ( -1\right ) ^{n}\frac{1}{n}e^{-k\left ( \frac{n\pi }{R}\right ) ^{2}t}\left ( \frac{n\pi }{R}\right ) ^{2}\sin \left ( \frac{n\pi }{R}r\right ) \tag{9} \end{equation} Comparing (8) and (9) shows they are the same expressions.
Verification that (7) satisfies the boundary condition.
When r=R, therefore (7) gives, when replacing r by R\begin{align*} u\left ( R,t\right ) & =u_{0}\left ( 1+\frac{2R}{R\pi }\sum _{n=1}^{\infty }\left ( -1\right ) ^{n}\frac{1}{n}e^{-k\left ( \frac{n\pi }{R}\right ) ^{2}t}\sin \left ( \frac{n\pi }{R}R\right ) \right ) \\ & =u_{0}\left ( 1+\frac{2R}{R\pi }\sum _{n=1}^{\infty }\left ( -1\right ) ^{n}\frac{1}{n}e^{-k\left ( \frac{n\pi }{R}\right ) ^{2}t}\sin \left ( n\pi \right ) \right ) \\ & =u_{0}\left ( 1+0\right ) \\ & =u_{0} \end{align*}
But n is integer. Hence \sin \left ( n\pi \right ) =0 for all n. And the above becomes\begin{align*} u\left ( R,t\right ) & =u_{0}\left ( 1+0\right ) \\ & =u_{0} \end{align*}
Verified.
Verification that (7) satisfies the initial conditions u\left ( r,0\right ) =0 for r<R.
At t=0 (7) becomes\begin{align*} u\left ( r,0\right ) & =u_{0}\left ( 1+\frac{2R}{r\pi }\sum _{n=1}^{\infty }\left ( -1\right ) ^{n}\frac{1}{n}\sin \left ( \frac{n\pi }{R}r\right ) \right ) \\ & =u_{0}+\frac{2R}{r\pi }u_{0}\sum _{n=1}^{\infty }\frac{\left ( -1\right ) ^{n}}{n}\sin \left ( \frac{n\pi }{R}r\right ) \\ & =u_{0}+\frac{2R}{r\pi }u_{0}\left ( -\sin \left ( \frac{\pi }{R}r\right ) +\frac{1}{2}\sin \left ( \frac{2\pi }{R}r\right ) -\frac{1}{3}\sin \left ( \frac{3\pi }{R}r\right ) +\frac{1}{4}\sin \left ( \frac{4\pi }{R}r\right ) -\cdots \right ) \end{align*}
I could not simplify the above by hand, but using the computer, I verified numerically it is zero for 0<r<R for a given R and given u_{0}.
Consider the Helmholtz equation
\begin{equation} \nabla ^{2}u\left ( r,\theta \right ) +k^{2}u\left ( r,\theta \right ) =0 \tag{1} \end{equation}
inside the circle r=r_{0} with the boundary condition u\left ( r_{0},\theta \right ) =f\left ( \theta \right ) . The solution can be written in the form u\left ( r,\theta \right ) =\int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) G\left ( r,\theta ;\theta ^{\prime }\right ) d\theta ^{\prime }. Find the Green function G.
solution
I will solve (1) directly and then compare the solution obtain to u\left ( r,\theta \right ) =\int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) G\left ( r,\theta ;\theta ^{\prime }\right ) d\theta ^{\prime } in order to read off the Green function expression. (1) in polar coordinates becomes
u_{rr}+\frac{1}{r}u_{r}+\frac{1}{r^{2}}u_{\theta \theta }+k^{2}u=0
Writing u\left ( r,\theta \right ) =R\left ( r\right ) \Theta \left ( \theta \right ) , the above PDE becomes\begin{align*} R^{\prime \prime }\Theta +\frac{1}{r}R^{\prime }\Theta +\frac{1}{r^{2}}\Theta ^{\prime \prime }R+k^{2}R\Theta & =0\\ \frac{R^{\prime \prime }}{R}+\frac{1}{r}\frac{R^{\prime }}{R}+\frac{1}{r^{2}}\frac{\Theta ^{\prime \prime }}{\Theta }+k^{2} & =0\\ r^{2}\frac{R^{\prime \prime }}{R}+r\frac{R^{\prime }}{R}+r^{2}k^{2} & =-\frac{\Theta ^{\prime \prime }}{\Theta }=m \end{align*}
Where m is the separation constant. The eigenvalue problem is taken as \Theta ^{\prime \prime }+m\Theta =0 Due to periodicity of the solution on the disk, then \Theta \left ( -\pi \right ) =\Theta \left ( \pi \right ) and \Theta ^{\prime }\left ( -\pi \right ) =\Theta ^{\prime }\left ( \pi \right ) . These boundary conditions restrict m to only positive integer values. Hence let m=n^{2} and the solution to the above becomes \Theta _{\alpha }\left ( \theta \right ) =A_{n}\cos \left ( n\theta \right ) +B_{n}\sin \left ( n\theta \right ) Now the radial ODE is \begin{align*} r^{2}\frac{R^{\prime \prime }}{R}+r\frac{R^{\prime }}{R}+r^{2}k^{2} & =\alpha ^{2}\\ r^{2}R^{\prime \prime }+rR^{\prime }+\left ( r^{2}k^{2}-n^{2}\right ) R & =0\\ R^{\prime \prime }+\frac{1}{r}R^{\prime }+\left ( k^{2}-\frac{n^{2}}{r^{2}}\right ) R & =0 \end{align*}
This is Bessel ODE whose solutions are (since n are integers) is R_{\alpha }\left ( r\right ) =C_{n}J_{n}\left ( kr\right ) +E_{n}Y_{n}\left ( kr\right ) But Y_{n}\left ( kr\right ) blows up at r=0, hence it is rejected leaving solution R_{n}\left ( r\right ) =C_{n}J_{n}\left ( kr\right ) . Hence the final solution is\begin{equation} u\left ( r,\theta \right ) =\sum _{m=1}^{\infty }\left ( A_{n}\cos \left ( n\theta \right ) +B_{n}\sin \left ( n\theta \right ) \right ) J_{n}\left ( kr\right ) \tag{2} \end{equation} Where the constant C_{n} is merged with the other two constants. Now, at r=r_{0} we are told that u\left ( r_{0},\theta \right ) =f\left ( \theta \right ) . Hence the above becomes f\left ( \theta \right ) =\sum _{m=1}^{\infty }\left ( A_{n}\cos \left ( n\theta \right ) +B_{n}\sin \left ( n\theta \right ) \right ) J_{n}\left ( kr_{0}\right ) By orthogonality of \cos \left ( n\theta \right ) ,\sin \left ( n\theta \right ) we find the Fourier cosine and Fourier sine coefficients A_{n},B_{n} as\begin{align*} A_{n}J_{n}\left ( kr_{0}\right ) \frac{1}{\pi } & =\int _{0}^{2\pi }f\left ( \theta \right ) \cos \left ( n\theta \right ) d\theta \\ B_{n}J_{n}\left ( kr_{0}\right ) \frac{1}{\pi } & =\int _{0}^{2\pi }f\left ( \theta \right ) \sin \left ( n\theta \right ) d\theta \end{align*}
Substituting the above back into the solution found in (2) results in\begin{align} u\left ( r,\theta \right ) & =\sum _{m=1}^{\infty }\left [ \left ( \frac{\pi }{J_{n}\left ( kr_{0}\right ) }\int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) \cos \left ( n\theta ^{\prime }\right ) d\theta ^{\prime }\right ) \cos \left ( n\theta \right ) +\left ( \frac{\pi }{J_{n}\left ( kr_{0}\right ) }\int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) \sin \left ( n\theta ^{\prime }\right ) d\theta ^{\prime }\right ) \sin \left ( n\theta \right ) \right ] J_{n}\left ( kr\right ) \nonumber \\ & =\sum _{m=1}^{\infty }\frac{\pi }{J_{n}\left ( kr_{0}\right ) }\left ( \int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) \cos \left ( n\theta ^{\prime }\right ) \cos \left ( n\theta \right ) d\theta ^{\prime }+\int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) \sin \left ( n\theta ^{\prime }\right ) \sin \left ( n\theta \right ) d\theta ^{\prime }\right ) J_{n}\left ( kr\right ) \tag{3} \end{align}
Using trig relations\begin{align*} \cos A\cos B & =\frac{1}{2}\left ( \cos \left ( A+B\right ) +\cos \left ( A-B\right ) \right ) \\ \sin A\sin B & =\frac{1}{2}\left ( \cos \left ( A-B\right ) -\cos \left ( A+B\right ) \right ) \end{align*}
Then (3) becomes u\left ( r,\theta \right ) =\sum _{m=1}^{\infty }\frac{\pi }{2J_{n}\left ( kr_{0}\right ) }\left ( \int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) \left ( \cos \left ( n\left ( \theta ^{\prime }+\theta \right ) \right ) +\cos \left ( n\left ( \theta ^{\prime }-\theta \right ) \right ) \right ) d\theta ^{\prime }+\int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) \left ( \cos \left ( n\left ( \theta ^{\prime }-\theta \right ) \right ) -\cos \left ( n\left ( \theta ^{\prime }+\theta \right ) \right ) \right ) d\theta ^{\prime }\right ) J_{n}\left ( kr\right )
Which is simplified to, after combining both integrals to one\begin{align*} u\left ( r,\theta \right ) & =\sum _{m=1}^{\infty }\frac{\pi }{2J_{n}\left ( kr_{0}\right ) }\left ( \int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) \left ( \cos \left ( n\left ( \theta ^{\prime }+\theta \right ) \right ) +\cos \left ( n\left ( \theta ^{\prime }-\theta \right ) \right ) +\cos \left ( n\left ( \theta ^{\prime }-\theta \right ) \right ) -\cos n\left ( \theta ^{\prime }+\theta \right ) \right ) d\theta ^{\prime }\right ) J_{n}\left ( kr\right ) \\ & =\sum _{m=1}^{\infty }\frac{\pi }{2J_{n}\left ( kr_{0}\right ) }\left [ \int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) 2\cos \left ( \theta ^{\prime }-\theta \right ) d\theta ^{\prime }\right ] J_{n}\left ( kr\right ) \\ & =\sum _{m=1}^{\infty }\left [ \int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) \frac{\pi }{J_{n}\left ( kr_{0}\right ) }\cos \left ( \theta ^{\prime }-\theta \right ) d\theta ^{\prime }\right ] J_{n}\left ( kr\right ) \end{align*}
Exchanging integration with summation gives u\left ( r,\theta \right ) =\int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) \left ( \sum _{m=1}^{\infty }\frac{\pi }{J_{n}\left ( kr_{0}\right ) }\cos \left ( \theta ^{\prime }-\theta \right ) J_{n}\left ( kr\right ) \right ) d\theta ^{\prime } Comparing the above to u\left ( r,\theta \right ) =\int _{0}^{2\pi }f\left ( \theta ^{\prime }\right ) G\left ( r,\theta ;\theta ^{\prime }\right ) d\theta ^{\prime } Shows that Green function is G\left ( r,\theta ;\theta ^{\prime }\right ) =\sum _{m=1}^{\infty }\frac{\pi }{J_{n}\left ( kr_{0}\right ) }\cos \left ( \theta ^{\prime }-\theta \right ) J_{n}\left ( kr\right ) Where r_{0} is radius of disk. It is symmetric in \theta as expected.