Many problems in mathematical physics and astrophysics can be modeled by the so-called Lane-Emden type equation defined in the form
With the initial conditions:
where “ ” is the differentiation with respect to x, is a nonlinear function of x and y, is the non-homogeneous term and are constants. The closed form solution of the Lane-Emden type Equation (1) is always enabled  in the neighborhood of the point , which is called a singular point, for the above initial conditions.
Taking and in Equations (1) and (2) respectively, we get
subject to initial conditions:
Equations (3) and (4) are known as the classical Lane-Emden equation.
Choosing and in Equations (1) and (2) respectively, we get the so-called white-dwarf equation that expresses the gravitation potential of the degenerate white-dwarf stars. Isothermal gas spheres  are similarly modeled by
with the initial conditions:
The parameter k in Equation (3) has physical significance in the range so Equation (3) has a closed form solutions for  but may be resorted to numerical solutions for other values of k but the singularity at is still a challenge to the numerical solutions. The semi-analytic solutions can be found by Adomian’s Decomposition Method (ADM) and perturbation techniques and mostly convergent in restricted regions so another approach such as Pade’s method was required to expand these regions    . Many semi-analytical methods and numerical methods like Homotopy analysis method    , Picard method  and spectral methods    were used to solve different types of differential equations. Many researches solved the Lane-Emden equation; like in  the authors applied Variational Iteration Method to Lane-Emden equation using some transformation. Singh et al.  suggested an effective analytic algorithm using Modified Homotopy Analysis Method (MHAM), at which convergence regions could be adjusted without using Pade’s technique. Also another algorithm was proposed in  , using an iterative method which is a hybrid of ADM and Variational Iteration Method (VIM). Lastly, Parand et al.  proposed an approximation algorithm using Hermite collocation method at which the solution is reduced to the solution of a system of algebraic equations.
In the present paper, we aim effectively to employ Variation of Parameters Method (VPM) coupling with an unknown auxiliary parameter h to solve the homogeneous as well as non-homogeneous Lane-Emden equation. VPM is free from calculation of the so-called Adomian’s polynomials and is widely applicable because of the reliability of it and the reduction in the size of the computational domain. A comparison with the standard VPM is made to show the reliability and efficiency of the applied algorithm which is simple and is accurately approximates the solution over a large domain.
2. Variation of Parameters Method with an Auxiliary Parameter
Variation of Parameters Method (VPM) was first proposed by Ma    and it has been thoroughly used by many researchers to solve a wide class of initial and boundary value problems     . Ghaneai and Hosseini  inserted an auxiliary parameter in variational iteration algorithm to obtain solution of wave-like and heat-like equations and concluded that an auxiliary parameter provides a simple way to control and adjust the convergence region of approximate solution in a large domain. In  VIM coupled with an auxiliary parameter are used to predict the multiple solutions of nonlinear boundary value problems with great success. Also an auxiliary parameter is inserted in VPM to solve some equations arising in physics   . To convey the basic step of VPM for solving Lane-Emden equation, consider the general nonlinear ordinary differential equation in operator form as follows:
with the initial conditions:
According to the Variation of Parameters Method  ,  ,  and  , the solution of the homogeneous Equation (11):
is given as
For the particular solution, the constants and in Equation (13) are replaced by the functions and respectively, thus we have
where and . The functions and can be determined as
where W is the Wronskian of and i.e.
Substituting the values of and in Equation (15), yields
Hence, the general solution of Equation (12) is
Applying the initial conditions (10) in Equation (20) and solve for and , we get and .
Substituting the values of and in Equation (20) yields
which can be solved iteratively as
Consequently, an exact solution can be obtained when
The valid values of h can be determined by means of the so-called h-curves    . According to this h-curve, it can be easily determine the valid region of h, which corresponds to the line segment nearly parallel to the horizontal axis. Note that the Variation of Parameters Method     provides the following iterative scheme for Equation (9):
3. Approximate Solution of the Lane-Emden Equation
Now we apply Variation of Parameters Method (VPM) coupled with an unknown auxiliary parameter developed in Section 2 for solving Equation (1). The recursive formulas (22) and (24) for this problem is obtained respectively as
Consider the following form of Lane-Emden equation    :
where k is a constant. With the initial conditions:
Note that the exact solutions for cases of are known, VPM developed in Section 2 is applied to solve Equations (26) and (27) for some cases of the constant k.
Case (1) for (k = 0):
The closed form solution is: . According to standard VPM, we have the following from iterative scheme (24):
The result obtained by standard VPM is not valid for large values of x. Now, using the iterative scheme (22) we have:
First few terms of solution are
It can be seen in Figure 1 that the admissible range of h is . The h curve of for 20th order approximation when , absolute error for15th order approximation by standard VPM and absolute error when are shown in Figure 1.
Figure 1. For (k = 0) in Example 1 (a) the h curve for the 20th order approximation; (b) Absolute error by the standard VPM; (c) Absolute error when .
Case (2) for (k = 1):
The closed form solution is . According to standard VPM, the absolute error of is shown in Figure 2(a) which also shows that the result is not valid for large value of x. Now, using the iterative scheme (22) we have the first few terms of solution:
Figure 2. For (k = 1) in Example 1 (a) Absolute error by the standard VPM; (b) Absolute error when ; (c) the h curve for the 30th order approximation; (d) Series and exact solutions.
Other terms of the solution can also be obtained in a similar way.
From Figure 2(c) the admissible range is and after plotting a number of h curves we get the best value of h which is equal to −0.35, therefore the approximate solution is given by
An approximate solution (not exact) was obtained, using VIM, by  and is valid only for . The solution is convergent in the interval as shown in Figure 2(d). The accuracy of absolute error is remarkably improved for 30th order approximation, as illustrated in Figure 2(b).
Case (3) for (k = 2):
In this case we analyze the dependence of the convergence regions on the value of , The h curve for the ninth-order approximation at is shown in Figure 3, evaluating . Figure 4 shows the solution at respectively, it is obvious that when the value of increases monotonically, the region of convergence monotonically decreases.
Case (4) for (k = 3, k = 4):
For these cases, a comparison is made through Table 1 and Table 2, respectively between the solution of Equation (26) and that obtained by  . It is shown that the obtained solution is accurate. h-curves at are shown in Figure 5(a) and Figure 5(b), respectively. The obtained results by the used method are compared with standard VPM through Table 3 and Table 4 respectively. It is clear from tables that the standard VPM is not efficient for large values of x and insertion of h to the method is very effective step in the solution.
Case (5) for (k = 5):
The closed form solution is . Adopting the similar
Figure 3. The h curve for (k = 2) in Example 1 for the 9th order approximation.
Figure 4. For (k = 2) in Example 1 (a) solution with ; (b) solution with ; (c) solution with ; (d) solution with .
Figure 5. F or Example 1 (a) the h curve for the seventh-order approximation for (k = 3); (b) the h curve for the fifth-order approximation for (k = 4).
Table 1. Exact value, , the absolute error and the error obtained by  for .
Table 2. Exact value, , the absolute error and the error obtained by  for .
Table 3. Solutions obtained by the standard VPM and VPM with an auxiliary parameter for .
Table 4. Solutions obtained by the standard VPM and VPM with an auxiliary parameter for .
procedure as in the previous examples, Figure 6(a) shows the absolute error of which confirm that the obtained result by standard VPM is not valid for large value of x. The first few terms of solution are
Other terms of the solution can also be obtained in a similar way, h can be chosen in the range
, from Figure 6(b), the h curve for
, the value of h equal to
For the Isothermal gas spheres equation:
with the initial conditions ( ), we can use Taylor series expansion of
Figure 6. For (k = 5) in Example 1 (a) Absolute error by the standard VPM; (b) The h curve; (c) Absolute error when ; (d) Series and exact solutions.
Table 5 shows the results obtained with standard VPM and VPM with h. Now applying the iterative scheme (22) we have:
Other terms of the solution can also be obtained in a similar way to get:
(with h equal to −0.395) which is the same as the solution obtained by  using ADM, and  using MHAM. Figure 7 shows the h curve for . Table 5 shows the comparison of obtained here and those obtained by  and  . The comparison between the solutions obtained by the standard VPM and VPM with an auxiliary parameter for isothermal gas sphere equation is tabulated in Table 6.
Figure 7. The h curve for isothermal gas sphere equation in Example 2.
Table 5. Comparison of between present method and series solution given by  and  for isothermal gas sphere equation in Example 2..
Table 6. Comparison between the solutions obtained by the standard VPM and VPM with an auxiliary parameter for isothermal gas sphere equation in Example 2.
Consider the following problem  :
With the initial conditions
having as exact solution. Now, using the iterative scheme (22) we have
Hence, which is the exact solution of Equation (45) with as the range of h is as shown in Figure 8(a). The exact solution and series solution for ( ) absolute error by standard VPM and by VPM with h are shown in Figures 8(b)-(e), respectively. Figure 8(c) shows that VPM with an auxiliary parameter gives better approximation than MHAM  . The obtained solution, in comparison to that obtained by  , is very accurate. Table 7 shows the comparison of obtained here and those obtained by  .
In the next examples, the nonhomogeneous Lane-Emden equations are considered.
With the initial conditions
Which has the following exact solution: . Now using the iterative scheme (22) we have:
Figure 8. For Example 3 (a) The h curve; (b) Exact solution; (c) Series solution with ; (d) Absolute error by standard VPM; (e) Absolute error by VPM with .
Table 7. Comparison between the solutions obtained by the present method and series solution given by Parand et al.  for Example 3.
After number of iterations we get an approximate solution for h in the range of at which by choosing an auxiliary parameter . Figure 9(a) shows the h curve at . Exact solution and series solution are shown in Figure 9(b) and Figure 9(c) respectively. Like all previous examples, the present method is more accurate than standard VPM especially for large value of x. The absolute error by standard VPM is shown in Figure 9(d) and the absolute error by VPM with is shown in Figure 9(e) which shows that the solution obtained here is more accurate in comparison to that obtained with MHAM  and also to that obtained by  . Table 8 shows the comparison of obtained after 20 approximations and those obtained by  .
Consider the equation:
with the initial conditions
Having as exact solution, using the iterative scheme (18) we have:
Figure 9. For Example 4 (a) the h curve; (b) Exact solution; (c) series solution when ; (d) absolute error by standard VPM; (e) absolute error by VPM when .
The exact solution ( ) can be obtained when n approaches to infinity by choosing a suitable value for h in the range which can be obtained from the h curve at shown in Figure 10(a), exact and series solutions with are shown in Figure 10(b).
Figure 10. For Example 5 (a) The h curve; (b) Exact and series solution.
Table 8. Comparison of given by present method at ( ) and series solution given by  for Example 4.
In this article, Variation of Parameter Method is used coupled with an auxiliary parameter for the Lane-Emden equation. Moreover, the iterative scheme (22) reduces to standard VPM for . The numerical results obtained in this research are very effective, completely reliable and powerful in obtaining analytical solutions of many problems. The approach used in the solution provides a high degree of accuracy only with a few iterations and reduces the size of calculations without using any restrictive assumptions.
 Semary, M.S. and Hassan, H.N. (2012) Series Solution for Continuous Population Models for Single and Interacting Species by the Homotopy Analysis Method. Communications in Numerical Analysis, 2012, Article ID: cna-00106.
 Semary, M.S. and Hassan, H.N. (2016) An Effective Approach for Solving MHD Viscous Flow Due to a Shrinking Sheet. Applied Mathematics and Information Sciences, 10, 1425-1432.
 Semary, M.S., Hassan, H.N. and Radwan, A.G. (2017) Single and Dual Solutions of Fractional Order Differential Equations Based on Controlled Picard’s Method with Simpson Rule. Journal of the Association of Arab Universities for Basic and Applied Sciences, 24, 247-253.
 Hassan, H.N. (2016) An Accurate Numerical Solution for the Modified Equal Width Wave Equation Using the Fourier Pseudo-Spectral Method. Journal of Applied Mathematics and Physics, 4, 1054-1067.
 Hassan, H.N. (2017) An Efficient Numerical Method for the Modified Regularized Long Wave Equation Using Fourier Spectral Method. Journal of the Association of Arab Universities for Basic and Applied Sciences, 24, 198-205.
 Hassan, H.N. and Saleh, H.K. (2010) The Solution of the Regularized Long Wave Equation Using the Fourier Leap-Frog Method. Zeitschrift für Naturforschung A, 65, 268-276.
 Dehghan, M. and Shakeri, F. (2008) Approximate Solution of a Differential Equation Arising in Astrophysics Using Variational Iteration Method. New Astronomy, 13, 53-59.
 Singh, O.P., Pandey, R.K. and Singh, V.K. (2009) An Analytic Algorithm of Lane-Emden Type Equations Arising in Astrophysics Using Modified Homotopy Analysis Method. Computer Physics Communications, 180, 1116-1124.
 Baranwal, V.K., Pandey, R.K., Tripathi, M.P. and Singh, O.P. (2012) An Analytic Algorithm of Lane-Emden-Type Equations Arising in Astrophysics—A Hybrid Approach. Journal of Theoretical and Applied Physics, 6, 2251-7235.
 Parand, K., Dehghan, M., Rezaei, A.R. and Ghaderi, S.M. (2010) 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, 1096-1108.
 Ma, W.-X. and You, Y.-C. (2005) Solving the Korteweg-de Vries Equation by Its Bilinear Form: Wronskian Solutions. Transactions of the American Mathematical Society, 357, 1753-1778.
 Ma, W.-X., Wu, H.-Y. and He, J.-S. (2007) Partial Differential Equations Possessing Frobenius Integrable Decompositions. Physics Letters A, 364, 29-32.
 Mohyud-Din, S.T., Noor, M.A. and Noor, K.I. (2009) Modified Variation of Parameters Method for Second-Order Integro-Differential Equations and Coupled Systems. World Applied Sciences Journal, 6, 1139-1146.
 Ghaneai, H. and Hosseini, M.M. (2015) Variational Iteration Method with an Auxiliary Parameter for Solving Wave-Like and Heat-Like Equations in Large Domains. Computers & Mathematics with Applications, 69, 363-373.
 Semary, M.S. and Hassan, H.N. (2015) A New Approach for a Class of Nonlinear Boundary Value Problems with Multiple Solutions. Journal of the Association of Arab Universities for Basic and Applied Sciences, 17, 27-35.
 Sikandar, W., Khan, U., Ahmed, N. and Mohyud-Din, S.T. (2018) Variation of Parameters Method with an Auxiliary Parameter for Initial Value Problems. Ain Shams Engineering Journal, 9, 1959-1963.
 Sikandar, W., Khan, U., Ahmed, N. and Mohyud-Din, S.T. (2017) Optimal Solutions for Homogeneous and Non-Homogeneous Equations Arising in Physics. Results in Physics, 7, 216-224.
 Ramos, J.I. (2008) On the Variational Iteration Method and Other Iterative Techniques for Nonlinear Differential Equations. Applied Mathematics and Computation, 199, 39-69.