# L-Stable Derivative-Free Error-Corrected Trapezoidal Rule for Burgers' Equation with Inconsistent Initial and Boundary Conditions

Authors: Lajja Verma
Journal: International Journal of Mathematics and Mathematical Sciences
 (2012)

---

## Abstract

We combine suitable arithmetic average approximations, with explicit backward Euler formula, and derive a third-orderL-stable derivative-free error-corrected trapezoidal rule (LSDFECT). Then, we apply LSDFECT rule to the linearized Burgers' equation with inconsistent initial and boundary conditions and test its stability and exactness. We use Mathematica 7.0 for computation.


---

## Body

## 1. Introduction

Nonlinear phenomena play a crucial role in applied mathematics and physics and in mechanics and biology. Burgers’ equation has been found to describe various kinds of phenomena such as a mathematical model of turbulence and the approximate theory of flow through a shock wave traveling in a viscous fluid. The equation has found applications in field as diverse as number theory, gas dynamics, heat conduction, and elasticity. In past years, a growing interest has been given to the propagation of nonlinear waves in nonlinear dynamical systems. A nonlinear wave is one of the fundamental objects of nature. These waves appear in a great array of contexts such as hydrodynamics, nonlinear optics, plasmas, solid-state physics, nuclear physics, and many other nonlinear phenomena. The systems are often described by nonlinear partial differential equations with constant or variable coefficients. In this paper, we are proposing a method which can deal efficiently with numerical solutions of Burgers’ equation which contains inconsistencies in terms of initial and boundary conditions.L-stable numerical methods are sought whenever inconsistencies are present in the initial and boundary conditions or in some other cases. So, in the present paper, we use derivative-free error-corrected trapezoidal rule [1], arithmetic average approximation, and explicit Euler’s backward formula on y'=f(t,y) and develop L-stable derivative-free error-corrected trapezoidal type quadrature formula.We use these results to study the numerical solution of one-dimensional quasilinear partial differential equation, that is, Burgers’ equation(1.1)∂w∂t+w∂w∂x-νd2∂2w∂x2=0,(x,t)∈[0,1]×[0,T],
with initial condition
(1.2)w(x,0)=g(x),0<x<1,
and Dirichlet boundary conditions
(1.3)w(0,t)=0=w(1,t),0≤t≤T,
where νd>0 is a coefficient of viscous diffusion and g is sufficiently smooth given function.We consider some examples of Burgers’ equation with Dirichlet boundary condition and transform it into linear diffusion equation with Neumann boundary condition for which CN produces unacceptable oscillations while the presentL-stable derivative-free error-corrected trapezoidal type rule provides both stable and accurate result.The Crank and Nicolson [2] scheme is A-stable not L-stable and produces undesirable oscillations in the case of inconsistencies in the initial and boundary conditions [3] or when time step taken is large [4]. As an alternative to CN, Chawla et al. [5] have proposed generalized trapezoidal formula (GTF(α)) which is second-order L-stable schemes for α>0 and gives quite stable and accurate approximations for a judicious choice of the parameter α.For higher-order accuracy, Chawla et al. [6] proposed a modified A-stable Simpson’s 1/3 rule (ASIMP) but, since it lacks L-stability, ASIMP also produces unwanted oscillations. To rectify this, Chawla and Evans [7] and Pandey et al. [8] presented a L-stable version of the Simpson’s 1/3 rule and Simpson’s 3/8 rule, respectively, and employed it to derive a third-order time integration scheme for the diffusion equation. These L-stable schemes provide stable and accurate approximations.This paper is organized in 4 sections. In Section2, we discuss derivative-free error-corrected trapezoidal rule. In Section 3 we derive L-stable version of derivative free error corrected trapezoidal rule and discuss its local truncation error and stability and apply it to linearized Burgers’ equation. Finally, in Section 4, we illustrate our findings numerically.


## 2. Derivative-Free Error-Corrected Trapezoidal Rule

Letf(x) be defined over the interval [a,b]. Let us divide the interval [a,b] into n subintervals [xi-1,xi] for i=1(1)n, where xi=x0+i*h, x0=a and xn=b. Then
(2.1)∫ab‍f(x)dx=∑i=1n‍Trap[xi-1,xi],∫ab‍f(x)dx=12{∑i=1n-1‍Simp1/3[xi-1,xi+1]+Simp1/3[x0,x1]+Simp1/3[xn-1,xn]},
where
(2.2)Trap[xi-1,xi]=h2[f(xi-1)+f(xi)],Simp1/3[xi-1,xi+1]=h3[f(xi-1)+4f(xi)+f(xi+1)].Now derivative-free corrector is the difference between (2.1) and is given by (cf. [1])
(2.3)-h4(f0+fn)+h3(f1/2+fn-(1/2))-h12(f1+fn-1).Hence derivative-free error-corrected trapezoidal rule is given by(2.4)h[12f0+f1+f2+⋯+fn-1+12fn]-h4(f0+fn)+h3(f1/2+fn-(1/2))-h12(f1+fn-1).The above quadrature formula will have error proportional to(b-a)h4, and it does not contain any derivative term.


## 3.L-Stable Derivative-Free Error-Corrected Trapezoidal-Type Rule

Definition 3.1.
A method is said to beL-stable if it is A-stable, and, in addition, when applied to the scalar test equation y'=λy, λ∈ℂ, Re(λ)<0, it yields yn+1=Ψ(z)yn, z=λh (h is step size), where |Ψ(z)|→0 as Re(z)→-∞.Consider the first-order initial value problem(3.1)y′(t)=f(t,y),y(t0)=η.The derivative-free error-corrected trapezoidal-type rule is given by(3.2)yn+1=yn+h36[3(fn+fn+1)+11(fn+(1/3)+fn+(2/3))+4(fn+(1/6)+fn+(5/6))].We start with following suitable weighted average approximations(3.3)yn+1/3=2yn+yn+13,yn+2/3=yn+2yn+13,yn+1/6=5yn+yn+16,yn+5/6=yn+5yn+16.The explicit backward Euler formula (Chawla et al. [9]) is
(3.4)yn=yn+1-hfn+1.Now we use (3.3) to obtain the following approximations:
(3.5)y¯n+1/3=4yn+5yn+19-29hfn+1,y¯n+2/3=yn+8yn+19-29hfn+1,y¯n+1/6=23yn+13yn+136-736hfn+1,y¯n+5/6=yn+11yn+112-112hfn+1.Definef¯i=f(xi,y¯i) and use it to derive our new LSDFECT rule which is given by
(3.6)yn+1=yn+h36[3(fn+fn+1)+11(f¯n+1/3+f¯n+2/3)+4(f¯n+1/6+f¯n+5/6)].


### 3.1. Local Truncation Error

Since(3.7)yn+1/3=2yn+yn+13-19h2yn′′-481h3yn′′′-13324h4yn(4)+O(h5),yn+2/3=yn+2yn+13-19h2yn′′-581h3yn′′′-19972h4yn(4)+O(h5),yn+1/6=5yn+yn+16-572h2yn′′-351296h3yn′′′-21531104h4yn(4)+O(h5),yn+5/6=yn+5yn+16-572h2yn′′-551296h3yn′′′-45531104h4yn(4)+O(h5),yn=yn+1-hfn+1+12h2yn′′+13h3yn′′′+18h4yn(4)+O(h5).We deduce that(3.8)y¯n+1/3=yn+1/3-281h3yn′′′+181h4yn(4)+O(h5),y¯n+2/3=yn+2/3-181h3yn′′′-2243h4yn(4)+O(h5),y¯n+1/6=yn+1/6-136h2yn′′+231296h3yn(3)+O(h4),y¯n+5/6=yn+5/6+136h2yn′′′-191296h3yn(3)+O(h4).
Again since
(3.9)yn+1=yn+h36[3(fn+fn+1)+11(fn+1/3+fn+2/3)+4(fn+1/6+fn+5/6)]+O(h5),
it can be shown that
(3.10)yn+1=yn+h36[3(fn+fn+1)+11(f¯n+1/3+f¯n+2/3)+4(f¯n+1/6+f¯n+5/6)]+tn(h),
where
(3.11)tn(h)=O(h4).Thus scheme is of third order.


### 3.2. Stability

Applying scheme (3.6) to the test equation
(3.12)y'=-λy,λ>0,
and setting z=λh, we have
(3.13)yn+1=Ψ(z)yn,
where Ψ(z)=((1-(1/3)z)/(1+(2/3)z+(1/6)z2)). Since Ψ(z)→0 as z→∞, the scheme (3.6) is L-stable.Therefore, it follows that time integration scheme (3.6) is of third order and L-stable. We call it LSDFECT rule for the first-order initial value problem.


## 3.1. Local Truncation Error

Since(3.7)yn+1/3=2yn+yn+13-19h2yn′′-481h3yn′′′-13324h4yn(4)+O(h5),yn+2/3=yn+2yn+13-19h2yn′′-581h3yn′′′-19972h4yn(4)+O(h5),yn+1/6=5yn+yn+16-572h2yn′′-351296h3yn′′′-21531104h4yn(4)+O(h5),yn+5/6=yn+5yn+16-572h2yn′′-551296h3yn′′′-45531104h4yn(4)+O(h5),yn=yn+1-hfn+1+12h2yn′′+13h3yn′′′+18h4yn(4)+O(h5).We deduce that(3.8)y¯n+1/3=yn+1/3-281h3yn′′′+181h4yn(4)+O(h5),y¯n+2/3=yn+2/3-181h3yn′′′-2243h4yn(4)+O(h5),y¯n+1/6=yn+1/6-136h2yn′′+231296h3yn(3)+O(h4),y¯n+5/6=yn+5/6+136h2yn′′′-191296h3yn(3)+O(h4).
Again since
(3.9)yn+1=yn+h36[3(fn+fn+1)+11(fn+1/3+fn+2/3)+4(fn+1/6+fn+5/6)]+O(h5),
it can be shown that
(3.10)yn+1=yn+h36[3(fn+fn+1)+11(f¯n+1/3+f¯n+2/3)+4(f¯n+1/6+f¯n+5/6)]+tn(h),
where
(3.11)tn(h)=O(h4).Thus scheme is of third order.


## 3.2. Stability

Applying scheme (3.6) to the test equation
(3.12)y'=-λy,λ>0,
and setting z=λh, we have
(3.13)yn+1=Ψ(z)yn,
where Ψ(z)=((1-(1/3)z)/(1+(2/3)z+(1/6)z2)). Since Ψ(z)→0 as z→∞, the scheme (3.6) is L-stable.Therefore, it follows that time integration scheme (3.6) is of third order and L-stable. We call it LSDFECT rule for the first-order initial value problem.


## 4. LSDFECT Rule for Partial Differential Equation

### 4.1. The Hopf and Cole Transformation

Hopf [10] and Cole [11] suggested that (1.1) can be reduced to linear heat equation by a nonlinear transformation given by
(4.1)ψ=-νd(logϕ),w=ψx.Using (4.1) in (1.1), we have
(4.2)∂ϕ∂t=νd2∂2ϕ∂x2,(4.3)ϕx(0,t)=ϕx(1,t)=0,(4.4)ϕ(x,0)=f(x).The Fourier series solution to the linearized heat equation (4.2) is
(4.5)ϕ(x,t)=A0+∑n=1∞‍Anexp(-νdn2π2t2)cosnπx
with Fourier coefficients at t=0 as
(4.6)A0=∫01‍exp(-1νd∫0x‍w0(ξ)dξ)dx,An=2∫01‍exp(-1νd∫0x‍w0(ξ)dξ)cos(nπx)dx,
where w0(ξ)=w(ξ,0).Hence, the exact solution of (1.1) is given by
(4.7)w(x,t)=πνd(∑n=1∞‍Anexp(-νdn2π2t/2)nsinnπxA0+∑n=1∞‍Anexp(-νdn2π2t/2)cosnπx).


### 4.2. LSDFECT Rule for Linear Parabolic Equation with Neumann Boundary Condition

We consider the application of LSDFECT rule for the time integration of diffusion equation (4.2).As a first step, replacing the spatial derivative by the central second-order finite difference formula, we obtain the resulting semidiscretization as follows.For a positive integerN, define the spatial grid xi=ih, i=0,1,…,N, where h=1/N. For a positive integer M, define the temporal grid tj=jk, j=0,1,…,M. Now set ϕi(t)=ϕ(xi,t) and so forth. Let
(4.8)Φ(t)=[ϕ0(t),ϕ1(t),…,ϕN(t)]T,
and J is a (N+1)×(N+1) matrix given by
(4.9)J=(2-200⋯00-12-10⋯000-12-1⋯00⋮⋮⋱⋱⋱⋱⋮000-12-100000-12-100000-22).
The semidiscretization of (4.2) with the Neumann boundary conditions
(4.10)ϕx(0,t)=0=ϕx(1,t)
can be written as follows:
(4.11)∂Φ(t)∂t=-νd2h2JΦ(t),
where Φ(0)=[f(x0),f(x1),…,f(xN)]. Now setting ρ=νdk/2h2 and applying LSDFECT rule for the time integration of (4.11), we get
(4.12)Φj+1=Φj-136ρJ[3(Φj+Φj+1)+11(Φ¯j+1/3+Φ¯j+2/3)+4(Φ¯j+1/6+Φ¯j+5/6)],
where
(4.13)Φ¯j+1/3=19[4Φj+(5I+2ρJ)Φj+1],Φ¯j+2/3=19[Φj+(8I+2ρJ)Φj+1],Φ¯j+1/6=136[23Φj+(13I+7ρJ)Φj+1],Φ¯j+5/6=112[Φj+(11I+ρJ)Φj+1],
and I is a (N+1)×(N+1) identity matrix. Using (4.12) to (4.13), we deduce
(4.14)(I+23ρJ+16ρ2J2)Φ¯j+1=(I-13ρJ)Φ¯j,
which is the LSDFECT rule for the linearized Burgers’ equation, which is the same as the L-stable scheme developed by Chawla and Evans [7] and later by Pandey et al. [8].


### 4.3. Unconditional Stability

We write (4.14) in the following form:
(4.15)Φ¯j+1=QΦ¯j,j=0,1,2…,
where
(4.16)Q=(I+23ρJ+16ρ2J2)-1(I-13ρJ).
Using the Gershgorin Circle Theorem, one can prove that J has positive eigenvalues. Then, it is easy to see that all the eigenvalues of Q are inside unit circle. Thus, LSDFECT rule is an unconditionally stable method for all ρ>0.


## 4.1. The Hopf and Cole Transformation

Hopf [10] and Cole [11] suggested that (1.1) can be reduced to linear heat equation by a nonlinear transformation given by
(4.1)ψ=-νd(logϕ),w=ψx.Using (4.1) in (1.1), we have
(4.2)∂ϕ∂t=νd2∂2ϕ∂x2,(4.3)ϕx(0,t)=ϕx(1,t)=0,(4.4)ϕ(x,0)=f(x).The Fourier series solution to the linearized heat equation (4.2) is
(4.5)ϕ(x,t)=A0+∑n=1∞‍Anexp(-νdn2π2t2)cosnπx
with Fourier coefficients at t=0 as
(4.6)A0=∫01‍exp(-1νd∫0x‍w0(ξ)dξ)dx,An=2∫01‍exp(-1νd∫0x‍w0(ξ)dξ)cos(nπx)dx,
where w0(ξ)=w(ξ,0).Hence, the exact solution of (1.1) is given by
(4.7)w(x,t)=πνd(∑n=1∞‍Anexp(-νdn2π2t/2)nsinnπxA0+∑n=1∞‍Anexp(-νdn2π2t/2)cosnπx).


## 4.2. LSDFECT Rule for Linear Parabolic Equation with Neumann Boundary Condition

We consider the application of LSDFECT rule for the time integration of diffusion equation (4.2).As a first step, replacing the spatial derivative by the central second-order finite difference formula, we obtain the resulting semidiscretization as follows.For a positive integerN, define the spatial grid xi=ih, i=0,1,…,N, where h=1/N. For a positive integer M, define the temporal grid tj=jk, j=0,1,…,M. Now set ϕi(t)=ϕ(xi,t) and so forth. Let
(4.8)Φ(t)=[ϕ0(t),ϕ1(t),…,ϕN(t)]T,
and J is a (N+1)×(N+1) matrix given by
(4.9)J=(2-200⋯00-12-10⋯000-12-1⋯00⋮⋮⋱⋱⋱⋱⋮000-12-100000-12-100000-22).
The semidiscretization of (4.2) with the Neumann boundary conditions
(4.10)ϕx(0,t)=0=ϕx(1,t)
can be written as follows:
(4.11)∂Φ(t)∂t=-νd2h2JΦ(t),
where Φ(0)=[f(x0),f(x1),…,f(xN)]. Now setting ρ=νdk/2h2 and applying LSDFECT rule for the time integration of (4.11), we get
(4.12)Φj+1=Φj-136ρJ[3(Φj+Φj+1)+11(Φ¯j+1/3+Φ¯j+2/3)+4(Φ¯j+1/6+Φ¯j+5/6)],
where
(4.13)Φ¯j+1/3=19[4Φj+(5I+2ρJ)Φj+1],Φ¯j+2/3=19[Φj+(8I+2ρJ)Φj+1],Φ¯j+1/6=136[23Φj+(13I+7ρJ)Φj+1],Φ¯j+5/6=112[Φj+(11I+ρJ)Φj+1],
and I is a (N+1)×(N+1) identity matrix. Using (4.12) to (4.13), we deduce
(4.14)(I+23ρJ+16ρ2J2)Φ¯j+1=(I-13ρJ)Φ¯j,
which is the LSDFECT rule for the linearized Burgers’ equation, which is the same as the L-stable scheme developed by Chawla and Evans [7] and later by Pandey et al. [8].


## 4.3. Unconditional Stability

We write (4.14) in the following form:
(4.15)Φ¯j+1=QΦ¯j,j=0,1,2…,
where
(4.16)Q=(I+23ρJ+16ρ2J2)-1(I-13ρJ).
Using the Gershgorin Circle Theorem, one can prove that J has positive eigenvalues. Then, it is easy to see that all the eigenvalues of Q are inside unit circle. Thus, LSDFECT rule is an unconditionally stable method for all ρ>0.


## 5. Numerical Illustrations

To illustrate the computational performance of the LSDFECT rule and to compare its performance with CN, we consider five examples. First four examples contain inconsistent initial and boundary conditions. It is easily seen from Figures1, 2, 3, and 4 that LSDFECT rule gives us stable and accurate result whereas CN method has large oscillations near boundary points.Figure 1
Comparison of the numerical solution by our method (LSDFECT), CN method, and exact solution of Problem 1 at different space points forN=80, νd=2, and k=0.01 at T=0.1.Figure 2
Comparison of the numerical solution by our method (LSDFECT), CN method, and exact solution of Problem 2 at different space points forN=80, νd=2, and k=0.01 at T=0.1.Figure 3
Comparison of the numerical solution by our method (LSDFECT), CN method, and exact solution of Problem 3 at different space points forN=80, νd=2, and k=0.01 at T=0.1.Figure 4
Comparison of the numerical solution by our method (LSDFECT), CN method, and exact solution of Problem 4 at different space points forN=80, νd=2, and k=0.01 at T=0.1.In Problem 5, we have considered consistent initial and boundary conditions. We have calculated solutions forh=0.1 and k=4.0 for νd=0.02 for time T=80. We observe that, since time step size k>1/π, CN method gives large oscillations which is unacceptable, while LSDFECT rule gives quite stable and accurate solutions (see Figure 5).Figure 5
Comparison of the numerical solution by our method (LSDFECT), CN method, and exact solution of Problem 5 at different space points forh=0.1, νd=0.02, and k=4 at T=80.


### 5.1. Problem 1

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.1)w(x,0)=2,0<x<1.
The exact solution is given by (4.7), where
(5.2)A0=∫01‍exp(-2xνd)dx,An=2∫01‍exp(-2xνd)cos(nπx)dx.


### 5.2. Problem 2

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.3)w(x,0)=2x,0<x<1.
The exact solution is given by (4.7), where
(5.4)A0=∫01‍exp(-x2νd)dx,An=2∫01‍exp(-x2νd)cos(nπx)dx.


### 5.3. Problem 3

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.5)w(x,0)=cosπ4x,0<x<1.
The exact solution is given by (4.7), where
(5.6)A0=∫01‍exp(-4πνdsinπ4x)dx,An=2∫01‍exp(-4πνdsinπ4x)cos(nπx)dx.


### 5.4. Problem 4

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.7)w(x,0)=sinπ4x,0<x<1.
The exact solution is given by (4.7), where
(5.8)A0=∫01‍exp(4πνd(cosπ4x-1))dx,An=2∫01‍exp(4πνd(cosπ4x-1))cos(nπx)dx.


### 5.5. Problem 5

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.9)w(x,0)=sinπx,0<x<1.The exact solution is given by (4.7), where
(5.10)A0=∫01‍exp(-1πνd(1-cosπx))dx,An=2∫01‍exp(-1πνd(1-cosπx))cosnπxdx.


## 5.1. Problem 1

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.1)w(x,0)=2,0<x<1.
The exact solution is given by (4.7), where
(5.2)A0=∫01‍exp(-2xνd)dx,An=2∫01‍exp(-2xνd)cos(nπx)dx.


## 5.2. Problem 2

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.3)w(x,0)=2x,0<x<1.
The exact solution is given by (4.7), where
(5.4)A0=∫01‍exp(-x2νd)dx,An=2∫01‍exp(-x2νd)cos(nπx)dx.


## 5.3. Problem 3

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.5)w(x,0)=cosπ4x,0<x<1.
The exact solution is given by (4.7), where
(5.6)A0=∫01‍exp(-4πνdsinπ4x)dx,An=2∫01‍exp(-4πνdsinπ4x)cos(nπx)dx.


## 5.4. Problem 4

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.7)w(x,0)=sinπ4x,0<x<1.
The exact solution is given by (4.7), where
(5.8)A0=∫01‍exp(4πνd(cosπ4x-1))dx,An=2∫01‍exp(4πνd(cosπ4x-1))cos(nπx)dx.


## 5.5. Problem 5

Consider (1.1) with the boundary conditions (1.3) and initial condition
(5.9)w(x,0)=sinπx,0<x<1.The exact solution is given by (4.7), where
(5.10)A0=∫01‍exp(-1πνd(1-cosπx))dx,An=2∫01‍exp(-1πνd(1-cosπx))cosnπxdx.


## 6. Conclusion

Using derivative-free error-corrected trapezoidal rule, suitable average approximation, and explicit backward Euler formula, we derive third-order integration scheme in time which isL-stable. The present method provides accurate and stable results whenever inconsistencies are present.



---

*Source: 821907-2012-09-27.xml*