International Journal of Scientific & Engineering Research, Volume 5, Issue 3, March2014 1348
ISSN 22295518
ChebyshevCollocation Spectral Method for
Class of ParabolicType Volterra Partial
IntegroDifferential Equations
G. I. ElBaghdady, M. S. ElAzab, W. S. ElBeshbeshy.
Abstract— The main purpose of this paper is to introduce a novel numerical method for parabolic Volterra partial integrodifferential equations based on Chebyshevcollocation spectral scheme. In the present work, the parabolic Volterra integrodifferential equation is converted to two coupled equivalent Volterra equations of the second kind. Then, we approximate the integration by replacing the integral function by its interpolating polynomials with Lagrange basis functions in terms of the Chebyshev polynomials instead of using Gauss quadrature approximation to obtain a linear algebraic system. Finally, some numerical examples are presented to illustrate the efficiency and accuracy of the proposed method.
Index Terms— Chebyshev spectral collocation method, Differentiation matrix, Lagrange basis function, Parabolic Volterra integro differential equations (PVIDE).
Cintegrodifferential equation
—————————— ——————————
Izadi and M. Dehghan apply Legendre spectral method to
ONSIDER the one dimensional parabolic partial Volterra
cation in [8], [9] and the references therein. In [10] F. Fakhar
t
¶t u  D u = k (x ,t 
0
with the initial condition
u (x , 0) = u 0 (x ),
s )Du (x , s )ds + f (x ,t ),
x Î Ω,
(1)
(2)
(PVIDE).
Spectral methods have a considerable attention in the last few
years; see [11], [12], [13], [14], and [15]. Spectral methods are
nice and powerful approach for the numerical solution for
and Dirichlet boundary conditions are assigned on the bound ary
ordinary or partial differential equations to high accuracy on a simple domain and if the data defining the problem are
u = 0, (x ,t ) Î ¶Ω´ I ,
where D is the Laplace operator in (x , t ) Î Q º
(3)
Ω´ I , Ω is the
smooth, see [16]. Also Volterra integral and ordinary Volterra
integrodifferential equations have a wide interest by using
bounded interval [−1, 1], with boundary ∂Ω ∈{−1,1}
and
many methods. In [17] H. Brunner used Collocation methods
for secondorder Volterra integrodifferential equations. Mul
I º (0, T ) for a given fixed number T > 0. For implementa
tion of highorder methods such as spectral methods, the known functions f, k and u 0 (x ) are assumed to be sufficiently smooth; real valued functions.
Problems of type (1)(3) arise in many applications; for ex
ample, it describes compression of poroviscoelastic media [1],
[2], the nonlocal reactive flows in porous media [4], [5], [6] and
heat conduction through materials with memory term [2], [3].
Many authors have been considered the numerical solution
of partial Volterra integrodifferential equations by many methods; for example, in [7] (Amiya K. Pani, et al.) use Ga lerkin finite element method and ADI orthogonal spline collo
————————————————
• G. I. ElBaghdady is currently Ph. student degree in Engineering Physics & Mathematics Dept., Faculty of Engineering, Mansoura Univ., Egypt, PH+201005208927. Email: amoun1973@ya hoo.co m
• M. S. ElAzab is currently Prof. in Engineering Physics & Mathemat
ics Dept., Faculty of Engineering, Mansoura University, Egypt, PH
+201227379809. Email: ms_elazab@hotmail.com
• W. S. ElBeshbeshy is currently Dr. in the same Dept., Email:
welbeshbeshy@yahoo.com
tistep collocation method is also used for Volterra integral equations in [18]. Chebyshev spectral collocation method for the solution of Volterra integral and ordinary Volterra integro differential equations are discussed in [19]. In [20] Tang intro duces Legendrespectral method with its error analysis for ordinary Volterra integrodifferential equation of the second kind. Another spectral method using Legendre spectral Ga lerkin method was introduced for secondkind Volterra inte gral equations in [21]. Most papers that mentioned before
were devoted to VIEs and ordinary VIDEs, but in this article
we considered the PVIDEs.
Our goal in this article is to apply a Chebyshevcollocation
method for both space and time variables that are an extension
of the method presented in [22], [23].
The organization of this paper is as follows. In Section 2, we
present the Chebyshevcollocation spectral scheme for discre tizing the introduced problem. As a result a set of algebraic linear equations are formed and a solution of the considered problem is discussed. In Section 3, we present some nu
IJSER © 2014 http://www.ijser.org
International Journal of Scientific & Engineering Research, Volume 5, Issue 3, March2014 1349
ISSN 22295518
merical examples to demonstrate the effectiveness of the proposed method. Section 4 gives some concluding remarks.
we convert the integration interval from [0, T (τ + 1)] to [1, τ]
2
in (6), so that (6) will be
2 ¶U
¶ 2U τ
In this section we derive the Chebyshevcollocation method to
T ¶τ
 =
¶x 2
K (x , τ 
 1
ρ)DU (x , ρ)d ρ + g (x , τ),
(7)
problem (1)(3). Before this we will introduce some basic properties of the most commonly used set of orthogonal poly
where K(x, τ 
ρ) = k(x, T (τ 
2
ρ)). Define the auxiliary function
nomials; Chebyshev polynomials.
Φ(x , τ) = T
τ
K (x , τ 
ρ)U
(x ,ρ)d ρ + g (x , τ).
2 ò 1 xx
(8)
The Chebyshev polynomials {T n (x )}, n = 0, 1, ..., are the Eig en functions of the singular SturmLiouville problem
In order to approximate problem (1)(3) by spectral methods, we restart (7) as two equivalent Volterra Integrodifferential equations by using (8) as the following
d dT
(x ) n 2
ìï 2 ¶U
¶ 2U
( 1
x 2 n ) +
T n (x ) = 0,
ï  = Φ(x , τ),
dx dx
1 x 2
ïí T ¶τ ¶x
ï τ
(9)
they are also orthogonal with respect to L2 inner product on
ï Φ(x , τ) = T
K (x , τ 
ρ)U
(x ,ρ)d ρ + g (x , τ),
ïî ò xx
 1 2  1
the interval [1, 1] with the weight function ω(x) = (1
x2 ) 2
by integration of both sides of the first part in (9) over the in
1
ò T (x)T (x)ω(x)dx =
c j π δ ,
terval [1, τ], we get
 1 j k
2 jk
ìï T
ï U (x , τ) = u (x ) +
τ
[Φ(x , x) + U
(x , x)]d x,
where δ jk is the Kronecker delta, c0 = 2 and c j = 1 for j ³ 1. ï
0 2 ò 1 xx
(10)
The Chebyshev polynomials satisfy the following threeterm
ï Φ(x , τ) = T
τ
K (x , τ 
ρ)U
(x ,ρ)d ρ + g (x , τ).
recurrence relation
îï 2 ò 1 xx
Tn+ 1(x) = 2xTn (x)  Tn 1(x),
n ³ 1,
Let the collocation points be the set of (N+1)(M+1) points
T0 (x) = 1,
T1(x) = x,
{(x , τ )} in which {x =  cos( π i)}N
are the Chebyshev
and
i j i N i= 0
2Tn (x) =
1
n + 1
'
n+ 1
(x) 
1
n  1
'
n 1
(x),
n ³ 2,
(4)
GaussLobatto nodes (CGL nodes) and τ j are the Chebyshev
2 j + 1
T0 (x) = T1 (x),
2T1(x) = 0.5T2 (x).
Gauss nodes (CG nodes) defined as {τ j =  cos(
2M + 2
π)}M .
' '
A unique feature of the Chebyshev polynomials is their explic
Equation (10) holds at {(xi , τ j )}
it relation with a trigonometric function:
ìï
ï U (x i , τ j ) = u 0 (x i ) +
T τ j
ò
[Φ(x i , x) + U xx (x i , x)]d x,
Tn (x) = cos(n arccos x).
(5) ï
í
2  1
1 £ i £ N  1, 0 £
j £ M ,
(11)
ï Φ x τ = T
τ j
K x τ
 ρ U x
ρ d ρ + g x τ
Now, consider problem (1)(3) on domain Ω. Because of the
ï ( i , j )
ï
2 ò 1
( i , j
) xx (
i , ) ( i ,
j ),
orthogonally property of the Chebyshev polynomials on the ïî
interval [1, 1], we transfer (1) from [0, T] to an equivalent
0 £ i £ N , 0 £
j £ M ,
problem defined in [1, 1], by using the substitution
For approximating the integral terms in (11) the integral inter
val will transfer from [1, τj ] to a fixed one [1, 1] by using a
t = T (τ + 1),
2
τ Î [ 1,1],
simple linear transformation
τ + 1
τ  1
then (1) will convert to
ρ(τ j , θ) = x(τ j , θ) =
j θ + j ,
θ Î [ 1,1],
2 ¶U
¶ 2U
T ( τ+ 1) T 2 2
T ¶τ
 = 2
¶x 2 0
k (x , (τ + 1) 
2
s )Du (x ,s )ds + g (x , τ), (6)
where {θk }k = 0
are roots of the (p+1)th Chebyshev polynomials.
in which
U := U (x , τ) = u (x ,T (τ + 1)),
g (x , τ) = f (x ,T (τ + 1)),
Then (11) becomes:
ì U (x , τ ) = u (x ) +
ï +
2 2 ï
T (τ j
1) 1
[Φ(x , x(τ , θ)) + U
(x , x(τ , θ))]d θ,
for all x, τ Î Ω.
.
4  1
í
i j xx i j
(12)
By using the following change of variable
ï Φ(x i , τ j ) = g (x i , τ j ) +
ï
s = T (ρ + 1),
ρ Î [ 1,1],
ï T (τ j + 1)
1
K (x , τ
 ρ(τ , θ))U
(x ,ρ(τ , θ))d θ,
2 î . 4  1
i j j xx i j
Now, we use
IJSER © 2014 http://www.ijser.org
International Journal of Scientific & Engineering Research, Volume 5, Issue 3, March2014 1350
ISSN 22295518
N M
M (x ,σ ) := l
(x )F
(σ )φ (x , τ ),
ïì
U = u (x
T (τ + 1)
) + j
ΦN
∑∑ n m n m
n = 0m = 0
N M
(13)
ï i , j 0 i
ï 4
ï 1
U M (x ,σ ) := ∑∑l
(x )F
(σ )U (x , τ ),
ï .ò
I M [ΦM (x , x(τ , θ)) + (U M ) (x , x(τ , θ))]d θ,
N n m n m
n = 0m = 0
íï  1
ï
T (τ j
+ 1)
(16)
to approximate the functions Φ and U, where 𝑙𝑛 (𝑥) and
ï Φi , j = g i , j +
ï 4
𝐹𝑚 (𝜎) are the nth and mth Lagrange basis functions corre ï 1
ï M  M
sponding to nonuniform meshes of {xi } and {τj } respectively.
After enforcing the homogeneous boundary conditions at
ïî
where 𝐼𝑁
. I N [K (x i , τ j
 1
ρ(τ j , θ))(U N )xx (x i , ρ(τ j , θ))]d θ,
𝑥0 = −1 and 𝑥𝑁 = 1 the first and the last terms in the interpo
lation polynomial of U are omitted. Therefore, we have
𝑀 is the interpolation operator associated with the
Chebyshev mesh points ��𝑥𝑖 , 𝜏𝑗 ��, defined as the following;
N  1 M
N å å
n m n m
N  1 M
N å å
n m n m
I M Q :=
l (x )F
(σ)Q (x , τ ).
U M (x , σ) :=
n = 1 m = 0
l (x )F (σ)U (x
, τ ).
(14)
n = 1 m = 0
Now each equation in (16) can be reformulated respectively as
Now we can approximate 𝑈𝑥𝑥 in (12) by using the interpola tion polynomial of 𝑈𝑁 from the previous equation as follows
N  1 M
Φi , j = g i , j +
T (τ
+ 1) N  1 M
å å
n = 1 m = 0
l n (x i )U (x n , τm )
(17)
N xx
å å n m n m
. l (x )
K (x , τ
 ρ(τ , θ ))F
(ρ(τ , θ ))
F (θ)d θ,
(U M ) (x , σ) :=
n = 1 m = 0
l '' (x )F (σ)U (x
, τ ).
(15)
p  1 p
å z i å
i j j k m j k
1
 1 k
Where 𝑙′′ (𝑥) = 𝐷2
is the second derivative of the Lagrange
z = 1
k = 0
T (τ
+ 1) N M
interpolation function 𝑙𝑛 (𝑥) which is a polynomial of degree
N2, which can be defined as introduced in [24], [25], and [26]
the second derivative of the differentiation matrix 𝐷𝑁+1. Now,
U i , j = u 0 (x i ) +
N  1 M
. l '' (x
j [å å
4 n = 0 m = 0
)U (x , τ )]
l n (x i )Φ(x n , τm ) +
we write the entries of 𝐷2 = [𝐷2 ] for {𝑥 } as the following
å å
n = 1 m = 0
p  1
n i n m
p 1
ïì æ
÷ö
. l (x )
F (x(τ , θ ))
F (θ)d θ,
ï 2D D 2 = í
i ,k çD
i ,i 
1 ÷,
x i  x k ÷ø
i ¹ k ,
å z i
z = 1
å
k = 0
m j k
 1 k
1
i ,k
ï N
ï 2
ï i ,l
i = k ,
Now we discuss an efficient way to find dk
= ∫−1
Fk (θ)dθ. First
îï l = 0 ,i¹ l
we express 𝐹𝑗 (𝑠) in terms of Chebyshev functions as in [20]:
N
where 𝐷𝑖,𝑘 the entries of the socalled differentiation matrix,
F (s) = ωC ∑ (T (x ) / γ
)T (s),
(19)
which has a dimension of (N + 1). The entries of the differenti
ation matrix can be defined in [27] for (CG) points as the fol
where ωC
j j p j p p
p = 0
lowing
ïì 2N 2 + 1
j is the Chebyshev weight corresponding to Cheby
shev points {x }N and
N
ï  ,
i = k = 0,
π
γ = ∑T 2 (x )ωC =
p = 0,
(20)
6
ï i + k
p p i i
i =1
π / 2,
1 ≤ p < N ,
c ( 1)
ï i ,
i ¹ k ,
and γ
= π / 2 if {x }N
is the Chebyshev Gauss or the Cheby
Di ,k
= ï c k
x i  x k
x
shev Gauss Radau points, where if we use {x }N as the Cheby
ï

ï 2 1 ï
i ,
x 2 )
1 £ i = k £ N  1,
shev Gauss Lobatto then γ N
1
= π. From (19) we can now calcu
ï 2N 2 + 1
ï ,
i = k = N ,
late dk = ∫−1
Fk (θ)dθ as following
îï 6
N 1
d = ωC ∑ (T (x ) / γ )∫
T (θ)dθ.
with 𝑐𝑖 = 2 for i = 0, N, and 𝑐𝑖 = 1 otherwise. Now we approx
k k p k p
p = 0
−1 p
imate the integration in (12) by replacing the integral function
by its interpolation polynomial approximation of Φ and ap
proximation of 𝑈𝑥𝑥 from (13) and (15) respectively in (12), and
writing, Φ(x i , τ j ) = Φi , j ,U (x i , τ j ) = U i , j , g (x i , τ j ) = g i , j . Then
1
To compute ∫−1Tp (θ)dθ, we use the recurrence relation (4) for
Chebyshev polynomials yields
2
1 , p is even number,
our goal is to find U i , j so we obtain the following
∫−1
T (θ)dθ = 1 − p2
0, othrwise.
Now rewrite (17) and (18) in the matrix form as
IJSER © 2014 http://www.ijser.org
International Journal of Scientific & Engineering Research, Volume 5, Issue 3, March2014 1351
ISSN 22295518
Φ( N +1)(M +1) = G( N +1)( M + ) + +U( N −1)( M +1) ,
To solve the coupled equations system in (17) and (18), we
U
( N Φ−1)( M +1)
= U−1
+UA
( N +1)(M +1) + B ( N −1)( M +1) ,
convert them to a linear algebraic system as follows
(I − B − A+)U( N −1)( M +1) = U−1 + AG( N +1)( M +1) .
the following
After solving the previous system, we obtain an approxima tion to U( N −1)( M +1) , then the approximation to the original prob
( N +1)(M +1) ,
G( N +1)( M + ) = vec g , ,
0 ≤ i ≤ N ,
i j
i j
lem for all x ∈ (1, 1), and t ∈ [0, T] can found by
U( N −1)( M +1) = vec U , ,
1 ≤ i ≤ N − 1,
0 ≤ j ≤ M ,
N −1 M 2
1
∑∑ n m n , m
i j
u0 (x2 )
u0 (x3 ) 3
u0 (xN −1 )
u(x,t) ≈
n =1 m = 0
l (x)F (
T
t − )U .
U = vec
−1 ,
u0 (x2 )
u0 (x3 ) 3
u0 (xN −1 )
in which the vec operator reshapes any matrix into a vector by placing columns of the matrix below each other from the first to the last. For the other matrices each one can defined as block ones as the following:
In order to test the utility of the proposed new method, we devoted this section to some numerical examples to view the efficiency and accuracy of the method in the previous one. In
our implementation , we set T = 1, p = N = M, and let {θ }p be
• A = ( Ai ), A is a matrix with dimension of
the ChebyshevGauss points with the corresponding weights
(N + 1)(N − 1) × (N + 1)2 in which the first and last (N+1)
𝑖
ωC =
π / (N
+ 1).
To show the efficiency of the previous method
columns are zeros and the other blocks �𝐴𝑗 � forms a di
agonal matrix in which its entries are given by
p
for our problems in comparison with the exact solution, we
calculate for different values of N the maximum error defined
by
( Ai )
T(τ + 1)
= j ∑
d F (ξ(τ ,θ )),
i ,i +1
4 k = 0 ax
T(τ
,
+ 1)
.
j k m j k
E
= m U
− u x j
∞ i N
with (0 ≤ m ≤ M, 0 ≤ j ≤ M, 1 ≤ i ≤ N − 1, and 1 ≤ n ≤ N – 1).
1≤ ≤ −1
i , j i 2
0≤ j ≤ M
𝑗
For each block in (𝐴𝑖 ) we obtain a matrix with dimension
of (N + 1) × (N + 1). The shape of the global matrix A will
be
Example 3.1: Consider the linear problem (1)(3) with the ker nel k(x, t) = exp (x2t) and we choose the forcing function f(x, t) so that u(x, t) = (1 x2) exp (t) is the exact solution.
1 0 3 0
j
A =
0 3 0
N −1
j
• + = (+i ), + is a matrix with dimension
(N + 1)2 × (N + 1)(N − 1) in which the first and last (N+1)
rows are zeros, and each block matrix (𝐿𝑖 ) has a dimen sion (N + 1) × (N + 1). The entries of (𝐿𝑖 ) can give by the
following
( i )
T(τ + 1) p
= j ( ,
− ( , )) ( ( , )) 2 ,
+
i +1,n
∑dk K xi τ j
k = 0
ρ τ θ
Fm ρ τ j θk
Di ,n
with (0 ≤ m ≤ M, 0 ≤ j ≤ M, 1 ≤ i ≤ N − 1, and 1 ≤ n ≤ N – 1)
respectively.
• Finally the matrix B = (Bi ), B has the dimension
(N + 1)(N − 1) × (N + 1)(N − 1). The entries of matrices (𝐵𝑖 ) is
obtained as the following
(Bi )
i ,n
T(τ j + 1) p
= ∑dk Fm (ξ(τ j ,θk ))
2
i ,n
4 k = 0
with (0 ≤ m ≤ M, 0 ≤ j ≤ M, 1 ≤ i ≤ N − 1, and 1 ≤ n ≤ N – 1) respectively. Each matrix in the previous equation of di mension (N + 1) × (N + 1).
For each previous matrix we can use kronecker product to find them well.
Fig. 1. Approximated and Exact solutions respectively for x
∈ (1, 1) and t ∈ [0, 1] at N = 12.
IJSER © 2014 http://www.ijser.org
International Journal of Scientific & Engineering Research, Volume 5, Issue 3, March2014 1352
ISSN 22295518
TABLE 1
MAXIMUM L∞ ERRORS FOR EXAMPLE (3.1)
CPU time(s) indicates the time for all calculations for the program from start to the end.
Example 3.2: Consider the linear problem (1)(3) with the ker nel k(x, t) = exp (x2t) and we choose the forcing function f(x, t) so that u(x, t) = (1 x4) exp (x+t) is the exact solution.
TABLE 2
MAXIMUM L∞ ERRORS FOR EXAMPLE (3.2)
CPU time(s) indicates the time for all calculations for the program from start to the end.
Fig. 2. Approximated and Exact solutions respectively for x
∈ (1, 1) and t ∈ [0, 1] at N = 16.
It can be seen that the errors in Tables (1)(2) decay rapidly, which is confirmed by spectral accuracy.
All the computations are carried out in double precision arithmetic using Matlab 7.9.0 (R2009b). To obtain sufficient accurate calculations, variable arithmetic precision (vpa) is employed with digit being assigned to be 32. The code was executed on a second generation Intel Core i52410M, 2.3 Ghz Laptop.
This article used a new competitive numerical scheme based on developing Chebyshev spectral collocation method to find the approximate solution of parabolic Volterra integro differential equations. The main advantage of using Cheby shev scheme instead of using Legendre one is that its quadra ture points have explicit and simple expressions as well as the corresponding weights. This enables us to avoid the complex computation of the Legendre quadrature points and the corre sponding weights. Moreover, we made a minor modification on approximating the integration by replacing the integral function by its interpolating polynomials instead of using Gauss quadrature approximation and this increases the accu racy of the suggested method. The numerical examples given in this work have demonstrated the potential of the newly proposed numerical scheme in solving parabolic Volterra in tegrodifferential and similar equations even with using a small number of collocation points.
[1] E. G. Yanik, and G. Fairweather, “Finite element methods for parabolic and hyperbolic partial integrodifferential equations”, Nonlinear Anal. (1988), 12, pp. 785–809.
[2] M. Renardy, W. Hrusa, and J. Nohel, “Mathematical Problems in Viscoelasticity”, Pitman monographs and Surveys in pure Appl. Math., New York: wiley. (1987), 35.
[3] T. A. Burton, “Volterra Integral and Differential Equations“, (second edition), Mathematics in Science and Engineering. Elsevier, vol. 202,
2005.
[4] G. Dagan, “The significance of heterogeneity of evolving scales to transport in porous formations“, Water Resour. Res., (1994), 30, pp.
33273337.
[5] J. Cushman, and T. Glinn, “Nonlocal dispersion in media with con tinuously evolving scales of heterogeneity“, Transport in Porous Media, (1993), 13, pp. 123138.
IJSER © 2014 http://www.ijser.org
International Journal of Scientific & Engineering Research, Volume 5, Issue 3, March2014 1353
ISSN 22295518
[6] R. Ewing, Y. Lin, and J. Wang, “A numerical approximation of non Fickian flows with mixing length growth in porous media“, Acta. Math. Univ. comenian. (2001), LXX, pp. 7584.
[7] Amiya K. Pani, and G. Fairweather, “H1Galerkin mixed finite element methods for parabolic partial integrodifferential equations”, IMA Journal of Numerical Analysis. (2002), 22, pp. 231252.
[8] Amiya K. Pani, G. Fairweather, and Ryan I. Fernandes, “ALTERNATING DIRECTION IMPLICIT ORTHOGONAL SPLINE COLLOCATION METHODS FOR AN EVOLUTION EQUATION WITH A POSITIVE TYPE MEMORY TERM”, SIAM J. NUMER. ANAL. (2008), 46, pp.
344364.
[9] Amiya K. Pani, G. Fairweather, and Ryan I. Fernandes, “ADI orthogonal spline collocation methods for parabolic partial integrodifferential equa tions”, IMA Journal of Numerical Analysis. (2010), 30, pp. 248276.
[10] F. FakharIzadi, and M. Dehghan, “The spectral methods for parabol ic Volterra integrodifferential equations”, Journal of Computational and Applied Mathematics. (2011), 235, pp. 40324046.
[11] C. Canuto, M. Y. Hussaini, A. Quarteroni, and T. A. Zang, “Spectral Methods: Fundamentals in Single Domains”, SpringerVerlag, 2006.
[12] L. N. Trefethen, “Spectral Methods in MAT+AB”, Society for Industrial and Applied Mathematics (SIAM), Philadelphia, PA, 2000.
[13] C. Canuto, M. Y. Hussaini, A. Quarteroni, and T. A. Zang, “Spectral
Methods in Fluid Dynamics”, SpringerVerlag, 1988.
[14] W. Guo, G. Labrosse, and R. Narayanan, “The Application of the Cheby shevSpectral Method in Transport Phenomena”, Notes in Applied and Computational Mechanics 68, SpringerVerlag Berlin Heidelberg,
2012.
[15] J. P. Boyd, “Chebyshev and Fourier Spectral Methods”, Dover Publica tions, Mineola, NY, 2001.
[16] J. Shen, and T. Tang, “Spectral and HighOrder Methods with applica tions”, Science Press, Beijing, China. 2006.
[17] M. Aguilar, and H. Brunner,” Collocation methods for secondorder
Volterra integrodifferential equations”, Appl. Numer. Math. (1988),
4, pp. 455470.
[18] D. Conte, and B. Paternoster, “Multistep collocation methods for Volterra Integral Equations”, Applied Numerical Mathematics. (2009), 59, pp. 17211736.
[19] Tobin A. Driscoll, “Automatic spectral collocation for integral, in tegrodifferential, and integrally reformulated differential equa tions”, Journal of Computational Physics. (2010), 229, pp. 59805998.
[20] T. Tang, “On spectral methods for Volterra integral equations and the convergence analysis”, J. Comput. Math. (2008), 26, pp. 825837.
[21] Z. Wan, Y. Chen, and Y. Huang, “Legendre spectral Galerkin method for secondkind Volterra integral equations”, Front. Math. China. (2009), 4, pp. 181193.
[22] Zhendong Cu, and Yanping Chen, “Chebyshev spectralcollocation method for Volterra integral equations”, Recent Advances in Scien tific Computing and Applications, Contemporary Mathematics. (2013), 586, pp. 163170.
[23] Wu Hua, and Zhang Jue, “ChebyshevCollocation Spectral Method for Volterra IntegroDifferential equations”, Journal of Shanghhai University (Natural Scince). (2011), 17, pp. 182188.
[24] B. Costa, and W. S. Don, “On the computation of high order pseudo
spectral derivatives”, Appl. Numer. Math. (2000), 33, pp. 151159. [25] R. Baltensperger, and M. R. Trummer, “Spectral Differencing with a
twist”, J. Indus. and Appl. Math. (2003), 24, pp. 14651487.
[26] Elsayed M. E. Elbarbary, and Salah M. ElSayed, “Higher order pseudospectral differentiation matrices”, Appl. Numer. Math. (2005),
55, pp. 425438.
[27] J. S. Hesthaven, S. Gottlieb, and D. Gottlieb, “Spectral Methods for
TimeDependent Problems”, Cambridge University Press. 21, pp. 96
99, 2007.
IJSER © 2014 http://www.ijser.org