Helmholtz equation has attracted much attention in many fields such as electromagnetic cavity scattering problems  , wave propagation  and acoustic problems  . Many methods have been proposed to solve the Helmholtz equation in general Cartesian coordinates, such as finite difference method  , finite element method  and other methods  . This equation is important for both theory and applications. Its theoretical significance has a variable coefficient under the first derivative, which renders the existing fourth-order compact finite difference methods inapplicably. From the standpoint of applications, the Helmholtz equation in polar coordinates appears in many scattering problems.
In general Cartesian coordinates, high-order methods for solving the Helmholtz equation have been well developed since the accuracy is related to the amount of the grid points with the wave number increases. Manohar et al.  proposed second- and sixth-order finite difference schemes for solving the Helmholtz equation which requires the derivatives and much computational cost. Nabavi et al.  further developed a compact nine-point sixth-order finite difference scheme and obtained a sixth-order approximation for the Neumann boundary condition. Sutmann  developed a new compact sixth-order finite difference scheme of sixth-order for three-dimensional Helmholtz equations.
In polar coordinates, the symmetric problem can be described more concisely. A finite difference method was proposed to solve the parabolic equation in polar coordinates in  . Britt et al.  constructed a high-order compact difference scheme for the Helmholtz equation. Su et al.  proposed a fourth-order method for solving Helmholtz equation with discontinuous wave number and Dirichlet boundary condition. A novel compact scheme based on finite difference discretizations and geometric grid has been developed to solve two-dimensional mildly non-linear elliptic equations in polar co-ordinate constituting singular terms  .
The existed works didn’t discuss the Neumann boundary condition and the sparsity of the coefficient matrix. In this paper, we develop a fourth-order scheme for the Helmholtz equation in polar coordinates with Dirichlet and Neumann boundary conditions. The sparse matrix form of the linear system is derived. With the help of the sparsity of the coefficient matrix of the linear system, the computational cost is remarkably reduced. Moreover, we give the approximation for the Neumann boundary condition. The feasibility and the order of method are verified by the Helmholtz equation with Dirichlet and Neumann boundary condition.
The paper is outlined as follows. In Section 2, a fourth-order finite difference scheme and the sparse matrix form for the Helmholtz equation in polar coordinates are derived. In Section 3, the approximation for the boundary condition is obtained. Two numerical experiments of the high-order algorithm are presented in Section 4. The paper is concluded in Section 5.
2. Fourth-Order Finite Difference Scheme
2.1. Fourth-Order Approximation
We consider the following two-dimensional Helmholtz equation in polar coordinates
where k is the wave number, D is a given region and f is a known function. For convenience, we rewrite the above equation as
Define a uniform mesh on the region . , . denotes the fourth-order solution of u at point .
First, we employ a fourth-order approximation on the left side of Equation (2) and obtain
Then approximating the derivatives of u with the second-order accuracy by central differences, we have
By substituting Equation (5) into Equation (4), we can obtain a nine-point stencil approximation for .
We now turn to the fourth-order approximation for the left side of Equation (3). Applying a standard central difference operator to the left item of Equation (3), we have
and the error
Differentiating both sides of (3) with respect to and combing (7), we have
Similarly, all derivatives of with respect to can be approximated as follows
By substituting Equation (9) into Equation (8), the fourth-order finite difference scheme for is obtained.
Therefore, combing Equation (4) and Equation (8), we have the overall fourth-order finite difference form for the Helmholtz equation in polar coordinates
Replacing the derivatives of and by Equations (5) and (9), we can derive
For simplicity, we rewrite the above Equation (11) in the following form
2.2. The Sparse Matrix Form of the Scheme
The relationship between the nine points can be illustrated in three parts in matrix form , which represent the horizontal, vertical, diagonal relationships respectively,
and denotes the Kronecker product, is the identity matrix. F is the right side item containing and the boundary conditions which will be discussed in the following section.
3. Boundary Condition
Implementation of Dirichlet boundary conditions at , is straightforward. And the right side item of Equation (13) can be written as
and are vectors in dimensions.
Implementation of the Dirichlet boundary conditions at , is more complicated. We consider the following Neumann boundary condition and it can be extended to the general cases
where is a constant, is a given function.
The ghost points are utilized to derive the difference scheme on the boundary can be written as
Therefore, the Neumann boundary condition can be approximated as
We can obtain
Substituting j for in Equation (12) gives
Combining Equations (18) and (20), we can eliminate and obtain
Therefore, collaborating Equations (13) and (21), the global system can be written as follows
We can observe the sparsity of the global system of the Helmholtz equation with the Neumann boundary condition in Figure 1.
4. Numerical Experiments
4.1. Example 1
We first consider the problem with Dirichlet boundary condition in polar coordinates as follows
and the exact solution is .
As we can see from Figure 2 and Figure 3, the numerical solution derived by the proposed method is highly consistent with the exact solution. Moreover, we depict the numerical solution in Cartesian coordinates in the right side of Figure 3. Furthermore, in order to test the computational order of the proposed method, we give the error between the numerical solution and the exact solution with different grid points and k in Table 1 and Table 2. The order is calculated by the following equation
It can be clearly seen from Table 1 and Table 2 that the finite difference scheme can reach the fourth-order when the wavenumber k is relatively small. As the mesh is refined and the number of grid points increases, the error becomes smaller and the accuracy tends to be fourth-order and gradually stabilized. When the wavenumber k increases, the error becomes oscillatory. Table 3
Figure 1. The sparsity of the global system with Neumann boundary condition.
Figure 2. The numerical solution (left) and the exact solution with and .
Figure 3. The error inside the region (left) and the numerical solution in Cartesian coordinates (right) with and .
Table 1. The errors and order of the proposed method with and .
Table 2. The errors and order of the proposed method with and .
Table 3. Computational time (s) for solving the Helmholtz equation with .
gives the comparison of the computational time(s) for solving the Helmholtz equation in polar coordinates, where Method I and Method II denote methods with and without the utilization of the sparsity of the coefficient matrix of the linear system. It can be observed that with the help of the sparsity of the coefficient matrix of the linear system, the computational cost is remarkably reduced.
4.2. Example 2
This example is a Helmholtz equation in polar coordinates with Neumann boundary condition
and the exact solution is . As we can see from Figures 4-6 that the
Figure 4. The numerical solution (left) and the exact solution with and .
Figure 5. The error inside the region (left) and the error on the top boundary (right) with and .
Figure 6. The numerical solution in Cartesian coordinates.
numerical solution is well agreed with the exact solution with and .
In this paper, we propose a high-order fast algorithm for solving the two-dimensional Helmholtz equation with Dirichlet and Neumann boundary conditions in polar coordinates. We develop a fourth-order accurate compact finite difference approximation to the Helmholtz equation. The sparse matrix form for the Helmholtz equation in polar coordinates is obtained which improves the efficiency for the computation process. Two numerical experiments have demonstrated the validity of the fourth-order algorithm.
This research was supported by the Fundamental Research Funds for the Central Universities (No. 2018MS129).
 Zhao, M.L. (2013) A Fast High Order Iterative Solver for the Electromagnetic Scattering by Open Cavities Filled with Inhomogeneous Media. Advances in Applied Mathematics and Mechanics, 5, 235-257. https://doi.org/10.4208/aamm.12-m12119
 Zhao, C. and Liu, T. (2003) Nonlinear Artificial Boundaries for Transient Scalar Wave Propagation in a Two-Dimensional Infinite Homogeneous Layer. International Journal for Numerical Methods in Engineering, 58, 1435-1456. https://doi.org/10.1002/nme.703
 Singer, I. and Turkel, E. (1998) High-Order Finite Difference Methods for the Helmholtz Equation. Computer Methods in Applied Mechanics and Engineering, 163, 343-358. https://doi.org/10.1016/S0045-7825(98)00023-1
 Harari, I. and Hughes, T. (1991) Finite Element Methods for the Helmholtz Equation in an Exterior Domain: Model Problems. Computer Methods in Applied Mechanics and Engineering, 87, 59-96. https://doi.org/10.1016/0045-7825(91)90146-W
 Kashirin, A., Smagin, S. and Taltykina, M. (2016) Osaic-Skeleton Method as
Applied to the Numerical Solution of Three-Dimensional Dirichlet Problems for the Helmholtz Equation in Integral Form. Computational Mathematics
and Mathematical Physics, 56, 612-625.
 Nabavi, M., Siddiqui, M. and Dargahi, J. (2007) A New 9-Point Sixth-Order Accurate Compact Finite Difference Method for the Helmholtz Equation. Journal of Sound and Vibration, 37, 972-982. https://doi.org/10.1016/j.jsv.2007.06.070
 Sutmann, G.
(2007) Compact Finite Difference Schemes of Sixth-Order for the Helmholtz Equation. Journal of Computational and Applied Mathematics, 203, 15-31.
 Zhang, R.S., Lu, G.Z. and Abomakhleb, G. (2017) Finite-Difference Solution of the Parabolic Equation under Horizontal Polar Coordinates. IEEE Antennas and Wireless Propagation Letters, 16, 2931-2934. https://doi.org/10.1109/LAWP.2017.2753220
 Britt, S., Tsynkov, S. and Turkel, E. (2001) Numerical Simulation of Time-Harmonic Waves in Inhomogeneous Media Using a Compact High Order Schemes. Communications in Computational Physics, 9, 520-541. https://doi.org/10.4208/cicp.091209.080410s
 Su, X.L., Feng, X.F. and Li, Z.L. (2016) Fourth-Order Compact Schemes for Helmholtz Equation with Piecewise Wave Numbers in the Polar Coordinates. Journal of Computational Mathematics, 34, 499-510. https://doi.org/10.4208/jcm.1604-m2015-0290
 Jha, N., Monhanty, R.K. and Kumar, N. (2017) Compact-FDM for Midly Nonlinear Two-Space Dimensional
Elliptic BVPs in Polar Coordinate System and Its Convergence Theory. Journal of Computational and Applied Mathematics, 3, 255-270.