# A Legendre Wavelet Spectral Collocation Method for Solving Oscillatory Initial Value Problems

Authors: A. Karimi Dizicheh; F. Ismail; M. Tavassoli Kajani; Mohammad Maleki
Journal: Journal of Applied Mathematics
 (2013)

---

## Abstract

In this paper, we propose an iterative spectral method for solving differential equations with initial values on large intervals. In the proposed method, we first extend the Legendre
wavelet suitable for large intervals, and then the Legendre-Guass collocation points of the Legendre wavelet are derived. Using this strategy, the iterative spectral method converts the differential equation to a set of algebraic equations. Solving these algebraic equations yields an approximate solution for the differential equation. The proposed method is illustrated by some numerical examples, and the result is compared with the exponentially fitted Runge-Kutta method. Our proposed method is simple and highly accurate.


---

## Body

## 1. Introduction

In this paper, we are concerned with the numerical solution of the initial value problem:(1)F(x,y(x),y′(x),…,y(l-1)(x),y(l)(x))=0,x∈[0,T),
with the initial conditions
(2)y(x0)=y0,…,y(l-1)(x0)=y0(l-1),
on a large domain [0,T). This kind of initial value problems appear in many practical life models and are fundamentally important in both theory and applications. It has applications in many branches of science including celestial mechanics, fluid mechanics, heat wave equations, astrophysics, quantum chemistry, and electronics; for example, see [1–4]. Therefore, it is worthwhile to develop numerical methods to approximate their solutions.Spectral methods are very powerful tools for solving many kinds of differential equations (and recently integral equations) arising in various fields of science and engineering [5, 6]. There are two main advantages for spectral methods. One is associated with high accuracy, that is also called “exponential convergence" meaning that the errors are exponentially small. The second advantage is that they can be easily implemented. These are two effective properties which have encouraged many active researchers to use them for different equations. Among many types of spectral methods, collocation methods are specific approaches that are more applicable and widely used. Spectral collocation methods have been used by many authors; for examples, see [7–12].In recent years, wavelets have found their way into many different fields of theoretical and practical sciences. Many researchers started using various wavelets [13–15] for analyzing problems of high computational complexity. It is proved that wavelets are powerful tools for exploring new problems and solving differential equations. The Haar wavelets are alternative tools that have been studied for similar purposes; see [16–19].In this paper, we introduce a reliable algorithm based on the Legendre wavelets spectral method to obtain numerical solutions of nonlinear ordinary differential equations on large intervals. There is a difficulty involved with a direct usage of spectral methods and their numerical integration for this type of equations. By introducing a new Legendre wavelet interpolant approximation, these difficulties rest mainly with the improved Legendre waveletspectral method. The method of this paper solves the problem iteratively. The validation of this technique is especially tested for solutions having oscillatory behavior. Numerical examples are included to demonstrate the efficiency and high accuracy of the proposed method.The remainder of the paper is organized as follows. In Section2, we explain the basic properties of the Legendre wavelets and Legendre polynomials. The family of the Legendre wavelet is an important example; see Kajani and Vencheh [20]. Then we introduce an interpolation with the Legendre wavelet on shifted Legendre-gauss points and its coefficients. We apply the Legendre wavelets on the spectral method for solving differential equations in large intervals in Section 3. Section 4 is devoted to the numerical results, compared to the exponentially fitted Runge-Kutta method; see Berghe et al. [21, 22]. Finally, in Section 5, the report ends with a brief conclusion.


## 2. The Legendre Wavelets Spectral Method

In this section, we describe spectral method on the Legendre wavelets. First, the Legendre wavelets on large domain and then spectral method are briefly introduced. Next, we state how the interpolation can be performed using the Legendre wavelets.


### 2.1. Review of Legendre Wavelets and Legendre Polynomials

Continuous wavelets transform (see [20]) is defined by
(3)ψa,b(t)=|a|-1/2ψ(t-ba),a,b∈ℝ,a≠0,
for any continuous function ψ, where a is called the scale and b is referred by the shift. A discrete wavelet transform is given by
(4)ψm,n(t)=a-m/2ψ(a-mt-nb),
in which the scale and shift are replaced by (am,namb) for any m,n∈ℤ. Here, we intend to consider a family of discrete wavelets ψm,n(t) when ψ is derived from a Legendre polynomial, and it depends on five arguments, namely, (k,T,n,m,t); T is a large integer, m is the order of the Legendre polynomial, t stands for the normalized time,
(5)k=1,2,…,T,n=1,2,…,2k-1T.
Therefore, we define the Legendre wavelets on the interval [0,T) by
(6)ψn,m(t)={(m+12)1/22k/2Lm(2kt-2n+1)n-12k-1≤t<n2k-1,0otherwise,
where Lm(t) is the famous mth-order Legendre polynomial with the orthogonal weight function w(t)=1. An orthonormal set of the Legendre wavelets is presented in [20]. The Legendre polynomials are obtained by the recursive formulas:
(7)L0(t)=1,L1(t)=t,Lm+1(t)=2m+1m+1tLm(t)-mm+1Lm-1(t),m=1,2,….
Further, the Legendre-Gauss quadrature formula is defined as follows
(8)∫-11‍g(x)dx≃∑j=0M-1‍wjg(xj),
where the Legendre-Gauss collocation points x0<x1<⋯<xM-1 are the roots of LM(x) in (-1,1) and {wj}j=0M-1 are corresponding weights. While explicit formulas for the quadrature nodes are not known, the quadrature weights can be expressed by the following relation:
(9)wj=2(1-xj2)(LM′(xj))2,j=0,1,…,M-1.
Thanks to the property of the Legendre-Gauss quadrature formula, it follows that
(10)∫ab‍f(x)dx=b-a2∫-11‍f((b-a)x+a+b2)dx≃b-a2∑j=0M-1wjf((b-a)xj+a+b2)=∑j=0M-1‍w-jf((b-a)xj+a+b2),
where
(11)w-j=(b-a)wj2,j=0,…,M-1.
For more details about the Legendre polynomials, see [6].


### 2.2. Interpolation by the Legendre Wavelets

A functionf∈L2[0,T) can be expanded in terms of the Legendre wavelets as
(12)f(x)=∑n=1∞‍∑m=0∞‍yn,mψn,m(x),
where ψn,m(x) are the Legendre wavelets defined by (6) and yn,m are given by
(13)yn,m=∫(n-1)/2k-1n/2k-1‍f(t)ψn,m(t)dt.
Equation (13) can be approximated using the quadrature rule in (8) as follows:
(14)yn,m≃∑j=0M-1‍w-jf(xnj)ψn,m(xnj),
where
(15)w-j=wj2k,xnj=xj2k+2n-12k,j=0,…,M-1,n=1,…,2k-1T.
Substituting (13) into (12) gives rise to
(16)f(x)=∫(n-1)/2k-1n/2k-1‍f(t)∑m=0∞‍∑n=0∞‍ψn,m(x)ψn,m(t)dt.
Equation (16) gives the “completeness relation” for the basis set, namely,
(17)∑m=0∞‍∑n=0∞‍ψn,m(x)ψn,m(t)=δ(x-t),
and is therefore an identity. If (14) is substituted into (12) (truncated at M-1 and 2k-1T) we have
(18)f(xni)=∑j=0M-1‍f(xnj)w-j∑m=0M-1‍∑n=02k-1T‍ψn,m(xni)ψn,m(xnj).
The discrete version of (17) is given by
(19)w-j∑m=0M-1‍∑n=02k-1T‍ψn,m(xni)ψn,m(xnj)=δij.
Equations (12) (truncated at M-1 and 2k-1T), and (14) also supply interpolation functions Inj(x) governed by
(20)Inj(x)=w-j∑m=0M-1‍ψn,m(xnj)ψn,m(x),j=0,1,…,M-1,n=1,…,2k-1T.
In particular, Inj(xni)=δij. This implies that
(21)f(x)≃∑n=12k-1T‍∑j=0M-1‍Inj(x)f(xnj),
since the following holds:
(22)∑n=12k-1T‍∑m=0M-1‍∑j=0M-1‍w-jf(xnj)ψn,m(xnj)ψn,m(x)=∑n=12k-1T‍∑j=0M-1‍∑m=0M-1‍w-jf(xnj)ψn,m(xnj)ψn,m(x)=∑n=12k-1T‍∑j=0M-1‍w-j∑m=0M-1‍f(xnj)ψn,m(xnj)ψn,m(x).


## 2.1. Review of Legendre Wavelets and Legendre Polynomials

Continuous wavelets transform (see [20]) is defined by
(3)ψa,b(t)=|a|-1/2ψ(t-ba),a,b∈ℝ,a≠0,
for any continuous function ψ, where a is called the scale and b is referred by the shift. A discrete wavelet transform is given by
(4)ψm,n(t)=a-m/2ψ(a-mt-nb),
in which the scale and shift are replaced by (am,namb) for any m,n∈ℤ. Here, we intend to consider a family of discrete wavelets ψm,n(t) when ψ is derived from a Legendre polynomial, and it depends on five arguments, namely, (k,T,n,m,t); T is a large integer, m is the order of the Legendre polynomial, t stands for the normalized time,
(5)k=1,2,…,T,n=1,2,…,2k-1T.
Therefore, we define the Legendre wavelets on the interval [0,T) by
(6)ψn,m(t)={(m+12)1/22k/2Lm(2kt-2n+1)n-12k-1≤t<n2k-1,0otherwise,
where Lm(t) is the famous mth-order Legendre polynomial with the orthogonal weight function w(t)=1. An orthonormal set of the Legendre wavelets is presented in [20]. The Legendre polynomials are obtained by the recursive formulas:
(7)L0(t)=1,L1(t)=t,Lm+1(t)=2m+1m+1tLm(t)-mm+1Lm-1(t),m=1,2,….
Further, the Legendre-Gauss quadrature formula is defined as follows
(8)∫-11‍g(x)dx≃∑j=0M-1‍wjg(xj),
where the Legendre-Gauss collocation points x0<x1<⋯<xM-1 are the roots of LM(x) in (-1,1) and {wj}j=0M-1 are corresponding weights. While explicit formulas for the quadrature nodes are not known, the quadrature weights can be expressed by the following relation:
(9)wj=2(1-xj2)(LM′(xj))2,j=0,1,…,M-1.
Thanks to the property of the Legendre-Gauss quadrature formula, it follows that
(10)∫ab‍f(x)dx=b-a2∫-11‍f((b-a)x+a+b2)dx≃b-a2∑j=0M-1wjf((b-a)xj+a+b2)=∑j=0M-1‍w-jf((b-a)xj+a+b2),
where
(11)w-j=(b-a)wj2,j=0,…,M-1.
For more details about the Legendre polynomials, see [6].


## 2.2. Interpolation by the Legendre Wavelets

A functionf∈L2[0,T) can be expanded in terms of the Legendre wavelets as
(12)f(x)=∑n=1∞‍∑m=0∞‍yn,mψn,m(x),
where ψn,m(x) are the Legendre wavelets defined by (6) and yn,m are given by
(13)yn,m=∫(n-1)/2k-1n/2k-1‍f(t)ψn,m(t)dt.
Equation (13) can be approximated using the quadrature rule in (8) as follows:
(14)yn,m≃∑j=0M-1‍w-jf(xnj)ψn,m(xnj),
where
(15)w-j=wj2k,xnj=xj2k+2n-12k,j=0,…,M-1,n=1,…,2k-1T.
Substituting (13) into (12) gives rise to
(16)f(x)=∫(n-1)/2k-1n/2k-1‍f(t)∑m=0∞‍∑n=0∞‍ψn,m(x)ψn,m(t)dt.
Equation (16) gives the “completeness relation” for the basis set, namely,
(17)∑m=0∞‍∑n=0∞‍ψn,m(x)ψn,m(t)=δ(x-t),
and is therefore an identity. If (14) is substituted into (12) (truncated at M-1 and 2k-1T) we have
(18)f(xni)=∑j=0M-1‍f(xnj)w-j∑m=0M-1‍∑n=02k-1T‍ψn,m(xni)ψn,m(xnj).
The discrete version of (17) is given by
(19)w-j∑m=0M-1‍∑n=02k-1T‍ψn,m(xni)ψn,m(xnj)=δij.
Equations (12) (truncated at M-1 and 2k-1T), and (14) also supply interpolation functions Inj(x) governed by
(20)Inj(x)=w-j∑m=0M-1‍ψn,m(xnj)ψn,m(x),j=0,1,…,M-1,n=1,…,2k-1T.
In particular, Inj(xni)=δij. This implies that
(21)f(x)≃∑n=12k-1T‍∑j=0M-1‍Inj(x)f(xnj),
since the following holds:
(22)∑n=12k-1T‍∑m=0M-1‍∑j=0M-1‍w-jf(xnj)ψn,m(xnj)ψn,m(x)=∑n=12k-1T‍∑j=0M-1‍∑m=0M-1‍w-jf(xnj)ψn,m(xnj)ψn,m(x)=∑n=12k-1T‍∑j=0M-1‍w-j∑m=0M-1‍f(xnj)ψn,m(xnj)ψn,m(x).


## 3. Solving IVPs on a Large Domain

Consider the initial value problem governed by (1) and (2). In order to solve this problem, we divide the interval [0,T) into some subintervals given by [(n-1)/2k-1,n/2k-1) for n=1,…,2k-1T. Consequently, x∈[(n-1)/2k-1,n/2k-1) implies that
(23)ψpj(x)=0,foranyp≠n.
Thus, Ipj(x)=0 for all p≠n. As a result, the Legendre wavelet interpolant approximation to the function y on the nth subinterval follows (21) and is given by
(24)y(x)≃Yn(x)=∑p=12k-1T‍∑j=0M-1‍Ipj(x)ypj=∑j=0M-1‍Inj(x)ynj,forx∈[n-12k-1,n2k-1).
Inductively, we define that
(25)Yn(s)(x):=∑j=0M-1‍Inj(s)(x)ynj,s=1,…,m.
Applying the points {xnj∣n=1,…,2k-1T,j=l,…,M-1} into (1) leads to
(26)F(xnj,Yn(xnj),Yn′(xnj),…,Yn(l-1)(xnj),Yn(l)(xnj))=0.
Besides, the initial conditions for the nth subinterval can be approximated as follows:
(27)limx→((n-1)/2k-1)-Yn-1(s)(x)-Yn(s)(n-12k-1)=0,foranys=0,…,l-1.Equations (26) and (27) give a system of M algebraic equations. By solving this algebraic system, we obtain {ynj∣n=1,…,2k-1T,j=0,…,M-1}. Then, the approximate solution on the nth subinterval is obtained by substituting ynj into (24).Our proposed procedure starts with settingn=1. The initial conditions in (2) are used for the approximate solution for the first subinterval [0,1/2k-1). This allows to obtain the initial conditions for the second subinterval from (27). Consequently, the approximate solution of the second subinterval is computed. The procedure is repeated till the last subinterval. Finally, the approximate solution of the IVP in (1) and (2) on the whole domain [0,T) is obtained by
(28)y(x)≃∑n=12k-1T‍Yn(x)=∑n=12k-1T‍∑j=0M-1‍Inj(x)ynj.


## 4. Numerical Experiments

We apply the method presented in the last section and numerically solve five problems to demonstrate the efficiency of the proposed method. We also compare the results with exponentially fitted Runge-Kutta method [21, 22]. Further, the global errors evaluated at the end point of the interval are tabulated. The following notations are used in the tables.LWSM stands for “the Legendre wavelets on spectral method,”EFRKMB denotes “the exponentially fitted Runge-Kutta method using Vanden Berghe's technique withh=0.0625”; see [21, 22].Example 1.
Consider the initial value problem(29)d2udx2=-u(u2+v2)-3/2,d2vdx2=-v(u2+v2)3/2,
for 0≤x≤7, where u(0)=1,u′(0)=0,v(0)=0,v′(0)=1.The exact solutions of this problem areu(x)=cos(x),v(x)=sin(x). In order to obtain a solid criterion for a global error assessment, we estimate the error at the endpoints. The Euclidean norms of the end-point global errors are provided in Table 1. Table 1 demonstrates that LWSM gives a better approximation than that of EFRKMB.Table 1
Numerical comparison for Example1.
Methods
M
k
Errors
LWSM
10
4
3.30
×
1
0
-
15
12
4
5.30
×
1
0
-
20
10
5
7.36
×
1
0
-
18
12
5
2.85
×
1
0
-
23
EFRKMB
2.03
×
1
0
-
14Example 2.
Let(30)y′′=-30sin(30x),
for 0≤x≤10, where y(0)=0,y′(0)=1.The exact solution of this initial value problem is given byy(x)=sin(30x)/30. Table 2 represents the Euclidean norms of the numerical errors at the endpoints associated with LWSM and EFRKMB. It indicates that our approach is highly more accurate than EFRKMB in this example. Furthermore, the error decreases for any increase in M and k as it is also the case for Example 1.Table 2
Numerical comparison for Example2.
Methods
M
k
Errors
LWSM
26
3
5.69
×
1
0
-
14
27
3
8.24
×
1
0
-
16
25
4
2.11
×
1
0
-
21
27
4
3.48
×
1
0
-
24
20
5
1.96
×
1
0
-
19
25
5
8.32
×
1
0
-
29
EFRKMB
1.70
×
1
0
-
14Example 3.
Let(31)u′′=-u+.001cos(x),
for 0≤x≤1000, where u(0)=1,u′(0)=0. The exact solution of this initial value problem is u(x)=cos(x)+0.0005xsin(x). The numerical data associated with the end-point global errors are provided in Table 3 which give us similar conclusions to Examples 1 and 2 as expected.Table 3
Numerical comparison for Example3.
Methods
M
k
Errors
LWSM
12
4
1.06
×
1
0
-
19
12
3
1.08
×
1
0
-
16
15
3
1.95
×
1
0
-
24
EFRKMB
2.07
×
1
0
-
7Example 4.
Let(32)y′′+y3+y=Bcos(Ωx),y(0)=0.20042672806900,y′(0)=0,
where 0≤x≤300,  B=0.002, and Ω=1.01.The exact solution is governed by(33)y(x)=0.200179477536cos(Ωx)+0.246946143×10-3cos(3Ωx)+0.304016×10-6cos(5Ωx)+0.374×10-9cos(7Ωx).
The numerical comparison of the two methods using the Euclidean norms of the end-point global errors is provided in Table 4. Similar conclusions to those of previous examples are drawn in this table.Table 4
Numerical comparison for Example4.
Methods
M
k
Errors
LWSM
12
4
2.06
×
1
0
-
18
12
3
2.08
×
1
0
-
15
15
3
2.95
×
1
0
-
21
EFRKMB
2.40
×
1
0
-
9


## 5. Conclusions

We define an interpolation using the Legendre wavelets for the shifted Legendre-Guass collocation points. Then, this gives rise to a set of algebraic equations, and finally its solution constitutes the desired approximate solutions. Our approach is based on a partitioning of the domain interval and solving a given IVP iteratively on each interval. This approach is called an iterative spectral method and is suitable for solving IVP on large intervals.The main characteristic behind this technique is that it iteratively solves the IVP and the defined interpolation substantially reduces the computational effort. This leads to highly accurate numerical results as depicted by several examples. The proposed method is simple and easy to implement in practical problems. It can be applied for solving other mathematical equations, such as integral equations, integro-differential equations and partial differential equations.



---

*Source: 591636-2013-05-07.xml*