Skip to content
#edge computing Preprint

A Fully Matrix-Free Three-Grid Preconditioner for the Time-Harmonic Maxwell Equations at Extreme Scale

Aug 2026 · 0 citations · 30 references
Mathematics Computer Science

Abstract

Three-dimensional time-harmonic Maxwell simulations generate massive complex indefinite systems whose mesh coarsening is strictly limited by phase accuracy. Although matrix-free finite element kernels utilize GPU throughput efficiently, standard multilevel solvers are ultimately bottlenecked by the memory and communication costs of exact coarse-grid factorizations. We present a fully matrix-free, factorization-free three-grid preconditioner for curl-conforming N{\'e}delec discretizations with perfectly matched layers (PML) and optimally blended quadrature. The method employs an outer FGMRES to solve the unshifted fine-grid equation, while an intermediate-grid correction is computed by a fixed-work FGMRES preconditioned with a complex-shifted $2h$--$4h$ cycle. This strategically confines the complex shift to an auxiliary preconditioner, preserving the physical Maxwell operator. A local Fourier analysis derives the blended Maxwell branches and compatible edge transfers, identifying robust shift and Jacobi damping parameters. Validated against the analytical Maxwell Green tensor, our approach demonstrates extreme scalability: using a single solver configuration, both homogeneous and highly heterogeneous systems with approximately 10.89 billion complex edge unknowns are solved in 42.0--72.0 seconds on just 64 NVIDIA A100 GPUs.

View source

Similar papers

Preprint Aug 2026

A Massively Parallel Three-Grid Preconditioner for the High-Frequency Helmholtz Equation

Accurate simulation of three-dimensional time-harmonic wave propagation over many wavelengths requires control of phase error and efficient solution of large indefinite systems. We develop a three-grid solver based on the compact 27-point interpolated optimized finite-difference (IOFD) discretization. Its wavenumber-dependent stencil supports a fine-grid resolution of six points per shortest wavelength and an unshifted physical correction on the \(2h\) grid at only three points per shortest wavelength. The method retains unshifted IOFD operators on the \(h\) and \(2h\) grids, while a complex-shifted \(2h\)--\(4h\) auxiliary cycle preconditions a factorization-free iterative approximation of the coarse inverse. Restricting the shift to this auxiliary cycle preserves the propagative character of the coarse correction. Comparison with the outgoing Green function confirms phase and relative-amplitude accuracy on a sequence of meshes up to \(6144^3\), with the largest problem spanning approximately 1024 wavelengths per coordinate. The same fixed solver configuration retains robust convergence across smooth, discontinuous, high-contrast, and geophysical velocity models and exhibits scalable parallel performance. In particular, a problem spanning approximately 340 wavelengths in each coordinate direction is solved in 18.1 seconds on just 64 NVIDIA A100 GPUs.

Shubin Fu, Yitong Wang, Zixiao Zhao · 0 citations
Preprint Jul 2026

R3MG-C: a high-order algebraic-geometric multilevel preconditioner for continuous finite element discretizations

Algebraic multigrid (AMG) methods are robust and efficient black-box preconditioners for linear systems arising from low-order discretizations of elliptic partial differential equations, but their performance often deteriorates for high-order methods. We introduce an algebraic-geometric multilevel preconditioner for continuous lagrangian finite element discretizations. The method automatically constructs prolongation operators and a Galerkin hierarchy using only the coordinates of finite element support points, without requiring a prescribed mesh hierarchy or domain decomposition. The support-point are recursively partitioned by an R-tree algorithm based on axis-aligned bounding boxes, producing a hierarchy of agglomerates. In contrast to classical AMG, which typically uses piecewise-constant aggregate coarse spaces, our method embeds local discontinuous polynomial spaces of degree $p'>0$, defined on the agglomerate boxes, through continuous nodal interpolation. This yields a high-order extension of the Nicolaides coarse space. A two-level analysis quantifies how the coarse polynomial degree offsets the effects of large agglomerates and high-order fine discretizations, under uniform box-regularity, interpolation-stability, and stable-decomposition assumptions. Numerical experiments in two and three dimensions show that the resulting V-cycle, used as a conjugate-gradient preconditioner, maintains stable iteration counts and remains effective in regimes where standard AMG deteriorates.

Davide Polverino, Marco Feder, Luca Heltai · 0 citations
Preprint Jul 2026

An unfitted boundary algebraic equation method with Calder\'on preconditioning for 2D Stokes flow in irregular geometry

We present an unfitted boundary algebraic equation method for the two-dimensional exterior/interior Stokes equations on a staggered MAC grid. By constructing an explicit free-space pair of velocity and pressure lattice Green's functions (LGFs) from free-space Laplace LGFs, we represent homogeneous fields using sources supported exclusively on thin staggered boundary layers. This formulation imposes physical Dirichlet data at cut points via local interpolation, while sampled-normal rank updates remove hydrostatic null modes associated with single or multiple obstacles. The workflow parallels that of classical boundary integral formulations and requires no artificial boundary conditions for exterior flows, but follows a discretize-then-represent route and does not require singular/near-singular quadrature. The resulting dense boundary system is solved via GMRES, utilizing a componentwise discrete Calder\'on preconditioner built from the scalar Laplace kernel and padded FFTs for fast volume convolutions. Extensive numerical validation, including multiply connected domains, narrow gaps, and Moffatt eddies, confirms discrete incompressibility to solver accuracy and recovers the expected Moffatt eddy scaling. We achieve second-order velocity and pressure convergence and bound maximum discrete divergence within numerical accuracy. The discrete Calder\'on preconditioner reduces the condition number by orders of magnitude and yields nearly mesh-independent conditioning in exterior configurations, while remaining effective---though more demanding---for narrow-gap and fine-grid interior problems.

Wenjun Ying, Qing Xia · 0 citations
Preprint Aug 2026

Analysis of Block Jacobi/Gauss-Seidel and additive/multiplicative Schwarz preconditioning through the theory of GLT sequences, with applications to domain decomposition discretizations

When a linear differential problem is discretized by a linear numerical method characterized by a mesh fineness parameter $n$, the computation of the numerical solution reduces to solving a linear discrete problem identified by a matrix $A_n$ whose size grows with $n$. The sequence of discretization matrices $\{A_n\}_n$ often falls within the class of generalized locally Toeplitz (GLT) sequences, even when the numerical method belongs to the family of domain decomposition methods (DDMs), as illustrated herein through examples. Four widely used preconditioners for DDM discretization matrices are the block Jacobi (BJ), block Gauss--Seidel (BGS), additive Schwarz (AS), and multiplicative Schwarz (MS) preconditioners. In this paper, we provide formal definitions of the BJ/BGS/AS/MS preconditioners for arbitrary multilevel block matrices. These definitions and the associated notations are inspired by the theory of GLT sequences and are proposed as alternatives to those commonly used by the DDM community. We analyze the structure of the BJ/BGS/AS/MS preconditioners when applied to multilevel block matrices $A_n$ belonging to a GLT sequence $\{A_n\}_n$. Every GLT sequence $\{A_n\}_n$ is uniquely associated with a special function $\kappa$ called symbol. We prove that, if $\{A_n\}_n$ is a GLT sequence with symbol $\kappa$, then the sequences of the BJ, BGS, and MS preconditioners are GLT sequences with symbol $\kappa$. For the AS preconditioner, we prove that $\{P_n^{AS}(A_n)\}_n$ is a GLT sequence with symbol $\kappa^{AS}\approx\kappa$, and $\kappa^{AS}=\kappa$ whenever the overlaps in the subdomains used for the construction of $P_n^{AS}(A_n)$ vanish as $n\to\infty$. A numerical validation of these results in the context of isogeometric DDMs is presented.

C. Garoni, Abdessadek Rifqui, S. Serra-Capizzano · 0 citations
Preprint Aug 2026

Fast high-order solvers for the Lippmann--Schwinger equation in piecewise-smooth heterogeneous media

This article presents a fast, high-order Nystr\"om solver for the two-dimensional Lippmann--Schwinger equation arising from time-harmonic scattering by penetrable, piecewise-smooth heterogeneous media. Relying on high-order evaluation of the Newtonian potential on unstructured grids adapted to interfaces of discontinuity, the methodology achieves high-order accuracy using existing fast algorithms such as the fast multipole method. As an iterative method the solver exhibits rapid convergence when coupled to a preconditioning strategy that exploits a class of structured-grid solvers---fast solution methods offering quasi-linear time and memory complexity and nearly-constant iteration counts, but long limited in accuracy. The preconditioning strategy couples the Nystr\"om discretization---given by high-order quadrature nodes over a (curved) unstructured mesh conforming to the support of the spatially varying contrast---to a uniform Cartesian grid underlying the fast preconditioner via a pair of transfer operators. The resulting preconditioner inherits the frequency-robust behavior of its Cartesian counterpart without sacrificing the geometric flexibility and high-order accuracy of the unstructured discretization. We prove that invertibility of the proposed preconditioner holds under explicit conditions on the mesh sizes and on the Cartesian preconditioner. Numerical experiments demonstrate that the preconditioned system requires significantly fewer GMRES iterations than its unpreconditioned counterpart, with iteration counts almost independent of mesh size and wavenumber, and illustrate the method's robustness for inhomogeneities with piecewise-smooth refractive indices and jump discontinuities across interfaces.

T. Anderson, Juan Burbano-Gallegos, Luiz M. Faria et al. · 0 citations
Open access Jul 2026

MeshFEM: A Block-accelerated Solver for Nonlinear Finite Elements

We introduce a high-performance framework for solving nonconvex optimization problems arising in simulation and geometry processing applications. Our framework especially benefits high-dimensional problems whose unknowns are spatial coordinates in ℝd and whose objectives are expressed as a sum of element energies over small stencils. The Hessians of these problems have a block-sparse structure, where the nonzero entries are grouped into dense d × d blocks. Many works have exploited this structure to speed up the Hessian-vector products and nodal smoothers of iterative solvers. In contrast, accelerating direct sparse linear solvers, preferred for applications requiring high accuracy on unstructured grids, has been less explored. Our framework leverages block structure to accelerate all phases of a Newton-type minimization algorithm, from sparsity pattern construction to system assembly, matrix factorization, and solves. A fundamental piece of the framework is BlockCatamari, our highly tuned, block-accelerated multifrontal sparse Cholesky code that achieves significantly faster factorizations than existing direct solvers on shared-memory systems. We combine this solver with our fast parallel Hessian assembly routine, heuristics for conditionally enabling per-element Hessian projections, and strategies for minimizing symbolic factorization re-computations. We show that this leads to dramatic speedups on a large suite of benchmark problems involving injective surface parametrization and elasticity simulation with contact, and we further demonstrate the ease with which new energies can be defined at the different levels of abstraction offered by our framework.

Haleh Mohammadian, Xinzhuo Hu, Roi Poranne et al. · 1 citation

Related blog posts

Microsoft Research Blog Aug 31, 2026

GigaPath-Flash and GigaTIME-Flash: Toward population-scale discovery with efficient pathology foundation models

What if pathology foundation models could do more with less? GigaPath-Flash and GigaTIME-Flash cut computational demands while maintaining strong performance, opening the door to larger studies and broader exploration. The post GigaPath-Flash and GigaTIME-Flash: Toward population-scale discovery with efficient pathology foundation models appeared first on Microsoft Research.

MIT News · Artificial Intelligence Aug 27, 2026

Looking beyond natural sequences

A new machine-learning framework aims to improve the success rate of computational protein design while moving away from results that reproduce sequences found in nature.