On a direct method for solving Helmholtz's type equations in 3-d rectangular regions

A direct solution method for solving elliptic pde's of the type kx(z) · ∂2ϕ/∂x2 + ky(z) · ∂2ϕ/∂y2 + kz · ∂2ϕ/∂z2 + σ(z) · ϕ = f(x, y, z) in 3D parallelepipeds with kz = const and kx(z), ky(z), σ(z) continuous functions of z, is presented. The spatial derivatives are approximated using the Hermite approach (Mchrstellen-verfahren) with O(h6) truncation error for Dirichlet boundary conditions or for periodic solutions of the problem. For Neumann conditions, it seems that in order to retain the direct character of the numerical algorithm employed, one should approximate the first spatial derivatives on the boundary by means of conventional schemes having a truncation error of O(h3) type rather than O(h6) which accordingly reduce the overall accuracy of the results. Despite the substantial reduction of the overall accuracy for Neumann conditions, this case has not been exluded, because the structure of the difference equations remains invariant for problems in which instead of known values of first-order normal derivatives at the boundaries, these very boundaries constitute symmetry planes of the solution. This feature allows a direct solution method to be used for such a problem, whereas the O(h6) truncation error of the difference schemes employed is retained. The given pde is discretised on a three-dimensional grid and the set of difference equations is formulated as a linear system of matrix equations whose solution is found by a suitable decomposition of unknows based on knowledge of the eigenvalues and eigenvectors of simple tridiagonal matrices. A hint for extending the applicability of the method—by means of a coordinate transformation—in cylindrical domains with an annular cross section, is also given.