Symmetric Boundary Condition for Laplacian on Net of Regular Hexagons

Hexagonal grid methods are found useful in many research works, including numerical modeling in spherical coordinates, in atmospheric and ocean models, and simulation of electrical wave phenomena in cardiac tissues. Almost all of these used standard Laplacian and mostly on one configuration of regular hexagons. In this work, discrete symmetric boundary condition and energy product for anisotropic Laplacian are investigated firstly on general net of regular hexagons, and then generalized to its most extent in twoor three-dimensional cell-center finite difference applications up to the usage of symmetric stencil in central differences. For analysis of Laplacian related applications, this provides with an approach in addition to the M-matrix theory, series method, functional interpolations and Fourier vectors.


Introduction
Hexagonal (Hex) grid methods are of interest in many research studies: (Pickering,1986) on direct method, (Makarov, Mararov & Moskal'kov, 1993) giving a formula without proof, (Bystrytskyi & Mosklkov, 2001) on seven-point method on rectangular grid with explicit form of eigenpairs, (Zhou & Fulton, 2009) with periodic boundary condition (BC), (Heikes & Randall, 1995, part I,II) and (Heikes, Randall & Konor, 2013) on numerical modeling in spherical coordinates, (van Eck & Kors, 2005) on action potential in heart modeling via algebraic method without using diffusion in form of differential equation, (Nickovic,Gavrilov & Tosic, 2002) showing advantages of Hex grids over commonly used square grids for use in atmospheric and ocean models.In the article by (Lee,Tien,Luo and Luk,2014), Hex grid finite difference (FD) methods are derived in a finite volume (FV) approach involving standard Laplacian, and used in the simulation of electrical wave phenomena propagated in two-dimensional reversed-C type cardiac tissues, exhibiting both linear and spiral waves more efficiently than similar computation carried on rectangular FVs.We note these cited works all used standard Laplacian and mostly on one configuration of regular hexagons.
In two-dimensional applications of configurations consisting of (subset of) Cartesian type regular hexagons, we denote the radius of hexagons by r, the height by h(= √ 3 2 r), and the center-to-center distance by d(= 2h).Near a typical center node, P 0 = (x 0 , y 0 ), the six neighbor (center) nodes are P j = (x j , y j ) = (x 0 , y 0 ) + d(cos ξ j , sin ξ j ), ξ j = φ + jπ 3 + π 6 , 1 ≤ j ≤ 6. (1) Here the phase angle, φ, is the configurarion parameter.Two particular instances are called type I (φ = 0) and type II (φ = −π 6 ) for convenience.Hexagon centers in lattices of these two types are indexed as for an orthogonal Cartesian mesh as shown in Table 1, while the geometry and neighborhood of a general Hex FV shown in Table 2. Indexing rules are illustrated in Figs. 1,2,3 and 4. For convenience, we abuse the notations and denote FV centers in a neighborhood (Figs. 3 and 4) by an ordered list, Type I : {P j } 6 j=0 = {P, P N , P NW , P S W , P S , P S E , P NE }, Type II : {P j } 6 j=0 = {P, P NE , P NW , P W , P S W , P NW , P E }. (2) We note for applications that a two-dimensional irregular domain may be approximated by a sequence of (not necessarily Cartesian) nets of hexagons.Actually, our work in numerical modeling of ECG depends on this (Algorithm 1 in (Lee,Tien,Luo & Luk, 2014)).
Table 2. Local geometry at a regular hexagon : six vertices and six neighbor centers Concerning the negative anisotropic Laplacian with positive constant diffusivities D 1 and D 2 , we observe the following.
Lemma 1 (Reflection principle for anisotropic Laplacian.)The two configurations, type I and II regular hexagons centered at the origin together with the anisotropic Laplacian, are convertible from each other by applying reflection with respect to the main diagonal in the xy-plane, and therefore interchanging the two symbol lists (Figs. ) of the indices.
The focus in subsequent discussion is on net of type I hexagons.
We note that spectral analysis of iterative methods solving the discrete anisotropic Laplacian on a net of hexagons seems not as easy as the analysis on square grids (Suli,1993) and (Karaa & Zhang, 2003), since finite trigonometric series is incomplete for the error analysis (even) on a single regular hexagon (McCartin,2002(McCartin, ,2003)).
In the world of differential equations, for example (Strauss, 2008), the term symmetric boundary condition is defined so as to make a (real) operator symmetric.On the other hand, as a practice for long time in the engineering literatures, the phrase symmetric boundary condition may mean, differently, that the computational domain is reduced by halving and the numerical BC on the virtual separating edge is of homogeneous Neumann type : (Kim & Huh, 2000), (Xu & Soares, 2013), and (Pal, Lan, Li, Hirleman & Ma, 2015).Same as such with spherical (reflexive) symmetric boundary condition in similar situations.
We are with the operator-theoretic view (Eqs. (5,20 and 24) in current work).
As for the remaining sections, symmetric boundary condition for the Laplacian is introduced for smooth scalar functions in Section 2, with detailed discussion in Section 3. The theory is generalized and simplified in Section 4 with simple assumption and argument.Discussed in Section 5 are many examples of symmetric boundary condition for Laplacian, including non-product type pairs (generators) for an invariant subspace of some operator on type I Hex grid.

Function Symmetric Boundary Condition for Laplacian
We focus here on standard Laplacian (Eq.( 3) with D 1 = D 2 = 1) and recall the Green's first and second identities, in which Ω is a domain with ∂Ω its piecewise smooth boundary such that these formulas are valid.The outward normal on the boundary is denoted by ∂⃗ n.
Definition 1.A pair of (distinct) functions f, g ∈ C 2 (Ω) satisfies symmetric boundary condition (for the Laplacian), if so that the Laplacian is symmetric on f and g, ⟨ f, ∇ 2 g⟩ = ⟨∇ 2 f, g⟩.
If every pair in a family of functions satisfies the symmetric boundary condition, we say the Laplacian is symmetric on the family.
Notice that, with or without satisfication of the symmetric boundary condition, we may consider symmetrization of the energy product, in which the boundary functional is Assuming the symmetric boundary condition, there may exist simplification of the boundary functional and the energy product.
Example 1.Some particular cases assume Accordingly, the energy product is Two classical examples include pairs both satisfying homogeneous Dirichlet or homogeneous Neumann BC.
Example 2. A pair satisfying Robin BC that provides another example which satisfies the symmetric boundary condition, owing to The corresponding energy product is then positive definite, under an additional assumption that c(x) > 0 on ∂Ω.
Very general results in discrete symmetric boundary condition are presented in Section 4. For the motivation, we discuss next the discretization (Fig. 1) of a Cartesian net of type I regular hexagons.

Discrete Symmetric Boundary Condition on Net of Regular Hexagons
The type I neighborhood topology (Fig. 3) and the geometry (Eq.( 2)) are detailed below.exterior (ghost) boundary interior and exterior with even i For the symmetry of the negative hexagonal seven-point Laplacian ( L 7 f we use backward (in vertical direction) differences to derive Similarly, by interchanging f and g, Setting the goal ⟨L 7 f, g⟩ − ⟨ f, L 7 g⟩ = 0 leads to backward difference version of the symmetric boundary condition, Alternatively, we can make use of forward (in vertical direction) differences.
Also, by interchanging f and g, Taking average of Eqs.(10 and 13), (11 and 14), respectively, leads to To summarize in compact form.Let P ′ represent a member in the neighborhood of P ∈ Ω, that is, symbolically, P ′ ∈ { P N , P NW , P S W , P S , P S E , P NE } Then, in terms of backward differences (using three lower neighbors, Eq. ( 10)), and, in terms of forward differences (using three upper neighbors, Eq. ( 13)), Taking average of these two leads to a compact form of Eq. ( 16), Similarly, with f and g interchanged, we obtain compact forms of Eqs.(11 and 14), and their average, now as a compact form of Eq. ( 17), In summary.
Theorem 2 (Symmetric boundary condition and energy product for negative seven-point Laplacian on Cartesian net of type I regular hexagons.) If the above (implied) symmetric boundary condition is satisfied, then the bilinear form is well-defined and symmetric.
Proof.Taking difference of Eqs.(19 and 18) yields the symmetric boundary condition, while taking average leads to the discrete product.
Remark 1.The Theorem is valid on a general (not necessarily Cartesian) net of type I hexagons, because there is no usage of integral indices in the relevant discussion.Even more general case is applicable with central differencing in a proper setup, as presented in Section 4.
Example 3. (Stencil-truncation.) Assume cell-average Dirichlet BC that the grid data vanish at ghost nodes, f P ′ = g P ′ = 0, ∀P ′ ∈ ∂ e , the symmetric boundary condition (Eq.( 20)) is then satisfied.The simplified energy product is certainly positive-definite.We note the boundary functional is defined with multiplicities at interior boundary nodes.
Example 4. (Torus.)A Cartesian net of hexagons (Fig. 1), as a computational domain with periodic BC, corresponds to a torus with no exterior boundary (ghost) node.The symmetric boundary condition for Laplacian is satisfied.The boundary functional vanishes and the discrete energy product is positive semi-definite.⟨ f, f ⟩ L 7 = 0 if and only if f is a (single) constant with the discrete topology being path-connected.
There are more examples along this line, with the periodic BC replaced by twist BC or by mixing of periodic and twist, resulting in applications of dynamical systems for the real projective plane or a Klein bottle.Open field problems subject to homogeneous Dirichlet or Neumann BC are discussed as corollaries to the general result in the next section.

Symmetric Boundary Condition for Discrete Laplacian by Cell-center Finite Difference
For applications involving discrete (negative) Laplacian in the form with Ω and ∂ e denoting respectively the (disjoint) set of interior grid nodes and ghost nodes, we assume very general assumptions that (i) the discrete neighborhood topology is reflexive so that being-neighbor-to is a symmetric relation among interior nodes, (ii) the central difference (CD) stencil is symmetric, assuming proper orientation (ordering) consistently in all (local) neighbor lists, such that As an example, consider standard Laplacian on a net of type I regular hexagons (Figs. 1 and 3, and Eq. ( 3)), the homogeneous discrete neighborhood is the relation which is symmetric among interior nodes.We abuse the notation slightly and use

N(P) ≡ { N(P), NW(P), S W(P), S (P), S E(P), NE(P)
} or, with implied dependence, The consistently ordered homogeneous symmetric stencil (Lee,2015) is, and, in order, We add in passing, for anisotropic negative Laplacian (Eq.( 3) on type I Hex net, that (Lee,2015) The following relations are helpful.
Lemma 3 We establish a very general result, as follows.
Theorem 4 The symmetric boundary condition for and up to satisfication of which, the resulting energy product is By taking difference of the last two equations and using Eq. ( 22), we obtain which yields the symmetric boundary condition (Eq.( 24)), while taking the average and using Eq. ( 23), with boundary functional Thus ends the proof.
Corollary 5 Suppose every exterior (ghost) boundary node, P ′ ∈ ∂ e , is of homogeneous Dirichlet, homogeneous Neumann or Robin type such that, respectively, either f P + f P ′ = g P + g P ′ = 0, or f P − f P ′ = g P − g P ′ = 0, or = const P ′ , at each pair (P, P ′ ), then the symmetric boundary condition (Eq.( 24)) is satisfied.
Corollary 6 (i) If a mixed type homogeneous BC consists of at least one Dirichlet and some Neumann node(s) with vanishing boundary functional (owing to f P ′ / f P = f P ′ / f P = ±1 at boundary), then the energy product simplifies to and is (symmetric and) positive-definite.
(ii) If pure homogeneous Neumann BC ( f P ′ − f P = g P ′ − g P = 0) holds, then the symmetric boundary condition is satisfied, the boundary functional vanishes, the energy product is reduced (to Eq. ( 28)) and is positive semi-definite such that ⟨ f, f ⟩ L = 0 if and only if f is constant in each (path-)connected component of Ω.
(iii) In case of (function) Robin BC (Eq.( 9)) with c(x) > 0 on ∂Ω, the central difference approximation with 2h being the node-to-node distance, implies Therefore, the energy product (Eq.( 25)) is positive definite, up to c := c (P+P ′ )/2 > 0 in practice.We note the boundary functional is defined with various multiplicities at ghost nodes.

Symmetric Boundary Condition for Five-point Laplacian
In a setup of rectangular grid, the scaled negative five-point Laplacian reads in which 1 ≤ i ≤ n x , 1 ≤ j ≤ n y at interior nodes, and i = 0, n x + 1, j = 0, n y + 1 at ghost nodes.The general theory (Eqs. (24 and 25)) specializes as follows.
Theorem 7 Let the indices run through all interior nodes, 1 ≤ i ≤ n x , 1 ≤ j ≤ n y .
(i) The discrete 2D symmetric boundary condition for five-point Laplacian is (ii) Suppose the above symmetric boundary condition is satisfied, then the following expression ) defines a symmetric bilinear form.
The above theory were actually studied firstly in details, and motivated the discussion of Hex grid case (Section 3) and then the very general case (Section 4).We omit any meta-analysis in deriving Eqs.(30 and 31).Instead, we specialize for a commonly encountered application on rectangular grid.
Theorem 8 (Product-form symmetric boundary condition.)Suppose that the application data are (separable) in product forms, then the symmetric boundary condition (Eq.( 30)) simplifies to which is satisfied if any of the following four conditions holds

Discussion
We give here examples of sets of (basis) vectors which all satisfy pairwisely the symmetric boundary condition for discrete Laplacian.The private BC in each case is indicated when it is convenient.
We consider 2D half-integral nodes with extensions to ghost nodes.Depending on applications, several Cartesian type bases exist as summarized in Table 3, in which the one-dimensional components are defined next.

One-dimensional Fourier Vectors
With (n, k, t k ) denoting (n x , i, x i ) or (n y , j, y j ), the Fourier half-wave Sine and Cosine, and quarter-wave Sine and Cosine vectors are, respectively, The implied boundary values (k = 0, n) satisfy private BCs that We allow the extension that v 0 i = u n+1 i ≡ 0 and refer to Figs. 5 and 6 for several low degree instances of these vectors.
The bases defined by Eq. ( 38) satisfy many properties, among which we mention the following.
(1) (Central-)Even-odd symmetry in half-wave vectors. With In particular, and therefore as twist and periodic BCs, respectively, that (2) Symmetry in quarter-wave vectors.
Orthogonal bases exist in various forms.We state a few.
Lemma 9 Every set in the following generates an orthogonal basis of R n . (1)

Figure 1 .Figure 2 .
Figure 1.Lattice of type I regular hexagons in natural order by columns

Figure 5 .Figure 6 .
Figure 5. Fourier half-wave Sine and Cosine vectors of degrees up to eight

Table 1 .
Lattices of type I and II regular hexagons