# Numerically Solving a Second Order Nonlinear ODE

Okay, I have this not so pretty 2nd order non-linear ODE I should be able to solve numerically.

$$f''(R) + \frac{2}{R} f'(R)=\frac{0.7}{R} \left( \frac{1}{\sqrt{f(R)}} - \frac{0.3}{\sqrt{1-f(R)}} \right),$$

$$f(1)=1.$$

The function around the origin is behaving very wildly.

I was thinking of breaking this guy up into a system of two first order ODE's and then solve, but I have no idea how to set this up. What method should I use to set up the system of ODE's?

If there is some other method rather than numerically solving a system of differential equations, please feel welcome to share. Thanks.

-
Is $f(1)=1$ really what you want? This will give division by zero on the right-hand side! Also, do you really have a boundary value problem (conditions at two different points $R=0$ and $R=1$), or did you mean an initial value problem with two conditions at one point (say $R=0$)? An IVP is relatively straightforward to integrate numerically, but a BVP is harder. –  Hans Lundmark Oct 16 '10 at 13:39
I am trying to replicate a graph of the solution from an old paper. I can show a picture of what the graph looks like but I am not sure how. Anyway the solution is for R between 0 and 1. At f(1)=1 at least from what the graphs shows. My other boundary condition is actually wrong (I will edit the main post), but f seems to be approaching positive infinity as R approached 0. –  Shagster_84 Oct 16 '10 at 21:43
Well I am not able to post the graph yet. I don't enough reputation points (i need 10). I might post a link later. –  Shagster_84 Oct 16 '10 at 21:54
+1 to give you the reputation for posting the graph. –  whuber Oct 17 '10 at 15:04
More than posting the graph, where did this graph come from, and how did you construct your DE from said graph? (+1 so you can post the graph). –  Ｊ. Ｍ. Oct 17 '10 at 16:20

The general method to reinterpret a higher-order ODE as a system of first order ODEs is to regard the derivatives of the unknown function as additional unknown functions. In your case, regard $f'$ as a new function $g$. Then the system of first order ODEs consists of two equations, the first being the original equation with $f'$ replaced by $g$ and $f''$ replaced by $g'$, and the second equation being $f'=g$.

-
Okay in that case I will end up with: g'(R)+(2/R)g(R)=(.7/R)((1/sqrt(h))-((0.3)/sqrt(1-h))), g=h' by letting f=h and f'=g. How will this effect the one boundary condition f(1)=1 I knew before the substitution? –  Shagster_84 Oct 16 '10 at 22:45
You don't need h. You just say f'=g and write the equation you have in terms of g',g, and f. It looks like yours with f instead of h. You still use your f(1)=1 (which divides by zero as stated before), but need a second boundary condition. –  Ross Millikan Oct 17 '10 at 0:52
Okay, could I use an approximation? –  Shagster_84 Oct 19 '10 at 15:39
For second-order ODEs you need two boundary conditions, otherwise you end up with a family of solutions (unless you're trying to pick a specific member using, say, optimization). Typically, you'll need a b.c. on f', which, after the transformation mentioned above, gives you a b.c. on g. To solve numerically the (stiff) ODE, the first thing to try is a Runge-Kutta method such as ode45 in Matlab. –  Dominique Oct 25 '11 at 22:51