Pseudospectral implementation of the Einstein-Maxwell system

TL;DR

This work delivers a first construction of a pseudospectral, first-order symmetric-hyperbolic Maxwell evolution embedded in the Einstein equations within the BAMPS code. By introducing constraint-damping via auxiliary fields and boundary-preserving conditions, it achieves boundary stability for the electrovacuum system and enables robust, high-accuracy simulations. The authors validate the approach through convergence tests, Reissner–Nordström and Kerr quasinormal mode benchmarks, and a strong-field, near-threshold collapse–like scenario, demonstrating reliable EM–GR coupling in electrovacuum. The methodology and verification suite establish a solid foundation for exploring electrovacuum dynamics, critical phenomena, and EM–gravitational wave interactions in strong gravity environments.

Abstract

Electromagnetism plays an important role in a variety of applications in gravity that we wish to investigate. To that end, in this work, we present an implementation of the Maxwell equations within the adaptive-mesh pseudospectral numerical relativity code BAMPS. We perform a thorough analysis of the evolution equations as a first order symmetric hyperbolic system of PDEs. This includes both the construction of the characteristic variables for use in our penalty boundary communication scheme, as well as radiation controlling, constraint preserving outer boundary conditions which, for the first time in a numerical context, are shown to be boundary-stable. After choosing a formulation of the Maxwell constraints that we may solve for initial data, we move on to show a suite of numerical tests. Our simulations, both within the Cowling approximation, and in full non-linear evolution, demonstrate rapid convergence of error with resolution, as well as consistency with known quasinormal decay rates on the Kerr background. Finally we evolve the electrovacuum equations of motion with strong data, a good representation of typical critical collapse runs.

Figures (10)

• Figure 1: In bamps, the domain is divided into three subdomains of different geometry. Each one of them is divided in patches. Inside each patch, the solution is represented as a linear combination of Chebyshev polynomials. The number of patches in the inner region is given to match the angular resolution outside.
• Figure 2: Snapshots of $\sigma E^\phi$ in the $xz$ plane at different times. We can see the expected behavior of propagation of a dipole wave. This is only a subset of the simulation domain; the true boundary is placed at double the radius of the plot.
• Figure 3: Convergence of various error norms as function of time $t$ with the number of collocation points $n$ for a simulation of an electromagnetic dipole wave on a flat background. Left: Norm of the difference to the analytic solution to the Cauchy problem. Observe that here it is expected, and observed, that once the outer boundary has a significant effect on the numerical solution we can no longer see convergence to a solution of the Cauchy problem. Middle: Maxwell constraint norm. We observe here that the constraints continue to converge even beyond $t/\sigma\simeq20$, in line with our interpretation of the left panel. Right: Norm of the difference between subsequent resolutions. The norms of the differences are calculated in the $x$ axis, as in equation \ref{['eq:norm_x_axis']}, whereas the constraint monitor is calculated with a norm over the whole domain as defined in Eq. \ref{['eq:maxwell_constraint_monitor']}.
• Figure 4: Comparison of the constraint monitor for different number of patches, scaled by $N^2$ in all sub-domains. The result is consistent with high-order polynomial convergence, and the most accurate run shows less error than the ones in Fig. \ref{['fig:dipole_results']}.
• Figure 5: Constraint violation at time $t/\sigma = 25$, for both convergence tests, increasing the number of collocation points (p refinement) and increasing the number of patches (h refinement). Notice that the bottom axis scale is linear in the number of collocation points $n$ (left) and logarithmic in the number of patches $N$ (right).