Comparison of numerical values of with solution (26).
In this chapter, we present a modified version of the spectral relaxation method for solving singular initial value problems for some Emden-Fowler equations. This study was motivated by the several applications that these equations have in Science. The first step of the method of solution makes use of linearisation to solve the model problem on a small subinterval of the problem domain. This subinterval contains a singularity at the initial instant. The first step is combined with using the spectral relaxation method to recursively solve the model problem on the rest of the problem domain. We make use of examples to demonstrate that the method is reliable, accurate and computationally efficient. The numerical solutions that are obtained in this chapter are in good agreement with other solutions in the literature.
- Emden-Fowler equations
- Lane-Emden equations
- singular initial value problem
- spectral relaxation method
- numerical method
The singular initial value problem
for the Lane-Emden Eq. (1) models several phenomena such as the thermal behaviour of a spherical cloud of gas acting under the mutual attraction of its molecules , the temperature variation of a self gravitating star, the kinetics of combustion , thermal explosion in a rectangular slab  and the density distribution in isothermal gas spheres . Moreover, Eq. (1) has been used many a time as a benchmark for new methods.
A particular case of Eq. (1) is the Emden-Fowler equation of the first kind:
As mentioned in , Eq. (2) represents the dimensionless form of the governing equation for the gravitational potential of a Newtonian self-gravitating, spherically symmetric, polytropic fluid. The equation provides a useful approximation for stars.
A more general form for Eq. (2) is the Emden-Fowler equation
which can be written as
an equation which was discussed in . An existence result for the solution is given therein under certain conditions on and .
Exact solutions are available for particular cases of Eq. (3) , but not for the general case according to the best of our knowledge. This is motivation enough for seeking approximate solutions. To this end several approximate analytical methods were used by other researchers to solve Eq. (3). Van Gorder  made use of the Homotopy analysis method (HAM) and its variant, the Optimal homotopy analysis method, to solve a boundary value problem for the Lane-Emden equation of the second kind. The two respective analytical solutions that they obtained were in strong agreement. The Homotopy perturbation method (HPM) is another variant of the HAM that was used by Chowdhury and Hashim  to solve an initial value problem for Eq. (3). Their analytical solutions were the same as those that were obtained by Wazwaz  using the Adomian decomposition method (ADM). Chowdhury and Hashim observed that the HPM was less computationally expensive than the ADM. Wazwaz  made use of the variational iteration method (VIM) to solve both initial value problems and boundary value problems for Eq. (3) and for some inhomogeneous Emden-Fowler equations. The results that they obtained demonstrated the reliability and effectiveness of the VIM.
Some numerical methods have been used by other researchers to approximate solutions to Eq. (3). Many of these numerical methods fall in the class of
on the problem domain , and the method uses integrals of Haar wavelets as basis functions. The TWM uses roots of shifted Legendre polynomials as collocation points, and as basis functions the method uses Taylor wavelets which are special functions that defined in terms of Taylor polynomials. Convergence results for the Taylor wavelet solution were presented in . The RDF-DQM uses collocation points
on the problem domain and the method makes use of Radial basis functions. Unlike making use of collocation methods for solving Eq. (3) Van Gorder and Vajravelu  used the Runge–Kutta-Felhberg 4-5 (RKF45) method to validate the analytical solutions that they obtained from using the HAM and from using the traditional power series method. The RKF45 method is an embedded Runge–Kutta-pair which makes use of an adaptive stepsize to control the method and to ensure stability properties such as -stability. See  for more details on the RKF45 method.
In this chapter we make use of a modified version of the spectral relation method (SRM) to solve an initial value problem for Eq. (3). We denote our method by MSRM. The SRM was successfully used to solve fluid flow problems by for example Motsa , Motsa
In Section 2 we describe the MSRM for the model problem. In Section 3 we make use of examples to demonstrate the accuracy and computational efficiency of the MSRM. Section 4 concludes this chapter.
2. The MSRM for the model problem
We seek an approximate solution to
where are given constants and and are given functions.
We follow the idea behind the solution method by Ramos , where the singularity at is isolated in a sufficiently small subinterval of where . The point splits interval into two subintervals: and . A linearisation method is used to solve Eq. (5) restricted to , i.e.,
The method by Ramos proceeds by using the same linearisation method on the subintervals of , i.e.
2.1 Near the singularity
Let be a sufficiently small number. Restrict problem (5) to and re-arrange to get
If we Taylor expand about and neglect higher order terms we get
where denotes . Consequently, Eq. (7) can be replaced by
We take and as initial values for problem (7) restricted to in the next section.
2.2 Away from the singularity
We seek satisfying
In this section we begin by describing the SRM for problem (11). In practical applications it is usually important to obtain a solution to (11) which possesses a prescibed degree of accuracy. To this end we make use of the SRM on a partition of the problem domain . This is our last task in this section.
The first step of the SRM for (11) is to let and so that we obtain the equivalent problem
for which upon making use of the change of variable
for where As described in  the next step of the SRM mimicks the Gauss–Seidel method for linear systems and it yields the iteration
for where and on we formed a grid consisting of the Chebyshev-Gauss-Lobatto collocation points
If the initial approximations and to and , respectively, are prescribed then Eqs. (16) and (17) generate sequences and of consecutive approximations. To this end we assume that and are known at the end of the
where is the Chebyshev differentiation matrix,
is a diagonal matrix, and with . Moreover,
and denotes matrix transpose.
We prescribe and by requiring that
Moreover, we assume that
The initial conditions
The larger is the less reliable is the SRM. As a workaround to this problem is we subdivide interval into a disjoint union of non-overlapping subintervals as detailed in Eq. (6). Given the the model problem on , we use the SRM to compute estimate to on . We make use of to compute initial values for the problem on . We repeat this procedure for the problem on and continue in a similar manner until we exhaust . Shown in Algorithm 2.2 is an outline of the MSRM for problem (5).
In this section we make use of some examples to investigate the accuracy and computational effeciency of the MSRM.
Example 1 We look for a numerical solution to the problem
which Wazwaz  solved using the VIM and obtained the approximate.
When we apply the MSRM on (25) we get the following components for constructing the numerical solution.
for problem (8).
2. Initial values
for problem (11).
The MSRM generates a numerical solution to problem (25) which is plotted together with the analytical solution (26) in Figure 1(a). Figure 1(a) shows a good agreement between the numerical and analytical solutions. For a more detailed comparion of the two solutions see Table 1. Table 1 and Figure 1(b) show that the absolute error of the numerical solution is no more than . Thus the numerical solution agrees with the analytical solution in the first decimal places. This degree of accuracy was achieved by the MSRM upon choosing , setting the maximum number of iterations as and imposing the stopping criterion
|MSRM||Solution (26)||Absolute error|
for the iterative method, where is the tolerance and is the infinity norm. It took only iterations for the MSRM to achive the given degree of accuracy.
Example 2 We consider the initial value problem
that was solved by Wazwaz  using the ADM to obtain the approximate analytical solution
Applying the MSRM on (29) produces the following building blocks for constructing the numerical solution.
for problem (8).
2. Initial values
for problem (11).
A comparison of the MSRM solution to problem (29) with the analytical solution (30) is shown in Figure 2(a). The graph suggests that the two solutions are exactly the same. However, a closer look at the two solutions is provided in Table 2 and we observe that the analytical and numerical solutions are slightly different. A plot of the absolute error of the MSRM solution over a grid on the problem domain is shown in Figure 2(b). We observe that the absolute error does not exceed . Hence the MSRM solution and the analytical solution agree to within decimal places. For the MSRM to achieve this degree of accuracy we chose and . We observed that the MSRM stopped at iteration for problem (29).
Example 3 We seek a numerical solution to
The MSRM for (34) gives the following components for constructing the numerical solution.
for problem (8).
2. Initial values
for problem (11).
3. Right hand sides
Figure 3(a) shows that the numerical solution agrees well with the analytical solution (34). For a closer look at how the numerical and analytical solutions compare see Table 3. We observe that the MSRM solution agrees with the analytical solution (34) on the problem domain to within at least decimal places. A plot of the absolute error at these grid points on is shown in Figure 3(b). We observe that the absolute error in the MSRM increases as me move away from the singular point , but the absolute error never exceeds . This degree of accuracy is achieved with and . Moreover, only iterations were required to achieve this degree of accuracy.
In this chapter we presented a modified spectral relaxation method (denoted by MSRM) for solving singular initial value problems for some Emden-Fowler equations. We made use of some examples of the model problem to demonstrate that the MSRM is reliable, accurate and computationally efficient. The method provided a reliable treatment of the singular point. The MSRM solutions were compared with analytical solutions that were obtained using other methods, i.e., the Variational iteration method and the Adomian decomposition method. There was agreement between the solutions that were compared in the first decimal places. A possible way of increasing the degree of accuracy of the MSRM would be to increase the tolerance for the method. This and other ways for optimizing the method could consistute future work. In all the examples that were considered, it took at most iterations for the MSRM to converge. Hence the method exhibited rapid convergence.