Hybrid Spectral Difference Methods for an Elliptic Equation

1 Introduction

Over the last two decades, high-order methods have received considerable attention because of their potential in delivering higher accuracy with lower cost than low-order methods. Many types of high-order methods have been developed to deal with a diverse range of problems [26]. See also conference proceedings on spectral and high-order methods [18].

Finite difference (FD) methods are amongst the most popular methods solving differential equations in science and engineering. In an effort to develop high-order methods, compact finite difference schemes are proposed for various problems [19, 9, 24, 23]. FD methods have low geometric flexibility, in particular, for high-order schemes.

It is well known that finite element formulations allow high-order approximations with high geometric flexibility. But, their implementation is not so straightforward. Recently, nodal discontinuous Galerkin (DG) methods and hybridized DG methods have been proposed to enhance efficiency of the DG schemes. Nodal DG methods [8] presented by Hesthaven and Warburton provide an efficient way of implementing high-order DG methods for various problems including Poisson equations and Maxwell equations. Hybridized DG methods [3, 12, 13, 14, 15] exploit Lagrange multipliers for hybridization to relax constraints imposed on the cell-interfaces, and static condensation via Schur complement allows efficient implementation while keeping most of the advantages of DG methods. Recently, hybrid high-order methods [5, 4, 2] have been presented in the context of HDG methods.

On the other hand, in order to develop high-order conservative schemes, Wang and his collaborators proposed and studied so-called spectral volume (SV) methods [25, 22] and spectral difference (SD) methods [21]. These methods have been successfully applied to, but limited to, mostly hyperbolic conservation laws. The SV method is designed as a high-order finite volume scheme by considering some local reconstruction idea in the DG schemes. The SD method utilizes a local pseudospectral representation of the solution to obtain spectral-like resolutions. The SD method first defines two sets of points, i.e., the solution points and flux points in each cell. Then, the conservative variables defined at the solution points provide the DOFs. To evolve these DOFs, the differential form of the governing equation is applied at solution points and the divergence of flux is evaluated in terms of values at flux points, cf. [27]. Note that placements of these points are staggered to enhance stability of the approximation to hyperbolic conservation laws. Stability of SD methods is investigated in [1, 10].

In this article, we develop and analyze a locally conservative hybrid spectral difference method (HSD) for the Poisson equation. Similar to the standard SD method, we use two sets of points, i.e., the solution points and flux points in each cell. But, placements of these points are chosen to be axiparallel unlikely to the standard SD method and no grid staggering is needed. Moreover, in our scheme, the differential form of the governing equation is applied at interior solution points to provide the cell finite difference, while flux points are related to the interface finite difference and provide global DOFs. After static condensation built in the procedure, a global stiffness matrix system is obtained in the interface unknowns only. The HSD aims to provide an alternative for the finite difference method and the finite element method for partial differential equations of elliptic type.

Figure 1

A Q2-mesh: |R1|=h1×k1, |R2|=h2×k1.

For simplicity, we present our scheme for the Poisson model problem with the Dirichlet boundary condition:

where Ω is an axiparallel polygonal domain. Let 𝒯h be a partition of Ω by rectangular cells. The simplest HSD is briefly explained, based on Figure 1. The hybrid difference is composed of two kinds of finite difference approximations as follows. The 5-point stencil FD approximation of the Poisson equation, for example, at the cell centered point η4 is the cell finite difference:

and the FD approximation of the flux continuity at η5 is the interface finite difference:

Static condensation on cell-centered unknowns reduces the global couplings, resulting in the system of equations in the cell interface unknowns only.

Numerical experiments on high-order hybrid difference methods (hybrid spectral difference) for elliptic and flow problems, including the Poisson equation, Stokes and Navier–Stokes equations, can be found in [11, 16]. The HSD can be viewed as a finite difference version of the hybridized discontinuous Galerkin method [3, 12, 13, 15]. The HSD is comparable with the finite difference method (FDM). The main difference between the FDM and HSD is that the FD formula of a single type is deployed for all interior nodes in the FDM, while the cell finite difference and the interface finite difference are combined in the HSD. The finite difference method is simple to implement and it can solve many physical problems efficiently [6, 20]. The HSD is as easy to implement as the FDM, and it apparently seems to possess several advantages over the FDM [11]. Those advantages are listed below.

1. The method can be applied to nonuniform grids, retaining the optimal order of convergence, and numerical methods with an arbitrarily high-order convergence can be obtained easily.

2. Problems on a complicated geometry can be treated reasonably well, and the boundary condition can be imposed exactly on the exact boundary.

3. The inf-sup stability is obtained without introducing staggered grids [16].

4. The mass conservation property holds in each cell, and flux continuity holds across inter-cell boundaries.

5. The embedded static condensation property considerably reduces global degrees of freedom.

The aim of the current paper is to propose and analyze the HSD for the Poisson equation. Numerical analysis of the Poisson equation is a very first step toward theoretical understanding of the more complicated PDEs that arise from flow problems or elasticity problems. There can be a couple of different approaches in deriving the FD formulas; one based on the Lagrange polynomial interpolation and the other based on the polynomial degrees of precision. Both approaches yield the same FD formulas for the same derivatives. The latter approach was previously used in the study of pseudospectral (PS) methods. Indeed, the PS method can be seen either as the limit of increasing order FD methods, or as approximations by basis functions, such as Fourier or Chebyshev [7].

The paper is organized as follows. In Section 2 we derive one-dimensional finite difference formulas based on the polynomial degrees of precision. In Section 3, a one-dimensional coercivity result is derived in a semi-discrete energy norm. In Section 4 we extend the coercivity result in Section 3 to the two-dimensional case. In Section 5 numerical experiments are performed for the Poisson equations and numerical results confirm our numerical analysis. Moreover, the Poisson equation on a quarter disk is considered for numerical experiments, and the HSDs provide quite good numerical results. Section 6 concludes the paper with some remarks.

2 Finite Difference Formulas

We derive one-dimensional finite difference formulas in terms of degrees of precision.

Let us introduce an increasing sequence of abscissas {η0,η1,…,ηm} in [a,a+h]. For μ=1,2,…,m, consider the problem to find (wj⁢k(μ))j=0,…,m;k=0,…,m such that

for k=0,…,m. For fixed μ and k, equation (2.1) has a unique solution (wj⁢k(μ))j=0,…,m as system (2.1) forms a Vandermonde matrix system. In particular, it follows from (2.1) that

for k=0,…,m.

Throughout the paper, the notation L≲R means that there exists a constant C>0, independent of h, such that L≤C⁢R.

3 One-Dimensional Coercivity Analysis

Consider m+1 abscissas {ξ0,ξ1,…,ξm} with ξ0=-1 and ξm=1 on the reference interval [-1,1] and the corresponding abscissas {η0,η1,…,ηm} on an interval I=[a,a+h] that are given as ηj=ϕ⁢(ξj), where ϕ⁢(x)=h⁢(x+1)/2+a.

Let us introduce a quadrature on the function space C⁢([a,a+h]):

Here, f~=f⁢(ϕ-1). From here and on, the weights {σj}j=1m-1 and the interior abscissas {ξj}j=1m-1 are the Gaussian weights and abscissas on the reference interval [-1,1], respectively. Therefore, the Gaussian quadrature has the polynomial degree of precision 2⁢m-3.

Let us introduce a bilinear form:

where

Since the bilinear form (3.1) is defined by using the function values of u and v at the given m+1 abscissas on an interval I, we can regard u and v as interpolating polynomials of degree m, that is, u,v∈Pm⁢(I).

As a conclusion the operator 𝒜h is symmetric and elliptic for each order m as long as the interior abscissas are chosen so that the corresponding Gaussian quadrature has the degree of precision 2⁢m-3.

4 Numerical Analysis in Two Dimensions

Let Ω be a rectangular domain with the boundary Γ. In the hybrid difference method we need a rectangular partition 𝒯h of the domain Ω. The partition 𝒯h is shape regular and quasi uniform, and the parameter h represents the maximum diameter of the partition. The notation Kh represents the skeleton of a mesh generation 𝒯h. To define the hybrid difference method we need to generate nodes in each cell, and Figure 2 represents the reference cell configuration (a two-dimensional version of the five abscissa case). In Figure 2 the abscissas {ηi⁢j}i,j=0,…,m (here, m=4) are represented and the nodes {ηi⁢j}i,j=0,m are introduced as virtual abscissas (for convenience of presentation). Let 𝒩⁢(Ω) denote the set of nodes in the closure of a domain and 𝒩⁢(Kh) the set of nodes on its skeleton.

Let C⁢(Ω) be the space of continuous functions in Ω and introduce an equivalence relation u≡v for u,v∈C⁢(Ω) if u⁢(η)=v⁢(η) for all η∈𝒩h⁢(Ω)∪𝒩h⁢(Γ). Denote by Ch⁢(Ω) the space of its equivalent classes. By using the equivalences of the function spaces we regard C~h⁢(R)=Qm⁢(R), where Qm⁢(R) is the space of polynomials in R of degree less than or equal to m in each variable. Then we define

Figure 2

A Q4-mesh: |R|=h1×h2. The point values at the four vertices {ηj⁢k:i,j=0,4} are not used in computation.

Let us introduce the FD approximations for the Laplacian and normal derivatives:

Here, the superscript h in Δh and ∂νh represents the level of discretization as well. The Gaussian quadratures on the reference cell R and its boundary ∂⁡R are defined as follows:

There exists a natural composite quadrature defined on the whole domain Ω as follows:

Let us consider an elliptic problem with the Dirichlet boundary condition:

Then, the solution u satisfies the cell problem

for each R∈𝒯h and it satisfies the interface condition

Here, ∂ν⁡u:=∂⁡u∂⁡ν and ν and ν′ are the outward unit normal vectors from R and R′, respectively.

The hybrid spectral difference method (HSD) is composed of two kinds of finite difference approximations. The cell finite difference solves the cell problem:

and the interface finite difference solves the interface condition

Of course, we require uh⁢(p)=0 for p∈𝒩⁢(Kh)∩Γ. The HSD, (4.2) and (4.3) can be unified in a variational form:

It is easy to see that the exact solution satisfies

The continuous function space with zero trace is represented by C0⁢(Ω). Introducing eh=uh-u, we have

The following theorem is essential for unique solvability of the HSD (4.4).

Let us introduce the norms

and the energy semi-norm

where

Figure 3

N-stacks of rectangles Sj and lines Lk for Q3.