Efficient solution of Poisson’s equation in linear combination of atomic orbitals (LCAO) electronic structure calculations

An accurate and efficient method is described for the evaluation of electrostatic contributions in LCAO electronic structure calculations. The charge density ρ(r) is decomposed into ρ(1)(r), a component whose rapid variation near any nucleus reproduces that of ρ(r) to a very good approximation, and a remainder density δρ(r)≡ρ(r)−ρ(1)(r), which is thereby guaranteed to be slowly varying in space. The power of the decomposition resides in the fact that ρ(1)(r) can be expressed exactly as a sum of one-center densities, without the use of any fit procedure. Because ρ(1)(r) is a sum of one-center multipolar densities, the Hartree potential is a function with a simple one-dimensional integral representation, and its matrix elements can be obtained by performing one-dimensional integrals over it. Since δρ(r) is spatially slowly varying, the Hartree potential to which it corresponds and the matrix elements of this potential can accurately be evaluated on a relatively coarse coordinate space mesh, using fast Fourier transforms. The method is illustrated via molecular structure calculations for N2 and NH3. The calculations are accurate to a few percent when the required integrals over δρ(r) and δV(r) are performed on a mesh of spacing 0.4 a.u. The N–N bond length and stretch frequency are found to equal 2.10 a.u. and 2.3×103 cm−1, respectively. The equilibrium N–H bond length and H–N–H angle are calculated to be 1.93 a.u. and 105°, respectively, while the NH3 inversion barrier turns out to equal 0.25 eV. These results are in good agreement with earlier calculations.