Grid operators

anterpolate!(grid_values, positions, charges, grid, basis) -> grid_values

Anterpolation (paper eq. 7): scatter particle charges to the grid using basis function Φ evaluated at offsets (r_i − x_m) / h. For each particle, contributions are added to all grid points within the basis support window.

Periodic axes wrap; open axes silently drop out-of-range contributions.

grid_values is cleared before scattering. Operates in-place; returns grid_values.

interpolate!(potentials, positions, grid_values, grid, basis) -> potentials

Interpolation (paper eq. 12): gather grid potentials onto particle sites using basis function Φ. The transpose of anterpolate!.

potentials is overwritten with the gathered values. Returns potentials.

interpolate_grad!(grads, positions, grid_values, grid, basis) -> grads

Gradient of interpolate with respect to particle position: returns the field gradient at each particle site, ∇e_i = Σ_m ∇Φ((r_i−x_m)/h)·e_m / h.

Used by force evaluation (∇ on the dual basis side).

restrict!(dst, src, dst_grid, src_grid, basis) -> dst

Grid-to-grid restriction (paper eq. 8). For each destination grid point, gather contributions from the surrounding source grid points, weighted by Φ((r_src − r_dst)/h_dst).

The canonical MSM restriction is the special case where dst_grid has exactly twice the spacing of src_grid (a "coarsen by 2" step); the implementation is generic in the two grids' spacings and origins.

This is the gather form of the operator — each output cell dst[n_dst] is computed independently from src. Threaded over destination points via OhMyThreads. Periodic axes wrap via wrap_index; open axes drop out-of-range contributions.

dst is overwritten.

prolong!(dst, src, dst_grid, src_grid, basis) -> dst

Grid-to-grid prolongation (paper eq. 11). For each destination grid point, gather from source grid points within the basis support window, weighted by Φ((r_dst − x_src)/h_src).

The transpose of restrict!. Canonical MSM prolongation is the special case where src_grid has exactly twice the spacing of dst_grid.

dst is overwritten (not accumulated). Callers wanting to add to an existing field should do dst .+= prolong(...) themselves.

grid_cutoff!(e, q, grid, splitting, l) -> e

Compute the level-l "grid cutoff" potential (paper eq. 9):

e[m] = Σ_n k_l(r_m − r_n) q[n]

where the sum is over destination grid points n within the compact support of k_l. Periodic axes wrap; open axes drop out-of-bounds contributions.

e is overwritten.

build_stencil(splitting, l, h_grid) -> (stencil, smax)

Precompute the level-l grid-cutoff stencil for the given splitting at grid spacing h_grid. Returns a (2·smax+1)-shaped tensor of stencil values together with the per-axis half-width smax.

The level-l "long-range" kernel k_l(r) (1 ≤ l ≤ L-1) has compact support |r| ≤ a_{l+1} = 2^l · a. The stencil radius is smax_α = ceil(a_{l+1} / h_grid_α).

top_level!(e, q, grid, splitting) -> e

Top-level grid potential (paper eq. 10):

e[m] = Σ_n k_L(r_m − r_n) q[n]

evaluated as a direct sum over all grid pairs (no cutoff — k_L has no compact support). Periodic axes apply minimum-image displacement.

When the grid is fully periodic the top grid is normally reduced to a single point along each periodic axis (see top_grid_size); the charge at that point is then set to zero by [apply_neutralising_background!] if the splitting requires it.

e is overwritten.

apply_neutralising_background!(q, splitting, grid) -> q

For splittings flagged with requires_neutralising_background == true and a fully periodic grid, zero out the (single-point) top-level grid charge to implement the neutralising-background trick (paper §2.1).

For partial periodicity or splittings that don't require it, this is a no-op. Returns q for chaining.

top_grid_size(g_finest, L) -> NTuple{D,Int}

The MSM top-level grid extent, derived from the finest grid by halving each axis L−1 times. Periodic axes are clipped to a minimum of 1; open axes are simply halved.