ANALYSIS OF A FOURTH ORDER SCHEME AND APPLICATION OF LOCAL DEFECT CORRECTION METHOD
ANALYSIS OF A FOURTH ORDER SCHEME AND APPLICATION OF LOCAL DEFECT CORRECTION METHOD
Journal of Applied Mathematics & Informatics. 2014. Sep, 32(3_4): 511-527
• Received : July 20, 2013
• Accepted : September 24, 2013
• Published : September 28, 2014
PDF
e-PUB
PPT
Export by style
Share
Article
Author
Metrics
Cited by
TagCloud
ALI ABBAS

Abstract
This paper provides a new application similar to the Local Defect Correction (LDC) technique to solve Poisson problem − u ′′( x ) = ƒ ( x ) with Dirichlet boundary conditions. The exact solution is supposed to have high activity in some region of the domain. LDC is combined with a fourth order compact scheme which is recently developed in Abbas (Num. Meth. Partial differential equations, 2013). Numerical tests illustrate the interest of this application. AMS Mathematics Subject Classification : 65N55, 65N55.
Keywords
1. Introduction
The Local Defect Correction (LDC) method was introduced byW. Hackbusch [8] for solving elliptic boundary value problem. LDC is a domain decomposition technique in which the local domain fully overlaps the global one. It is a generic iterative algorithm of multiscale type for the resolution of discrete problem. Two or more grids of calculation can be used. We refer to [6] for applications to problems in combustion and numerical simulations of the flow and heat transfer in a glass tank. An analysis of the LDC technique in combination with finite difference discretizations is presented in [4] . In [7] the method is extended to include adaptivity, multilevel refinement, domain decomposition and regredding. The LDC method is combined with finite volume discretizations in [5] and finite elements discretizations in [1 , 2] .
Let’s briefly outline the basic version of the LDC technique. Consider the elliptic boundary value problem
PPT Slide
Lager Image
where L is a linear elliptic differential operator, ƒ and g are the source term and Dirichlet boundary condition, respectively. To discretize (1), we first choose a global coarse grid (grid spacing H ), which we denote by Ω H . Let LH be a discrete operator approximating the continuous operator L . An initial approximation
PPT Slide
Lager Image
i = 0, on Ω H is obtained by solving the system
PPT Slide
Lager Image
In (2), the vector ƒ H contains the source term ƒ and the Dirichlet boundary condition g as well. Suppose that LH is invertible. Moreover, suppose that the exact solution u ( x, y ) of (1) has a high activity region in some (small) part of the domain. We select a subdomain Ω l ⊂ Ω such that the high activity region of u is contained in Ω l . The subdomain Ω l is discretized by a local fine grid (grid spac-ing h ). We denote it by
PPT Slide
Lager Image
. The fine grid
PPT Slide
Lager Image
is built such that
PPT Slide
Lager Image
This means that the coarse grid points that lie in the region of refinement are also points of
PPT Slide
Lager Image
. In order to formulate a discrete problem on
PPT Slide
Lager Image
, we should define artificial boundary conditions on Ґ, the interface between Ω l and Ω \ Ω l , ( Figure 1 ). Actually, we use an interpolation operator Ph,H to obtain artificial boundary conditions on Ґ. The operator Ph,H gives the values of the fine grid points on Ґ h by interpolating the values of the coarse grid points on the interface Ґ H .
PPT Slide
Lager Image
Discretization in one dimension. The big points represent the coarse grid ΩH, the small points represent the fine domain .
On Ω l \ Ґ, we obtain boundary conditions using the Dirichlet boundary condi-tions g in (1) b . On the fine grid
PPT Slide
Lager Image
we consider the following discrete problem, at iteration i = 0,
PPT Slide
Lager Image
where
PPT Slide
Lager Image
In (3), i represents an index of iteration and
PPT Slide
Lager Image
is a square matrix representing the dependence of the fine grid solution on the artificial boundary conditions on Ґ. In (3), the matrix
PPT Slide
Lager Image
(assumed to be invertible) is a discrete approximation of the operator L on the subdomain Ω l . In the right term in (3),
PPT Slide
Lager Image
contains the source term ƒ and the Dirichlet boundary condition g on Ω l \ Ґ given in (1). Notice that Ґ ∩ Ω can be empty as in the figure 1 . If we know the values of the defect dH = LH ( u | ΩH )− ƒH , we can use it to find more accurate approximation on the coarse grid. This can be obtained by replacing dH in the right hand side of (2). However, we do not know the exact solution of the exact continuous problem then we cannot calculate dH . We will use
PPT Slide
Lager Image
calculated on the fine grid to approximate dH . Indeed, for all i , index of iteration, the function
PPT Slide
Lager Image
is the function on the coarse grid defined by:
PPT Slide
Lager Image
For each index of iteration i we define
PPT Slide
Lager Image
by
PPT Slide
Lager Image
Actually,
PPT Slide
Lager Image
provides an estimate of the local discrestization error at all points of
PPT Slide
Lager Image
. We observed numerically it is better to use the approximation (6) on a proper subset of
PPT Slide
Lager Image
only. In particular, the points of the coarse grid near Ω l should be excluded. This leads to introduce a “ Safety region ” denoted by
PPT Slide
Lager Image
( Figure 2 ). The estimation of the local error of discretization of the coarse grid is placed in the right hand side of the equations corresponding to coarse grid points belonging to
PPT Slide
Lager Image
only. Then we apply the correction step on the coarse
PPT Slide
Lager Image
PPT Slide
Lager Image
The small points consist of points of Ωh, the big ones are the points of the safety region .
PPT Slide
Lager Image
As (7) contains estimations of the local error of the discretization on the coarse grid, we expect
PPT Slide
Lager Image
to be a more accurate approximation than
PPT Slide
Lager Image
Then we obtain better boundary conditions on Ґ and a solution on the fine local grid by (3) with i = 1. By performing the iterations on the index i , we obtain the following algorithm:
Algorithm 1.1. Initialization
- Solve the problem(2) on the coarse grid Ω H with i = 0. We obtain the vector
PPT Slide
Lager Image
- Solve the problem(3) on the fine local grid
PPT Slide
Lager Image
with i = 0 . We obtain the vector
PPT Slide
Lager Image
Iteration i = 1, 2...
- Compute
PPT Slide
Lager Image
by (5).
- Compute
PPT Slide
Lager Image
by (6).
- Re-solve the problem (7) on the coarse grid Ω H . We obtain the vector
PPT Slide
Lager Image
- Solve the problem (3) on the local fine grid
PPT Slide
Lager Image
We obtain the vector
PPT Slide
Lager Image
In practice, one iteration is enough to obtain a satisafactory approximation on the composite grid Ω H,h which is the grid formed by the union of the coarse grid and the fine grid,
PPT Slide
Lager Image
The algorithm 1.1 is an elementary version of the LDC method. Several generalizations are possible:
• Use of several fine grids. Notice that the local problems are independent each other and can be solved simultaneously. Recursive refinement where a local fine grid can be a coarse grid for another local fine grid.
• Use of different discretizations. We refer to [6] for a detailed analysis of the behavior of the convergence for the Poisson problem solved with classic five points finite difference scheme.
2. Main results
2.1. Analysis of a fourth order scheme. We recall and analyse a new scheme called hermitian Box-scheme (HB-scheme), for the Poisson problem in one dimension based on the previous work [12] . This scheme appears to be fourth order accurate in practice. It has been successfully extended to two dimensions in [9] and three dimensions in [10] . Let’s recall the principle of this scheme. Consider the one-dimensional Poisson problem on the interval Ω = ( a, b ) with length L = b a ,
PPT Slide
Lager Image
Equation (8) is recast in mixed form:
PPT Slide
Lager Image
We lay out on Ω a regular grid xj = a + jh , 0 ≤ j N with stepsize h = L/N , ( Figure 3 ). The unknowns are denoted by uj u ( xj ) and ux,j u ′( xj ). The vectors U , Ux ∈ ℝ N−1 stand for the unknowns at internal points,
PPT Slide
Lager Image
Grid in dimension 1. The two boundary points x0, xN (denoted by “◦”). The points xi, i = 1, · · · ,N − 1 are the interior points (denoted by “•”).
PPT Slide
Lager Image
In analogy to the original box-scheme [13] , the HB-scheme is deduced from the integration of the two equations (9) a,b on a box Kj =] x j-1 , x j+1 [ of length 2 h .
Suppose given the averaged values of the source term ƒ ( x ) on the box Kj ,
PPT Slide
Lager Image
In practice, (11) can be approximated using Simpson formula
PPT Slide
Lager Image
The conservation equation (9) a becomes
PPT Slide
Lager Image
Second, the equation (9) b is integrated on the box Kj . This yields
PPT Slide
Lager Image
Approximating the integral in the left-hand side of (14) by Simpson formula suggests the following fourth-order hermitian approximation
PPT Slide
Lager Image
Equations (15) require approximations of the derivatives on the boundaries. Here, we consider the following third-order approximations
PPT Slide
Lager Image
Relations (16) are obtained using Taylor expansions. In summary, the equations (13), (15), (16) with Dirichlet boundary conditions translate to the following HB-scheme: Find u = ( ui ) 0≤iN and v = ( vi ) 0≤iN solution of:
PPT Slide
Lager Image
This HB-scheme is found numerically to be fourth order accurate as it is shown in the numerical tables.
Lemma 2.1. Suppose that Π 0 ƒj is approximated by Simpson formula:
PPT Slide
Lager Image
and ( uj ) j , ( ux,j ) j are two sequences verifying
PPT Slide
Lager Image
Then ( uj ) j verifies
PPT Slide
Lager Image
Proof . Denote by Ej , j , the j-th equation of (19) a . Using the sum
PPT Slide
Lager Image
and (19) b , we find that the left term of (21) verifies
PPT Slide
Lager Image
The right term of (21) verifies
PPT Slide
Lager Image
where the equality (20) for j ∈ ℤ.
We refer to [11] for a detailed proof of the following lemma concerning the maximum principles of the standard three point Laplacian scheme.
Lemma 2.2. Let u = ( ui ) 0≤iN . Suppose that
PPT Slide
Lager Image
then
PPT Slide
Lager Image
Similarly, if
PPT Slide
Lager Image
then
PPT Slide
Lager Image
The scheme (17) verifies the following discrete maximum and minimum prin-ciples.
Proposition 2.1 (Principles of maximum and minimum). Let u = ( ui ) 0≤iN be the solution of HB-scheme (17). Suppose that ƒi ≤ 0, ∀ 1 ≤ i N − 1, then
PPT Slide
Lager Image
Similarly, if ƒi ≥ 0, ∀ 1 ≤ i N − 1, u verifies
PPT Slide
Lager Image
This means that u attains its maxima and minima on the boundary points x 0 and xN .
Proof . Suppose that ƒi ≤ 0, ∀ 1 ≤ i N − 1 then
PPT Slide
Lager Image
By (20),
PPT Slide
Lager Image
If N is even, using the maximum principle of Lemma 2.2 we obtain
PPT Slide
Lager Image
In the same manner if N is odd, we obtain
PPT Slide
Lager Image
It results for all N ,
PPT Slide
Lager Image
Then to prove (25) it is enough to prove that
PPT Slide
Lager Image
Suppose that
PPT Slide
Lager Image
We have
PPT Slide
Lager Image
This gives u x,0 u x,2 and u x,1 u x,3 . The HB-scheme (17) verifies
PPT Slide
Lager Image
Moreover,
PPT Slide
Lager Image
Therefore,
PPT Slide
Lager Image
But we have u 2 u 0 , then
PPT Slide
Lager Image
therefore
PPT Slide
Lager Image
where (32). Suppose that
PPT Slide
Lager Image
The equations (17) give
PPT Slide
Lager Image
In addition, using (34), u x,N−2 ux,N , therefore
PPT Slide
Lager Image
as u N−1 u N−3 , and uN u N−2 , we get
PPT Slide
Lager Image
We deduce
PPT Slide
Lager Image
where
PPT Slide
Lager Image
Finally,
PPT Slide
Lager Image
where (25). In the same manner we can prove the minimum principle (26).
Corollary 2.1 (existence and uniqueness of solution). The HB-scheme (17) has unique solution
PPT Slide
Lager Image
Proof . We have (17) a linear system with 2 N equations and 2 N unknowns. To prove the existence and uniqueness of the solution, it is enough to prove that Π 0 ƒ i = 0, 1 ≤ i N − 1, implies u = v = 0 . Let u = ( ui ) 0≤iN and v = ( vi ) 0≤iN such that
PPT Slide
Lager Image
Using the discrete maximum principle of proposition 2.1, we obtain
PPT Slide
Lager Image
then
PPT Slide
Lager Image
which proves that the matrix associated to the HB-scheme (17) has nul kernel. As it is square matrix then it is invertible, where the existence and uniqueness of the solution of HB-scheme.
2.2. Matrix form of hermitian Box-scheme. We summarize the finite-difference and matrix notations used in the sequel.
• The tridiagonal matrix
PPT Slide
Lager Image
is
PPT Slide
Lager Image
• The Simpson matrix Ps
PPT Slide
Lager Image
is
PPT Slide
Lager Image
where I is the identity matrix of order N − 1.
• The antisymmetric matrix K
PPT Slide
Lager Image
given by
PPT Slide
Lager Image
• Denoting ( ei ) 1≤iN−1 the canonical basis of ℝ N−1 , the matrices F 1 , F 2
PPT Slide
Lager Image
are defined by
PPT Slide
Lager Image
Let us turn now to the matrix form of the scheme. We use the notation UL = u 0 , UR = uN , ux,L = u x,0 , ux,R = ux,N for the boundary values.
Proposition 2.2. For all linear approximation of the derivatives on the bound-ary Ux,L, Ux,R in terms of U, Ux and the values on the boundary UL, UR, there exist matrices
PPT Slide
Lager Image
such that
PPT Slide
Lager Image
Proof . Any linear approximation of u x,0 in terms of u = ( ui ) 0≤iN , ux = ( ux,i )>) 0≤iN and u 0 , uN can be written in general form as:
PPT Slide
Lager Image
with
PPT Slide
Lager Image
Similarly,
PPT Slide
Lager Image
with
PPT Slide
Lager Image
Multiplying (53) by e 1 and (54) by e N−1 we get
PPT Slide
Lager Image
PPT Slide
Lager Image
Using the fact that
PPT Slide
Lager Image
and
PPT Slide
Lager Image
we can write (56) as
PPT Slide
Lager Image
Finally to complete the proof, let
PPT Slide
Lager Image
be the matrices
PPT Slide
Lager Image
We claim that (16) translates into matrix form as [9]
PPT Slide
Lager Image
where the matrices A, B, C are given by
PPT Slide
Lager Image
and the vector of derivative is given in terms of U and the boundary conditions by
PPT Slide
Lager Image
where the matrices D and ε are
PPT Slide
Lager Image
Using (59) and (61), Ux can be eliminated which gives the expression of the HB-scheme in the sole unknown U as
PPT Slide
Lager Image
where F = [Π 0 ƒ 1 , · · · , Π 0 ƒ N-1 ] T 0 ƒj is given in (11)). The matrices
PPT Slide
Lager Image
are
PPT Slide
Lager Image
If needed, the gradient is recovered as a postprocessing by (61).
2.3. Two-grid refinement algorithm. Here, a multiscale technique is com-bined with the hermitian Box-scheme. We consider the problem
PPT Slide
Lager Image
The domain Ω is discretized with stepsize H = 1/ N such that the discrete coarse points are
PPT Slide
Lager Image
The coarse grid is Ω H . The unknowns vectors are UH and
PPT Slide
Lager Image
such that
PPT Slide
Lager Image
To avoid any confusion between variables, we perform a change of notations for the sequel
PPT Slide
Lager Image
The hermitian Box-scheme corresponding to the problem (65) has the following matrix form
PPT Slide
Lager Image
with F = F = [Π 0 ƒ 1 , Π 0 ƒ 2 , Π 0 ƒ N-1 ] T . By comparing with (2), let
PPT Slide
Lager Image
By resolving (67) we obtain the solution vector
PPT Slide
Lager Image
i = 0, on Ω hc and the global vector solution
PPT Slide
Lager Image
PPT Slide
Lager Image
The domain Ω l has been choosen such that Ω l = ( a 1 , b 1 ), with a 1 , b 1 two fixed scalars between a and b, a a 1 < b 1 b . The choice of a 1 and b 1 is adapted such that for all N there exist i 0 and j 0 between 1 and N − 1 such that
PPT Slide
Lager Image
that is to say a 1 and b 1 are two points of the coarse grid Ω hc . In order that
PPT Slide
Lager Image
has the same number of points as Ω H , we choose a 1 and b 1 such that
PPT Slide
Lager Image
Suppose that we know the values g a1 = u ( a 1 ) and g b1 = u ( b 1 ) where u ( x ) is the exact solution. Hence, on the domain Ω l we can write
PPT Slide
Lager Image
with ul = u | Ωl and ƒl = ƒ | Ωl . The domain Ω l is discretized with stepsize
PPT Slide
Lager Image
such that the discrete fine points are
PPT Slide
Lager Image
with
PPT Slide
Lager Image
Notice the equality between N and Nƒ when
PPT Slide
Lager Image
The fine grid is
PPT Slide
Lager Image
The unknown vectors are
PPT Slide
Lager Image
such that
PPT Slide
Lager Image
and
PPT Slide
Lager Image
1 ≤ i N − 1. The interface between Ω l and Ω is Ω l = { a 1 , b 1 }. Based on the matrix form (63), the hermtian Box-scheme corresponding to (71) is
PPT Slide
Lager Image
with Fl = [Π 0 ƒ 1 , Π 0 ƒ 2 , ..., Π 0 ƒ N-1 ] T ∈ ℝ N-1 . The problem is that we do not know the values g a1 and g b1 on Ω l in (71). To approach these values, we replace it by artificial boundary conditions
PPT Slide
Lager Image
with i 0 and j 0 are the indices of the points a 1 and b 1 in the coarse grid defined in (70). Here, in particular, the artificial boundary Ω l = { a 1 , b 1 } is composed of two points of the coarse grid Ω hc then there is no need to interpolate therefore the matrix
PPT Slide
Lager Image
is zero. By comparing (72) and (3), let
PPT Slide
Lager Image
By solving (72) with i = 0 we get
PPT Slide
Lager Image
then we obtain the approximation
PPT Slide
Lager Image
Following (5) the vector
PPT Slide
Lager Image
verifies
PPT Slide
Lager Image
The defect vector
PPT Slide
Lager Image
is calculated by
PPT Slide
Lager Image
Finally, the algorithm is the following:
Algorithm 2.1. Initialization
- Solve the problem (67) on the coarse grid Ω hc . The vector solution is denoted by
PPT Slide
Lager Image
- Solve the problem (75) on the local grid
PPT Slide
Lager Image
The vector solution is denoted by
PPT Slide
Lager Image
Iteration i = 1, 2, ...
- Compute the vector
PPT Slide
Lager Image
by (76).
- Compute the defect vector
PPT Slide
Lager Image
by (77).
- Solve on the coarse grid
PPT Slide
Lager Image
with ε > 0. The solution vector is denoted by
PPT Slide
Lager Image
- Solve the problem (75) on the fine grid. The vector solution is denoted by
PPT Slide
Lager Image
2.4. Numerical results. In this part, we display some numerical results prov-ing the interest of the technique presented in the last section. The calculation is performed using Matlab. The CPU time is calculated by the functions tic and toc. Notice that the CPU time can be extremely reduced using FFT or another programming language [3] . The algorithm used so far is the algorithm 2.1 after only one iteration. We intend to study in another paper the behavior of the convergence of this algorithm. In the numerical tables, uex , ux,ex are the exact solution and derivative, u , ux are the computed solution and derivative on the uniform grid and
PPT Slide
Lager Image
are the computed solution and derivative on the composite domain. We use the discrete norm L to estimate the error:
PPT Slide
Lager Image
The convergence rate is calculated using the formula
PPT Slide
Lager Image
where eN is the error obtained on a mesh of size N .
Error and convergence rate for Test 1 on uniform grid
PPT Slide
Lager Image
Error and convergence rate for Test 1 on uniform grid
Error and convergence rate for Test 1 on composite grid
PPT Slide
Lager Image
Error and convergence rate for Test 1 on composite grid
Test 1: We consider the Gaussian function
PPT Slide
Lager Image
This function has a high activity region around
PPT Slide
Lager Image
The results obtained on the composite grid are compared with those computed by the HB-scheme on the uniform grid with fine stepsize hƒ . Observe in tables 1 and 2 the same accuracy obtained by the two approaches with less points on the composite grid.
Test 2: Suppose that the exact solution is
PPT Slide
Lager Image
on Ω = (0, 1). This function has a deep gradient around
PPT Slide
Lager Image
We choose
PPT Slide
Lager Image
The numerical results are given in tables 3 and 4.
3. Conclusion
This paper provides a new application of the multiscale technique combined with a new high order scheme called HB-scheme. The high order accuracy of this scheme is particularly interesting in applications where physical fields have to be calculated as an outcome of a potential solution of a Poisson problem (electromagnetism, gravitation). The main features of this work are the high order of precision on the boundaries and the calculation of the derivatives. This application illustrates the gain in memory and CPU time using this technique as it is shown in the numerical tables. For instance, we have obtained for Test 1 the same precision using 3072 points instead of 4096. This gain is expected to be more interesting in higher dimensions. A generalization of this application to higher dimensions using cubic spline interpolations, based on the previous results in [9 , 10] will be a part of future work.
Error and convergence rate for Test 2 on uniform grid
PPT Slide
Lager Image
Error and convergence rate for Test 2 on uniform grid
Error and convergence rate of Test 2 on composite grid
PPT Slide
Lager Image
Error and convergence rate of Test 2 on composite grid
BIO
Ali Abbas received Ph.D at University of Paul Verlaine, France. Since 2012, he has been at Lebanese International University. His research interests include finite differences and finite volumes methods.
Department of Mathematics, Lebanese International University, Rayak-Bekaa, Lebanon.
e-mail: ali.abbas@liu.edu.lb
References
Huang P. , Feng X. , Su H. (2012) A new defect-correction method for the stationary Navier-Stokes equations based on local Gauss integration Mathematical Methods in the Applied Sciences 35 1033 - 1046    DOI : 10.1002/mma.1618
Huang P. , Feng X. , He Y. (2013) Two-level defect correction Oseen iterative stabilized finite element methods for the stationary Navier-Stokes equations Applied Mathematical Modelling 37 728 - 741    DOI : 10.1016/j.apm.2012.02.051
Abbas A. 2011 Schémas compacts hermitiens- Algorithmes rapides pour la discrétisation des équations aux dérivées partielles Univ. Paul Verlaine - Metz
Ferket Peter J. J. , Reusken Arnold A. (1996) A finite difference discretization method on composite grids Computing 56 343 - 369    DOI : 10.1007/BF02253460
Hof Bas van ’t 1998 Numerical aspects of laminar flame simulation Eindhoven University of Technology - Eindhoven
Anthonissen M.G.H. 2001 Local Defect correction Techniques : Analysis and Application to combustion Eindhoven University of Technology
Ferket P.J.J. , Reusken A.A. (1996) Further Analysis of the Local Defect Correction Method Computing 56 117 - 139    DOI : 10.1007/BF02309341
Hackbush W. (1984) Local Defect Correction Method and Domain Decomposition Techniques Computing 5 89 - 113
Abbas A. , Croisille J-P. (2011) A fourth order Hermitian Box-Scheme with fast solver for the Poisson problem in a square J. Sci. Comput 49 239 - 367    DOI : 10.1007/s10915-010-9458-y
Abbas A. 2013 A fourth order Hermitian Box-Scheme with fast solver for the Poisson problem in a cube, Numer. Methods. Partial differential equations
Ben-Artzi M. , Croisille J-P. , Fishelov D. 2013 Navier-Stokes equations in planar do-mains Imperial College Press
Croisille J-P. (2006) A Hermitian Box-Scheme for one-dimensional elliptic equations - Applica-tion to problems with high contrasts in the ellipticity Computing 78 329 - 353    DOI : 10.1007/s00607-006-0181-3
Keller H. B. 1971 A new difference scheme for parabolic problems Academic Press