A novel semi-analytical approach for solving nonlinear Volterra integro-differential equations

Applied Mathematics and Computation 263 (2015) 25–35

Contents lists available at ScienceDirect

Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc

A novel semi-analytical approach for solving nonlinear Volterra integro-differential equations Kyunghoon Kim, Bongsoo Jang∗ Department of Mathematical Sciences, Ulsan National Institute of Science and Technology (UNIST), Ulsan, 689-798, Republic of Korea

a r t i c l e

i n f o

a b s t r a c t

Keywords: Volterra integro-differential equation Taylor series Differential transform method Gronwall inequality

In this work, we present an efficient semi-analytical method based on the Taylor series for solving nonlinear Volterra integro-differential equations, namely the differential transform method (DTM). The DTM provides a recursive relation for the coefficients of the Taylor series that is derived from the given equations. We provide a new recursive relation for the nonlinear Volterra integro-differential equations with complex nonlinear kernels. Since the DTM is based on the Taylor series, it is difficult to obtain accurate approximate solutions in a large domain. To overcome this difficulty, the standard DTM is applied in each subdomain, called the multistage differential transform method (MsDTM). We also present an convergence analysis for the proposed method. To demonstrate the efficiency of the proposed method, several numerical examples are performed and support the results in our analysis. © 2015 Elsevier Inc. All rights reserved.

1. Introduction This paper concerns the following nonlinear Volterra integro-differential equation:

u (x) = f (x) + λ



x 0

k(x, t)g(u(t)) dt,

u(0) = u0 ,

(1)

where the kernel k(x, t), f(x) and g(u(x) are sufficiently smooth functions. The nonlinear integro-differential equations play an important role to describe many phenomena arising in science and engineering fields such as elasticity, population dynamics, fluid mechanics and so on [1]. There exist many numerical methods for solving integral equations such as collocation method [2], Taylor polynomial [3], Tau method [4], wavelet-Galekin method [5], Adomain decomposition method [6], Homotopy analysis method [7], Homotopy analysis transform method [8–11], Homotopy perturbation method [12], Sinc collocation [13], etc. In the present work, we propose an efficient numerical method for solving the nonlinear Volterra integro-differential equation by using the Taylor series expansion which is called the differential transform method (DTM) [14–18]. In DTM it is the key to derive a recursive relation for the coefficient of the Taylor series from the given equation. Since the DTM is based on the Taylor series, the convergence can be taken into account in radius of convergence. The main feature of the proposed method is that we apply the DTM in each sub domain and construct the recursive relation of the coefficient of the Taylor series. The article is organized as follows. In Section 2, a brief description of the DTM and its several properties are presented. In Section 3, we give the derivation of the recursive relation. The convergence analysis of the proposed method is described in Section 4. Several numerical results are demonstrated and all numerical approximations support the theoretical results in Section 5. Finally a conclusion is presented in Section 6. ∗

Corresponding author. Tel.: +82 52 217 3136; fax: +82 52 217 3101. E-mail address: [email protected] (B. Jang).

http://dx.doi.org/10.1016/j.amc.2015.04.011 0096-3003/© 2015 Elsevier Inc. All rights reserved.

26

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35 Table 1 Fundamental operations for the one dimensional DTM. Original function w(x)

Differential transform W(k)

u(x) ± v(x) xr u(x)v(x)

U(k) ± V(k) δ (k − r) k r=0 U (r )V (k − r ) (k + 1)U(k + 1)

du(x) dx

U (k−1) k

u(x)dx

2. Description of differential transform method Let us describe the basic definition and some properties of the standard DTM. The differential transform U(m) of a given function u(x) at x = xi is defined as

U (m) =

  1 dm u(x) . m! dxm x=xi

Correspondingly, the inverse differential transform u(x) of U(m) at x = xi is defined by

u(x) =

∞ 

U (m)(x − xi )m .

m=0

Several fundamental properties of the differential transform are listed in Table 1. In DTM it is most important to derive a recursive relation of the differential transforms from the given equation. Since the nonlinear function g(x) is included in the integral form, we firstly discuss to find the differential transform G(m) of g(x). The followings are recursive algorithms to calculate the differential transforms of several nonlinear functions [19].  ∞ m m Theorem 2.1. Let g(u(t)) = eau(t) , where a is a constant. Suppose that g(u(t)) = ∞ m=0 G(m)t and u(t) = m=0 U (m)t . Then we have

⎧ ⎪ eaU(m), m = 0, ⎪ ⎨ m−1  G(m) = a ⎪ (m − r)G(r)U(m − r), m ≥ 1. ⎪ ⎩m r=0

Theorem 2.2. Let g(u(t)) = ln (a + bu(t)), where a and b are constants. Suppose that g(u(t)) = Then we have

⎧ ln(a + bU (m)), ⎪ ⎪ ⎪ ⎪ ⎪ b ⎪ ⎨ U (m), G(m) = a + bU (0) ⎪ m−2 ⎪  r+1 ⎪ b ⎪ ⎪ G(r + 1)U (m − 1 − r)], [U ( m ) − ⎪ ⎩ a + bU (0) m

∞

m=0

G(m)tm and u(t) =

∞

m=0

U (m)tm .

m = 0, m = 1,

(2)

m ≥ 1,

r=0

Theorem 2.3. Let us consider g(u(t)) = sin (au(t)) and h(u(t)) = cos (au(t)), where a is a constant. Suppose that g(u(t)) =  ∞ m m h(u(t)) = ∞ m=0 H(m)t and u(t) = m=0 U (m)t . Then we have

⎧ sin(aU (m)), ⎪ ⎪ ⎨  G(m) = a m−1 ⎪ (m − r)H(r)U(m − r), ⎪ ⎩m

∞

m=0

G(m)tm ,

m = 0, m ≥ 1,

(3)

r=0

⎧ cos(aU (m)), ⎪ ⎪ ⎨ m−1 H(m) = a  ⎪ − (m − r)G(r)U(m − r), ⎪ ⎩ m

m = 0, m ≥ 1.

(4)

r=0

3. Differential transform method for the nonlinear Volterra integro-differential equation In this section, we describe how to apply the DTM to the nonlinear Volterra integro-differential equation in (1). To do this let us consider the model problem (1) as the following simplified form

u (x) = f (x) + λv(x),

(5)

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35

where

v(x) =



x 0

k(x, t)g(u(t)) dt.

Suppose that u(x) =

∞

m=0

27

(6)

U (m)xm . Applying the basic properties of DTM gives the following recursive relation

(m + 1)U(m + 1) = F (m) + λV (m), U(0) = u0 , m = 0, 1, 2, . . ..

(7)

where F(m) and V(m) are the differential transforms of f(x) and v(x), respectively. To solve the recursive Eq. (7) it is essential to find the differential transforms F(m) and V(m). It is easy to obtain F(k) from the given function f(x). Thus, it is a key to find the differential transform V(m) in solving the recursive relation (7). Now, let us obtain the differential transform V(m) of the Volterra term v(x). To do this we consider the Taylor series of k(x, t) at x = 0

k(x, t) =

∞ 

Kˆt (j)xj .

j=0

Then we have

v(x) =



x 0

k(x, t)g(u(t)) dt =

∞ 



x

xj 0

j=0

Kˆt (j)g(u(t)) dt.

(8)

It is noted that the differential transform Kˆt (j) is a function of t. Thus ,we can consider the Taylor series of Kˆt (j) and g(u(t)) at t=0

Kˆt (j) =

∞ 

K (j, p)tp ,

g(u(t)) =

p=0

∞ 

G(p)tp .

p=0

Applying the fundamental operation of multiplication in Table 1, we have

v(x) =

∞ 

 x

j

0 p=0

j=0

=

∞ x

∞  ∞  j=0 p=0



p 



K (j, r)G(p − r) tp dt

r=0

 p  1 K (j, r)G(p − r) xj+p+1 . p+1

(9)

r=0

Expanding Eq. (9) with each term separate, we can find the Taylor of v(x). Set m = j + p + 1. Since j  0 and p  0, the index m is m  1. Thus, the Taylor series of v(x) can be written by

v(x) =

∞ 

V (m)xm ,

m=1

where

V (m) =

m−1  j=0

⎡ ⎤ m−j−1 1 ⎣  K (j, r)G(m − j − r − 1)⎦. m−j r=0

Combined with (7) we have the following result. Theorem 3.1. The recursive relation of differential transform U(m), m = 0, 1, . . ., of the nonlinear Volterra integro-differential Eq. (1) is written by

U (m + 1) =

where

V (m) =

m−1  j=0

⎧ 1 ⎪ ⎪ ⎨ m + 1 F (m), ⎪ ⎪ ⎩

m = 0,

1 {F (m) + λV (m)}, m ≥ 1, m+1

(10)

⎡ ⎤ m−j−1 1 ⎣  K (j, r)G(m − j − r − 1)⎦. m−j r=0

One can obtain the U(m) by solving the recursive relation (10). Then u(x) can be approximated by the finite-term Taylor series. In other words,

u(x) ≈

n  m=0

U (m)tk ≡ sn (t).

28

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35

It is well known that the approximated error is given by

|u(x) − sn (x)| ≤

M

|x|n+1 ,

(n + 1)!

where |u(n + 1) (x)|  M. Thus, it is easy to see that, for the fixed term n, the error is getting bigger if |x| is increasing. To reduce the error, the more terms are added. But it requires cumbersome computational work. Moreover, it is possible that the error blows up even with many terms if the value x is the outside radius of convergence for the Taylor series of u(x). To overcome this difficulty arising in the standard DTM, the multistage differential transform method(MsDTM) is introduced [17,18]. The basic idea of the MsDTM is to apply the standard DTM to each sub-domain. To describe the MsDTM we consider the equally spaced nodal points xl ; 0 = x0 < x1 <  < xN − 1 < xN = T and xl + 1 − xl = T/N. On lth sub-domain [xl , xl + 1 ]  Dl , we define u(x)|Dl ≡ ul (x) and f (x)|Dl = fl (x). The differential transform Ul (m) of ul (x) at x = xl is defined by

  1 dm ul (x) . m! dxm x=xl

Ul (m) =

The inverse differential transform of Ul (m) is defined by ∞ 

ul (x) =

Ul (m)(x − xl )m .

(11)

m=0

Remark. In the standard MsDTM, the differential transform Ul (m) can be obtained by solving the same recursive relation with (10). That is,

(m + 1)Ul (m + 1) = Fl (m) + λVl (m), m = 0, 1, 2, · · · .

(12)

where Fl (m) and Vl (m) are the differential transforms of fl (x) and vl (x) on Dl , respectively. However, it is not enough to find the correct Ul (m) by solving the above recursive relation because the Volterra term v(x) is the global operator. In what follow we describe how to derive the correct recursive relation of the differential transforms in MsDTM that consists of the characteristic in the global operator. For x  Dl , l > 0, the model problem (1) can be rewritten by

ul (x) = fl (x) + λ

l−1   i=0

xi+1 xi

k(x, t)g(ui (t)) dt + λ



x xl

k(x, t)g(ul (t)) dt

(13)

ni Suppose that ui (x), i = 0, 1, , l − 1, are approximated by ui (x) ≈ si,ni (x), where si,ni (x) = m=0 Ui (m)(x − xi )m . To find the differential transform Ul (m) on Dl which satisfies Eq. (13), let us consider the Taylor series expansion of k(x, t) at x = xl in the first integral term as ∞ 

k(x, t) =

Kˆt,l (j)(x − xl )j .

(14)

j=0

Then we have



xi+1

xi

k(x, t)g(ui (t)) dt =

∞ 

Vˆi (j)(x − xl )j ,

(15)

j=0

where

Vˆi (j) =



xi+1 xi

Kˆt,l (j)g(ui (t)) dt.

The second integral term of (13) is exactly same if the interval (0, x) in (8) is substituted by (xl , x). Thus, if we consider the Taylor series expansions of g(u(t)) and Kˆt,l (j) at t = xl as

Kˆt,l (j) =

∞  p=0

Kl (j, p)(t − xl )p ,

g(u(t)) =

∞ 

Gl (p)(t − xl )p ,

p=0

where Kˆt,l (j) is obtained by (14), then Theorem 3.1 with (15) and (16) gives the explicit formulae of Ul (m) as follows

(16)

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35

29

Theorem 3.2. The differential transforms Ul (m), l > 0, m = 0, 1, , can be formulated by

⎫ ⎧ ⎧ l−1 ⎬ ⎨ ⎪  ⎪ 1 ⎪ ⎪ Vˆi (m) , Fl (m) + λ ⎪ ⎪ ⎭ ⎨m + 1 ⎩ i=0 ⎫ ⎧ Ul (m + 1) = ⎪ l−1 ⎬ ⎨ ⎪  1 ⎪ ⎪ ⎪ Vˆi (m) + λVl (m) , Fl (m) + λ ⎪ ⎩m + 1 ⎩ ⎭

m = 0, (17)

m ≥ 1,

i=0

where

Vˆi (m) = Vl (m) =



xi+1

xi m−1  j=0

Kˆt,l (j)g(ui (t)) dt.

⎤ ⎡ m−j−1 1 ⎣  Kl (j, r)Gl (m − j − r − 1)⎦. m−j

(18)

r=0

Remark. In order to compute the value of Vˆi (m) in (18) it is necessary to have the exact value of ui (t). However, the MsDTM gives the ni th partial sum of the Taylor series of ui (t). Thus, we can approximate the value of Vˆi (m) by substituting si,ni (t) into ui (t) with an appropriate numerical scheme. That is,

Vˆi (m) ≈



xi+1 xi

Kˆt,l (j)g(si,ni (t)) dt.

(19)

Once the differential transforms Ul (m) are obtained, the solution ul (x) can be approximated as follows nl 

ul (x) ≈

Ul (m)(x − xl )m ≡ sl,nl (x).

m=0

Thus, the solution u(t) can be approximated by

u(x) ≈

N−1 

sl,nl (x)χl (x),

(20)

l=0

where

χl (x) =

 1, x ∈ Dl 0, otherwise

4. Analysis of convergence In this section, we discuss about the convergence of the approximated solution obtained by (20). In order to analyze the convergence we assume that the solution u(x), the kernel function k(x, t) and the nonlinear function g(u(x)) of (1) in the domain D = [0, 1] satisfy the following conditions: i. u(n + 1) (x)  C[0, 1] with ||u(n + 1) ||  μ. ii. |k(x, t)|  M for all 0  x, t  1. iii. g(u(x)) is satisfied in Lipschitz condition such that

|g(u(x)) − g(v(x))| ≤ L|u(x) − v(x)|.

(21)

We define the equally spaced nodal points xi : 0 = x0 < x1 <  < xN = 1 with the length h = 1/N and define the subintervals Dl = [xl , xl + 1 ]. Let us remind that the solution ul (x) in Dl satisfies the following equation

ul (x) = fl (x) + λ

l−1   i=0

xi+1

xi

k(x, t)g(ui (t))dt + λ



x xl

k(x, t)g(ul (t)) dt.

(22)

As stated in (19), the approximate solution by the proposed MsDTM in Dl is the Taylor series of uˆ l (t) that solves the following problem

uˆ l (x) = fl (x) + λ

l−1   i=0

xi+1

xi

k(x, t)g(sˆi,ni (t)) dt + λ



x xl

k(x, t)g(uˆ l (t)) dt,

(23)

where sˆi,ni (x) is the ni th partial sum of the Taylor series of uˆ i (x) and the initial condition is uˆ l (xl ) = sˆl−1,nl−1 (xl ). To simplify the analysis we assume that the number of terms for each partial sum sˆi,ni (x) by ni = n in Di = (xi , xi + 1 ). The following two theorems will be used to analyze the convergence [20].

30

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35

Theorem 4.1 (Gronwall’s inequality). Let u(x), v(x) and w(x) be continuous and nonnegative functions on [a, b], and let w(x) be a nondecreasing. If

u(x) ≤ w(x) +



x a

v(t)u(t) dt,

x ∈ [a, b],

then u(x) satisfies

u(x) ≤ w(x) exp



x a

 v(t) dt ,

x ∈ [a, b].

Theorem 4.2 (Discrete Gronwall’s inequality). Let ε i  0 and ηi  0, i = 0, 1, , n and a  ε 0 . If

n ≤ a +

n−1 

ηi i , n ≥ 1,

i=0

then ε n satisfies

n ≤ a exp

n−1 



ηi , n ≥ 1.

i=0

Theorem 4.3. Let ul (x) and uˆ l (x) be the solutions of (22) and (23), respectively. Assume that ul (x), uˆ l (x), k(x, t) and the nonlinear function g(x) satisfy the assumptions (21). Then we have

||ul − uˆ l ||L∞ (Dl ) → 0 as h → 0. Moreover,

||ul − uˆ l ||L∞ (Dl ) ≤ c

hn+1

(n + 1)!

,

where c is a constant and independent of h. Proof. Since uˆ l (x) is the solution of (23) with the initial condition sˆl−1,n (xl ), we have



uˆ l (x) = sˆl−1,n (xl ) +

x xl

f (s) ds + λ



⎡

x

⎣

xl

⎤

l−1   i=0

xi+1 xi

k(s, t)g(sˆi,n (t)) dt⎦ ds + λ

From (13), ul (x) also satisfies the following:

ul (x) = ul (xl ) +



x

f (s) ds + λ

xl



x xl

⎡ l−1   ⎣ i=0

⎤ xi+1 xi

k(s, t)g(ui (t)) dt⎦ ds + λ



x xl





x



xl

s xl

s xl

k(s, t)g(uˆ l (t)) dtds.

k(s, t)g(ul (t)) dtds.

Let us define the error el (x) = ul (x) − uˆ l (x). Using the assumption (21) we have

|el (x)| ≤ |ul (xl ) − sˆl−1,n (xl )| + λML



l−1  xi+1 x

xl

i=0

xi

|ui (t) − sˆi,n (t)| dt + λML

≤ |ul (xl ) − sˆl−1,n (xl )| + λh2 ML

l−1 



x xl



s

xl

|el (t)| dtds

||ui − sˆi,n ||L∞ (Di ) + λhML



i=0

x xl

|el (t)| dt

Since

vi (t) = sˆi,n (t) +

v(i n+1)(ξ ) (t − xi )n+1 , for some ξ ∈ (xi , t), (n + 1)!

we have

||ui − sˆi,n ||L∞ (Dl ) ≤ ||ui − vi ||L∞ (Dl ) + ||vi − sˆi,n ||L∞ (Dl ) ≤ ||ei ||L∞ (Dl ) +

μ hn+1 . (n + 1)!

(24)

Since ul (xl ) = ul − 1 (xl ), we have

|ul (xl ) − sˆl−1,n (xl )| ≤ ||ul−1 − vl ||L∞ (Dl−1 ) + ||vl−1 − sˆl−1,n ||L∞ (Dl−1 ) μ hn+1 . ≤ ||el−1 ||L∞ (Dl−1 ) + (n + 1)!

(25)

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35

31

Combined with (24) and (25) we have

|el (x)| ≤ (μ + λh2 μlML)

 hn+1 + (λh2 ML + δi,l−1 )||ei ||L∞ (Di ) + λhML (n + 1)! l−1

i=0



x xl

|el (t)| dt,

where δ i, j = 1 if i = j and δ i, j = 0 if i j. Now, applying the Gronwall’s inequality in Theorem 4.1, we get

⎫ l−1 ⎬  hn+1 2 + ||el ||L∞ (Dl ) ≤ (μ + α h μl) (α h + δi,l−1 )||ei ||L∞ (Di ) exp (α h2 ), ⎭ ⎩ (n + 1)! ⎧ ⎨

2

i=0

where α = λML. The discrete Gronwall’s inequality in Theorem 4.2 gives

||el ||L∞ (Dl ) ≤ [(μ + α h2 μl) exp (α h2 ) exp {(α h2 l + 1) exp (α h2 )}]

hn+1 . (n + 1)!

Then it is clear that ||el ||L∞ (Dl ) → 0 as h → 0. Moreover, since l < N = 1/h, lh < 1. The proof is completed by taking the constant c as

c = (μ + αμ) exp (α) exp {(α + 1) exp (α)}.



Theorem 4.4. Let ul (x) and uˆ l (x) be the solutions of (22) and (23), respectively. Assume that ul (x), uˆ l (x), k(x, t) and the nonlinear function  g(x) satisfy the assumptions (21), and that sˆl,n (x) is the approximation of uˆ l (x) by using the MsDTM. Let uˆ h (x) = N−1 sˆ (x)χl (x). l=0 l,n Then we have

||u − uˆ h ||L∞ (D) → 0 as h → 0. Moreover,

||u − uˆ h ||L∞ (D) ≤ c

hn+1 , (n + 1)!

where c is a constant and independent of h. Proof. Since the sˆl,n (x) is the nth partial sum of the Taylor series of uˆ l (x),

||uˆ l − sˆl,n ||L∞ (Dl ) ≤ μ

hn+1 , (n + 1)!

l = 0, 1, . . ., N − 1.

(26)

For l = 0, u0 (x) = uˆ 0 (x). Then we have

||u0 − sˆ0,n ||L∞ (D0 ) = ||uˆ 0 − sˆ0,n ||L∞ (D0 ) ≤ μ

hn+1

(n + 1)!

.

(27)

For l  1, Theorem 4.3 with (26) gives

||ul − sˆl,n ||L∞ (Dl ) ≤ ||ul − uˆ l ||L∞ (Dl ) + ||uˆ l − sˆl,n ||L∞ (Dl ) ≤ (d + μ)

hn+1 , (n + 1)!

(28)

where d is a constant and independent of h. From (27) and (28), we have

||u − uˆ h ||L∞ (D) = max ||ul − sˆl,n ||L∞ (Dl ) ≤ (d + μ) 0≤l≤N−1

hn+1 . (n + 1)!

It completes the proof.  5. Numerical illustrations In this section, we show that the proposed scheme(MsDTM) is very efficient and accurate in obtaining the approximate solution for the model problem (1) with various nonlinear functions g(u(t)) such as power, exponential, logarithmic and trigonometric functions. To test the proposed scheme, we define the equally spaced nodal points xi : 0 = x0 < x1 <  < xN = 1 with the length h = 1/N. We define the approximation uˆ h (x) by

uˆ h,n (x) =

N−1 

sˆl,n (x)χl (x),

l=0

where sˆl,n is the nth partial sum of the Taylor series of uˆ l (x) that is the solution of (23). Here, sˆl,n (x) can be obtained by using the Theorem 3.2.

32

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35 Table 2 a , n = 1, 3, 5 in Example 1. Maximum errors and convergence rates for Eh,n h

a ||Eh,1 ||∞

rate

a ||Eh,3 ||∞

rate

a ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

2.375(−1) 1.298(−1) 4.911(−2) 1.483(−2) 4.021(−3) 1.042(−3)

0.872 1.402 1.728 1.883 1.948

2.177(−2) 3.405(−3) 3.101(−4) 2.300(−5) 1.562(−6) 1.017(−7)

2.676 3.457 3.753 3.880 3.941

1.099(−3) 4.251(−5) 9.851(−7) 1.861(−8) 3.193(−10) 5.216(−12)

4.693 5.431 5.726 5.864 5.935

Table 3 b , n = 1, 3, 5 in Example 1. Maximum errors and convergence rates for Eh,n h

b ||Eh,1 ||∞

rate

b ||Eh,3 ||∞

rate

b ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

1.832(−1) 7.328(−2) 2.023(−2) 4.724(−3) 1.076(−3) 2.510(−4)

1.322 1.857 2.098 2.134 2.099

1.518(−2) 1.499(−3) 7.837(−5) 3.312(−6) 1.348(−7) 5.760(−9)

3.340 4.258 4.564 4.619 4.548

8.181(−4) 1.881(−5) 2.332(−7) 2.284(−9) 2.028(−11) 1.648(−13)

5.443 6.334 6.674 6.815 6.943

Remark. From the Theorems 2.1–2.3, it is easy to see that the differential transforms Gi (m) and Uˆ i (m) are obtained sequentially.  Thus, if uˆ i (t) is approximated by uˆ i (t) ≈ nm=0 Uˆ i (m)(x − xi )m ≡ sˆi,n (t), then the nonlinear function gi (u(t)) is also approximated n m by gi (u(t)) ≈ m=0 Gi (m)(x − xi ) ≡ s˜i,n (t). Then we can consider the two cases in evaluating Vˆi (m) in (18) as follows

(a). Vˆi (m) ≈ (b). Vˆi (m) ≈



xi+1 xi



xi+1 xi

Kˆt,l (j)s˜i,n (t) dt, Kˆt,l (j)g(sˆi,n (t)) dt,

(29)

where the case (b) is suggested in the proposed scheme in (19). For each case, any numerical integration formula can be applied to approximate the value of Vˆi (m). However, it is clear that the computational cost in the case (b) is cheap. a (x) and Eb are defined by The error functions Eh,n h,n a Eh,n (x) = uˆ ah,n (x) − u(x),

b Eh,n (x) = uˆ bh,n (x) − u(x),

where uˆ ah,n (x) and uˆ bh,n (x) represent the approximations obtained by each case in (29), respectively. And we define the rate of convergence by

 i  E  h,n rate = log2  i E

h/2,n

   , 

i = a, b

(30)

For all numerical examples we present the maximum error and the rate of convergence (30) as the length h varies from 1/2 to 1/64. Example 1. Consider the following nonlinear integro-differential equation:

u (x) =

3 x 1 3x e − e + 2 2



0

x

ex−t u3 (t) dt,

(31)

where u(0) = 1. It is easy to see that the exact solution u(x) = ex . In Tables 2 and 3, the numerical results are presented. As h is getting smaller, the maximum errors is also decreasing and the rate of convergence follows from the fact that it is approaching n + 1, where n is the number of partial sum of the Taylor series. For h = 1/64 and n = 5, the rate of convergence with the case (b) in (29) is one order higher than expected in Theorem 4.4, is called superconvergence. Example 2. Consider the following nonlinear integro-differential equation:

u (x) = −2 sin x −

1 2 cos x − cos(2x) + 3 3



x 0

cos(x − t)u2 (t) dt,

(32)

where u(0) = 1 and the exact solution u(x) = cos x − sin x. The numerical results are demonstrated in Tables 4 and 5. For the approximate solution with the case (a) in (19), the rate of convergence is approaching the value supported by Theorem 4.4 as the length h is decreasing. In case (b), it gives a more accurate approximate solution than expected. But, it is not as accurate as in the previous example.

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35

33

Table 4 a Maximum errors and convergence rates for Eh,n , n = 1, 3, 5 in Example 2. h

a ||Eh,1 ||∞

rate

a ||Eh,3 ||∞

rate

a ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

8.166(−2) 2.239(−2) 5.724(−3) 1.442(−3) 3.619(−4) 9.063(−5)

1.867 1.968 1.988 1.995 1.998

1.958(−3) 1.529(−4) 1.035(−5) 6.675(−7) 4.229(−8) 2.660(−9)

3.679 3.885 3.954 3.980 3.991

3.341(−5) 7.200(−7) 1.258(−8) 2.054(−10) 3.271(−12) 5.129(−14)

5.536 5.839 5.936 5.972 5.995

Table 5 b , n = 1, 3, 5 in Example 2. Maximum errors and convergence rates for Eh,n h

b ||Eh,1 ||∞

rate

b ||Eh,3 ||∞

rate

b ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

7.335(−2) 1.724(−2) 4.037(−3) 9.736(−4) 2.390(−4) 5.922(−5)

2.089 2.095 2.052 2.026 2.013

1.648(−3) 8.445(−5) 4.204(−6) 2.241(−7) 1.275(−8) 7.567(−10)

4.286 4.328 4.229 4.136 4.074

2.765(−5) 3.155(−7) 3.177(−9) 3.381(−11) 3.983(−13) 4.885(−15)

6.454 6.634 6.554 6.407 6.350

Table 6 a , n = 1, 3, 5 in Example 3. Maximum errors and convergence rates for Eh,n h

a ||Eh,1 ||∞

rate

a ||Eh,3 ||∞

rate

a ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

3.997(−2) 8.560(−3) 1.968(−3) 4.723(−4) 1.157(−4) 2.864(−5)

2.223 2.121 2.059 2.029 2.014

4.841(−3) 2.336(−4) 1.243(−5) 7.122(−7) 4.256(−8) 2.600(−9)

4.373 4.232 4.125 4.065 4.033

7.828(−4) 8.842(−6) 1.096(−7) 1.496(−9) 2.178(−11) 3.281(−13)

6.468 6.334 6.195 6.102 6.053

Table 7 b , n = 1, 3, 5 in Example 3. Maximum errors and convergence rates for Eh,n h

b ||Eh,1 ||∞

rate

b ||Eh,3 ||∞

rate

b ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

4.026(−2) 8.850(−3) 2.057(−3) 4.953(−4) 1.215(−4) 3.009(−5)

2.186 2.105 2.054 2.027 2.014

4.942(−3) 2.450(−4) 1.316(−5) 7.566(−7) 4.528(−8) 2.769(−9)

4.334 4.218 4.120 4.062 4.032

8.006(−4) 9.103(−6) 1.139(−7) 1.562(−9) 2.279(−11) 3.437(−13)

6.459 6.321 6.188 6.099 6.051

Example 3. Consider the following nonlinear integro-differential equation:

u (x) = x +

3 2 1 x + − (1 + x)2 ln(1 + x) + 2 1+x



x 0

(x − t)2 e−u(t) dt,

(33)

where u(0) = 0 and the exact solution u(x) = ln (1 + x). In Tables 6 and 7, the maximum error and the rate of convergence for the approximate solutions are listed as the length h is decreasing. Here, for both cases in (19), the similar numerical results are shown. Example 4. Consider the following nonlinear integro-differential equation:

u (x) = f (x) + where

f (x) =

 0

x

(x − t) ln(1 + u(t))dt,

(34)

     1 12 − 8 1 + x − 4x(−6 + 5 1 + x) − 12x2 ln(1 + 1 + x) , 8 + 9x2 + √ 24 1+x

with the initial condition u(0) = 1, and the exact solution u(x) =

√ 1 + x.

In Tables 8 and 9, for both cases in (19), the similar numerical results are presented as the length h is decreasing.

34

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35 Table 8 a , n = 1, 3, 5 in Example 4. Maximum errors and convergence rates for Eh,n h

a ||Eh,1 ||∞

rate

a ||Eh,3 ||∞

rate

a ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

8.464(−3) 1.955(−3) 4.666(−4) 1.137(−4) 2.806(−5) 6.969(−6)

2.114 2.067 2.036 2.019 2.010

7.536(−4) 3.868(−5) 2.141(−6) 1.252(−7) 7.558(−9) 4.641(−10)

4.284 4.175 4.096 4.050 4.026

9.937(−5) 1.166(−6) 1.502(−8) 2.101(−10) 3.095(−12) 4.730(−14)

6.413 6.279 6.160 6.085 6.032

Table 9 b , n = 1, 3, 5 in Example 4. Maximum errors and convergence rates for Eh,n h

b ||Eh,1 ||∞

rate

b ||Eh,3 ||∞

rate

b ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

8.414(−3) 1.925(−3) 4.579(−4) 1.115(−4) 2.749(−5) 6.824(−6)

2.128 2.072 2.038 2.020 2.010

7.498(−4) 3.787(−5) 2.082(−6) 1.214(−7) 7.322(−9) 4.494(−10)

4.307 4.185 4.100 4.052 4.026

9.895(−5) 1.138(−6) 1.451(−8) 2.022(−10) 2.973(−12) 4.508(−14)

6.443 6.292 6.166 6.087 6.044

Table 10 a , n = 1, 3, 5 in Example 5. Maximum errors and convergence rates for Eh,n h

a ||Eh,1 ||∞

rate

a ||Eh,3 ||∞

rate

a ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

7.519(−2) 1.347(−2) 2.809(−3) 6.398(−4) 1.525(−4) 3.726(−5)

2.481 2.262 2.134 2.068 2.033

2.334(−2) 1.864(−3) 1.147(−4) 6.845(−6) 4.136(−7) 2.534(−8)

3.646 4.023 4.066 4.049 4.029

4.569(−4) 7.396(−6) 3.520(−7) 7.498(−9) 1.322(−10) 2.178(−12)

5.949 4.393 5.553 5.825 5.924

Table 11 b , m = 1, 3, 5 in Example 5. Maximum errors and convergence rates for Eh,n h

b ||Eh,1 ||∞

rate

b ||Eh,3 ||∞

rate

b ||Eh,5 ||∞

rate

1/2 1/4 1/8 1/16 1/32 1/64

7.488(−2) 1.416(−2) 3.650(−3) 9.728(−4) 2.534(−4) 6.477(−5)

2.402 1.956 1.908 1.941 1.968

1.686(−2) 5.120(−4) 1.485(−5) 5.123(−7) 2.156(−8) 1.064(−9)

5.041 5.108 4.857 4.571 4.341

1.776(−4) 7.017(−6) 8.454(−8) 7.468(−10) 6.190(−12) 5.063(−14)

4.661 6.375 6.823 6.915 6.934

Example 5. Consider the following nonlinear integro-differential equation:

u (x) = f (x) +



0

x

ex−t sin(u(t)) dt,

(35)

where f(x) = −e1 − x + e−(1 − x) (cos (e) − cos (e1 − x )) with the initial condition u(0) = e, and the exact solution u(x) = e1 − x . The numerical results are listed in Tables 10 and 11. As stated in Example 1, the rate of convergence for the approximate solution with case (b) in (19) is approaching the superconvergence as the length h is decreasing and the number of partial sum n in the Taylor series is increasing. 6. Conclusion In this paper, we propose an efficient semi-analytical approach, namely the multistage differential transform method(MsDTM) for solving nonlinear Volterra integro-differential equations. The basic idea of the proposed method is to obtain an approximate solution by using the Taylor series in each subdomain. In DTM, it is the key to find a recursive relation for the differential transforms from the given equation. We derive the recursive relation of the differential transforms for the nonlinear Volterra integro-differential equation with complex nonlinear kernels. The maximum error and the rate of convergence of the approximate solutions by the MsDTM are demonstrated. It is shown that all numerical illustrations support the results in the theory of convegence. For some examples whose solutions are smooth enough, the rate of convergence is one order faster than expected in the result of theory, called a superconvergence. One of the strong advantage in MsDTM it is very simple to implement. All we

K. Kim, B. Jang / Applied Mathematics and Computation 263 (2015) 25–35

35

have to do is to construct the recursive relation for the differential transforms. Moreover, as shown in convergence analysis the scheme is highly accurate. Thus, it can be concluded that the proposed scheme is very efficient in finding accurate approximate solutions. Here, all computations are performed by using Mathematica 8.0. Acknowledgment B. Jang was supported by Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education (NRF-2014R1A1A2A16051147). References [1] I. Abdul, Introduction to Integral Equations with Application, Wiley, New York, 1999. [2] H. Brunner, Collocation Methods for Volterra Integral and Related Functional Equation Method, Cambridge University Press, 2004. [3] K. Maleknejad, Y. Mahmoudi, Taylor polynomial solution of high-order nonlinear Volterra–Fredholm integro-differential equations, Appl. Math. Comput. 145 (2003) 641–653. [4] S.M. Hosseini, S. Shahmorad, Numerical solution of a class of integro-differential equations by the tau method with an error estimation, Appl. Math. Comput. 136 (2–3) (2003) 559–570. [5] K. Maleknejad, F. Mirzaee, Using rationalized Haar wavelet for solving linear integral equations, Appl. Math. Comput. 160 (2005) 579–587. [6] S.M. El-Sayed, M.R. Abdel-Aziz, A comparison of Adomian’s decomposition method and wavelet-Galerkin method for solving integro-differential equations, Appl. Math. Comput. 136 (2003) 151–159. [7] A. Shidfar, A. Molabahrami, A. Babaei, A. Yazdanian, A series solution of the nonlinear Volterra and Fredholm integro-differential equations, Commun. Nonlinear. Sci. Numer. Simul. 15 (2010) 205–215. [8] S. Kumar, A new analytical modelling for fractional telegraph equation via Laplace transform, Appl. Math. Model. 38 (13) (2014) 3154–3163. [9] S. Kumar, M.M. Rashidi, New analytical method for gas dynamics equation arising in shock fronts, Comput. Phys. Commun. 185 (7) (2014) 1947–1954. [10] S. Kumar, A numerical study for solution of time fractional nonlinear shallow-water equation in oceans, Z. Naturforsch A 68a (2013) 1–7. [11] S. Kumar, Numerical computation of time-fractional Fokker Planck equation arising in solid state physics and circuit theory, Z. Naturforsch. A 68a (2013) 777–784. [12] S. Kumar, O.P. Singh, Numerical inversion of the Abel integral equation using homotopy perturbation method, Z. Naturforsch. A 65a (2009) 677–682. [13] M. Zarebnia, Sinc numerical solution for Volterra integro-differential equation, Commun. Nonlinear Sci. Numer. Simul. 15 (3) (2010) 700–706. [14] M.J. Jang, C.L. Chen, Y.C. Liu, Two-dimensional differential transform for partial differential equations, Appl. Math. Comput. 121 (2001) 261–270. [15] B. Jang, Comments on “solving a class of two-dimensional linear and nonlinear Volterra integral equations by the differential transform method”, J. Comput. Appl. Math. 233 (2) (2009) 224–230. [16] B. Jang, Solving linear and nonlinear initial value problems by the projected differential transform method, Comput. Phys. Commun. 181 (2010) 848–854. [17] A. Gökdo˘gan, M. Merdan, A. Yildirim, A multistage differential transformation method for approximate solution of Hantavirus infection model, Commun. Nonlinear. Sci. Numer. Simul. 17 (2012) 1–8. [18] Z.M. Odibat, C. Bertelle, M.A. Aziz-Alaoui, G.H.E. Duchamp, A multi-step differential transform method and application to non-chaotic or chaotic systems, Comput. Math. Appl. 59 (4) (2010) 1462–1472. [19] S.-H. Chang, I.-L. Chang, A new algorithm for calculating one-dimensional differential transform of nonlinear functions, Appl. Math. Comput. 195 (2008) 799–808. [20] R. Li, Numerical Solution of Differential Equation, High Education Press, Beijing, 1996.