An effective numerical method to solve a class of nonlinear singular boundary value problems using improved differential transform method

In this work, an effective numerical method is developed to solve a class of singular boundary value problems arising in various physical models by using the improved differential transform method (IDTM). The IDTM applies the Adomian polynomials to handle the differential transforms of the nonlinearities arising in the given differential equation. The relation between the Adomian polynomials of those nonlinear functions and the coefficients of unknown truncated series solution is given by a simple formula, through which one can easily deduce the approximate solution which takes the form of a convergent series. An upper bound for the estimation of approximate error is presented. Several physical problems are discussed as illustrative examples to testify the validity and applicability of the proposed method. Comparisons are made between the present method and the other existing methods.

emerges from the modeling of steady state oxygen diffusion in a spherical cell with Michaelis-Menten uptake kinetics (Lin 1976;McElwain 1978). In this case, u(x) represents the oxygen tension; δ and µ are positive constants involving the reaction rate and the Michaelis constant. Hiltmann and Lory (1983) proposed the existence and uniqueness of the solution for b = 1 and a = c. Analytical bounding functions were given in Anderson and Arthurs (1985). The numerical methods to solve the SBVPs for this case have attracted a reasonable amount of research works, such as the finite difference method (FDM) (Pandey 1997), the cubic spline method (CSM) (Rashidinia et al. 2007;Ravi and Bhattacharya 2006), the Sinc-Galerkin method (SGM) (Babolian et al. 2015), the Adomian decomposition method (ADM) and its modified methods (Khuri and Sayfy 2010;Wazwaz et al. 2013;Singh and Kumar 2014), the variational iteration method (VIM) (Ravi and Aruna 2010; Wazwaz 2011), the series expansion technique (SEM) (Turkyilmazoglu 2013) and the B-spline method (BSM) (Çağlar et al. 2009).
The second case arises in the study of the distribution of heat sources in the human head (Flesch 1975;Gray 1980;Duggan and Goodman 1986), in which α = 2 and In Duggan and Goodman (1986), point-wise bounds and uniqueness results were presented for the SBVPs with the nonlinear function f(x, u) of the forms given by (4) and (5). Quite a little amount of works by using different approaches, including the FDM (Pandey 1997), the CSM (Rashidinia et al. 2007;Ravi and Bhattacharya 2006) and the SGM (Babolian et al. 2015), have been proposed to obtain the approximate solutions of this case.
Besides, Chandrasekhar (1939) derived another case for α = 2, b = 0 and which γ is a physical constant. This case is in connection with the equilibrium of thermal gas thermal (Ames 1968). The numerical solution of this kind of equation for γ = 5 was considered by using various methods, such as the FFDM (Chawla et al. 1988), the VIM (Ravi and Aruna 2010), the SEM (Turkyilmazoglu 2013) and the modified Adomian decomposition method (Singh and Kumar 2014). All the aforementioned methods can yield a satisfied result. However, each of these methods has its own weaknesses. For example, the VIM (Ravi and Aruna 2010; Wazwaz 2011) has an inherent inaccuracy in identifying the Lagrange multiplier, and fails to solve the equation when the nonlinear function f(x, u) is of the forms (5) and (6). Those methods such as the FDM (Pandey 1997;Chawla et al. 1988), the SEM (Turkyilmazoglu 2013), the SGM (Babolian et al. 2015) and the spline method (Rashidinia et al. 2007;Ravi and Bhattacharya 2006;Çağlar et al. 2009) require a tedious process and huge volume of computations in dealing with the linearization or discretization of variables. The ADM (Wazwaz et al. 2013) needs to obtain the corresponding Volterra integral form of the given equation, via which one can overcome the difficulty of singular behavior at x = 0. The modified ADM (Khuri and Sayfy 2010;Kumar and Singh 2010) needs to introduce a twofold indefinite integral operator to give better and accurate results; moreover, the success of method in (Singh and Kumar 2014) relies on constructing Green's function before establishing the recursive relation for applying the ADM to derive the solution components. All those manners are at the expense of computation budgets. Besides, none of above methods is applied to handle the equations with all forms of nonlinearities (4-7).
In recent years, a lot of attentions have been devoted to the applications of differential transform method (DTM) and its modifications. The DTM proposed by Pukhov (1980Pukhov ( , 1982Pukhov ( , 1986 at the beginning of 1980s. However, his work passed unnoticed. , Zhou (1986 reintroduced the DTM to solve the linear and nonlinear equations in electrical circuit problems. The DTM is a semi-numerical-analytic method that generates a Taylor series solution in the different manner. In the past forty years, the DTM has been successfully applied to solve a wide variety of functional equations; see Xie et al. (2016) and the references therein.
Although being powerful, there still exist some difficulties in solving various of equations by the classical DTM. Some researchers have devoted to deal with these obstacles so as to extend the applications of the DTM. For example, in view of the DTM numerical solution cannot exhibit the real behaviors of the problem, Odibat et al. (2010) proposed a multi-step DTM to accelerate the convergence of the series solution over a large region and applied successfully to handle the Lotka-Volterra, Chen and Lorenz systems. In Gökdoğan et al. (2012), Momani and Ertürk (2008) suggested an alternative scheme to overcome the difficulty of capturing the periodic behavior of the solution by combining the DTM, Laplace transform and Padé approximants. Another difficulty is to compute the differential transforms of the nonlinear components in a simple and effective way. By using the traditional approach of the DTM, the computational difficulties will inevitably arise in determining the transformed function of an infinity series. Compared to the traditional method, Chang and Chang (2008) proposed a relatively effective algorithm for calculating the differential transform through a derived recursive relation. Yet, by using their method, it is inevitable to increase the computational budget, especially in dealing with those differential equations which have two or more nonlinear terms being investigated. Recently, the authors Elsaid (2012), Fatoorehchi and Abolghasemi (2013) disclosed the relation between the Adomian polynomials and the differential transform of nonlinearities, and developed an inspiring approach to handle the nonlinear functions in the given functional equation. Meanwhile, the problem of tedious calculations in dealing with nonlinear problems by using the ADM has also been improved considerably by Duan (2010aDuan ( , b, 2011. All of these effective works make it possible to broaden the applicability and popularity of the DTM considerably. The aim of this work is to develop an efficient approach to solve the SBVPs (1-3) with those nonlinear terms (4-7). This scheme is mainly based on the improved differential transform method (IDTM), which is the improved version of the classical DTM by using the Adomian polynomials to handle the differential transforms of those nonlinear functions (4-7). No specific technique is required in dealing with the singular behavior at the origin. Meanwhile, unlike some existing approaches, the proposed method tackles the problem in a straightforward manner without any discretization, linearization or perturbation. The numerical solution obtained by the proposed method takes the form of a convergent series with those easily computable coefficients through the Adomian polynomials of those nonlinear functions as the forms of (4-7).
The rest of the paper is organized as follows. In the next section, the concepts of DTM and Adomian polynomials are introduced. Algorithm for solving the problem (1-3) and an upper bound for the estimation of approximate error are presented in Sect. 3. Sect.4 shows some numerical examples to testify the validity and applicability of the proposed method. In Sect. 5, we end this paper with a brief conclusion.

Adomian polynomial
In the Adomian decomposition method (ADM), a key notion is the Adomian polynomials, which are tailored to the particular nonlinearity to easily and systematically solve nonlinear differential equations. The interested readers are referred to Adomian (1990Adomian ( , 1994 for the details of the ADM.
For the applications of decomposition method, the solution of the given equation in a series form is usually expressed by and the infinite series of polynomials for the nonlinear term f(u), where A m is called the Adomian polynomials, and depends on the solution components u 0 , u 1 , . . . , u m . The traditional algorithm for evaluating the Adomoan polynomials A n was first provided in Adomian and Rach (1983) by the formula A large amount of works (Duan 2010b(Duan , b, 2011Adomian and Rach 1983;Rach 2008Rach , 1984Wazwaz 2000;Abbaoui et al. 1995;Abdelwahid 2003; Azreg-AÏnou 2009) have been applied to give the more effective computational method for the Adomian polynomials. For fast computer generation, we favor Duan's Corollary 3 algorithm (Duan 2011) among all of these methods, as it merely involves the analytic operations of addition and multiplication without the differentiation operator, which is eminently convenient for symbolic implementation by computer algebraic systems such as Maple and Mathematics. The method to generate the Adomian polynomials in Duan (2011) is described as follows: such that It is worth mentioning that Duan's algorithm involving (11) and (12) has been testified to be one of the fastest subroutines on record (Duan 2011), including the fast generation method given by Adomian and Rach (1983).

Differential transform
The differential transform of the kth differentiable function u(x) at x = 0 is defined by and the differential inverse transform of U(k) is described as where u(x) is the original function and U(k) is the transformed function.
For the practical applications, the function u(x) is expressed by a truncated series and Eq. (14) can be written as It is not difficult to deduce the transformed functions of the fundamental operations listed in Table 1.

Method of solution of SBVPs (1-3)
We want to find the approximate solution of the problem (1-3) with the type: where the coefficients U (0), U (1), . . . , U (N ) are determined using the following steps: • According to the definition (13) of the differential transform and the boundary value condition (2), we have Suppose that where β is a real parameter to be determined. • Multiplying both sides of Eq. (1) by variable x, we have Applying the differential transform (13) to Eq. (19), we get the following recurrence relation: where F(k) is the differential transform of the nonlinear function f (x, u) = f (u). • Using Lemma 3.1 in Fatoorehchi and Abolghasemi (2013), we compute F(k) through the Adomian polynomials A k : Remark 1 Lemma 3.1 in Fatoorehchi and Abolghasemi (2013) indicates that the differential transforms and the Adomian polynomials of nonlinear functions have the same mathematical structure such that we can derive the differential transforms of any

Original function
Transformed function nonlinear functions by merely calculating the relevant Adomian polynomials but with constants instead of variable components.
• Substituting (21) into (20), and then combining the relations (16-18), we obtain the truncated series solution of the problem (1-3) as follows: • Imposing the truncated series solution (22) on the boundary condition (3), we obtain a nonlinear algebraic equation with unknown parameter β: Solving Eq. (23), and substituting the value of β into (22), we obtain the final result. An upper bound for the estimation of approximate error is presented in the following lemma.
k! x k is the Taylor polynomial of the unknown function u(x) at

We then have
Combining the relations (25-27), it follows that Thus, the proof is completed.

Numerical examples
In this section, based on the discussion in Sect. 3, we report numerical tests of five classical examples discussed frequently to testify the validity and applicability of the proposed method. All the numerical computations were performed using Maple and Matlab on personal computer. For comparison, we computed the absolute error defined by and the maximal absolute error by where u(x) is the exact solution and u N (x) is the truncated series solution with degree N.
Example 1 Consider the following nonlinear SBVP in the study of isothermal gas sphere (Singh and Kumar 2014;Ravi and Aruna 2010;Chawla et al. 1988): subject to the boundary conditions The exact solution of this problem is given by u(x) = 3 3+x 2 . It is also known as the Emden-Fowler equation of the first kind. In what follows, we shall solve it with the proposed algorithm.
Firstly, we set The Adomian polynomials of nonlinear term f (x, u) = −u 5 (x) in this problem are computed as Furthermore, according to the relations (20) and (21), we obtain the differential transforms U(k) of the unknown function u(x) By using Eq. (22), we obtain the truncated series solution for N = 10 as follows: Secondly, imposing the truncated series solution (33) on the boundary conditions u(1) = √ 3/2, we get a nonlinear algebraic equation. By solving it, the unknown parameter β is computed as Finally, substituting (34) into (33), we get the approximate solution with degree 10 In Table 2, we compare the absolute errors (29) of numerical results obtained by the present method, the VIM (Ravi and Aruna 2010) and the modified ADM using Green functions (GIDM) (Singh and Kumar 2014) for N = 12. Table 3 lists the theoretical estimate errors (24) and the maximal absolute errors (30) of the approximate solutions for changing approximation levels, and shows a comparison of the maximal absolute errors with the GIDM (Singh and Kumar 2014) and the FFDM (Chawla et al. 1988). We can see from Table 3 that the accuracy of our computational results is getting better as the approximation level is increasing. Moreover, our numerical solution u 10 (x) has an accuracy of O(10 −4 ), whereas the GIDM (Singh and Kumar 2014) needs to employ 14 terms to archive this goal as shown in Table 1 of Singh and Kumar (2014); numerical solution with even 64 terms obtained by the FFDM (Chawla et al. 1988) still hovers at this level.
(34) β = 1.000553890.   A comparison of the absolute errors (29) of the numerical solutions for N = 10, 20, 40 obtained by the present method and the modified decomposition method (BSDM) (Khuri and Sayfy 2010) is described in Table 4. Table 5 lists the maximal absolute errors (30) of those numerical results derived from the proposed method, the BSM (Çağlar et al. 2009) and the FFDM (Chawla et al. 1988). And also, we list the theoretical estimate errors (24) in Table 5 for comparison. It can be seen from Tables 4 and 5 that one can obtain the better approximate solution by using the present method compared to the other mentioned methods, even if we take the relative smaller N. Moreover, the theoretical estimate errors, the absolute errors and the maximal absolute errors all decrease as the increase of N. Therefore, evaluation of more components of the numerical solution will reasonably improve the accuracy.
Example 3 Consider the following nonlinear SBVP in the study of steady-state oxygen diffusion in a spherical cell (Babolian et al. 2015;Khuri and Sayfy 2010;Wazwaz 2011;Çağlar et al. 2009): subject to the boundary conditions  where δ and µ are often taken as 0.76129 and 0.03119, respectively. We take the value of α as 1, 2 and 3. The Adomian polynomials of nonlinear term f (x, u) = δu(x) u(x)+µ in this problem are computes as Proceeding as before, we compute the approximate solution u 12,2 (x) for N = 12 and α = 2, and show a comparison of the numerical results compared to the other existing methods in Table 6, from which one can see that the results of our computations are in good agreement with those ones obtained by the SGM (Babolian et al. 2015), the BSDM (Khuri and Sayfy 2010), the VIM (Wazwaz 2011) and the BSM (Çağlar et al. 2009).
Moreover, since there is no exact solution of this problem, we instead investigate the absolute residual error functions and the maximal error remainder parameters, which are the measures of how well the numerical solution satisfies the original problem (37-38). The absolute residual error functions are and the maximal error remainder parameters are In Fig. 1, we plot the absolute residual error functions |ER N ,2 (x)| for N = 2 through 12 by step 2. Besides, the maximal error remainder parameters MER N ,α for the same N and α = 1, 2, 3 are listed in Table 7, from which it is interesting to point out that for a given N the accuracy of our approximate solutions increases with the increase of α. Moreover, Fig. 1 and Table 7 show clearly that the accuracy of our method is getting better as the approximation level is increasing for a fixed α. The logarithm plots of the value of MER 2,α through MER 12,α for α = 1, 2, 3 are displayed in Fig. 2, which demonstrates an approximately exponential rate of convergence for the obtained truncated series solutions and thus the presented method converges rapidly to the exact solution.
Example 4 Consider the following nonlinear SBVP which arises in the study of the distribution of heat sources in the human head (Pandey 1997;Rashidinia et al. 2007;Ravi and Bhattacharya 2006;Babolian et al. 2015;Khuri and Sayfy 2010;Singh and Kumar 2014;Çağlar et al. 2009;Duggan and Goodman 1986): subject to the boundary conditions We consider the following two cases: Case one: a = b = 1. Case two: a = 0.1, b = 1.   The Adomian polynomials of nonlinear term f (x, u) = −e −u(x) in this problem are computed as Again no exact solution exists for this equation, hence it was handled numerically. Table 8 describes the numerical results of the first case obtained by the proposed method at the order of approximation N = 12 and the other existing methods, including the FDM (Pandey 1997), the non-polynomial cubic spline method (NPCSM) (Rashidinia et al. 2007), the CSM (Ravi and Bhattacharya 2006) Table 9. One can seen from two Tables that our computations are in good line with the results obtained by the other approaches compared. In fact, at the approximation level for N = 12, the maximal absolute error is found to be order of magnitude O(10 −7 ) for the first case, and O(10 −9 ) for the second case.
Example 5 Consider the following SBVP with nonlinear term different from the forms (4-7) which arises in the radial stress on a rotationally symmetric shallow membrane cap (Singh and Kumar 2014;Ravi and Aruna 2010): The Adomian polynomials of nonlinear term f (x, u) = 1 2 − 1 8u 2 (x) in this problem are computed as (42) u ′ (0) = 0, u(1) = 1.  Fig. 3, we plot the absolute residual error functions |ER N (x)| for N = 4 through 14 by step 2. The logarithm plot for the maximal error remainder parameters MER N for the same N is shown in Fig. 4, which demonstrates an approximately exponential rate of convergence of the obtained truncated series solutions and thus the presented method converges rapidly to the exact solution. Even though there is no exact solution for this problem, the following 10th order approximation has an accuracy of O(10 −8 ) and can be used for practical applications

Conclusion
In this work, a reliable approach based on the IDTM is presented to handle the numerical solutions of a class of nonlinear SBVPs arising in various physical models. This scheme takes the form of a truncated series with easily computable coefficients via the Adomian polynomials of those nonlinearities in the given problem. With the proposed algorithm, there is no need of discretization of the variables, linearization or small perturbation. Numerical results show that the proposed method works well for the SBVPs (1-3) with a satisfying low error. Besides, it is obvious that evaluation of more components of the approximate solution will reasonably improve the accuracy of truncated series solution by using the proposed method. Comparisons of the results reveal that the present method is very effective and accurate. Moreover, we are convinced that the A 0 = 1 2 − 1 8U 2 (0) , IDTM can be extended to solve the other type of functional equations involving nonlinear terms more easily as the Adomian polynomials are applicable for any analytic nonlinearity and can be generated quickly with the aid of the algorithm proposed by Duan.
It is necessary to point out that algebraic Eq. (23) is a nonlinear one, and we shall inevitably encounter the bad roots while solving it. The criterion to separate the good root from a swarm of bad ones is convergence because it represents the value of unknown function at the origin and will not change for the different N.