International Journal of Differential Equations

Volume 2018, Article ID 9216703, 8 pages

https://doi.org/10.1155/2018/9216703

## On Solving the Poisson Equation with Discontinuities on Irregular Interfaces: GFM and VIM

^{1}Faculty of Industrial Engineering, Mechanical Engineering and Computer Science Department, University of Iceland, Reykjavik 107, Iceland^{2}Mechanical Engineering Department, University of California, Santa Barbara, CA 93106, USA^{3}Computer Science Department, University of California, Santa Barbara, CA 93106, USA

Correspondence should be addressed to Ásdís Helgadóttir; si.ih@ehsidsa

Received 15 May 2018; Accepted 3 October 2018; Published 17 October 2018

Guest Editor: Manoj Kumar

Copyright © 2018 Ásdís Helgadóttir et al. 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

We analyze the accuracy of two numerical methods for the variable coefficient Poisson equation with discontinuities at an irregular interface. Solving the Poisson equation with discontinuities at an irregular interface is an essential part of solving many physical phenomena such as multiphase flows with and without phase change, in heat transfer, in electrokinetics, and in the modeling of biomolecules’ electrostatics. The first method, considered for the problem, is the widely known Ghost-Fluid Method (GFM) and the second method is the recently introduced Voronoi Interface Method (VIM). The VIM method uses Voronoi partitions near the interface to construct local configurations that enable the use of the Ghost-Fluid philosophy in one dimension. Both methods lead to symmetric positive definite linear systems. The Ghost-Fluid Method is generally first-order accurate, except in the case of both a constant discontinuity in the solution and a constant diffusion coefficient, while the Voronoi Interface Method is second-order accurate in the -norm. Therefore, the Voronoi Interface Method generally outweighs the Ghost-Fluid Method except in special case of both a constant discontinuity in the solution and a constant diffusion coefficient, where the Ghost-Fluid Method performs better than the Voronoi Interface Method. The paper includes numerical examples displaying this fact clearly and its findings can be used to determine which approach to choose based on the properties of the real life problem in hand.

#### 1. Introduction

The Poisson equations with discontinuities across irregular interfaces emerge in applications such as multiphase flows with and without phase change, in heat transfer, in electrokinetics, or in the modeling of biomolecules’ electrostatics. Several numerical methods have been proposed to solve this system, each with their own advantages and disadvantages. One approach is in the context of Discontinuous Galerkin methods, an extension of the finite element method (FEM); e.g., see [1–12] and the references therein. Finite element methods lead to symmetric positive definite linear systems that can be efficiently solved with fast iterative solvers [13]. In addition, FEM-type approaches can derive and use* a priori* error estimates to refine the mesh where higher resolution is needed. However, FEM-type methods rely on the quality of the underlying mesh, which is often difficult to obtain in cases where the irregular domain undergoes large deformations. In this case, it is challenging to generate a mesh with elements that pass a quality measure needed to ensure accurate solutions.

Differential quadrature method is worth mentioning as presented in [14, 15]. There Lagrange interpolation and/or modified cubic B-splines are used, depending on boundary conditions, to approximate the solution to two-dimensional nonlinear hyperbolic partial differential equations. Irregular interfaces and discontinuities were, however, not studied there.

Another approach is within the context of finite difference methods (FDM); e.g., see [16–25] and the references therein. For finite difference methods the grid is Cartesian (uniform or adaptive), which leads to a straightforward grid generation process. However, interfaces must be represented by other means and the treatment of boundary conditions requires additional considerations.

In order to capture the interface and enforce the correct boundary conditions, several approaches have been explored. The immersed boundary method (for example, [26–30]) uses the -formulation that smears out the solution profile across the interface. This produces algorithms that are straightforward to implement since they are similar to solving the same equations on a regular domain. However, the smearing of the solution introduces errors near the interface. The immerse interface method (IIM) (see, for example, [31–37]) is a sharp interface method that leads to second-order accurate solution, albeit it is not a robust second-order method since it reaches its order by minimizing the truncation error. The method leads to neither symmetric nor positive definite linear systems and the application of the immerse interface method may be difficult in three spatial dimensions.

In [20], Liu* et al. * introduced a Ghost-Fluid Methodology (a finite difference approach) that treats the jump condition in a dimension-by-dimension framework. This leads to a linear system that is symmetric positive definite and the jump conditions only affect the right-hand side of the linear system, leading to an easy-to-implement method. However, the GFM suffers a loss in accuracy due to smearing of the tangential derivative of the discontinuity at the interface caused by this dimension-by-dimension framework, as indicated in [20].

The Voronoi Interface Method (VIM) [38] avoids the loss of accuracy of the GFM of [20] by constructing local Voronoi partitioning near the interface. The cells adjacent to the interface therefore have their faces orthogonal to the fluxes of the solution, hence providing a configuration that can leverage the Ghost-Fluid Methodology to its fullest. The advantage is therefore that the solution is second-order accurate. A drawback is that the solution is computed at the cells’ center of the Voronoi partition and an additional interpolation step is required, if the solution is needed on the original Cartesian mesh. As is the case of GFM of [20], the linear system is symmetric positive definite and only its right-hand side is modified by the jump conditions. Finally, we note that even though the construction of a global Voronoi partition may be difficult and costly, the construction of a local partition, i.e., only for cells that are adjacent to the interface, is straightforward. In [38], the authors used the Voro++ library of [39].

In this paper we highlight the performance of both GFM and VIM. It is worth noting that both methods are unconditionally stable since they are implicit. Both methods have also shown that they can handle high ratios of discontinuity in their coefficients.

#### 2. Governing Equations and Numerical Methods

We consider the Poisson equation with variable coefficient and discontinuities (jumps) across an irregular interface, , which splits the computational domain into two domains and , both in , . The governing equation is where and are given. Here, is bounded from below by a positive constant and is in . This equation is supplemented by jump conditions on the irregular interface, where denotes a jump in across . The functions and are given. Either Dirichlet, Neumann, or Robin boundary conditions can be imposed at the boundaries of the computational domain, . In order to represent the irregular interface, we use the signed distance level-set function, , which is positive inside , negative inside , and zero on . The outward normal to the irregular interface can be computed from the level-set function by

##### 2.1. The Ghost-Fluid Method

A detailed description of the Ghost-Fluid Method is given in [20], we will explore here the main aspects in two spatial dimensions. The discontinuity in , where refers to the derivative in the normal direction to the interface, is . The discontinuity in , where refers to the derivative in the tangential direction to the interface, is . These equations lead to and . However, in order to devise a method in a dimension-by-dimension framework, the Ghost-Fluid Method smears out the discontinuity in the tangential derivative leading to the simplification and . The discretization at each grid point is then given by where if and are on the same side of the interface or by otherwise, where is the value of of the node in and is the value of of the node in . Here refers to the value of adjacent to the interface in the domain and and are the cells’ sizes in the - and -directions, respectively. The left-hand side thus gives the same symmetric positive definite matrix as the one generated by the standard five-point stencil discretization of the Poisson equation on regular domains and only the right-hand side is altered when a discontinuity occurs. Furthermore, and are only activated if there is a discontinuity present in the local five point stencil. Here, and and the contribution from the left and right grid points to the current grid point and and are the contributions from the bottom and top grid points.

We give the details for and in the -direction, referring the reader to the original paper [20] for the description for and in the -direction, since they follow the same procedure. If and have opposite signs, define , , where is the jump at the grid node and . is then defined as Similarly, if and have opposite signs, define , , and . is then defined as

The Ghost-Fluid Method leads to a symmetric positive definite linear systems that captures the discontinuity in the normal derivative while smearing out the discontinuity in the tangential direction to the interface.

##### 2.2. The Voronoi Interface Method

We present a summary of the method and refer the reader to [38] for a detailed description. For a given set of seeds, the Voronoi cell of a seed is defined as the points of space that are closer to that seed than any other. The union of the Voronoi cells is a tessellation of space, and the first step for the Voronoi Interface Method is to generate the Voronoi mesh associated with the background mesh chosen, which in our case is a uniform Cartesian grid. The seeds are of two types,(i)if a cell of the background mesh is not crossed by the irregular interface , its center is a seed of the Voronoi mesh,(ii)if a cell of the background mesh is crossed by the irregular interface , we locate the projection of its center onto the interface and generate two seeds located on either side of the projected point at a distance .

The Voronoi partition associated with those seeds is constructed using a simple geometric algorithm. We refer the reader to the Voro++ library [39] for an efficient tool to compute Voronoi partitions in both two and three spatial dimensions.

Equation (1) is then discretized on the Voronoi mesh using a finite volume approach. The complete derivation is presented in [38] and leads to the discretization of the interaction between point and point where any variable of the form would be the quantity at point , is the distance between and , is the length of the face between and , , and is the volume of the Voronoi cell associated with point . The discontinuities at the interface only affect the right-hand side of the linear system where the irregular interface is located and the system is symmetric positive definite.

We remark that the solution is provided at the center of the Voronoi cells. If the solution is needed on the original background mesh, then an interpolation step is required. This can be done, for example, with least square interpolations and it does not impact the order of accuracy of the solution as long as the order of the polynomial interpolants is high enough. In this article, we work with second-order polynomial interpolants.

#### 3. Numerical Experiments: Results

We present a pair of two-dimensional numerical examples where the two methods are compared. Both examples have a star shaped irregular interface. The first example has a coefficient that is constant in the whole domain. This example has a constant discontinuity in the solution across the irregular interface but neither a discontinuity in the normal nor tangential derivatives at the interface. The second example has a coefficient that is not constant and has a discontinuity across the irregular interface. This example has a nonconstant discontinuity in the solution and nonconstant discontinuities in the normal and tangential components of the gradient of the solution at the irregular interface.

##### 3.1. Constant -Coefficient

Let us consider in two spatial dimensions in with the level set function . We take when and an exact solution of . For the region , we take and the exact solution to be . A representation of the solution and the Voronoi mesh are given in Figure 1. The comparison of the Ghost-Fluid Method and the Voronoi Interface Method is shown in Table 1 for the maximum error in the solution, in Table 2 for the average error in the solution, in Table 3 for the maximum error of the gradient, and in Table 4 for the average error of the gradient. For the Voronoi Interface Method, we present the errors on both the Voronoi mesh and the Cartesian mesh. For the Voronoi mesh, the gradients are computed on the faces of the Voronoi cells. For example, the solution only experiences a constant discontinuity in its solution and no discontinuities in its normal or tangential derivatives nor is there a discontinuity in the coefficient, which is constant in the entire domain. Both methods give second-order accuracy in the solution in the - and -norms and the solution’s gradient in the -norm. However, the Voronoi method gives first-order accuracy for the gradient of the solution in the -norm while the GFM gives second-order accuracy for such simplified problems. A likely explanation for this difference is that the fluxes in the volume of fluid derivation for the Voronoi Interface Method are not necessarily computed at the center of the faces between two points. We also observe that the interpolation from the Voronoi mesh to the Cartesian mesh does not impact the solution itself but affects its gradient. However, the order of accuracy is conserved. We conclude that for such simple problems the GFM is preferable.