# A new analytical technique for solving Lane - Emden type equations arising in astrophysics.

1 Introduction

Many problems of science and engineering lead to different types of differential equations and it is still very hard to solve them in the presence of strong nonlinearity. Many numerical methods have been dealt with in order to solve these equations. Alternatively, there is great interest in discovering methods for analytic approximate solutions. Recently, there has been much attention devoted to investigate better and more efficient analytical techniques such as homotopy decomposition method [1], auxiliary equation method [2], homotopy perturbation method [3, 4], Taylor collocation method [5], Sumudu transform method [6] and Adomian decomposition method [7, 8].

Emden-Fowler equation is one of the most important differential equations of mathematical physics [9, 10]. It distinctively characterizes many scientific phenomena. The generalized Emden-Fowler equation is defined in the following form

y" + [alpha](x)y' + [beta](x)[gamma](y) = 0 (1.1)

subject to conditions

y(0)= A, y (0)= B

where A, B are constants and [alpha](x), [beta](x), [gamma](y) are some arbitrary functions. For different [gamma](y), the Eq. (1.1) has been subject of many studies in the literature such as the theory of stellar structure, thermionic currents and isothermal gas spheres [11]. When [alpha](x) = [k/x], [beta]{x) = [[beta].sub.0][x.sup.r],[gamma](y) = [y.sup.s] (k and [[beta].sub.0] are constants, s and r are real numbers), Eq.(1.1) reduces to the classic Emden-Fowler equation:

y" + [k/x]y' + [[beta].sub.0][x.sup.r][y.sup.s] = 0; y(0) = A,y'(Q) = B (1.2)

Furthermore, by choosing r = 0 and k = 2, we get the standard Lane-Emden equation

y" + [2/x]y' + [[beta].sub.0][y.sup.s] = 0 (1.3)

with supplementary conditions

y(0)= A, y' (0)= 0

which arises in astrophysics. Eq. (1.2) is also used to model the thermal behavior of a spherical cloud of gas acting under the mutual attraction of its molecules [12]. Many analytical techniques have been considered by various researchers to obtain the approximate solutions for these types of equations [13-16].

In this study, we derive a new effective technique namely optimal perturbation iteration method (OPIM) to get a new approximate solutions for nonlinear problems. Main idea of this method is essentially based on optimal homotopy asymptotic method [17, 18] and perturbation iteration technique [19]. Both of these methods have been recently developed and they have been successfully implemented to some strongly nonlinear systems [20-24]. We apply OPIM to obtain more reliable approximate solutions to the Lane-Emden equations

y" + [2/x]y' + [beta](x)[gamma](y) = 0; y(0) = A, y'(0) = 0 (1.4)

with different choices of [beta](x), [gamma](y). It is also indicated that this method enables us to control the convergence of solution series for the given illustrations.

2 Optimal Perturbation Iteration Algorithm

In this section, the following formulation is given to explain the basic concept of OPIM for second order differential equations.

(a) Write the governing differential equation as:

F(y", y', y, [epsilon]) = A + g(x) = 0 (2.1)

where [epsilon] is the auxiliary perturbation parameter, y = y(x) is the unknown function and g(x) is the source term. Furthermore, (2.1) can be decomposed into A = L + N where L is the linear simpler part, which can be easily managed and N is the remaining part which is more crucial for algorithms of OPIM. Here we have a great freedom to choose linear part L.

(b) Approximate solution is taken as

[y.sub.n+1] = [y.sub.n] + [[epsilon]([y.sub.c]).sub.n] (2.2)

with one correction term in the perturbation expansion. Inserting (2.2) into the (2.1) and expanding the remaining part (N) in a Taylor series with first derivatives yields optimal perturbation iteration algorithm (OPIA):

N([y''.sub.n], [y'.sub.n], [y.sub.n], 0) + [N.sub.y] ([y.sub.c])[n.sup.[epsilon]] + [N'.sub.y], ([y'.sub.c])[n.sup.[epsilon]] + [N.sub.y"],([y".sub.c])[n.sup.[epsilon]] + [N.sup.[epsilon]][epsilon] = -L - g(x). (2.3)

It should be noted that all derivatives and functions are calculated at [epsilon] = 0. To describe the iterative scheme, first correction term [([y.sub.c]).sub.0] can be computed from the algorithm (2.3) by using a first guess [y.sub.0] and initial condition(s). (2.3) may seem complicated at first, but it should not be forgotten that we use the general form of the differential equations of second order to illustrate the proposed method. Actually, most differential equations in literature contain only some of the nonlinear terms y, y', y". So, the algorithm (2.3) reduces to some simple mathematical expressions in many cases.

c) Use the following equation

[y.sub.n+1] = [y.sub.n] + [[S.sub.n] ([y.sub.c]).sub.n] (2.4)

to increase the accuracy of the results and effectiveness of the method. Here [S.sub.0], [S.sub.1], [S.sub.2],... are convergence control parameters which provide us to adjust the convergence.

Proceeding for n = 0,1,..approximate solutions are found as:

[y.sub.1] = y(x, [S.sub.0]) = [y.sub.0] + [S.sub.0] ['([y.sub.c]).sub.0] [y.sub.2] (x, [S.sub.0], [S.sub.1])= [y.sub.1] + [S1[S.sub.1]([y.sub.c]).sub.1] (2.5) [y.sub.m](x, [S.sub.0], ..., [S.sub.m-1]) = [y.sub.m-1] + [S.sub.m-1] [([y.sub.c]).sub.m-1]

d) Substitute the approximate solution [y.sub.m] into the Eq.(2.1) and the general problem results in the following residual:

R(x,[S.sub.0],...,[S.sub.m-1]) = A ([y.sub.m](x,[S.sub.0],...,[S.sub.m-1])) + g(x). (2.6)

Obviously, when R(x, [S.sub.0],..., [S.sub.m-1]) = 0 then the approximation [y.sub.m](x, [S.sub.0],..., [S.sub.m-1]) will be the exact solution. Generally it doesn't happen in nonlinear equations, but the functional can be minimized as:

[mathematical expression not reproducible] (2.7)

where a and b are selected from the domain of the problem. Optimum values of [S.sub.0], [S.sub.1],... can be obtained from the conditions

[[partial derivative]j/[[partial derivative]S.sub.0]] = [[partial derivative]j/[[partial derivative]S.sub.1]] = ... = [[partial derivative]j/[[partial derivative]S.sub.m-1]] = 0. (2.8)

The constants [S.sub.0], [S.sub.1],... can also be defined from

R([x.sub.0], [S.sub.i]) = R([x.sub.1], [S.sub.i]) = ... = R([x.sub.m-1], [S.sub.i])= 0, i = 0, 1,...,m - 1 (2.9)

where [x.sub.i] [member of] (a, b). Putting these constants into the last one of the Eqs. (2.5), the approximate solution of order m is determined. For much more information about finding these constants please see [25].

3 Applications

In this section, we apply OPIA to solve the Lane-Emden type problems. Obtained results show that the new algorithm gives better results when compared with many other methods in literature.

Example 3.1. Consider the following homogeneous Lane-Emden equation of the first kind:

y" + [2/x]y' + [y.sup.5] = 0, y(0) = 1, y'(0) =0, x [greater than or equal to] 0. (3.1)

is a basic equation in the theory of stellar structure [26, 27]. In [28], exact solution of (3.1) is given by

[mathematical expression not reproducible] (3.2)

Perturbation parameter [epsilon] can be inserted into (3.1) as:

F(y", y', y, [epsilon]) = y" + [epsilon][2/x]-y' + [[epsilon]y.sup.5] = 0 = L + N = 0 (3.3)

where

L = y" and N = [epsilon] ([2/x]y' + [y.sup.5]). (3.4)

For this problem, we do not have the term y" in the remaining part N. Therefore, Eq. (2.3) is simplified to:

N + [N.sub.y] ([y.sub.c])[n.sup.[epsilon]] + [N.sub.y'], ([y.sub.c'])[n.sup.[epsilon]] + [N.sub.[epsilon]][epsilon] = - L (3.5)

Using the Eqs. (2.2), (3.4), (3.5) and setting [epsilon] = 1 yields the following algorithm

[([y.sub.c])".sub.n] = -{[y.sub.n])"-[2/x]{[y.sub.n])' - [([y.sub.n]).sup.5]. (3.6)

One may select [y.sub.0] = 1 as a starting guess which satisfies the given initial conditions. Substituting [y.sub.0] into the Eq. (3.6), yields first order problem:

[([y.sub.c])".sub.n] = -1, y(0)= y'(0)= 0. (3.7)

Solving the (3.7) gives first correction term:

[([y.sub.c]).sub.0] = -[[x.sup.2]/2] (3.8)

Now, Eq.(3.8) is inserted into Eq.(2.4) to obtain first approximate solution:

[y.sub.1] = 1 - [[[S.sub.0][x.sup.2]/2] (3.9)

By using the procedure mentioned in Section 2, one obtains the following approximate solutions:

[mathematical expression not reproducible] (3.10)

[mathematical expression not reproducible] (3.11)

It should be emphasized that [y.sub.3] does not represent the third correction term; rather it is the approximate solution after the third iteration. Unknown constants can be determined from the residual

R(x, [S.sub.0], [S.sub.1], [S.sub.2]) = L ([y.sub.3](x, [S.sub.0], [S.sub.1], [S.sub.2])) + N ([y.sub.3](x, [S.sub.0], [S.sub.1], [S.sub.2])) (3.12)

for the third order approximation. Using the Eq. (2.9) with x = 1, 2, 3 yields

[S.sub.0] = 1.53406577, [S.sub.1] = 0.98031243, [S.sub.2] = -0.10588147 (3.13)

Substituting these constants into the Eq.(3.11), the approximate solution of the third order is obtained as:

[y.sub.3] (x) = 1 - 0.166667[x.sup.2] + 0.0416667[x.sup.4] - 0.0115741[x.sup.6] + 0.00337577[x.sup.8]0.00101275[x.sup.10] + 0.000309553[x.sup.12] - 0.0000961702[x.sup.14] + 0.0000305113[x.sup.16]7.757926108 x [10.sup.-6][x.sup.18] - 8.2105715765 x [10.sup.-6][x.sup.20] + 0.0000259462[x.sup.22] - 0.00002895[x.sup.24] - 3.5887084471241 x [10.sup.-7][x.sup.26] + 0.0000302234[x.sup.28] - 7.866052351076 x [10.sup.-6][x.sup.30] - 0.00003015[x.sup.32] + 4.4502399674 x [10.sup.-6][x.sup.34] + 0.00003140[x.sup.36] + 6.5466232556 x [10.sup.-6][x.sup.38] - 0.0000288707[x.sup.40] - 0.000020677[x.sup.42] + 0.000018640[x.sup.44] + 0.000031354[x.sup.46] + 2.724316366046 x [10.sup.-6][x.sup.48] - 0.0000345355[x.sup.50] - 0.000011951[x.sup.52] + 0.000034677[x.sup.54] + 0.000017377[x.sup.56] - 0.0000476036[x.sup.58 + 0.0000265864[x.sup.60]- 4.989939036859478 x [10.sup.-6][x.sup.62] (3.14)

It is obvious that as the number of iterations increase, the approximate solution becomes more and more complicated which requires symbolic computer programs. Mathematica 9.0 is used to deal with the complex calculations in this work. Repeating the same steps by using Mathematica, one can easily get higher order approximate solutions. Dehghan et al.[29] and Singh et al. [30] obtained approximate solutions for this problem using variational iteration method(VIM) and homotopy analysis method(HAM). Figure 3.1 shows that OPIM yields better results than those obtained by VIM and HAM.

Example 3.2. Consider the Lane-Emden type equation :

y" + [2/x]y' + [e.sup.y] = 0, y(0) = y'(0) = 0 (3.15)

which represents the isothermal gas spheres equation in the case that the temperature stays constant [15, 29, 31].

Reconsider the Eq. (3.15) as:

F(y", y',y, [epsilon]) = y" + ([2[epsilon]/x]y' +[e.sup.[epsilon]y]) =L+N=O (3.16)

where

L = y" and N = ([2[epsilon]/x]y' +[e.sup.[epsilon]y]). (3.17)

There is not y" term in part N. Thus, (2.3) reduces to:

N + [N.sub.y] ([y.sub.c])[n.sup.[epsilon]] + [N.sub.y'], ([y.sup.c'])[n.sup.[epsilon]] + [N.sub.[epsilon]][epsilon] = - L. (3.18)

With the help of the Eqs. (2.2), (3.17), (3.18) and setting [epsilon] = 1, OPIA is constructed as:

([y.sub.c])"n = -([y.sub.n])" - [2/x]([y.sub.n])' - [y.sub.n] -1, y(0) = y'(0) = 0. (3.19)

By choosing [y.sub.0] = 0 as a starting guess, (3.19) turns into first order problem:

([y.sub.c])"0 = -([y.sub.o])"- [2/x]([y.sub.n])'- [y.sub.n] -1, y(0) = y'(0) = 0. (3.20)

Solving (3.20) and using (2.4), one gets the following approximate solutions:

[y.sub.1] = -[1/2][S.sub.0][x.sup.2] (3.21)

[y.sub.2] = -[1/2][S.sub.0][x.sup.2] + [1/24] ([S.sub.1]) [x.sup.2] ([S.sub.0][x.sup.2] + 36[S.sub.0] - 12) (3.22)

[mathematical expression not reproducible] (3.23)

[mathematical expression not reproducible] (3.24)

Using the Eq. (2.9) with the residual

R(x, [S.sub.0], [S.sub.1], [S.sub.2], [S.sub.3]) = L ([y.sub.4](x, [S.sub.0], [S.sub.1], [S.sub.2], [S.sub.3])) + N ([y.sub.4](x, [S.sub.0], [S.sub.1], [S.sub.2], [S.sub.3])) (3.25)

the unknown constants are obtained as

[S.sub.0] = 2.0203622551, [S.sub.1] = -1.0201147822, [S.sub.2] = -0.9963202221, [S.sub.3] = 0.020789994 (3.26)

for the fourth order approximation. Substituting the Eq.(3.26) into (3.24) yields:

[y.sub.4] (x) = -0.15989962328[x.sup.2] + 0.007920058473[x.sup.4]+ - 0.005608897215631[x.sup.6] + 0.000039055217456[x.sup.8]. (3.27)

In [26, 31], the authors have used differential transform method (DTM) and Adomian decomposition method(ADM) to solve the (3.15) and they obtained the series solution:

[mathematical expression not reproducible] (3.28)

Figure 3.2 and Figure 3.3 are sketched for comparison of the higher order approximate solutions of OPIM and ADM-DTM solutions. It is also clear from the figures that OPIM solutions are valid in larger region.

Example 3.3. Consider the following homogeneous Lane-Emden type equation [32]:

y" + [2/x]-y' - [4[x.sup.2] + 6)y = 0, y (0) = 1, y'(0) = 0, 0 [less than or equal to] x [less than or equal to] 1. (3.29)

Exact solution ofthis problem is given as

[mathematical expression not reproducible] (3.30)

Rewrite the Eq. (3.29) as:

F{y",y', y,[epsilon]) = y" + [epsilon] ([2/x]y' - ([4x.sup.2] + 6)y) = L + N = 0 (3.31)

where

L = y" and N = [epsilon] ([2/x]y' - ([4x.sup.2] + 6)y). (3.32)

y" does not appear in the remaining part N. Therefore, Eq. (2.3) becomes:

N + [N.sub.y] ([y.sub.c])[n.sub.[epsilon]] + [N.sub.y'], ([y.sub.c'])[n.sup.[epsilon]] + [N.sub.[epsilon]][epsilon] = - L (3.33)

In the light of previous examples, OPIA can be established as

[([y.sub.c])".sub.n] = -{[y.sub.n])" -[2/x]([y.sub.n])' + ([4x.sup.2] + 6)[y.sub.n]. (3.34)

by using the Eqs. (2.2), (3.32)and (3.33). Accordingly, first order problem is obtained as:

[([y.sub.c])".sub.n] = ([4x.sup.2] + 6), y(0) = y'(0) = 0 (3.35)

Starting with the initial condition [y.sub.0] = 1, iterations can be reached as:

[y.sub.1] = 1 + [1/3][S.sub.0] ([x.sup.4] + [9x.sup.2]) (3.36)

[mathematical expression not reproducible] (3.37)

[mathematical expression not reproducible] (3.38)

Unknown constants can be determined from the residual

R(x, [S.sub.0], [S.sub.1], [S.sub.2]) = L ([y.sub.3](x, [S.sub.0], [S.sub.1], [S.sub.2])) + N ([y.sub.3](x, [S.sub.0], [S.sub.1], [S.sub.2])) (3.39)

for the third order approximation. Using the Eq. (2.9) with x = 0.3,0.6,0.9 yields

[S.sub.0] = 0.33423439, [S.sub.1] = 0.31859877, [S.sub.2] = 0.20764389 (3.40)

Substituting these constants into the Eq.(3.38), the approximate solution of the third order is obtained as:

[y.sub.3] (x) = 0.000135466[x.sup.12] + 0.00419221[x.sup.10] + 0.041724[x.sup.8] + 0.185681[x.sup.6]+ 0.483647[x.sup.4] + 1.0041[x.sup.2] + 1. (3.41)

Repeating the same steps one can get the following approximations:

[y.sub.4] (x) = 0.0000407289[x.sup.16] + 0.000166717[x.sup.14] + 0.00142143[x.sup.12] + 0.00831423[x.sup.10] + 0.0416732[x.sup.8] + 0.166665[x.sup.6] + 0.5[x.sup.4] + [x.sup.2] + 1. (3.42)

[y.sub.5] (x) = 4.529285401255387 x [10.sup.-7][x.sup.20] + 2.306750146551151 x [10.sup.-6][x.sup.18]+ 0.0000254168[x.sup.16] + 0.000197899[x.sup.14] + 0.00138916[x.sup.12]+ 0.00833324[x.sup.10] + 0.0416667[x.sup.8] + 0.166667[x.sup.6] + 0.5[x.sup.4] + [x.sup.2] + 1. (3.43)

This problem has been also investigated by Ozis et al using variational iteration method(VIM) and homotopy perturbation method (HPM) [32, 33]. Figures 3.4, 3.5 and Table 1 give important information on the convergence and the absolute errors for OPIA and other approximate solutions in literature. It is clear that the results obtained by OPIM are more accurate than those in [32, 33].

4 Conclusion

In this paper, OPIM is applied for the first time to investigate the new approximate solutions for Lane-Emden type differential equations. This new technique provides us to optimally control the convergence of solution series. Also, it gives a very good approximation even in a few terms to these kinds of nonlinear equations. The results obtained in this paper proves that the OPIM is a very effective technique for differential equations. It is worth mentioning that, a symbolic program is necessary for successive calculations after a few iterations. Mathematica 9 has been used to overcome the complicated calculations for this present research.

ACKNOWLEDGEMENT : The authors would like to thank The Scientific and Technological Research Council of Turkey (TUBITAK) for their financial support.

References

[1] Abdon Atangana and Adem Kilicman. Analytical solutions of boundary values problem of 2d and 3d Poisson and biharmonic equations by homotopy decomposition method. In Abstract and Applied Analysis, volume 2013. Hindawi Publishing Corporation, 2013.

[2] Zehra Pinar, Abhishek Dutta, Guido Beny, and Turgut Ozis. Analytical solution of population balance equation involving growth, nucleation and aggregation in terms of auxiliary equation method. Applied Mathematics & Information Sciences, 9(5):2467, 2015.

[3] Mehmet Giyas Sakar, Fatih Uludag, and Fevzi Erdogan. Numerical solution of time-fractional nonlinear PDEs with proportional delays by homotopy perturbation method. Applied Mathematical Modelling, 2016.

[4] Asma Ali Elbeleze, Adem Kilicman, and Bachok M Taib. Homotopy perturbation method for fractional Black-Scholes european option pricing equations using Sumudu transform. Mathematical Problems in Engineering, 2013, 2013.

[5] N Bildik and S Deniz. Comparison of solutions of systems of delay differential equations using Taylor collocation method, Lambert W function and variational iteration method. Scientia Iranica. Transaction D, Computer Science & Engineering, Electrical, 22(3):1052, 2015.

[6] Hasan Bulut, Haci Mehmet Baskonus, and Fethi Bin Muhammad Belgacem. The analytical solution of some fractional ordinary differential equations by the Sumudu transform method. In Abstract and Applied Analysis, volume 2013. Hindawi Publishing Corporation, 2013.

[7] Saeid Abbasbandy. Improving Newton--Raphson method for nonlinear equations by modified Adomian decomposition method. Applied Mathematics and Computation, 145(2):887-893, 2003.

[8] Sinan Deniz and Necdet Bildik. Application of Adomian decomposition method for singularly perturbed fourth order boundary value problems. In International Conference of Numerical Analysis and Applied Mathematics 2015 (ICNAAM 2015), volume 1738, page 290017. AIP Publishing, 2016.

[9] A Sami Bataineh, Mohd Salmi Md Noorani, and Ishak Hashim. Homotopy analysis method for singular IVPs of Emden--Fowler type. Communications in Nonlinear Science and Numerical Simulation, 14(4):1121-1131, 2009.

[10] K Parand, AR Rezaei, and A Taghavi. Lagrangian method for solving Lane--Emden type equation arising in astrophysics on semi-infinite domains. Acta Astronautica, 67(7):673-680, 2010.

[11] Afgan Aslanov. Determination of convergence intervals of the series solutions of Emden--Fowler equations using polytropes and isothermal spheres. Physics Letters A, 372(20):3555-3561, 2008.

[12] HR Marzban, HR Tabrizidooz, and M Razzaghi. Hybrid functions for nonlinear initial-value problems with applications to Lane--Emden type equations. Physics Letters A, 372(37):5883-5886, 2008.

[13] Shijun Liao. A new analytic algorithm of Lane--Emden type equations. Applied Mathematics and Computation, 142(1):1-16, 2003.

[14] Ji-Huan He. Variational approach to the Lane--Emden equation. Applied Mathematics and Computation, 143(2):539-541, 2003.

[15] K Parand, Mehdi Dehghan, AR Rezaei, and SM Ghaderi. An approximation algorithm for the solution of the nonlinear Lane--Emden type equations arising in astrophysics using Hermite functions collocation method. Computer Physics Communications, 181(6):1096-1108, 2010.

[16] James SW Wong. On the generalized Emden-Fowler equation. Siam Review, 17(2):339-360, 1975.

[17] Vasile Marinca, Nicolae Herisanu, Constantin Bota, and Bogdan Marinca. An optimal homotopy asymptotic method applied to the steady flow of a fourth-grade fluid past a porous plate. Applied Mathematics Letters, 22(2): 245-251, 2009.

[18] Vasile Marinca and Nicolae Herisanu. Application of optimal homotopy asymptotic method for solving nonlinear equations arising in heat transfer. International Communications in Heat and Mass Transfer, 35(6):710-715, 2008.

[19] Yigit Aksoy and Mehmet Pakdemirli. New perturbation-iteration solutions for Bratu-type equations. Computers & Mathematics with Applications, 59(8): 2802-2808, 2010.

[20] Vasile Marinca, Nicolae Herisanu, and Iacob Nemes. Optimal homotopy asymptotic method with application to thin film flow. Open Physics, 6(3): 648-653, 2008.

[21] Vasile Marinca and Nicolae Herisanu. Determination of periodic solutions for the motion of a particle on a rotating parabola by means of the optimal homotopy asymptotic method. Journal ofSound and Vibration, 329(9):1450-1459, 2010.

[22] Yigit Aksoy, Mehmet Pakdemirli, Saeid Abbasbandy, and Hakan Boyaci. New perturbation-iteration solutions for nonlinear heat transfer equations. International Journal of Numerical Methods for Heat & Fluid Flow, 22(7):814-828, 2012.

[23] Mehmet Senol, Ihsan Timucin Dolapci, Yigit Aksoy, and Mehmet Pakdemirli. Perturbation-iteration method for first-order differential equations and systems. In Abstract and Applied Analysis, volume 2013. Hindawi Publishing Corporation, 2013.

[24] Necdet Bildik and Sinan Deniz. Comparative study between optimal homotopy asymptotic method and perturbation-iteration technique for different types of nonlinear equations. Iranian Journal of Science and Technology, Transactions A: Science, pages 1-8.

[25] N Herisanu and Vasile Marinca. Accurate analytical solutions to oscillators with discontinuities and fractional-power restoring force by means of the optimal homotopy asymptotic method. Computers & Mathematics with Applications, 60(6):1607-1615, 2010.

[26] Abdul-Majid Wazwaz and Randolph Rach. Comparison of the Adomian decomposition method and the variational iteration method for solving the Lane-Emden equations of the first and second kinds. Kybernetes, 40(9/10): 1305-1318, 2011.

[27] Subrahmanyan Chandrasekhar and Subrahmanyan Chandrasekhar. An introduction to the study of stellar structure, volume 2. Courier Corporation, 1958.

[28] Harold Thayer Davis. Introduction to nonlinear differential and integral equations. Courier Corporation, 1962.

[29] Mehdi Dehghan and Fatemeh Shakeri. Approximate solution of a differential equation arising in astrophysics using the variational iteration method. New Astronomy, 13(1):53-59, 2008.

[30] Om P Singh, Rajesh K Pandey, and Vineet K Singh. An analytic algorithm of Lane--Emden type equations arising in astrophysics using modified homotopy analysis method. Computer Physics Communications, 180(7):1116-1124, 2009.

[31] Yasir Khan, Zdenek Svoboda, and Zdenek Smarda. Solving certain classes of Lane-Emden type equations using the differential transformation method. Advances in Difference Equations, 2012(1):1-11, 2012.

[32] Ahmet Yildirim and Turgut Ozis. Solutions of singular IVPs of Lane--Emden type by the variational iteration method. Nonlinear Analysis: Theory, Methods & Applications, 70(6):2480-2484, 2009.

[33] Ahmet Yildirim and Turgut Ozis. Solutions of singular IVPs of Lane--Emden type by homotopy perturbation method. Physics Letters A, 369(1):70-76, 2007.

Department of Mathematics, Faculty of Art and Sciences, Manisa Celal Bayar University, 45140 Manisa, Turkey.

e-mails: sinandeniz01@gmail.com, n.bildik@cbu.edu.tr

Sinan Deniz

Necdet Bildik (*)

(*) Corresponding author

Received by the editors in June 2016 - In revised form in December 2016.

Communicated by K. In't Hout.

2010 Mathematics Subject Classification : 34E20, 34A12, 85A04.

Key words and phrases : Optimal perturbation iteration method, singular initial value problems, Lane-Emden equation.
Table 1: Comparison of absolute errors for Example 3 at different
orders of approximations.

x                             Errors for OPIA
|[Y.sub.exact]-[Y.sub.5]]   |[Y.sub.exact]-[Y.sub.6]]

0.1   3.08426x[10.sup.-13]        2.22045x[10.sup.-16]
0.2   3.85914x[10.sup.-13]        2.22045x[10.sup.-16]
0.3   5.62883x[10.sup.-13]        1.00025x[10.sup.-17]
0.4   9.64347x[10.sup.-13]        2.44249x[10.sup.-15]
0.5   1.95532x[10.sup.-12]        1.50998x[10.sup.-14]
0.6   4.77649x[10.sup.-12]        1.01037x[10.sup.-13]
0.7   5.45375x[10.sup.-11]        5.02709x[10.sup.-13]
0.8   2.78031x[10.sup.-11]        2.05902x[10.sup.-12]
0.9   3.08426x[10.sup.-10]        6.38223x[10.sup.-12]
1     2.10201x[10.sup.-9]         2.90235x[10.sup.-12]

x                   Errors for VIM-HPM
|[Y.sub.exact]-[Y.sub.5]]   |[Y.sub.exact]-[Y.sub.6]]

0.1   1.11022x[10.sup.-15]        1.00128x[10.sup.-16]
0.2   5.72165x[10.sup.-12]        3.26406x[10.sup.-14]
0.3   7.47710x[10.sup.-10]        9.59788x[10.sup.-12]
0.4   2.38451x[10.sup.-8]         5.43454x[10.sup.-10]
0.5   3.51584x[10.sup.-7]         1.24994x[10.sup.-8]
0.6   3.18608x[10.sup.-6]         1.62772x[10.sup.-7]
0.7   0.0000206568                1.43282x[10.sup.-6]
0.8   0.000104921                 9.47740x[10.sup.-6]
0.9   0.000442699                 0.000050436
1     0.00161516                  0.000226273