Research Article  Open Access
Jun Liu, Yan Wang, "A Numerical Method of High Accuracy for Linear Parabolic Partial Differential Equations", Mathematical Problems in Engineering, vol. 2012, Article ID 497365, 14 pages, 2012. https://doi.org/10.1155/2012/497365
A Numerical Method of High Accuracy for Linear Parabolic Partial Differential Equations
Abstract
We report a new numerical algorithm for solving onedimensional linear parabolic partial differential equations (PDEs). The algorithm employs optimal quadratic spline collocation (QSC) for the space discretization and twostage Gauss method for the time discretization. The new algorithm results in errors of fourth order at the gridpoints of both the space partition and the time partition, and large time steps are allowed to save computational cost. The stability of the new algorithm is analyzed for a model problem. Numerical experiments are carried out to confirm the theoretical order of accuracy and demonstrate the effectiveness of the new algorithm.
1. Introduction
Quadratic spline collocation (QSC) is a kind of numerical methods for solving systems of differential equations, which gives rise to an approximate solution in the quadratic spline space. An interesting property of the QSC method is that the optimal order of convergence can be obtained by adding appropriate perturbations to the spatial differential operators. Such an idea is used in the smooth cubic spline collocation for a special linear initial value problem [1]. Then a modified version of the cubic spline collocation is formed in [2] for more general problems. On this basis, the cubic Bspline scaling function has been studied and applied widely for a variety of problems [3–8]. The optimal QSC method is applied for the solution of boundary value problems (BVPs) in [9–11] and extended to solve elliptic partial differential equations (PDEs) in [12], as well as parabolic PDEs in [13, 14].
In [13, 14], the optimal QSC is combined with the CrankNicolson (CN) method to formulate several revised onestep QSCCN algorithms for linear parabolic PDEs. These algorithms give fourthorder accuracy at the midpoints and gridpoints of the space partition and secondorder accuracy at the gridpoints of the time partition, while they require solving a tridiagonal linear system at each time step. However, the time step of the QSCCN algorithms should be small enough to achieve the overall high performance. In addition, the oscillations are likely to be introduced by the application of the CN method. In fact, we can employ the optimal QSC method directly for the parabolic PDE and use highorder numerical methods for the resulting collocation equations, such as the twostage Gauss method [15, 16].
Based on QSC and the twostage Gauss method, we can propose a QSCTG algorithm in this paper for linear parabolic PDEs. The new algorithm gives highorder accuracy at the gridpoints of both the space partition and the time partition. Because large time steps are allowed by the twostage Gauss method, we do not need too much computations to reach a comparable order of accuracy. Moreover, the QSCTG algorithm has competitive stability properties superior to the QSCCN algorithms presented in [13, 14] and immunes to oscillations.
The remainder of this paper is organized as follows. In Section 2, we formulate the QSCTG algorithm for onedimensional linear parabolic PDEs and present the corresponding order of accuracy. In Section 3, we analyze numerically the stability properties of the new algorithm. Numerical experiments are given in Section 4 to illustrate the effectiveness of the new algorithm.
2. The QSCTG Algorithm for Linear Parabolic PDEs
We consider the numerical solution of the onedimensional linear parabolic PDE subjecting to the initial condition and the boundary conditions For simplicity, we consider in this paper the homogeneous Dirichlet boundary conditions as follows: where are given functions, is the space domain, is the time interval, and the function is to be computed. Here, we assume that .
For convenience, we denote the spatial differential operator by Then (2.1) is equivalent to the equation
2.1. QSC for Parabolic PDEs
Following [13, 14], we consider a uniform partition of the space domain with mesh size . Denote as the set of the quadratic polynomials on , and let where is the set of the functions with continuous firstorder derivatives on .
Because of the three coefficients to be determined in each quadratic polynomial and the continuity of the firstorder derivative at the gridpoints, the dimension of the space is . Furthermore, let be a set of piecewise polynomial basis functions of the space .
The approximate solution in space to system (2.1), (2.2), and (2.4) can be written as where , are degrees of freedom. For any fixed value of , relations are needed to specify the approximate solution . Obviously, two relations of them can be obtained from the boundary conditions. Therefore, we have to choose other points on , from which the relations can be found. The points together with the two boundary points are called collocation points, which can be denoted by With the collocation points, we can obtain the following relations: where is the value at the point on the interpolation of in . For some specially chosen basis functions for , the relations (2.11) have simple forms.
Let We choose a set of quadratic Bsplines as the basis functions of the space . Following [13], we denote by the space of quadratic splines satisfying homogeneous Dirichlet boundary conditions. The dimension of the space is , and a set of its basis functions is Then the quadratic spline approximate solution of system (2.1), (2.2), and (2.4) can be written as The values of the quadratic spline basis functions and their derivatives at the collocation points are respectively.
Denote as a diagonal matrix with the diagonal elements listed in the brackets. Let The relations (2.11) lead to the following system of ODEs: where , The vector satisfies , where is the interpolation of at collocation points.
The spline collocation approximate solution from system (2.18) has secondorder accuracy. Similar to the way to get optimal spline collocation approximation for BVPs in [9, 13], the optimalorder approximation to the system (2.1), (2.2), and (2.4) can be obtained by the following system of ODEs with extra perturbations: where
By the system (2.20) from the optimal QSC, the discretization error at the midpoints and gridpoints of the uniform space partition can be given.
2.2. The TwoStage Gauss Method for the Collocation Equation
Denote the matrices and system (2.20) can be rewritten as We employ a kind of twostage implicit RungeKutta methods for system (2.23), with the following scheme: where , and is an approximation to . The RungeKutta method (2.24), which is based on GaussLegendre quadrature, is also called the twostage Gauss method, with the fourth order of accuracy.
Based on the results from scheme (2.24), we can obtain the approximate solution of system (2.1), (2.2), and (2.4) by . The hybrid algorithm by QSC and the twostage Gauss method (2.24) is called QSCTG algorithm in this paper, which gives rise to discretization errors of .
3. Stability of QSCTG
In this section, we consider the stability of QSCTG for a model problem as follows: where is a positive constant, is a given function, is the space domain, is the time interval, and is the unknown function.
For system (3.1), the optimal QSC gives rise to the following collocation equation: where is the interpolation of at the collocation points. The twostage Gauss method for system (3.2) is where , , and are defined in (2.20). Denote and . Substitute and of (3.3) into the first equation of (3.3), then we have where the matrix is of size , which can be regarded as the iteration matrix for the scheme described by (3.4). We denote by the th power of and follow the stability analysis presented in [14]. The stability of (3.4) is guaranteed if is bounded independently of .
For (3.4), we consider the quantities , with , , and , for several values of and . Figure 1 shows how behaves as increases, with , for several values of . It can be observed that the quantities are bounded by a constant which is independent of , and thus is bounded independently of . Therefore, the scheme of QSCTG for the model problem (3.1) is stable without any restriction on the time step size.
(a)
(b)
(c)
(d)
To illustrate the advantages on the stability of QSCTG, we recall the QSCCN0 algorithm presented in [13, 14] for the model problem (3.1). Applying the QSCCN0 algorithm to the model problem (3.1) gives rise to the linear equations of the form with the matrices where Let , which is the iteration matrix of the scheme described by (3.5). It has been analyzed in [13, 14] that the QSCCN0 algorithm for the model problem (3.1) is stable without any restriction on the time step size. The results by the comparison between the behaviors of and , with the same values of and , are also shown in Figure 1. It can be seen clearly that, for fixed value of , the upper bound of the quantities is smaller than the upper bound of the quantities .
In Figure 2, we compare the spectral radii of the iteration matrix with the spectral radii of the iteration matrix , versus , for different values of . We can observe that the spectral radii of the matrix are smaller than those of for all values of . It seems that the QSCTG algorithm is better than QSCCN0 as far as stability is concerned.
(a)
(b)
(c)
(d)
4. Numerical Experiments
We compute a linear parabolic PDE as follows: where is a constant. The exact solution of system (4.1) is
To perform the QSCTG algorithms, we employ a uniform partition for the space domain with mesh size and choose the collocation points The resulting collocation equation can be written as where the vector and the matrices , , and , of size , are the same as the matrices in system (2.20).
By the QSCTG algorithm, we can obtain an approximate solution to system (4.1). The resulting error is measured by where denotes the th entry of , and is the true solution at .
Case 1 (when ). We first investigate the contributions of the optimal quadratic spline collocation and the twostage Gauss method, respectively. If we choose the time step , which is small enough to ensure the errors are introduced by QSC, the observed errors with different values of and the corresponding orders of accuracy are shown in Table 1. Similarly, if we choose , which is big enough to ensure the errors are introduced by the time integration, the observed errors with different values of and the corresponding orders of accuracy are shown in Table 2. The results in Tables 1 and 2 confirm the discretization errors of for the QSCTG algorithm.
To show the advantages of the QSCTG algorithm, we employ the QSCCN0 algorithm presented in [13, 14] for comparison, which has been proved to be much efficient and unconditionally stable. In Table 3, we present the errors and the time cost (measured in seconds) by the QSCCN0 algorithm and the QSCTG algorithm for system (4.1). We can see that the QSCTG algorithm needs much less running time when achieving a desired high accuracy.
Furthermore, we notice that the QSCCN0 algorithm with behaves worse than that with . In fact, the approximate solutions by QSCCN0 contain spurious oscillations if the value of is large [17]. The QSCTG algorithm has no such restriction and immunes to oscillations.



Case 2 (when ). To show the advantages of QSCTG more clearly, we implement the QSCTG and QSCCN0 again for system (4.1) with . The observed errors and the running time are compared in Table 4. We can see that the QSCTG algorithm costs much less running time than QSCCN0, while it reaches much better accuracy.

5. Conclusions
We have proposed a QSCTG algorithm for solving linear onedimensional parabolic PDEs. The space discretization is dealt with by the optimal quadratic spline collocation, and the time discretization is treated by the twostage Gauss method. High order of accuracy both in space and time discretizations can be achieved. The QSCTG algorithm has been confirmed numerically to be unconditional stable. The results of the numerical experiments show that the QSCTG algorithm costs much less running time than the very efficient QSCCN0 algorithm, which is presented in [13, 14], when solving the same system and achieving the same high accuracy.
Acknowledgments
This work was supported by the Natural Science Foundation of China (NSFC) under Grant 11071192 and the International Science and Technology Cooperation Program of China under Grant 2010DFA14700. This work was partly supported by NSFC under Grant 11101127 (Y. Wang).
References
 J. C. Cavendish and C. A. Hall, “${L}_{\infty}$convergence of collocation and Galerkin approximations to linear twopoint parabolic problems,” Aequationes Mathematicae, vol. 11, pp. 230–249, 1974. View at: Google Scholar
 D. Archer, “An $O({h}^{4})$ cubic spline collocation method for quasilinear parabolic equations,” SIAM Journal on Numerical Analysis, vol. 14, no. 4, pp. 620–637, 1977. View at: Google Scholar
 M. Lakestani and M. Dehghan, “Numerical solution of Riccati equation using the cubic Bspline scaling functions and Chebyshev cardinal functions,” Computer Physics Communications, vol. 181, no. 5, pp. 957–966, 2010. View at: Publisher Site  Google Scholar
 M. Dehghan and M. Lakestani, “Numerical solution of nonlinear system of secondorder boundary value problems using cubic Bspline scaling functions,” International Journal of Computer Mathematics, vol. 85, no. 9, pp. 1455–1461, 2008. View at: Publisher Site  Google Scholar
 M. Dehghan and M. Lakestani, “The use of cubic Bspline scaling functions for solving the onedimensional hyperbolic equation with a nonlocal conservation condition,” Numerical Methods for Partial Differential Equations, vol. 23, no. 6, pp. 1277–1289, 2007. View at: Publisher Site  Google Scholar
 M. Lakestani and M. Dehghan, “Numerical solution of FokkerPlanck equation using the cubic Bspline scaling functions,” Numerical Methods for Partial Differential Equations, vol. 25, no. 2, pp. 418–429, 2009. View at: Publisher Site  Google Scholar
 M. Lakestani and M. Dehghan, “Numerical solutions of the generalized KuramotoSivashinsky equation using Bspline functions,” Applied Mathematical Modelling, vol. 36, no. 2, pp. 605–617, 2012. View at: Publisher Site  Google Scholar
 M. Lakestani, M. Dehghan, and S. Irandoustpakchin, “The construction of operational matrix of fractional derivatives using Bspline functions,” Communications in Nonlinear Science and Numerical Simulation, vol. 17, no. 3, pp. 1149–1162, 2012. View at: Publisher Site  Google Scholar
 E. N. Houstis, C. C. Christara, and J. R. Rice, “Quadraticspline collocation methods for twopoint boundary value problems,” International Journal for Numerical Methods in Engineering, vol. 26, no. 4, pp. 935–952, 1988. View at: Publisher Site  Google Scholar
 B. Bialecki, G. Fairweather, A. Karageorghis, and Q. N. Nguyen, “Optimal superconvergent one step quadratic spline collocation methods,” BIT. Numerical Mathematics, vol. 48, no. 3, pp. 449–472, 2008. View at: Publisher Site  Google Scholar
 G. Fairweather, A. Karageorghis, and J. Maack, “Compact optimal quadratic spline collocation methods for the Helmholtz equation,” Journal of Computational Physics, vol. 230, no. 8, pp. 2880–2895, 2011. View at: Publisher Site  Google Scholar
 C. C. Christara, “Quadratic spline collocation methods for elliptic partial differential equations,” BIT. Numerical Mathematics, vol. 34, no. 1, pp. 33–61, 1994. View at: Publisher Site  Google Scholar
 T. Chen, An efficient algorithm based on quadratic spline collocation and finite difference methods for parabolic partial differential equations [M.S. thesis], University of Toronto, Ontario, Canada, 2005.
 C. C. Christara, T. Chen, and D. M. Dang, “Quadratic spline collocation for onedimensional linear parabolic partial differential equations,” Numerical Algorithms, vol. 53, no. 4, pp. 511–553, 2010. View at: Publisher Site  Google Scholar
 A. M. Stuart and A. R. Humphries, Dynamical Systems and Numerical Analysis, vol. 2 of Cambridge Monographs on Applied and Computational Mathematics, Cambridge University Press, New York, NY, USA, 1996.
 J. C. Butcher, Numerical Methods for Ordinary Differential Equations, John Wiley& Sons, Chichester, UK, 2nd edition, 2008. View at: Publisher Site
 J. W. Thomas, Numerical Partial Differential Equations: Finite Difference Methods, vol. 22 of Texts in Applied Mathematics, Springer, New York, NY, USA, 1995.
Copyright
Copyright © 2012 Jun Liu and Yan Wang. 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.