On the calculation of the interaction potential in multiatomic systems

A numerical method for finding the potential of a multiatomic system in the real space is proposed. A distinctive feature of this method is the decomposition of the electron density $\rho$ and the potential $\varphi $ into two parts $\rho={{\rho }_{0}}+\hat{\rho }$ and $\varphi = {{\varphi }_{0}}+\hat {\varphi }$, where ${{\rho }_{0}}$ is the sum of the spherical atom densities and the potential ${{\varphi }_{0}}$ is generated by the density ${{\rho }_{0}}$. The potential $\hat\varphi$ is found by solving Poisson's equation. The boundary conditions are obtained by expanding the reciprocal distance between two points in a series in Legendre polynomials. To improve the accuracy of the method, the computation domain is decomposed into Voronoi polyhedra, and asymptotic estimates of iterations are used when the characteristic function is replaced by its smooth approximations. Poisson's equation is numerically solved using the two-grid method and the Fourier transform. An estimate $O({{h}^{{\gamma-1}}})$, where $h$ is the grid size and $\gamma$ is a fixed number greater than one, is obtained for the accuracy of the method. The error of the method is analyzed using a two-atom problem as an example.