# The Spectral Method for Solving Sine-Gordon Equation Using a New Orthogonal Polynomial

Authors: Zoleikha Soori; Azim Aminataei
Journal: ISRN Applied Mathematics
 (2012)

---

## Abstract

We have presented a numerical scheme for solving one-dimensional nonlinear sine-Gordon equation. We apply the
spectral method with a basis of a new orthogonal polynomial which is orthogonal over the interval[0,1] with weighting function
one. The results show the accuracy and efficiency of the proposed method.


---

## Body

## 1. Introduction

The sine-Gordon equation is a nonlinear hyperbolic partial differential equation (PDE). It was originally considered in the nineteenth century in the course of surface of constant negative curvature. This equation attracted a lot of attention in the 1970s due to the presence of soliton solutions [1, 2]. The sine-Gordon equation appears in a number of physical applications [3–5], including applications in the chain of coupled pendulums and modelling the propagation of transverse electromagnetic (TEM) wave on a superconductor transmission system. Consider the one-dimensional nonlinear sine-Gordon equation(1.1)∂2u(x,t)∂t2=∂2u(x,t)∂x2-sin⁡⁡(u(x,t)),(x,t)∈[0,1]×[t0,T],
with the initial conditions(1.2)u(x,t0)=f0(x),∂u∂t(x,t0)=f1(x),x∈[0,1],
and the boundary conditions(1.3)u(0,t)=g0(t),u(1,t)=g1(t),t≥t0.
In the last decade, several numerical schemes have been developed for solving (1.1), for instance, high-order solution of one-dimensional sine-Gordon equation using compact finite difference and DIRKN methods [3]. The authors of [6] proposed a numerical scheme for solving (1.1) using collocation and radial basis functions. Also, the boundary integral equation approach is used in [7]. Bratsos has proposed a numerical method for solving one-dimensional sine-Gordon equation and a third-order numerical scheme for the two-dimensional sine-Gordon equation in [8, 9], respectively. Also, in [10], a numerical scheme using radial basis function for the solution of two-dimensional sine-Gordon equation is used. In addition, several authors proposed spectral methods and Fourier pseudospectral method for solving nonlinear wave equation using a discrete Fourier series and Chebyshev orthogonal polynomials [11–13].In this paper, we have proposed a numerical scheme for solving (1.1) using spectral method with a basis of a new orthogonal polynomial. This polynomial is introduced by Chelyshkov in [14].The outline of this paper is as follows. In Section2, we introduce a new orthogonal polynomial. In Section 3, we apply spectral method for solving (1.1). In Section 4, we show the numerical results. Finally, a brief conclusion is given in Section 5.


## 2. The New Orthogonal Polynomial

In this section, we briefly introduce the new orthogonal polynomial, that is introduced in [14].Letn be a fixed whole number, and Pn is a sequence of polynomial, that is,(2.1)Pn={Pnk}k=n0,Pnk≡Pnk(x)=∑l=knτnklxl,
that satisfy the orthogonality relationship(2.2)∫01PnkPnldx={0,k≠l,1(k+l+1),k=l,
and standardization(2.3)sign⁡(τnkn)=(-1)n-k.
The coefficient τnkl of the polynomial Pnk is defined by requirements (2.2), (2.3), and the Gram-Schmidt orthogonalization procedure (without normalization).The explicit definition of the polynomial Pnk is as follows:(2.4)Pnk(x)=∑j=0n-k(-1)j(n-kj)(n+k+1+jn-k)xk+j,k=0,1,…,n.
This yields Rodrigue’s type representation as follows:(2.5)Pnk(x)=1(n-k)!1xk+1dn-kdxn-k(xn+k+1(1-x)n-k),k=0,1,…,n,and the orthogonality relationship (2.2) is confirmed by applying (2.5). It also follows from (2.5) that(2.6)∫01Pnk(x)dx=∫01xndx=1n+1.This polynomial can be connected to a fixed set of the Jacobi polynomials Pn(α,β)(ξ) [15], that is,(2.7)Pnk=xnPk-n(2n,0)(1-2x).


## 3. The Proposed Method

In this section, the introduced orthogonal polynomial is applied as the basis for spectral method to solve sine-Gordon equation.The approximate solutionuN(x,t) to the exact solution u(x,t) can be written in the form(3.1)u(x,t)≅UN(x,t)=∑n=0Nan(t)Pn0(x),wherean is the time-dependent quantities which must be determined. In this paper, we apply orthogonal polynomial Pnk for k=0. From (3.1), we get(3.2)UttN=∑n=0Nän(t)Pn0(x),UxxN=∑n=0Nan(t)P̈n0(x).Using  (3.2) and substituting in (1.1), we obtain(3.3)∑n=0Nän(t)Pn0(xi)=∑n=0Nan(t)P̈n0(xi)-sin⁡⁡(UN(xi,t)),i=1,2,…,N-1,whereän and P̈n0 denote the second-order derivatives of an and Pn0 with respect to t and x, respectively, wherein the collocation points are as xi=i/N,i=0,1,…,N.Equation (13) is rewritten as follows:(3.4)MÄ(t)=KA(t)-S1.For displaying this equation in matrix form, fori=1,2,…,N-1, and n=0,1,…,N, we define
(3.5)A(t)=[a0(t),a1(t),…,aN(t)]T,M=Min=Pn0(xi),K=Kin=P̈n0(xi),S1=[sin⁡⁡(U1),sin⁡⁡(U2),…,sin⁡⁡(UN-1)]T,Ui=UN(xi,t).The matricesM and K are positive definite matrices of order (N-1)×(N-1). The system given in (3.4) consists of (N-1) equations in (N+1) unknowns. To obtain a unique solution for this system, we need two additional equations. For this, we add two boundary conditions, which are given by (1.3), in the following forms:(3.6)UN(0,t)=∑n=0Nan(t)Pn0(0)=g0(t),t≥t0,UN(1,t)=∑n=0Nan(t)Pn0(1)=g1(t),t≥t0.From (3.6), we get(3.7)∑n=0Nän(t)Pn0(0)=d2g0(t)dt2,t≥t0,∑n=0Nän(t)Pn0(1)=d2g1(t)dt2,t≥t0.Therefore, by using (3.6), and (3.7), a new system is obtained as follows:(3.8)MÄ(t)=F(t,A(t)),F(t,A(t))=KA(t)-S+G(t),whereA(t), F(t,A(t)), S, and G(t) are vectors of order (N+1) with (N+1) components. The matrices M and K are of order (N+1)×(N+1), which are defined as follows:(3.9)M=(p00(0)p10(0)⋯pN0(0)p00(x1)p10(x1)⋯pN0(x1)⋮⋮⋱⋮p00(xN-1)p10(xN-1)⋯pN0(xN-1)p00(1)p10(1)⋯pN0(1)),(3.10)K=(00⋯0p̈00(x1)p̈10(x1)⋯p̈N0(x1)⋮⋮⋱⋮p̈00(xN-1)p̈10(xN-1)⋯p̈N0(xN-1)00⋯0).
Also, the vectors S and G(t) are defined as follows:(3.11)S=[0,sin⁡⁡(U1),sin⁡⁡(U2),…,sin⁡⁡(UN-1),0]T,G(t)=[d2g0(t)dt2,0,…,0,d2g1(t)dt2]T.
The nonlinear second-order system of ODEs (3.8) can be solved numerically using the fourth-order A-stable DIRKN method [3, 16]. From the initial conditions (1.2), we determine the initial vectors A0 and Ȧ0. Using (3.1), we get(3.12)UN(xi,0)=∑n=0Nan(t0)Pn0(xi),i=0,1,…,N,∂UN∂t(xi,0)=∑n=0Nȧn(t0)Pn0(xi),i=0,1,…,N.Equation (3.12) is system of (N+1) equations in (N+1) unknowns. These equations can be written in the matrix form MA0=b1 and MȦ0=b2, where
(3.13)b1=[f0(x0),…,f0(xN)]T,b2=[f1(x0),…,f1(xN)]T.


### 3.1. The DIRKN Method

For solving (3.8), we use DIRKN method [3, 16]:(3.14)MQi=MAn+ciΔtMȦn+Δt2∑j=1saijF(tn+cjΔt,Qj),i=1,…,s,MAn+1=MAn+ciΔtMȦn+Δt2∑j=1sb¯iF(tn+cjΔt,Qj),MȦn+1=MȦn+Δt∑j=1sbiF(tn+cjΔt,Qj),whereQi is the internal stages, An+1 and Ȧn+1 represent approximations to A(tn+1) and Ȧ(tn+1), respectively. Also we define(3.15)Δt=t0+jΔt,j=0,1,…,M,where Δt=(T-t0)/M. The RKN method can be denoted in Butcher,s notation by the table of coefficients:(3.16)c1a11⋯a1s⋮⋮⋱⋮csas1⋯assb̅1⋯b̅sb1⋯bsIn the DIRKN method,aij=0 as (i<j) and aii are equal. Also, in the above table of coefficients, s is the number of stages. We use the four-stage A-stable DIRKN method [16], with algebraic fourth order, which is defined by the following table of coefficients for solving all experiments.Using Table1, we can calculate Qi, i=1,…,4, separately from (3.14) with Newton-Raphson iteration method. We iterate until the following criterion is satisfied:(3.17)‖(Qik+1-Qik)‖∞<10-12,
where Qik is the value of Qi at the kth iteration of the Newton-Raphson method.Table 1
Four-stage A-stable DIRKN method.
0.65820.240.93000.1864849679520.240.0700−3.1486640861952.9456280221830.240.34181.491936101409−1.8163048528350.1030070016890.240.1085722379060.0127646089720.1695869477700.2090762053520.3176484432580.1823515567420.1823515567420.317648443258


## 3.1. The DIRKN Method

For solving (3.8), we use DIRKN method [3, 16]:(3.14)MQi=MAn+ciΔtMȦn+Δt2∑j=1saijF(tn+cjΔt,Qj),i=1,…,s,MAn+1=MAn+ciΔtMȦn+Δt2∑j=1sb¯iF(tn+cjΔt,Qj),MȦn+1=MȦn+Δt∑j=1sbiF(tn+cjΔt,Qj),whereQi is the internal stages, An+1 and Ȧn+1 represent approximations to A(tn+1) and Ȧ(tn+1), respectively. Also we define(3.15)Δt=t0+jΔt,j=0,1,…,M,where Δt=(T-t0)/M. The RKN method can be denoted in Butcher,s notation by the table of coefficients:(3.16)c1a11⋯a1s⋮⋮⋱⋮csas1⋯assb̅1⋯b̅sb1⋯bsIn the DIRKN method,aij=0 as (i<j) and aii are equal. Also, in the above table of coefficients, s is the number of stages. We use the four-stage A-stable DIRKN method [16], with algebraic fourth order, which is defined by the following table of coefficients for solving all experiments.Using Table1, we can calculate Qi, i=1,…,4, separately from (3.14) with Newton-Raphson iteration method. We iterate until the following criterion is satisfied:(3.17)‖(Qik+1-Qik)‖∞<10-12,
where Qik is the value of Qi at the kth iteration of the Newton-Raphson method.Table 1
Four-stage A-stable DIRKN method.
0.65820.240.93000.1864849679520.240.0700−3.1486640861952.9456280221830.240.34181.491936101409−1.8163048528350.1030070016890.240.1085722379060.0127646089720.1695869477700.2090762053520.3176484432580.1823515567420.1823515567420.317648443258


## 4. Numerical Experiments

In this section, we present the results of numerical experiments using the method introduced in Section3. The L2 and L∞ error norms are defined as follows:(4.1)L2=‖u-UN‖=[∑i=0N(ui-UiN)2]1/2,L∞=‖u-UN‖=max⁡0≤i≤N|ui-UiN|.We chooseΔt=1/1000, for solving all experiments.The computations associated with the experiments discussed above were performed in Maple 13 on a PC with a CPU of 2.4 GHZ.Experiment 1.
In this experiment [3, 6, 10], we consider (1.1) with t0=0 and the initial conditions
(4.2)f0(x)=0,f1(x)=4sech⁡(x),
with the boundary conditions (1.3), which can be obtained from the following exact solution:
(4.3)u(x,t)=4arctan(sech⁡(x)t).
In Table 2, the approximate solution for several final times (T) with different locations is shown. Figure 1 shows two type plots at t=0.5 of Experiment 1.Table 2
The approximate solution forN=16 of Experiment 1.
xT=0.1T=0.2T=0.3T=0.4T=0.50.00.39867461000.78958223941.16582717791.52202550841.85459043600.20.39088221000.77443854421.14413490521.49481830701.82299156440.40.36895306920.73117328291.08276411171.41753016381.73281295130.60.33662333890.66854503230.99144590401.30169712441.59654466890.80.29852448870.59376009410.88263028041.16245005681.43104813691.00.25885973290.51556922810.76808352301.01455564391.2534069046(a) The exact and approximate solutions forN=16 are shown. (b) The L∞ error norm for N=16 is shown.
(a)(b)Experiment 2.
In this experiment [3, 8], we consider (1.1) with t0=0 and the initial conditions
(4.4)f0(x)=0,(4.5)f1(x)=41+c2sech⁡(x1+c2),
with the boundary conditions (1.3), which can be obtained from the following exact solution:
(4.6)u(x,t)=4arctan(c-1sin⁡(ct1+c2)sech(x1+c2)),c=1.
In Tables3 and 4, we show the L2 and L∞ error norms, for several final times (T). As we see from Tables 3 and 4, for the given values of N, the L∞ error norm is less than the L2 error norm. Figure 2 shows four type plots at t=0.5 of Experiment 2.Table 3
TheL2 error norm for different values of N of Experiment 2.
TN=4N=8N=160.15.02×10-69.11×10-99.56×10-150.25.37×10-51.00×10-81.35×10-130.32.05×10-43.67×10-72.22×10-130.45.00×10-48.37×10-69.10×10-130.59.46×10-41.29×10-63.43×10-120.61.50×10-31.52×10-61.75×10-110.72.10×10-31.58×10-61.29×10-100.82.65×10-31.65×10-68.08×10-100.93.08×10-31.94×10-64.94×10-91.03.37×10-32.13×10-62.26×10-8Table 4
TheL∞ error norm for different values of N of Experiment 2.
TN=4N=8N=160.11.30×10-62.37×10-92.68×10-150.29.79×10-61.83×10-82.23×10-140.32.98×10-55.36×10-82.97×10-140.46.10×10-59.14×10-82.31×10-130.59.82×10-51.00×10-77.94×10-130.61.33×10-41.00×10-72.38×10-120.71.57×10-41.00×10-71.29×10-110.81.63×10-41.00×10-71.14×10-100.91.63×10-41.16×10-78.03×10-101.01.80×10-41.16×10-74.59×10-9(a) The exact and approximate solutions forN=16 are shown. (b) The L∞ error norm for N=16 is shown. (c) The L2 error norm for different values of N is shown. (d) The L∞ error norm for different values of N is shown.
(a)(b)(c)(d)Experiment 3.
In experiment [3, 9], we consider  (1.1) with t0=0 and the initial conditions
(4.7)f0(x)=4arctan(csinh⁡(x1-c2)),(4.8)f1(x)=0,
with the boundary conditions (1.3), which can be obtained from the following exact solution:
(4.9)u(x,t)=4arctan(csinh⁡(x/1-c2)cosh⁡(ct/1-c2)),c=0.5.
In Table5, we represent the L2 and L∞ error norms for different values of N at t=0.5. With increasing N, we obtain better results, and these error norms decrease five orders in magnitude. Figure 3 shows four-type plots at t=0.5 of Experiment 3.Table 5
Error norms for different values ofN of Experiment 3.
NL2error normL∞error norm49.13×10-37.36×10-465.29×10-43.77×10-584.13×10-52.94×10-6109.31×10-77.30×10-8122.00×10-73.82×10-8143.23×10-87.21×10-9161.87×10-83.94×10-9(a) The exact and approximate solutions forN=16 are shown. (b) The L∞ error norm for N=16 is shown. (c) The L2 error norm for different values of N is shown. (d) The L∞ error norm for different values of N is shown.
(a)(b)(c)(d)


## 5. Conclusion

In this paper, we have presented the spectral method for solving one-dimensional nonlinear sine-Gordon equation using a new orthogonal polynomial. Numerical experiments show that the spectral method is an efficient one, and the results for value ofN=16 is more accurate than the other values of N. The proposed method has higher order of accuracy even for small values of N (i.e., N=4). We obtain better results when we choose small value of Δt, in comparison to those obtained in the literature.



---

*Source: 462731-2012-02-02.xml*