# Analysis of IVPs and BVPs on Semi-Infinite Domains via Collocation Methods

Authors: Mohammad Maleki; Ishak Hashim; Saeid Abbasbandy
Journal: Journal of Applied Mathematics
 (2012)

---

## Abstract

We study the numerical solutions to semi-infinite-domain two-point boundary value problems and initial value problems. A smooth, strictly monotonic transformation
is used to map the semi-infinite domainx∈[0,∞) onto a half-open interval t∈[−1,1). The resulting finite-domain two-point boundary value problem is transcribed to a system of algebraic equations using Chebyshev-Gauss (CG) collocation, while the resulting initial
value problem over a finite domain is transcribed to a system of algebraic equations using Chebyshev-Gauss-Radau (CGR) collocation. In numerical experiments, the tuning of the map ϕ:[−1,+1)→[0,+∞) and its effects on the quality of the discrete approximation are analyzed.


---

## Body

## 1. Introduction

Boundary value problems (BVPs) on infinite intervals and initial value problems (IVPs) frequently occur in mathematical modelling of various applied problems. Typically, these problems arise very frequently in fluid dynamics, aerodynamics, quantum mechanics, electronics, astrophysics, and other domains of science, as examples, in the study of unsteady flow of a gas through a semi-infinite medium [1–3], heat transfer in the radial flow between parallel circular disks [4], draining flows [5], circular membranes [6–8], plasma physics [9, 10], radially symmetric solutions of semilinear elliptic equations [10–12], non-Newtonian fluid flows [13], study of stellar structure, thermal behavior of a spherical cloud of gas, isothermal gas sphere, theory of thermionic currents [14–16], and modeling of vortex solitons [17, 18].The analytical solutions to BVPs and IVPs are not readily attainable and thus the need for finding efficient computational algorithms for obtaining numerical solutions arises. During the last few years much progress has been made in the numerical treatment of boundary value problems over infinite intervals (see e.g., [19–23] and references therein). Moreover, different spectral methods have been proposed for solving problems on unbounded domains. One of the methods is through the use of polynomials that are orthogonal over unbounded domains, such as the Hermite spectral and the Laguerre spectral methods [24–29]. However, all of these algorithms need certain quadratures on unbounded domains, which introduce errors and so weaken the merit of spectral approximations.Another approach is replacing the infinite domain with[-L,L] and the semi-infinite interval with [0,L] by choosing L sufficiently large. This method is named as the domain truncation [30]. When we use finite difference methods or finite element methods to solve unbounded problems numerically, we often utilize domain truncation and impose certain conditions on artificial boundaries, which destroy the accuracy usually [22].Another direct approach for solving such problems is based on rational approximations. Christov [31] and Boyd [32, 33] developed some spectral methods on unbounded intervals by using mutually orthogonal systems of rational functions. Boyd [33] defined a new spectral basis, named rational Chebyshev functions on the semi-infinite interval, by mapping it to the Chebyshev polynomials. Guo et al. [34] introduced a new set of rational Legendre functions that are mutually orthogonal in L2(0,+∞). They applied a spectral scheme using the rational Legendre functions for solving the Korteweg-de Vries equation on the half-line. Boyd et al. [35] applied pseudospectral methods on a semi-infinite interval and compared the rational Chebyshev, the Laguerre and the mapped Fourier sine methods. The authors of [36–40] solved some IVPs and BVPs on semi-infinite intervals.The use of a suitable mapping to transfer infinite domains to the finite domains and then applying the standard spectral methods for the transformed problems in finite domains is another approach that is frequently used [41–45] and is used in this paper too.In this paper we describe two collocation methods for the numerical solution of semi-infinite-domain two-point boundary value problems and initial value problems. For mapping the semi-infinite domainx∈[0,∞) onto a finite half-open interval t∈[-1,+1), the algebraic and logarithmic transformations, which are smooth and strictly monotonic, are used. The Chebyshev-Gauss (CG), Chebyshev-Gauss-Radau (CGR), and Chebyshev-Gauss-Lobatto (CGL) points are then introduced. Then, the resulting finite-domain two-point boundary value problem is discretized using CG collocation, and the resulting initial value problem over a finite domain is discretized using CGR collocation in order to handle the initial slope condition based on the collocation points. Our collocation schemes avoid the singularity at t=+1 introduced by the change of variables because the final CG or CGR point is strictly smaller than t=+1. Furthermore, in the CG scheme for two-point boundary value problems, the function at x=+∞ is obtained directly from the final condition, while, in the CGR scheme for initial value problems, the function at x=+∞ is estimated after solving the problem. In addition, we find that if the exact solution to the original problem decays exponentially fast as x→∞, then better approximations to the problem can be attained by using a transformation x=ϕ(t) that grows more slowly near the point t=+1, while if the exact solution decays algebraically as x→∞, then better approximations to the problem can be attained by using a transformation that grows more rapidly near the point t=+1. Thus, by tuning the choice of the transformation, the accuracy of the discretization can be improved.This paper is organized as follows: in Section2, using algebraic and logarithmic transformations we transform the semi-infinite-domain BVPs and IVPs to finite-domain problems. In Section 3, we discuss the choices for collocation points and explain our method. Section 4 is devoted to some numerical examples and analysis of tuning the transformations. Finally, in Section 5 we provide conclusions.


## 2. Problem Statement and Domain Transformation

Consider the semi-infinite-domain two-point boundary value problem as(2.1)f(x,y(x),y′(x),y′′(x))=0,0⩽x<∞,(2.2)y(0)=c0,y(∞)=c1,
where f is a continuous function and c0 and c1 are known constants. Also, consider the initial value problem over a semi-infinite domain as defined by (2.1) with the conditions(2.3)y(0)=c0,y′(0)=c0′,
where c0 and c0′ are known constants and in this case we assume that limx→∞y(x) exists.We make the change of variablesx=ϕ(t), where ϕ is a differentiable, strictly monotonic function of t that maps the interval [-1,1) onto [0,∞). There are several typical mappings that relate infinite and finite domains to each other, but three specific types of mappings, as algebraic, logarithmic, and exponential, are more practical. Typical examples of such functions are the algebraic transformation(2.4)x=ϕa(t)=L1+t(1-t)K
and the logarithmic transformation(2.5)x=ϕl(t)=Lln⁡⁡(21-t),
where L is a positive parameter that is called mapping parameter, and in this paper we set K=1/2, 1, and 3/2. Note that other positive values for K are possible, but these values seem to be proper for a wide range of problems. Algebraic transformation with K=1 and exponential transformation were proposed in [41, 42], while other algebraic transformations are introduced in this paper.Boyd in [30, 42] offered guidelines for optimizing the mapping parameter L where L>0. Boyd used these guidelines to obtain approximations to known functions and not to find solutions to differential equations. In general, there is no way to avoid a small amount of trial and error in choosing L when solving problems on an infinite domain. Note that the optimum value of L varies with the number of collocation points, and, as Boyd [30] stated, for smooth functions the error is relatively insensitive to the precise choice of L in the neighborhood of the optimum value of L.As we will see in the numerical experiments, better discretization can be achieved by tuning the transformations to the problem. For a fixedL, the logarithmic transformation produces slower growth in x=ϕ(t) as t approaches +1 than that of algebraic transformations introduced here. Moreover, in algebraic transformations by increasing K faster growth in x=ϕ(t) is produced. However, the logarithmic transformation and the algebraic transformation with K=1/2 have almost the same behavior. In particular, we will see that if the exact solution to the original problem decays exponentially fast as x→∞, then both transformations (algebraic and logarithmic) work, but logarithmic transformation and algebraic transformation with K=1/2 give better results in comparison with other choices for algebraic transformations. On the other hand, if the exact solution decays algebraically as x→∞, then algebraic transformations with K=1 and K=3/2 give much better results. Furthermore, the smoothness of the exact solution in the transformed coordinate plays an important role in tuning the transformations to the problem.The transformations in (2.4)-(2.5) are especially convenient because they yield simple expressions for derivatives. With a general transformation x=ϕ(t) with ϕ′(t)>0, -1⩽t⩽1, the first and second derivatives of y(x) can be expressed in terms of t as follows:(2.6)y′(x)=1ϕ′(t)Y′(t),y′′(x)=1ϕ′2(t)Y′′(t)-ϕ′′(t)ϕ′3(t)Y′(t).
After changing variables from x to t, the nonlinear differential equation in (2.1) becomes(2.7)f(ϕ(t),Y(t),1ϕ′(t)Y′(t),1ϕ′2(t)Y′′(t)-ϕ′′(t)ϕ′3(t)Y′(t))=0.
Moreover, the boundary conditions in (2.2) become(2.8)Y(-1)=c0,Y(1)=c1,
while the initial conditions in (2.3) become(2.9)Y(-1)=c0,1ϕ′(-1)Y′(-1)=c0′.


## 3. Collocation Methods for Semi-Infinite-Domain BVPs and IVPs

In this section we formulate discrete approximations to the problems described in Section2. These schemes are based on global collocation using CG and CGR collocation points. As will be seen, the schemes for semi-infinite-domain two-point boundary value problems and for initial value problems are different in respect of their treatment of the origin and the infinity. For semi-infinite-domain two-point boundary value problems, the function value at the infinity is obtained directly from the boundary condition at the infinity, while for initial value problems the derivative at the origin is expanded based on collocation points to handle the initial slope condition, and the function value at the infinity is estimated after obtaining the approximate solution.


### 3.1. CG, CGR, and CGL Collocation Points

The CG, CGR, and CGL collocation points [46] lie on the open interval (-1,1), the half-open interval [-1,1) or (-1,1], and the closed interval [-1,1], respectively. The CG points contain neither −1 nor 1, the CGR points contain only one of the points −1 or 1, and the CGL points contain both −1 and 1. The CG points are the roots of TN(t) (Nth-degree Chebyshev polynomial) as(3.1)ti=-cos⁡⁡(2i-12Nπ),i=1,2,…,N,
the CGR points are defined by(3.2)ti=±cos⁡⁡(2i2N+1π),i=0,1,…,N,
and the CGL points are the roots of (1-t2)TN′(t) as(3.3)ti=-cos⁡⁡(iNπ),i=0,1,…,N.
It is known that the CG and CGL points are symmetric about the origin whereas the CGR points are asymmetric. In addition, the CGR points are not unique in that two sets of points exist (one including the point −1 and the other including the point +1). The CGR points that include the initial point are called the standard CGR points and the CGR points that include the terminal end point are often called the flipped CGR points.


### 3.2. Gauss Collocation Scheme for Infinite-Domain Two-Point Boundary Value Problems

Consider a set of CG collocation points(t1,t2,…,tN) defined by (3.1) on the interval (-1,1) and two additional noncollocated points t0=-1 and tN+1=1 (the terminal value, corresponding to x=+∞). The function Y∈L2[-1,1] is approximated by a polynomial of degree at most N+1 as(3.4)Y(t)≈IN+1(Y)(t)=∑j=0N+1YjLj(t),
where Yj=Y(tj) and(3.5)Lj(t)=∏k=0,k≠jN+1t-tktj-tk,j=0,1,…,N+1,
is a basis of (N+1)th-degree Lagrange polynomials satisfying Lj(ti)=δij. Notice that the basis includes the functions L0 and LN+1 corresponding to the initial value t0=-1 (x=0) and the terminal value tN+1=1 (x=+∞).The derivative matrixD(1)=(dij(1)) at CG collocation points ti, i=1,2,…,N is defined by(3.6)dij(1)=Lj′(ti),j=0,1,…,N+1,
and the second derivative D(2)=(dij(2)) is defined by(3.7)dij(2)=Lj′′(ti),j=0,1,…,N+1.
For computing the entries of D(1) and D(2), we follow the method proposed in [47]. According to [47] the entries of D(1) are computed by(3.8)dij(1)={cicj(xi-xj)-1,i≠j,-∑l=0,l≠iN+1dil(1),i=j,
where i=1,2,…,N, j=0,1,…,N+1, and ck=∏l=0,l≠kN+1(xk-xl). For large value of N, the computation of ratio ci/cj in dij(1) may cause round-off error, and therefore to avoid this problem as mentioned in [47] we can compute them as(3.9)bk=∑l=0,l≠kN+1ln⁡⁡(|xk-xl|),cicj=(-1)i+jebi-bj.
The entries of matrix D(2) can be computed recursively by the entries of D(1) as follows:(3.10)dij(2)={2(dii(1)dij(1)-(xi-xj)-1dij(1)),i≠j,-∑l=0,l≠iN+1dil(2),i=j.
Similarly, higher-order derivative matrices can be computed recursively by lower-order ones. For more details, see [47].Therefore, the approximations of the first and second derivatives of the series of (3.4) at CG collocation points ti, i=1,2,…,N, are obtained as follows:(3.11)Y′(ti)≈Di(1)⋅Ỹ,(3.12)Y′′(ti)≈Di(2)⋅Ỹ,
where Ỹ=[Y0,Y1,…,YN+1]T and Di(1) and Di(2) are the ith row of D(1) and D(2), respectively. The rectangular N×(N+2) matrices D(1) and D(2) formed by the coefficients dij(1) and dij(2), (i=1,…,N, j=0,…,N+1) are the first- and second-order Gauss differentiation matrices since they transform the approximation of Y(t) at t0,…,tN+1 to the first and second derivatives of Y(t) at the gauss collocation points t1,…,tN.Our discrete approximation to the semi-infinite-domain two-point boundary value problem in (2.1)-(2.2) is obtained as follows. After transformation of the domain, approximating the function Y using the series of (3.4) and evaluating the boundary conditions in (2.8) at noncollocated (but interpolation) points t0=-1 and tN+1=1 imply(3.13)Y(t0)≈Y0=c0,Y(tN+1)≈YN+1=c1.
As a result, N unknowns Y1,…,YN remain to be determined. By evaluating (2.7) at each collocation point ti, i=1,…,N, and replacing Y′(ti) and Y′′(ti) by their discrete approximations in (3.11)-(3.12) we obtain(3.14)f(ϕ(ti),Yi,1ϕ′(ti)Di(1)⋅Ỹ,1ϕ′2(ti)Di(2)⋅Ỹ-ϕ′′(ti)ϕ′3(ti)Di(1)⋅Ỹ)=0.
By using (3.14) we generate a set of N nonlinear algebraic equations for the unknowns Y1,…,YN that can be solved using Newton’s iterative method. It is well known that the initial guess for Newton’s iterative method is important especially for complicated problems. To choose the initial guess for our problem, in the first stage we set N=5 and apply Newton’s iterative method for solving N nonlinear equations by choosing the values of the linear function between the boundary conditions in (2.8) as our initial guesses. We then increase N by 5 or 10 and use the approximate solution at stage one as our initial guess for the next stage. We continue this approach until the results are similar up to a required number of decimal places for two consecutive stages.Note that, although the transformed two-point boundary value problem in (2.7) may have a singularity at t=-1 and the change of variable x=ϕ(t) must have a singularity at t=+1, we never evaluate the problem at the singularities, rather we evaluate it at the collocation points that are all strictly more than −1 and less than 1.With substituting the obtained values ofYj, j=0,1,…,N+1, in (3.4) the approximation of Y(t) is obtained. Finally, the approximate solutions on [0,∞) using the algebraic transformations in (2.4) with K=1/2, 1, and 3/2 are found as(3.15)y(x)≈∑j=0N+1YjLj(-2L2-x2+x8L2+x22L2),y(x)≈∑j=0N+1YjLj(x-Lx+L),y(x)≈∑j=0N+1YjLj(γ(-L2+3x2+γ)+L2-12L2x23x2γ),
respectively, where(3.16)γ=(-L6+18L4x2-54L2x4+6-6L6x6+81L4x8)1/3,
and the approximate solution on [0,∞) using the logarithmic transformation in (2.5) is found as(3.17)y(x)≈∑j=0N+1YjLj(1-2e-x/L).


### 3.3. Radau Collocation Scheme for Initial Value Problems

Consider thestandard CGR collocation points -1=t0<⋯<tN<+1 defined by (3.2) and the additional noncollocated point tN+1=1. As the CG scheme, the function Y is approximated using (3.4). It is worth mentioning that for the Radau scheme, the point t0 is a collocation point, while, for the Gauss scheme in the previous subsection, the point t0 is just an interpolation node. For IVPs, the CGR points are more appropriate because they allow approximation of the initial slope condition to be obtained based on collocation points.Again, the derivative matrixD(1) at CGR collocation points ti, i=0,1,…,N, is defined by(3.18)dij(1)=Lj′(ti),j=0,1,…,N+1,
and the second derivative D(2) is defined by(3.19)dij(2)=Lj′′(ti),j=0,1,…,N+1.
The elements of D(1) and D(2) are computed similarly to the previous subsection with the difference that here i=0,1,…,N. In addition, the approximations of the first and second derivatives of the series of (3.4) at CGR collocation points ti, i=0,1,…,N, are obtained using (3.11)-(3.12). Notice that, in the CGR scheme, unlike the CG scheme, the approximation of Y′(t0) is obtained based on collocation points. The rectangular (N+1)×(N+2) matrices D(1) and D(2) formed by the coefficients dij(1) and dij(2), (i=0,…,N, j=0,…,N+1) are the first- and second-order Radau differentiation matrices, respectively. These matrices transform the approximation of Y(t) at t0,…,tN+1 to the first and second derivatives of Y(t) at the Radau collocation points t0,…,tN.The discrete approximation to the initial value problem in (2.1) and (2.3) is obtained as follows. After transformation of the domain, by approximating the function Y using the series of (3.4) and evaluating the boundary conditions in (2.9) at the collocation (also interpolation) point t0=-1 and using (3.11) we obtain(3.20)Y(t0)≈Y0=c0,(3.21)1ϕ′(t0)Y′(t0)≈1ϕ′(t0)D0(1)⋅Ỹ=c0′.
Consequently, N+1 unknowns Y1,Y2,…,YN+1 remain to be determined. Unlike the CG scheme for infinite-domain BVPs, in the CGR scheme for IVPs the function value at the infinity, YN+1, appears in the unknowns. As the CG scheme, we evaluate (2.7) at collocation points ti for i=1,…,N and replace Y′(ti) and Y′′(ti) by their discrete approximations in (3.11)-(3.12) to generate N algebraic equations, by which in combination with (3.21) we obtain a set of N+1 algebraic equations for the unknowns Y1,Y2,…,YN+1. In this case, to choose the initial guess for Newton’s iterative method, in the first stage we set N=5 and apply Newton’s iterative method for solving N+1 nonlinear equations by choosing c0 in (2.9) as our initial guess. We then proceed the procedure described for the CG scheme.Note that the CGR scheme of this subsection can also be applied for infinite-domain two-point boundary values problems in Section3.2. Nevertheless, the errors using either the CG or CGR schemes are almost the same.


### 3.4. Assessment of Spectral Accuracy

The convergence properties of the approximation toy(x) can be determined from the behavior of the function Y(t)=y(ϕ(t)). Infinite-order accuracy is expected when Y(t) is infinitely smooth on [-1,1]. Assuming that y(x) itself is infinitely smooth on [0,+∞), the critical issue is the behavior of the derivatives of Y(t) at t=1. Moreover, unlike expansions on a finite domain, spectral approximations on a semi-infinite domain have two discretization parameters, the mapping parameter L in addition to the number of collocation points. Furthermore, the decay rate of the solution y(x) plays an important role in investigating the convergence properties. More precisely, if y(x) decays exponentially fast as x→∞, then better approximations to the problem can be attained by using a transformation x=ϕ(t) that grows more slowly near the point t=+1, while if y(x) decays algebraically as x→∞, then better approximations can be attained by using a transformation that grows more rapidly near the point t=+1. In the next section, we analyze the effects of the mentioned parameters via test examples.


## 3.1. CG, CGR, and CGL Collocation Points

The CG, CGR, and CGL collocation points [46] lie on the open interval (-1,1), the half-open interval [-1,1) or (-1,1], and the closed interval [-1,1], respectively. The CG points contain neither −1 nor 1, the CGR points contain only one of the points −1 or 1, and the CGL points contain both −1 and 1. The CG points are the roots of TN(t) (Nth-degree Chebyshev polynomial) as(3.1)ti=-cos⁡⁡(2i-12Nπ),i=1,2,…,N,
the CGR points are defined by(3.2)ti=±cos⁡⁡(2i2N+1π),i=0,1,…,N,
and the CGL points are the roots of (1-t2)TN′(t) as(3.3)ti=-cos⁡⁡(iNπ),i=0,1,…,N.
It is known that the CG and CGL points are symmetric about the origin whereas the CGR points are asymmetric. In addition, the CGR points are not unique in that two sets of points exist (one including the point −1 and the other including the point +1). The CGR points that include the initial point are called the standard CGR points and the CGR points that include the terminal end point are often called the flipped CGR points.


## 3.2. Gauss Collocation Scheme for Infinite-Domain Two-Point Boundary Value Problems

Consider a set of CG collocation points(t1,t2,…,tN) defined by (3.1) on the interval (-1,1) and two additional noncollocated points t0=-1 and tN+1=1 (the terminal value, corresponding to x=+∞). The function Y∈L2[-1,1] is approximated by a polynomial of degree at most N+1 as(3.4)Y(t)≈IN+1(Y)(t)=∑j=0N+1YjLj(t),
where Yj=Y(tj) and(3.5)Lj(t)=∏k=0,k≠jN+1t-tktj-tk,j=0,1,…,N+1,
is a basis of (N+1)th-degree Lagrange polynomials satisfying Lj(ti)=δij. Notice that the basis includes the functions L0 and LN+1 corresponding to the initial value t0=-1 (x=0) and the terminal value tN+1=1 (x=+∞).The derivative matrixD(1)=(dij(1)) at CG collocation points ti, i=1,2,…,N is defined by(3.6)dij(1)=Lj′(ti),j=0,1,…,N+1,
and the second derivative D(2)=(dij(2)) is defined by(3.7)dij(2)=Lj′′(ti),j=0,1,…,N+1.
For computing the entries of D(1) and D(2), we follow the method proposed in [47]. According to [47] the entries of D(1) are computed by(3.8)dij(1)={cicj(xi-xj)-1,i≠j,-∑l=0,l≠iN+1dil(1),i=j,
where i=1,2,…,N, j=0,1,…,N+1, and ck=∏l=0,l≠kN+1(xk-xl). For large value of N, the computation of ratio ci/cj in dij(1) may cause round-off error, and therefore to avoid this problem as mentioned in [47] we can compute them as(3.9)bk=∑l=0,l≠kN+1ln⁡⁡(|xk-xl|),cicj=(-1)i+jebi-bj.
The entries of matrix D(2) can be computed recursively by the entries of D(1) as follows:(3.10)dij(2)={2(dii(1)dij(1)-(xi-xj)-1dij(1)),i≠j,-∑l=0,l≠iN+1dil(2),i=j.
Similarly, higher-order derivative matrices can be computed recursively by lower-order ones. For more details, see [47].Therefore, the approximations of the first and second derivatives of the series of (3.4) at CG collocation points ti, i=1,2,…,N, are obtained as follows:(3.11)Y′(ti)≈Di(1)⋅Ỹ,(3.12)Y′′(ti)≈Di(2)⋅Ỹ,
where Ỹ=[Y0,Y1,…,YN+1]T and Di(1) and Di(2) are the ith row of D(1) and D(2), respectively. The rectangular N×(N+2) matrices D(1) and D(2) formed by the coefficients dij(1) and dij(2), (i=1,…,N, j=0,…,N+1) are the first- and second-order Gauss differentiation matrices since they transform the approximation of Y(t) at t0,…,tN+1 to the first and second derivatives of Y(t) at the gauss collocation points t1,…,tN.Our discrete approximation to the semi-infinite-domain two-point boundary value problem in (2.1)-(2.2) is obtained as follows. After transformation of the domain, approximating the function Y using the series of (3.4) and evaluating the boundary conditions in (2.8) at noncollocated (but interpolation) points t0=-1 and tN+1=1 imply(3.13)Y(t0)≈Y0=c0,Y(tN+1)≈YN+1=c1.
As a result, N unknowns Y1,…,YN remain to be determined. By evaluating (2.7) at each collocation point ti, i=1,…,N, and replacing Y′(ti) and Y′′(ti) by their discrete approximations in (3.11)-(3.12) we obtain(3.14)f(ϕ(ti),Yi,1ϕ′(ti)Di(1)⋅Ỹ,1ϕ′2(ti)Di(2)⋅Ỹ-ϕ′′(ti)ϕ′3(ti)Di(1)⋅Ỹ)=0.
By using (3.14) we generate a set of N nonlinear algebraic equations for the unknowns Y1,…,YN that can be solved using Newton’s iterative method. It is well known that the initial guess for Newton’s iterative method is important especially for complicated problems. To choose the initial guess for our problem, in the first stage we set N=5 and apply Newton’s iterative method for solving N nonlinear equations by choosing the values of the linear function between the boundary conditions in (2.8) as our initial guesses. We then increase N by 5 or 10 and use the approximate solution at stage one as our initial guess for the next stage. We continue this approach until the results are similar up to a required number of decimal places for two consecutive stages.Note that, although the transformed two-point boundary value problem in (2.7) may have a singularity at t=-1 and the change of variable x=ϕ(t) must have a singularity at t=+1, we never evaluate the problem at the singularities, rather we evaluate it at the collocation points that are all strictly more than −1 and less than 1.With substituting the obtained values ofYj, j=0,1,…,N+1, in (3.4) the approximation of Y(t) is obtained. Finally, the approximate solutions on [0,∞) using the algebraic transformations in (2.4) with K=1/2, 1, and 3/2 are found as(3.15)y(x)≈∑j=0N+1YjLj(-2L2-x2+x8L2+x22L2),y(x)≈∑j=0N+1YjLj(x-Lx+L),y(x)≈∑j=0N+1YjLj(γ(-L2+3x2+γ)+L2-12L2x23x2γ),
respectively, where(3.16)γ=(-L6+18L4x2-54L2x4+6-6L6x6+81L4x8)1/3,
and the approximate solution on [0,∞) using the logarithmic transformation in (2.5) is found as(3.17)y(x)≈∑j=0N+1YjLj(1-2e-x/L).


## 3.3. Radau Collocation Scheme for Initial Value Problems

Consider thestandard CGR collocation points -1=t0<⋯<tN<+1 defined by (3.2) and the additional noncollocated point tN+1=1. As the CG scheme, the function Y is approximated using (3.4). It is worth mentioning that for the Radau scheme, the point t0 is a collocation point, while, for the Gauss scheme in the previous subsection, the point t0 is just an interpolation node. For IVPs, the CGR points are more appropriate because they allow approximation of the initial slope condition to be obtained based on collocation points.Again, the derivative matrixD(1) at CGR collocation points ti, i=0,1,…,N, is defined by(3.18)dij(1)=Lj′(ti),j=0,1,…,N+1,
and the second derivative D(2) is defined by(3.19)dij(2)=Lj′′(ti),j=0,1,…,N+1.
The elements of D(1) and D(2) are computed similarly to the previous subsection with the difference that here i=0,1,…,N. In addition, the approximations of the first and second derivatives of the series of (3.4) at CGR collocation points ti, i=0,1,…,N, are obtained using (3.11)-(3.12). Notice that, in the CGR scheme, unlike the CG scheme, the approximation of Y′(t0) is obtained based on collocation points. The rectangular (N+1)×(N+2) matrices D(1) and D(2) formed by the coefficients dij(1) and dij(2), (i=0,…,N, j=0,…,N+1) are the first- and second-order Radau differentiation matrices, respectively. These matrices transform the approximation of Y(t) at t0,…,tN+1 to the first and second derivatives of Y(t) at the Radau collocation points t0,…,tN.The discrete approximation to the initial value problem in (2.1) and (2.3) is obtained as follows. After transformation of the domain, by approximating the function Y using the series of (3.4) and evaluating the boundary conditions in (2.9) at the collocation (also interpolation) point t0=-1 and using (3.11) we obtain(3.20)Y(t0)≈Y0=c0,(3.21)1ϕ′(t0)Y′(t0)≈1ϕ′(t0)D0(1)⋅Ỹ=c0′.
Consequently, N+1 unknowns Y1,Y2,…,YN+1 remain to be determined. Unlike the CG scheme for infinite-domain BVPs, in the CGR scheme for IVPs the function value at the infinity, YN+1, appears in the unknowns. As the CG scheme, we evaluate (2.7) at collocation points ti for i=1,…,N and replace Y′(ti) and Y′′(ti) by their discrete approximations in (3.11)-(3.12) to generate N algebraic equations, by which in combination with (3.21) we obtain a set of N+1 algebraic equations for the unknowns Y1,Y2,…,YN+1. In this case, to choose the initial guess for Newton’s iterative method, in the first stage we set N=5 and apply Newton’s iterative method for solving N+1 nonlinear equations by choosing c0 in (2.9) as our initial guess. We then proceed the procedure described for the CG scheme.Note that the CGR scheme of this subsection can also be applied for infinite-domain two-point boundary values problems in Section3.2. Nevertheless, the errors using either the CG or CGR schemes are almost the same.


## 3.4. Assessment of Spectral Accuracy

The convergence properties of the approximation toy(x) can be determined from the behavior of the function Y(t)=y(ϕ(t)). Infinite-order accuracy is expected when Y(t) is infinitely smooth on [-1,1]. Assuming that y(x) itself is infinitely smooth on [0,+∞), the critical issue is the behavior of the derivatives of Y(t) at t=1. Moreover, unlike expansions on a finite domain, spectral approximations on a semi-infinite domain have two discretization parameters, the mapping parameter L in addition to the number of collocation points. Furthermore, the decay rate of the solution y(x) plays an important role in investigating the convergence properties. More precisely, if y(x) decays exponentially fast as x→∞, then better approximations to the problem can be attained by using a transformation x=ϕ(t) that grows more slowly near the point t=+1, while if y(x) decays algebraically as x→∞, then better approximations can be attained by using a transformation that grows more rapidly near the point t=+1. In the next section, we analyze the effects of the mentioned parameters via test examples.


## 4. Illustrative Examples

We tested the methods presented in this paper on several problems such as the problem of unsteady flow of gas through a semi-infinite porous medium and the Flierl-Petviashvili (FP) equation. For the FP equation, a modified version of the Radau scheme of this paper is introduced in order to approximate the FP monopole. In this section we also analyze tuning the transformations to the problems. All solutions were obtained using the well-known symbolic software “Mathematica 8."


### 4.1. Semi-Infinite-Domain Two-Point Boundary Value Problems

Example 4.1.
We begin with a simple example to illustrate and analyze tuning the transformations. Consider the linear semi-infinite-domain two-point boundary value problem [22]
(4.1)-y′′(x)-2y′(x)+2y(x)=e-2x,y(0)=1,y(∞)=0.
The exact analytical solution is ye(x)=(1/2)(e-(1+3)x+e-2x). The problem in (4.1) was solved for N=(10,20,30) using the collocation method described in Section 3.2 with the CG collocation points described in Section 3.1 and with the strictly monotonic transformations of the domain t∈[-1,1) given in (2.4)-(2.5).
The maximum absolute error is defined by(4.2)EN=max⁡{|Ye(t)-IN+1(Y)(t)|:-1⩽t⩽1},
where IN+1(Y)(t) and Ye(t) are the approximated and exact solutions in the transformed coordinate, respectively. The maximum absolute errors obtained using the CG scheme of this paper, are given in Table 1 with the logarithmic transformation and the algebraic transformations with K=1/2, 1,  and 3/2 given in (2.4)-(2.5). As can be seen, for all four transformations, the errors decrease until N=30, demonstrating an approximately exponential convergence rate. It is seen that, for each value of N, there exists an optimal value for the mapping parameter L and this optimal value varies with the transformation. In addition, it is observed that changing the values of mapping parameter L improves the quality of the discretization, significantly. For example, for N=20, the maximum absolute error using the logarithmic transformation ϕl decreases from 2.72×10-8 for L=1 to 1.80×10-15 for L=5. Moreover, the logarithmic transformation ϕl allows accurate results to be obtained much more efficiently. For example, to achieve an error of less than 10-14 we require less than 20 collocation points using the logarithmic transformation with L=3, while for other transformations more collocation points are needed (see Table 1). However, the algebraic transformation ϕa with K=1/2 produces very satisfactory results too. Finally, it is seen that for each value of N, the best maximum absolute error for logarithmic transformation is less than that for algebraic transformations.
The different behavior of the transformations given in (2.4)-(2.5) is understood if we apply the change of variables to the continuous solution. The exact solution in the transformed coordinates is as follows:
(4.3)Ya(t)=12(e-(1+3)L((1+t)/(1-t)K)+e-2L((1+t)/(1-t)K)),Yl(t)=12((1-t2)(1+3)L+(1-t2)2L).
Here the subscripts a and l correspond to the algebraic and logarithmic transformations given in (2.4)-(2.5). An advantage of using the logarithmic change of variable and also the algebraic change of variable with K=1/2, as compared to other algebraic change of variables, is that for proper values of L they essentially move collocation points associated with large values of x to the left. Because the exact solution changes slowly when x is large, this leftward movement of the collocation points is beneficial since more collocation points are situated where the solution is changing most rapidly. The disadvantage of a logarithmic change of variables is seen in the function ln(1-t) where the growth is so slow near t=+1 that the transformed solution possesses a singularity in a derivative at t=+1. In other words, the jth derivative of a function of the form (1-t)α, where α>0 is not an integer, is singular at t=+1 for j>α. To achieve exponential convergence, Y(t) should be infinitely smooth. For this particular problem, it is seen that for K=1/2,1, and 3/2, Ya(t) and all its derivatives with respect to t approach 0 as t→1-, thus, algebraic transformations ϕa(t) would produce satisfactory results (see Table 1). On the other hand, Yl(t) is relatively smooth with [(1+3)L] derivatives at t=+1 if L is an integer and with [2L] derivatives if L is not an integer ([·] denotes the greatest integer value), although not infinitely smooth, and with a proper value for L, collocation points corresponding to large values of x, where the solution changes slowly, are moved to the left (when compared to the algebraic transformations) where the solution changes more rapidly. As a result, for 10⩽N⩽30, the logarithmic transformation yields a solution that is often more accurate than the other choices for ϕ.
The problem in (4.1) has been solved in [22] using domain truncation and a forth-order finite difference method (FFDM) and later in [23] using the precise time integration (PTI) method. In Table 2, a comparison is made between the values obtained using the present CG collocation method with the logarithmic transformation ϕl for L=5 and N=20, method in [22] with truncated domain [0,10] and h=1/64, and method in [23] with matrix exponential technique. It is seen that our method is more accurate than the FFDM and PTI methods, and even though the exact results become extremely small when x>6, our method is still capable of giving results with high precision.Table 1
Maximum absolute errors for Example4.1.
MappingNL123456ϕa(t)  K=1/2103.31×10-62.22×10-75.96×10-74.65×10-54.22×10-42.08×10-3208.48×10-83.42×10-109.20×10-122.95×10-131.16×10-134.12×10-11303.10×10-99.60×10-138.51×10-158.40×10-158.00×10-158.00×10-15ϕa(t)  K=1101.32×10-41.14×10-54.05×10-64.46×10-62.32×10-55.81×10-5204.18×10-74.00×10-81.34×10-87.20×10-96.43×10-94.77×10-9306.58×10-91.71×10-101.45×10-111.05×10-112.53×10-129.52×10-13ϕa(t)  K=3/2105.75×10-42.78×10-42.21×10-43.20×10-42.38×10-41.75×10-4205.47×10-68.18×10-73.55×10-76.47×10-89.80×10-84.63×10-8305.42×10-86.15×10-92.10×10-91.29×10-97.22×10-101.88×10-10ϕl(t)101.06×10-65.41×10-96.18×10-102.83×10-95.55×10-68.42×10-5202.72×10-83.11×10-124.25×10-152.95×10-151.80×10-156.10×10-15303.12×10-94.61×10-143.52×10-151.65×10-151.05×10-155.07×10-15Table 2
Approximate and exact values ofy(x) for Example 4.1.
xFFDM [22]PTI [23]Present methodExact11.00220×10-11.0021048×10-11.002104788741×10-11.002104788741×10-121.12776×10-21.1275892×10-21.127589195768×10-21.127589195768×10-241.76737×10-41.7670378×10-41.767037762946×10-41.767037762946×10-463.11068×10-63.1101148×10-63.110114828388×10-63.110114828396×10-685.63886×10-85.6428565×10-85.642859756603×10-85.642859752063×10-8108.13692×10-10—1.031258865734×10-91.031258873612×10-912——1.887856605151×10-111.887856203661×10-1114——3.457320525861×10-133.457322449053×10-1315——4.679133135965×10-144.678891146411×10-14Example 4.2.
Consider the following nonlinear semi-infinite-domain two-point boundary value problem:(4.4)1x2y′′(x)+92y2(x)y′(x)+3xy3(x)=0,y(0)=1,y(∞)=0.
The exact analytical solution is ye(x)=1/1+x3. This example was solved for N=(10,20,30,40) using the CG collocation method described in Section 3.2 with the strictly monotonic transformations of the domain t∈[-1,1) given in (2.4)-(2.5).
The maximum absolute errors obtained using the CG scheme of this paper are given in Table3. In this example, the exact solution decays algebraically; thus, algebraic transformations yield more accurate results in comparison with the logarithmic transformation. It is seen from Table 3 that for this problem the algebraic transformations are relatively not sensitive to the mapping parameter L. In addition, the algebraic transformation ϕa with K=3/2 allows accurate results to be obtained more efficiently, and as N increases the errors obtained using this transformation decay more rapidly than the errors obtained using other choices for ϕ. For example, for N=40, using ϕa with K=1 and L=5 the maximum absolute error is obtained as 2.18×10-7, while using ϕa with K=3/2 and L=5 the maximum absolute error is obtained as 1.05×10-9 (see Table 3). An advantage of using the algebraic change of variable ϕa with K=3/2, as compared to other changes of variables, is that it moves collocation points to the right. Because the exact solution decays slowly, this rightward movement of the collocation points is beneficial since more collocation points are situated where x is large.
Now, consider the exact continuous solution in the transformed coordinates using the transformations given in (2.4)-(2.5) as
(4.5)Ya(t)=11+(L(1+t)/(1-t)K)3,Yl(t)=11+(Lln⁡(2/(1-t)))3.
The subscripts a and l correspond to the algebraic and logarithmic transformations. The second reason that the transformation ϕa with K=3/2 yields more accurate results is related to the smoothness of Ya(t) and Yl(t). It can be easily seen that Ya(t) with K=3/2 has two derivatives as t→1-, while Ya(t) with K=1 has one derivative. In addition, the first derivatives of functions Ya(t) with K=1/2 and Yl(t) approach -∞ as t→1-, which induces a relatively large error when these transformations are used.Table 3
Maximum absolute errors for Example4.2.
MappingNL123456ϕa(t)  K=1/2105.84×10-32.23×10-31.31×10-21.79×10-23.62×10-27.23×10-2201.77×10-36.45×10-43.56×10-43.62×10-42.00×10-41.36×10-3307.63×10-42.40×10-41.39×10-49.40×10-56.12×10-54.95×10-5404.10×10-41.38×10-47.75×10-55.53×10-53.60×10-52.74×10-5ϕa(t)  K=1101.85×10-32.71×10-35.81×10-31.22×10-24.83×10-22.28×10-2202.62×10-51.47×10-53.45×10-59.53×10-51.23×10-43.61×10-4307.20×10-62.45×10-61.31×10-61.00×10-62.89×10-66.57×10-6402.84×10-61.00×10-65.55×10-73.40×10-72.18×10-72.60×10-7ϕa(t)  K=3/2107.24×10-37.19×10-31.48×10-21.12×10-21.92×10-26.81×10-2206.95×10-55.47×10-55.75×10-51.07×10-41.81×10-42.24×10-4307.52×10-75.63×10-78.65×10-71.57×10-62.98×10-63.82×10-6401.12×10-83.64×10-91.94×10-91.22×10-91.05×10-91.14×10-9L1510152025ϕl(t)106.78×10-22.31×10-27.22×10-23.75×10-1——206.40×10-21.11×10-23.85×10-35.58×10-31.77×10-21.22×10-2307.20×10-28.00×10-32.95×10-31.65×10-39.96×10-41.56×10-3403.82×10-27.15×10-32.55×10-31.43×10-39.00×10-46.95×10-4Example 4.3.
In this example we consider the nonlinear semi-infinite-domain two-point boundary value problem due to Kidder [2] as
(4.6)y′′(x)+2x1-αyy′(x)=0,0<α<1,
which appears in the problem of the transient flow of gas [2] within a one-dimensional semi-infinite porous medium. The typical boundary conditions imposed by the physical properties are
(4.7)y(0)=1,y(∞)=0.
The problem in (4.6)-(4.7) was handled by Kidder [2] where a perturbation technique is carried out to include terms of the second order. The principal limitation of his work is that the complexity of the calculations increases rapidly with increasing order of the terms. Moreover, Countryman and Kannan [48], without solving the problem in (4.6)-(4.7), showed that the lower and upper solutions β(x) and γ(x) for this problem are obtained from the linear infinite-domain two-point boundary value problems
(4.8)β′′(x)+2x1-αβ′(x)=0,β(0)=1,β(∞)=0,γ′′(x)+2xγ′(x)=0,γ(0)=1,γ(∞)=0,
respectively. It can be easily seen that the solutions β(x) and γ(x) are given by
(4.9)β(x)=1-erf⁡((1-α)-1/4x),γ(x)=1-erf⁡(x),
where the error function is defined by
(4.10)erf⁡(x)=2π∫0xe-s2ds.
Further,
(4.11)β(x)≤y(x)≤γ(x),
where y(x) is the solution to (4.6)-(4.7).
In this example the solution decays exponentially fast asx increases; thus, according to the results obtained in Examples 4.1 and 4.2, we solve this problem using the logarithmic transformation ϕl and also using the algebraic transformation ϕa with K=1/2. It is worth mentioning that for problems whose exact solution is not available, we first solve the problem using the algebraic transformation ϕa with K=1 and N=30 to determine the decay rate of the solution. We then tune the transformations, appropriately, to improve the quality of our discretization technique.
For comparison purposes, the estimated values ofy(x) for specific value of α=0.5 obtained using the perturbation method in [2], the CG scheme of this paper with ϕl for L=4 and N=50 together with the lower and upper solutions [48] given in (4.9), are listed in Table 4. In addition, in Table 5 we give the approximate values of y(x) obtained using the present CG method for α=0.5 and different transformations. It is seen that both transformations ϕa with K=1/2 and ϕl produce almost the same results. Moreover, the values obtained using the present method highly agree with the upper solutions given in Table 4, while the value of y(x) obtained by the perturbation method [2] deviates from the upper solution at x=1. In addition, as expected, changing the values of the mapping parameter L improves the precision of our collocation method. For example, at x=5 the value of y(x) obtained using both transformations for L=1, a bit deviate from the upper solution given in Table 4, while increasing L to L=4, improves the precision (see Table 5). The computational results for y(x) using the present CG scheme with ϕl, L=4, and N=50 together with the lower and upper solutions in (4.9) are shown in Figure 1. Finally, we note that unlike the perturbation technique used by Kidder [2], where the complexity of the calculations increases rapidly with the increase of order of terms, the present CG scheme leads to approximation through a reduced size of calculations.Table 4
Approximate, lower, and upper values ofy(x) for Example 4.3.
xKidder [2]Lower solution [48]Upper solution [48]Present method0.10.88165882830.86644223850.88753708400.8813646513920.20.76630767810.73660049050.77729741080.7658288130750.30.65653799950.61388347180.67137324050.6560006849620.40.55440240320.50112726610.57160764490.5538989417640.50.46136502950.40040596680.47950012220.4609427627610.60.37831093150.31293833170.39614390920.3779815845880.70.30559765460.23909374960.32219880620.3053523384920.80.24313254730.17848503630.25789903530.2429543915860.90.19046236810.13012388310.20309178760.1903342334771.00.15876898260.09260901950.15729920710.1467733125492.0—0.00076936960.00467773500.0042991828513.0—4.5262×10-70.00002209050.0000202927634.0—1.730×10-111.5417×10-81.416256×10-85.0—1.388×10-171.537×10-121.41205×10-12Table 5
Approximate values ofy(x) for Example 4.3 obtained using the present method with N=50.
xϕa  (K=1/2)ϕlL=1L=4L=1L=40.10.8813646513920.8813646513920.8813646513930.8813646513920.20.7658288130740.7658288130750.7658288130760.7658288130750.30.6560006849620.6560006849620.6560006849600.6560006849620.40.5538989417640.5538989417640.5538989417640.5538989417640.50.4609427627600.4609427627610.4609427627620.4609427627610.60.3779815845860.3779815845880.3779815845880.3779815845880.70.3053523384930.3053523384920.3053523384890.3053523384920.80.2429543915880.2429543915860.2429543915850.2429543915860.90.1903342334740.1903342334770.1903342334800.1903342334771.00.1467733125520.1467733125490.1467733125460.1467733125492.00.0042991828540.0042991828510.0042991828560.0042991828513.00.0000202927680.0000202927630.0000202927590.0000202927634.01.415843×10-81.416261×10-81.416474×10-81.416256×10-85.02.30431×10-121.40494×10-121.68339×10-121.41205×10-12Figure 1
Computational results ofy(x), lower and upper solutions for Example 4.3.


### 4.2. Initial Value Problems

Example 4.4.
Consider the Lane-Emden equation given by(4.12)y′′(x)+2xy′(x)+y5(x)=0,y(0)=1,y′(0)=0,
where limx→∞y(x)=0. The exact analytical solution is given in [50] by ye(x)=1/1+x2/3. The problem in (4.12) was solved for N=(20,30,40) using the CGR collocation scheme described in Section 3.3 with the strictly monotonic algebraic transformations of the domain t∈[-1,1) given in (2.4). As stated in Section 3, in the CGR scheme of this paper, the function value at t=+1 (x=+∞) is an unknown variable in the discretization and is approximated after obtaining the approximate solution. The maximum absolute errors together with the approximate values of y(∞) obtained using the CGR scheme of this paper and the algebraic transformations ϕa with K=1/2, 1, and 3/2 are given in Table 6.
Similarly to Example4.2, the exact solution of this problem decays algebraically, and, hence, algebraic transformations allow accurate results to be obtained more efficiently and they are relatively not sensitive to the mapping parameter L. But, unlike Example 4.2, here, ϕa with K=1 is more accurate than ϕa with K=3/2. Again, the different behavior of these transformations is understood by considering the exact solution in the transformed coordinates using the algebraic transformations as
(4.13)Ya(t)=11+(L(1+t)/(1-t)K)2/3.
It can be easily seen that Ya(t) with K=1 is infinitely smooth, while Ya(t) with K=3/2 has just one derivative as t→1- and also Ya′(t) with K=1/2 approaches -∞ as t→1-. As a result, for this particular problem, the algebraic transformation with K=1 yields a solution that is by far more accurate than other choices for ϕ.Table 6
Maximum absolute errors and approximations ofy(∞) for Example 4.4.
MappingNL(L,y(∞))1234ϕa(t)  K=1/2201.87×10-29.32×10-36.21×10-34.65×10-3(4,4.6534×10-3)301.26×10-26.30×10-34.20×10-33.15×10-3(4,3.1477×10-3)409.53×10-34.76×10-33.17×10-32.38×10-3(4,2.3788×10-3)ϕa(t)  K=1206.84×10-93.36×10-98.00×10-93.32×10-8(2,8.1966×10-11)306.65×10-133.34×10-131.79×10-128.00×10-12(2,3.1491×10-15)403.50×10-154.87×10-155.25×10-155.76×10-15(1,5.6380×10-16)ϕa(t)  K=3/2207.10×10-53.55×10-52.34×10-51.76×10-5(4,1.0517×10-5)302.20×10-51.09×10-57.25×10-65.42×10-6(4,3.2707×10-6)405.43×10-63.14×10-64.22×10-73.15×10-7(4,1.2510×10-7)Example 4.5.
Finally, in this example we consider the Flierl-Petviashvili (FP) equation [17, 18]. The FP monopole y(x) is the radially symmetric solution to the equation
(4.14)y′′(x)+1xy′(x)-y(x)-y2(x)=0,y(0)=α,y′(0)=0,
where the boundary condition is y(∞)=0. Flierl [17] and Petviashvili [18] independently derived this equation as a model for vortex solitons in the ocean and in plasmas, respectively. Equation (4.14) was investigated by Boyd in [20] and later by Wazwaz in [49].
In this problem the FP monopoleα is an unknown value that should be determined. Boyd [20] and Wazwaz [49] employed the boundary condition y(∞)=0 in combination with the Pade approximants in order to approximate the FP monopole α. In this paper, however, we do not use the boundary condition y(∞)=0. In fact, here, we aim to solve the following problem:
(4.15)y′′(x)+1xy′(x)-y(x)-y2(x)=0,y′(0)=0,
and then we approximate the value of α=y(0).
Our procedure for solving this problem and approximating the FP monopoleα is obtained by a little change in our CGR scheme as follows. After expressing the FP equation in (4.15) in terms of t using (2.7), in the first stage, we set N=5 and evaluate the mentioned equation at collocation points (t1,t2,…,t5) (in general (t1,t2,…,tN)) and noncollocated point t6=1 (in general tN+1=1), by which in combination with the initial slope condition in (3.21) we obtain a system of 7 (in general N+2) algebraic equations for unknowns Y0,Y1,…,Y6 (in general Y0,Y1,…,YN+1). Note that evaluating (2.7) at noncollocated point tN+1=1 is not versus our CGR collocation scheme because at tN+1=1 the coefficients of Y′(tN+1) and Y′′(tN+1) in (2.7) using the transformations of this paper will be zero and there is no need to approximate these derivatives based on collocation points. For N=5, the obtained system of algebraic equations has 64 (in general 2N+1) solution sets. These solution sets can be obtained using any symbolic software. We obtained them using “Mathematica” with the command “NSolve.” Now we evaluate each of the solution sets. The complex solution sets should be discarded. Next, as Boyd [20] stated from physical grounds, the solution y(x) for FP equation is negative and increasing; thus, between the real solution sets, after discarding the solution sets that are oscillatory or contain positive elements, the acceptable solution set is obtained. We then increase N by 5, but in this stage we solve the obtained system of algebraic equations using Newton’s iterative method and utilize the approximate solution at stage one as initial guess for Newton’s method. We continue this approach until the results for α≈Y0 are similar up to a required number of decimal places for two consecutive stages.
Solving this problem using the above procedure and usingϕa with K=1 and N=5 reveals that the exact solution of the FP equation decays almost exponentially; thus, for this problem both algebraic and logarithmic changes of variables would work. In this example we set the mapping parameter L=2, which provides faster convergence. Table 7 exhibits the approximated and exact values of FP monopole α. It is observed that using both algebraic and logarithmic transformations of this paper, the approximations converge monotonically to the exact α=y(0) as N increases. However, the algebraic transformation ϕa with K=1/2 provides more rapid convergence in comparison with other transformations. Moreover, it is seen that our procedure produces more accurate results than those obtained in [20, 49] by using the Pade approximants. The computational result for y(x) using the present CGR scheme with N=40, L=2, ϕa, and K=1/2 is shown in Figure 2.Table 7
Approximate and exact values ofα=y(0) for Example 4.5.
Nϕa  (K=1/2)ϕa  (K=1)ϕa  (K=3/2)ϕl5−2.416296577−2.488733071−2.221866799−2.38563630810−2.391957269−2.391760175−2.387175116−2.39195708315−2.391956196−2.391955001−2.391806999−2.39195629020−2.391956398−2.391956374−2.391958288−2.39195642825−2.391956405−2.391956410−2.391956155−2.39195639730−2.391956403−2.391956404−2.391956426−2.39195640535−2.391956403−2.391956403−2.391956405−2.39195640240−2.391956403−2.391956403−2.391956403−2.391956404Wazwaz [49], Pade [8/8]−2.3922139Boyd [20], Pade [10/10]−2.3919746Exact−2.3919564Figure 2
Computational results ofy(x) for Example 4.5 using the CGR scheme for ϕa, K=1/2, L=2, and N=40.


## 4.1. Semi-Infinite-Domain Two-Point Boundary Value Problems

Example 4.1.
We begin with a simple example to illustrate and analyze tuning the transformations. Consider the linear semi-infinite-domain two-point boundary value problem [22]
(4.1)-y′′(x)-2y′(x)+2y(x)=e-2x,y(0)=1,y(∞)=0.
The exact analytical solution is ye(x)=(1/2)(e-(1+3)x+e-2x). The problem in (4.1) was solved for N=(10,20,30) using the collocation method described in Section 3.2 with the CG collocation points described in Section 3.1 and with the strictly monotonic transformations of the domain t∈[-1,1) given in (2.4)-(2.5).
The maximum absolute error is defined by(4.2)EN=max⁡{|Ye(t)-IN+1(Y)(t)|:-1⩽t⩽1},
where IN+1(Y)(t) and Ye(t) are the approximated and exact solutions in the transformed coordinate, respectively. The maximum absolute errors obtained using the CG scheme of this paper, are given in Table 1 with the logarithmic transformation and the algebraic transformations with K=1/2, 1,  and 3/2 given in (2.4)-(2.5). As can be seen, for all four transformations, the errors decrease until N=30, demonstrating an approximately exponential convergence rate. It is seen that, for each value of N, there exists an optimal value for the mapping parameter L and this optimal value varies with the transformation. In addition, it is observed that changing the values of mapping parameter L improves the quality of the discretization, significantly. For example, for N=20, the maximum absolute error using the logarithmic transformation ϕl decreases from 2.72×10-8 for L=1 to 1.80×10-15 for L=5. Moreover, the logarithmic transformation ϕl allows accurate results to be obtained much more efficiently. For example, to achieve an error of less than 10-14 we require less than 20 collocation points using the logarithmic transformation with L=3, while for other transformations more collocation points are needed (see Table 1). However, the algebraic transformation ϕa with K=1/2 produces very satisfactory results too. Finally, it is seen that for each value of N, the best maximum absolute error for logarithmic transformation is less than that for algebraic transformations.
The different behavior of the transformations given in (2.4)-(2.5) is understood if we apply the change of variables to the continuous solution. The exact solution in the transformed coordinates is as follows:
(4.3)Ya(t)=12(e-(1+3)L((1+t)/(1-t)K)+e-2L((1+t)/(1-t)K)),Yl(t)=12((1-t2)(1+3)L+(1-t2)2L).
Here the subscripts a and l correspond to the algebraic and logarithmic transformations given in (2.4)-(2.5). An advantage of using the logarithmic change of variable and also the algebraic change of variable with K=1/2, as compared to other algebraic change of variables, is that for proper values of L they essentially move collocation points associated with large values of x to the left. Because the exact solution changes slowly when x is large, this leftward movement of the collocation points is beneficial since more collocation points are situated where the solution is changing most rapidly. The disadvantage of a logarithmic change of variables is seen in the function ln(1-t) where the growth is so slow near t=+1 that the transformed solution possesses a singularity in a derivative at t=+1. In other words, the jth derivative of a function of the form (1-t)α, where α>0 is not an integer, is singular at t=+1 for j>α. To achieve exponential convergence, Y(t) should be infinitely smooth. For this particular problem, it is seen that for K=1/2,1, and 3/2, Ya(t) and all its derivatives with respect to t approach 0 as t→1-, thus, algebraic transformations ϕa(t) would produce satisfactory results (see Table 1). On the other hand, Yl(t) is relatively smooth with [(1+3)L] derivatives at t=+1 if L is an integer and with [2L] derivatives if L is not an integer ([·] denotes the greatest integer value), although not infinitely smooth, and with a proper value for L, collocation points corresponding to large values of x, where the solution changes slowly, are moved to the left (when compared to the algebraic transformations) where the solution changes more rapidly. As a result, for 10⩽N⩽30, the logarithmic transformation yields a solution that is often more accurate than the other choices for ϕ.
The problem in (4.1) has been solved in [22] using domain truncation and a forth-order finite difference method (FFDM) and later in [23] using the precise time integration (PTI) method. In Table 2, a comparison is made between the values obtained using the present CG collocation method with the logarithmic transformation ϕl for L=5 and N=20, method in [22] with truncated domain [0,10] and h=1/64, and method in [23] with matrix exponential technique. It is seen that our method is more accurate than the FFDM and PTI methods, and even though the exact results become extremely small when x>6, our method is still capable of giving results with high precision.Table 1
Maximum absolute errors for Example4.1.
MappingNL123456ϕa(t)  K=1/2103.31×10-62.22×10-75.96×10-74.65×10-54.22×10-42.08×10-3208.48×10-83.42×10-109.20×10-122.95×10-131.16×10-134.12×10-11303.10×10-99.60×10-138.51×10-158.40×10-158.00×10-158.00×10-15ϕa(t)  K=1101.32×10-41.14×10-54.05×10-64.46×10-62.32×10-55.81×10-5204.18×10-74.00×10-81.34×10-87.20×10-96.43×10-94.77×10-9306.58×10-91.71×10-101.45×10-111.05×10-112.53×10-129.52×10-13ϕa(t)  K=3/2105.75×10-42.78×10-42.21×10-43.20×10-42.38×10-41.75×10-4205.47×10-68.18×10-73.55×10-76.47×10-89.80×10-84.63×10-8305.42×10-86.15×10-92.10×10-91.29×10-97.22×10-101.88×10-10ϕl(t)101.06×10-65.41×10-96.18×10-102.83×10-95.55×10-68.42×10-5202.72×10-83.11×10-124.25×10-152.95×10-151.80×10-156.10×10-15303.12×10-94.61×10-143.52×10-151.65×10-151.05×10-155.07×10-15Table 2
Approximate and exact values ofy(x) for Example 4.1.
xFFDM [22]PTI [23]Present methodExact11.00220×10-11.0021048×10-11.002104788741×10-11.002104788741×10-121.12776×10-21.1275892×10-21.127589195768×10-21.127589195768×10-241.76737×10-41.7670378×10-41.767037762946×10-41.767037762946×10-463.11068×10-63.1101148×10-63.110114828388×10-63.110114828396×10-685.63886×10-85.6428565×10-85.642859756603×10-85.642859752063×10-8108.13692×10-10—1.031258865734×10-91.031258873612×10-912——1.887856605151×10-111.887856203661×10-1114——3.457320525861×10-133.457322449053×10-1315——4.679133135965×10-144.678891146411×10-14Example 4.2.
Consider the following nonlinear semi-infinite-domain two-point boundary value problem:(4.4)1x2y′′(x)+92y2(x)y′(x)+3xy3(x)=0,y(0)=1,y(∞)=0.
The exact analytical solution is ye(x)=1/1+x3. This example was solved for N=(10,20,30,40) using the CG collocation method described in Section 3.2 with the strictly monotonic transformations of the domain t∈[-1,1) given in (2.4)-(2.5).
The maximum absolute errors obtained using the CG scheme of this paper are given in Table3. In this example, the exact solution decays algebraically; thus, algebraic transformations yield more accurate results in comparison with the logarithmic transformation. It is seen from Table 3 that for this problem the algebraic transformations are relatively not sensitive to the mapping parameter L. In addition, the algebraic transformation ϕa with K=3/2 allows accurate results to be obtained more efficiently, and as N increases the errors obtained using this transformation decay more rapidly than the errors obtained using other choices for ϕ. For example, for N=40, using ϕa with K=1 and L=5 the maximum absolute error is obtained as 2.18×10-7, while using ϕa with K=3/2 and L=5 the maximum absolute error is obtained as 1.05×10-9 (see Table 3). An advantage of using the algebraic change of variable ϕa with K=3/2, as compared to other changes of variables, is that it moves collocation points to the right. Because the exact solution decays slowly, this rightward movement of the collocation points is beneficial since more collocation points are situated where x is large.
Now, consider the exact continuous solution in the transformed coordinates using the transformations given in (2.4)-(2.5) as
(4.5)Ya(t)=11+(L(1+t)/(1-t)K)3,Yl(t)=11+(Lln⁡(2/(1-t)))3.
The subscripts a and l correspond to the algebraic and logarithmic transformations. The second reason that the transformation ϕa with K=3/2 yields more accurate results is related to the smoothness of Ya(t) and Yl(t). It can be easily seen that Ya(t) with K=3/2 has two derivatives as t→1-, while Ya(t) with K=1 has one derivative. In addition, the first derivatives of functions Ya(t) with K=1/2 and Yl(t) approach -∞ as t→1-, which induces a relatively large error when these transformations are used.Table 3
Maximum absolute errors for Example4.2.
MappingNL123456ϕa(t)  K=1/2105.84×10-32.23×10-31.31×10-21.79×10-23.62×10-27.23×10-2201.77×10-36.45×10-43.56×10-43.62×10-42.00×10-41.36×10-3307.63×10-42.40×10-41.39×10-49.40×10-56.12×10-54.95×10-5404.10×10-41.38×10-47.75×10-55.53×10-53.60×10-52.74×10-5ϕa(t)  K=1101.85×10-32.71×10-35.81×10-31.22×10-24.83×10-22.28×10-2202.62×10-51.47×10-53.45×10-59.53×10-51.23×10-43.61×10-4307.20×10-62.45×10-61.31×10-61.00×10-62.89×10-66.57×10-6402.84×10-61.00×10-65.55×10-73.40×10-72.18×10-72.60×10-7ϕa(t)  K=3/2107.24×10-37.19×10-31.48×10-21.12×10-21.92×10-26.81×10-2206.95×10-55.47×10-55.75×10-51.07×10-41.81×10-42.24×10-4307.52×10-75.63×10-78.65×10-71.57×10-62.98×10-63.82×10-6401.12×10-83.64×10-91.94×10-91.22×10-91.05×10-91.14×10-9L1510152025ϕl(t)106.78×10-22.31×10-27.22×10-23.75×10-1——206.40×10-21.11×10-23.85×10-35.58×10-31.77×10-21.22×10-2307.20×10-28.00×10-32.95×10-31.65×10-39.96×10-41.56×10-3403.82×10-27.15×10-32.55×10-31.43×10-39.00×10-46.95×10-4Example 4.3.
In this example we consider the nonlinear semi-infinite-domain two-point boundary value problem due to Kidder [2] as
(4.6)y′′(x)+2x1-αyy′(x)=0,0<α<1,
which appears in the problem of the transient flow of gas [2] within a one-dimensional semi-infinite porous medium. The typical boundary conditions imposed by the physical properties are
(4.7)y(0)=1,y(∞)=0.
The problem in (4.6)-(4.7) was handled by Kidder [2] where a perturbation technique is carried out to include terms of the second order. The principal limitation of his work is that the complexity of the calculations increases rapidly with increasing order of the terms. Moreover, Countryman and Kannan [48], without solving the problem in (4.6)-(4.7), showed that the lower and upper solutions β(x) and γ(x) for this problem are obtained from the linear infinite-domain two-point boundary value problems
(4.8)β′′(x)+2x1-αβ′(x)=0,β(0)=1,β(∞)=0,γ′′(x)+2xγ′(x)=0,γ(0)=1,γ(∞)=0,
respectively. It can be easily seen that the solutions β(x) and γ(x) are given by
(4.9)β(x)=1-erf⁡((1-α)-1/4x),γ(x)=1-erf⁡(x),
where the error function is defined by
(4.10)erf⁡(x)=2π∫0xe-s2ds.
Further,
(4.11)β(x)≤y(x)≤γ(x),
where y(x) is the solution to (4.6)-(4.7).
In this example the solution decays exponentially fast asx increases; thus, according to the results obtained in Examples 4.1 and 4.2, we solve this problem using the logarithmic transformation ϕl and also using the algebraic transformation ϕa with K=1/2. It is worth mentioning that for problems whose exact solution is not available, we first solve the problem using the algebraic transformation ϕa with K=1 and N=30 to determine the decay rate of the solution. We then tune the transformations, appropriately, to improve the quality of our discretization technique.
For comparison purposes, the estimated values ofy(x) for specific value of α=0.5 obtained using the perturbation method in [2], the CG scheme of this paper with ϕl for L=4 and N=50 together with the lower and upper solutions [48] given in (4.9), are listed in Table 4. In addition, in Table 5 we give the approximate values of y(x) obtained using the present CG method for α=0.5 and different transformations. It is seen that both transformations ϕa with K=1/2 and ϕl produce almost the same results. Moreover, the values obtained using the present method highly agree with the upper solutions given in Table 4, while the value of y(x) obtained by the perturbation method [2] deviates from the upper solution at x=1. In addition, as expected, changing the values of the mapping parameter L improves the precision of our collocation method. For example, at x=5 the value of y(x) obtained using both transformations for L=1, a bit deviate from the upper solution given in Table 4, while increasing L to L=4, improves the precision (see Table 5). The computational results for y(x) using the present CG scheme with ϕl, L=4, and N=50 together with the lower and upper solutions in (4.9) are shown in Figure 1. Finally, we note that unlike the perturbation technique used by Kidder [2], where the complexity of the calculations increases rapidly with the increase of order of terms, the present CG scheme leads to approximation through a reduced size of calculations.Table 4
Approximate, lower, and upper values ofy(x) for Example 4.3.
xKidder [2]Lower solution [48]Upper solution [48]Present method0.10.88165882830.86644223850.88753708400.8813646513920.20.76630767810.73660049050.77729741080.7658288130750.30.65653799950.61388347180.67137324050.6560006849620.40.55440240320.50112726610.57160764490.5538989417640.50.46136502950.40040596680.47950012220.4609427627610.60.37831093150.31293833170.39614390920.3779815845880.70.30559765460.23909374960.32219880620.3053523384920.80.24313254730.17848503630.25789903530.2429543915860.90.19046236810.13012388310.20309178760.1903342334771.00.15876898260.09260901950.15729920710.1467733125492.0—0.00076936960.00467773500.0042991828513.0—4.5262×10-70.00002209050.0000202927634.0—1.730×10-111.5417×10-81.416256×10-85.0—1.388×10-171.537×10-121.41205×10-12Table 5
Approximate values ofy(x) for Example 4.3 obtained using the present method with N=50.
xϕa  (K=1/2)ϕlL=1L=4L=1L=40.10.8813646513920.8813646513920.8813646513930.8813646513920.20.7658288130740.7658288130750.7658288130760.7658288130750.30.6560006849620.6560006849620.6560006849600.6560006849620.40.5538989417640.5538989417640.5538989417640.5538989417640.50.4609427627600.4609427627610.4609427627620.4609427627610.60.3779815845860.3779815845880.3779815845880.3779815845880.70.3053523384930.3053523384920.3053523384890.3053523384920.80.2429543915880.2429543915860.2429543915850.2429543915860.90.1903342334740.1903342334770.1903342334800.1903342334771.00.1467733125520.1467733125490.1467733125460.1467733125492.00.0042991828540.0042991828510.0042991828560.0042991828513.00.0000202927680.0000202927630.0000202927590.0000202927634.01.415843×10-81.416261×10-81.416474×10-81.416256×10-85.02.30431×10-121.40494×10-121.68339×10-121.41205×10-12Figure 1
Computational results ofy(x), lower and upper solutions for Example 4.3.


## 4.2. Initial Value Problems

Example 4.4.
Consider the Lane-Emden equation given by(4.12)y′′(x)+2xy′(x)+y5(x)=0,y(0)=1,y′(0)=0,
where limx→∞y(x)=0. The exact analytical solution is given in [50] by ye(x)=1/1+x2/3. The problem in (4.12) was solved for N=(20,30,40) using the CGR collocation scheme described in Section 3.3 with the strictly monotonic algebraic transformations of the domain t∈[-1,1) given in (2.4). As stated in Section 3, in the CGR scheme of this paper, the function value at t=+1 (x=+∞) is an unknown variable in the discretization and is approximated after obtaining the approximate solution. The maximum absolute errors together with the approximate values of y(∞) obtained using the CGR scheme of this paper and the algebraic transformations ϕa with K=1/2, 1, and 3/2 are given in Table 6.
Similarly to Example4.2, the exact solution of this problem decays algebraically, and, hence, algebraic transformations allow accurate results to be obtained more efficiently and they are relatively not sensitive to the mapping parameter L. But, unlike Example 4.2, here, ϕa with K=1 is more accurate than ϕa with K=3/2. Again, the different behavior of these transformations is understood by considering the exact solution in the transformed coordinates using the algebraic transformations as
(4.13)Ya(t)=11+(L(1+t)/(1-t)K)2/3.
It can be easily seen that Ya(t) with K=1 is infinitely smooth, while Ya(t) with K=3/2 has just one derivative as t→1- and also Ya′(t) with K=1/2 approaches -∞ as t→1-. As a result, for this particular problem, the algebraic transformation with K=1 yields a solution that is by far more accurate than other choices for ϕ.Table 6
Maximum absolute errors and approximations ofy(∞) for Example 4.4.
MappingNL(L,y(∞))1234ϕa(t)  K=1/2201.87×10-29.32×10-36.21×10-34.65×10-3(4,4.6534×10-3)301.26×10-26.30×10-34.20×10-33.15×10-3(4,3.1477×10-3)409.53×10-34.76×10-33.17×10-32.38×10-3(4,2.3788×10-3)ϕa(t)  K=1206.84×10-93.36×10-98.00×10-93.32×10-8(2,8.1966×10-11)306.65×10-133.34×10-131.79×10-128.00×10-12(2,3.1491×10-15)403.50×10-154.87×10-155.25×10-155.76×10-15(1,5.6380×10-16)ϕa(t)  K=3/2207.10×10-53.55×10-52.34×10-51.76×10-5(4,1.0517×10-5)302.20×10-51.09×10-57.25×10-65.42×10-6(4,3.2707×10-6)405.43×10-63.14×10-64.22×10-73.15×10-7(4,1.2510×10-7)Example 4.5.
Finally, in this example we consider the Flierl-Petviashvili (FP) equation [17, 18]. The FP monopole y(x) is the radially symmetric solution to the equation
(4.14)y′′(x)+1xy′(x)-y(x)-y2(x)=0,y(0)=α,y′(0)=0,
where the boundary condition is y(∞)=0. Flierl [17] and Petviashvili [18] independently derived this equation as a model for vortex solitons in the ocean and in plasmas, respectively. Equation (4.14) was investigated by Boyd in [20] and later by Wazwaz in [49].
In this problem the FP monopoleα is an unknown value that should be determined. Boyd [20] and Wazwaz [49] employed the boundary condition y(∞)=0 in combination with the Pade approximants in order to approximate the FP monopole α. In this paper, however, we do not use the boundary condition y(∞)=0. In fact, here, we aim to solve the following problem:
(4.15)y′′(x)+1xy′(x)-y(x)-y2(x)=0,y′(0)=0,
and then we approximate the value of α=y(0).
Our procedure for solving this problem and approximating the FP monopoleα is obtained by a little change in our CGR scheme as follows. After expressing the FP equation in (4.15) in terms of t using (2.7), in the first stage, we set N=5 and evaluate the mentioned equation at collocation points (t1,t2,…,t5) (in general (t1,t2,…,tN)) and noncollocated point t6=1 (in general tN+1=1), by which in combination with the initial slope condition in (3.21) we obtain a system of 7 (in general N+2) algebraic equations for unknowns Y0,Y1,…,Y6 (in general Y0,Y1,…,YN+1). Note that evaluating (2.7) at noncollocated point tN+1=1 is not versus our CGR collocation scheme because at tN+1=1 the coefficients of Y′(tN+1) and Y′′(tN+1) in (2.7) using the transformations of this paper will be zero and there is no need to approximate these derivatives based on collocation points. For N=5, the obtained system of algebraic equations has 64 (in general 2N+1) solution sets. These solution sets can be obtained using any symbolic software. We obtained them using “Mathematica” with the command “NSolve.” Now we evaluate each of the solution sets. The complex solution sets should be discarded. Next, as Boyd [20] stated from physical grounds, the solution y(x) for FP equation is negative and increasing; thus, between the real solution sets, after discarding the solution sets that are oscillatory or contain positive elements, the acceptable solution set is obtained. We then increase N by 5, but in this stage we solve the obtained system of algebraic equations using Newton’s iterative method and utilize the approximate solution at stage one as initial guess for Newton’s method. We continue this approach until the results for α≈Y0 are similar up to a required number of decimal places for two consecutive stages.
Solving this problem using the above procedure and usingϕa with K=1 and N=5 reveals that the exact solution of the FP equation decays almost exponentially; thus, for this problem both algebraic and logarithmic changes of variables would work. In this example we set the mapping parameter L=2, which provides faster convergence. Table 7 exhibits the approximated and exact values of FP monopole α. It is observed that using both algebraic and logarithmic transformations of this paper, the approximations converge monotonically to the exact α=y(0) as N increases. However, the algebraic transformation ϕa with K=1/2 provides more rapid convergence in comparison with other transformations. Moreover, it is seen that our procedure produces more accurate results than those obtained in [20, 49] by using the Pade approximants. The computational result for y(x) using the present CGR scheme with N=40, L=2, ϕa, and K=1/2 is shown in Figure 2.Table 7
Approximate and exact values ofα=y(0) for Example 4.5.
Nϕa  (K=1/2)ϕa  (K=1)ϕa  (K=3/2)ϕl5−2.416296577−2.488733071−2.221866799−2.38563630810−2.391957269−2.391760175−2.387175116−2.39195708315−2.391956196−2.391955001−2.391806999−2.39195629020−2.391956398−2.391956374−2.391958288−2.39195642825−2.391956405−2.391956410−2.391956155−2.39195639730−2.391956403−2.391956404−2.391956426−2.39195640535−2.391956403−2.391956403−2.391956405−2.39195640240−2.391956403−2.391956403−2.391956403−2.391956404Wazwaz [49], Pade [8/8]−2.3922139Boyd [20], Pade [10/10]−2.3919746Exact−2.3919564Figure 2
Computational results ofy(x) for Example 4.5 using the CGR scheme for ϕa, K=1/2, L=2, and N=40.


## 5. Conclusions

Two collocation methods have been presented for the numerical solution of nonlinear semi-infinite-domain two-point boundary value problems and initial value problems using global collocation at Chebyshev-Gauss and Chebyshev-Gauss-Radau points. It was shown that the transformed problems that arise from a transformation followed by Gauss or Radau collocation include an approximation to the solution atx=+∞. The results of this paper indicate that the use of Chebyshev-Gauss and Chebyshev-Gauss-Radau points lead to accurate approximations to continuous nonlinear semi-infinite-domain boundary and initial value problems. Algebraic and logarithmic transformations were utilized and a new category of algebraic transformations was introduced. For problems whose exact solution decays algebraically, the transformations were relatively insensitive to the mapping parameter L. Furthermore, by tuning the transformations used to map the semi-infinite domain [0,∞) to a finite interval, it was possible to improve the accuracy in the discrete approximation, significantly. However, the proper choice of transformation should be based on the criterions that the solution to the problem be smooth in the transformed coordinate and the transformation should adjust spacing of the collocation points based on the decay rate of the solution.



---

*Source: 696574-2012-05-14.xml*