Coupling all-electron full-potential density functional theory with grid-based continuum embeddings

TL;DR

This work provides a practical bridge between all-electron full-potential DFT and grid-based continuum embedding by introducing a smoothing scheme that preserves atomic multipole moments when mapping dense nuclear cusps onto regular grids. The modular interface supports self-consistent coupling, combining $G^{\text{solv}}$ with vacuum terms to yield accurate solvent-corrected energies and forces, while a low-pass FFT derivative filter mitigates high-frequency noise. Benchmark tests on water, NaF, and CO on Pt(111) demonstrate converged total energies to within $\sim$1 meV and forces to $\lesssim$10 meV/Å, with clear recommendations for grid densities, multipole orders, and filter parameters depending on the solvent model (SSCS vs SCCS). The approach enables interoperable coupling of Environ with all-electron codes like FHI-aims and paves the way for broader modular multiscale simulations that integrate density embedding and electrostatics across diverse electronic structure packages.

Abstract

Recent advances in continuum embedding models have enabled the incorporation of solvent and electrolyte effects into density functional theory (DFT) simulations of material surfaces, significantly benefiting electrochemistry, catalysis, and other applications. To extend the simulation of diverse systems and properties, the implementation of continuum embedding models into the Environ library adopts a modular programming paradigm, offering a flexible interface for communication with various DFT programs. The speed and scalability of the current implementation rely on a smooth definition of the key physical properties of the atomistic system, in particular of its electronic density. This has hindered the coupling of Environ with all-electron simulation packages, as the sharp electron density peaks near atomic nuclei are difficult to represent on regular grids. In this work, we introduce a novel smoothing scheme that transforms atom-centered electron densities into a regular grid representation while preserving the accuracy of electrostatic calculations. This approach enables a minimal and generic interface, facilitating seamless interoperability between Environ and all-electron DFT programs. We demonstrate this development through the coupling of Environ with the FHI-aims package and present benchmark simulations that validate the proposed method.

Figures (8)

• Figure 1: Communication workflow between FHI-aims (blue/top) and Environ (red/bottom). At each SCF step, the components of $\rho^\text{MP}$ are smoothened in a way that conserves the long range multipole moments of the core regions, cf. \ref{['main_revised:sec:smoothing']}. Shown as inserts are examples of original $l=0$ and $l=1$ components (blue) and their smoothened counterparts (red). The smoothened density is evaluated on Environ's regular grid. Environ defines the cavity and solves the specified electrostatic problem. It computes any additional contributions and corrections specified by the user. To avoid artifacts in $\nabla s$, Environ's FFT core subroutines are modified to include a low-pass filter in the derivatives, cf. \ref{['main_revised:sec:lowpass']}. The solvent contributions to the effective potential $V_\text{eff}^\text{solv}$ in the KS operator, the free energy in solution $G^\text{solv}$, and the forces $\mathbf{F}^\text{solv}_{at}$ (symbols shortened in figure for clarity) are communicated to FHI-aims. The potential is interpolated back onto the overlapping atom-centered grids and added to the vacuum contribution computed by FHI-aims. An additional energy correction is computed, as described in \ref{['main_revised:sec:full_E']}. The KS equations are solved, the density is updated, and the process is repeated to self-consistency.
• Figure 2: Water molecule in implicit solvent. Difference (left) and overlap (right) between $f^\text{switch}_{\text{H}^1,l\neq0}$ of one hydrogen atom H$^1$ and inverted cavity function $1-s$ of entire molecule for SCCS (top) and SSCS (bottom), cut through the molecular plane. $s$ computed from converged $\rho^\text{el}$ for SCCS. Contrast in right images chosen to capture maximum value across both images. Radius $r^\text{solv}_{\text{H}^1}$ of dashed circle estimates shortest distance of H$^1$ from any point outside. Inside of a sphere with half of that radius (hatched area), $\delta \widetilde{\rho}_{\text{H}^1,lm}(r)$ is replaced by a polynomial $P_{\text{H}^1,lm}(r)$. In transition region of $f^\text{switch}_{\text{H}^1,l\neq0}$, original $\delta \widetilde{\rho}_{\text{H}^1,lm}(r)$ is used but integration errors in this region are still partially compensated in inner (hatched) region. For $l=0$, integration errors over entire $\mathbb{R}^3$ are compensated. Corresponding plots for the O atom, different cavity definitions, and ionic systems shown in \ref{['si-sec:separation_smoothing_solvent']} in the SI.
• Figure 3: Gradient component $\partial s / \partial x$ for the same system as in \ref{['main_revised:fig:regions']}, SCCS. Red and blue show positive and negative values, respectively, at an arbitrary scale. $\nabla s$ computed from $\nabla \Bar{\rho}^\text{el}$ via chain rule after SCF convergence. Titles refer to computation of $\nabla\Bar{\rho}^\text{el}$. Top left: FFT with low-pass. Top right: FFT without low-pass. Bottom left: computed in FHI-aims using the analytical derivative of \ref{['main_revised:eq:rho_MP']}. Filtered FFT version correctly reproduces this reference. Bottom right: for reference, FFT without low-pass, the DFT program being Quantum Espressogiannozzi2017 instead of FHI-aims. While the overall less smooth shape of the cavity in FHI-aims is a result of its atom-centered basis set and is consistently present, cf. also \ref{['main_revised:fig:regions']}, the ripples in the unfiltered FFT derivative are an artifact that is absent in the analytical derivatives and can be removed by a low-pass filter. Corresponding plots for $\nabla^2 s$ and differences of filtered to unfiltered and direct gradients shown in \ref{['si-fig:gradients_diff', 'si-fig:laplacians']} in the SI.
• Figure 4: Convergence of total energy with Environ grid density parameter $e_\text{cut}$ for different test systems and solvent models. Dashed lines only as visual guideline.
• Figure 5: Convergence of total forces of selected atoms in different test systems and solvent models with Environ grid density parameter $e_\text{cut}$. Dashed lines only as visual guideline.