Magnetic reconnection on a Klein bottle Open Access

During reconnection, oppositely directed components of a magnetic field come together at an X-point and undergo a topological reordering that enables the transfer of magnetic energy to the surrounding plasma. In the collisionless plasmas of space and astrophysics, reconnection has particular importance due to its role in catalyzing energy release over a wide range of scales in flares, jets, and other phenomena.

The form of the magnetic field in the Harris sheet imposes constraints on the boundary conditions. Periodic boundaries are typically used in the horizontal (x) direction, but the anti-symmetric nature of $B$⁠, $Bx(y=−Ly/2)=−Bx(y=Ly/2)$⁠, precludes the use of periodic boundaries in the vertical (y) direction. Two options are usually considered. The first is the inclusion of hard conducting walls at the top and bottom boundaries of the domain. Appropriate boundary conditions are imposed on the electromagnetic fields (e.g., Dirichlet for $Ex$⁠, $By$⁠, and $Ez$ and von Neumann for the other components) while particles specularly reflect from the walls. The influential GEM Challenge² used these conditions and so hereinafter we will use its name to denote them. For the GEM Challenge boundary conditions, the computational domain is topologically equivalent (homeomorphic) to a cylinder.

Either choice comes with limitations. Simulations using the GEM Challenge boundary conditions will stagnate when the reconnected flux accumulating in the magnetic islands presses against the hard boundaries. Fully periodic boundary conditions employing the field described in Eq. (2) can mitigate this stagnation by offsetting the X-points on the two current layers by $Lx/2$⁠, thus allowing the flux accumulating in the island of one layer to drive inflow toward the X-line on the other. However, keeping the amount of flux available to reconnect the same for both types of boundary conditions (equivalently, keeping the aspect ratios of the current sheets the same) requires a doubling of the size of the computational domain, and hence the computational cost, in the fully periodic case. Less-common alternatives, e.g., open boundaries or driven reconnection,^3–6 increase the computational complexity and can only partially mitigate these issues (often while raising new ones).

Here we present a novel type of periodic boundary condition in which the computational domain is homeomorphic to a Klein bottle. Making this choice allows the simulation to enjoy the benefits of the two major options discussed above: the computational expense of the GEM Challenge conditions coupled with the self-driving of the fully periodic domain. We show that these boundary conditions preserve the important properties of magnetic reconnection. Section II discusses our computational methods, including the implementation of the new boundary condition. Results from the simulations are presented in Sec. III, while Sec. IV offers our conclusions.

The simulations are performed with p3d, a massively parallel particle-in-cell (PIC) code with a long history of simulating magnetic reconnection.⁷ It employs units in which a reference density $n0$ and magnetic field strength $B0$ define the units of length, the proton inertial length $dp0=c/ωp0$ (where $ωp0=4πn0e2/mp$ is the proton plasma frequency), and time, the proton cyclotron time $Ωp0−1=mpc/eB0$⁠. Velocities are measured in units of a reference Alfvén speed $vA0=B02/4πmpn0$⁠, while electric fields and temperatures are normalized to $vA0B0/c$ and $mpvA02$⁠, respectively. The simulations presented here follow particles in a 2x3v (also known as 2.5D) phase space in which no variations are allowed in the third spatial dimension (i.e., $∂/∂z=0$⁠).

With the exception of a slight alteration in the size of the computational domain, the parameters are chosen to match those of the GEM Challenge.² In order to reduce the separation between the spatial scales of the protons and the electrons, and hence the computational expense, the electron mass is taken to be $me/mi=0.04$⁠, which implies the electron inertial length $de=0.2dp$⁠. For similar computational reasons, the ratio of the speed of light to the Alfvén speed (equivalently, the ratio of the proton plasma and cyclotron frequencies) is taken to be 20. The initial electron and proton temperatures are $1/12$ and $5/12$⁠, respectively. The timestep is chosen to satisfy the CFL condition and, as a consequence, particles can travel at most one grid cell per timestep.

For later reference, it is useful to discuss the layout of the computational domain in more detail than usual. To parallelize a computation, p3d breaks the domain into subdomains in configuration space, each of which is assigned to a single processor. The processors are organized into a rectangular grid of size $P=Px×Py×Pz$⁠, and each processor is further subdivided into $Nx×Ny×Nz$ grid points. For notational simplicity, we will assume $Nz=Pz=1$ in what follows, but the extension to multiple processors and grid points in the third dimension is straightforward. Each processor has an overall label that can range from 0 to $P−1$ as well as a label in each dimension: $px∈{0,1,2,…,Px−1}$⁠, $py∈{0,1,2,…,Py−1}$⁠, and $pz∈{0,1,2,…,Pz−1}$⁠. The overall label and the triplet $(px,py,pz)$ are unique to each processor and it is straightforward, given the dimensions of the computational domain, to convert between the two. For example, the processor at the bottom-left of the domain, $P=0$⁠, has triplet $(0,0,0)$⁠. The processor to its immediate right, $P=1$⁠, has triplet $(1,0,0)$⁠, and the processor immediately above is $P=Px$ with triplet $(0,1,0)$⁠. For communication purpose, every processor needs to know the identities of the other processors bordering it, and while this calculation is simple for those on the interior of the computational domain, it can be more complicated for those on the boundaries.

The novel boundary condition presented here is based on the topology of a Klein bottle, a two-dimensional non-orientable manifold. (The perhaps-more-familiar Möbius strip is another example of a non-orientable manifold.) A Klein bottle is a surface with no distinct front or back side that can be obtained by identifying and gluing together specific edges of a rectangle in a non-trivial way, resulting in a structure that cannot be embedded smoothly in three-dimensional space without self-intersections or singularities (see Fig. 1 for a representation). Higher-dimensional analogs of Klein bottles can be similarly constructed. For example, mapping the edges of a rectangular prism in a particular way produces a system homeomorphic to a Klein geometry in $ℝ4$⁠.

In numerical simulations, the identification and gluing-together process is equivalent to describing how processors abut one another when crossing a boundary. This process is given a pictorial representation in Fig. 2. In this visualization, small squares represent individual processors; squares painted the same color neighbor one another in the vertical direction. (Only the top and bottom boundaries are of interest, as the horizontal direction is assumed to be periodic in each panel.) In the configuration shown, with $P=Px×Py×Pz=16×8×1=128$⁠, the processors are numbered from 0 to 127, with those on the bottom row having $P∈{0,1,…,15}$⁠. The fact that the doubly periodic boundary conditions described in Sec. I are topologically equivalent to a torus are demonstrated in the first panel. Traveling downward from the bottom row brings one to the top row, where the processors have $P∈{112,113,…,127}$⁠. In other words, processors next to the boundary (on either side) with triplets of general form $(i,j,0)$ border processors with triplets of the form $[i,Py−(j+1),0]$⁠. The overplotted arrow illustrates how a vector transforms when moving between the top and bottom rows of processors. Note that interpreting the arrow as a magnetic field vector makes it clear why a system with a single current sheet—across which the field must switch direction—cannot use these boundary conditions.

To map a two-dimensional sheet to a Klein bottle, one connects the left and right sides together as usual. However, processors on the top and bottom boundaries now adjoin counterparts at mirrored positions on the x-axis, as shown in Fig. 2(ii). A processor on the top or bottom boundary with triplet of the form $(i,j,0)$ borders a processor with triplet $[Px−(i+1),Py−(j+1),0]$⁠. (Note the set of neighboring grid points within processors are also mirrored, so that the leftmost grid point on, for example, the bottom row of processor 1 is considered to be directly above the rightmost grid point in the top row of processor 126.) As a consequence, vectors flip orientation across the boundary, as shown by the overplotted arrow. Implementation of this boundary condition requires several algorithmic changes. Processors on the top and bottom rows have a new set of neighboring processors, and the x and z components of any vector quantity, e.g., particle velocities and electromagnetic fields, must be flipped across the boundary. (In practice, processors in p3d interact with neighbors through guard cells. Any changes to the signs of gridded quantities required by the boundary conditions occur when these cells are updated.) This means that only one current sheet needs to be included in the domain, as the reversal in $Bx$⁠, is naturally accommodated.

However, the basic Klein boundary condition shown in Fig. 2(ii) is not ideal for reconnection simulations. The reason becomes clear after consideration of the magnetic configuration at late times. Reconnection proceeding at the X-line (assumed to be at the center of the domain) produces an island of magnetic flux centered on the left/right (they are equivalent) boundary. This island grows, eventually becoming large enough to approach the top and bottom boundaries. In runs employing the GEM Challenge boundary condition, this island presses against the conducting walls, producing a back-pressure that resists further reconnection and eventually stifles the simulation. The same issue bedevils the Klein boundary condition, although in that case the back-pressure is not provided by conducting walls but by the island itself—growth of the island in, for instance, the bottom-left of the domain presses against the island in the upper-right.

Fortunately, a straightforward adjustment that mimics the behavior of doubly periodic boundary conditions can alleviate this issue by, in essence, creating self-driving reconnection. Arranging for the island growing on the left and right sides of the domain to emerge from the top and bottom boundaries in the middle of the domain extends the time that reconnection will occur. The phase shift is produced by displacing the set of neighboring processors for either the top or bottom boundary by $Lx/2$⁠, forming the configuration displayed in Fig. 2(iii). Returning to the example of a $16×8$-processor domain, each processor with $P∈{0,1,…,7}$ now neighbors the corresponding processor $P∈{119,118,…,112}$⁠, and vice versa, while each processor $P∈{8,9,…,15}$ abuts the corresponding processor $P∈{127,126,…,120}$⁠. More generally, a processor on the border with triplet $(i,j,0)$ directly borders a processor with triplet $([Px/2−(i+1)] modPx,Py−(j+1),0)$⁠. As with the unshifted Klein boundaries, the set of neighboring grid points are again mirrored across the boundary so, for example, the leftmost grid point in the bottom row of processor 1 is considered to be directly above the rightmost grid point of processor 118.

We now present results from simulations using the various types of boundary conditions discussed above. The computational domain measures $Lx×Ly=51.2di×12.8di$ for the GEM Challenge and shifted Klein cases and is doubled in $Ly$ for the doubly periodic case. (In other words, each reconnecting current sheet has the same aspect ratio in all runs.) The reconnecting component of the magnetic field has the initial form shown in either Eq. (1) or (2), while the density has the asymptotic value of $0.2n0$ and varies as necessary to ensure force balance in the y direction. Section III A considers the case where the asymptotic fields are anti-parallel, as in the original GEM Challenge, while Sec. III B adds a $Bz$ component that reduces the shear angle between the asymptotic fields.

Figure 3 shows a comparison between runs with the GEM Challenge and the shifted Klein boundary conditions with the three panels in each column depicting the simulations at $tΩcp=20$⁠, 40, and 60, respectively. In each, the out-of-plane current density is shown in color and is overplotted with magnetic field lines. For the GEM Challenge boundaries (left column), the islands of reconnected magnetic flux grow freely at first (top panel) until they begin to interact with the walls at the top and bottom of the domain (middle panel). The increasing pressure in the island exerts a growing force on plasma ejected from the X-point and suppresses further reconnection (bottom panel).

In contrast, the panels in the right column of Fig. 3 show results from an otherwise-identical simulation that employs the shifted Klein boundary conditions described above. By the time of the middle panel the effects of the boundary conditions are already noticeable in the magnetic field lines since, unlike in the GEM Challenge simulation, they are free to pass through the top and bottom boundaries. By late time (bottom panel), the islands of reconnected flux have grown substantially and, due to the shift at the boundaries, continue to drive reconnection at the X-line.

The lack of hard walls at the top and bottom boundaries makes the Klein boundary conditions similar to the doubly periodic case. Figure 4 presents a comparison of several quantities at $tΩcp=60$⁠, a time of robust reconnection, from runs with both boundary conditions. The panels from the Klein boundaries (left column) show the entire domain, while those from the doubly periodic case (right column) show only one of the two current sheets. The panels depict, from top to bottom, $Jez$⁠, $vex$⁠, $vix$⁠, $Bz$⁠, and $Ey$⁠. Each horizontal pair of images uses the same color scale and so may be directly compared. The two columns are quite similar, suggesting that the Klein boundary conditions do not significantly affect the evolution of the system.

We next compare all three runs in Fig. 5, which plots the reconnected flux, $Δψ$⁠, and the reconnection rate, $d(Δψ)/dt$⁠, vs time for the different cases: GEM Challenge in blue dashed, doubly periodic in red dot-dashed, and shifted Klein in black. Contours of $ψ$ are the magnetic field lines shown in Figs. 1 and 4, and in the first and third cases, $Δψ$ is merely the difference in $ψ$ between the X- and O-points. In the doubly periodic case, the curve plots the average of the differences observed in the two current sheets. Each of the simulations exhibit an initial overshoot where the reconnection rate spikes to $≈0.25$⁠, but in the GEM Challenge case that spike is followed by a steady decline, with reconnection essentially ending after $tΩcp≈50$ when $Δψ$ plateaus. In contrast, the doubly periodic and shifted Klein runs both settle into a state with reasonably constant reconnection between $tΩcp=40$ and $tΩcp=80$ with $d(Δψ)/dt≈0.15$⁠. The large excursions near $tΩcp=60$ correspond to the time when the once-reconnected field lines surrounding the magnetic islands interact with an X-point a second time, either that on the opposite current sheet for the doubly periodic case or the original one for the shifted Klein (see Fig. 1). It is interesting to note that the reconnection rate after these fluctuations roughly remains unchanged.

The elapsed time (measured on a clock) for the three computations differed by less than $3%$⁠. However, the GEM Challenge and shifted Klein cases were run on the same number of processors, while the doubly periodic case, for which the computational domain was doubled in size, used twice as many. Hence, while the doubly periodic and shifted Klein cases found very similar physical results, the former incurred twice the computational expense of the latter. (For the relatively small simulations considered here, p3d exhibits near-perfect weak scaling and so any extra overhead incurred by doubling the number of processors is minimal.)

Adding a spatially constant guide field, $Bg$⁠, to the initial conditions of the GEM Challenge and doubly periodic cases is straightforward. The Klein case is more subtle. Recall that the x and z components of vectors are flipped when crossing the Klein boundary. As a consequence, a spatially constant guide field would exhibit a discontinuous jump—a grid point at the top or bottom boundary with $Bz=Bg$ would be adjacent to one with $Bz=−Bg$⁠. Such a configuration would generate a large current and lead to instability. This behavior arises because Klein bottles are non-orientable, which means there is ambiguity in the definition of “out of the page” and “into the page.” (This is perhaps most easily seen with the Möbius strip. If one takes an arrow pointing normal to the surface and translates it once around the strip to the starting location, the head will point in the opposite direction.)

To sidestep this obstacle, we can slightly modify the boundary condition. Write $Bz$ as the sum of a constant background $Bz0$ and a perturbation: $Bz=Bz0+δBz$⁠. Then, rather than apply the Klein boundary to the full magnitude, only change the sign of the perturbation $δBz$ when it crosses the boundary. Note that in the run discussed in Sec. III A, $Bz0=0$ so this procedure is identical to modifying the full $Bz$⁠. It might be argued that, for consistency, such a modification should be applied to any z component (e.g., $Ez$ or $Jz$⁠), but in the initial states considered here all other quantities have no uniform component and so making such an adjustment would have no effect.

We now compare results from three simulations with $Bg0=1$⁠. Figure 6 shows the evolution of the GEM Challenge boundary conditions (left) and the Klein boundaries with the modification discussed above. As is typical, the addition of a guide field breaks a symmetry of the system and produces differences between the separatrices, while the reconnection layer becomes more likely to form plasmoids, as can be seen in the middle panels. As with the anti-parallel case shown in Fig. 1, the GEM Challenge system is restricted by the walls and stagnates by the last panel. The Klein boundaries again allow flux to pass through them as can be seen by the field lines intersecting the top and bottom of the domain.

A comparison between the Klein and doubly periodic runs shows more significant differences than in the anti-parallel case. Figure 7 compares the two in the same format as Fig. 4, but at an earlier time, $tΩcp=45$⁠. The reason for this choice will be discussed below. The plasmoids in the Klein case distort the reconnecting layer, but the overall structure is generally similar. (Note that the color bar is the same for each pair of images so direct comparisons can be made.) The largest differences come in panels G and H, which show $Bz$ in the two runs. While the fields in the islands of reconnected flux are roughly comparable, bands of low-strength (blue) field are apparent just upstream of the current sheet in the doubly periodic case. These features represent the last vestiges of unreconnected plasma being pushed toward the X-line by an island expanding on the other (unseen) sheet. These features do not appear in the Klein boundary case, although hints are visible at the top and bottom of the domain. Still, the boundary conditions are clearly beginning to retard the inflow.

The effects of the various boundary conditions can also be seen in Fig. 8, which shows the reconnected fluxes and reconnection rates for the three runs in the same format as Fig. 5. After the peak at $tΩcp≈25$⁠, the GEM Challenge run essentially stagnates, with the reconnection rate falling to zero. The other two cases exhibit a similar earlier peak and decline, but both have a later period of increasingly fast reconnection with the Klein case peaking at $tΩcp≈45$⁠, the time shown in Fig. 7. After this point, however, the Klein system stagnates, while the doubly periodic system continues to enjoy robust reconnection. It is worth noting that simulations are usually halted once current sheets begin interacting with each other, so the lack of agreement after $tΩcp≈50−55$ is not significant from a practical standpoint.

The elapsed times for the three computations were again similar, differing by less than $4%$⁠. However, as previously, the doubly periodic case incurred twice the computational expense due to its larger domain.

We have shown that a novel boundary condition based on the topology of a Klein bottle can be used in numerical simulations of magnetic reconnection. The new boundary condition combines the advantages of other common options by being both self-reinforcing (as in systems with doubly periodic boundary conditions) and computationally compact (as seen in the single current sheet of the GEM Challenge boundaries).

The extension of Klein boundary conditions to full three-dimensional (3 × 3v) simulations should be straightforward. The additional dimension can be physically pictured as adding a “thickness” to the Klein bottle in Fig. 2. From a computational perspective, the most significant change is the identification of neighboring processors across the Klein boundary. Neighboring processors in z must change in the manner shown in panel (ii) of Fig. 3—recall that the z component of vectors flips across the Klein boundary—but do not have to undergo the shift of panel (iii).

However, regardless of the dimensionality, Klein boundary conditions do have limitations. Symmetric reconnection, in which the asymptotic plasmas on the two sides of the current sheet share the same characteristics, is an idealization, albeit a common one, that sidesteps many of the complexities of actual systems. In reality, asymmetries can exist in any or all of the plasma density, plasma temperature, strength of the reconnecting component of the field, strength of the guide field, etc., and at some locations (e.g., the terrestrial magnetopause), reconnection is nearly always asymmetric. Such complications are straightforward to include with GEM Challenge or doubly periodic boundary conditions because the two sides of a reconnecting current sheet represent distinct plasmas. For Klein boundaries, there is no distinction between the two sides, and so simulating asymmetric reconnection while maintaining the associated benefits (low computational cost and non-stagnation) is not straightforward.

Finally, the topological equivalence of the doubly periodic system to a torus implies the existence of certain conservation laws. Specifically, consider the surface integral of Faraday's Law, $∂B/∂t=K∇×E$⁠, where the proportionality constant K depends on the choice of units. Performing a surface integral and invoking Stokes' theorem transform the right-hand side into a boundary integral, which vanishes because a torus has no boundary. Hence, the magnetic flux through the plane of the simulation must be constant, which can be verified numerically. A Klein bottle also has no boundary, but it is non-orientable, and so Stokes' theorem does not apply. As a consequence, the magnetic flux through the plane of the simulation can, and does, change when Klein boundary conditions are employed. While this difference does not appear to affect the dynamics of reconnection, further research will be necessary to explore whether it has physical significance.

The authors acknowledge the support of NASA Grant Nos. 80NSSC22K0352 and 18-DRIVE18_2-0029, “Our Heliospheric Shield (Grant No. 80NSSC22M0164). The authors also acknowledge UMD's TREND program, funded through No. NSF PHY-2150399. The authors would like to thank P. R. Colarco and E. A. Colarco for providing visual inspiration. The simulations were carried out at the National Energy Research Scientific Computing Center (NERSC). The data used to perform the analysis and construct the figures for this paper are preserved at https://doi.org/10.5281/zenodo.12785961.

The authors have no conflicts to disclose.

Luke Xia: Investigation (equal); Methodology (equal); Software (equal); Writing – original draft (equal). M. Swisdak: Funding acquisition (lead); Investigation (equal); Methodology (equal); Project administration (lead); Software (equal); Writing – original draft (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.