Neutron Diffusion Equation/Spherical Geometry Source Problem

[deebrwnfn11]

Homework Statement

Solve for the flux distribution using the 1D neutron diffusion equation in a finite sphere for a uniformly distributed source emitting S₀ neutrons/cc-sec.

My problem right now is that I can't figure out the boundary conditions for this problem. We usually work with point sources in infinite domains when working with spherical geometry, so I am unfamiliar with setting up boundary conditions for finite spheres.

Homework Equations

The governing differential equation is:
<br /> \frac{1}{r^2} \frac{d}{dr} r^2 \frac{d\phi}{dr} - \frac{1}{L^2}\phi(r) = \frac{-S_0}{D}<br />

with a boundary condition:
<br /> \phi(r=R) = 0<br />
This problem is going to be turned into a computer code, so we were told not to use extrapolated boundary conditions. I just need help finding the boundary condition for r = 0

As I mentioned above, we've never worked with finite spheres before so I am not 100% certain that the solution involves sin and cos. Maybe sinh and cosh?
<br /> \phi(r) = \frac{C_1}{r}sin(\frac{r}{L}) + \frac{C_2}{r}cos(\frac{r}{L}) + \frac{S_0L^2}{D}<br />

The Attempt at a Solution

Physically, I know that the flux profile will be flat at the center of the sphere. I can't impose that the derivative equals zero because that would lead to a zero as the denominator in the coefficient terms. I've considered using a limit as r → 0, but that would cause the coefficient terms to shoot up to infinity. I don't think I can say that the neutron current at the center is equal to some value because there is a distributed source.

I am completely lost in finding this boundary condition. I think what I have written so far for the equations is right. If someone could point me in the right direction, or enlighten me of a mistake I have made I would be forever grateful.

 

[SteamKing]

Why do you need a boundary condition for r = 0? Isn't r = 0 within the sphere? Why isn't the BC Phi = 0 at r = R sufficient?

Obviously, you can substitute your trial function Phi(r) into the differential equation and see if the DE is satisfied. Have you tried this?

 

[deebrwnfn11]

I did check my trial function by plugging it into the original ODE and it didn't take me long to figure out that it wasn't going to satisfy it. So I dusted off the old Differential equations textbook and found out what I was doing wrong. The actual solution is:

<br /> \phi(r) = \frac{C_1}{r}sinh(\frac{r}{L}) + \frac{C_2}{r}cosh(\frac{r}{L}) + \frac{S_0L^2}{D}<br />

with the boundary conditions \phi(r=R) = 0 and \lim_{r \to 0} \phi(r) = finite

I reconsidered using the limit boundary condition after re-remembering L'hopital's rule. So my flux (with coefficients solved for) is:

<br /> \phi(r) = \frac{-S_0L^2R}{rDsinh(\frac{R}{L})}sinh(\frac{r}{L}) + \frac{S_0L^2}{D}<br />

Thank you for suggesting that I check my solution, I often forget to do that when I'm solving differential equations.