Efficient Boundary Value Problem Solution for a Lane-emden Equation

An efficient method for determining an initial guess to an iterative solution to the boundary value problem y ′′ + (k/x)y ′ + δe y = 0 solved subject to y ′ (0) = 0 and y(1) = 0 is proposed. This initial guess overcomes the instability that occurs at the boundary y ′ (0) = 0 for k ≤ 1. When the iterative method becomes singular we can conclude that the maximum value of the critical parameter δ has been determined.


INTRODUCTION
In this paper we propose an efficient method for determining an initial guess for an iterative solution to a boundary value problems for Lane-Emden type equations.A generalized Lane-Emden equation of the second-kind is given by [3] y ′′ + α x y ′ + βx ν−1 e ny = 0, where ′ = d/dx.α, β, ν and n are constants determined by the physics of the problem under investigation.Of particular interest is the Lane-Emden equation solved subject to the boundary conditions that models the steady state temperature distribution in an infinite rectangular slab (k = 0), infinite cylinder (k = 1) or a spherical vessel (k = 2) [8].Geometries corresponding to non-integer values of k have been investigated by Balakrishnan et al. [2].The constant δ is known as the Frank-Kamanetskii parameter [8].The value of δ is related to the critical temperature at which ignition in a thermal explosion takes place.Harley and Momoniat [7] have shown that for values of the shape factor k > 1 the boundary condition y ′ = 0 is stable on the line x = 0 while it is unstable for k ≤ 1. Harley and Momoniat [7] have indicated that an appropriate choice of initial guess is sufficient to overcome this instability.In this paper we show how such a suitable initial guess can be determined.The results obtained here can easily be applied to other boundary value problems.
To show where the discrepancy comes from we consider the analytical solution to (2) when k = 1.When k = 1 (2) is written as It can easily be shown that (4) admits the solutions [8,4,6] both of which satisfy y ′ (0) = 0. c 1 and c 2 are constants to be determined.Imposing the boundary condition y(1) = 0 on (5) we find that and We plot the four analytical solutions given by (5) where the constants c 1 and c 2 are given by ( 6) and ( 7) respectively in Figure 1 below.The solutions overlap and result in two solutions which satisfy the boundary conditions (3).
In Figure 2 we plot the temperature at the centre of the vessel y(0) against the critical parameter δ given by the solutions (5).We label the solutions y 1 and y 2 .We note that the solutions give a bifurcation which yields a critical value at which δ is maximum.For the case k = 1 this value is given by δ = 2.The corresponding value for y(0) when δ = 2 is given by y(0) = log 4 ≈ 1.38629.
From the results for the analytical solution it is clear that only one branch of the solution needs to be found to obtain the critical value of δ.In the next section we propose initial guesses that can be used to determine the lower branch of the solutions and hence the critical values of δ for the cases k ≤ 1 and k > 1.
The rest of the paper is divided up as follows: in section 2 we derive approximate solutions that make useful initial guesses.Concluding remarks are made in section 3.

APPROXIMATE SOLUTIONS
The first step is to determine a suitable initial guess that can be used to determine a solution to (2) solved subject to (3).We note that (2) has two parameters δ and k.To make analytical progress in determining a suitable initial guess we assume that The model equation ( 2) can then be written as We investigate an approximate solution admitted by ( 9) of the form Substituting (10) into ( 9) and separating by coefficients of k we obtain the the system of equations Solving the system (11) and imposing the boundary conditions obtained from (3) we find that We now show how (13) can be used as an initial guess in a numerical scheme to solve (2) for k ≪ 1.We approximate the singularity in (2) at x = 0 by where we have used the fact that (1/x)y ′ ≈ y ′′ at x = 0. Hence we solve the system We divide the interval [0, 1] up into n + 1 equidistant points We use the notational convention y i = y(x i ).Derivatives are approximated at x i by the central differences The boundary conditions (3) imply that The resulting nonlinear system of equations we have to solve is given by  . . where The system ( 19) is reduced to the linear system by writing it as where we assume an initial y vector with components given by The system (21) was presented in Harley and Momoniat [7] and is solved iteratively.
In Figure 3 we plot the changing values for y(0) for increasing δ obtained from the numerical solution for different values of k.The numerical iteration for the initial guess (13) becomes singular for k = 0.1 when δ = 0.98 and y(0) = 1.26032.For k = 0.5 the numerical iteration becomes singular when δ = 1.4 with the corresponding value of y(0) = 1.48926.It is at these critical values where the upper solution takes over.We can conclude that we have achieved the maximum values of y(0) when the iteration becomes singular.
For the case k ≫ 1 we divide (9) by k to obtain For large k we determine an approximate solution of the form Substituting (24) into (23) and separating by coefficients of 1/k we obtain the system Solving the system (25) subject to we find that We plot the varying values of y(0) for increasing values of δ in Figure 4.As indicated above we find that the iteration becomes singular for k = 1.5 when δ = 3.32 and y(0) = 3.74816.For k = 2.0 the iteration becomes singular when δ = 2.64 and y(0) = 3.98378.

CONCLUDING REMARKS
In this paper we have shown how an efficient guess for a Lane-Emden can be determined by matching the order of the parameters and then finding a series solution of the simplified equation.The initial guesses obtained here can easily be used in the package bvp4c in MATLAB [9] to overcome the instability in the boundary conditions indicated by Harley and Momoniat [7].We showed how these initial guesses can be used in an iterative solution to the problem to determine the variation in y(0) with increasing δ.When compared to an analytical solution for the case k = 1 we observed that the initial guesses obtained here tend to a lower solution.When the iteration in the lower solution becomes singular we find the maximum value of δ.The values at which the iteration becomes singular can be refined by using a different method of solving the nonlinear system like Newton's method.Alternatively, the values for the tolerances can be set lower increasing the computing time.

Figure 3 :
Figure3: Plot of the temperature at the centre of the vessel, y(0), against the critical parameter δ obtained from the numerical solution for k = 0.1 and k = 0.5 using (13) as an initial guess.

Figure 4 :
Figure 4: Plot of the temperature at the centre of the vessel, y(0), against the critical parameter δ obtained from the numerical solution for k = 1.5 and k = 2.0 using (27) as an initial guess.