A boundary value approach for solving three-dimensional elliptic and hyperbolic partial differential equations

In this article, the boundary value method is applied to solve three dimensional elliptic and hyperbolic partial differential equations. The partial derivatives with respect to two of the spatial variables (y, z) are discretized using finite difference approximations to obtain a large system of ordinary differential equations (ODEs) in the third spatial variable (x). Using interpolation and collocation techniques, a continuous scheme is developed and used to obtain discrete methods which are applied via the Block unification approach to obtain approximations to the resulting large system of ODEs. Several test problems are investigated to elucidate the solution process.


Background
This paper is devoted to the numerical computation of the three dimensional elliptic and hyperbolic PDEs. The general three dimensional problem (without cross derivatives) is given as defined in the domian � = {(x, y, z)/(x, y, z) ∈ (L 1 , L 2 , L 3 ) × (L 4 , L 5 , L 6 )} with boundary ∂�. Associated with (1) are the initial conditions or the boundary conditions and the Dirichlet boundary conditions (1) a ∂ 2 u ∂x 2 + b ∂ 2 u ∂y 2 + c ∂ 2 u ∂z 2 + d ∂u ∂x + e ∂u ∂y + f ∂u ∂z = G(x, y, z, u), (x, y, z) ∈ � ⊂ R 3 u(L 1 , y, z) = g 1 (y, z), u x (L 1 , y, z) = g 2 (y, z) ( y, z) ∈ (L 2 , L 3 ) × (L 5 , L 6 ) u(L 1 , y, z) = g 1 (y, z), u(L 4 , y, z) = g 3 (y, z) ( y, z) ∈ (L 2 , L 3 ) × (L 5 , L 6 ) u(x, L 2 , z) = h 1 (x, z), u(x, L 5 , z) = h 2 (x, z) u(x, y, L 3 ) = h 3 (x, y), u(x, y, L 6 ) = h 4 (x, y) The unknown function u, the variable coefficients a, b, c, d, e, f and the forcing term G are assumed to be continuously differentiable and have the required partial derivatives on . The numerical approximation of (1) has received great attention from researchers in the past few decades. This is because most modeled physical and engineering processes results in PDEs of the form (1). The standard convection diffusion equation is obtained if a = b = c � = 0 and G ≡ G(x, y, z) for all x, y, z ∈ (see Ge and Zhang 2002). Many transport processes including fluid flows and heat transfer can be modeled by a convective diffusive equation which describes the convection and diffusion of various physical quantities such as heat and momentum (see Roache 1976). Several numerical schemes have been developed for (1) which include the compact difference scheme (Ge and Zhang 2002;Spotz 1995;Spotz and Carey 1996;Zhang 1998a, b, c, d;Zhang et al. 2000), compact alternating direct implicit scheme (Cui 2010;Liao and Sun 2010), finite difference methods (Jain et al. 1992;Mohanty and Singh 2006), method of lines (Brastos 2007;Dehghan and Shokri 2008) , Jacobi elliptic function method (Bhrawy and Abdelkawy 2013), the use of collocation and radial basis functions (Shakeri and Dehghan 2008), Riccati-Bernoulli sub ODE method (Yang et al. 2015), the use of the expansion methods (He 1998;Alam et al. 2014;Roshid et al. 2013Roshid et al. , 2014a among others. The method of lines approach is commonly used for solving PDEs whereby the PDEs is converted into a system of ODEs by replacing the appropriate derivatives by finite difference approximations. Our objective is to convert the 3D PDEs into a system of ODEs by replacing two of the spatial derivatives by central difference approximations. The resulting system is then solved using the BVM. Specifically, we discretize the y and z variable with mesh spacings We then define the vectors and where u m,n (x) ≈ u(x, y m , z n ) and G m,n (x) =≈ G(x, y m , z n , u m,n ).
Also, we replace the partial derivatives ∂ 2 u(x, y, z) ∂y 2 , ∂ 2 u(x, y, z) ∂z 2 , ∂u(x, y, z) ∂y and ∂u(x, y, z) ∂z occuring in (1) by central difference approximations Problem (1) leads to the semidiscretized system which can be written in the form subject to the initial conditions or the boundary conditions ) arising from the semidiscretized system (2) which is expressed in the form (3) and whose solution is sought using the BVMs. We note that G is a vector of constants. The BVMs are a class of linear multistep methods (LMM) with step number k and whose k additional conditions are not only imposed at the beginning of the integration process but also at the end so that they form a discrete analog of the continuous boundary value problems. Thus, they are used for the numerical approximation of both initial and boundary value problems. They have been used for the solution of first order 1D initial and boundary value problems and their convergence and linear stability properties have been fully discussed in Brugnano and Trigiante (1998). Recently, Biala and Jator (2015) developed BVMs for the direct solution of systems of the general second order ODEs. One main feature of the BVMs is that they can be used in the same way for solving both initial and boundary value problems. Therefore, such methods are the best candidate for solving the semididscretized PDEs in (3). The BVMs simultaneously generates approximate solutions (u 1,m,n , u 2,m,n , . . . , u l,m,n ) T to the exact solution (u m,n (x 1 ), u m,n (x 2 ), . . . , u m,n (x l )) T on the entire interval of integration. This approach has the advantage of producing smaller global errors than those produced by the step-by-step methods due to the presence of accumulated errors at each step in the step-by-step method. The paper is organized as follows: in "Description of the BVM", we derive a continuous scheme which is used to formulate the BVM as well as investigate the properties of the BVM. The computational complexities associated with the method is addressed in "Computational procedures". Several numerical test examples are given in "Numerical examples" to show the accuracy of the method. We give some concluding remarks in "Conclusion".

Development of the continuous BVM
In this section, we discuss the construction of the continuous scheme, via the interpolation and collocation approach (Jator and Li 2012), from which (6) is derived.

are basis functions and W is a vector given by
Proof We define the polynomial basis functions where α i+1,j and h 2 β i+1,j are coefficients to be uniquely determined. Substituting (10) into (7), we have which may be written as and expressed in the form where Imposing conditions (8) on (11), we obtain a system of seven equations, which can be expressed as V = LW where L = (ℓ 0 , ℓ 1 , . . . , ℓ 6 ) T is a vector of seven undetermined coefficients. Using the Crammer's rule, the elements of L can be obtained and are given as where V j is obtained by replacing the jth column of V by W. We rewrite (11) as (9) using the newly found elements of L.

Convergence analysis
We discuss the convergence of the BVMs in the following theorem Theorem 2 Let U be an approximation of the solution vector U for the system obtained on a partition π L : We also define the vectors where L(h) is the local truncation error vector of the formulas in (12). The exact form of the system formed by (12) and for sufficiently small h, A + M is a monotone matrix and thus nonsingular. Therefore and which shows that the methods are convergent and the global error is of order O(h 6 ).

Computational procedures
The method (12) can also be expressed in block form as where the positive integer Ŵ = N /4 is the number of blocks, and A 0 , A 1 , B 0 , and B 1 are matrices each of dimension 8 whose entries are given by the coefficients of (12). Equation (1) is converted to (3) by discretizing the partition π M and π N , given by where are constant stepsizes of the partition π M and π N respectively, m = 1(1)M, n = 1(1)N , M and N are positive integers and m, n are the grid index in the y and z direction respectively.
The resulting system of ODEs (3) is then solved on the partition π L . The block unification of (17) lead to a large system of finite difference equations which is then solved to provide all solutions of (3) on .
Step 3 The system obtained in Step 2 is solved using the feature NSolve in Mathematica 10.0 for linear problems and FindRoot (which incorporates the Newton's method) for nonlinear problems.
Step 4 The solution of (1) is approximated by the solution in Step 3 as

Numerical examples
In this section, some examples are investigated to show the reliability and efficiency of the proposed scheme in this paper.

Test 1
We consider the Laplace equation with non zero forcing term G(x, y, z) and with zero boundary values on the entire ∂� given in Zhang (1998c) and whose solution is given as where = [0, 1] 2 × [0, L 4 ]. Table 1 shows the errors in the l ∞ norm with different mesh sizes. Figure 1 shows the plot of the exact, approximate and error function when x = 0.5.

Test 4
We also consider the singularly perturbed convection diffusion equation given in Mohanty and Singh (2006) ǫ ∂ 2 u ∂x 2 + ∂ 2 u ∂y 2 + ∂ 2 u ∂z 2 = ∂u ∂x  where σ 2 = 2π 2 + 1 4ǫ 2 . We have solved the problem using ǫ = 0.1 to 1.0. The l ∞ norms are given in Table 4 with different meshsizes. Figure 4 also shows the plot of the exact, approximate and error function when x = 0.5.

Test 5
The singularly perturbed elliptic boundary value problem given in Mohanty and Singh (2006) defined in the domain = [0, 1] 3 with boundary ∂� and subject to the Dirichlet boundary conditions and where the forcing term G is set to satisfy the exact solution given as u(x, y, z) = x 2 cosh(y) sinh(z). We have solved the problem using ǫ = 0.001. The l ∞ norms are given in Table 5 with different meshsizes and for different values of α. Figure 5 also shows the plot of the exact, approximate and error function when x = 0.5. ǫ ∂ 2 u ∂x 2 + ∂ 2 u ∂y 2 + ∂ 2 u ∂z 2 + α x ∂u ∂x = G(x, y, z)

Conclusion
In this paper, we have developed a highly accurate 3D problem solver. This has been achieved by the discretization of two of the spatial variables and the construction of a continuous BVM via the interpolation and collocation approach for solving the resulting semidiscretized system. The results given in "Numerical examples" show that the approach is highly efficient and highly accurate.
Authors' contributions SJ proposed the algorithm. TB developed, analyzed and implement the methods. SJ also contributed to the implementation of the methods. Both authors read and approved the final manuscript.