# An Unconditionally Stable Parallel Difference Scheme for Telegraph Equation

Authors: A. Borhanifar; Reza  Abazari
Journal: Mathematical Problems in Engineering
 (2009)

---

## Abstract

We use an unconditionally stable parallel difference scheme to solve
telegraph equation. This  method is based on domain decomposition concept and using
asymmetric Saul'yev schemes for  internal nodes of each sub-domain and alternating
group implicit method for sub-domain's interfacial nodes. This new method has several
advantages such as: good parallelism, unconditional stability and better accuracy than
original Saul'yev schemes. The details of implementation and proving stability are
briefly discussed. Numerical experiments on stability and accuracy are also presented.


---

## Body

## 1. Introduction

The large-scale scientific and engineering computations need the parallel computer with massively parallel processors of higher speed and large memory and also need effective parallel numerical methods and parallel algorithms. Some numerical methods directly have parallel character and have great efficiency on using in the parallel computer, but much numerical methods usually need to be reconstructed to be more appropriate for the parallel computing [1].Consider the Telegraph equation of the form(1.1)∂2u∂t2+c∂u∂t=a2∂2u∂x2+bu+F(x,t),
where c>0 and b<0. Numerical solution of this equation based on finite difference has been successfully proposed in [2–4].There are two types of schemes for numerical solution of time-dependent partial differential equations, implicit and explicit schemes. The former has no restriction on its time stepping. But in each time step one has to solve a global system of equations. The implementation on parallel computer is not straightforward due to its global nature. The latter is easy to program and implement on parallel computer. However, it suffers the severely restricted time step size from stability requirement [5]. Recently, a so-called alternating explicit-implicit scheme has been studied by many authors, which is based on the concept of domain decomposition and a combination of implicit schemes and explicit schemes [6, 7]. In these approaches, a physical domain is divided into several sub-domains. At each time step, the problem is solved in some sub-domain implicitly and others explicitly. For the next time step the solution scheme is changed (implicit is changed to explicit and vice versa). The main advantage of alternating explicit-implicit schemes is their parallelism. Although some families of these schemes are unconditionally stable, some families of these schemes still give rise to restriction involving the time step, however, this restriction is much less severe than for a fully explicit scheme. In spite of intrinsic parallel nature and good stability of alternating explicit-implicit schemes, these schemes still have implicit phase, that is, solution of problem implicitly over implicit subdomains. Implicit phase leads to solution of linear or nonlinear algebraic system of equations at each time level that is not desirable due to high computational cost, especially for large-scale problems. Therefore an unconditionally stable fully explicit scheme with good parallelism is desirable. In [7], Saul'yev has introduced two unconditionally stable fully explicit asymmetric schemes for solution of diffusion equation. Unfortunately the Saul'yev schemes are not intrinsically parallel.In this paper we developed a new family of group explicit scheme for solution of Telegraph equations. The presented scheme is based on domain decomposition concept and using asymmetric Saul'yev schemes for internal nodes of each sub-domain and Evan's group explicit method [8] for sub-domain’s boundary nodes. The rest of this paper is organized as follows. In Section 2, we construct the scheme and describe details of implementation. In Section 3 the stability analysis is presented and in Section 4, some numerical examples are given to show the stability and accuracy of presented scheme.


## 2. Construction of New Group Explicit Scheme

In this section the numerical solution of the following problem is considered:(2.1)∂2u∂t2+c∂u∂t=a2∂2u∂x2+bu+F(x,t),(x,t)∈Ω,u(x,0)=f(x),x∈[0,1],ut(x,0)=g(x),x∈[0,1],ux(0,t)=α(t),t∈[0,T],ux(1,t)=β(t),t∈[0,T],
where c>0, b<0, f(x) and g(x)
are initial conditions, and α(t) and β(t) are prescribed boundary conditions. The domain Ω=[0,1]×[0,T] will  be divided into a M×N mesh with spatial step size h=1/M in x direction and the time step size k=T/N, respectively. Grid points (xi,tj) are given by xi=ih,i=0,1,2,…,M, and tj=jk,j=0,1,2,…,N, in which M and N are positive integers. We use uij to denote the finite difference approximations of u(ih,jk). In [7], Saul’yev has introduced two asymmetric approximations for the second derivative uxx as follows:
(2.2)A=ui+1j-uij-uij+1+ui-1j+1h2,B=ui+1j+1-uij+1-uij+ui-1jh2
By combining with the usual form for utt+cut-bu, one gets the Saul’yev formula,
(2.3)uij-1-2uij+uij+1k2+cuij+1-uij-12k-buij={a2A+Fij,a2B+Fij.
So for the Saul'yev scheme A, we have
(2.4)-ra2ui-1j+1+(1+ck2+ra2)uij+1=(2+bk2-ra2)uij+ra2ui+1j-(1-ck2)uij-1+k2Fij,
for i=0,1,2,…,M, where r=k2/h2. Although the above approximation does not appear explicit, because uij+1 and ui-1j+1 are on the left-hand side, a suitable use of the equation makes it explicit. Hence we write (2.4) in the following form:
(2.5)uij+1=11+ck/2+ra2[ra2ui-1j+1+(2+bk2-ra2)uij+ra2ui+1j-(1-ck2)uij-1+k2Fij],
for i=0,1,2,…,M. It is important to note that the Saul’yev scheme A is explicit if the calculation begins at the left boundary, and moves to the right, so that only the single value uij+1 is unknown. It can be easily seen that this explicit formula is unconditionally stable. Replacing Saul’yev B in (2.3), then
(2.6)(1+ck2+ra2)uij+1-ra2ui+1j+1=ra2ui-1j+(2+bk2-ra2)uij-(1-ck2)uij-1+k2Fij,
for i=M,M-1,…,1,0. The computational molecule of formulae (2.4) and (2.6) is shown in Figure 1. The Saul’yev formula B is in the following form:
(2.7)uij+1=11+ck/2+ra2[ra2ui+1j+1+ra2ui-1j+(2+bk2-ra2)uij-(1-ck2)uij-1+k2Fij],
for i=M,M-1,…,1,0. The Saul’yev scheme B has a similar explicit nature if the calculations proceed to the left from the right-hand boundary. Application of Saul’yev formulae A and B alternatively give better accuracy and error distribution. Obviously the Saul’yev scheme A or B is not intrinsically parallel.The computational molecule for the Saul’yev formulaeA and B.For parallel implementation of Saul’yev scheme, we use concept of domain decomposition. For convenience, we only decompose the domain[0,1] into two sub-domains, which are [0,xξ] and [xξ+1,1]. If the formula A is used for the sub-domain 1, ([0,xξ]), and the calculation begins at the left boundary, and moves to the right and in similar manner, the formula B is used for the sub-domain 2, ([xξ+1,1]), and the calculations proceed to the left from the right-hand boundary, then the calculation of each sub-domain is independent from another. For the next time step the formula B is used for the sub-domain 1 and formula A is used for the sub-domain 2 and the calculation is started from points xξ and xξ+1 for sub-domains 1 and 2, respectively. But we apply the initial condition ut(x,0)=g(x) only in the j=1, therefore, we consider two cases j=0,1 and j≥2. In case 1, and j=0, for i=0,1,…,ξ the solution is approximated by (2.4)
(2.8)-ra2ui-11+(1+ck2+ra2)ui1=(2+bk2-ra2)ui0+ra2ui+10-(1-ck2)ui-1+k2Fi0.
By initial condition ut(x,0)=g(x), we have ui-1=ui1-2kg(xi) , then
(2.9)-ra2ui-11+(2+ra2)ui1=(2+bk2-ra2)ui0+ra2ui+10+2k(1-ck2)g(xi)+k2Fi0,
while for i=ξ+1,ξ+2,…,M the solution is approximated by (2.6)
(2.10)(1+ck2+ra2)ui1-ra2ui+11=ra2ui-10+(2+bk2-ra2)ui0-(1-ck2)ui-1+k2Fi0.
Similarly, by initial condition ut(x,0)=g(x), we have
(2.11)(2+ra2)ui1-ra2ui+11=ra2ui-10+(2+bk2-ra2)ui0+2k(1-ck2)g(xi)+k2Fi0.
Also, in case 1, and j=1, for sub-domain 1, we use formula (2.6) and for sub-domain 2, we use formula (2.4), then at node (xξ,t2), the solution is approximated by (2.6)
(2.12)(1+ck2+ra2)uξ2-ra2uξ+12=ra2uξ-11+(2+bk2-ra2)uξ1-(1-ck2)uξ0+k2Fξ1,
and at node (xξ+1,t2), the solution is approximated by (2.4)
(2.13)-ra2uξ2+(1+ck2+ra2)uξ+12=(2+bk2-ra2)uξ+11+ra2uξ+21-(1-ck2)uξ+10+k2Fξ+11.
Assum that μ=1+ck/2+ra2, and ω=2+bk2-ra2, then the matrix form of (2.12) and (2.13) will be
(2.14)[μ-ra2-ra2μ][uξ2uξ+12]=[ω00ω][uξ1uξ+11]+[ra2uξ-11-(1-ck2)uξ0+k2Fξ1ra2uξ+21-(1-ck2)uξ+10+k2Fξ+11],
then, we have
(2.15)[uξ2uξ+12]=1μ2-r2a4{[μωra2ωra2ωμω][uξ1uξ+11]+[μra2ra2μ][ra2uξ-11-(1-ck2)uξ0+k2Fξ1ra2uξ+21-(1-ck2)uξ+10+k2Fξ+11]},
or, we have
(2.16)[uξ2uξ+12]=1μ2-r2a4×[ra2μuξ-11+μωuξ1+ra2ωuξ+11+r2a4uξ+21-μ(1-ck2)uξ0-ra2(1-ck2)uξ+10+μk2Fξ1+ra2k2Fξ+11r2a4uξ-11+ra2ωuξ1+μωuξ+11+ra2μuξ+21-ra2(1-ck2)uξ0-μ(1-ck2)uξ+10+ra2k2Fξ1+μk2Fξ+11]
which can be diagrammatically represented by the computational molecule given in Figure 2.The computational molecule for node (xξ,t2) (a) and node (xξ+1,t2) (b).
(a)(b)But, in case 2, (j≥2), for j=2,4,…, the solution at the nodes i=0,1,…,ξ is approximated by (2.4), and the solution at the nodes i=ξ+1,ξ+2,…,M is approximated by (2.6), where j=2. Similary, at the node (xξ,tj+1), the solution is approximated by (2.6),
(2.17)(1+ck2+ra2)uξj+1-ra2uξ+1j+1=ra2uξ-1j+(2+bk2-ra2)uξj-(1-ck2)uξj-1+k2Fξj,
while at node (xξ+1,tj+1) the solution is approximated by (2.4),
(2.18)-ra2uξj+1+(1+ck2+ra2)uξ+1j+1=(2+bk2-ra2)uξ+1j+ra2uξ+2j-(1-ck2)uξ+1j-1+k2Fξ+1j.
If (2.17) and (2.18) are rewritten in matrix form:
(2.19)[μ-ra2-ra2μ][uξj+1uξ+1j+1]=[ω00ω][uξjuξ+1j]+[ra2uξ-1j-(1-ck2)uξj-1+k2Fξjra2uξ+2j-(1-ck2)uξ+1j-1+k2Fξ+1j],
whose
(2×2) matrix of coefficient can easily be inverted, the equation can be written in explicit form as
(2.20)[uξj+1uξ+1j+1]=1μ2-r2a4{[μωra2ωra2ωμω][uξjuξ+1j]+[μra2ra2μ][ra2uξ-1j-(1-ck2)uξj-1+k2Fξjra2uξ+2j-(1-ck2)uξ+1j-1+k2Fξ+1j]}.
This simplifies to
(2.21)[uξj+1uξ+1j+1]=1μ2-r2a4×[ra2μuξ-1j+μωuξj+ra2ωuξ+1j+r2a4uξ+2j-μ(1-ck2)uξj-1-ra2(1-ck2)uξ+1j-1+μk2Fξj+ra2k2Fξ+1jr2a4uξ-1j+ra2ωuξj+μωuξ+1j+ra2μuξ+2j-ra2(1-ck2)uξj-1-μ(1-ck2)uξ+1j-1+ra2k2Fξj+μk2Fξ+1j]
which similarity can be diagrammatically represented by the computational molecule given in Figure 2. The flow chart of this method is displayed in Figure 3. We use ▸ to denote the scheme (2.4), ◂ to denote the scheme (2.6) and ⊳ and ⊲ to denote the scheme (2.16).Figure 3
The diagram of the parallel implementation of asymmetric Saul’yev schemes for two sub-domains.Implementation of this algorithm for more than two sub-domains is straightforward as follows (see Figure4). For more than two sub-domains, we have :Figure 4
The diagram of the parallel implementation of asymmetric Saul’yev schemes.(i)
Division of spatial computational domain to desired number of sub-domains, based on load balancing constraint.(ii)
Indication of sweeping direction for each sub-domain. Sweeping direction of each sub-domain must be in opposite direction of its neighbors. For example, we use LR (left to right) direction for odd sub-domains and RL (right to left) direction for even sub-domains. This sweeping direction is inversed after each time step.(iii)
Updating start node of each sub-domain with (2.13) and remained nodes with (2.5) or (2.6). If start node is located at physical boundaries, the prescribed boundary value is used (without using (2.13).


## 3. Stability Analysis and Convergence

### 3.1. Stability Analysis

In this section we prove unconditional stability of presented scheme by the matrix method. Pay attention that in the best case, if the number of the sub-domains is reduced to one   sub-domain, the presented scheme converted to asymmetric Saul'yev schemes that are unconditionally stable (see [9]) and in the worst case, if the number of sub-domains increased such that each sub-domain only has two nodes the presented scheme is converted to alternating group explicit method that is unconditionally stable (see [6]).To prove the stability of presented scheme, suppose that the spatial computational domain is divided to number ofκ sub-domains (κ is even for simplicity) and each sub-domain contains number of λ nodes. Also the first sub-domain contains the left Neumann boundary node and the last sub-domain contains the right Neumann boundary node as additional nodes. Similarly, the presented scheme can be expressed in two cases. In the first case, considering j=0,1, then we have :
(3.1)(2I+ra2G1)U1=((2+bk2)I-ra2G2)U0+B1,((1+ck2)I+ra2G2)U2=((2+bk2)I-ra2G1)U1+B2,
where Uj=(u0j,u1j,…,uM-1j,uMj)T, j=0,1,2,…, and
(3.2)B1=(2k(1-ck2)g(x0)+k2F(x0,t0)-2ra2hα(t1)2k(1-ck2)g(x1)+k2F(x1,t0)2k(1-ck2)g(x2)+k2F(x2,t0)⋮2k(1-ck2)g(xM-2)+k2F(xM-2,t0)2k(1-ck2)g(xM-1)+k2F(xM-1,t0)2k(1-ck2)g(xM)+k2F(xM,t0)+2ra2hβ(t1)),B2=(-(1-ck2)u00+k2F(x0,t1)-2ra2hα(t1)-(1-ck2)u10+k2F(x1,t1)-(1-ck2)u20+k2F(x2,t1)⋮-(1-ck2)uM-20+k2F(xM-2,t1)-(1-ck2)uM-10+k2F(xM-1,t1)-(1-ck2)uM0+k2F(xM,t1)+2ra2hβ(t1)),
and in the second case, j=2,4,…,
(3.3)((1+ck2)I+ra2G1)Uj+1=((2+bk2)I-ra2G2)Uj+Bj+1,((1+ck2)I+ra2G2)Uj+2=((2+bk2)I-ra2G1)Uj+1+Bj+2,
where
(3.4)Bj+1=(-(1-ck2)u0j-2+k2F(x0,tj)-2ra2hα(tj+1)-(1-ck2)u1j-2+k2F(x1,tj)-(1-ck2)u2j-2+k2F(x2,tj)⋮-(1-ck2)uM-2j-2+k2F(xM-2,tj)-(1-ck2)uM-1j-2+k2F(xM-1,tj)-(1-ck2)uMj-2+k2F(xM,tj)+2ra2hβ(tj+1)),Bj+2=(-(1-ck2)u0j-1+k2F(x0,tj+1)-2ra2hα(tj+1)-(1-ck2)u1j-1+k2F(x1,tj+1)-(1-ck2)u2j-1+k2F(x2,tj+1)⋮-(1-ck2)uM-2j-1+k2F(xM-2,tj+1)-(1-ck2)uM-1j-1+k2F(xM-1,tj+1)-(1-ck2)uMj-1+k2F(xM,tj+1)+2ra2hβ(tj+1)),
In both cases, G1 and G2 are(3.5)G1=(1-1-110⋱⋱⋱-11001-1⋱⋱⋱01-101-1-110⋱⋱⋱⋱⋱⋱⋱⋱⋱-11001-1⋱⋱⋱01-1-11)(M+1)×(M+1)G2=(1-101-1⋱⋱⋱01-1-110⋱⋱⋱-110-11001-101-1⋱⋱⋱⋱⋱⋱⋱⋱⋱01-1-110⋱⋱⋱-110-11)(M+1)×(M+1)
assume that
(3.6)Pκ=(-11-11⋱-11-11)λ×(λ+1),Qκ=(1-11-1⋱1-11-1)λ×(λ+1),
then we have
(3.7)G1=(1-1P1Q2⋱Pκ-1Qκ-11),G2=(1-10Q1P2⋱Qκ-1Pκ0-11),
By eliminating U1 from (3.1) we obtain
(3.8)U2=Φ̂U0+Υ̂,
where
(3.9)Φ̂=((1+ck2)I+ra2G2)-1((2+bk2)I-ra2G1)(2I+ra2G1)-1((2+bk2)I-ra2G2),Υ̂=((1+ck2)I+ra2G2)-1((2+bk2)I-ra2G1)(2I+ra2G1)-1B1+((1+ck2)I+ra2G2)-1B2.
similarly, by eliminating Uj+1 from (3.3) we obtain
(3.10)Uj+2=ΦUj+Υ,j=2,4,…,
where
(3.11)Φ=((1+ck2)I+ra2G2)-1((2+bk2)I-ra2G1)((1+ck2)I+ra2G1)-1×((2+bk2)I-ra2G2),Υ=((1+ck2)I+ra2G2)-1((2+bk2)I-ra2G1)((1+ck2)I+ra2G1)-1Bj+1+((1+ck2)I+ra2G2)-1Bj+2.
By Lax-Richtmyer condition of stability, the below inequality must be satisfy
(3.12)∥Φ̂∥≤1,∥Φ∥≤1.
then we have
(3.13)∥Φ̂∥≤∥((1+ck2)I+ra2G2)-1∥·∥(2+bk2)I-ra2G1∥·∥(2I+ra2G1)-1∥·∥(2+bk2)I-ra2G2∥,∥Φ∥≤∥((1+ck2)I+ra2G2)-1∥·∥(2+bk2)I-ra2G1∥·∥((1+ck2)I+ra2G1)-1∥·∥(2+bk2)I-ra2G2∥.
Note that the matrixs G1 and G2 are singular, also must be M∣κ, then for M and κ, the eignvalues of G1 and G2 are
(3.14)eig(G1)={1,1,1,…,1︸M-κ-2,2,0,2,0,…,2,0},eig(G2)={1,1,1,…,1︸M-κ,2,0,2,0,…,2,0}.
It easy to see that [3, 4], the Φ̂ and Φ are nonsingular and for all r>0(3.15)∥Φ̂∥≤1,∥Φ∥≤1.
Then the following theorem satisfies for Telegraph equation.Theorem 3.1.
The parallel implementation of asymmetric Saul'yev schemes for the solution of Telegraph equation (2.1) is unconditionally stable for all r>0.


### 3.2. Convergence

The real important of the concept of consistency lies in a theorem by Lax [10], which state that If a linear finite-difference equation is consistent with a properly posed linear initial value problem then stability guarantees convergence ofutoUas the mesh lengths tend to zero, whereuandUare numerical and exactsolutionsof problem. Consistency can be defined in either of two equivalent but slightly different ways.The more general definition is as follows [11, pages 40–45]. Let L(U)=0, represent the partial differential equation in the independent variables
x and t, with exact solution U, and F(u)=0 represent the approximation finite difference equation with exact solution u. Let v be a continuous function of x and t with a sufficient number of continues derivatives to enable L(v) to be evaluated at the point (ih,jk).Then the truncation errorTij(v) at the point (ih,jk) is defined by
(3.16)Tij(v)=Fij(v)-L(vij).IfTij(v)→0ash→0,k→0,then the difference equation is said to be consistent or compatible with the partial differential equation. With this definition Tij gives an indication of the error resulting from the replacement of L(vij) by Fij(v).Most authors putv=U, because L(U)=0. It then follows that
(3.17)Tij(U)=Fij(U),
and the truncation error coincides with the local truncation error. The difference equation is then consistent if the limiting value of the local truncation error is zero ash→0,k→0. We now illustrate the convergence of covering difference equation (3.10). Substituting ϵmn=Umn-umn into (3.10), we get
(3.18)ϵn+2=Φϵn+Tn,n=2,4,…,
where ϵn=[ϵ1nϵ2n…ϵM-1n]T, and Tn=[T1nT2n…TM-1n]T, is the truncation error of (2.3), and for m=1,2,…,M-1, we have Tmn=O(k2+h2). Then
(3.19)∥ϵn+2∥≤∥Φ∥∥ϵn∥+∥Tn∥,n=2,4,…≤∥Φ∥n/2+1∥ϵ0∥+∑s=1n/2‍∥Φ∥s-1∥T2s-2∥,≤∥Φ∥n/2+1∥ϵ0∥+I-∥Φ∥n/2I-∥Φ∥max1≤s≤n∥T2s-2∥.
Using Theorem 3.1, we have ∥Φ∥≤1, so we have
(3.20)∥ϵn+2∥≤∥ϵ0∥+11-∥Φ∥max1≤s≤n∥T2s-2∥,
where
(3.21)max1≤s≤n∥T2s-2∥=M(k2+h2),
and M is a constant. Therefore
(3.22)∥ϵn+2∥≤∥ϵ0∥+M(k2+h2).
From this result we deduce unconditional convergence as k→0 and h→0.


## 3.1. Stability Analysis

In this section we prove unconditional stability of presented scheme by the matrix method. Pay attention that in the best case, if the number of the sub-domains is reduced to one   sub-domain, the presented scheme converted to asymmetric Saul'yev schemes that are unconditionally stable (see [9]) and in the worst case, if the number of sub-domains increased such that each sub-domain only has two nodes the presented scheme is converted to alternating group explicit method that is unconditionally stable (see [6]).To prove the stability of presented scheme, suppose that the spatial computational domain is divided to number ofκ sub-domains (κ is even for simplicity) and each sub-domain contains number of λ nodes. Also the first sub-domain contains the left Neumann boundary node and the last sub-domain contains the right Neumann boundary node as additional nodes. Similarly, the presented scheme can be expressed in two cases. In the first case, considering j=0,1, then we have :
(3.1)(2I+ra2G1)U1=((2+bk2)I-ra2G2)U0+B1,((1+ck2)I+ra2G2)U2=((2+bk2)I-ra2G1)U1+B2,
where Uj=(u0j,u1j,…,uM-1j,uMj)T, j=0,1,2,…, and
(3.2)B1=(2k(1-ck2)g(x0)+k2F(x0,t0)-2ra2hα(t1)2k(1-ck2)g(x1)+k2F(x1,t0)2k(1-ck2)g(x2)+k2F(x2,t0)⋮2k(1-ck2)g(xM-2)+k2F(xM-2,t0)2k(1-ck2)g(xM-1)+k2F(xM-1,t0)2k(1-ck2)g(xM)+k2F(xM,t0)+2ra2hβ(t1)),B2=(-(1-ck2)u00+k2F(x0,t1)-2ra2hα(t1)-(1-ck2)u10+k2F(x1,t1)-(1-ck2)u20+k2F(x2,t1)⋮-(1-ck2)uM-20+k2F(xM-2,t1)-(1-ck2)uM-10+k2F(xM-1,t1)-(1-ck2)uM0+k2F(xM,t1)+2ra2hβ(t1)),
and in the second case, j=2,4,…,
(3.3)((1+ck2)I+ra2G1)Uj+1=((2+bk2)I-ra2G2)Uj+Bj+1,((1+ck2)I+ra2G2)Uj+2=((2+bk2)I-ra2G1)Uj+1+Bj+2,
where
(3.4)Bj+1=(-(1-ck2)u0j-2+k2F(x0,tj)-2ra2hα(tj+1)-(1-ck2)u1j-2+k2F(x1,tj)-(1-ck2)u2j-2+k2F(x2,tj)⋮-(1-ck2)uM-2j-2+k2F(xM-2,tj)-(1-ck2)uM-1j-2+k2F(xM-1,tj)-(1-ck2)uMj-2+k2F(xM,tj)+2ra2hβ(tj+1)),Bj+2=(-(1-ck2)u0j-1+k2F(x0,tj+1)-2ra2hα(tj+1)-(1-ck2)u1j-1+k2F(x1,tj+1)-(1-ck2)u2j-1+k2F(x2,tj+1)⋮-(1-ck2)uM-2j-1+k2F(xM-2,tj+1)-(1-ck2)uM-1j-1+k2F(xM-1,tj+1)-(1-ck2)uMj-1+k2F(xM,tj+1)+2ra2hβ(tj+1)),
In both cases, G1 and G2 are(3.5)G1=(1-1-110⋱⋱⋱-11001-1⋱⋱⋱01-101-1-110⋱⋱⋱⋱⋱⋱⋱⋱⋱-11001-1⋱⋱⋱01-1-11)(M+1)×(M+1)G2=(1-101-1⋱⋱⋱01-1-110⋱⋱⋱-110-11001-101-1⋱⋱⋱⋱⋱⋱⋱⋱⋱01-1-110⋱⋱⋱-110-11)(M+1)×(M+1)
assume that
(3.6)Pκ=(-11-11⋱-11-11)λ×(λ+1),Qκ=(1-11-1⋱1-11-1)λ×(λ+1),
then we have
(3.7)G1=(1-1P1Q2⋱Pκ-1Qκ-11),G2=(1-10Q1P2⋱Qκ-1Pκ0-11),
By eliminating U1 from (3.1) we obtain
(3.8)U2=Φ̂U0+Υ̂,
where
(3.9)Φ̂=((1+ck2)I+ra2G2)-1((2+bk2)I-ra2G1)(2I+ra2G1)-1((2+bk2)I-ra2G2),Υ̂=((1+ck2)I+ra2G2)-1((2+bk2)I-ra2G1)(2I+ra2G1)-1B1+((1+ck2)I+ra2G2)-1B2.
similarly, by eliminating Uj+1 from (3.3) we obtain
(3.10)Uj+2=ΦUj+Υ,j=2,4,…,
where
(3.11)Φ=((1+ck2)I+ra2G2)-1((2+bk2)I-ra2G1)((1+ck2)I+ra2G1)-1×((2+bk2)I-ra2G2),Υ=((1+ck2)I+ra2G2)-1((2+bk2)I-ra2G1)((1+ck2)I+ra2G1)-1Bj+1+((1+ck2)I+ra2G2)-1Bj+2.
By Lax-Richtmyer condition of stability, the below inequality must be satisfy
(3.12)∥Φ̂∥≤1,∥Φ∥≤1.
then we have
(3.13)∥Φ̂∥≤∥((1+ck2)I+ra2G2)-1∥·∥(2+bk2)I-ra2G1∥·∥(2I+ra2G1)-1∥·∥(2+bk2)I-ra2G2∥,∥Φ∥≤∥((1+ck2)I+ra2G2)-1∥·∥(2+bk2)I-ra2G1∥·∥((1+ck2)I+ra2G1)-1∥·∥(2+bk2)I-ra2G2∥.
Note that the matrixs G1 and G2 are singular, also must be M∣κ, then for M and κ, the eignvalues of G1 and G2 are
(3.14)eig(G1)={1,1,1,…,1︸M-κ-2,2,0,2,0,…,2,0},eig(G2)={1,1,1,…,1︸M-κ,2,0,2,0,…,2,0}.
It easy to see that [3, 4], the Φ̂ and Φ are nonsingular and for all r>0(3.15)∥Φ̂∥≤1,∥Φ∥≤1.
Then the following theorem satisfies for Telegraph equation.Theorem 3.1.
The parallel implementation of asymmetric Saul'yev schemes for the solution of Telegraph equation (2.1) is unconditionally stable for all r>0.


## 3.2. Convergence

The real important of the concept of consistency lies in a theorem by Lax [10], which state that If a linear finite-difference equation is consistent with a properly posed linear initial value problem then stability guarantees convergence ofutoUas the mesh lengths tend to zero, whereuandUare numerical and exactsolutionsof problem. Consistency can be defined in either of two equivalent but slightly different ways.The more general definition is as follows [11, pages 40–45]. Let L(U)=0, represent the partial differential equation in the independent variables
x and t, with exact solution U, and F(u)=0 represent the approximation finite difference equation with exact solution u. Let v be a continuous function of x and t with a sufficient number of continues derivatives to enable L(v) to be evaluated at the point (ih,jk).Then the truncation errorTij(v) at the point (ih,jk) is defined by
(3.16)Tij(v)=Fij(v)-L(vij).IfTij(v)→0ash→0,k→0,then the difference equation is said to be consistent or compatible with the partial differential equation. With this definition Tij gives an indication of the error resulting from the replacement of L(vij) by Fij(v).Most authors putv=U, because L(U)=0. It then follows that
(3.17)Tij(U)=Fij(U),
and the truncation error coincides with the local truncation error. The difference equation is then consistent if the limiting value of the local truncation error is zero ash→0,k→0. We now illustrate the convergence of covering difference equation (3.10). Substituting ϵmn=Umn-umn into (3.10), we get
(3.18)ϵn+2=Φϵn+Tn,n=2,4,…,
where ϵn=[ϵ1nϵ2n…ϵM-1n]T, and Tn=[T1nT2n…TM-1n]T, is the truncation error of (2.3), and for m=1,2,…,M-1, we have Tmn=O(k2+h2). Then
(3.19)∥ϵn+2∥≤∥Φ∥∥ϵn∥+∥Tn∥,n=2,4,…≤∥Φ∥n/2+1∥ϵ0∥+∑s=1n/2‍∥Φ∥s-1∥T2s-2∥,≤∥Φ∥n/2+1∥ϵ0∥+I-∥Φ∥n/2I-∥Φ∥max1≤s≤n∥T2s-2∥.
Using Theorem 3.1, we have ∥Φ∥≤1, so we have
(3.20)∥ϵn+2∥≤∥ϵ0∥+11-∥Φ∥max1≤s≤n∥T2s-2∥,
where
(3.21)max1≤s≤n∥T2s-2∥=M(k2+h2),
and M is a constant. Therefore
(3.22)∥ϵn+2∥≤∥ϵ0∥+M(k2+h2).
From this result we deduce unconditional convergence as k→0 and h→0.


## 4. Numerical Experiments

Example 4.1.
In this section the accuracy and stability of presented parallel scheme is numerically investigated by solution of a test problem. Consider the following problem:(4.1)∂2u∂t2+3∂u∂t=∂2u∂x2-u,(x,t)∈[0,1]×[0,1],u(x,0)=-sin(x),x∈[0,1],ut(x,0)=sin(x),x∈[0,1],ux(0,t)=-e-t,t∈[0,1],ux(1,t)=-e-tcos1,t∈[0,1].
The exact solution is given by
(4.2)u(x,t)=-e-tsin(x).
The average of absolute errors of this problem for M=N=96, and κ=2,4,6,8, is given in Table 1. The figures of this data, that is, |∑j=1N‍uexact(x,tj)-uapprox(x,tj)|/N, are shown in Figure 5.Table 1
The average of absolute errors of test Example4.1.
xNumber of sub-domain24682191.88767059e-0042.01585278e-0042.20992690e-0042.42050716e-0044192.17086802e-0042.30713785e-0042.17202537e-0042.39552996e-0046192.57232399e-0042.39418130e-0042.62471306e-0042.54550100e-0048193.03646182e-0042.90480449e-0042.85129788e-0042.77156840e-00410193.21312054e-0043.15940605e-0043.10202033e-0043.06178032e-00412193.68769542e-0043.69093725e-0043.64261355e-0043.62145246e-00414194.10539234e-0044.13602591e-0043.86283671e-0043.83301138e-00416194.42670434e-0044.22685287e-0044.22098694e-0043.97233016e-00418194.61241511e-0044.42709395e-0044.20277475e-0043.97016774e-00414.64209163e-0044.45913196e-0044.23724691e-0044.00709711e-004(a)The number of sub−domain=2. (b) The number of sub−domain =4. (c) The number of sub−domain =6. (d) The number of sub−domain =8.
(a)(b)(c)(d)


## 5. Conclusion

In this paper, we have discussed a new three-level implicit parallel difference scheme ofO(k2+h2) for the different  solution of the Telegraph equation (2.1). The scheme is stable for r>0, and applicable to singular equations. It is hoped that the new idea presented in this paper will be useful for the development of an unconditionally stable parallel difference formula of O(k4+h4) for the numerical solution of the Telegraph equation (2.1).



---

*Source: 969610-2009-06-14.xml*