This paper investigates the linear stability of a force-free equilibrium in a plasma cylinder of finite aspect ratio, against the two-fluid resistive tearing mode. An analytic dispersion relation is derived by extending to cylindrical geometry the slab geometry boundary layer theory applicable to two-fluid tearing modes for high beta and general ion skin depths [E. Ahedo and J. J. Ramos, Plasma Phys. Controlled Fusion 51, 055018 (2009)]. The cylindrical dispersion relation shows the dependence of the mode growth rate and real frequency on the ion skin depth, through different regimes that range from the single-fluid MHD limit to the electron MHD limit. It also shows that the non-zero real frequency of the mode arises due to the combination of two-fluid and cylindrical effects. A numerical solution of the complete set of normal-mode equations that resolves the fine-scale singular layer is carried out, for a wide range of resistivity and ion skin depth values. The numerically obtained eigenvalues agree very well with the analytic dispersion relation and the agreement improves the smaller the resistivity and the larger the ion skin depth are. Comparison between the numerical eigenfunctions and the inner solutions of the boundary layer theory shows that the eigenfunctions develop imaginary parts within the resonant layer, also due to the combination of two-fluid and cylindrical effects.
I. INTRODUCTION
The resistive tearing instability has been studied extensively as one of the physical mechanisms that drive magnetic reconnection, an important effect that is crucial to many space and laboratory plasma phenomena. The classical tearing mode is described in the framework of single-fluid resistive MHD.1 However, when the length scales of interest become comparable to the ion skin depth or to the ion Larmor radius, the small scale physics that includes the two-fluid and finite-Larmor-radius effects cannot be neglected. Two-fluid effects are considered as a mechanism for enhancement of the tearing mode instability in force-free equilibria.2–7 This can be studied in the framework of two-fluid (or Hall) MHD, which is one of the extended MHD models and includes the scale length of the ion skin depth. Extended MHD models are used to investigate small scale effects on macroscopic phenomena,8 and nonlinear numerical simulations with extended MHD models are presently a very active field of research.9–12 As the numerical models become more sophisticated, corresponding theories are also desirable for verification.7,13,14 The work of Ref. 7 derived an analytic dispersion relation for two-fluid tearing modes in the slab geometry, for a force-free equilibrium with constant density and temperature, allowing for general values of the ion skin depths and the plasma beta. This dispersion relation was benchmarked successfully against the extended-MHD simulation code NIMROD.10
In order to provide a more advanced theoretical model, this paper investigates the two-fluid tearing mode instability in a plasma cylinder of finite aspect ratio. Tearing modes in cylindrical geometry1,15–18 add several complications due to the curvature, the radial dependence, and the periodicity in the poloidal direction. Cylinder-like structures are observed in space plasmas such as flux ropes in the heliosphere as well as in laboratory plasmas. We derive an analytic dispersion based on the boundary layer theory with asymptotic matching, assuming a singular perturbation near a resonant surface as in previous studies.1,6,7,15,19 This new dispersion relation adds the cylindrical effects to the dispersion relation of Ref. 7 in the high beta regime. A fully numerical solution of the complete set of normal-mode equations that resolves the fine-scale singular layer is then carried out. These numerical results are compared with the predictions of the analytic boundary layer theory for a wide range of resistivity and ion skin depth values.
The paper is organized as follows. In Sec. II, we introduce the model equations for the two-fluid tearing instability in cylindrical geometry. In Sec. III, we derive the analytic dispersion relation using the boundary layer theory. In Sec. IV, the equilibrium profiles adopted for our analysis are introduced. In Sec. V, the solution of the analytic dispersion relation is obtained and compared with approximations. In Sec. VI, the fully numerical analysis of the complete system of linear eigenmode equations is carried out. A summary and discussion are given in Sec. VII.
II. THE MODEL
The goal of this work is to investigate the effects of two-fluid physics and cylindrical geometry on the linear resistive tearing mode, neglecting the finite-Larmor-radius, electron inertia, and equilibrium pressure gradient effects. To this end, we consider the following system of resistive, cold-ion, Hall-MHD equations:
The resistivity η, which will be assumed to be small and will matter only in the close vicinity of the mode resonant surface, will be treated as a constant evaluated at such a resonant surface. In order to close this system, a constant electron temperature, pe/n = Te0 = constant, will also be assumed. This system will be linearized about a force-free, cylindrical equilibrium with vanishing density gradient and flow velocity: n0 = constant, u0 = 0. Accordingly, the equilibrium magnetic field is
and the equilibrium current density is j0 = λ0B0, where
with the prime denoting differentiation with respect to r. The squared Alfvén and sound velocities are defined as and , respectively.
The linear stability analysis will consider normal-mode perturbations of the form
with wavevector . It is convenient to introduce the parallel wavevector function
such that F(rs) = 0 at the mode resonant surface rs, and the perpendicular wavevector function
The main variables of interest are the radial component of the perturbed magnetic field B1r, the radial component of the fluid displacement ξ ≡ u1r/γ, and the non-ideal part of the z-component of the perturbed magnetic field
All the other components of the perturbation can be eliminated algebraically in favor of these three and, without further approximations, the resulting normal-mode system for (B1r, ξ, Q) is
III. ANALYTIC DISPERSION RELATION
An analytic, perturbative solution of the normal-mode system (12)–(14) can be obtained in the limit of small resistivity, characterized by large values of the dimensionless Lundquist number, S ≡ cA/(ηk) ≫ 1. The normalized growth rate will be proportional to a negative fractional power of S and much less than unity. Therefore, resistive and inertial effects are only important within a narrow layer near the resonant surface rs. Outside that layer, the solution can be approximated by its marginally stable, ideal-MHD limit defined by the conditions γ = η = 0, which yield Q = 0, B1r = iFξ and
This radial equation for B1r is to be solved subject to a regularity condition at r = 0 and an appropriate boundary condition at the plasma edge r = rW. In general, one such solution can be found that is continuous but has a discontinuous derivative at the mode resonant surface. This discontinuity defines the tearing mode stability index
The inner solution in a narrow layer about r = rs takes into account the resistive and inertial effects and will match asymptotically the ideal-MHD outer solution, regularizing its singularity and allowing for a smooth global representation of the mode. Within this layer, the radial variable x ≡ (r – rs) ≪ rs is to be used as the argument of the perturbed quantities, while the equilibrium quantities can be approximated by their constant values at r = rs. The exception is the parallel wavenumber function F (9) that vanishes at r = rs and must be represented by the leading term of its Taylor expansion
where LB is the magnetic shear length
which, for a periodic cylinder of length 2πR and equilibrium inverse rotational transform , is
The macroscopic lengths, r, LB, k−1, and are formally ordered as comparable, whereas the microscopic length x scales as with . Then, the radial derivative of a perturbed quantity within the layer will be taken as , and only the leading terms as δ → 0 will be retained. Three additional approximations will be made. The first one is the “constant-rB1r” (or “constant-ψ”) approximation that neglects but retains in the inner layer equations and is valid provided . The second approximation is the neglect of ηk2 when compared to γ, which is justified by the anticipated result that the growth rate will be proportional to a fractional (less than one) power of the resistivity. The third approximation is the restriction to a deeply subsonic regime, where γ ≪ kcS and the fluid perturbation becomes incompressible with a growth rate independent of the ratio . This is guaranteed for an “effectively high-β” characterized by , which is a well satisfied condition in most situations of interest. With these assumptions, the general normal-mode system (12)–(14) reduces to the singular layer system
Here, B1r is to be considered as a normalizing constant and the integral across the singular layer is to be equated to the macroscopic jump . We have introduced and and will use these quantities for normalization such as and in the rest of this section. Then, (16) and (20) yield the matching condition
Only the parts of ξ and Q that are odd with respect to x contribute to this matching condition that will determine the mode dispersion relation, so the terms in (21) and (22) that drive the even parts of ξ and Q can be ignored here. Thus, introducing the ion skin depth and the scaled variables
where
and ρ and σ are two dimensionless parameters that take, respectively, into account the effects of the cylindrical geometry and the two-fluid physics
The coupled system (27) and (28) for and can be diagonalized and its solution can be expressed as the linear combination of two uncoupled eigenfunctions W±, with appropriately scaled arguments. This solution is
where
and the diagonalizing eigenfunctions W± satisfy the uncoupled equations
where
It is noted that, even if y is complex due to a complex γ, as will be seen in Secs. V and VI, the integration with respect to the imaginary part of y does not contribute to Eq. (38) as long as the integrand is a regular function. Equation (36) has the Hermite series solution
where H2n+1 are the odd Hermite polynomials. After integrating with respect to y, the ensuing series in (38) can be summed to arrive at the compact result
In summary, we have obtained the following analytic dispersion relation for two-fluid, Hall-resistive tearing modes in cylindrical geometry, in the “effectively high-β” incompressible regime:
where
Several known results can be recovered as special limits of this dispersion relation. The single-fluid limit corresponds to σ → 0, hence
which is the dispersion relation for single-fluid resistive tearing modes in cylindrical geometry, derived in Ref. 15. For σ → ∞, we obtain
which is the generalization to cylindrical geometry (i.e., finite ρ) of the electron-MHD result derived for the slab geometry in Ref. 4 and coincides also with Eq. (A10) of Ref. 17 in the high beta limit without equilibrium pressure gradient. Finally, the slab geometry limit corresponds to ρ → 0, hence
which is the slab geometry dispersion relation for Hall-MHD tearing modes in their high beta regime derived in Ref. 7 and used in the NIMROD code benchmark of Ref. 10. Equations (43) and (45) include also the dispersion relation for single-fluid resistive MHD in the slab geometry limit1,20
IV. EQUILIBRIUM PROFILES
In order to specify a force-free cylindrical equilibrium upon which we shall base our numerical calculations, we consider the following choice of profiles:
For numerical purposes, we normalize lengths to a, inverse lengths to a−1, magnetic fields to Bc and time to a/cAc, where . The thus normalized quantities are denoted by overbars. We choose the parameters ε = 0.25, q0 = 0.25, m = 2, and that satisfy the assumptions of the analytic dispersion relation. Then, the other parameters result in , and . The profiles of and are shown in Fig. 1. It is noted that we have chosen the poloidal mode number m = 2 since, for m = 1, can be too large to satisfy the assumption for the analytic dispersion relation and, for large m, can become negative making the MHD tearing mode stable.16,20
V. SOLUTIONS OF THE ANALYTIC DISPERSION RELATION
where
The growth rate for single-fluid MHD in the slab geometry, , is obtained from D(0, 0) as
where
In the present equilibrium, is calculated from the numerical solution of the outer solution. The equation for in the outer regions is written in the normalized form as
We consider the region with a fixed boundary at rW = 7.5a. Equation (55) is solved numerically with a matrix solver, taking the sixth order finite differences with 100 000 intervals, subject to the boundary conditions
where . The numerical solution for is shown in Fig. 2 and is calculated as
After substituting the values of (58) and the other parameters, the dispersion relation (51) is solved numerically for by root finding methods. Although ρ = 0.09851 is given by our choice of cylindrical equilibrium, we also consider the case of ρ = 0 (the slab limit) with the other parameters fixed in order to examine the cylindrical effect. In the case of single-fluid MHD (σ = 0), the tearing mode is purely growing and, taking , the growth rate is for ρ = 0.09851, which shows a slight stabilizing effect of the cylindrical geometry compared to for the slab geometry. Figure 3 shows the dependence of the real and imaginary parts of on the ion skin depth for ρ = 0 and ρ = 0.09851. The real part of (i.e., the growth rate) is enhanced by the two-fluid effect. The effect of cylindrical geometry is slightly stabilizing for each value of . The imaginary part of (i.e., the real frequency of the mode) arises when both of the two-fluid and cylindrical effects exist and increase with . In the slab geometry, a real frequency associated with the diamagnetic drifts would appear if the equilibrium had a pressure gradient.21,22
Since in the present case and ρ ≪ 1, we derive approximations of the growth rate and the real frequency by using the expansions of D(ρ, σ) shown in Appendix A. The dispersion relation becomes algebraic in two limits of σ, σ ≪ 1, and σ ≫ 1. In the case of σ ≪ 1, the dispersion relation (51) is expanded with the ordering σ ∼ ϵ, ρ ∼ ϵ, and , where and . The real and imaginary parts of the dispersion relation are, respectively, given by
where
Equation (59) shows small corrections to due to σ and ρ. The second term in the parentheses of the RHS of (59) shows the enhancement of the growth due to the two-fluid effect. The third term, on the other hand, shows stabilization due to the cylindrical effect since C2 > 0. Equation (60) shows that the real frequency appears due to the combination of the two-fluid and cylindrical effects. In the limit of σ ≫ 1, the dispersion relation (51) (with ρ > 0 and σ ≫ 1) is expanded with the ordering ρ ∼ ϵ and . The growth rate is then given by
This shows the small correction due to the cylindrical effect, which is also stabilizing compared to the dispersion relation for electron MHD in the slab geometry (ρ = 0 and σ ≫ 1)4
The real frequency is given by
Figure 4 shows the dependence on the ion skin depth of the real and imaginary parts of (now in logarithmic scales), compared with the approximations for σ ≪ 1 and σ ≫ 1. The dispersion relation (51) is well approximated by Eqs. (59) and (60) for and by (63) and (65) for . Figure 5 shows the dependence of the real and imaginary parts of on the resistivity , compared with the approximations for σ ≪ 1 and σ ≫ 1. For and , the growth rate is well approximated by (59) that yields the single-fluid MHD scaling , whereas for it is well approximated by (63) that yields the electron MHD scaling . For , the real frequency is zero, for it is well approximated by (59) that yields a scaling proportional to , and for it is well approximated by (63) that yields a scaling proportional to . Thus, we conclude that covers the range of ion skin depths of interest to show the two-fluid effects on the tearing instability in cylindrical geometry, for the present set of equilibrium parameters.
VI. NUMERICAL EIGENMODE ANALYSIS
In order to compare the results derived from the analytic dispersion relation with those obtained with a fully numerical approach, we solve numerically the complete set of normal-mode equations on a sufficiently fine mesh that resolves the fine-scale singular layer and yields the global eigenfunctions. For convenience in the numerical analysis, we solve the eigenmode equations for instead of Eqs. (12)–(14) and calculate and using the relations
The boundary conditions at and are
The condition (69) is the fixed boundary condition for radial perturbations. We use the condition (70) on only as a boundary condition at . The eigenmode equations for are numerically solved to obtain eigenvalues of and eigenfunctions. We discretize the variables by taking the sixth order finite differences with 100 000 intervals in the region with the boundary conditions (68)–(70) to obtain a matrix equation. We find the eigenvalue of that satisfies the equation that the determinant of the matrix of the coefficients is zero. Substituting this solution of , the matrix equation is solved to obtain the eigenfunctions. We choose β = 0.1 as a high beta value.
Figures 6–9 show the eigenfunctions normalized to for different values of and . In (a) of Figs. 6–9, the eigenfunctions agree well with the outer solution of Fig. 2 in the outer region. For the single-fluid MHD cases (, Figs. 6 and 8), and are real, and and are imaginary. They all become complex due to the two-fluid effect (see Figs. 7 and 9). Figure 10 shows the variation of the real and imaginary parts of with the ion skin depth for , as obtained from the fully numerical solution of the eigenmode equations and compared with the analytic dispersion relation. There is excellent agreement in the range . Figure 11 shows the variation of the real and imaginary parts of with the resistivity , compared with the analytic dispersion relation. As expected, the agreement deteriorates with increasing since the analytic result is based on the asymptotic matching method of a singular perturbation theory for and requires that the diffusive part of the perturbation be localized within an infinitesimally narrow layer around the resonant surface. This is illustrated in the single-fluid MHD case () by Figs. 6(c) and 8(c) which show how the eigenfunction for is much narrower than for .
In order to further clarify the structure of the two-fluid tearing eigenfunctions and their parametric dependences on the resistivity and the ion skin depth, as well as the parametric dependence of the growth rate eigenvalue, we calculate the eigenfunctions and that determine the analytic dispersion relation of Sec. III from the eigenfunctions and and the eigenvalue of of the numerical solution. The two representations are related by
where and are extracted using
In (a)–(d) of Figs. 12 and 13, is dominant when is small while, in (a)–(d) of Figs. 14 and 15, is dominant when is large. As increases, the width of decreases, while the width of increases. These behaviors are the same as in the slab geometry.7 As pointed out in Refs. 2 and 7, the inner solutions of two-fluid tearing mode vary on two length scales. These scale lengths are for σ ≫ 1 as in Ref. 7
Equations (74) and (75) give, respectively, the scale lengths associated with the magnetic diffusion and two-fluid effect, such that d1 < d0 < d2. The diffusion in the narrower scale width enhances the growth of the tearing mode.2 Since the diffusive length scale of for is narrower than that of for or , the diffusion is dominant at large . This explains the better agreement between the analytical and numerical growth rates at large and , compared to the same but smaller . Equation (74) also shows that cases with smaller and larger require higher spatial resolution of calculation. The numerical code used in this study had good resolution up to and .
Figures 12–15 also show that and have imaginary parts that do not exist in the slab geometry.7 Approximations of and for σ ∼ ρ ∼ ϵ ≪ 1 can be obtained by expanding Eqs. (33) and (34) with respect to σ and ρ as
Equations (76) and (77) show that the imaginary parts of both and appear due to the combination of two-fluid and cylindrical effects. Equation (77) also shows that the imaginary part of (i.e., ) appears due to the cylindrical effect even in the single-fluid MHD case. Equations (76) and (77) are compared with the numerical eigenfunctions in Figs. 12 and 16 for small values of and σ: σ ≃ 0, 0.0013 and 0.13, and , and 10−2. In (a) of Figs. 12 and 16, is dominant and does not much depend on , in agreement with Eq. (76) that tells that is of order unity and the two-fluid effect is second-order. In (b) of Figs. 12 and 16, is proportional to but much smaller than the first-order of , in agreement with Eq. (76) that tells that is proportional to ρσ. Figures 12(e) and 16(c) show that is proportional to and comparable to the first order of , in agreement with Eq. (77) that tells that is proportional to only σ. Figures 12(f) and 16(d) show that is finite even for and does not much depend on , in agreement with Eq. (77) that tells that is proportional to only ρ.
VII. SUMMARY AND DISCUSSION
Using singular perturbation boundary layer techniques, we have derived analytically a two-fluid tearing mode dispersion relation for force-free equilibria in a plasma cylinder of the finite aspect ratio. This dispersion relation shows the variation of the growth rate and the real frequency of the mode with the ion skin depth and includes previous theories as special limits. We have found that the growth rate is slightly reduced relative to the slab geometry result due to the cylindrical curvature effect and that the real frequency appears due to the combination of the two-fluid and cylindrical effects. Scaling laws for the real frequency and the growth rate as functions of the resistivity have been obtained by expansion of the dispersion relation in different asymptotic limits of the ion skin depth. We have also solved numerically the complete set of normal-mode equations and we have compared the results with those from the analytic dispersion relation for a wide range of resistivity and ion skin depth values. The agreement is uniformly good for the dependence on the ion skin depth and the agreement with regard to the dependence on the resistivity improves for smaller resistivity and larger ion skin depth. We have analyzed the two-fluid and cylindrical effects on the numerical eigenfunctions in comparison with the inner solutions of the boundary layer theory and shown that the imaginary parts of the solutions arise due to the combination of two-fluid and cylindrical effects.
In the present results, the real frequency is small compared with the growth rate. Equations (59), (60), (63), and (65) suggest that this is due to the small value of ρ in the considered equilibrium and that the real frequency could be comparable to the growth rate and large enough to be observed in other simulations12 when ρ ∼ 1. A parametric study of the dependence on ρ, consistent with a set of different equilibrium profiles, will be shown elsewhere.
ACKNOWLEDGMENTS
This work was partially supported by Japan/U.S. Cooperation in Fusion Research and Development, Grant-in-Aid for Scientific Research, KAKENHI (24740375), the NIFS Collaborative Research Program (NIFS15KNSS057 and NIFS15KNXN306) and by U.S. D.O.E Grant Nos. DEFG02-91-ER54109 and DEFC02-08ER54969.
APPENDIX A: EXPANSIONS OF D (ρ, σ)
For arbitrary ρ and σ ≪ 1, the expansion of D(ρ, σ), Eq. (42), with respect to σ up to the second order is
For ρ ≪ 1 and arbitrary σ, the expansion of D(ρ, σ), Eq. (42), with respect to ρ up to the second order is
For ρ ∼ σ ≪ 1, the expansion of D(ρ, σ) with respect to ρ and σ up to the second order is obtained from (A1) or (A2) as
For ρ ≪ 1 and the limit of σ ≫ 1, the expansion of D(ρ, σ) with respect to ρ is obtained from (44) or (A2) as