We investigate the influence of background shear flow on linear resistive tearing instabilities with Joule heating for two compressible plasma slab configurations: a Harris current sheet and a force-free, shearing magnetic field that varies its direction periodically throughout the slab, possibly resulting in multiple magnetic nullplanes. To do so, we exploit the latest version of the open-source, magnetohydrodynamic spectroscopy tool Legolas. Shear flow is shown to dramatically alter tearing behavior in the presence of multiple magnetic nullplanes, where the modes become propagating due to the flow. Finally, the tearing growth rate is studied as a function of resistivity, showing where it deviates from analytic scaling laws, as well as the Alfvén speed, the plasma-β, and the velocity parameters, revealing surprising nuance in whether the velocity acts stabilizing or destabilizing. We show how both slab setups can produce growth rate regimes, which deviate from analytic scaling laws, such that systematic numerical spectroscopic studies are truly necessary, for a complete understanding of linear tearing behavior in flowing plasmas.

In many plasma configurations, from astrophysical to laboratory settings, the most dramatic and violent events are often driven by magnetic reconnection. During the reconnection process, magnetic field lines break and reconnect, drastically altering the magnetic topology, e.g., by unraveling a knotted field structure. Consequently, stored magnetic energy is converted into thermal and kinetic energy, accelerating particles and producing matter outflows. It is frequently observed in various phenomena, like coronal mass ejections,1 the heliospheric current sheet,2 and Earth's magnetosphere.3 

Due to its importance and ubiquity, reconnection has been a topic of great debate ever since Sweet4 and Parker5 first proposed their model, which was later followed by Petschek's model6 to achieve faster reconnection rates in accordance with observations. Irrespective of the specific reconnection process, its initialization in current sheets is intrinsically linked to tearing instabilities, which occur due to magnetic shear, and create the necessary reconnection points (X-points). These tearing instabilities come in collisional7 (resistive) and collisionless8 (mediated by the Hall effect and/or electron inertia) varieties, both of which have been studied extensively to understand their role in triggering fast magnetic reconnection.

To reach the fast reconnection regime, the influence of several effects on the resistive tearing mode has been studied for decades. For thin current sheets with a width on the order of the ion inertial length, the importance of the Hall term in the generalized Ohm's law was already shown forty years ago by Terasawa,9 who showed that it enhances the growth rate of the resistive tearing mode. Following this work, many more studies have honed in on the growth rate modification of resistive tearing in the Hall regime.10–15 However, the role of the Hall term in magnetic reconnection is not limited to the enhancement of the resistive tearing growth rate. It has also been linked to finite Larmor radius effects;16–18 it couples to the anisotropic electron pressure tensor, which enhances reconnection rates;19,20 and it can induce a transition to whistler-mediated reconnection.8,21,22 Finally, Liu et al.23 recently showed that the Hall term plays an important role in attaining the proper geometry for fast reconnection.

However, kinetic simulations suggest that the collisionless terms in the generalized Ohm's law do not dominate in collisional plasmas until the current sheet thins to the ion inertial scale,24,25 though collisionless reconnection does occur at larger scales in plasmas with negligible collisionality due to electron inertia.26 Hence, for thicker, sufficiently collisional current sheets, the initial unstable perturbation is expected to be of the resistive tearing variety, with a negligible contribution due to collisionless terms. How these sheets then reduce to the ion inertial scale to establish fast reconnection is not fully understood, though processes like fractal reconnection27 have been suggested.

This then begs the question which other factors might influence the growth rate outside of the Hall regime. It is known that the growth rate is affected by various environmental factors, notably by the background flow.28,29 In this endeavor, the influence of equilibrium flow on the resistive tearing mode has already been studied extensively using both analytic30–34 and numerical28,29,35–37 methods. Flow's influence is also not limited to the resistive tearing type, as its impact on collisionless tearing has been shown for flow parallel to the reconnecting field38 and to a guide field.39 Furthermore, the role of flow, and particularly flow shear, extends well beyond the modification of tearing growth rates, as it may introduce the Kelvin–Helmholtz instability (KHI) into the system. This instability regularly becomes the dominant instability,30,40 though it is also observed to co-exist with the plasmoid instability (secondary tearing) in current sheets,41 at plasma interfaces,42 and in turbulent reconnection.43 

A profound understanding of instabilities and how to suppress them is also of vital importance in laboratory plasmas, and particularly for fusion research. In fusion devices like tokamaks, tearing instabilities lead to the formation of magnetic islands, which in turn disrupt plasma confinement. There have been many experiments,44,45 linear studies,46–48 and non-linear simulations49–51 studying the effects of both toroidal and poloidal flow. In these torus-like geometries, it is now generally accepted that the toroidal flow shear has a stabilizing influence on the system by suppressing island formation.

In this work, we revisit the linear analysis of the resistive tearing mode in the presence of background flow using numerical means. Unlike earlier literature, we eliminate the need for approximations, like incompressibility, by employing the magnetohydrodynamic (MHD) spectroscopic code Legolas (Refs. 15, 52, and 53, https://legolas.science). With this code, we explore plasma stability parametrically for a selection of configurations and show that equilibrium flow can both enhance and suppress the growth rate of the resistive tearing mode depending on the parameter regime. The effect of equilibrium flow has been studied analytically by Chen and Morrison,34 and we connect our results to the power-law predictions from their study.

To study the linear properties of the resistive tearing instability, we linearize the dimensionless, compressible, resistive MHD equations
ρ t = · ( ρ v ) ,
(1)
ρ v t = p ρ v · v + J × B ,
(2)
ρ T t = ρ v · T ( γ 1 ) p · v + ( γ 1 ) η J 2 ,
(3)
B t = × ( v × B ) × ( η J )
(4)
around a one-dimensionally varying equilibrium (x-dependence only) and assume Fourier solutions
f 1 ( r , t ) = f ̂ 1 ( x ) exp [ i ( k y y + k z z ω t ) ]
(5)
for perturbed quantities f 1 { ρ 1 , v 1 , T 1 , B 1 }. Here, ρ, v, T, and B are the usual quantities density, bulk velocity, temperature, and magnetic field, respectively, with subscripts 0 and 1 differentiating between equilibrium quantities and perturbations. Furthermore, p = ρ T denotes the pressure, J = × B the current density, and η = 10 4, unless specified otherwise, and γ = 5 / 3 are the resistivity and adiabatic index, respectively. Finally, k = k y e ̂ y + k z e ̂ z and ω are the wave vector and frequency. The resulting algebraic problem is solved for the frequency and corresponding perturbed quantities by Legolas.

Note that our choice of energy equation, Eq. (3), including Joule heating, along with compressibility and a constant resistivity, differs from the isothermal closure used in the original derivation of the tearing mode.7 

Now consider a semi-infinite plasma confined in the x direction between two perfectly conducting plates, described in Cartesian coordinates. In such a plasma, we examine the resistive tearing instability for two separate equilibrium configurations: a Harris current sheet and a force-free magnetic field.

1. Harris current sheet

For the first configuration, we consider the popular Harris current sheet, as presented in Ref. 28. In their simulation setup, they consider a typical Harris current sheet
B 0 ( x ) = B c tanh ( x a B ) e ̂ y ,
(6)
which they supplement with a similar velocity profile
v 0 ( x ) = v c tanh ( x a v ) e ̂ y
(7)
and a uniform density ρ 0 = 1. The temperature profile is simply obtained by demanding that the total equilibrium pressure (the sum of plasma pressure and magnetic pressure) is constant,54 i.e.,
x ( ρ 0 T 0 ( x ) + 1 2 B 0 2 ( x ) ) = 0
(8)
with a solution
T 0 ( x ) = ( p c 1 2 B 0 2 ( x ) ) / ρ 0
(9)
for any constant p c, set to p c = 1 throughout this work, unless specified. The profiles are visualized in Fig. 1(a).
FIG. 1.

(a) The Harris current sheet equilibrium, Eqs. (6)–(9). (b) Equilibrium with the force-free magnetic field, Eqs. (11). Parameter values were chosen for visual clarity.

FIG. 1.

(a) The Harris current sheet equilibrium, Eqs. (6)–(9). (b) Equilibrium with the force-free magnetic field, Eqs. (11). Parameter values were chosen for visual clarity.

Close modal
Since B 0 ( x ) vanishes at x = 0, so does
F ( x ) = k | k | · B 0 ( x )
(10)
for any wave vector k, which we choose to be along the y axis, k = k y e ̂ y. Hence, the magnetic nullplane (or resonant plane), defined as the location where F vanishes,7 is always at x = 0 in this configuration, and magnetic reconnection will occur here.

The locations of the perfectly conducting walls at x w = ± 15 a B are chosen such that their effect on the tearing instability is negligible (according to Ref. 55, the effect of the conducting walls is negligible if they occur at a position | x w | 10 a B). Due to the sharp transitions in the equilibrium profiles near the origin and approximately constant behavior away from the center, the problem is solved using a non-uniform grid concentrated near the origin, as described in the  Appendix.

2. Force-free magnetic field

In the second configuration, the magnetic field strength is fixed while its direction varies continuously throughout the plasma slab, similar to the configuration in Refs. 56 and 57. This force-free magnetic field profile is complemented with a constant density and temperature, and a linear velocity profile, as defined in Ref. 58, Sec. 14.3.3
ρ 0 ( x ) = ρ c , B 0 ( x ) = sin ( α x ) e ̂ y + cos ( α x ) e ̂ z , T 0 ( x ) = β 0 B 0 2 2 ρ c , v 0 ( x ) = v c x e ̂ y ,
(11)
where ρ c , v c, and α are constant parameters. As the notation suggests, β0 is the equilibrium plasma-β. Note that | B 0 | = 1 such that the dimensionless Alfvén speed c A equals c A = 1 / ρ c. The equilibrium is visualized in Fig. 1(b).

For this equilibrium, we choose k proportional to e ̂ y (and thus parallel to v 0), such that there is a magnetic nullplane at x 0 = 0, where F ( x 0 ) = 0. The perfectly conducting boundaries are placed symmetrically around the nullplane, at x w = ± 0.5. Note that additional magnetic nullplanes appear in this interval for sufficiently large α.

1. ψ-regimes

In the analytic literature surrounding the resistive tearing instability (notably Refs. 7 and 34), the mode is usually classified by the behavior of the normalized, magnetic x-perturbation amplitude B ̂ 1 x, called ψ in Ref. 7 and subsequent literature, in a resistive layer [ x 0 δ , x 0 + δ ] around the magnetic nullplane at x = x0, where F ( x 0 ) = 0. If ψ is approximately constant across this resistive layer, the mode is called a constant-ψ mode. In general, this corresponds to short wavelengths (large k). On the other hand, for longer wavelengths (small k), the variation in ψ throughout the resistive layer is not negligible. In this case, the tearing mode is called a nonconstant-ψ mode. This distinction is important because, analytically, the growth rate of the tearing mode scales differently with resistivity η for the two regimes. In the static (and small v 0) case, the growth rate scales as η 3 / 5 for constant-ψ modes, whereas for nonconstant-ψ modes, the growth rate scales as η 1 / 3.7,34,40

To visualize how ψ changes with k, ψ is shown for the tearing mode of the static Harris sheet from Sec. II A 1 ( B c = 1, aB = 1, v c = 0) for a selection of wavenumbers in Fig. 2(a).

FIG. 2.

(a) ψ for different wavenumbers in a Harris sheet. (b) ψ / ψ and (c) J ̂ z for the k = 0.5 case from panel (a). Δ and δ indicate the numerical matching quantity and resistive layer halfwidth, respectively.

FIG. 2.

(a) ψ for different wavenumbers in a Harris sheet. (b) ψ / ψ and (c) J ̂ z for the k = 0.5 case from panel (a). Δ and δ indicate the numerical matching quantity and resistive layer halfwidth, respectively.

Close modal
Alongside ψ and the resistive layer halfwidth δ, the literature quantifies the matching quantity7,
Δ = [ 1 ψ ψ x ] x = x 0 x = x 0 + .
(12)
In the analytic approach, ψ is obtained by solving the linearized, ideal MHD equations outside of the resistive layer, because resistivity is negligible there, and matching them with the resistive solution inside the layer. However, as an artifact of this approach, the resulting solution has a discontinuity in ψ . To eliminate discontinuities from the calculations, the matching quantity Δ is introduced in the analytic approach. Consequently, this quantity appears in the analytic scaling laws and imposes the condition δ | Δ | < 1 for the constant-ψ approximation to hold.7,34,40
Numerically, ψ has no true discontinuity. Instead, there is a steep but smooth reversal across the magnetic nullplane. Hence, we here define Δ as the difference in ψ / ψ between the extrema on either side of the nullplane, as seen for the static Harris sheet with k = 0.5 in Fig. 2(b). As shown in Ref. 59, defining δ as the distance from the nullplane to the nearest inflexion point of J ̂ z, i.e., J ̂ z ( δ ) = 0, as shown in Fig. 2(c), is consistent with boundary layer theory. Though J ̂ z is symmetric around the nullplane, we take
δ = | x M x m | 2
(13)
with J ̂ z ( x M ) = max ( J ̂ z ) and J ̂ z ( x m ) = min ( J ̂ z ), to limit the impact of the discretized grid and the numerical differentiation of J ̂ z. However, if background flow is included, all perturbations become complex. To maintain consistency with the static configurations, the complex factor is chosen in such a way that Im ( B ̂ 1 x ) is positive and symmetric around the nullplane. For ψ / ψ, this factor is eliminated and the real part is used to define Δ , but this allows us to define δ based on the now odd function Re ( J ̂ z ).

To distinguish between analytic and numerical matching quantities and resistive layer halfwidths, we will refer to them with subscripts A and N, respectively.

2. Shear ratio

Similarly to F in Eq. (10), for the equilibrium flow, we define the angle-modulated Alfvén Mach number (like Ref. 34)
G ( x ) = k · v 0 ( x ) | k | c A .
(14)
Here, c A indicates the dimensionless Alfvén speed c A = | B 0 | / ρ 0. Following Ref. 34, the expression,
R 0 = | G ( x 0 ) F ( x 0 ) | ,
(15)
where the prime denotes the derivative with respect to x, acts as a diagnostic parameter to quantify the relative strength of the flow shear compared to the magnetic shear. We will refer to R0 as the shear ratio.

3. Scaling laws

Now that we have introduced the ψ-regimes and shear ratio R0, we can differentiate between the various conditions under which the growth rate scaling laws were derived as a function of the resistivity η. The scaling laws are summarized in Table I as they were obtained by Chen and Morrison.34 

TABLE I.

Analytic scaling laws of the tearing growth rate with the resistivity η.34 

Constant-ψ Nonconstant-ψ
R 0 1  γ η 3 / 5  γ η 1 / 3 
R 0 1  γ η 1 / 2  γ η 1 / 3 
R 0 > 1  Stabilized  Stabilized 
Constant-ψ Nonconstant-ψ
R 0 1  γ η 3 / 5  γ η 1 / 3 
R 0 1  γ η 1 / 2  γ η 1 / 3 
R 0 > 1  Stabilized  Stabilized 

4. Relative quantities

Throughout this article, we compare the tearing growth rate under the influence of equilibrium flow to the equivalent configuration without flow. In these cases, we opt to use the relative growth rate γ, which we define as
γ = Im ( ω flow ) Im ( ω no flow ) Im ( ω no flow ) .
(16)
Hence, γ ranges from –1, which means the tearing instability is completely stabilized, to + . If γ = 0, background flow does not alter the growth rate. Similarly, we also define the relative numerical matching quantity
Δ rel = Δ N , flow Δ N , n o flow Δ N , n o flow .
(17)

5. Solvers

Throughout this work, two main solvers are used in Legolas, as described in Ref. 53. The first one is the default solver, QR-cholesky, which results in the full spectrum, but is quite slow. The second one is the shift-invert Arnoldi solver, which only computes a selection of modes and is thus faster, but requires a target to converge around. This latter method is preferred for parameter studies once we have an approximation of the growth rate to act as the target.

In this work, we numerically investigate the flow-sheared resistive tearing mode in two different configurations. In Sec. III A, we first present visualizations of the perturbed magnetic field and flow in a Harris sheet (see Sec. II A 1) due to the linear tearing mode, both in the absence and presence of equilibrium flow. Then, we systematically vary the parameters appearing in the shear ratio R0 to evaluate how the growth rate is affected. In Sec. III B, we repeat this parameter study for the force-free magnetic field configuration (see Sec. II A 2), where we also vary the plasma-β.

For the Harris current sheet setup (Sec. II A 1), the effect of shear flow on the resistive tearing mode was probed in Ref. 28 using non-linear, incompressible MHD simulations by computing the reconnection rate for a selection of test cases. Here, we compute the linear growth rate using the compressible equations for various parameter combinations, initially to look at how the growth rate scales with resistivity, and later density. Assuming k = k y e ̂ y, the shear ratio reduces to
R 0 = a B v 0 a v B 0 .
(18)
Hence, we also verify that when the flow shear exceeds the magnetic shear, i.e., R 0 > 1, the tearing instability is fully stabilized,34 by varying the velocity parameters.

1. Linear tearing of the Harris sheet

To begin, we look at the influence of the resistive tearing mode on the magnetic field and flow. To do so, we consider the Harris sheet presented in Sec. II A 1, both with and without the background flow, Eq. (7). In both cases, we set our parameters to k = 0.12 e ̂ y , ρ 0 = 1 , B c = 1, and aB = 1, and use v c = 0.25 and a v = 0.75 when equilibrium flow is included. This is solved on a symmetric, accumulated grid as described in the  Appendix with parameters p 1 = 0.2 , p 2 = 0 , p 3 = 0.01, and p 4 = 5, resulting in 327 grid points. In this configuration, the sheet is expected to tear up into magnetic islands, as represented schematically in Fig. 3.

FIG. 3.

Schematic representation of magnetic field lines altered by the tearing instability.

FIG. 3.

Schematic representation of magnetic field lines altered by the tearing instability.

Close modal

For the flowless case, Figs. 4(a) and 4(b) show Legolas's solutions for the B ̂ x and B ̂ y perturbation amplitudes from Eq. (5), respectively, where the arbitrary factor was chosen in such a way that the maximal density perturbation equals 1%, i.e., max ( ρ 1 ) = 0.01 ρ 0, at t = 0. Similarly, Figs. 4(d) and 4(e) show the amplitudes in the case with background flow. In either case, Bz is not perturbed. Note that in the flowless case, one component is purely real and one is purely imaginary, while the presence of the background shear flow makes both components fully complex.

FIG. 4.

Magnetic field during linear Harris sheet tearing. No flow: (a) B ̂ x perturbation amplitude; (b) B ̂ y perturbation amplitude; (c) magnetic field lines and B 1 y-visualization (color) for max ( ρ 1 ) = 0.01 ρ 0. With flow: (d) B ̂ x perturbation amplitude; (e) B ̂ y perturbation amplitude; (f) magnetic field lines and B 1 y-visualization (color) for max ( ρ 1 ) = 0.01 ρ 0.

FIG. 4.

Magnetic field during linear Harris sheet tearing. No flow: (a) B ̂ x perturbation amplitude; (b) B ̂ y perturbation amplitude; (c) magnetic field lines and B 1 y-visualization (color) for max ( ρ 1 ) = 0.01 ρ 0. With flow: (d) B ̂ x perturbation amplitude; (e) B ̂ y perturbation amplitude; (f) magnetic field lines and B 1 y-visualization (color) for max ( ρ 1 ) = 0.01 ρ 0.

Close modal

Substituting these perturbation amplitudes in Eq. (5) for t = 0 and y ranging from 0.5 λ to 0.5 λ ( 26.18 ), with λ = 2 π / | k | the wavelength, results in the linear perturbation of the magnetic field, presented as field lines and with a color map of the B 1 y component, in Figs. 4(c) and 4(f), for the case without and with flow, respectively. Note that since the amplitudes are normalized and can be multiplied with any complex factor, in principle, the (dimensionless) time t presented here cannot be linked to a physical time without context. Further note that the range of x-coordinates was limited to focus on the behavior near the sheet. As expected from the schematic representation in Fig. 3, the visualization in Fig. 4(c) reveals a magnetic island in the center, with dipping in the previously straight magnetic field lines on either side. If the system were to evolve linearly, however, due to the sharp peaks in the By perturbation amplitude on either side of the magnetic nullplane [see panel (b)], the magnetic field lines inside the island would dip progressively deeper inwards at the nullplane (because they cannot cross), until a magnetic field line meets itself again at the origin, dividing the island into two smaller islands (closed magnetic fieldlines) on either side of the nullplane. Since this is not observed in non-linear simulations, the time when this behavior starts to develop in the linear solution marks an upper bound on the transition time from the linear to the non-linear regime. The inclusion of background shear flow does not dramatically alter this magnetic field structure, aside from introducing a slight shear deformation of the field lines across the nullplane. For the chosen parameters, this cannot be seen clearly at t = 0, shown in Fig. 4(f), but the effect becomes more apparent as the perturbation grows.

Similarly, Fig. 5 presents the v ̂ x and v ̂ y perturbation amplitudes from Eq. (5) for the flowless case in panels (a) and (b) and for the case with flow in panels (d) and (e). Again, the z-component, vz, is not perturbed, and the flowless case has a purely real and purely imaginary v 1-component. Identically to the magnetic field perturbation, the velocity perturbation amplitude also becomes complex when background flow is added. Contrary to the magnetic field perturbation amplitude though, where the newly introduced real/imaginary parts are an order of magnitude smaller than the original part, the real and imaginary parts of the velocity perturbation amplitudes are of a comparable order of magnitude.

FIG. 5.

Flow during linear Harris sheet tearing. No flow: (a) v ̂ x perturbation amplitude; (b) v ̂ y perturbation amplitude; (c) magnetic field lines, stream vectors for | v | < 0.5 × 10 3, and v 1 y-visualization (color) for max ( ρ 1 ) = 0.01 ρ 0. With flow: (d) v ̂ x perturbation amplitude; (e) v ̂ y perturbation amplitude; (f) magnetic field lines, (fixed-length) stream vectors for | v | < 0.24, and v 1 y-visualization (color) for max ( ρ 1 ) = 0.01 ρ 0.

FIG. 5.

Flow during linear Harris sheet tearing. No flow: (a) v ̂ x perturbation amplitude; (b) v ̂ y perturbation amplitude; (c) magnetic field lines, stream vectors for | v | < 0.5 × 10 3, and v 1 y-visualization (color) for max ( ρ 1 ) = 0.01 ρ 0. With flow: (d) v ̂ x perturbation amplitude; (e) v ̂ y perturbation amplitude; (f) magnetic field lines, (fixed-length) stream vectors for | v | < 0.24, and v 1 y-visualization (color) for max ( ρ 1 ) = 0.01 ρ 0.

Close modal

Once more, Figs. 5(c) and 5(f) display a 2D visualization of the magnetic field lines for the same y-interval and time as Figs. 4(c) and 4(f). Here, however, the panels are colored by the magnitude of the y-component of the velocity perturbation, v 1 y. In addition, the flow patterns are further highlighted with stream vectors up to a certain magnitude | v |, to show variations in the regions of smaller speed. From the color map in Fig. 5(c), it is immediately clear that the velocity in the flowless case is much higher at the magnetic nullplane than further away from it. In fact, the velocity indicates an inflow toward the center of the magnetic island along the nullplane. Away from the sheet's center, the stream vectors cross the magnetic field lines, moving outward from the island's center, though with a much smaller speed than the inflow speed. The addition of background flow in Fig. 5(f) changes the picture significantly. While the flow speed is still much higher in the center of the sheet, the velocity now roughly follows the magnetic field at the island edges, resulting in a flow within the island. Again, the islands are too small at t = 0 to clearly distinguish the flow pattern inside, but it becomes apparent as the perturbation grows.

2. Matching quantity and resistive layer halfwidth

Before turning to the growth rate scaling, we evaluate the numerical matching quantity and resistive layer halfwidth, how they compare to analytic values, and how they are affected by introducing flow. Once again, we first consider the static Harris sheet with k = 0.5 e ̂ y , ρ 0 = 1 , B c = 1, and aB = 1 ( v 0 = 0) in this section. For this static Harris sheet configuration, the analytic matching quantity is given exactly by7,40
Δ A = 2 a B ( 1 k a B k a B ) .
(19)
Hence, for the chosen parameters, we find Δ A = 3. Additionally, for constant-ψ modes, the analytic resistive layer halfwidth is given by40 
δ A η 2 / 5 ( Δ A ) 1 / 5 ( k B 0 ) 2 / 5 ,
(20)
which evaluates to δ A 4.13 × 10 2 here, where B 0 ( x ) is evaluated at the nullplane (x = 0). Note that this satisfies the constant-ψ condition δ A | Δ A | 0.124 < 1.

Numerically (grid parameters p 1 = 0.2 , p 2 = 0 , p 3 = 0.01 , p 4 = 5), as explained in Sec. II B, we find Δ N 2.56 and δ N 3.5 × 10 2 (note the use of due to the x-discretization in Legolas). Again, the constant-ψ condition is satisfied, δ N | Δ N | 0.090 < 1. Hence, though the analytic and numerical values for Δ and δ are not in perfect agreement, they are of the same order and reasonably close. Additionally, though the wavenumber k = 0.5 is neither small nor large considering that the configuration is tearing unstable for k 1, the constant-ψ approximation holds, and thus, we expect to recover a growth rate scaling proportional to η 3 / 5 in the static case.

Analytically, Hofmann30 showed that if the velocity v 0 is proportional to the Alfvén velocity (and thus the magnetic field B 0) everywhere, the matching quantity is unaltered from the static case. For the Harris sheet with velocity profile Eq. (7), v 0 B 0 is only true if av = aB. However, as shown in Fig. 6(a), the ψ / ψ ratio is not identical for the static and stationary ( v c = 0.5, av = 1) case. This is further highlighted by the difference between both cases, shown in Fig. 6(b). Consequently, neither is the numerical matching quantity. Now, we find a numerical matching quantity of Δ N 3.77 and a resistive layer halfwidth of δ N 5.5 × 10 2. Hence, both the numerical matching quantity and resistive layer halfwidth depend on the velocity, even if v 0 B 0. However, δ N | Δ N | 0.21 < 1 still holds. Since ψ / ψ is identical in our incompressible approximation, which also eliminates the Joule heating,15 this deviation from Hofmann's result is presumably due to the use of a constant resistivity rather than a convectively perturbed resistivity.

FIG. 6.

(a) ψ / ψ for k = 0.5 in a static and stationary Harris sheet. (b) Stationary minus static difference in ψ / ψ from panel (a). (c) Static Δ and δ as a function of k = k y e ̂ y. (d) Static Δ and δ as a function of aB. (e) Δ N and (f) δ N for k = 0.5 in a Harris sheet with flow of various transition halfwidths.

FIG. 6.

(a) ψ / ψ for k = 0.5 in a static and stationary Harris sheet. (b) Stationary minus static difference in ψ / ψ from panel (a). (c) Static Δ and δ as a function of k = k y e ̂ y. (d) Static Δ and δ as a function of aB. (e) Δ N and (f) δ N for k = 0.5 in a Harris sheet with flow of various transition halfwidths.

Close modal

The discrepancy between analytic and numerical Δ and δ, as well as the difference between the static and stationary Harris sheet, raises the question how the numerical matching quantity and resistive layer halfwidth depend on the various parameters. In Fig. 6(c), both the analytic and numerical matching quantities are shown side by side as a function of k for the static Harris sheet. For large wavenumbers, i.e., in the constant-ψ regime, the analytic and numerical Δ are in good agreement, but for small wavenumbers, the analytic matching quantity diverges to + , whereas Δ N is observed to decrease again. The maximal value of Δ N in our set of discrete k-values was achieved for k = 0.09, with δ N | Δ N | 0.62, indicating a strong numerical deviation from the analytic results in the traditional nonconstant-ψ regime. Simultaneously, no maximum is observed in δ N, which continues to increase as k decreases. Note that the steplike behavior of δ N is due to the discretization of the x-coordinate in Legolas [also in panel (f)]. This deviation from analytic results in the nonconstant-ψ regime is not surprising. As k decreases, ψ ( 0 ) approaches zero and Δ N approaches an x 1 scaling. Due to the definition of Δ N, the evaluation of ψ / ψ occurs near the edge of the resistive layer, and thus, Δ N δ N 1. Hence, Δ N decreases as δ N increases for k 0. Consequently, Δ N is not a proper numerical equivalent of Δ A in the nonconstant-ψ regime, and thus, δ N | Δ N | cannot be used to distinguish between constant-ψ and nonconstant-ψ modes.

Similarly to Fig. 6(c), Fig. 6(d) presents the dependence of Δ N and δ N on the magnetic reversal halfwidth aB. Here, the numerical matching quantity appears to follow the analytic scaling reasonably well, though it is consistently larger than Δ A for our choice of k = 0.5. The numerical resistive layer halfwidth δ N, on the other hand, is observed to be constant as a function of aB, which is surprising considering that Furth et al.7 derived an a B 1 dependence for the width of the region of discontinuity.

Finally, again adding the velocity profile Eq. (7), Figs. 6(e) and 6(f) show the dependence of Δ N and δ N on the transition halfwidth av of the velocity profile for various flow speeds v c. As expected, the effect of flow on the matching quantity increases with the flow speed v c. For sufficiently high v c and small av, the numerical matching quantity is observed to increase initially with av, whereas the resistive layer halfwidth decreases. As av increases further, Δ N reaches a maximum and decreases again steadily. For the resistive layer halfwidth, we observe a monotone decrease with av, which is more pronounced for larger v c.

3. Growth rate scaling with resistivity

In the seminal work by Furth, Killeen, and Rosenbluth,7 the authors derive power laws for the scaling of the incompressible tearing growth rate as a function of the resistivity η, for small values of η. Here, we introduce compressibility and consider a wide range of resistivity values. However, contrary to their work, we assume a constant resistivity without convective perturbations, and include Joule heating. From their derivations, they conclude that the growth rate scales as a power law in η, Im ( ω ) η p, and that p depends on whether or not ψ is approximately constant across the magnetic nullplane. This ψ-classification remains equally important when shear flow is added.34 

In Sec. III A 2, we have established that the tearing mode of the Harris sheet from Sec. II A 1 falls into the constant-ψ regime for the parameters k = 0.5 e ̂ y , ρ 0 = 1 , B c = 1, aB = 1, av = 1, and v c sufficiently small. Now, we vary the resistivity η to compare to the literature's corresponding scaling laws. In Fig. 7(a), the flowless case is compared to the case with flow profile Eq. (7) for a selection of R0 values, which are expected to scale differently based on the analytic power laws presented in Ref. 34. (Note that the no flow case is not clearly visible since it almost coincides with the R 0 = 0.1 case.) All Legolas runs were performed using grid parameters p 1 = 0.2 , p 2 = 0 , p 3 = 0.01, and p 4 = 5 (327 grid points). The desired value of R0 was obtained by setting v c = a v B c R 0 / a B.

FIG. 7.

Resistive tearing growth rate of a Harris sheet, Eq. (6), as a function of (a) the resistivity η; (b) the density ρ0; for parameters k = 0.5 e ̂ y , ρ 0 = 1 , B c = 1, aB = 1, and (a) av = 1; (b) v c = 0.1, if flow was included. In (b), the density where the maximal velocity of the equilibrium configuration coincides with the Alfvén speed is indicated with a dotted line. (c) Growth rate scaling with ρ0 for constant Alfvén speed ( B c = ρ 0) and p c = 10 6.

FIG. 7.

Resistive tearing growth rate of a Harris sheet, Eq. (6), as a function of (a) the resistivity η; (b) the density ρ0; for parameters k = 0.5 e ̂ y , ρ 0 = 1 , B c = 1, aB = 1, and (a) av = 1; (b) v c = 0.1, if flow was included. In (b), the density where the maximal velocity of the equilibrium configuration coincides with the Alfvén speed is indicated with a dotted line. (c) Growth rate scaling with ρ0 for constant Alfvén speed ( B c = ρ 0) and p c = 10 6.

Close modal

Since the analytic scaling laws are derived under the assumption that there is a clear temporal separation between the resistive diffusion time τ R = η 1 (in our dimensionless units) and the Alfvén time τ A = a B ρ 0 B c 1,7,59 comparison to these power laws is only meaningful for η < 10 2 in this case. Therefore, to be on the safe side, all fitted power laws from Table I were limited to η < 10 3 in Fig. 7(a). In accordance with the work of Chen and Morrison,34 the scaling of the static Harris sheet's growth rate is found to be close to η 3 / 5, as expected in the constant-ψ regime. For R 0 = 0.1 1, the growth rate is slightly larger than the static case, but the same scaling law seems to hold. As R0 increases to R 0 = 0.5, the scaling law changes to η 1 / 2, which is also in line with the literature for R 0 1. However, for R 0 = 0.8, none of the analytic power laws are a good fit. At small η, this case is observed to scale as η 0.44. Note that this power law lies between the constant-ψ η 1 / 2 and nonconstant-ψ η 1 / 3 scaling laws, though the product δ N | Δ N | remains significantly smaller than 1, once again confirming it is not a good metric to classify ψ. In addition, the growth rate deviates strongly from a power law well before η reaches 10 2.

As expected, aside from the growth rate dropoff in the R 0 = 0.8 case, the other cases also deviate from the power laws as η crosses the 10 2 threshold. Note that growth rates above 10 2 should be interpreted with care, since diffusion of the equilibrium was neglected and might affect the dynamics significantly.

4. Density variation

To demonstrate that the growth rate does not care about the specific density value, but only about the Alfvén speed, first consider the static Harris sheet with k = 0.5 e ̂ y , B c = 1, aB = 1, and varying density (grid parameters: p 1 = 0.2 , p 2 = 0 , p 3 = 0.01 , p 4 = 5). This is shown in Fig. 7(b), where a fit to the data shows that the growth rate scales approximately as ρ 0 0.232. Comparing this to the analytic (incompressible) growth rate scaling ρ 0 1 / 5 [see, e.g., Ref. 7, Eq. (54)], while reasonably close, the deviation is significant. This may be due to compressibility, especially considering that the plasma-β goes to at the nullplane in this configuration.

Since both the Alfvén speed c A and sound speed c s are proportional to ρ 1 / 2 (due to our choice of T 0 ρ 0 1 profile), the question becomes how the growth rate scales with either speed. To answer this, we now consider a second static Harris sheet with the same parameters as before, except we set B c = ρ 0 and p c = 10 6. These choices ensure that the Alfvén speed is constant when varying the density while T 0 p c / ρ 0, and thus, c s ρ 0 1 / 2. In this case, we find a constant growth rate, shown in Fig. 7(c). Hence, we conclude that the growth rate does not depend on the specific density or sound speed, but scales with the Alfvén speed.

The effect of the inclusion of a non-zero flow of fixed size v c = 0.1 (and various transition halfwidths) depends on the case. In the latter case, the v c / c A ratio is constant, and the growth rate modification does not vary with ρ0. In the first case however, the flow modifies the tearing growth rate variation as a function of the density ρ0, or better, because the ratio v c / c A changes. In fact, the flow introduces a critical density above which the tearing mode is fully damped as v c / c A 1. This is highlighted in Fig. 7(b) by the dotted line representing where v c = c A. However, from the figure it is clear that the tearing mode is not always damped exactly when v c reaches c A, but depends on the flow transition halfwidth av. Therefore, the Alfvén speed does not necessarily act as a transition value in the equilibrium speed with respect to tearing suppression. Though the tearing instability may vanish at a density lower than where the equilibrium speed equals the Alfvén speed, it appears the Alfvén speed still imposes an upper limit on the density above which the tearing mode vanishes, from the sharp dropoff there in the a v = 1.25 case.

5. Velocity variation

Since both parameters of the velocity profile (the maximal speed v c and halfwidth av) appear to play an important role, we vary both parameters simultaneously to identify the regions of stabilization and further destabilization of the tearing instability. Here, the maximal speed v c is kept sub-Alfvénic ( v c < 1) because the flow-induced KHI dominates in the super-Alfvénic regime.30 The result is shown in Fig. 8(a) for fixed parameters k = 0.5 e ̂ y , ρ 0 = 1 , B 0 = 1, and aB = 1 and was obtained with the shift-invert method on an accumulated grid with parameters p 1 = 0.2 , p 2 = 0 , p 3 = 0.01, and p 4 = 5 (327 grid points). Since this method requires an initial guess, the first guess was obtained from a run with the QR-cholesky solver for v c = 10 2 and av = 1, which was used to compute the growth rate for all values of av and v c = 10 2 with the shift-invert method. This array of growth rates was subsequently used as the initial guesses for the next value of v c, and so on. Visually speaking, each value (except those in the left column) in Fig. 8(a) was computed via shift-invert by providing the value to its immediate left as the initial guess. Note that while the tearing mode is fully damped in the top right corner of panel (a), the system is still unstable because the KHI appears here before v c reaches the Alfvén speed. This was again checked with the QR-cholesky solver in this region of the parameter space. (For the few parameter combinations where the shift-invert method failed to converge, the growth rate was calculated using the QR-cholesky solver.) Since the analytic scaling laws typically feature a scaling with the matching quantity, Fig. 8(b) shows the relative numerical matching quantity Δ rel for the parameter combinations where the tearing growth rate exceeds 10 4, and gray elsewhere.

FIG. 8.

(a) Relative tearing growth rate γ and (b) relative numerical matching quantity Δ rel with respect to the static case for combinations of the maximal speed v c and flow transition halfwidth av. The yellow dotted line indicates the magnetic field transition halfwidth aB and the orange dashed line represents R 0 = 1.

FIG. 8.

(a) Relative tearing growth rate γ and (b) relative numerical matching quantity Δ rel with respect to the static case for combinations of the maximal speed v c and flow transition halfwidth av. The yellow dotted line indicates the magnetic field transition halfwidth aB and the orange dashed line represents R 0 = 1.

Close modal

In this figure, a few things stand out. First of all, the resistive tearing instability is clearly stabilized at the R 0 = 1 (orange dashed) line for all considered parameter combinations, but is already fully suppressed before this line is reached. Since the analytic growth rate scales with Δ A34 and Δ A is only identical to the static case if v 0 B 0,30 the instability is only expected to vanish exactly at R 0 = 1 if av = aB, so this is not surprising. However, even for av = aB (yellow dotted line), the instability is already fully suppressed before R 0 = 1 is reached. This may be due to the earlier observation (see Sec. III A 2) that Δ N does differ between static and stationary Harris sheets here, even if v 0 B 0.

Second, contrary to the non-linear observation in Ref. 28 that there exists a single critical av value 0.35 where the transition from stabilizing to destabilizing occurs, we here observe that this critical av depends on the maximal speed v c. Furthermore, their critical value lies in our stabilizing region of the parameter regime across all velocities. Note though that the simulations in Ref. 28 are incompressible, include a non-zero viscosity, and lack Joule heating, any of which may affect this result.

Finally, comparing Figs. 8(a) and 8(b), there is a discrepancy between the variations in the growth rate and the numerical matching quantity. At smaller av, the growth rate is slightly stabilized while the matching quantity is slightly increased. Therefore, since we only varied the flow parameters, the influence of shear flow on the growth rate is not fully encapsulated in the numerical matching quantity, and the growth rate does not simply scale with Δ N.

Now, we turn to the setup from Sec. II A 2, where a magnetic field of fixed magnitude varies its direction periodically throughout the plasma slab. In this particular case, the shear ratio becomes
R 0 = | G ( 0 ) F ( 0 ) | = v c ρ c α .
(21)

While there are three parameters in this expression, our parametric study will only focus on the variation of the equilibrium density ρ c and velocity coefficient v c, for reasons which will become clear after a demonstration of the role of α.

1. Multiple tearing modes

Since the parameter α regulates how fast the magnetic field's direction varies along x, it also determines how many magnetic nullplanes the system has in a certain x-interval for a given wave vector. Consequently, the number of tearing modes supported by the configuration depends on α. Additionally, if the equilibrium velocity is described by an odd function, like the linear profile in Eqs. (11), the spectrum is symmetric with respect to the imaginary axis. This is illustrated in Figs. 9(a)–9(d), where we varied the parameter α for parameters ρ c = 1 , β 0 = 0.15 , v c = 0.15, and k = 1.5 e ̂ y at 251 grid points. For (a) α = π / 2, the magnetic shear is insufficient to induce a tearing instability. Increasing α without introducing an additional nullplane results in one non-propagating tearing mode (i.e., purely imaginary), visualized in (b) for α = 4.738 84 (slightly more than 3 π / 2). In the presence of three nullplanes for α = 5 π / 2, (c) shows a pair of forward–backward-propagating instabilities and one non-propagating one. Finally, (d) contains only two pairs of forward–backward-propagating tearing pairs for α = 4.1 π, despite the presence of five nullplanes in the domain. If the equilibrium flow is removed, all tearing modes become non-propagating, i.e., purely imaginary, as demonstrated in Fig. 10(a) for the case with α = 4.1 π.

FIG. 9.

[(a)–(d)] Parts of the spectra of a plasma slab with force-free magnetic field, Eqs. (11), for ρ c = 1 , β 0 = 0.15 , v c = 0.15, and k = 1.5 e ̂ y. The angular parameter α in the magnetic field profile determines the number of magnetic nullplanes and unstable modes, and takes a different value in each panel: (a) α = π / 2 (b) α = 4.738 84 (c) α = 5 π / 2, and (d) α = 4.1 π. [(e)–(j)] B ̂ 1-eigenfunctions of panel (d)'s [left: (e), (g), and (i)] dominant, forward-traveling tearing mode; [right: (f), (h), and (j)] less unstable, forward-traveling tearing mode. Magnetic nullplanes are indicated by dash-dotted lines, with blue lines marking where tearing occurs (lighter blue indicates that multiplication with a complex factor is required to highlight tearing behavior there).

FIG. 9.

[(a)–(d)] Parts of the spectra of a plasma slab with force-free magnetic field, Eqs. (11), for ρ c = 1 , β 0 = 0.15 , v c = 0.15, and k = 1.5 e ̂ y. The angular parameter α in the magnetic field profile determines the number of magnetic nullplanes and unstable modes, and takes a different value in each panel: (a) α = π / 2 (b) α = 4.738 84 (c) α = 5 π / 2, and (d) α = 4.1 π. [(e)–(j)] B ̂ 1-eigenfunctions of panel (d)'s [left: (e), (g), and (i)] dominant, forward-traveling tearing mode; [right: (f), (h), and (j)] less unstable, forward-traveling tearing mode. Magnetic nullplanes are indicated by dash-dotted lines, with blue lines marking where tearing occurs (lighter blue indicates that multiplication with a complex factor is required to highlight tearing behavior there).

Close modal
FIG. 10.

(a) Part of the spectrum of a plasma slab with force-free magnetic field, Eqs. (11), for ρ c = 1 , β 0 = 0.15 , α = 4.1 π, and k = 1.5 e ̂ y, without background flow ( v c = 0). [(b)–(e)] B ̂ 1-eigenfunctions for (b) ω 0.3100 i; (c) ω 0.2961 i; (d) ω 0.2553 i; (e) ω 0.1654 i. Magnetic nullplanes are indicated by dash-dotted lines, with darker lines marking where tearing occurs in that panel.

FIG. 10.

(a) Part of the spectrum of a plasma slab with force-free magnetic field, Eqs. (11), for ρ c = 1 , β 0 = 0.15 , α = 4.1 π, and k = 1.5 e ̂ y, without background flow ( v c = 0). [(b)–(e)] B ̂ 1-eigenfunctions for (b) ω 0.3100 i; (c) ω 0.2961 i; (d) ω 0.2553 i; (e) ω 0.1654 i. Magnetic nullplanes are indicated by dash-dotted lines, with darker lines marking where tearing occurs in that panel.

Close modal

For the flowless case in Fig. 10(a), the magnetic field perturbation amplitudes of the unstable modes are shown in panels (b) through (e). All positions of magnetic nullplanes are marked with a dash-dotted line. At the darker-colored nullplanes, we observe a dip in B ̂ 1 x and a sharp transition in B ̂ 1 y, similar to Figs. 4(a) and 4(b), indicative of tearing at this nullplane. Note that the dominant instability features tearing behavior at the three central nullplanes, whereas the secondary instability only tears up the central nullplane, while the tertiary mode is tearing up the same nullplanes as the dominant instability except for the central nullplane. Surprisingly, the least unstable mode does not show any signs of tearing. Additionally, the outer nullplanes do not show strong signs of tearing, presumably due to their proximity to the perfectly conducting boundaries, which exert a stabilizing influence.

With the addition of flow, however, the situation changes. Now, the spectrum in Fig. 9(d) no longer has a single dominant instability, but a dominant pair and less unstable pair of forward–backward-propagating instabilities. Furthermore, their perturbation amplitudes are now fully complex. For the dominant, forward-propagating [ Re ( ω ) > 0] instability, the B ̂ 1 perturbation amplitudes are shown in Figs. 9(e), 9(g), and 9(i), and similarly, those of the less unstable, forward-propagating instability are shown in Figs. 9(f), 9(h), and 9(j). The backward-propagating counterparts are not shown, but are the mirror image with respect to x = 0 of the forward-propagating perturbation amplitudes. The magnetic nullplanes are again indicated by dash-dotted lines.

Contrary to the static cases, all four modes appear to tear all three central nullplanes, colored in blue, though lighter blue lines indicate multiplication with a complex factor is needed to highlight the tearing behavior. Additionally, the dominant, forward-propagating mode's largest amplitudes variations (i.e., strongest tearing) occur at the positive intermediate and central nullplanes, with smaller variations at the negative intermediate nullplane. The less unstable, forward-propagating mode, on the other hand, has its largest amplitudes variations at the central nullplane, with smaller variations at both intermediate nullplanes. Hence, by introducing flow, the tearing behavior of all modes is altered significantly.

From now on, the value of α is set to α = 4.738 84 (i.e., the value used in Fig. 9(b) and Ref. 58). This ensures that the magnetic field makes between one-half and a full rotation in the considered domain, x [ 0.5 , 0.5 ], resulting in a single nullplane at x = 0 and a single tearing mode. Since there is only one nullplane, no multitearing occurs.60,61 A detailed look at multitearing is beyond the scope of this paper.

2. Matching quantity and resistive layer halfwidth

Like in the previous case of the Harris sheet, we evaluate the matching quantity and resistive layer halfwidth for the force-free magnetic field configuration before moving on to the scaling of the growth rate (for now, we set aside our conclusion from that section regarding Δ N). Again, we choose the parameters k = 1.5 e ̂ y , ρ c = 1 , α = 4.738 84, and β 0 = 0.15. When velocity is included, we set v c = 0.25. All runs were performed at 251 grid points. Numerically, the static case yields Δ N 5.67 and δ N 1.4 × 10 2, whereas the stationary case yields Δ N 5.84 and δ N 1.4 × 10 2. Hence, the constant-ψ condition is satisfied in both cases, with δ N | Δ N | 0.079 < 1 and δ N | Δ N | 0.082 < 1 for the static and stationary case, respectively. Therefore, we again expect a growth rate scaling proportional to η 3 / 5 in Sec. III B 3.

Again, we look at the influence of the flow profile on the ψ / ψ ratio. In Fig. 11(a), the ψ / ψ ratio is shown for the static and stationary case. The difference between both cases is shown in Fig. 11(b). Clearly, the ψ / ψ ratio does not change wildly with the addition of this flow profile. However, the difference between both cases is not negligible either, explaining the difference in numerical matching quantity calculated above.

FIG. 11.

(a) ψ / ψ for k = 1.5 in the static and stationary force-free B 0-configuration. (b) Stationary minus static difference in ψ / ψ from panel (a). (c) Static δ N , Δ N, and tearing growth rate as a function of k = k y e ̂ y. (d) Stationary δ N , Δ N, and tearing growth rate as a function of v c.

FIG. 11.

(a) ψ / ψ for k = 1.5 in the static and stationary force-free B 0-configuration. (b) Stationary minus static difference in ψ / ψ from panel (a). (c) Static δ N , Δ N, and tearing growth rate as a function of k = k y e ̂ y. (d) Stationary δ N , Δ N, and tearing growth rate as a function of v c.

Close modal

Now the question arises how strongly the numerical matching quantity is impacted by a variation in wavenumber or speed, and whether the growth rate scales with Δ N. In Fig. 11(c), the numerical matching quantity and resistive layer halfwidth are shown as a function of k for the static case. Similarly to the Harris sheet, the numerical matching quantity initially increases with k before decreasing again. Interestingly though, the resistive layer halfwidth is decreasing as k increases, but the product δ N | Δ N | never exceeds 0.20, which would imply that the constant-ψ approximation is valid across all k-values. However, as we concluded in Sec. III A 2, the definition of Δ N makes it unfit to compare to Δ A in the nonconstant-ψ regime, and ensures that δ N | Δ N | < 1 is an insufficient criterion to classify a perturbation as a constant-ψ mode. Additionally, though the growth rate variation with the wavenumber appears to follow a similar trend as the numerical matching quantity, it appears they are not directly proportional. This is further highlighted by Fig. 11(d), where the numerical matching quantity and resistive layer halfwidth are shown as a function of v c for the stationary case. Here, the numerical matching quantity is observed to increase with v c initially, whereas the growth rate starts to decline sooner as v c increases. Here, the resistive layer halfwidth is observed to be more or less constant. Consequently, we henceforth abandon the study of the numerical matching quantity and its relation to the growth rate for this configuration, and solely focus on the modification of the growth rate by the flow profile.

3. Growth rate scaling with resistivity

Once again, we now turn to a comparison of the growth rate scaling with analytic predictions. The parameters are the same as in Sec. III B 2, except that β takes on various values and v c = 0.15 when flow is included. As demonstrated, for this equilibrium we are dealing with a constant-ψ mode and thus expect a scaling of Im ( ω ) η 3 / 5 in the absence of flow. For the case with flow, we expect an identical scaling since R 0 0.03 1.34 This scaling is expected to hold for a / η > 10 2, since the boundary layer approach used in analytic works is no longer appropriate as the resistivity increases,59 and thus, substituting a quarter period for the transition halfwidth a = π / 2 α, for η < 3.3 × 10 3. Hence, all comparisons were limited to η < 10 3. Indeed, initially, the growth rate scales as η 3 / 5, as shown in Fig. 12(a) for 301 grid points, independently of β. However, as η approaches 10 3, the growth rate starts to deviate from this scaling, with smaller β-values deviating sooner, and eventually decreases again. Adding velocity to this variation in resistivity steepens the growth rate dropoff for lower β-values, as evidenced by Fig. 12(b), going as far as eliminating the instability entirely. Of course, care should again be taken in the interpretation of growth rates above η 10 2, since diffusion of the equilibrium was neglected.

FIG. 12.

Resistive tearing growth rate of a force-free magnetic field, Eqs. (11), for ρ c = 1 and k = 1.5 e ̂ y, as a function of η for various plasma-β (a) without flow and (b) for v c = 0.15. Both panels include three power-law fits ( η 3 / 5 , η 1 / 2, and η 1 / 3) for the case with β = 0.1.

FIG. 12.

Resistive tearing growth rate of a force-free magnetic field, Eqs. (11), for ρ c = 1 and k = 1.5 e ̂ y, as a function of η for various plasma-β (a) without flow and (b) for v c = 0.15. Both panels include three power-law fits ( η 3 / 5 , η 1 / 2, and η 1 / 3) for the case with β = 0.1.

Close modal

4. Density variation

Adopting the same approach as for the Harris sheet configuration, we here show that the growth rate of the force-free magnetic field configuration also varies with the Alfvén speed. To do so, we first consider the flowless case with parameters k = 1.5 e ̂ y and α = 4.738 84, for various values of β and varying the density ρ0 (301 grid points). As can be seen in Fig. 13(a), a similar growth rate scaling, ρ 0 0.254, to the Harris sheet is recovered for high plasma-β. For this configuration, the deviation from the analytic ρ 0 1 / 5 scaling is even greater, though an explanation might be that this is due to the incompressible approximation breaking down at high β. However, as β decreases, the growth rate deviates even more from a power law, especially for larger densities (smaller Alfvén speeds). This may be due to the inclusion of Joule heating, whose contribution to the energy equation is more significant for small β ( T 0 β).

FIG. 13.

Resistive tearing growth rate of a force-free magnetic field, Eqs. (11), for k = 1.5 e ̂ y, as a function of ρ0 for various plasma-β (a) without flow and (b) for v c = 0.15. The dotted vertical line indicates where the maximal equilibrium speed equals the Alfvén speed. (c) and (d) show the growth rate for the same parameters as in (a) and (b), respectively, but with B 0 scaled with ρ 0 such that the Alfvén speed is constant.

FIG. 13.

Resistive tearing growth rate of a force-free magnetic field, Eqs. (11), for k = 1.5 e ̂ y, as a function of ρ0 for various plasma-β (a) without flow and (b) for v c = 0.15. The dotted vertical line indicates where the maximal equilibrium speed equals the Alfvén speed. (c) and (d) show the growth rate for the same parameters as in (a) and (b), respectively, but with B 0 scaled with ρ 0 such that the Alfvén speed is constant.

Close modal

Now adding a velocity with v c = 0.15, behavior similar to the Harris sheet is observed in Fig. 13(b). Again, deviations from the flowless scaling are observed as the density approaches the critical value where the maximal flow in the domain equals the Alfvén speed (indicated with a dotted line). Above this threshold, the tearing mode is heavily damped and the growth rate goes to zero.

To show once again that this scaling is due to the change in Alfvén velocity, consider the above case but with B 0 the profile in Eqs. (11) multiplied with ρ 0, such that the Alfvén speed is constant. The resulting growth rates are shown in Figs. 13(c) and 13(d), for the static and stationary cases, respectively. In both cases, the growth rate is constant for all densities, with higher β resulting in a higher growth rate. Hence, the growth rate does not depend on the specific density, but on the Alfvén speed. Since β appears to affect the growth rate, its role is further investigated in conjunction with the velocity.

5. Velocity variation

Despite the simple velocity profile, the ( v c , β )-parameter space reveals surprising complexity. After assuming a constant density ρ 0 = 1, the present parametric survey varied v c and β simultaneously for wave vector k = 1.5 e ̂ y, where the v c parameter was limited to the interval [ 0 , 2 ], such that the equilibrium velocity remains sub-Alfvénic ( | v 0 | c A) on the entire domain, and β-values from 10 2 to 102 were studied. The results are shown in Fig. 14, where all runs in panel (a) were performed at 301 grid points, whereas 201 grid points were used in panel (b).

FIG. 14.

Resistive tearing mode growth rate of a force-free magnetic field, Eqs. (11), for ρ c = 1. (a) Absolute growth rate as a function of v c for various plasma-β values. (b) Relative growth rate for varying v c and β with respect to the flowless growth rate.

FIG. 14.

Resistive tearing mode growth rate of a force-free magnetic field, Eqs. (11), for ρ c = 1. (a) Absolute growth rate as a function of v c for various plasma-β values. (b) Relative growth rate for varying v c and β with respect to the flowless growth rate.

Close modal

As pointed out in Ref. 30, the introduction of flow in a system that is unstable to the resistive tearing mode can either stabilize or further destabilize the plasma. This is also immediately clear from Fig. 14(b), where blue indicates a stabilized system and red a strong increase in tearing growth rate. While the plasma is mostly destabilized further by the presence of flow for large β, as clearly evidenced by Fig. 14, the destabilizing effect does not scale monotonically with the velocity coefficient. Rather, the maximal destabilization appears at some intermediate value between small speeds and the Alfvén speed. For small to intermediate β ( 1), on the other hand, both stabilizing and destabilizing influences are observed in significant fractions of the velocity space, with the strongest stabilizing effect occurring at the Alfvén speed. Additionally, more than one stabilizing–destabilizing transition is observed in panel (b) along the speed axis for small β ( 1).

This dependence on the plasma-β, which already appeared in a less pronounced way in Sec. III B 3, is not surprising. The ion sound Larmor radius ρ s is known to affect the growth rate,59,62 and β enters in the radius through our definition of T0. Nevertheless, any dependencies on ρ s were derived using an isothermal closure, contrary to the equations used here. Since it is clear that the plasma-β affects the role of flow significantly, a more in-depth analysis of the influence of β on tearing modes in general offers interesting perspectives for future work.

In this work, we studied the linear regime of the resistive tearing mode in a compressible plasma with Joule heating in two different configurations featuring shear flow: a Harris current sheet and a force-free magnetic field varying its direction periodically throughout a plasma slab.

First, we visualized the magnetic field lines and flow patterns of the linearly perturbed Harris sheet in 2D, both in the absence and presence of a background flow. In either case, the magnetic field lines were pinched together periodically along the sheet to reconnect and form magnetic islands, as observed in simulations. If this perturbation is allowed to evolve linearly, the formation of two smaller sub-islands is observed inside the island. Since this behavior does not occur in non-linear simulations, the time where this formation is initiated in the linear evolution places an upper limit on the transition time from the linear to the non-linear regime.

For the velocity, the largest perturbation occurs at the magnetic nullplane for both cases with and without equilibrium flow, with plasma leaving the pinched regions and streaming toward the magnetic islands. The plasma further away from the nullplane has an almost negligible velocity compared to the plasma at the nullplane. However, for the flowless equilibrium, the plasma was observed to move away from the magnetic islands here, whereas the plasma rotates inside the inner islands if a background flow is included.

Next, we introduced a numerical equivalent of the matching quantity Δ and resistive layer halfwidth δ. However, for both configurations, the numerical quantities deviated from the analytic predictions. In the nonconstant-ψ regime in particular, the numerical matching quantity fails to capture the behavior of the analytic matching quantity, due to the lack of sharp transition in ψ / ψ near the nullplane for wider resistive layers. In addition, no clear scaling of the growth rate with the numerical matching quantity was observed. Therefore, we advocate against the use of the matching quantity as the sole quantifier of flow's influence on the tearing growth rate.

As a consequence of this deviation from analytic predictions in the matching quantity in the nonconstant-ψ regime, the product δ N | Δ N | is not a good indicator of the validity of the constant-ψ approximation. Nevertheless, the literature's scaling laws hold in the constant-ψ regime, as far as we have observed. Indeed, the growth rates of both configurations were found to scale as η 3 / 5 in the absence of flow, as expected for a constant-ψ mode. Subsequently, when flow was introduced, the transition to the constant-ψ η 1 / 2 scaling was also observed for considerable flow speeds. However, in the case of the Harris sheet, a different scaling was found for a case with even larger flow speed. Of course, the growth rate scaling with resistivity deviated from the analytic scaling for stronger resistivities, as is to be expected, especially in the presence of flow. In both cases, flow was observed to introduce a cutoff resistivity above which the tearing mode is damped. In the case of the Harris sheet, this cutoff even lay inside the regime where the analytic scaling law is still expected to hold, though this was for the case that did not follow the analytic scaling law anyway.

Afterward, we showed that the growth rate scales with the Alfvén speed, though the scaling was found to be closer to c A 1 / 2 than to the analytic c A 2 / 5. Additionally, a significant deviation from this scaling was observed for larger densities, i.e., smaller Alfvén speeds, in the force-free field case with small plasma-β. The importance of the plasma-β was then further highlighted in this configuration, where an intricate interplay between the plasma-β and flow speed was observed. For large plasma-β, the flow was observed to destabilize the plasma further, with the strongest destabilization occurring at intermediate flow speeds. For small plasma-β though, both stabilizing and destabilizing influences were observed, with the strongest stabilization occurring at the Alfvén speed. Additionally, more than one stabilizing–destabilizing transition was observed along the flow speed axis for small plasma-β. Presumably, the importance of the plasma-β in this study is due to the inclusion of Joule heating, which is absent in preceding work.

Finally, it was shown that the addition of flow to an equilibrium with multiple magnetic nullplanes and tearing modes can modify the tearing behavior significantly. The tearing behavior of all modes was altered, with all modes tearing all nullplanes, albeit to different degrees. This is in stark contrast to the static case, where only the dominant mode tore all nullplanes, and the secondary and tertiary modes only tore the central and non-central nullplanes, respectively.

Looking ahead, future work could focus more on the role of the plasma-β, especially at small values. Additionally, Legolas could be used to incorporate viscosity,15 or to investigate transitions in instability dominance, from resistive tearing to the Kelvin–Helmholtz instability, at near-Alfvénic speeds. Also, the Hall field14 and electron inertia effects15 are prime candidates for further investigation. Of course, a similar study to this one can be performed for the cylindrical tearing mode,63 to compare the effects of axial and azimuthal flow. Finally, due to the growth rate depending on the specific flow profile, Legolas could also be employed as a computationally inexpensive, diagnostic tool for concrete configurations, particularly for experiments and for comparison of linear theory to non-linear simulations.

This work was supported by funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program, Grant Agreement No. 833251 PROMINENT ERC-ADG 2018. J.D.J. acknowledges further funding by the UK's Science and Technology Facilities Council (STFC) Consolidated Grant No. ST/W001195/1. R.K. is further supported by Internal Funds KU Leuven through the Project No. C14/19/089 TRACESpace and an FWO Project No. G0B4521N.

We thank the referees for their constructive and thought-provoking comments, which helped to improve the manuscript significantly.

The authors have no conflicts to disclose.

Jordi De Jonghe: Formal analysis (lead); Investigation (lead); Methodology (equal); Visualization (lead); Writing – original draft (lead); Writing – review & editing (lead). Rony Keppens: Funding acquisition (lead); Methodology (equal); Writing – original draft (supporting); Writing – review & editing (supporting).

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

All functionalities to reproduce the results presented here are available in Legolas v2.1.1. The Legolas code is freely available under the GNU General Public License. For more information, visit https://legolas.science/.

Due to the heavily localized transitions in the equilibrium profiles of the Harris current sheet, an equidistant grid would require many more grid points than a centrally accumulated grid to properly resolve the region of steepest change in the middle. Therefore, this study opted for a grid constructed using the algorithm below64 for an interval x [ a , b ] and function f.

Declare array auxGrid 
Declare array finalGrid 
auxGrid ( 1 ) a 
i 1 
while auxGrid ( i ) < ( a + b ) / 2 do 
   auxGrid ( i + 1 ) auxGrid ( i ) + f ( auxGrid ( i ) ) 
   i i + 1 
end while 
finalGrid ( 1 ) a 
finalGrid ( 2 i 1 ) b 
κ ( ( a + b ) / 2 auxGrid ( i 1 ) ) / ( auxGrid ( i ) auxGrid ( i 1 ) ) 
for j from 1 to i − 1 do 
   finalGrid ( j + 1 ) auxGrid ( j ) + κ f ( auxGrid ( j ) ) 
   finalGrid ( 2 i j 1 ) a + b finalGrid ( j + 1 ) 
end for 
finalGrid ( i + 1 ) ( a + b ) / 2 
Declare array auxGrid 
Declare array finalGrid 
auxGrid ( 1 ) a 
i 1 
while auxGrid ( i ) < ( a + b ) / 2 do 
   auxGrid ( i + 1 ) auxGrid ( i ) + f ( auxGrid ( i ) ) 
   i i + 1 
end while 
finalGrid ( 1 ) a 
finalGrid ( 2 i 1 ) b 
κ ( ( a + b ) / 2 auxGrid ( i 1 ) ) / ( auxGrid ( i ) auxGrid ( i 1 ) ) 
for j from 1 to i − 1 do 
   finalGrid ( j + 1 ) auxGrid ( j ) + κ f ( auxGrid ( j ) ) 
   finalGrid ( 2 i j 1 ) a + b finalGrid ( j + 1 ) 
end for 
finalGrid ( i + 1 ) ( a + b ) / 2 
This results in a symmetric grid of 2 i 1 grid points around the center of the interval, ( a + b ) / 2. A symmetric grid is desired because it reduces the likelihood that the spectrum's symmetry, in the case of an odd flow profile, is broken by numerical errors. In our specific case, we used the Gaussian function
f ( x ) = p 1 ( p 1 p 3 ) exp ( ( x p 2 ) 2 2 p 4 )
(A1)
with p 2 = 0 to obtain a grid strongly accumulated around the interval's center 0.
1.
J.
Lörinčík
,
J.
Dudík
,
G.
Aulanier
,
B.
Schmieder
, and
L.
Golub
, “
Imaging evidence for solar wind outflows originating from a coronal mass ejection footpoint
,”
Astrophys. J.
906
,
62
(
2021
).
2.
T. D.
Phan
,
J. L.
Verniero
,
D.
Larson
,
B.
Lavraud
,
J. F.
Drake
,
M.
Øieroset
,
J. P.
Eastwood
,
S. D.
Bale
,
R.
Livi
,
J. S.
Halekas
,
P. L.
Whittlesey
,
A.
Rahmati
,
D.
Stansby
,
M.
Pulupa
,
R. J.
MacDowall
,
P. A.
Szabo
,
A.
Koval
,
M.
Desai
,
S. A.
Fuselier
,
M.
Velli
,
M.
Hesse
,
P. S.
Pyakurel
,
K.
Maheshwari
,
J. C.
Kasper
,
J. M.
Stevens
,
A. W.
Case
, and
N. E.
Raouafi
, “
Parker solar probe observations of solar wind energetic proton beams produced by magnetic reconnection in the near-sun heliospheric current sheet
,”
Geophys. Res. Lett.
49
,
e2021
GL
096986
, https://doi.org/10.1029/2021GL096986 (
2022
).
3.
Y.
Qi
,
T. C.
Li
,
C. T.
Russell
,
R. E.
Ergun
,
Y.-D.
Jia
, and
M.
Hubbert
, “
Magnetic flux transport identification of active reconnection: MMS observations in Earth's magnetosphere
,”
Astrophys. J. Lett.
926
,
L34
(
2022
).
4.
P. A.
Sweet
, “
The neutral point theory of solar flares
,” in
Electromagnetic Phenomena in Cosmical Physics
, edited by
B.
Lehnert
(
Cambridge University Press
,
1958
), Vol.
6
, p.
123
.
5.
E. N.
Parker
, “
Sweet's mechanism for merging magnetic fields in conducting fluids
,”
J. Geophys. Res.
62
,
509
520
, https://doi.org/10.1029/JZ062i004p00509 (
1957
).
6.
H. E.
Petschek
, “
Magnetic field annihilation
,” in AAS-NASA Symposium on the Physics of Solar Flares (NASA Special Publication,
1964
), Vol.
50
, p.
425
.
7.
H. P.
Furth
,
J.
Killeen
, and
M. N.
Rosenbluth
, “
Finite-resistivity instabilities of a sheet pinch
,”
Phys. Fluids
6
,
459
484
(
1963
).
8.
J.
Birn
,
J. F.
Drake
,
M. A.
Shay
,
B. N.
Rogers
,
R. E.
Denton
,
M.
Hesse
,
M.
Kuznetsova
,
Z. W.
Ma
,
A.
Bhattacharjee
,
A.
Otto
, and
P. L.
Pritchett
, “
Geospace environmental modeling (GEM) magnetic reconnection challenge
,”
J. Geophys. Res.: Space Phys.
106
,
3715
3719
, https://doi.org/10.1029/1999JA900449 (
2001
).
9.
T.
Terasawa
, “
Hall current effect on tearing mode instability
,”
Geophys. Res. Lett.
10
,
475
478
, https://doi.org/10.1029/GL010i006p00475 (
1983
).
10.
A.
Fruchtman
and
H. R.
Strauss
, “
Modification of short scale-length tearing modes by the Hall field
,”
Phys. Fluids B
5
,
1408
1412
(
1993
).
11.
J. D.
Huba
and
L. I.
Rudakov
, “
Hall magnetic reconnection rate
,”
Phys. Rev. Lett.
93
,
175003
(
2004
).
12.
F.
Pucci
,
M.
Velli
, and
A.
Tenerani
, “
Fast magnetic reconnection: “ideal” tearing and the Hall effect
,”
Astrophys. J.
845
,
25
(
2017
).
13.
E.
Papini
,
S.
Landi
, and
L.
Del Zanna
, “
Fast magnetic reconnection: Secondary tearing instability and role of the Hall term
,”
Astrophys. J.
885
,
56
(
2019
).
14.
C.
Shi
,
M.
Velli
,
F.
Pucci
,
A.
Tenerani
, and
M. E.
Innocenti
, “
Oblique tearing mode instability: Guide field and Hall effect
,”
Astrophys. J.
902
,
142
(
2020
).
15.
J.
De Jonghe
,
N.
Claes
, and
R.
Keppens
, “
Legolas: Magnetohydrodynamic spectroscopy with viscosity and Hall current
,”
J. Plasma Phys.
88
,
905880321
(
2022
).
16.
X.
Wang
and
A.
Bhattacharjee
, “
Nonlinear dynamics of the m = 1 instability and fast sawtooth collapse in high-temperature plasmas
,”
Phys. Rev. Lett.
70
,
1627
1630
(
1993
).
17.
R. G.
Kleva
,
J. F.
Drake
, and
F. L.
Waelbroeck
, “
Fast reconnection in high temperature plasmas
,”
Phys. Plasmas
2
,
23
34
(
1995
).
18.
X.
Wang
,
A.
Bhattacharjee
, and
Z. W.
Ma
, “
Collisionless reconnection: Effects of Hall current and electron pressure gradient
,”
J. Geophys. Res.
105
,
27633
27648
, https://doi.org/10.1029/1999JA000357 (
2000
).
19.
H. J.
Cai
and
L. C.
Lee
, “
The generalized Ohm's law in collisionless magnetic reconnection
,”
Phys. Plasmas
4
,
509
520
(
1997
).
20.
L.
Yin
,
D.
Winske
,
S. P.
Gary
, and
J.
Birn
, “
Hybrid and Hall-MHD simulations of collisionless reconnection: Dynamics of the electron pressure tensor
,”
J. Geophys. Res.
106
,
10761
10776
, https://doi.org/10.1029/2000JA000398 (
2001
).
21.
M. E.
Mandt
,
R. E.
Denton
, and
J. F.
Drake
, “
Transition to whistler mediated magnetic reconnection
,”
Geophys. Res. Lett.
21
,
73
76
, https://doi.org/10.1029/93GL03382 (
1994
).
22.
M. A.
Shay
,
J. F.
Drake
,
B. N.
Rogers
, and
R. E.
Denton
, “
Alfvénic collisionless magnetic reconnection and the Hall term
,”
J. Geophys. Res.
106
,
3759
3772
, https://doi.org/10.1029/1999JA001007 (
2001
).
23.
Y.-H.
Liu
,
P.
Cassak
,
X.
Li
,
M.
Hesse
,
S.-C.
Lin
, and
K.
Genestreti
, “
First-principles theory of the rate of magnetic reconnection in magnetospheric and solar plasmas
,”
Commun. Phys.
5
,
97
(
2022
).
24.
A.
Bhattacharjee
, “
Impulsive magnetic reconnection in the Earth's magnetotail and the solar corona
,”
Annu. Rev. Astron. Astrophys.
42
,
365
384
(
2004
).
25.
W.
Daughton
,
V.
Roytershteyn
,
B. J.
Albright
,
H.
Karimabadi
,
L.
Yin
, and
K. J.
Bowers
, “
Transition from collisional to kinetic regimes in large-scale reconnection layers
,”
Phys. Rev. Lett.
103
,
065004
(
2009
).
26.
B.
Coppi
, ““
Inertial” instabilities in plasmas
,”
Phys. Lett.
11
,
226
228
(
1964
).
27.
K.
Shibata
and
S.
Tanuma
, “
Plasmoid-induced-reconnection and fractal reconnection
,”
Earth Planets Space
53
,
473
482
(
2001
).
28.
J. H.
Li
and
Z. W.
Ma
, “
Nonlinear evolution of resistive tearing mode with sub-Alfvenic shear flow
,”
J. Geophys. Res.: Space Phys.
115
,
A09216
(
2010
).
29.
J. H.
Li
and
Z. W.
Ma
, “
Roles of super-Alfvénic shear flows on Kelvin-Helmholtz and tearing instability in compressible plasma
,”
Phys. Scr.
86
,
045503
(
2012
).
30.
I.
Hofmann
, “
Resistive tearing modes in a sheet pinch with shear flow
,”
Plasma Phys.
17
,
143
157
(
1975
).
31.
R. K.
Pollard
and
J. B.
Taylor
, “
Influence of equilibrium flows on tearing modes
,”
Phys. Fluids
22
,
126
131
(
1979
).
32.
R. B.
Paris
and
W. N.
Sy
, “
Influence of equilibrium shear flow along the magnetic field on the resistive tearing instability
,”
Phys. Fluids
26
,
2966
2975
(
1983
).
33.
G.
Einaudi
and
F.
Rubini
, “
Resistive instabilities in a flowing plasma. I. Inviscid case
,”
Phys. Fluids
29
,
2563
(
1986
).
34.
X. L.
Chen
and
P. J.
Morrison
, “
Resistive tearing instability with equilibrium shear flow
,”
Phys. Fluids B
2
,
495
507
(
1990
).
35.
X.
Zhang
,
L. J.
Li
,
L. C.
Wang
,
J. H.
Li
, and
Z. W.
Ma
, “
Influences of sub-Alfvénic shear flows on nonlinear evolution of magnetic reconnection in compressible plasmas
,”
Phys. Plasmas
18
,
092112
(
2011
).
36.
L. N.
Wu
and
Z. W.
Ma
, “
Linear growth rates of resistive tearing modes with sub-Alfv'enic streaming flow
,”
Phys. Plasmas
21
,
072105
(
2014
).
37.
C.
Shi
, “
Instabilities in a current sheet with plasma jet
,”
J. Plasma Phys.
88
,
555880401
(
2022
).
38.
M.
Faganello
,
F.
Pegoraro
,
F.
Califano
, and
L.
Marradi
, “
Collisionless magnetic reconnection in the presence of a sheared velocity field
,”
Phys. Plasmas
17
,
062102
(
2010
).
39.
E.
Tassi
,
D.
Grasso
, and
L.
Comisso
, “
Linear stability analysis of collisionless reconnection in the presence of an equilibrium flow aligned with the guide field
,”
Eur. Phys. J. D
68
,
88
(
2014
).
40.
D.
Biskamp
,
Magnetic Reconnection in Plasmas
(
Cambridge University Press
,
2000
), Vol.
3
.
41.
N. F.
Loureiro
,
A. A.
Schekochihin
, and
D. A.
Uzdensky
, “
Plasmoid and Kelvin-Helmholtz instabilities in Sweet-Parker current sheets
,”
Phys. Rev. E
87
,
013102
(
2013
).
42.
R.
Keppens
,
G.
Tóth
,
R. H. J.
Westermann
, and
J. P.
Goedbloed
, “
Growth and saturation of the Kelvin-Helmholtz instability with parallel and antiparallel magnetic fields
,”
J. Plasma Phys.
61
,
1
19
(
1999
).
43.
D.
Borgogno
,
D.
Grasso
,
B.
Achilli
,
M.
Romé
, and
L.
Comisso
, “
Coexistence of plasmoid and Kelvin-Helmholtz instabilities in collisionless plasma turbulence
,”
Astrophys. J.
929
,
62
(
2022
).
44.
Y.
Park
,
S.
Sabbagh
,
J.
Bialek
,
J.
Berkery
,
S.
Lee
,
W.
Ko
,
J.
Bak
,
Y.
Jeon
,
J.
Park
,
J.
Kim
,
S.
Hahn
,
J.-W.
Ahn
,
S.
Yoon
,
K.
Lee
,
M.
Choi
,
G.
Yun
,
H.
Park
,
K.-I.
You
,
Y.
Bae
,
Y.
Oh
,
W.-C.
Kim
, and
J.
Kwak
, “
Investigation of MHD instabilities and control in KSTAR preparing for high beta operation
,”
Nucl. Fusion
53
,
083029
(
2013
).
45.
J.
Shao
,
H.
Liu
,
Y.
Xu
,
Z.
Chen
,
T.
Wang
,
J.
Cheng
,
X.
Wang
,
J.
Huang
,
H.
Liu
,
X.
Zhang
,
K.
Xu
,
C.
Tang
, and
T. J.-T.
Team
, “
Effect of the toroidal flow and flow shear on the m/n = 2/1 tearing mode in J-TEXT tokamak
,”
Plasma Phys. Controlled Fusion
63
,
065017
(
2021
).
46.
M. S.
Chu
,
J. M.
Greene
,
T. H.
Jensen
,
R. L.
Miller
,
A.
Bondeson
,
R. W.
Johnson
, and
M. E.
Mauel
, “
Effect of toroidal plasma flow and flow shear on global magnetohydrodynamic MHD modes
,”
Phys. Plasmas
2
,
2236
2241
(
1995
).
47.
R. L.
White
and
R.
Fitzpatrick
, “
Effect of rotation and velocity shear on tearing layer stability in tokamak plasmas
,”
Phys. Plasmas
22
,
102507
(
2015
).
48.
H.
Cai
and
J.
Cao
, “
Influence of toroidal rotation on magnetic islands in tokamaks
,”
Nucl. Fusion
58
,
036008
(
2018
).
49.
X. L.
Chen
and
P. J.
Morrison
, “
Nonlinear interactions of tearing modes in the presence of shear flow
,”
Phys. Fluids B
4
,
845
854
(
1992
).
50.
A. I.
Smolyakov
,
E.
Lazzaro
,
M.
Azumi
, and
Y.
Kishimoto
, “
Stabilization of magnetic islands due to the sheared plasma flow and viscosity
,”
Plasma Phys. Controlled Fusion
43
,
1661
(
2001
).
51.
Z.
Ren
,
F.
Wang
,
H.
Cai
, and
J.
Liu
, “
Influence of toroidal rotation on nonlinear evolution of tearing mode in tokamak plasmas
,”
Plasma Phys. Controlled Fusion
65
,
015007
(
2023
).
52.
N.
Claes
,
J.
De Jonghe
, and
R.
Keppens
, “
Legolas: A modern tool for magnetohydrodynamic spectroscopy
,”
Astrophys. J. Suppl. Ser.
251
,
25
(
2020
).
53.
N.
Claes
and
R.
Keppens
, “
Legolas 2.0: Improvements and extensions to an MHD spectroscopic framework
,”
Comput. Phys. Commun.
291
,
108856
(
2023
).
54.
Sometimes, the temperature is chosen as the constant quantity and subsequently, the density profile is determined by Eq. (8) such as in, e.g., Ref. 58, Sec. 14.4.1.
55.
L.
Ofman
,
P. J.
Morrison
, and
R. S.
Steinolfson
, “
Nonlinear evolution of resistive tearing mode instability with shear flow and viscosity
,”
Phys. Fluids B
5
,
376
387
(
1993
).
56.
M. A.
Cross
and
G.
Van Hoven
, “
Magnetic and gravitational energy release by resistive instabilities
,”
Phys. Rev. A
4
,
2347
2353
(
1971
).
57.
G.
Van Hoven
and
M. A.
Cross
, “
Energy release by magnetic tearing: The nonlinear limit
,”
Phys. Rev. A
7
,
1347
1352
(
1973
).
58.
H.
Goedbloed
,
R.
Keppens
, and
S.
Poedts
,
Magnetohydrodynamics of Laboratory and Astrophysical Plasmas
(
Cambridge University Press
,
2019
).
59.
H.
Betar
,
D.
Del Sarto
,
M.
Ottaviani
, and
A.
Ghizzo
, “
Microscopic scales of linear tearing modes: A tutorial on boundary layer theory for magnetic reconnection
,”
J. Plasma Phys.
88
,
925880601
(
2022
).
60.
H. P.
Furth
,
P. H.
Rutherford
, and
H.
Selberg
, “
Tearing mode in the cylindrical tokamak
,”
Phys. Fluids
16
,
1054
1063
(
1973
).
61.
P. L.
Pritchett
,
Y. C.
Lee
, and
J. F.
Drake
, “
Linear analysis of the double-tearing mode
,”
Phys. Fluids
23
,
1368
1374
(
1980
).
62.
H.
Betar
,
D.
Del Sarto
,
M.
Ottaviani
, and
A.
Ghizzo
, “
Multiparametric study of tearing modes in thin current sheets
,”
Phys. Plasmas
27
,
102106
(
2020
).
63.
B.
Coppi
,
J. M.
Greene
, and
J. L.
Johnson
, “
Resistive instabilities in a diffuse linear pinch
,”
Nucl. Fusion
6
,
101
117
(
1966
).
64.
A.
Quarteroni
,
Numerical Models for Differential Problems
(
Springer-Verlag
,
2009
), Vol.
2
.