Deep inelastic vortex scattering: A third outcome for head-on collisions

Results are presented from numerical simulations of the flat-space nonlinear Maxwell-Klein-Gordon equations demonstrating deep inelastic scattering of $m=1$ vortices for a range of Ginzburg-Landau (or Abelian-Higgs) parameters ($\kappa$), impact parameters ($b$), and initial velocities ($v_0$). The threshold ($v_0^*$) of right-angle scattering is explored for head-on ($b=0$) collisions by varying $v_0$. Solutions obey time-scaling laws, $T\propto \alpha\ln(v_0-v_0^*) $, with $\kappa$-dependent scaling exponents, $\alpha$, and have $v_0^*$ that appear not to have the previously reported upper bound. The arbitrarily long-lived static intermediate attractor at criticality ($v_0=v_0^*$) is observed to be the $\kappa$-specific $m=2$ vortex solution. Scattering angles are observed for off-axis ($b\neq 0$) collisions for a wide range of $b$, $v_0$, and $\kappa$. It is shown that for arbitrarily small impact parameters ($b\rightarrow 0$), the unstable %but arbitrarily long-lived $\kappa$-dependent $m=2$"critical"vortex is an intermediate attractor and decays with a $\kappa$-\emph{independent} scattering angle of $135^{\circ}$, as opposed to either of the well-known values of $180^{\circ}$ or $90^{\circ}$ for $b=0$.


I. INTRODUCTION
Topological defects are known to appear in a wide variety of physical models that involve spontaneous symmetry breaking of gauge field theories. In three spatial dimensions, they can appear in the form of monopoles, strings/vortices, or domain walls and are often lumped together with their slightly more exotic cousins, textures [1]. This paper discusses vortices that are solutions to the axisymmetric Maxwell-Klein-Gordon equations of motion. They have a conserved topological charge and are characterized by a "winding number" that counts the number of times the phase of the complex scalar field rotates when propagated along a closed loop around the vortex axis. Such vortices are present in a range of contexts, including the Abelian-Higgs model of particle physics and the Ginzburg-Landau model of condensed matter physics. While these models have very different physical interpretations, they reduce to the same dimensionless action and have the same dimensionless Maxwell-Klein-Gordon equations of motion, thus giving rise to common phenomenology of great interest across many parts of the physics community.
Nielsen and Olesen explored vortices in the Abelian-Higgs model to better understand the possible string nature of fundamental particles, noting the approximate correspondence between a vortex and a Nambu string [2]. A great number of authors have studied Abelian-Higgs vortices in the context of early universe cosmology, where they can appear as topological defects (cosmic strings) when a grand unified theory has undergone spontaneous symmetry breaking [3].
While not yet experimentally observed, cosmic strings are believed by many to have been present in the early universe, where they may have influenced the large-scale structure currently observed in the universe [4][5][6]. * ehonda@alum.mit.edu Abrikosov vortices, on the other hand, are solutions to the Ginzburg-Landau equations of motion that describe magnetic flux tubes within superconductors. They were first predicted by Abrikosov in 1957 [7,8] and later observed in the laboratory by Cribier et al. by means of neutron diffraction in 1964 [9]. This paper adopts the conventions typically used in the context of Ginzburg-Landau vortices, where vortices and the "bulk" supporting them are characterized by a number of parameters. The Ginzburg-Landau parameter is given by κ = Λ/ξ and is the ratio of the London penetration depth of the magnetic field (Λ) to the coherence length of the complex scalar field (ξ). Superconductors with κ < 1/ √ 2 are referred to as Type-I superconductors; vortices in Type-I superconductors attract one another and typically coalesce, breaking down the superconductivity when in abundance. Superconductors with κ > 1/ √ 2 are referred to as Type-II superconductors; vortices in Type-II superconductors repel one another (do not coalesce) and form flux tube (Abrikosov) lattices when in abundance. Superconductors with κ = 1/ √ 2 are referred to as having critical coupling, and they support vortices that are noninteracting.
As topological solitons, vortices have a conserved topological charge and are characterized by their winding number, m. The equations of motion do not allow an m = 1 vortex to "unwind" itself, which ensures stability independent of κ. However, the stability of vortices with m > 1 is a κ-dependent phenomenon. Bogomol'nyi demonstrated using a trial function approach that vortices with κ ≤ 1/ √ 2 can be stable for m > 1, while vortices with κ > 1/ √ 2 are unstable for m > 1 [10]. The repulsive interaction and instability of vortices with m > 1 in a Type-II bulk gives rise to very interesting and nontrivial scattering.
In 1988, it was numerically demonstrated that the head-on collisions of two vortices could lead to deep inelastic right-angle scattering [11,12]. This remarkable phenomenon was predicted for a number of models for the noninteracting (critical coupling) cases where one can obtain Bogomol'nyi equations that reduce secondorder partial differential equations (PDEs) to simpler, first-order PDEs [13]. While this approach had the benefit of greatly reducing computational demands at a time when such calculations were prohibitive, the results were only relevant at or near critical coupling. The noninteracting (κ = 1/ √ 2) vortices would scatter at right angles, independent of velocity, while repulsive (κ > 1/ √ 2) vortices would transition from direct backscattering to right-angle scattering at a κ-dependent velocity. In 1992, Myers et al. performed a numerical study [14] of vortex scattering, clearly demonstrating right-angle scattering for repulsive vortices with coupling as strong as κ = 2 √ 2 (λ MRS = 32 in their Abelian-Higgs model). While the authors pushed the limits of computation at the time, the breadth of (b, v 0 ) parameter space, the resolution of simulations, and the range of κ were limited. This paper discusses new results and enhances previous results of vortex-vortex scattering within this model. The previously understood bound on the critical velocity, v * 0 (κ), that marks the transition from backscattering to right-angle scattering is shown not to exist. Further, it is shown that there exist intermediate attractor solutions that are arbitrarily long-lived and demonstrate a third type of scattering in the limit where the impact parameter goes to zero. While the results are applicable to both the Abelian-Higgs or Ginzburg-Landau models, conventions are adopted that are most consistent with the Ginzburg-Landau model.
The remainder of the paper is organized as follows. In Sec. II, the formalism is defined and the fully general equations of motion are presented. In Sec. III, the initial data and evolution equations are discussed. Section IV discusses head-on collisions and observed scaling laws and shows that in the center of mass frame, a static m = 2 vortex is the intermediate attractor for the threshold of backscattering and right-angle scattering. Section V presents results from off-axis collisions and discusses a new κ-independent scattering phenomenon. Appendix A shows how to obtain the dimensionless Maxwell-Klein-Gordon action from the dimensionful Ginzburg-Landau action and shows how to relate that dimensionless action to the actions in [12] and [14]. Appendix B discusses the computational methods used.

II. GENERAL FORMALISM, DEFINITIONS, AND CONVENTIONS
The dimensionless Maxwell-Klein-Gordon Lagrangian being studied here is easily obtained from the Abelian-Higgs or Ginzburg-Landau models by transforming to dimensionless variables (Appendix A), where the gauge covariant derivative and the Maxwell field tensor are given by respectively. The Maxwell equations are given by where is the Maxwell four-current. The Klein-Gordon equations for the components of the complex scalar field are given by where is the covariant d'Alembertian. Finally, it is worth noting that the signature of the metric is (−+++) and all quantities represent dimensionless variables.

III. INITIAL DATA AND EQUATIONS OF MOTION
Any analysis of vortex scattering begins with a single static vortex solution. The ansatz for an axisymmetric static vortex with winding number m is given by This work discusses scattering of m = 1 vortices with κ > 1/ √ 2 and only considers solutions with zero charge (A t = E R = 0) at t = 0. Assuming the Lorentz gauge, Eqs. (4), (5), (7), and (8) in cylindrical symmetry reduce to the following equations: The inner (R = 0) boundary conditions are given by since regularity at R = 0 implies that φ(R) ≈ cR m for R ≈ 0. The initial conditions χ 0 and (B z ) 0 are varied while "shooting" on the large-R boundary conditions, As is typical for shooting methods, the fields diverge before reaching the outer boundary, so solutions were smoothed and matched to their known large-R behavior.
To numerically simulate a collision of two vortices, individual static vortices must be boosted and combined. The vortices are assumed to be positioned symmetrically about the origin and boosted at one another in the x direction with equal magnitude velocities, v 0 , and with a given impact parameter. The Maxwell fields of the individually boosted vortices are added together, while the complex scalar fields are combined by multiplying the magnitudes and adding the phases of the individually boosted vortex solutions. The Π i = ∂ t φ i are obtained by computing the time derivative of the Lorentz-boosted static-vortex scalar field solutions in the boosted frame. These initial data are then time evolved with the hyperbolic equations of motion in Cartesian coordinates. The equations of motion for the electric and magnetic fields are obtained from (4) and (5), giving The evolution equations for the components of the gauge potential are obtained from the definition of the Maxwell tensor, and the Lorentz gauge condition, The complex scalar field evolution equations are obtained from (7) and (8), giving while the definition of Π i and commutation of partial derivatives give For all simulations discussed in this paper, these data are evolved with vortices starting at positions (−x 0 , y 0 ) and (x 0 , −y 0 ) with velocities (v 0 , 0) and (−v 0 , 0), respectively. The impact parameter b = 2y 0 is defined to be the distance between the vortices perpendicular to their collision velocity.

IV. CRITICAL HEAD-ON SCATTERING
This section builds upon, and in one case corrects, the current understanding of head-on (b = 0) scattering of repulsive vortices. It was clearly shown by [14] that repulsive vortices backscatter, θ s = 180 • , for low-velocity collisions, v 0 < v * 0 , and right-angle scatter, θ s = 90 • , for high-velocity collisions, v 0 > v * 0 , for some critical velocity, v * 0 . However, based on the simulations performed at the time, it appeared and was conjectured that the critical velocity approached an asymptotic value of v 0 ≈ 0.56 for large κ. Those investigations used couplings of 1 ≤ λ MRS ≤ 32 in the Abelian-Higgs model 1 , which corresponds to 1/2 ≤ κ ≤ 2 √ 2 in the Ginzburg-Landau model, and used a coherence length to the lattice spacing ratio (ξ/dx) of 5:1 at κ = 2 √ 2. This work uses a ξ/dx ratio of greater than 10:1 and maintains this ratio for increasing κ (decreasing ξ) while keeping the domain large enough to also capture the magnetic field dynamics that remain on the length scale Λ ≈ 1. This approach maintains energy conservation through the life of the simulations to a few parts in 10 3 for high values of κ and a few parts in 10 4 for smaller κ (Appendix B).
The threshold of backscattering and right-angle scattering is first explored for κ = 1. Initial values of v 0 are selected that are observed to be above and below v * 0 and are referred to as v + 0 and v − 0 , respectively. Vortices colliding with velocity v + 0 exhibit right-angle scattering, while vortices with velocity v − 0 backscatter. A bisecting procedure was performed to bracket v * 0 with v ± 0 values that become successively closer to v * 0 until machine (double) precision is reached. As v 0 approaches v * 0 from either side, one observes longer interaction; the two m = 1 vortices become a single m = 2 vortex with a lifetime that grows as v 0 gets closer to v * 0 . Figures 1 and 2 demonstrate the time evolution of head-on vortex scattering simulations for v 0 just below and above v * 0 , referred to as subcritical and supercritical scattering, respectively. observed to collide, overcoming the mutual repulsion, at t ≈ 17. For times 17 t 92, the m = 2 vortex is approximately static, and at t ≈ 92, the m = 2 vortex decays into two m = 1 vortices with θ s = 180. Similarly, Fig. 2 demonstrates supercritical scattering with v 0 tuned to one part in 10 13 above v * 0 . The vortices are observed to collide, to overcome mutual repulsion, to coalesce into an m = 2 vortex, and to eventually decay into two m = 1 vortices, only this time with θ s = 90. Based on conservation of topological charge, the presence of an approximate m = 2 vortex was well understood to exist previously [12]. However, the ability to fine-tune the collision to allow a static m = 2 vortex to be created and exist for an arbitrary amount of time is new to this work.
To confirm the presence of an m = 2 static vortex as the intermediate attractor, it is useful to consider the dynamic energy of the scalar field, which is defined here simply to be contribution of the kinetic energy of the scalar field from only the time derivative terms, i.e., without the "shape" energy from the spatial derivative terms. Figure 3 demonstrates the dynamic energy in the scalar field as a function of time for the same subcritical collision shown in Fig. 1. shows the components of the complex scalar field (left) and the phase of the scalar field (right) along a circular path centered on the vortex with radius equal to the coherence length, ξ = 1, for before, during, and after the collision. At times t = 0 and t = 110, the phase of the scalar field rotates once (m = 1) around the left vortex, whereas at t = 62.25, the phase rotates twice (m = 2) around the single vortex at the origin. While this result is not surprising, it definitively shows that the intermediate state during the deep inelastic scattering is indeed a vortex with winding number m = 2. Figures 1, 2, and 3 also begin to show the longevity of the m = 2 vortex and its role as an intermediate attrac-tor. The upper-right graph of Fig. 4 demonstrates the "lifetime" of each collision as a function of ln |v 0 − v * 0 | for κ = 1. Since all of the evolutions for a given κ start out with vortices in the same positions and end based on common exit criteria, the relative increase in lifetime is a good measure of the lifetime of the m = 2 intermediate vortex. It is clear that a time-scaling relationship is observed, where the slope α governs the scaling and the y intercept is ignored because it depends on the initial positions and exit criteria. This type of time scaling has been observed in many contexts, including the gravitational collapse of Yang-Mills and boson fields [15][16][17] and the selfinteraction of nonlinear oscillons in flat space [18][19][20]. In each of these cases, the time scaling describes the life of an unstable intermediate attractor solution that can be observed by tuning a parameter that results in one of two possible outcomes. In the case of gravitational collapse, the threshold is the creation (or not) of a black hole. In the case of oscillons, the threshold appears at resonances that separate solutions with n and n + 1 "shape mode" modulations. In this work, the discrete outcomes are right-angle and backward scattering resulting from the head-on (b = 0) collision of two m = 1 vortices. While the above discussion is for κ = 1 vortex collisions, there exist similar time-scaling relations for all of the repulsive (κ > 1/ √ 2) vortices explored. Figure  4 demonstrates the scaling laws for a wide range of κ, from κ = 0.8 to κ = 80. This dynamic range in κ corresponds to 2.56 ≤ λ MRS ≤ 25, 600 in the Abelian-Higgs model studied in [14], where the authors explored 1 ≤ λ MRS ≤ 32. Table I captures the results for the critical velocities, time-scaling exponents, grid size, and dimensions for each κ explored.
One observation that can immediately be made from Table I is that the previously understood upper bound on critical velocity of v * 0 ≈ 0.56 does not exist.
where the best fit gives m γ = 0.184, c γ = 0.97, and where κ 0 ≡ 1/ √ 2. While the factor of κ 0 is redundant and could be absorbed into c γ , this form nicely results in the relevant domain of ln(κ/κ 0 ) being [0, +∞]. The observed value of c γ ≈ 1 also suggests that which would imply that as κ approaches κ 0 , v * 0 approaches zero. This clear relationship not only provides new insight into vortex scattering, but it also provides convincing evidence that the only upper bound to v * 0 (κ) would be the speed of light, and this would occur in the limit of infinite κ. It is this author's belief that the work of [14] was groundbreaking and incredibly thorough for its time, but it is likely that the simulations performed therein were underresolved. The higher ratio of vortex size to lattice spacing used here for the entire range of κ being explored appears to be necessary to fully resolve vortex-vortex interactions.
Another observation one can make from Table I (and Fig. 6) is that there appears to be a relationship between the time-scaling exponent and the Ginzburg-Landau parameter, where the best fit gives m α = −0.85, c α = 0.02, and again, where κ 0 ≡ 1/ √ 2. This suggests that the timescaling exponent, α, is (nearly) inversely proportional to the Ginzburg-Landau parameter, κ, In the limit where κ approaches κ 0 , the vortices become less and less repulsive, the critical velocity goes to zero, and the lifetime of the intermediate attractor goes to infinity. This is actually just the trivial case of two noninteracting vortices located at the same position at rest. In the limit where κ approaches infinity, the time-scaling exponent goes to zero. This results in shorter lifetimes for larger κ, as seen in Fig. 4. This should also be intuitively expected in that with larger κ, the vortices are more repulsive, and it becomes harder to suppress (through fine-tuning v 0 ) the increasingly dominant unstable mode, where the scaling exponent α is inversely proportional to the Lyapunov exponent of the unstable mode, similar to [21,22]. These results not only adjust the previous understanding of head-on (b = 0) vortex collisions (the upper bound of v * 0 ), but they also provide additional insight into headon collision dynamics that were not previously understood. It is definitively shown that the m = 2 vortex is the intermediate attractor on the threshold of rightangle scattering, and more interestingly, the new scaling relations observed here strongly suggest that with appropriate fine-tuning of v 0 , one can predictably create an arbitrarily long-lived m = 2 vortex from the head-on collision of two m = 1 vortices.

V. OFF-AXIS SCATTERING
This section presents new results for off-axis scattering of m = 1 vortices. Initial data are again prepared as described in Sec. III, only now for b = 0. Figure 7 shows results of scattering angle as a function of impact parameter, θ s (b), for κ = 1/ √ 2 (top) and κ = 1 (bottom) with linearly sampled b, and where each θ s (b) curve is for a given v 0 . The θ s (b) for noninteracting vortices converge to a single point, θ s = 90 • , at b = 0, whereas the θ s (b) for repulsive vortices converge to either θ s = 90 • or θ s = 180 • at b = 0, depending on whether v 0 is above or below the critical velocity (v * 0 ≈ 0.271 for κ = 1), respectively. 5 different v 0 . As such, these data enhance the results of [12] for noninteracting (λ SR = 1) collisions and [14] for repulsive (λ MRS = 4) vortices and serve as a starting point to better understand the off-axis deep inelastic scattering of vortices.
By observing θ s (b) curves with v 0 approaching v * 0 , one can explore the threshold of right-angle scattering by using two dynamical degrees of freedom, b and v 0 . to more easily observe the behavior in the b → 0 limit. To allow Fig. 7 (bottom) and previous work [14] to be easily compared to Fig. 8, the top graph in Fig. 8 displays data as a function of linear b while the bottom graph displays the same data as a function of ln(b). While the same v 0 used in Fig. 7 are used in Fig. 8, an additional 26 values of v 0 are also used that logarithmically approach v * 0 from above and below in constant intervals of ln |v where n ranges from 2 to 14 with ∆n = 1, and the b are sampled at 10 samples per decade for a total of 141 points for each θ s (b), resulting in data from 4371 vortex-vortex collisions. It is immediately apparent that there exists more structure in the θ s (b) curves than was originally observed in Fig. 7 and [12,14]. Figure 8 clearly demonstrates that as v 0 approaches v * 0 from above or below, the θ s (b) "latch on" to an intermediate attractor solution that plateaus at θ s ≈ 135 • . Figure 9 shows the time  Figure  10 plots ln(b ± ) as a function of ln |v 0 − v * 0 |, which clearly shows that as v 0 approaches v * 0 , the b ± go to zero. Fitting the data to results in a good fit to a straight line with best-fit values of m b = 1.00 and c b = 2.56.   Figure 11 shows θ s (b) curves for κ = 4 (top) and κ = 10 (bottom), for v 0 = 0.2, 0.4, 0.6, 0.8, and (1 ± 10 −n ) · v * 0 , for n = 2 through n = 14 at an increment of ∆n = 2, and with b sampled at 10 samples per decade, for a total of 2538 collisions for each κ. Figure 12 shows θ s (b) curves for κ = 40 (top) and κ = 80 (bottom), for v 0 = 0.8, 0.85, and (1 ± 10 −n ) · v * 0 , for n = 3 through n = 12 at an increment of ∆n = 3, and with b sampled at five samples per decade, for a total of 710 collisions for each κ. While the density of b samples and the number of v 0 explored are decreased with increasing κ due to significantly increasing computational demands (while maintaining ξ/dx > 10), one begins to see common phenomena across the entire range of κ explored.
For all κ observed, as v 0 approaches v * 0 , the θ s (b) curves approach an intermediate attractor solution that plateaus at θ s = 135 • with b ± that approach zero. Figure  13 displays plots of ln(b ± ) as a function of ln |v 0 − v * 0 | for κ = 4, 10, 40, and 80. It is again clear that there is a linear relationship and (46) appears to accurately describe b ± , especially as v 0 approaches v * 0 . Table II  best-fit values for the range of κ explored. The numerical closeness of m b to 1 suggests the possibility that m b = 1 exactly, which would give a simple relationship, where c b = e c b . Independent of whether m b = 1 exactly, it has been shown that at the threshold of right-angle scattering, θ s (b) plateaus to θ s = 135 • and b ± → 0. Since the b ± at which θ s (b) diverges from 135 • goes to zero, θ s (b) can be shown to be 135 • for an arbitrarily small b. This behavior is demonstrated for a range of κ and persists for arbitrarily small b, as long as v 0 is sufficiently fine-tuned to v * 0 . So, while precisely head-on (b = 0) scattering only allows θ s = 90 • or θ s = 180 • , at criticality (v 0 = v * 0 ) there exists an intermediate attractor m = 2 vortex solution that decays with a scattering angle of θ s = 135 • for arbitrarily small impact parameters.

VI. CONCLUSIONS
Results have been presented from numerical simulations of the flat-space nonlinear Maxwell-Klein-Gordon equations with a spontaneously broken symmetry demonstrating deep inelastic scattering of m = 1 vortices for a wide range of Ginzburg-Landau parameters (κ), impact parameters (b), and initial velocities (v 0 ).
Head-on (b = 0) collision simulations were conducted that explored the threshold of right-angle scattering by varying the initial velocity. While the fact that repulsive vortices exhibit right-angle scattering above a critical velocity (v * 0 ) has been known for decades, until this work, the threshold of the transition from backward scattering to right-angle scattering has not been thoroughly explored or understood. This work corrects previous work that appeared to demonstrate an upper bound to the critical velocity. In addition to this clarification, a new phenomenon was observed for b = 0 collisions. It was shown that solutions obey time-scaling laws, T ∝ α ln |v 0 − v * 0 |, with κ-dependent scaling exponents α that vary according to a simple power-law dependence on κ − κ 0 , where κ 0 = 1/ √ 2 is the Ginzburg-Landau parameter for noninteracting vortices. While expected, it was definitively shown that the intermediate attractor solution on the threshold is the m = 2 vortex solution for the κ being investigated.
Off-axis (b = 0) collision simulations were conducted to gain new insight into vortex-vortex scattering, especially for small impact parameters and with v 0 ≈ v * 0 . It is shown that for arbitrarily small (nonzero) impact parameters, the unstable intermediate attractor m = 2 vortex decays with a κ-independent scattering angle of 135 • . This was demonstrated through careful analysis of scattering angle curves θ s (b) and their divergence away from θ s = 135 • in the v 0 → v * 0 and b → 0 limits. Simply put, for an appropriately fine-tuned initial velocity, an arbitrarily small impact parameter can be shown to scatter at θ s = 135 • . As such, this work has demonstrated the existence of a third outcome for "head-on" (b → 0) vortex collisions.

Appendix A: Dimensionless Lagrangian
The Ginzburg-Landau system is described by the following dimensionful Lagrangian containing a Maxwell field coupled to a complex scalar field in a gauge covariant way with a spontaneously broken U (1) symmetry: where A µ is the Maxwell gauge potential, φ is the order parameter, and |φ 2 | is the number of superconducting particles per unit volume. Following standard practice, to convert (A1) to a dimensionless Lagrangian, one converts each coordinate and field variable to a dimensionless form, where one can use the standard definitions, to obtain the dimensionless Lagrangian, whereD µφ = 1 κ∂ µ − iÂ µ φ , and (A12) The Lagrangian (A11) is used throughout this paper but with the hats removed for clarity. The work by Shellard and Ruback [12] uses the Lagrangian while the work by Myers et al. [14] uses the Lagrangian These two Lagrangians are physically identical to (A11), and they can all be compared to one another with the conversions in Table III. A quantity that is commonly discussed is the Ginzburg-Landau (or Abelian-Higgs) coupling parameter for noninteracting vortices that separates attractive (Type I) vortices from repulsive (Type II) vortices. Table III (bottom row) shows that κ = 1/ √ 2 corresponds to λ SR = 1 and λ MRS = 2e 2 (where the convention e = 1 is usually used).   [12,14]. Conversion from one model to another is performed by substituting values from one column with their respective values from another column.

Appendix B: Computational Methods
This appendix briefly discusses the numerical methods used in this paper. All vortex collisions were performed using a second-order implicit Crank-Nicolson finite difference scheme coded in C++ and CUDA and were run on a single NVIDIA GeForce RTX 2060 Super graphics processing unit (GPU) with 8 GB of RAM. This GPU platform was ideally suited for these calculations and the grid sizes used; all of the grid functions could be stored directly on the GPU device simultaneously (roughly 60% memory utilization for the largest, κ = 80, grids), resulting in minimal inefficiencies from moving data over the PCI-e bus.
Given the large dynamic range of length scales, where the magnetic fields change on the order of Λ ≈ 1 while the scalar field changes (at a vortex) on the order of the coherence length (ξ = 1/κ), it is imperative that one use fine enough grid resolution to properly resolve the salient features of the system. This was done by maintaining a ratio of ξ/dx > 10 for all simulations, which led to simulations on large grids (Table I).
A common method used to determine solution quality is to monitor how well energy is conserved in the simulation over time. This is done by integrating the energy density in the spatial grid at any time and adding that to the integrated momentum flux along the boundary from the beginning of the simulation to the time of interest. The energy density of the system is given by the sum of the energy densities, ρ mass = 1 2 φ 2 1 + φ 2 2 A 2 t + A 2 x + A 2 y + A 2 z , and (B4) where ρ V 0 is set such that ρ V = 0 when φ is at the expectation value of the spontaneously broken symmetry. The momentum densities of the system are given by j cur j mass (B13) Figure 14 shows the final mass error, where M final is the integrated energy density at the end of the simulation, M initial is the integrated energy density at t = 0 , and M radiated is the integrated outgoing momentum flux over the outer boundary of the computational domain over the entire time of the simulation. Most of the simulations conducted had final M error between 10 −4 and 10 −3 . While the ratio ξ/dx > 10 was maintained for all simulations, the increase in v * 0 as κ increases leads to Lorentz contraction of the vortex in the x direction by a factor of γ * 0 that reduces the effective number of grid points per vortex in the direction of motion, leading to lower effective resolution and higher final M error values for larger κ. Despite this reduction in resolution, the simulations conducted here have typical mass errors that are still 1, if not 2, orders of magnitude lower than previous work.