Journal of Applied Mathematics

Volume 2012 (2012), Article ID 924765, 18 pages

doi:10.1155/2012/924765

Research Article

Differential Quadrature Solution of Hyperbolic Telegraph Equation

B. Pekmen1,2 and M. Tezer-Sezgin1,3

¹Institute of Applied Mathematics, Middle East Technical University, 06800 Ankara, Turkey

²Department of Mathematics, Atilim University, 06836 Ankara, Turkey

³Department of Mathematics, Middle East Technical University, 06800 Ankara, Turkey

Received 10 April 2012; Accepted 15 June 2012

Academic Editor: Roberto Barrio

Copyright © 2012 B. Pekmen and M. Tezer-Sezgin. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Differential quadrature method (DQM) is proposed for the numerical solution of one- and two-space dimensional hyperbolic telegraph equation subject to appropriate initial and boundary conditions. Both polynomial-based differential quadrature (PDQ) and Fourier-based differential quadrature (FDQ) are used in space directions while PDQ is made use of in time direction. Numerical solution is obtained by using Gauss-Chebyshev-Lobatto grid points in space intervals and equally spaced and/or GCL grid points for the time interval. DQM in time direction gives the solution directly at a required time level or steady state without the need of iteration. DQM also has the advantage of giving quite good accuracy with considerably small number of discretization points both in space and time direction.

1. Introduction

Consider the second-order linear hyperbolic telegraph equation in one-dimensional space:

and the second-order two-dimensional hyperbolic telegraph equation:

where and .

The hyperbolic differential equations are the basis for fundamental equations of atomic physics and commonly used in signal analysis for transmission and propagation of electrical signals. Equations (1.1)-(1.2) are damped wave equations for of which the solution is of great importance in wave phenomena and telegraph equation for . The initial conditions are assumed to be

with Dirichlet type boundary conditions

or Neumann type boundary conditions

for (1.1), and

with the boundary conditions

or

for (1.2).

The functions are continuous functions defined on , and similarly are defined and continuous on .

For one-dimensional telegraph equation there are several studies which give numerical solutions. Mohebbi and Dehghan [1] give a compact finite difference approximation of fourth order in space and use collocation method for the time direction. For large values of time level they are forced to increase the degree of polynomial in collocation approach to obtain accurate results. In [2], a scheme similar to finite difference method is proposed using collocation points and approximating the solution with thin plate splines radial basis functions. Dehghan and Lakestani [3] and Saadatmandi and Dehghan [4] make use of Chebyshev cardinal functions and shifted Chebyshev polynomials, respectively, for expanding the approximate solution of one-dimensional hyperbolic telegraph equation. In [4], the advantage is obtaining the closed form of the approximate solution. Dehghan and Ghesmati [5] have applied also dual reciprocity boundary element method (DRBEM) for solving second-order one-dimensional hyperbolic telegraph equation. For time discretization, Crank-Nicolson finite difference has been used, and the solution has been obtained iteratively.

Dehghan and Mohebbi [6] have used the same idea given in [1] for solving the two-dimensional linear hyperbolic telegraph equation which is the combination of finite difference in space and collocation in time directions. In [7], Mohanty and Jain introduced a new unconditionally stable alternating direction implicit (ADI) scheme of second order accurate for two-dimensional telegraph equation. The solution is progressed in time direction by splitting the systems in x- and y-directions and solving two systems for each time level. Mohanty [8] and Mohanty et al. [9] extended their studies for linear hyperbolic equations with variable coefficients in two-space dimensions. Mohanty [10] and Mohanty et al. [11] extended also ADI method solution procedure to two- and three-space dimensional telegraph equations. Ding and Zhang [12] proposed a three-level compact difference scheme of fourth order for the solution of two-dimensional, second-order, nonhomogenous linear hyperbolic equation for positive coefficients. Meshless method has also been used by Dehghan and Ghesmati [13] and Dehghan and Shokri [14] for solving two-dimensional telegraph equation. The conventional moving least squares approximation is exploited in order to interpolate the solution by using monomials from the Pascal triangle in [13], and thin plate splines radial basis functions are used for the approximation of the solution in [14]. In both of these studies, another time integration scheme has been used (finite difference) for time derivatives, and the solutions are obtained iteratively.

The usage of DQM in both space and time directions is encountered in [15]. The authors proposed to solve time-dependent problems (first order time derivatives) by a block-marching methodology in time direction. In each time block, DQM is applied both in space and time directions. The novelty of this approach is in the higher order of accuracy and less computational effort compared to 4-stage Runge-Kutta method.

In this paper, we employ the differential quadrature method in both time and space directions to obtain numerical solutions of one- and two-dimensional linear hyperbolic telegraph equations which contain second-order time derivatives. The use of DQM in time direction, which also discretizes the initial condition , automatically results in an overdetermined system. The numerical scheme then provides the solution at any time level without an iteration. This makes the main difference from the conventional time integration methods. The solution is obtained directly at all required time levels by solving one overdetermined system which contains the solution at the grid points in space and time directions. The Gauss-Chebyshev-Lobatto (GCL) points are used in space direction whereas either equally spaced or GCL grid points can be taken in time direction. The numerical procedure requires very small number of grid points in space directions and appropriate number of time grid points for reaching a certain time level.

The organization of the paper is as follows. In Section 2, we introduce the polynomial and Fourier-based DQM and discuss the formulation of the method both in time and space directions for (1.1)-(1.2). In Section 3, we discuss the accuracy and efficiency of the proposed method applying to some test problems in one- and two-dimensional cases. Section 4 gives the conclusion about the findings.

2. Method of Solution

2.1. Differential Quadrature Method

The differential quadrature method approximates the derivative of a function at a grid point by a linear summation of all the functional values in the whole computational domain.

When the function is approximated by a high-order polynomial, Shu [16] presented some explicit formulations to compute the weighting coefficients. The weighting coefficients in that weighted sum are determined using the Lagrange interpolation polynomial which has no limitation on the choice of the grid points. This leads to the polynomial-based differential quadrature method (PDQ). When Fourier series expansion is used, Fourier expansion-based differential quadrature (FDQ) is presented.

Suppose the degree of the polynomial is , then a single variable function and its first- and second-order derivatives can be approximated at a grid point by PDQ approach as

where and is the number of grid points in the whole (one-dimensional) domain, and are the Lagrange interpolated polynomials. The weighting coefficients for the first-order derivative are [16]

where

The weighting coefficients for the second-order derivative are given in terms of weighting coefficients for the first-order derivatives as

As it can be seen from these equations, there is a recurrence relationship for the second and so higher-order derivatives.

When the function on any interval is approximated by a Fourier series expansion in the form

the Lagrange interpolated trigonometric polynomials are taken as

where

Thus, we can write

In FDQ approach for nonperiodic problems, weighting coefficients for the first- and second-order derivatives are, respectively:

where .

Different grid distributions may be chosen. In this study, we choose the Gauss-Chebyshev-Lobatto (GCL) grid points in space direction which are clustered through the end points and equally spaced (ES) or GCL grid points in time direction. GCL points are the extremal values of Chebyshev polynomials and are defined as

in where is the degree of the Chebyshev polynomial.

2.2. Application of DQM To Time and Space Directions for Hyperbolic Telegraph Equation

Application to one- and two-dimensional hyperbolic telegraph equations is going to be shown using PDQ method. FDQ application will be similar, the only difference being in the weighting coefficients.

2.2.1. One-Dimensional Hyperbolic Telegraph Equation

PDQ approximations for the derivatives in (1.1) can be written as

where , and are the number of discretization points, and , and , are the weighting coefficients for the first- and second-order derivatives, in space and time directions, respectively. These coefficients are computed by (2.2a)–(2.4b) when PDQ is used. Therefore, (1.1) will be discretized as

where and .

GCL grids in -space on an interval are computed as

In time direction, equally spaced (ES) grid points on an interval as

are used.

The main system equation (2.12) may be denoted as an algebraic system

where is the unknown vector to be determined with the entries at the grid points , and the matrix consists of the known weighting coefficients. The known vector contains the function values as entries. The initial condition is inserted to the system (2.15) directly modifying the matrix and the vector . The other initial condition is expanded by DQM formulation as

Equation (2.16) will also be a system which can be described as

where the matrix and the vector contain coefficients and as entries.

If the Dirichlet type boundary conditions are given, these conditions are directly inserted to (2.15) and (2.17) with the given initial conditions . So, the coefficient matrices and will be of size and , respectively.

The systems (2.15) and (2.17) form an overdetermined system. Therefore, least square method or QR factorization will be made use of for obtaining the solution vector .

If the Neumann type boundary conditions are given, the sizes of the new coefficient matrix in (2.15) and in (2.17) will be and , respectively. Neumann boundary conditions are also discretized using PDQ as

which can be formed as system

where is a matrix of size containing the coefficients and the vector contains as entries. Then, the whole system (2.15) with , (2.17) with and (2.19) again is an overdetermined system which will be solved either by QR or least square method.

Also, the ordering of the unknown vector is important since the structure of or depends on this ordering. To get a well-conditioned matrix , the order of the unknown vector is arranged as

where if boundary conditions are Dirichlet type and if Neumann type boundary conditions are given.

2.2.2. Two-Dimensional Hyperbolic Telegraph Equation

The space and time derivatives in (1.2) can be discretized by using PDQ as

and the discretized form of (1.2) takes the form

where are the number of discretization points in directions and , respectively.

GCL grids for - and -spaces are taken as

when . Furthermore, time direction is divided equally as in (2.14) or by GCL grid points on an interval as

The initial condition and the Neumann type boundary conditions are added to the system equation (2.22) discretizing them using PDQ approximations as follows:

As in one-dimensional telegraph equation, the DQM discretized system (2.22) will be solved together with the initial and boundary conditions. If the boundary conditions are of Dirichlet type, they will be inserted to the overdetermined system (2.22) combined with the system resulting from initial conditions. For Neumann boundary conditions, the system (2.22) and (2.25) will be solved together.

One of the difficulties in two-dimensional hyperbolic telegraph equation is that the system to be solved becomes larger as are increased. This causes more memory and CPU usage. To overcome this problem, the system is reduced by removing the entries on the coefficient matrix of the system which correspond to known information (e.g., initial condition and the Dirichlet type boundary conditions). Meanwhile, the right hand side of the reduced system is also modified taking into account removed known entries in the coefficient matrix.

The order of the unknown vector in this case is organized to get a well-conditioned system. Consider the unknown vector as a matrix whose each row entry corresponds to a time level. Notationally,

where , and the matrix of size will be rewritten as a vector writing columns consecutively.

The solution will be obtained by solving only one system with the initial and boundary conditions being all inserted, and the aforementioned reduction of known entries is performed to reduce the size of the system. The solution vector contains all required time level values in it.

The solvability of the overdetermined system of equations depends on the column rank of the coefficient matrix which is in this case. When the initial and/or Neumann type boundary conditions are discretized using DQM, and added to the system, the row size is certainly greater than which makes the system overdetermined. The choice of the grid points in both space and time domains affects the stability of the system. As mentioned in the Shu’s book [16], the solution with GCL grid points becomes more stable than equally spaced grid points in both space and time directions. Moreover, appropriate choice of , and makes the final coefficient matrix full column rank.

3. Numerical Results

The discretization is performed by taking GCL grid points in x- and y-directions. When equally spaced points are taken in time direction on ,

in which is the required uptime level.

To measure the accuracy of the method, we use the following errors which can be defined as [3, 13]

where nod is the total number of all grid points in space direction, and the errors can be computed at any time level.

We present some numerical results of one- and two-dimensional hyperbolic telegraph equation for different and values to observe the accuracy and efficiency of the DQ method. Exact solutions of the test problems are available. Thus, the variations of RMS and relative errors can be obtained for discussing accuracy of DQM solutions.

3.1. Problem 1

We consider the one-dimensional hyperbolic telegraph equation (1.1) in the interval , with . The exact solution is taken as [2]

In this case, . The initial conditions are given as

We solve the equation first by taking Dirichlet type boundary conditions and then with the Neumann type boundary conditions. The Dirichlet type boundary conditions from the exact solution are

The Neumann type boundary conditions are

Table 1 shows RMS errors obtained by using Dirichlet type boundary conditions with different (equally spaced time grid is used) values. It is noticed that accuracy is achieved even with a coarse mesh () in time direction. It is improved with a finer mesh ().

Table 1: RMS errors for Problem 1 with using Dirichlet boundary conditions.

0.5

1

2

3

Table 2 reports RMS errors obtained by using Neumann boundary conditions with different number of discretization points in . As can be seen from Table 2, increasing from to improves the results in terms of accuracy, which was dropped with Neumann boundary conditions being used, with . Table 3 gives RMS errors when FDQ weighting coefficients are used in space direction and PDQ weighting coefficients in time direction for the same Neumann problem. For this particular problem, the accuracy of the results by using FDQ even with doubled from the accuracy of the results obtained by PDQ.

Table 2: RMS errors for Problem 1 with using PDQ with Neumann boundary conditions.

0.5

1

2

3

Table 3: RMS errors for Problem 1 with using FDQ with Neumann boundary conditions.

0.5

1

2

3

Figure 1 depicts the very well agreement of DQM and exact solutions at different time levels; even equally spaced time grid points are used. The number of GCL grid points in space is small (), and is considerably large compared to other time integration schemes.

Figure 1: Numerical and exact solutions at different times for Problem 1, with , using Neumann boundary conditions.

3.2. Problem 2

We consider now two-dimensional telegraph equation (1.2) in the region , with . The analytical solution given by [13] is

with the initial conditions

and the boundary conditions are Dirichlet type

We extract from the analytical solution as

Tables 4 and 5 present the comparison of the DQM solution and the exact solution in terms of relative errors for different values when equally spaced time grid points are used. In this two-dimensional problem, FDQ approximation for space derivatives, keeping PDQ approximation for time derivatives, is also studied. As can be seen from the Tables 4 and 5, both PDQ and FDQ approximations in space directions give almost the same accuracy even with coarse grids ( and ) taking a little more space grid points in FDQ. Furthermore, as is increased ( is decreased), the accuracy is increased.

Table 4: Relative errors of Problem 2 with using PDQ.

0.5

1

2

3

Table 5: Relative errors of Problem 2 using using FDQ.

0.5

1

2

3

3.3. Problem 3

In this case, (1.2) is taken with and in the region . The analytical solution is given by [6]

with the initial conditions

and the boundary conditions are of Dirichlet and Neumann type on the walls parallel to each other:

The inhomogeneity is extracted from the exact solution as

Relative errors between numerical solution (PDQ) and exact solution for different values are presented in Table 6. Even with a large , about accuracy is reached. As is decreased, accuracy is increased as expected. In Table 7, the error decreases with a higher value. Again, it is noticed that with considerably small number of grid points both in space direction () and time direction (  or  ), at least accuracy is obtained.

Table 6: Relative errors of Problem 3 with .

1

3

5

Table 7: Relative errors of Problem 3 with .

1

3

5

3.4. Problem 4

Consider the homogenous () (1.2) in the region , with . The exact solution is given by [12]

with the initial conditions

and Dirichlet and Neumann boundary conditions on parallel walls are

Table 8 compares the numerical solution (both PDQ and FDQ) with the exact solution for different time levels in terms of RMS errors. FDQ weighting coefficients are used in space directions remaining PDQ weighting coefficients in time direction (FDQ-PDQ) and PDQ approximation is used in both directions (PDQ-PDQ). In this problem, we use GCL grid points in both space and time directions. FDQ-PDQ approximation gives better accuracy than PDQ-PDQ approximation for this problem.

Table 8: RMS errors for of Problem 5 using .

0.75

2.25

3

To emphasize the importance of the GCL grid points in time direction, we present RMS errors versus the number of grid points in time direction in Figures 2 and 3. The number of GCL grid points in space direction () at is fixed, and only PDQ-PDQ approximation is considered. When the number of equally spaced points in time exceeds 28, the system becomes rank deficient. This means that large number of ES grid points in time direction causes unstable solution. On the other hand, GCL grid points in time direction still give very good accuracy for as can be seen in Figure 3. Moreover, RMS error using 11 to 23 ES grid points increases faster than RMS error using the same number of GCL grid points in time. The accuracy with grid points in time direction remains in a scale between and .

(a)

(b)

Figure 2: RMS errors at versus ES grid points in time direction for Problem 4.

Figure 3: RMS errors at versus GCL grid points in time direction for Problem 4.

This observation is important for physical problems which require the solution at a high time level. In that case, the number of grid points in time direction should be large and endure oscillations of the solution. Thus, GCL grid points in time are more preferable than ES grid points.

4. Conclusion

The polynomial-based or Fourier expansion-based DQM in space domain and polynomial-based DQM in time domain have been proposed to solve one- and two-dimensional hyperbolic telegraph equations. Differential quadrature method has the capability of producing highly accurate results using considerably small number of grid points, and thus resulting in minimal computational effort. We have used Gauss-Chebyshev-Lobatto grid in space directions which is a nonuniform grid distribution clustering through the end points. For the time direction, equally spaced points may be used for achieving the solution at a specific time level. Considering stability, GCL grids in time direction are more reliable than ES grids when the solution is required at a high time level. The accuracy and efficiency of the results have been demonstrated by root mean square (RMS) and relative errors. The use of FDQ increases the accuracy for one-dimensional telegraph equation especially when the exact solution contains trigonometric functions. For the two-dimensional hyperbolic telegraph equations, the combination of FDQ-PDQ in space and time directions is preferable.

The advantage of the use of DQM both in space and time directions lies in the fact that the solution can be obtained at a required time level by solving one system. There is no need to employ an iteration between the time levels, and large time steps can be used.

References

A. Mohebbi and M. Dehghan, “High order compact solution of the one-space-dimensional linear hyperbolic equation,” Numerical Methods for Partial Differential Equations, vol. 24, no. 5, pp. 1222–1235, 2008.

M. Dehghan and A. Shokri, “A numerical method for solving the hyperbolic telegraph equation,” Numerical Methods for Partial Differential Equations, vol. 24, no. 4, pp. 1080–1093, 2008.

M. Dehghan and M. Lakestani, “The use of Chebyshev cardinal functions for solution of the second-order one-dimensional telegraph equation,” Numerical Methods for Partial Differential Equations, vol. 25, no. 4, pp. 931–938, 2009.

A. Saadatmandi and M. Dehghan, “Numerical solution of hyperbolic telegraph equation using the Chebyshev tau method,” Numerical Methods for Partial Differential Equations, vol. 26, no. 1, pp. 239–252, 2010.

M. Dehghan and A. Ghesmati, “Solution of the second-order one-dimensional hyperbolic telegraph equation by using the dual reciprocity boundary integral equation (DRBIE) method,” Engineering Analysis with Boundary Elements, vol. 34, no. 1, pp. 51–59, 2010.

M. Dehghan and A. Mohebbi, “High order implicit collocation method for the solution of two-dimensional linear hyperbolic equation,” Numerical Methods for Partial Differential Equations, vol. 25, no. 1, pp. 232–243, 2009.

R. K. Mohanty and M. K. Jain, “An unconditionally stable alternating direction implicit scheme for the two space dimensional linear hyperbolic equation,” Numerical Methods for Partial Differential Equations, vol. 17, no. 6, pp. 684–688, 2001.

R. K. Mohanty, “An operator splitting method for an unconditionally stable difference scheme for a linear hyperbolic equation with variable coefficients in two space dimensions,” Applied Mathematics and Computation, vol. 152, no. 3, pp. 799–806, 2004.

R. K. Mohanty, M. K. Jain, and K. George, “On the use of high order difference methods for the system of one space second order nonlinear hyperbolic equations with variable coefficients,” Journal of Computational and Applied Mathematics, vol. 72, no. 2, pp. 421–431, 1996.

R. K. Mohanty, “New unconditionally stable difference schemes for the solution of multi-dimensional telegraphic equations,” International Journal of Computer Mathematics, vol. 86, no. 12, pp. 2061–2071, 2009.

R. K. Mohanty, M. K. Jain, and U. Arora, “An unconditionally stable ADI method for the linear hyperbolic equation in three space dimensions,” International Journal of Computer Mathematics, vol. 79, no. 1, pp. 133–142, 2002.

H. Ding and Y. Zhang, “A new fourth-order compact finite difference scheme for the two-dimensional second-order hyperbolic equation,” Journal of Computational and Applied Mathematics, vol. 230, no. 2, pp. 626–632, 2009.

M. Dehghan and A. Ghesmati, “Combination of meshless local weak and strong (MLWS) forms to solve the two dimensional hyperbolic telegraph equation,” Engineering Analysis with Boundary Elements, vol. 34, no. 4, pp. 324–336, 2010.

M. Dehghan and A. Shokri, “A meshless method for numerical solution of a linear hyperbolic equation with variable coefficients in two space dimensions,” Numerical Methods for Partial Differential Equations, vol. 25, no. 2, pp. 494–506, 2009.

C. Shu, Q. Yao, and K. S. Yeo, “Block-marching in time with DQ discretization: an efficient method for time-dependent problems,” Computer Methods in Applied Mechanics and Engineering, vol. 191, no. 41-42, pp. 4587–4597, 2002.

C. Shu, Differential Quadrature and Its Application in Engineering, Springer, 2000.