# Solving The Lane-Emden Equation

## Introduction

Most of theoretical astrophysics involves the numerical solution of various differential equations, both ordinary DEs (ODEs) and partial DEs (PDEs). Our aim here is to solve the Lane-Emden equation. This is a second order ODE.

## The Lane Emden Equation

(NOTE: Sometimes the equations does not render upon page load. You might want to refresh the page for the equations to display properly.)

The dimensionless Lane-Emden equation is:

$$\dfrac{1}{{\xi}^2}\dfrac{d}{d\xi}{\xi}^2\dfrac{d\theta}{d\xi} = - {\theta}^n$$

Here $\xi$ is the scaled radius, $\theta$ is the scaled temperature, and $n$ is the index.

Where the initial conditions at $\xi = 0$ are, $\theta = 1$ and $\dfrac{d\theta}{d\xi} = 0$.

Now we will be looking at how to approach a solution to such a second order ODE.

We first need to write the 2nd order Lane-Emden equation as two first order equations. In the code we will use $t$ for $\xi$, $y$ for $\theta$ and $z$ for $\dfrac{d\theta}{d\xi}$. The Lane-Emden equation is then:

$$\dfrac{d}{dt}({t}^2\dfrac{dy}{dt}) = - {t}^2{y}^n$$

We write this as two first order equations by setting $z = \dfrac{dy}{dt}$. The two equations thus become:

$\dfrac{dy}{dt} = z$
$\dfrac{dz}{dt} = -\dfrac{2}{t}z - {y}^n$

and the initial conditions become:

$y(0) = 1$ and $z(0) = 0$

## Solution using Euler Method

I will update the contents later.

## Solution using Runge-Kutta Method

### The Runge-Kutta Method

The Euler method is a particular a class of techniques known as "Runge-Kutta" method, which have the general form:
$$y_{i+1} = y_{i} + h\phi$$
where, $\phi$ is some approximation to the slope. In Euler Method $\phi = k_{1} = f(x_{i},y_{i})$.

The most used Runge-Kutta method is the 4th order Runge-Kutta Method (hereafter RK4), which has a global error of order $h^{4}$. In this method we take:

$k_{1} = f(x_{i},y_{i})$
$k_{2} = f(x_{i}+\dfrac{h}{2},y_{i}+\dfrac{hk_{1}}{2})$
$k_{3} = f(x_{i}+\dfrac{h}{2},y_{i}+\dfrac{hk_{1}}{2})$
$k_{4} = f(x_{i}+h,y_{i}+hk_{3})$

and

$y_{i+1} = y_{i} + h\phi$

where

$\phi = \dfrac{1}{6}(k_{1}+2k_{2}+2k_{3}+k_{4})$

### Using RK4 for 2nd order ODE

Any second order differential equation can be written as two coupled first order equations,

$\dfrac{dy}{dt} = z$   and   $\dfrac{dz}{dt} = F(t, y, z)$

with initial conditions $y(0) = y_{0}$   and   $z(0) = z_{0}$ and step size $h$;

These coupled equations can be solved numerically using a fourth order Runge-Kutta method as follows:

$k_{y1} = z(0) = z_{0}$     and     $k_{z1} = F(t, y_{0}, z_{0})$

$k_{y2} = z_{0} + \dfrac{hk_{z1}}{2}$     and     $k_{z2} = F(t+\dfrac{h}{2}, y_{0}+\dfrac{hk_{y1}}{2}, z_{0}+\dfrac{hk_{z1}}{2})$

$k_{y3} = z_{0} + \dfrac{hk_{z2}}{2}$     and     $k_{z3} = F(t+\dfrac{h}{2}, y_{0}+\dfrac{hk_{y2}}{2}, z_{0}+\dfrac{hk_{z2}}{2})$

$k_{y4} = z_{0} + hk_{z3}$     and     $k_{z4} = hF(t+h, y_{0}+hk_{y3}, z_{0}+hk_{z3})$

$dy = \dfrac{h}{6}(k_{y1}+2k_{y2}+2k_{y3}+k_{y4})$   and   $dz = \dfrac{h}{6}(k_{z1}+2k_{z2}+2k_{z3}+k_{z4})$

So now the RK4 approximation of next values are:

$y_{1} = y_{0}+dy$    and    $z_{1} = z_{0}+dz$

and so on.

#### Pseudo code

Declare function F(t,y,z)
Declare Tinitial
Declare Tfinal
Declare N = Number of steps
Declare h = (Tinitial - Tfinal)/N
Declare Yinitial
Declare Zinitial
Y = Yinitial
Z = Zinitial
Print Yinitial, Zinitial
LOOP for T from Tinitial to Tfinal with increment of h
ky1 = Z;
kz1 = F(t, Y, Z);

ky2 = (Z + h*kz1/2);
kz2 = F(t + h/2, Y + h*ky1/2, Z + h*kz1/2);

ky3 = (Z + h*kz2/2);
kz3 = F(t + h/2, Y + h*ky2/2, Z + h*kz2/2);

ky4 = (Z + h*kz3);
kz4 = F(t + h, Y + h*ky3, Z + h*kz3);

dy = h/6*(ky1 + ky2 + ky3 + ky4);
dz = h/6*(kz1 + kz2 + kz3 + kz4);

Y = Y + dy;
Z = Z + dz;

Print Y, Z
END LOOP

Define function F(t,y,z) {
equation for F(t,y,z) goes here
}

Show Comments