Open access peer-reviewed chapter

# A Multi-Domain Spectral Collocation Approach for Solving Lane-Emden Type Equations

Written By

Motsa Sandile Sydney, Magagula Vusi Mpendulo, Goqo Sicelo Praisegod, Oyelakin Ibukun Sarah and Sibanda Precious

Submitted: October 22nd, 2015 Reviewed: March 10th, 2016 Published: August 24th, 2016

DOI: 10.5772/63016

From the Edited Volume

## Numerical Simulation

Edited by Ricardo Lopez-Ruiz

Chapter metrics overview

View Full Metrics

## Abstract

In this work, we explore the application of a novel multi-domain spectral collocation method for solving general non-linear singular initial value differential equations of the Lane-Emden type. The proposed solution approach is a simple iterative approach that does not employ linearisation of the differential equations. Spectral collocation is used to discretise the iterative scheme to form matrix equations that are solved over a sequence of non-overlapping sub-intervals of the domain. Continuity conditions are used to advance the solution across the non-overlapping sub-intervals. Different Lane-Emden equations that have been reported in the literature have been used for numerical experimentation. The results indicate that the method is very effective in solving Lane-Emden type equations. Computational error analysis is presented to demonstrate the fast convergence and high accuracy of the method of solution.

### Keywords

• Multi-domain
• Lame-Emden
• Spectral Relaxation method
• Collocation
• Astrophysics

## 1. Introduction

In the most general form, Lane-Emden type equations are given as

y+γxy+f(x,y)=g(x),x[0,),y(0)=α0,y(0)=β0E1

where the prime denotes differentiation with respect to x, f(x, y) is a non-linear function, g(x) is a prescribed function and γ, α0, β0 are known constants. In recent years, problems described by this class of differential equations have been widely investigated by many researchers because of their applications in astronomy, mathematical biology, mathematical physics, non-Newtonian fluid mechanics, and other areas of science and engineering. From a solution method viewpoint, it has been observed that, owing to the singularity at x = 0, Lane-Emden type equations are not trivial to solve. For this reason, the equations are normally used as benchmark equations for testing the effectiveness and robustness of new analytical and numerical methods of solution.

Analytical approaches that have recently been used in solving the Lane-Emden equations are mostly based on truncated series expansions. Examples include the Adomian decomposition method [13], differential transformation method [4, 5], Laplace transform [6, 7], homotopy analysis method [810], power series expansions [1114] and variational iteration method [1517]. Being power series based, the above methods have a small region of convergence and are not suitable for generating solutions in very large values of x. For this reason, most analytical approaches have only reported solutions of Lane-Emden type equations on small interval [0,1] on x axis. Despite this limitation, analytical approaches have been found to be desirable because they easily overcome the difficulty caused by the singularity at x = 0.

To overcome the limitations of analytical solution methods, several numerical approaches have been proposed for the solution of Lane-Emden type equations. Numerical methods based on spectral collocation have been found to be particularly effective. Collocation methods that have been reported recently for the solution of Lane-Emden type equations include the Bessel collocation method [18, 19], Jacobi-Gauss collocation method [20], Legendre Tau method [21], Sinc-collocation method [22], Chebyshev spectral methods [23], quasi-linearisation based Chebyshev pseudo-spectral method [24, 25] and a collocation method based on radial basis functions [26]. The discretisation scheme of collocation based methods is only implemented on interior nodes of the discretised domain. This property makes it possible for these collocation methods to overcome the difficulty of dealing with the singular point.

In this work, we present a multi-domain spectral collocation method for solving Lane-Emden equations. The method is based on the innovative idea of reducing the governing non-linear differential equations to a system of first-order equations which are solved iteratively using a Gauss-Seidel–like relaxation approach. The domain of the problem is divided into smaller non-overlapping sub-intervals on which the Chebyshev spectral collocation method is used to solve the iteration scheme. The continuity condition is used to advance the solution across neighbouring sub-intervals. The advantage of the approach is that it does not use Taylor-series based linearisation methods to simplify the non-linear differential equations. The method is free of errors associated with series truncation. The algorithm is also very easy to develop and yields very accurate results using only a few discretisation nodes. The accuracy of the method is validated against known results from the literature. The aim of the study is to explore the applicability of the multi-domain spectral collocation method to Lane-Emden type equations over semi-infinite domains. The results confirm that the method is suitable for solving all types of Lane-Emden equations.

## 2. A multi-domain pseudospectral relaxation method

In this section, we describe the development of the multi-domain collocation algorithm for the solution of Lane-Emden type equations. This algorithm addresses some limitations of standard collocation methods. Without loss of generality, we express the second-order generalised Lane-Emden Eq. (1) as a system of first-order ordinary differential equations (ODEs). If we let hy0, the initial value problem (1) transforms to

y(x)=h(x),y(0)=α0,E2
h(x)+γxh+f(x,y)=g(x),h(0)=β0,E3

where y′(x) = h(x). The following iterative scheme for solving (2) and (3) is introduced:

yi+1(x)=hi(x),yi+1(0)=α0,E4
hi+1(x)+γxhi+1(x)+f(x,yi+1)=g(x),hi+1(0)=β0.E5

Assuming that an initial approximation h0(x) is given, Eq. (4) can be solved for yi+1(x) and the solution can be used immediately in Eq. (5) which is, in turn, solved for hi+1. Thus, at each iteration level i + 1, Eqs. (4) and (5) form a pair of linear decoupled first-order differential equations. The solution procedure is discussed below.

Below, we describe the development of the multi-domain approach for solving the system of first-order Eqs. (4) and (5). The multi-domain technique approach assumes that the main interval can be decomposed into p non-overlapping sub-intervals. Let xΛ, where Λ=[a,b] is the interval where the solution of Eq. (1) exists. The sub-intervals are defined as

Λk=(xk1,xk),k=1,2,,p,with,a=x0<x1<x2<<xp=b.E6

The solution procedure assumes that the solution at each sub-interval Λk, denoted by y(k)(x), can be approximated by a Lagrange interpolation polynomial of the form

y(k)(x)s=0Ny(k)(xs)Ls(x),E7

which interpolates y(k)(x) at selected points, chosen to be the Chebyshev-Gauss-Lobatto, defined by

{τs}s=0N={cos(πsN)}s=0N.E8

The choice of the interpolation function (7) with Lagrangre polynomial basis and grid points (8) enables the use of simple formulas for converting the continuous derivatives to discrete matrix-vector form at the collocation points τs. The function Ls(x) is the characteristic Lagrange cardinal polynomial defined as

Ls(x)=s=0skMxxkxsxk,E9

that obey the Kronecker delta equation:

Ls(xk)=δsk={0 if sk1 if s=k..E10

The solution method seeks to solve each of the system of first-order ODEs independently in each kth sub-interval Λk using the solutions obtained in the preceding interval Λk1 as initial conditions. In the first interval the solution is computed on [x0,x1] and is labelled y(1)(x). The value of the solution at the last node y(1)(x1) is used as an initial condition when seeking a solution in the second sub-interval Λ2. This procedure is repeated in each Λk interval with the solutions marched across each interval using the following continuity condition:

y(k)(xk1)=y(k1)(xk1),h(k)(xk1)=h(k1)(xk1).E11

Therefore, in each sub-interval Λk, the following system of ODEs is solved:

dyi+1(k)(x)dx=hi(k)(x),yi+1(k)(xk1)=αk1E12
dhi+1(k)(x)dx+γxhi+1(k)(x)=g(x)f(x,yi+1(k)),hi+1(k)(xk1)=βk1,E13

where

αk1=y(k1)(xk1),βk1=h(k1)(xk1),k=1,2,,p.

Eqs. (12) and (13) cannot be solved exactly because of the non-linear function f(x, y). Accordingly, we employ the spectral collocation at the grid points xj for j=0,1,2,,N for each sub-interval Λk. Before the spectral method is applied, each sub-interval [xk1,xk] is transformed to [1,1] using the linear transformation:

x=xkxk12τ+xk+xk12,τ[1,1].E14

The derivatives at the collocation points are evaluated as

dy(k)dx|x=xj=s=0Ny(k)(xs)dLs(xj)dx=2Δxks=0Ny(k)(τs)dLs(τj)dτ=s=0NDjsy(k)(τs),E15

where Δxk=xkxk1, Djs=2ΔxkDjs with Djs=dLs(τj)dτ being the jth and sth entry of the standard first derivative Chebyshev differentiation matrix of size (N+1)×(N+1) as defined in [27]. Evaluating Eqs. (12) and (13) at the grid points τj for j=0,1,2,,N gives

s=0NDjsyi+1(k)(τs)=hi(k)(τj),yi+1(k)(τN)=αk1,E16
s=0NDjshi+1(k)(τs)+γxjhi+1(k)(τj)=g(xj)f(xs,yi+1(k)(τj)),hi+1(k)(τN)=βk1.E17

Including the relevant initial conditions in Eqs. (16) and (17) yields

s=0N1Djsyi+1(k)(τs)=Vi(k)(τj),E18
s=0N1Djshi+1(k)(τs)+γxjhi+1(k)(τj)=Wi(k)(τj).E19

where

Vi(k)(τj)=hi(k)(τj)DjNyi+1(k)(τN),E20
Wi(k)(τj)=g(xj)f(xj,yi+1(k)(τj))DjNhi+1(k)(τN),E21

Eqs. (18) and (19) can be expressed in matrix form as follows:

DYi+1(k)=Vi(k),E22
(D+ζ)Hi+1(k)=Wi(k),E23

where the vectors Y(k),H(k),V(k) and W(k) are respectively given by

Y(k)=[y(k)(τ0),y(k)(τ1),y(k)(τ2),,y(k)(τN1)]T,E24
H(k)=[h(k)(τ0),h(k)(τ1),h(k)(τ2),,h(k)(τN1)]T,E25
V(k)=[v(k)(τ0),v(k)(τ1),v(k)(τ2),,v(k)(τN1)]T,E26
W(k)=[w(k)(τ0),w(k)(τ1),w(k)(τ2),,w(k)(τN1)]T,E27

and T denotes the transpose of the vector. The N × N diagonal matrix ζ is given by

ζ=[γx0γx1γxN1].E28

The approximate values of y(k)(x) and dy(k)(x)dx are obtained by iteratively solving Eqs. (22) and (23), respectively, for i=0,1,2,, starting from suitable initial approximation.

## 3. Numerical experiments

In this section, we discuss the numerical experiments to be used to demonstrate the accuracy and general performance of the multi-domain pseudospectral relaxation method. In these numerical experiments, we have selected equations with known exact solutions, and to determine the accuracy of the method, we find the relative error. The relative error is defined as

Ej=|ye(xj)ya(xj)||ye(xj)|E29

where Ej is the relative error at a grid point xj, ye(xj) and ya(xj) are the exact and approximate solutions at a grid point xj, respectively.

### Example 1

We first consider the linear, homogeneous Lane-Emden equation, with variable coefficients:

y+2xy2(2x2+3)y=0,x>0,subjecttoy(0)=1andy(0)=0,E30

which has the exact solution:

y(x)=ex2.

Eq. (30) has been solved by various researchers using different techniques such as the variational iteration method and the homotopy-pertubation method [16, 28, 29].

### Example 2

Secondly, the non-linear, homogeneous Lane-Emden equation, with variable coefficients is considered:

y+2xy+4(2ey+ey/2)=0,x>0,subjecttoy(0)=0andy(0)=0E31

which has the exact solution:

y(x)=2ln(1+x2).

Eq. (31) has been solved by [29, 30, 31] using the lie symmetry, the homotopy-pertubation method and the variational iteration method.

### Example 3

In this example, we consider the non-linear, variable coefficients, homogeneous Lane-Emden equation:

y+2xy6y=4yln(y),x>0,subjecttoy(0)=1andy(0)=0E32

which has the exact solution:

y(x)=ex2.

Ramos [28] solved Eq. (32) using the variation iteration method, while Yildirim [16] solved the same equation using the linearisation method.

### Example 4

We consider the non-linear, homogeneous Lane-Emden equation:

y+6xy+2y(7+ln(y2))=0,x>0,subjecttoy(0)=1andy(0)=0E33

which has the exact solution [32]:

y(x)=ex2.

This example was solved by Wazwaz [32] using the Adomian decomposition method.

### Example 5

We consider the non-linear, homogenous Lane-Emden equation, which represent an infinite circular cylinder in astrophysics:

y+2xy+y36x6=0,x>0,subjecttoy(0)=0andy(0)=0,E34

which has the exact solution [33]:

y(x)=x2.

This example was solved by Yin et al. [33] using the modified Laplace decomposition method.

### Example 6

Lastly, we consider the Lane-Emden equation with the general form:

y+2xy+ym=0,x>0,subjecttoy(0)=1andy(0)=0.E35

Eq. (35) is the standard Lane-Emden equation that models the thermal behaviour of a spherical cloud of glass acting under the mutual attraction of its molecules and subject to classical laws of thermodynamics [34]. The values of m in the interval [0,5] are most physically interesting to study. The equation is linear when m = 0 and m = 1 and non-linear for values of m>1. Wazwaz [35] gave the general solution of Eq. (35) in series form as

y(x)=116x2+m120x4m(8m5)3.7!x6+m(70183m+122m2)9.9!x8+m(31501080m+12642m25032m3)45.11!x10+.E36

Analytical solutions for m=0,1 and m=5 are given as [35]

y(x)=113!x2,y(x)=sinxx,andy(x)=(1+x23)12,E37

respectively. In this example, we consider the Lane-Emden Eq. (35) for m = 5. We therefore consider the equation:

y+2xy+y5=0,x>0,subjecttoy(0)=1andy(0)=0E38

with the exact solution:

y(x)=(1+x23)1/2.E39

## 4. Results and discussion

In this section, we discuss and present the results obtained using the proposed algorithm. We used the six examples in the previous section. The results were generated using MATLAB 2013. To validate the accuracy, computational time and general performance of the method, we computed relative errors and computational time of each of the numerical examples. The level of accuracy of the algorithm at a particular level is determined by the relative error Rk defined by

Rk=|ye(xk)ya(xk)||ye(xk)|,0kN,E40

where N is the number of grid points, ya(xk) is the approximate solution and ye(xk) is the exact solution at the grid point xk. The graphs were all generated using N = 4. For each numerical experiment, we present the relative error and the corresponding approximate solution for N = 4 and N = 6 in a tabular form. The central processing unit computational time is displayed. Graphs showing an excellent agreement between the analytical and approximate solutions are presented for each numerical experiment. These graphs validate the accuracy of the method. Error graphs showing the distribution of the relative errors are also presented. These error graphs are in excellent agreement with the results presented in the tables for all the numerical experiments used in this chapter.

N = 4 N = 6
x Exact Approximate Relative error Approximate Relative error
0.2 1.060313 1.060313 1.884728e−14 1.060313 4.251109e−014
0.4 1.140018 1.140018 3.505912e−14 1.140018 6.505415e−014
0.6 1.371416 1.371416 7.852585e−14 1.371416 1.066980e−013
0.8 1.787189 1.787189 1.413886e−13 1.787189 1.474757e−013
1.0 2.522988 2.522988 2.355112e−13 2.522988 1.864022e−013
1.2 3.858367 3.858367 3.776361e−13 3.858367 2.234047e−013
1.4 6.391979 6.391979 5.987445e−13 6.391979 2.591455e−013
1.6 11.471251 11.471251 9.227679e−13 11.471251 2.937560e−013
1.8 22.301278 22.301278 1.396995e−12 22.301278 3.316738e−013
2.0 46.966942 46.966942 2.066546e−12 46.966942 3.633883e−013
CPU Time (sec) 0.659414 0.629082 0.659414

### Table 1.

Analytical, approximate solutions and relative errors for Example 1.

Table 1 shows the exact, approximate solution and the relative error for Eq. (30). For N = 4, the multi-domain spectral relaxation method gives a relative error of approximately 10−12 and for N = 6, the relative error is on average 10−13. Increasing the number of grid points results in a more accurate solution. The results obtained from the multi-domain spectral relaxation method are remarkable since few grid points give accurate results in a large domain. Using a few grid points ensures that the numerical method converges within few seconds.

Figure 1 shows the relative error displayed in Table 1 for N = 4. The results in Figure 1 are in excellent agreement with those in Table 1. Figure 2 shows the analytical and approximate solutions. Since the approximate solution is superimposed on the exact solutions, this implies that the multi-domain pseudospectral relaxation method converged to the exact solution over the domain x[0,2]. Table 1 and Figures 1 and 2 validate the accuracy and computational efficiency of the multi-domain pseudospectral relaxation method for Eq. (30).

The results obtained from approximating the solution to Eq. (31) using the multi-domain pseudospectral relaxation method are shown in Table 1 and Figures 3 and 4. In Table 2, we display the exact solution, approximate solution and the relative error of Eq. (31) in Example 2. For N = 4, the multi-domain spectral relaxation method gives a relative error of approximately 10−11. For N = 4 the relative error is approximately 10−13. We observe that increasing the number of grid points decreases the relative error. The multi-domain pseudospectral relaxation method uses a few grid points to achieve accurate results in the domain x[0,20]. A maximum of N = 6 grid points ensured that the numerical method converged to an error of 10−13 within a fraction of a second.

Figure 3 shows the relative error displayed in Table 2 for N = 4. The results in Figure 3 are in excellent agreement with those in Table 2. Figure 4 shows the analytical and approximate solutions of Eq. (31). The approximate solution superimposed on the exact solutions shows that the multi-domain pseudospectral relaxation method converged to the exact solution over the domain x[0,20]. The match between the exact and approximate solutions in Table 2 and Figures 3 and 4 validates the accuracy and computational efficiency of the multi-domain pseudospectral relaxation method for Eq. (31).

N = 4 N = 6
x Exact Approximate Relative error Approximate Relative error
2 –3.271943 –3.271943 8.810376e−011 –3.271943 2.704972e−014
4 –5.697684 –5.697682 1.000287e−011 –5.697684 6.202217e−014
6 –7.243401 –7.243399 8.316853e−012 –7.243401 2.136148e−014
8 –8.365152 –8.365151 1.308867e−011 –8.365152 3.825343e−015
10 –9.243421 –9.243420 1.383787e−011 –9.243421 4.807193e−015
12 –9.964487 –9.964487 1.318604e−011 –9.964487 5.172115e−015
14 –10.575872 –10.575872 1.205994e−011 –10.575872 6.720993e−015
16 –11.106445 –11.106445 1.080384e−011 –11.106445 3.199792e−016
18 –11.575028 –11.575028 9.561315e−012 –11.575028 6.293749e−015
20 –11.927621 –11.927622 8.504839e−012 –11.927621 1.738910e−014
CPU Time (sec) 0.649135 0.602359 0.649135

### Table 2.

Analytical, approximate solutions and relative errors for Example 2.

The results obtained from approximating the solution to Eq. (32) are given in Table 3 and Figures 5 and 6. Table 3 shows the exact solution, the approximate solution and the relative error of Eq. (32). For N = 4, the multi-domain spectral relaxation method gives a relative error of approximately 10−12, while for N = 6, the relative error is approximately 10−13. We observe that increasing the number of grid points decreases the relative error and hence increases the accuracy of the method. This pseudospectral method uses a few grid points to achieve accurate results in the domain x[0,2]. N = 6 grid points ensured that the numerical method converged to an error of 10−13 within few seconds.

N = 4 N = 6
x Exact Approximate Relative error Approximate Relative error
0.2 1.029322 1.029322 2.502345e−014 1.029322 1.682611e−014
0.4 1.146713 1.146713 1.316722e−013 1.146713 4.124439e−014
0.6 1.383892 1.383892 3.732052e−013 1.383892 6.915367e−014
0.8 1.809228 1.809228 8.178657e−013 1.809228 1.078787e−013
1.0 2.562286 2.562286 1.581696e−012 2.562286 1.493997e−013
1.2 3.931024 3.931024 2.814656e−012 3.931024 1.883216e−013
1.4 6.53322 6.53322 4.761036e−012 6.53322 2.300241e−013
1.6 11.762306 11.762306 7.722919e−012 11.762306 2.750095e−013
1.8 22.940410 22.940410 1.211045e−011 22.940410 3.176323e−013
2.0 48.467816 48.467816 1.845750e−011 48.467816 3.667955e−013
CPU Time (sec) 1.069773 1.010628 1.069773

### Table 3.

Analytical, approximate solutions and relative errors for Example 3.

The relative error shown in Table 3 is displayed in Figure 5. The results in Figure 5 are in excellent agreement with those in Table 3. Figure 6 shows the analytical and approximate solutions of Eq. (32). The approximate solution being superimposed on the exact solutions implies that the multi-domain pseudospectral relaxation method converged to the exact solution over the domain x[0,2]. The match between the exact and approximate solutions in Table 3 and Figures 5 and 6 validates the accuracy and computational efficiency of the multi-domain pseudospectral relaxation method for Eq. (32).

The results obtained from approximating the solution to Eq. (33) are given in Table 4 and Figures 7 and 8. Table 4 shows the exact solution, the approximate solution and the relative error of Eq. (33). For N = 4, the multi-domain spectral relaxation method gives a relative error of approximately 10−13. For N = 6, the relative error is also approximately 10−13. Increasing the number of grid points decreases the relative error. Thus, a maximum of N = 6 grid points ensures convergence of the method. N = 6 grid points ensured that the numerical method converged to an error of 10−13 in a fraction of seconds.

N = 4 N = 6
x Exact Approximate Relative error Approximate Relative error
0.2 0.974097 0.974097 7.978218e−015 0.974097 2.803774e−014
0.4 0.877179 0.877179 1.215047e−014 0.877179 6.999178e−014
0.6 0.729173 0.729173 1.248514e−014 0.729173 1.079508e−013
0.8 0.559538 0.559538 1.805602e−014 0.559538 1.341305e−013
1.0 0.396355 0.396355 3.221242e−014 0.396355 1.606419e−013
1.2 0.259177 0.259177 6.232707e−014 0.259177 1.756296e−013
1.4 0.156446 0.156446 1.011254e−013 0.156446 1.754615e−013
1.6 0.087174 0.087174 1.525094e−013 0.087174 1.404105e−013
1.8 0.044840 0.044840 1.897188e−013 0.044840 3.528213e−014
2.0 0.021292 0.021292 1.404623e−013 0.021292 2.123230e−013
CPU Time (sec) 1.080828 1.046274 1.080828

### Table 4.

Analytical, approximate solutions and relative errors for Example 4.

Figure 7 shows the relative error displayed in Table 4 for N = 4. The results in Figure 7 are in excellent agreement with those in Table 4. Figure 8 shows the analytical and approximate solutions. Since the approximate solution is superimposed on the exact solutions, this implies that the multi-domain pseudospectral relaxation method converged to the exact solution over the domain x[0,2]. Table 4 and Figures 7 and 8 validate the accuracy and computational efficiency of the multi-domain pseudospectral relaxation method for Eq. (33).

The results obtained from approximating the solution to Eq. (34) are given in Table 5 and Figures 9 and 10. Table 5 shows the exact solution, the approximate solution and the relative error of Eq. (34). The multi-domain pseudospectral relaxation method gives a relative error of approximately 10−14. For N = 6, the relative error is approximately 10−14. Increasing the number of grid points decreases the relative error and thus implying that the numerical method converged to the exact solution. The pseudospectral collocation method uses a few grid points to achieve accurate results in the domain x[0,2]. The numerical method converged to an error of 10−14 in a fraction of a second.

N = 4 N = 6
x Exact Approximate Relative error Approximate Relative error
0.2 0.026244 0.026244 4.759186e−015 0.026244 1.692155e−014
0.4 0.131044 0.131044 1.397903e−014 0.131044 3.473577e−014
0.6 0.315844 0.315844 2.495721e−014 0.315844 5.518706e−014
0.8 0.580644 0.580644 3.537301e−014 0.580644 6.864276e−014
1.0 0.925444 0.925444 4.138845e−014 0.925444 8.445643e−014
1.2 1.350244 1.350244 3.880967e−014 1.350244 9.981979e−014
1.4 1.855044 1.855044 4.093663e−014 1.855044 1.098825e−013
1.6 2.439844 2.439844 3.749517e−014 2.439844 8.700337e−014
1.8 3.104644 3.104644 2.445989e−014 3.104644 4.706026e−014
2.0 3.849444 3.849444 4.614580e−016 3.849444 6.575777e−015
CPU Time (sec) 0.705910 0.638383 0.705910

### Table 5.

Analytical, approximate solutions and relative errors for Example 5.

Figure 9 shows the relative error displayed in Table 5 for N = 4. The results in Figure 9 are in excellent agreement with those in Table 5. Figure 10 shows the analytical and approximate solutions. Since the approximate solution is superimposed on the exact solutions, this implies that the multi-domain pseudospectral relaxation method converged to the exact solution over the domain x[0,2]. Table 5 and Figures 9 and 10 validate the accuracy and computational efficiency of the multi-domain pseudospectral relaxation method for Eq. (34).

The results obtained from approximating the solution to Eq. (38) are given in Table 6 and Figures 11, 12 and 13. Table 6 shows values obtained for the exact and approximate solution together with the respective relative error values of Eq. (38). For N = 4 the multi-domain pseudospectral collocation method gives a relative error of approximately 1011 and 10−13 for N = 6. An increase in the number of grid points results in a decrease in the relative error. This implies that the numerical method converges to the exact solution of Eq. (38). The multi-domain pseudospectral relaxation method used a few grid points to achieve accurate results in the domain x[0,20]. N = 6 grid points ensured that the numerical method converged to an error of 10−13 within a few seconds as shown in Table 6.

N = 4 N = 6
x Exact Approximate Relative error Approximate Relative error
2 0.650926 0.650926 1.991857e−011 0.650926 3.189482e−014
4 0.395693 0.395693 2.135316e−011 0.395693 3.745701e−014
6 0.276499 0.276499 6.204217e−012 0.276499 5.842239e−014
8 0.211100 0.211100 8.877991e−012 0.211100 1.519920e−013
10 0.170333 0.170333 2.323778e−011 0.170333 2.328536e−013
12 0.142624 0.142624 3.707619e−011 0.142624 3.096185e−013
14 0.122609 0.122609 5.055364e−011 0.122609 4.021536e−013
16 0.107492 0.107492 6.378059e−011 0.107492 4.983484e−013
18 0.095677 0.095607 7.683911e−011 0.095677 5.994828e−013
20 0.087057 0.087057 8.847457e−011 0.087057 6.926360e−013
CPU Time (sec) 1.199892 1.046417 1.199892

### Table 6.

Analytical, approximate solutions and relative errors for Example 6.

Figure 11 shows the plot for different values of m. The results are in good agreement with those obtained by [25, 26]. Figure 12 displays the relative error graph of Eq. (38). The results in Figure 12 are in excellent agreement with those obtained in Table 6. The comparison between the exact solution and approximate solution of Eq. (38) is shown in Figure 13. The superimposition of the approximate solution on the exact solution implies that the multi-domain pseudospectral relaxation method converged to the exact solution over the domain x[0,20]. Table 6 and Figures 11, 12 and 13 give a validation of the accuracy and computational efficiency of the multi-domain pseudospectral relaxation method for Eq. (38).

## 5. Conclusion

In this work, we presented a multi-domain spectral collocation method for solving Lane-Emden equations. This numerical method was used to solve six Lane-Emden equations. The results obtained were remarkable in the sense that using few grid points, we were able to achieve accurate results. We were able to compute the results using minimal computational time. The approximate solutions were in excellent agreement with the exact solutions in all the numerical experiments. We presented error graphs, approximate and exact solution graphs to show accuracy of the method. Tables showing relative errors were also generated to show accuracy and computational efficiency of the numerical method presented.

This approach is useful for solving other nonlinear, singular initial value problems. This approach is an alternative to the already existing list of numerical methods that can be used to solve such equations. This numerical approach can be extended to solve time-dependent Lane-Emden equations and other singular value type of equations.

## Acknowledgments

This work is based on the research supported in part by the National Research Foundation of South Africa (Grant No: 85596).

## References

1. 1. S.G. Hosseini, S. Abbasbandy, Solution of Lane-Emden type equations by combination of the spectral method and adomian decomposition method, Mathematical Problems in Engineering, vol. 2015, 10 p., 2015, Article ID 534754. doi:10.1155/2015/534754
2. 2. A.M. Wazwaz, A new method for solving singular initial value problems in the second-order ordinary differential equations, Applied Mathematics and Computation, vol. 128, no. 1, pp. 45–57, 2002.
3. 3. A.M. Wazwaz, Adomian decomposition method for a reliable treatment of the Emden-Fowler equation, Applied Mathematics and Computation, vol. 161, no. 2, pp. 543–560, 2005.
4. 4. V.S. Ertürk, Differential transformation method for solving differential equations of Lane-Emden type, Mathematical & Computational Applications, vol. 12, no. 3, pp. 135–139, 2007.
5. 5. Y. Khan, Z. Svoboda, Z. Šmarda, Solving certain classes of Lane-Emden type equations using the differential transformation method, Advances in Difference Equations, vol. 2012, p. 174, 2012. doi:10.1186/1687-1847-2012-174
6. 6. F. Yin, J. Song, F. Lu, H. Leng, A coupled method of Laplace transform and Legendre wavelets for Lane-Emden-type differential equations, Journal of Applied Mathematics, vol. 2012, 16 p., 2012, Article ID 163821. doi:10.1155/2012/163821
7. 7. F. Yin, W. Han, J. Song, Modified Laplace decomposition method for Lane-Emden type differential equations, International Journal of Applied Physics and Mathematics, vol. 3, no. 2, 2013. doi:10.7763/IJAPM.2013.V3.184
8. 8. A.S. Bataineh, M.S.M. Noorani, I. Hashim, Homotopy analysis method for singular IVPs of Emden-Fowler type, Communications in Nonlinear Science and Numerical Simulation, vol. 14, no. 4, pp. 1121–1131, 2009.
9. 9. S. Liao, A new analytic algorithm of Lane-Emden type equations, Applied Mathematics and Computation, vol. 142, no. 1, pp. 1–16, 2003.
10. 10. O.P. Singh, R.K. Pandey, V.K. Singh, An analytic algorithm of Lane-Emden type equations arising in astrophysics using modified Homotopy analysis method, Computer Physics Communications, vol. 180, no. 7, pp. 1116–1124, 2009.
11. 11. A. Aslanov, Determination of convergence intervals of the series solutions of Emden-Fowler equations using polytropes and isothermal spheres, Physics Letters A, vol. 372, no. 20, pp. 3555–3561, 2008.
12. 12. C. Hunter, Series solutions for polytropes and the isothermal sphere, Monthly Notices of the Royal Astronomical Society, vol. 328, no. 3, pp. 839–847, 2001.
13. 13. C. Mohan, A.R. Al-Bayaty, Power-series solutions of the Lane-Emden equation, Astrophysics and Space Science, vol. 73, no. 1, pp. 227–239, 1980.
14. 14. I.W. Roxburghm, L.M. Stockman, Power series solutions of the polytrope equations, Monthly Notices of the Royal Astronomical Society, vol. 303, no. 3, pp. 466–470, 1999.
15. 15. J.H. He, Variational approach to the Lane-Emden equation, Applied Mathematics and Computation, vol. 143, no. 2–3, pp. 539–541, 2003.
16. 16. A. Yildirim, T. Öziş, Solutions of singular IVPs of Lane-Emden type by the variational iteration method, Nonlinear Analysis, vol. 70, no. 6, pp. 2480–2484, 2009.
17. 17. A.M. Wazwaz, The variational iteration method for solving the Volterra integro-differential forms of the Lane-Emden equations of the first and the second kind, Journal of Mathematical Chemistry, vol. 52, no. 2, pp. 613–626, 2014.
18. 18. S. Yüzbasi, A numerical approach for solving the high-order linear singular differential-difference equations, Computers and Mathematics with Applications, vol. 62, pp. 2289–2303, 2011.
19. 19. S. Yüzbasi, M. Sezer, An improved Bessel collocation method with a residual error function to solve a class of Lane-Emden differential equations, Mathematical and Computer Modelling, vol. 57, pp. 1298–1311, 2013.
20. 20. A.H. Bharwy, A.S. Alofi, A Jacobi-Gauss collocation method for solving nonlinear Lane-Emden type equations, Communications in Nonlinear Science and Numerical Simulation, vol. 17, pp. 62–70, 2012.
21. 21. K. Parand, M. Razzaghi, Rational legendre approximation for solving some physical problems on semi-infinite intervals, Physica Scripta, vol. 69, no. 5, pp. 353–357, 2004.
22. 22. K. Parand, A. Pirkhedri, Sinc-collocation method for solving astrophysics equations, New Astronomy, vol. 15, no. 6, pp. 533–537, 2010.
23. 23. J.P. Boyd, Chebyshev spectral methods and the Lane-Emden problem, Numerical Mathematics: Theory, Methods and Applications, vol. 4, no. 2, pp. 142–157, 2011.
24. 24. S.S. Motsa, P. Sibanda, A new algorithm for solving singular IVPs of Lane-Emden type, in Proceedings of the 4th International Conference on Applied Mathematics, Simulation, Modelling (ASM’10), pp. 176–180, Corfu Island, Greece, July 2010.
25. 25. S.S. Motsa, S. Shateyi, A successive linearization method approach to solve Lane-Emden type of equations, Mathematical Problems in Engineering, vol. 2012, 14 p., 2012, Article ID 280702. doi:10.1155/2012/280702
26. 26. K. Parand, S. Abbasbandy, S. Kazem, A.R. Rezaei, An improved numerical method for a class of astrophysics problems based on radial basis functions, Physica Scripta, vol. 83, no. 1, 2011, Article ID 015011.
27. 27. L.N. Trefethen, Spectral methods in MATLAB, SIAM, Philadelphia, 2000.
28. 28. J.I. Ramos, Linearization techniques for singular initial-value problems of ordinary differential equations, Journal of Applied Mathematics and Computation, vol. 161, pp. 525–542, 2005.
29. 29. M.S.H. Chowdhury, I. Hashim, Solution of a class of singular second-order IVPs by homotopy-perturbation method, Physics Letters A, vol. 365, pp. 439–447, 2007.
30. 30. C.M. Khalique, P. Ntsime, Exact solutions of the Lane-Emden-type equation, New Astronomy, vol. 13, no. 7, pp. 476–480, 2008.
31. 31. K. Parand, M. Dehghan, A.R. Rezaei, S.M. Ghaderi, An approximation algorthim for the solution of the nonlinear Lane-Emden type equation arising in Astorphisics using Hermit functions Collocation method, Computer Physics Communications, vol. 181, pp. 1096–1108, 2010.
32. 32. A.M. Wazwaz, A new method for solving singular initial value problems in second-order ordinary differential equations, Applied Mathematics and Computation, vol. 128, pp. 45–57, 2002.
33. 33. F.K. Yin, W.Y. Han, J.Q. Song, Modified Laplace decomposition method for Lane-Emden type differential equations, International Journal of Applied Physics and Mathematics, vol. 3, pp. 98–102 no. 2, March 2013.
34. 34. N. T. Shawagfeh, Nonperturbative approximate solution of Lane-Emden equation, Journal of Mathematical Physics, vol. 34, pp. 4364–4369, 1993.
35. 35. A.M. Wazwaz, A new algorithm for solving differential equations of Lane-Emden type, Applied Mathematics and Computation, vol. 118, pp. 287–310, 2001.

Written By

Motsa Sandile Sydney, Magagula Vusi Mpendulo, Goqo Sicelo Praisegod, Oyelakin Ibukun Sarah and Sibanda Precious

Submitted: October 22nd, 2015 Reviewed: March 10th, 2016 Published: August 24th, 2016