Bernstein Series Solution of a Class of Lane-Emden Type Equations

The purpose of this study is to present an approximate solution that depends on collocation points and Bernstein polynomials for a class of Lane-Emden type equations with mixed conditions. The method is given with some priori error estimate. Even the exact solution is unknown, an upper bound based on the regularity of the exact solution will be obtained. By using the residual correction procedure, the absolute error can be estimated. Also, one can specify the optimal truncation limit n which gives a better result in any norm. Finally, the effectiveness of the method is illustrated by some numerical experiments. Numerical results are consistent with the theoretical results.


Introduction
Lane-Emden type equation that is presented in (1) models many phenomena in mathematical physics and astrophysics [1,2].Consider   () + 2    () +  () = 0,  > 0,   (0) = 0,  (0) = ,  is a constant. ( It describes the equilibrium density distribution in selfgravitating sphere of polytrophic isothermal gas [3].On the other hand [3], it plays an important role in various fields such as stellar structure [2], radiative cooling, and modeling of clusters of galaxies.It is a nonlinear ordinary differential equation that has a singularity at the origin.In the neighborhood of  = 0, it has an analytic solution [1].
It is labeled by the names of the astrophysicists Lane [4] and Robert Emden.

Recent Works.
Recently, a number of numerical methods are used for handling the Lane-Emden type problems based on perturbation techniques or series solutions.Adomian decomposition method [6,7] which provides a convergent series solution has been used to solve (1) [8][9][10].Wazwaz [8] gave an algorithm to overcome the difficulty of the singular point in using Adomian decomposition method [1].
The quasilinearization method [11][12][13] can be considered as an example for iteration methods.Its fast convergence, Mathematical Problems in Engineering monotonicity, and numerical stability were analyzed by Krivec and Mandelzweig [12].They verified this method on scattering length calculations in the variable phase approach to quantum mechanics.They also showed that the iterations converge uniformly and quadratically to the exact solution.The method gives accurate and stable answers for any coupling strengths, including super singular potentials for which each term of the perturbation theory diverges.
The Legendre wavelet method was given by Yousefi [14] to solve Lane-Emden equation.This method was used to convert Lane-Emden equations to integral equations and was expanded the solution by Legendre wavelets with unknown coefficients.
Ramos [15] applied a piecewise linearization method to solve the Lane-Emden equation.This method provided piecewise linear ordinary differential equations that can be easily integrated.Furthermore, it has given accurate results for hypersingular potentials, for which perturbation methods diverge.Homotopy analysis method (HAM) and modified HAM have also been used [16,17] to solve (1).Parand et al. [18] proposed a collocation method based on a Hermite function collocation (HFC) method for solving some classes of Lane-Emden type equations which are nonlinear ordinary differential equations on the semi-infinite domain.A matrix method was given by Yuzbasi for solving nonlinear Lane-Emden type equations.Moreover, Yuzbasi and Sezer [5] applied a matrix method that depends on Bessel polynomials to solve (2).They estimated the absolute errors by using the residual correction procedure.In this study, a similar method to [5] was constructed.In addition, error analysis of the matrix method was developed.
In 2012, Pandey and coworkers [19][20][21][22] studied five methods.First, Pandey et al. [19] gave a numerical method for solving linear and nonlinear Lane-Emden type equations using Legendre operational matrix of differentiation.Second, Pandey et al. [20] studied a numerical method to solve linear and nonlinear Lane-Emden type equations using Chebyshev wavelet operational matrix.Third, Kumar et al. [21] presented a method for linear and nonlinear Lane-Emden type equations using the Bernstein polynomial operational matrix of integration.Fourth, Pandey and Kumar [22] proposed a numerical method for solving Lane-Emden type equations arising in astrophysics using Bernstein polynomials.This method is similar to the method used in the present study.And finally, a shifted Jacobi-Gauss collocation spectral method was proposed by Bhrawy and Alofi [23] for solving the nonlinear Lane-Emden type equation.
This paper is organized as follows.In Section 2, some definitions and theorems are given.The method is presented in Section 3. First, a matrix form for each term in (2) is found.Substituting these matrix forms into (2) gives a matrix equation, fundamental matrix equation.Then, a linear system by using collocation points is obtained.For the error analysis, in Section 4, some theorems that give some upper bounds for the absolute errors are presented.One of them guarantees the convergence if the solution is polynomial.The second one gives an upper bound in the case of the exact solution being unknown under the regularity condition.The residual correction procedure to estimate the absolute errors is also given so that the optimal truncation limit  can be specified.On the other hand, this procedure gives a new approximate solution.Some numerical examples are given to illustrate the method.

Preliminaries
Bernstein polynomials of th-degree are defined by where  is the maximum range of the interval [0, ] over which the polynomials are defined to form a complete basis [24].We substitute the relation into ( 5) and obtain the relation Let us consider  + 1 pairs (  ,   ).The problem is to find a polynomial   , called interpolating polynomial, such that The points   are called interpolation nodes.If  ̸ = , the problem is over-or underdetermined.
Theorem 2 (see [25]).Let  0 ,  1 , . . .,   be  + 1 distinct nodes, and let  be a point belonging to the domain of a given function .Assume that  ∈  +1 (  ), where   is the smallest interval containing the nodes  0 ,  1 , . . .,   and .Then, the interpolation error at the point  is given by where  ∈   .
The interpolation error may be reduced by using the roots of Chebyshev polynomials

Fundamental Relations
Let   be Bernstein series solution of (2).Let us find the matrix forms of   and  ()  .  can be written as where Therefore,  ()  can be written as On the other hand, B ()  () can be written as [26][27][28]] where For X () (), the relation is obtained where Substituting ( 24) into (22) yields Putting ( 26) into (19) yields the matrix form for  ()  as By substituting (19) and ( 27) into (2), we obtain a matrix equation as By using the collocation points 0 <  0 <  1 < . . .<   ≤  in (28), one obtains the fundamental matrix equation ] , We can obtain the corresponding matrix form for conditions (3), by means of the relation (27), as follows: On the other hand, the matrix forms for the conditions can be written as where Replacing the condition matrices (31) by any two rows of [W, G], we get the augmented matrix as [ W, G].Let the collocation points be selected such that the rank of W is  + 1.Therefore, the unknown matrix A is obtained as

Error Analysis and Estimation of the Absolute Error
In this section, some upper bounds of the absolute error are given by using Lagrange and Hermite interpolation polynomials.Also, an estimation of the error based on residual correction is given.
Theorem 3 (see [29]).Let  be a nonsingular matrix and  ̸ = 0 a vector.If  and x =  +  are, respectively, the solutions of the systems  =  and x =  + , one has Let  be the exact solution of (2) and    the interpolation polynomial of it on the nodes { 0 ,  1 , . . .,   }.If  ∈  +1 [0, ], then we can write  as  =    +   , where If   is the Bernstein series solution of (2), then it satisfies (2) on the nodes.So,   and    are the solutions of WA = G and WÂ = G + ΔG, respectively, where Proof.Since the exact solution is polynomial, for  ≥ deg(),   () = 0; the right hand side of (37) is zero.
The following theorem can be used for the estimation of the absolute error when the exact solution is unknown.Hence, an upper bound depending on Theorem 6.Let   and  be Bernstein series solution and the exact solution of (2), respectively.Let the interpolation nodes contain 0 and .Let  ∈ The first term on the right hand side can be bounded by (17) since ∈  (3) [0, ].
If the exact solution is unknown, the following steps can be used to find an upper bund of the absolute error.First, we construct the differential equation of   .If  3−1 ∈  3−1 is the Hermite interpolation polynomial of  on the nodes {  1 ,   2 , . . .,   −2 } ∪ {0, } ⊂ { 0 ,  1 , . . .,   }, then   satisfies the following differential equation: with the conditions Since   is a polynomial, the method gives the exact solution by Corollary 5 under the condition deg(  ) ≤ .Thus,   is obtained by finding Bernstein series solution of (42) so that an upper bound of the error is obtained depending on  (3) .
The following procedure, residual correction (e.g., see, [30][31][32]), can be given for the estimation of the absolute error.Moreover, one can estimate the optimal  giving minimal absolute error using this procedure.The procedure is basic.
First, adding and subtracting the term to (2) yields the following differential equation, which admits   :=  −   as an exact solution: with the conditions Let  *  be Bernstein series solution to (45).If ‖  −  *  ‖ ≤  is sufficiently small, the absolute error can be estimated by  *  .Hence, the optimal  for the absolute error can be obtained measuring the error functions

Numerical Examples
In this section, some numerical examples are given to illustrate the method.Some examples are given with their error estimation by using Theorem 4.Moreover, for these examples, the ∞-norms of the error function   , the estimate error function  *  , and the absolute error of the corrected Bernstein series solution   +  *  given in Corollary 7 are calculated for some  and .The optimal truncation limit  is specified for each example.All calculations are done in Maple 15.Since  = 0 is a singular point, the equidistant nodes are selected as {( + 1)/( + 1) :  = 0 ⋅ ⋅ ⋅ }.

Example 8. Consider the Lane-Emden equation
Applying the method for  = 4 on the equidistant nodes, Bernstein series solution is obtained as which is the exact solution [14].
Example 9. Let us consider the equation with the boundary conditions (1) = 0 and   (0) = 0.The exact solution of ( 49) is [5] For different values , the norms and the upper bounds of the absolute errors are obtained on the equidistant nodes by using Theorem 4. Also, estimations of the absolute errors for  = 12 and the norms of the absolute errors for corrected Bernstein series solutions,   +  * 12 , are calculated on the Chebyshev nodes.All results are given in Table 1.The absolute error function for  = 10 and the estimation of the error function,  *   12 , are plotted in Figure 1.As seen from Table 1, the optimal truncation limit  is specified as  = 16, which gives us the best approximation from the set { 3 ,  4 , . . .,  18 }.Moreover, the expected upper bounds are consistent with the absolute errors.Adding  * 12 to   yields the better results in ∞-norm for 3 ≤  ≤ 12.  having () =   2 as exact solution [14,18,33].Assuming that the exact solution  is unknown and  ∈  () [0, ], an upper bound depending on  () is obtained by Theorem 6.The errors for  = 9 and  = 12 are given in Tables 2 and  3, respectively.To obtain   and   , the equidistant nodes and the Chebyshev collocation nodes are used, respectively.Here,  8 and  11 are the Hermite interpolation polynomials on the sets {0,  4 , 1} and {0,  4 ,  8 , 1}, respectively. 9 and  12 represent the values of  (9) and  (12) in ∞-norms, respectively.By calculating ‖ (12) ‖ ∞ and using Theorem 6, the upper bounds of the absolute errors on the equidistant nodes are given in Table 4 by comparison with the absolute error.As seen from the table, these upper bounds bound the absolute error on some reference points.In Table 5, a comparison between Bernstein series solutions for  = 10, 20 and the approximate solution obtained by the Hermite functions collocation (HFC) method [18] for  = 30,  = 6, and  = 2 is given.The results are as follows.
with the boundary conditions (1) = 5.5 and   (0) = 0.The exact solution of ( 49) is [14,33]  For different values  and  = 15, the norms of the absolute errors, the estimations of the absolute errors, and the corrected absolute errors are obtained on the equidistant nodes and given in Table 6.As seen from Table 6, for  ≤ 16, corrected absolute errors are better than the absolute errors.Moreover, residual correction procedure estimates the absolute errors accurately.
which has exact solution (sin )/.To show the effect of working with high accurate computations, Bernstein series solutions are obtained for digits 20 and digits 40.The results are given in Tables 7 and 8 for digits 20 and digits 40, respectively.Table 9 shows the comparison of the Bernstein series solution and the approximate solution given by Pandey et al. [19].
Clearly, norms of the absolute errors decrease to  = 12, and then they increase after that point.These results can be achieved by increasing digits number as in Table 8.Hence, working with high accuracy may yield more accurate results.

Conclusions
To solve Lane-Emden type equations numerically, we introduce a matrix method depending on Bernstein polynomials and collocation points.The method is given with their error analysis.By using Lagrange and Hermite interpolation polynomials, some upper bounds obtained in Section 4 whenever the exact solution is sufficiently smooth.Also the residual correction procedure is given to estimate the absolute error.Even if the exact solution is unknown, one can find an upper bound for the absolute error as in Example 10.Numerical results are consistent with the theoretical results.As in Example 11, increasing number of digits may decrease the round-off error; therefore, more accurate results can be obtained.On the other hand, for  ≤ , corrected Bernstein series solution,   +  *  , is a better approximation than   in ∞-norm in the tables.As a disadvantage of the method, even if Bernstein series solution for  ≫ 20 can be obtained, the results may not be reliable since cond( W) increases.
As a future work, we will shortly extend our study to nonlinear Lane-Emden type differential equation.The error analysis of the method can be improved.The conditions that guarantee the convergence of the method will be explored.

Corollary 7 .
If   is Bernstein series solution to (2), then   +  *  is also an approximate solution for (2).Moreover, its error function is   −  *  .Note that the approximate solution   +  *  is a better approximation than   in the norm for ‖  −  *  ‖ ≤ ‖ −   ‖.Let us call the approximate solution   +  *  as corrected Bernstein series solution.
*  for different  values in any norm.

Table 1 :
The ∞-norms of the absolute errors, estimations of the absolute errors, the ∞-norms of the corrected absolute errors, and upper bounds of the absolute errors for Example 9.

Table 2 :
The values of the absolute error at some points assuming that the exact solution is unknown for Example 10 ( = 9).

Table 3 :
The values of the absolute error at some points assuming that the exact solution is unknown for Example 10 ( = 12).

Table 4 :
Comparison with the absolute errors and their estimated upper bounds obtained by Theorem 6 for Example 10.

Table 6 :
The ∞-norms of the absolute errors, estimations of the absolute errors, and ∞-norms of the corrected absolute errors for Example 11.

Table 7 :
The ∞-norms of the absolute errors, estimations of the absolute errors, and ∞-norms of the corrected absolute errors for Example 12 (digits: 20).

Table 8 :
The ∞-norms of the absolute errors, estimations of the absolute errors, and ∞-norms of the corrected absolute errors for Example 11 (digits: 40).