Stability of exact force-free electrodynamic solutions and scattering from spacetime curvature

Recently, a family of exact force-free electrodynamic (FFE) solutions was given by Brennan, Gralla and Jacobson, which generalizes earlier solutions by Michel, Menon and Dermer, and other authors. These solutions have been proposed as useful models for describing the outer magnetosphere of conducting stars. As with any exact analytical solution that aspires to describe actual physical systems, it is vitally important that the solution possess the necessary stability. In this paper, we show via fully nonlinear numerical simulations that the aforementioned FFE solutions, despite being highly special in their properties, are nonetheless stable under small perturbations. Through this study, we also introduce a three-dimensional pseudospectral relativistic FFE code that achieves exponential convergence for smooth test cases, as well as two additional well-posed FFE evolution systems in the appendix that have desirable mathematical properties. Furthermore, we provide an explicit analysis that demonstrates how propagation along degenerate principal null directions of the spacetime curvature tensor simplifies scattering, thereby providing an intuitive understanding of why these exact solutions are tractable, i.e. why they are not backscattered by spacetime curvature.

Figure 1

The ratio $|B_{\text{null}}|/|B_{\text{mono}}|$ for the ${\text{null}}^{+}$ solution studied, where ${\mathbf{B}}_{\text{null}}$ and ${\mathbf{B}}_{\text{mono}}$ are the null and monopole contributions in Eq. (2.33), while the norm is defined as $|B_{\text{null}}|=\sqrt{{\mathbf{B}}_{\text{null}}·{\mathbf{B}}_{\text{null}}}$. Shown is a cross section in the $xz$ plane, which extends radially from $R_{-}=1.95\mathrm{M}$ to $R_{+}=195\mathrm{M}$.

Figure 2

The top left figure shows the electric field lines, and the top right figure shows the magnetic field lines for the same ${\text{null}}^{+}$ solution as shown in Fig. 1; the lines are colored by $E^2$ and $B^2$, respectively. The bottom figures plot the time variation of the $E_x$ component, showing the wave propagating inwards.

Figure 3

Half of the computational domain and the spectral grid (spectral collocation points are at the intersection of the black lines); the semitransparent sphere represents the event horizon of the Schwarzschild black hole.

Figure 4

Top: The $L_2$ norm of the error measure $\mathrm{\Delta }{B}^2$ for the constant $\mathbf{B}$ field simulation. Bottom: The $L_2$ norm of the error measure $\mathrm{\Delta }{E}^2$. The simulation time is shown in units of the light-crossing time $R_{+}$ from the origin to the outer boundary $R_{+}$. There is only one curve being plotted in each panel, which oscillates quickly giving the appearance of a band.

Figure 5

Top: The $L_2$ norm of the error measure $\mathrm{\Delta }{B}^2$ for the ${\text{null}}^{+}$ simulations. Bottom: The $L_2$ norm of the error measure $\mathrm{\Delta }{E}^2$. There are ten resolutions being plotted, with the number of radial collocation points given by $k+6$, and the maximum $l$ for the angular $Y^{lm}$ decomposition given by $2k+7$. In both panels, the error decays with increasing $k$, so that $k=0$ corresponds to the topmost (black) line, and the lower lines are, in turn, $k=1,2,3,\dots$. The errors decline approximately linearly in log scale, indicating that our pseudospectral code achieves exponential convergence as expected over a wide range of resolutions, with the highest resolution being limited in accuracy by machine precision.

Figure 6

Top row: Contours of $\Re {\varphi }_0$ that show the structure of the wave train. It has a $Y^{2,2}$ angular profile, and consists of fast radial oscillations with frequency $\mathrm{\Omega }$. It is confined into a radial top-hat-style envelope of width 200 M. The shape of the wave train shows that our test is fully three dimensional, without any axisymmetry. The figures on the left and right depict the wave train at earlier and later times, showing that it is initially moving inwards toward the black hole. The green semitransparent surface is the equatorial slice of the computational domain, included to indicate its extent. The grid structure on this slice is also plotted to provide a visual description of the density of collocation points at $K=7$. Bottom figure: The $\Re {\varphi }_0$ is depicted by warping half of the equatorial slice with a vertical displacement proportional to the $\Re {\varphi }_0$ values at collocation points, and then connecting them via “straight lines” into the gray surface (i.e. the gray surface is a trivial interpolation between the collocation points and is not the actual $\Re {\varphi }_0$ as represented by the spectral functions). The green surface is the top-hat-like radial envelope function. This figure also corresponds to $K=7$.

Figure 7

The log of ${\varphi }_2^{(2,2)}$ plotted as a time series. The curves corresponding to different wave train frequencies $\mathrm{\Omega }$ are shifted slightly relative to each other vertically, so that the curves will not lie right on top of each other, and we can see more clearly. The flat region for the curves correspond to the wave train itself, moving back outwards after having scattered around the black hole. The sloped region corresponds to the exponentially decaying QNM. We can see that the wave train has very different frequencies, but the QNMs excited in all cases share the same frequency and decay rate.

Figure 8

The dashed black lines indicate the fitting results (they extend only over the fitting interval) corresponding to the fitted parameters in Table 1, while the curves being fitted are the negatively sloped segments of the curves in Fig. 7.

Figure 9

The top figure shows the wave train and the QNM tail for the four resolutions of the $\mathrm{\Omega }=0.4$ simulation. The bottom figure shows the differences between resolutions. Because the differences are approximately equally displaced in the vertical direction, which is plotted in log scale, the simulation indeed achieves exponential convergence.

Figure 10

Left: The right-hand side of Eq. (3.1), or explicitly, the initial $-\nabla ·\mathbf{B}$ value on a vertical slice of the computational domain before we apply our initial data solver (the $z$ axis corresponds to $\theta =0$ and $\theta =\pi$). Right: The left-hand side of Eq. (3.1), which is the value of ${\nabla }^2\mathrm{\Phi }$ that emerges after solving for the initial data.

Figure 11

Top: The $L_2$ norm of the difference measure $\mathrm{\Delta }{B}^2$ for the ${\text{null}}^{+}$ simulations initially perturbed in the ${\varphi }_0$ values. Bottom: The $L_2$ norm of the difference measure $\mathrm{\Delta }{E}^2$. The arrangement of the collocation points is the same as for Fig. 5.

Figure 12

Top four panels: The relative difference $|\mathrm{\Delta }\mathbf{B}|/|\mathbf{B}|$ on a vertical slice of the computational domain. The magnitude of the relative difference is indicated by the height of the surface as well as the coloring. Different panels correspond to different times. The initial perturbation seen at $t=0$ propagates inwards, creating the pattern shown at $t=60$. After scattering around the black hole, the perturbations propagate outwards, forming the pattern seen at $t=150$, and eventually begin to exit the computational domain as seen at $t=270$. Bottom two panels: The initial and later perturbations on the horizontal equatorial slice of the computational domain; they show that the perturbations are three dimensional in nature, without axisymmetry.

Figure 13

The numerical error in ${\mathbf{B}}^P$ (as seen in $|\mathrm{\Delta }\mathbf{B}|/|\mathbf{B}|$) being generated at subdomain boundaries (signified by dense concentrations of black lines). The gray frames highlight the locations of the error on such boundaries. Note that the warping scale and the color map are different from Fig. 12, but the plane shown is the same as those appearing in Fig. 12.

Figure 14

Example congruence of geodesics in the equatorial plane of the Schwarzschild black hole, with a nonvanishing impact parameter $b$ translating into $\mathcal{P}=2/27M^2$ (not the value we use for the simulations, but instead chosen to accentuate the perturbation). The red circle is the location of the event horizon. Angle $Y$ is the angle at which the dashed geodesic strikes the origin.

Figure 15

Similar to Fig. 10, but for initial data with a perturbed propagation direction. Left: The right-hand side of Eq. (3.1). Right: The left-hand side of Eq. (3.1).

Figure 16

Top: The $L_2$ norm of the difference measure $\mathrm{\Delta }{B}^2$ for the ${\text{null}}^{+}$ simulations initially perturbed in the propagation direction. Bottom: The $L_2$ norm of the difference measure $\mathrm{\Delta }{E}^2$.

Figure 17

The relative difference $|\mathrm{\Delta }\mathbf{B}|/|\mathbf{B}|$ on a vertical slice of the computational domain. The initial perturbation seen at $t=0$ propagates mostly outwards, creating the patterns seen at $t=60$ and $t=100$, and exits the computational domain as seen at $t=200$.

Figure 18

Similar to Fig. 15, but for initial data with a perturbed propagation direction based on a time-dependent unperturbed solution. Left: The right-hand side of Eq. (3.1). Right: The left-hand side of Eq. (3.1).

Figure 19

Top: The $L_2$ norm of the difference measure $\mathrm{\Delta }{B}^2$ for the initially perturbed time-dependent ${\text{null}}^{+}$ simulation. Bottom: The $L_2$ norm of the difference measure $\mathrm{\Delta }{E}^2$.

Figure 20

The absolute difference $|\mathrm{\Delta }\mathbf{B}|$ on a vertical slice of the computational domain. The initial perturbation seen at $t=0$ propagates inwards creating the patterns seen at $t=70$ and $t=120$. By $t=200$, the perturbation has been almost entirely absorbed by the black hole. We have shown the absolute, rather than the relative difference, because $|\mathbf{B}|$ increases quickly when we approach the black hole, so it is more difficult to see what is going on in the inner regions with a relative difference plot.

Figure 21

The principal directions (black lines) and the asymptotes (red lines) on curved surfaces. The figure on the left depicts the situation with a more generic surface, while the figure on the right depicts a more degenerate surface where two asymptotes become coincident and the principal direction in between them becomes flat.