Fast Ewald summation efficiently evaluates Coulomb interactions and is widely used in molecular dynamics simulations. It is based on a split into a short-range and a long-range part, where evaluation of the latter is accelerated using the fast Fourier transform (FFT). The accuracy and computational cost depend critically on the mollifier in the Ewald split and the window function used in the spreading and interpolation steps that enable the use of the FFT. The first prolate spheroidal wavefunction (PSWF) has optimal concentration in real and Fourier space simultaneously, and is used when defining both a mollifier and a window function. We provide a complete description of the method and derive rigorous error estimates. In addition, we obtain closed-form approximations of the Fourier truncation and aliasing errors, yielding explicit parameter choices for the achieved error to closely match the prescribed tolerance. Numerical experiments confirm the analysis: PSWF-based Ewald summation achieves a given accuracy with significantly fewer Fourier modes and smaller window supports than Gaussian- and B-spline-based approaches, providing a superior alternative to existing Ewald methods for particle simulations.
Long-range interactions such as Coulomb and Stokes potentials arising in particle-based systems are central to molecular dynamics, fluid dynamics, and wave propagation. A major computational challenge is that direct evaluation of such interactions scales quadratically with the number of particles. Ewald summation, introduced in 1921 [12] reduces this cost and has since become a standard tool in large-scale simulations.
In Ewald’s original formulation, the three-dimensional Laplace kernel G(r) = 1/r (r = |x|, x ∈ R 3 ) is decomposed into a mollified kernel M (r) and a residual kernel R(r) according to
where erf(x) = 2 √ π
x 0 e -t 2 dt and erfc(x) = 1 -erf(x) are the error and complementary error functions. The parameter σ > 0 is a width-parameter of the mollification that controls the decay of both terms. In the original Ewald-split (1), the mollified kernel M is the convolution of G with a Gaussian 1 √ πσ e -r 2 /σ 2 , while the residual is R = G -M . Because a Gaussian is not compactly supported, R must be truncated numerically outside a cutoff radius r c > 0, with truncation error controlled by a numerical tolerance that is typically denoted by ε.
The Ewald split enables efficient evaluation of the potential ϕ : R 3 → R at the particle locations {x j } n j=1 in a triply periodic domain Ω = [0, L) 3 (cubic for clarity of presentation; the method extends directly to general periodic lattices). Using (1), it can be written as
where {ρ j } n j=1 are particle strengths (e.g., charges), r ∈ Z 3 indexes the periodic images, and the prime indicates that the self-interaction (i = j, r = 0) is omitted. The Ewald split (1) makes the otherwise conditionally convergent sum (2) convergent provided that the charge-neutrality condition n j=1 ρ j = 0 holds. In (3), the mollified part M is evaluated in Fourier space, where M (k) decays rapidly, while the residual R captures local interactions within a cutoff radius r c . This reduces the cost of evaluating ϕ from O(n 2 ) to O(n 3/2 ) for close to randomly uniform particle distributions [15].
Several fast Ewald variants accelerate the Fourier-space sum further to O(n log n) using the fast Fourier transform (FFT) [7]. These include the Particle Mesh Ewald (PME) [8], Smooth Particle Mesh Ewald (SPME) [11], Particle-Particle-Particle-Mesh Ewald (PPPM or P 3 M) [15], the Spectral Ewald (SE) [2,20,28], and the Particle-Particle NFFT Ewald (P 2 NFFT) method [16,22]. These schemes spread particle data onto a uniform grid using a localized window function, apply the FFT, scale in Fourier space, and interpolate back. The accuracy and efficiency therefore depend on two distinct choices of functions: (i) the mollifier used in the Ewald split, and (ii) the window function used in spreading and interpolation. Both influence the spectral decay and thus the number of Fourier modes and grid points required for a target accuracy. Classical methods differ mainly in these choices: PME, SPME, and P 3 M employ low-order Bsplines as window functions, while SE and P 2 NFFT use Gaussians or various types of approximations to the first prolate spheroidal wavefunction.
The first prolate spheroidal wavefunction of order zero (PSWF), studied extensively by Slepian and collaborators [30,31] and in [23,24], is optimally concentrated in both real and Fourier space. This makes it attractive both as a mollifier (for splitting) and as a window function (for spreading/interpolation), though in principle different PSWFs with different parameters may be used for the two roles. Historically, the absence of a closed form limited practical use, and approximations such as the Kaiser-Bessel and exponential of semicircle functions were employed instead [5,18], including in SE and P 2 NFFT implementations [22,28]. With modern evaluation algorithms [23] and polynomial approximations [5], PSWFs have become computationally practical. and Fourier space (right). For equal effective support, the PSWF mollifier requires roughly half as many Fourier modes to achieve the same accuracy.
Until recently, PSWFs have primarily been considered as window functions in fast Ewald summation. However, employing them also as mollifiers in the Ewald split yields further advantages, reducing the number of required Fourier modes by around a factor of two in each dimension. Figure 1 illustrates this effect by comparing the decay of the PSWF and Gaussian functions.
The idea of using PSWFs for kernel splitting was first introduced in the Dual-split Multilevel Kernel (DMK) framework [17]. The PSWF split and the DMK algorithm have since also been applied to Stokes potentials [1]. While the DMK framework established the basic concept, a complete mathematical description of PSWF-based Ewald methods is still lacking.
In this paper, we provide such a description, incorporating the first PSWF of order zero both as a mollifier and a window function, with distinct parameters for the two roles. Our main contributions are:
• Explicit definiti
This content is AI-processed based on open access ArXiv data.