skip to content
Aymen Hafeez

Recall from the first part of this post a diffusion equation was derived to approximate the movement of neutrons in a nuclear reactor, which in two dimensions was,

d2ϕdr2+1rdϕdr+B2ϕ=0. \frac{\text{d}^2 \phi}{\text{d} r^2} + \frac{1}{r} \frac{\text{d} \phi}{\text{d} r} + B^{2}\phi = 0.

The main reason I was interested in solving this was due to the presence of non-constant coefficients, which was not a general form I had come across. Solving second order ODE’s of this form can be done using the Frobenius method, which produces an infinite series solution. We begin by defining an expression in the form which we want the solution to be:

ϕ=∑n=0∞anrn+k\begin{equation*} \begin{aligned} \phi = \sum_{n=0}^{\infty}a_nr^{n+k} \end{aligned} \end{equation*}

The first and second derivatives are, therefore,

dϕdr=∑n=0∞(n+k)anrn+k−1\begin{equation*} \begin{aligned} \frac{\text{d} \phi}{\text{d} r} = \sum_{n=0}^{\infty}(n + k)a_n r^{n+k-1} \end{aligned} \end{equation*} d2ϕdr2=∑n=0∞(n+k−1)(n+k)anrn+k−2\begin{equation*} \begin{aligned} \frac{\text{d}^2 \phi}{\text{d} r^2} = \sum_{n=0}^{\infty}(n+k-1)(n + k)a_n r^{n+k-2} \end{aligned} \end{equation*}

Substituting the derivatives back into the ordinal expression gives

∑n=0∞(n+k−1)(n+k)anrn+k−2+1r∑n=0∞(n+k)anrn+k−1+B2∑n=0∞anrn+k=0\begin{equation*} \begin{aligned} \sum_{n=0}^{\infty}(n+k-1)(n+k)a_nr^{n+k-2} + \frac{1}{r}\sum_{n=0}^{\infty}(n + k)a_n r^{n+k-1} + B^2\sum_{n=0}^{\infty}a_nr^{n+k} = 0 \\ \end{aligned} \end{equation*}

Ideally we would want each rr term to be raised to same power. Expanding out the first two terms of the rn+kr^{n+k} expression, followed by some manipulation, gives:

∑n=0∞(n+k−1)(n+k)anrn+k−2+1r∑n=0∞(n+k)anrn+k−1+B2∑n=2∞anrn+k−2=0\begin{equation*} \begin{aligned} \sum_{n=0}^{\infty}(n+k-1)(n+k)a_nr^{n+k-2} + \frac{1}{r}\sum_{n=0}^{\infty}(n + k)a_n r^{n+k-1} + B^2\sum_{n=2}^{\infty}a_nr^{n+k-2} = 0 \\ \end{aligned} \end{equation*} ∑n=2∞(n+k−1)(n+k)anrn+k−2+1r∑n=2∞(n+k)anrn+k−1+B2∑n=2∞an−2rn+k−2+k(k−1)a0rk−2+k(k+1)a1rk−1+1r(ka0rk−1+(k+1)a1rk)=0\begin{equation*} \begin{aligned} \sum_{n=2}^{\infty} (n+k-1)(n+k)a_{n}r^{n+k-2} + \frac{1}{r}\sum_{n=2}^{\infty} (n+k)a_{n}r^{n+k-1} + B^{2}\sum_{n=2}^{\infty} a_{n-2}r^{n+k-2} \\ + k(k-1)a_{0}r^{k-2} + k(k+1)a_{1}r^{k-1} + \frac{1}{r}\left( ka_{0}r^{k-1}+(k+1)a_{1}r^{k} \right) = 0 \nonumber \end{aligned} \end{equation*} ∑n=2∞(n+k−1)(n+k)anrn+k−2+1r∑n=2∞(n+k)anrn+k−1+B2∑n=2∞an−2rn+k−2+k(k−1)a0rk−2+k(k−1)a1rk−1+ka0rk−2+(k+1)a1rk−1=0\begin{equation*} \begin{aligned} \sum_{n=2}^{\infty} (n+k-1)(n+k)a_{n}r^{n+k-2} + \frac{1}{r}\sum_{n=2}^{\infty} (n+k)a_{n}r^{n+k-1} + B^{2}\sum_{n=2}^{\infty} a_{n-2}r^{n+k-2} \\ + k(k-1)a_{0}r^{k-2} + k(k-1)a_{1}r^{k-1} + ka_{0}r^{k-2} + (k+1)a_{1}r^{k-1} = 0 \end{aligned} \end{equation*}

We can equate the coefficients on the right side of the equation to zero. For the rk−2r^{k-2} terms we get

k(k−1)a0+ka0=0a0[k(k−1)+k]=0  ⟹  k=0 ∴ a0=0\begin{equation*} \begin{aligned} k(k - 1)a_{0} + k a_{0} = 0 \\ a_{0}\left[ k(k - 1) + k \right] = 0 \\ \implies k = 0 \ \therefore \ a_{0} = 0 \end{aligned} \end{equation*}

Doing the same with the rk−1r^{k-1} terms:

k(k+1)a1+(k+1)a1=0a1[k(k+1)+(k+1)]=0a1(k+1)(k+1)=0a1=0\begin{equation*} \begin{aligned} k(k + 1)a_{1} + (k + 1)a_{1} = 0 \\ a_{1}\left[ k(k + 1) + (k + 1) \right] = 0 \\ a_{1}(k + 1)(k + 1) = 0 \\ a_{1} = 0 \\ \end{aligned} \end{equation*}

The previous expression, therefore, becomes

∑n=2∞(n+k)(n+k−1)anrn+i−2+∑n=2∞(n+k)anrn+k−2+B2∑n=2∞an−2rn+k−2=0 \sum_{n=2}^{\infty} (n + k)(n + k - 1)a_{n}r^{n+i-2} + \sum_{n=2}^{\infty} (n + k)a_{n}r^{n+k-2} + B^{2}\sum_{n=2}^{\infty} a_{n-2}r^{n+k-2} = 0 ∑n=2∞{[(n+k)(n+k−1)+(n+k)]an+B2an−2}rn+k−2=0 \sum_{n=2}^{\infty} \left\{ \left[ (n + k)(n + k - 1) + (n + k) \right] a_{n} + B^{2}a_{n-2}\right\} r^{n+k-2} = 0

We can see from this that the coefficients of rn+k−2r^{n+k-2} must equate to 0, and so,

[(n+k)(n+k−1)+(n+k)]an+B2an−2=0\begin{equation*} \begin{aligned} \left[ (n + k)(n + k - 1) + (n + k) \right] a_{n} + B^{2}a_{n-2} = 0 \end{aligned} \end{equation*} an(k)=−B2an−2(n+k)(n+k−1)+(n+k)=−B2an−2(n+k)(n+k−1+1)=−B2an−2(n+k)2\begin{equation*} \begin{aligned} a_{n}(k) &= -\frac{B^{2}a_{n-2}}{(n + k)(n + k - 1) + (n + k)} \\ &= - \frac{B^{2}a_{n-2}}{(n + k)(n + k - 1 + 1)} \\ &= - \frac{B^{2}a_{n-2}}{(n + k)^{2}} \\ \end{aligned} \end{equation*}

If we now let k=0k = 0,

an(0)=−B2an−2n2 a_{n}(0) = - \frac{B^{2}a_{n-2}}{n^{2}}

We know that a1=0a_{1} = 0, so we can say that

a1=a3=a5=a2n+1=0 a_{1} = a_{3} = a_{5} = a_{2n+1} = 0

For k=0k = 0, let n=2mn = 2m:

a2m(0=−B2a2m−2(2m2 a_{2m}(0 = - \frac{B^{2}a_{2m-2}}{(2m^{2}}

We can see that, for even values of mm,

a2=−B2a022 a_{2} = - \frac{B^{2}a_{0}}{2^{2}} a4=−B2a222⋅22=−B2a022⋅−B222⋅22=B4a024⋅22=B4a024(22⋅12) a_{4} = \frac{-B^{2}a_{2}}{2^{2} \cdot 2^{2}} = \frac{-B^{2}a_{0}}{2^{2}} \cdot \frac{-B^{2}}{2^{2}\cdot 2^{2}} = \frac{B^{4}a_{0}}{2^{4}\cdot 2^{2}} = \frac{B^{4}a_{0}}{2^{4}(2^{2} \cdot 1^{2})} a6=−B2a42232=−B22232⋅B4a024(2⋅1)2=−B6a026(3⋅2⋅1)2a_{6} = \frac{-B^{2}a_{4}}{2^{2}3^{2}} = \frac{-B^{2}}{2^{2}3^{2}} \cdot \frac{B^{4}a_{0}}{2^{4}(2 \cdot 1)^{2}} = \frac{-B^{6}a_{0}}{2^{6}(3 \cdot 2 \cdot 1)^{2}} a8=−B2a62242=−B22242⋅−B6a026(3⋅2⋅1)2=B8a028(4⋅3⋅2⋅1)2 a_{8} = \frac{-B^{2}a_{6}}{2^{2}4^{2}} = \frac{-B^{2}}{2^{2}4^{2}} \cdot \frac{-B^{6}a_{0}}{2^{6}(3 \cdot 2 \cdot 1)^{2}} = \frac{B^{8}a_{0}}{2^{8}(4 \cdot 3 \cdot 2 \cdot 1)^{2}}

And so, in general we can say that

a2m=(−1)mB2ma022m(m!)2 a_{2m} = \frac{(-1)^{m}B^{2m}a_{0}}{2^{2m}(m!)^{2}}

We can put this back into out original expression for ϕ\phi to get a solution

ϕ1=a0∑n=0∞(−1)mB2m22m(m!)2r2m \phi_{1} = a_{0} \sum_{n=0}^{\infty} \frac{(-1)^{m}B^{2m}}{2^{2m}(m!)^{2}}r^{2m}

This is known as the Bessel function of the first kind (i.e. for k=0k = 0 ).

ϕ1=J0(Br) \phi_{1} = \mathcal{J}_{0} (Br)