Abstract and Applied Analysis

Volume 2013 (2013), Article ID 686483, 10 pages

http://dx.doi.org/10.1155/2013/686483

Research Article

Numerical Solution and Simulation of Second-Order Parabolic PDEs with Sinc-Galerkin Method Using Maple

Aydin Secer

Department of Mathematical Engineering, Yildiz Technical University, Davutpasa, 34210 Istanbul, Turkey

Received 4 March 2013; Accepted 12 June 2013

Academic Editor: Adem Kılıçman

Copyright © 2013 Aydin Secer. 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.

An efficient solution algorithm for sinc-Galerkin method has been presented for obtaining numerical solution of PDEs with Dirichlet-type boundary conditions by using Maple Computer Algebra System. The method is based on Whittaker cardinal function and uses approximating basis functions and their appropriate derivatives. In this work, PDEs have been converted to algebraic equation systems with new accurate explicit approximations of inner products without the need to calculate any numeric integrals. The solution of this system of algebraic equations has been reduced to the solution of a matrix equation system via Maple. The accuracy of the solutions has been compared with the exact solutions of the test problem. Computational results indicate that the technique presented in this study is valid for linear partial differential equations with various types of boundary conditions.

1. Introduction

Sinc methods for differential equations were originally introduced by Stenger in [1–3]. The sinc functions were first analyzed in [4, 5] and a detailed research of the method for two-point boundary-value problems can be found in [6, 7]. In [8], parabolic and hyperbolic problems are presented in detail. To solve a problem arising from chemical reactor theory, the properties of the sinc-Galerkin method are used to reduce the computation of nonlinear two-point boundary-value problems to some algebraic equations in [9]. A computer algorithm for sinc method to solve numerically the linear and nonlinear ODEs and their simulations has been presented in [7, 10], respectively. The full sinc-Galerkin method is developed for a family of complex-valued partial differential equations with time-dependent boundary conditions [9]. A study of the performance of the Galerkin method using sinc basis functions for solving Bratu’s problem is presented in [11]. In [12] a numerical algorithm has been presented for recovering the unknown function and obtaining a solution to the inverse ill-posed problem. They have presented a Galerkin method with the sinc basis functions in both space and time domains for solving the direct problem. A sinc-collocation method has been developed for solving linear systems of integrodifferential equations of Fredholm and Volterra type with homogeneous boundary conditions in [13].

2. Sinc-Approximation Formula for PDEs

We use the sinc-Galerkin method as mentioned in [1] to derive an approximate solution of the following:

For the equation given above, the sinc-Galerkin scheme can be developed in both space and time directions as follows.

In general, approximations can be constructed for infinite, semi-infinite, and infinite intervals and both spatial and time spaces will be introduced. Define the function

which is a conformal mapping from , the eye-shaped domain in the -plane, onto the infinite strip, , where

A more general form of sinc basis according to intervals can be given as follows:

where

and the conformal maps for both directions

are used to define the basis functions on the intervals and , respectively. represents the mesh sizes in the space direction and the time direction, respectively. The sinc nodes and are chosen so that , .

Here the function is an inverse mapping of . We may define the range of on the real line as

For the evenly spaced nodes on the real line, the image which corresponds to these nodes is denoted by

where , for all .

The sinc basis functions in (4) do not have a derivative when tends to 0 or 1. We modify the sinc basis functions as

where

For the temporal space, we construct an approximation by defining the function

which is a conformal mapping from , the wedge-shaped domain in the -plane, onto the infinite strip, , where

derived from composite translated functions

for .

Here and . We may define on the real line as

For the evenly spaced nodes on the real line, the image which corresponds to these nodes is denoted by

where , for all .

A list of conformal mappings may be found in Table 1 [14].

Table 1: Conformal mappings and nodes for several subintervals of .

ϕ 

Definition 1. Let be the class of functions that are analytic in and satisfy

where

and on the boundary of satisfy

The proof of the following theorems can be found in [1].

Theorem 2. Let be , , and then for sufficiently small

where

For the sinc-Galerkin method, the infinite quadrature rule must be truncated to a finite sum; the following theorem indicates the conditions under which exponential convergence results.

Theorem 3. If there exist positive constants , and such that

then the error bound for the quadrature rule (19) is

The infinite sum in (19) is truncated with the use of (20) to arrive at (22).

Making the selections

where is integer part of statement,

Theorems 2 and 3 can be used to approximate the integrals that arise in the formulation of the discrete systems corresponding to two-point BVPs.

3. Discrete Solutions Scheme for Two-Point BVPs

In ordinary differential equations

on , sinc solution is assumed as an approximate solution in the form of series with terms

The coefficients are determined by orthogonalizing the residual with respect to the sinc basis functions where . An inner product for two continuous functions such as and can be given by the following formula

where is the weight function and is chosen depending on boundary conditions. If we implement the above inner product rule in orthogonalization, this yields the discrete sinc-Galerkin system:

Now, we are going to derive discrete sinc-Galerkin system for PDEs. Assume is the approximate solution of (1).Then, the discrete system takes the following form:

The coefficients are determined by orthogonalizing the residual with respect to the sinc basis functions where for . In this case the inner product takes the following form:

The choice of the weight function in the double integrand depends on the boundary conditions, the domain, and the partial differential equation. Therefore, the discrete Galerkin system is

4. Matrix Representation of the Derivatives of Sinc Basis Functions at Nodal Points

The sinc-Galerkin method actually requires the evaluated derivatives of sinc basis functions at the sinc nodes, . The th derivative of with respect to , evaluated at the nodal point , is denoted by

The expressions in (14) for each and can be stored in a matrix . For

where

The chain rule has been used for the -derivative of product sinc functions. For example, when ,

Now, we are going to develop discrete form for (1). We choose for special case the parameters as follows for the spatial dimension:

and for the temporal space as

The discrete form of (1) can be given the following form:

We solve this by taking our approximating basis functions to be

If we apply sinc-quadrature rules with the help of (32)–(37) on the definite integral given (38) by using (39), we can get the following matrix system.

Let denote a diagonal matrix, whose diagonal elements are and nondiagonal elements are zero. Then (38) reproduces the following matrixes accordingly.

Firstly we set the coefficient matrix as follows:

Finally, for the right side function given (1) can be written in the following matrix form:

Using (32)–(37) we arrive at a matrix system given in [1] as follows:

Finally, by using Maple Computer Algebra Software, the matrix system (42) can be solved by using LU or QR decomposition method and unknown coefficients can be found. After calculation of we get approximate solution as follows:

5. Numerical Simulation

The example in this section will illustrate the sinc method.

Example 4. This problem has been addressed in [1]. The following equation is given in Dirichlet-type boundary condition:

The particular solution of (44) can be calculated via separation variables rules and can be given as follows:

For (44) we choose sinc components here in the following:

According to the above parameters, the approximate solution simulation of (44) has been given in Figure 1 and numerical results also can be found in Table 2.

Table 2: Numerical results.

Exact solution
Sinc-Galerkin solution
Error

,
0.1
0.35446621870000
0.59102026517764600
0.23655404647764600

0.11
0.00001833412226
−0.08327827240548060
0.08329660652774060

0.21

−0.03068110846399060
0.03068110941228980

0.31

−0.01083286679129120
0.01083286679134030

0.41

−0.000420168677914081
0.000420168677914083

0.51

−0.00108069489771883
0.00108069489771883

0.61

−0.00492511086695295
0.00492511086695295

0.71

0.00331706758955819
0.00331706758955819

0.81

−0.00303031084399243
0.00303031084399243

0.91

0.00342519240062265
0.00342519240062265

(a) The sinc-Galerkin solutions according to the grid points size N

(b) Intersection of surfaces for . The left-side figure is normal perspective and the second one has medium perspective view

(c) The exact solution (or particular solution)

Figure 1: Simulation of approximate solution.

6. Conclusions

We have developed a Maple algorithm to solve and simulate second-order parabolic PDEs with Dirichlet-type boundary conditions based on sinc-Galerkin approximation on some closed real intervals and the method has been compared with the exact solutions. When compared with other computational approaches, this method turns out to be more efficient in the sense that selection parameters and changing boundary conditions and also giving different problems to the algorithms. The accuracy of the solutions improves by increasing the number of sinc grid points . The method presented here is simple and uses sinc-Galerkin method that gives a numerical solution, which is valid for various boundary conditions. Several PDEs have been solved by using our technique in less than 20 seconds. All computations and graphical representations have been prepared automatically by our algorithm.

Appendix

See Algorithm 1.

restart:

with(linalg):

with(LinearAlgebra):

N:=16:

d:=Pi/2:

h:=0.75/sqrt(N);

s:=.5/sqrt(N);

SIZE:=2*N+1;

delta0:=i,j->piecewise(j=i,1,j<>i,0):

delta1:=i,j->piecewise(i=j,0,i<>j,((-1)(i-j))/(i-j)):

delta2:=i,j->piecewise(i=j,(-Pi2)/3,i<>j,-2*(-1)(i-j)/(i-j)2):

I_0:=Matrix(SIZE,delta0):

I_1:=Matrix(SIZE,delta1):

I_2:=Matrix(SIZE,delta2):

xk:=1/2+1/2*tanh(k*h/2);

xk_func:=unapply(xk,k);

phi:=unapply(log((x)/(1-x)),x);

Dphi:=unapply(simplify(diff(phi(x),x)),x);

g:=unapply(simplify(1/diff(phi(x),x)),x);

Dg:=unapply(diff(g(x),x),x);

gk:=unapply(subs(x=xk,g(x)),k);

Dgk:=unapply(subs(x=xk,Dg(x)),k);

g_Div_Dphi:=unapply(g(x)/Dphi(x),x);

g_Div_Dphi_k:=unapply(subs(x=xk,g_Div_Dphi(x)),k);

#Temporal  Spaces;

tl:=exp(l*s);

tl_func:=unapply(tl,l);

gamm:=unapply(log(t),t);

Dgam:=unapply(diff(gamm(t),t),t);

g_gamm:=unapply(1/diff(gamm(t),t),t);

g_gamm_l:=unapply(subs(t=tl,g_gamm(t)),l);

gamm_Div_Dgam:=unapply(g_gamm(t)/Dgam(t),t);

gamm_Div_Dgam:=  t->  t2;

gamm_Div_Dgam_l:=unapply(subs(t=tl,gamm_Div_Dgam(t)),l);

GenerateDiagonalAm  :=  proc(  x  )

                                                local  i:=1:

                                                  local  A:=Matrix(SIZE):

                                                      for  i  from  1  by  1  to  SIZE

                                                       do

                                                          Ai,i:=evalf(x(-N+i-1)):

                                                       end  do:

                                                      return  A;

                                        end  proc;

B:=  -2*h*MatrixMatrixMultiply(I_0,GenerateDiagonalAm(gk))+

      MatrixMatrixMultiply(I_1,GenerateDiagonalAm(Dgk))+1/h*I_2:

DD:=h*GenerateDiagonalAm(g_Div_Dphi_k):

C:=MatrixMatrixMultiply(GenerateDiagonalAm(g_gamm_l),s*(I_0-I_1)):

E:=s*GenerateDiagonalAm(gamm_Div_Dgam_l):

Fkl:=unapply((Pi2-4)*sin(Pi*xk_func(k-N-1))*exp(-tl_func(l-N-1)),k,l);

F:=evalf(Matrix(SIZE,Fkl)):

V:=Matrix(SIZE,v):

EQN_SYS:=evalf(

    MatrixMatrixMultiply(

                          MatrixMatrixMultiply(

                                  Matrix(inverse(DD)),B

                                                                ),V

                            )

                    )

        +evalf(

                MatrixMatrixMultiply(

                        MatrixMatrixMultiply(

                                  V,C),Matrix(inverse(E)

                                                              )

                                                      )

                    ):

SYS:=:

            for  i  from  1  by  1  to  SIZE

                  do

                        for  j  from  1  by  1  to  SIZE

                                do

                                      SYS:=op(SYS),EQN_SYS(i,j)=F(i,j);

                                end  do:

                  end  do:

vars:=seq(seq(v(i,j),i=1..2*N+1),j=1..2*N+1):

A,b:LinearAlgebraGenerateMatrix(evalf(SYS),vars):

c:=linsolve(A,b):

CoeffMatrix=Matrix(SIZE):

cnt:=1;

for  i  to  SIZE  do

        for  j  to  SIZE  do

            CoeffMatrixj,i:=ccnt:

            cnt:=cnt+1

        end  do;

end  do;

CoeffMatrix:=Matrix(CoeffMatrix,SIZE):

ApproximateSol:=unapply(

                                          evalf(

                                                        sum(

                                                            sum("CoeffMatrix"[m+N+1,n+N+1]

                                  *sin(Pi*(phi(x)-m*h)/h)/(Pi*(phi(x)-m*h)/h)

                              *sin(Pi*(gamm(t)-n*s)/s)/(Pi*(gamm(t)-n*s)/s)

                                              ,m=-N..N),

                                          n=-N...N)

                                                    )

                                                      +exp(-4*t)*sin(Pi*x)

                                            ,x,t):

plot3d(ApproximateSol(x,t),x=0..1,t=0..1):

Exact:=unapply(exp(-Pi2*t)*sin(Pi*x),x,t);

plot3d(Exact(x,t),x=0..1,t=0..1);

plot3d({Exact(x,t),ApproximateSol(x,t)},x=0..1,t=0..1);

XX:=.6;

#Numerical  Comparision  EXACT

for  j  from  0.1  to  10  by  1

      do

        evalf(Exact(XX,j)):

      od;

#Numerical  Comparision  APPROX

for  j  from  0.1  to  10  by  1

      do

        evalf(ApproximateSol(XX,j)):

      od;

#Numerical  Comparision  ERROR

for  j  from  0.1  to  10  by  1

      do

        abs(evalf(evalf(ApproximateSol(XX,j)-Exact(XX,j)))):

      od;

Algorithm 1

Conflict of Interests

The author declares that he has no conflict of interests.

References

F. Stenger, “A Sinc-Galerkin method of solution of boundary value problems,” Mathematics of Computation, vol. 33, no. 145, pp. 85–109, 1979.

F. Stenger, “Summary of Sinc numerical methods,” Journal of Computational and Applied Mathematics, vol. 121, no. 1-2, pp. 379–420, 2000.

F. Stenger, “Approximations via Whittaker's cardinal function,” Journal of Approximation Theory, vol. 17, no. 3, pp. 222–240, 1976.

E. T. Whittaker, “On the functions which are represented by the expansions of the interpolation theory,” Proceedings of the Royal Society of Edinburg, vol. 35, pp. 181–194, 1915.

J. M. Whittaker, Interpolation Function Theory, Cambridge Tracts in Mathematics and Mathematical Physics, no. 33, Cambridge University Press, London, UK, 1935, by E. F. Beckenbach, McGraw-Hill, New York, NY, USA, 1961.

J. Lund, “Symmetrization of the sinc-Galerkin method for boundary value problems,” Mathematics of Computation, vol. 47, no. 176, pp. 571–588, 1986.

A. Secer and M. Kurulay, “The Sinc-Galerkin method and its applications on singular Dirichlet-type boundary value problems,” Boundary Value Problems, vol. 2012, article 126, 2012.

K. M. McArthur, K. L. Bowers, and J. Lund, “Numerical implementation of the Sinc-Galerkin method for second-order hyperbolic equations,” Numerical Methods for Partial Differential Equations, vol. 3, no. 3, pp. 169–185, 1987.

A. Saadatmandi, M. Razzaghi, and M. Dehghan, “Sinc-Galerkin solution for nonlinear two-point boundary value problems with applications to chemical reactor theory,” Mathematical and Computer Modelling, vol. 42, no. 11-12, pp. 1237–1244, 2005.

A. Secer, M. Kurulay, M. Bayram, and M. A. Akinlar, “An efficient computer application of the Sinc-Galerkin approximation for nonlinear boundary value problems,” Boundary Value Problems, vol. 2012, article 117, 2012.

J. Rashidinia, K. Maleknejad, and N. Taheri, “Sinc-Galerkin method for numerical solution of the Bratu's problems,” Numerical Algorithms, vol. 62, no. 1, pp. 1–11, 2013.

A. Shidfar and A. Babaei, “The sinc-Galerkin method for solving an inverse parabolic problem with unknown source term,” Numerical Methods for Partial Differential Equations, vol. 29, no. 1, pp. 64–78, 2013.

M. Zarebnia and M. G. A. Abadi, “A numerical sinc method for systems of integro-differential equations,” Physica Scripta, vol. 82, no. 5, Article ID 055011, 2010.

J. Lund and K. L. Bowers, Sinc Methods for Quadrature and Differential Equations, SIAM, Philadelphia, Pa, USA, 1992.