Accurate Error Estimates and Optimal Parameter Selection in Ewald Summation for Dielectrically Confined Coulomb Systems

TL;DR

This work develops a rigorous error framework for Ewald summation in dielectric-confined quasi-2D Coulomb systems, modeling polarization via an infinite image-charge series and truncating at level $M$. It derives sharp bounds for both image-truncation errors and the electrostatic layer correction (ELC) contributions, and reformulates the ICM-Ewald2D method into an efficient 3D Ewald formulation with padding length $L_z$, enabling $O(N\, ext{log}\,N)$ or $O(N)$ scaling. A key insight is the leading-order analysis of ELC errors, which reveals regimes where truncation can either degrade or improve accuracy, thereby explaining non-monotonic error behavior reported in prior work. The authors propose a practical, case-aware parameter-selection strategy that jointly selects $M$, $L_z$, and Ewald parameters to meet a prescribed tolerance with minimal cost, and validate the approach on prototypical dielectric-confined systems. The results provide actionable guidance for accurate and efficient MD simulations of dielectric-confined Coulomb systems.

Abstract

Dielectrically confined Coulomb systems are widely employed in molecular dynamics (MD) simulations. Despite extensive efforts in developing efficient and accurate algorithms for these systems, rigorous and accurate error estimates, which are crucial for optimal parameter selection for simulations, is still lacking. In this work, we present a rigorous error analysis in Ewald summation for electrostatic interactions in systems with two dielectric planar interfaces, where the polarization contribution is modeled by an infinitely reflected image charge series. Accurate error estimate is provided for the truncation error of image charge series, as well as decay rates of energy and force correction terms, as functions of system parameters such as vacuum layer thickness, dielectric contrasts, and image truncation levels. Extensive numerical tests conducted across several prototypical parameter settings validate our theoretical predictions. Additionally, our analysis elucidates the non-monotonic error convergence behavior observed in previous numerical studies. Finally, we provide an optimal parameter selection strategy derived from our theoretical insights, offering practical guidance for efficient and accurate MD simulations of dielectric-confined systems.

Figures (12)

• Figure 1: Relative errors in force ($\mathcal{E}_r$) as a function of the truncation parameter $M$ for the image charge series. The dashed lines represent the fitted curves with decay rates using our theoretical prediction Eq. \ref{['eq::Ferr1']}. In panel (a), we fix $\gamma_{\text{u}} = \gamma_{\text{d}} = 1$ and consider systems with varying heights $H = 0.5$, $1$, and $5$. In panel (b), we fix $H = 1$ while varying $\gamma_{\text{u}} = \gamma_{\text{d}} = \gamma$ with values of $0.3, 0.6,$ and $1$. In both panels, we fix $L_x=L_y=10$.
• Figure 2: Relative errors in force ($\mathcal{E}_r$) as a function of the padding ratio $P$ for systems without dielectric interfaces. We consider systems with heights $H = 0.5, 1, 5$ while fixing $L_x=L_y=10$. The padding ratio is defined as $P = (L_z - H) / L_x$. The dashed lines represent the fitted curves with decay rates using theoretical prediction Eq. \ref{['eq::ForceError']}.
• Figure 3: Relative errors in force ($\mathcal{E}_r$) for systems with dielectric interfaces. Here we fix $\gamma_{\text{u}}=\gamma_{\text{d}}=\gamma = 0.6$, $L_x=L_y=10$ and $H = 0.5$. Panel (a) illustrates errors as a function of padding ratio $P$ with fixed image charge layers ($M=5,~15,~25$); panel (b) illustrates errors as a function of $M$ with fixed padding ratios ($P = 1,~3,~5$). The dashed lines in (a) and (b) represent the fitted curves with decay rates using theoretical predictions Eq. \ref{['eq::lhs_less1']} and Eq. \ref{['eq::Ferr1']}, respectively.
• Figure 4: Relative errors in force ($\mathcal{E}_r$) for systems with dielectric interfaces. Here we fix $\gamma_{\text{u}}=\gamma_{\text{d}}=\gamma = 1$, $L_x=L_y=10$ and $H = 0.5$. Panel (a) illustrates errors as a function of padding ratio $P$ with fixed image charge layers ($M=25,~35,~45$); panel (b) illustrates errors as a function of $M$ with fixed padding ratios ($P = 1,~3,~5$). The dashed lines in (a) and (b) represent the fitted curves with decay/growth rates using theoretical predictions Eq. \ref{['eq::lhs_bigger1']} and Eq. \ref{['eq::Ferr1']}.
• Figure 5: Relative errors ($\mathcal{E}_r$) in electrostatic energy for systems of 2:1 electrolytes with dieletric interfaces. Here we consider the same system setup studied by Yuan et al.yuan2021particle using the ICM-PPPM method with $\gamma_{\text{u}}=\gamma_{\text{d}}=\gamma$ = 0.6, 0.95, and 1 (from left to right), respectively. In each panel, we fix $L_x = L_y = 15$, $H = 5$, and consider $L_z=$$45$, $75$, and $105$. Finally, the dashed lines represent the fitted curves using the sum of Eqs. \ref{['eq:Uerr']} and \ref{['eq::U_ELC']} (with coefficients determined by fitting).