Turbulence driven magnetic reconnection causing long-wavelength magnetic islands

Magnetic reconnection caused by turbulence in a current sheet is studied by means of numerical simulations of fluid equations. It is found that turbulence produces long-wavelength magnetic islands even if the current sheet is so thick that spontaneous magnetic reconnection does not occur. Thus, turbulence modifies the threshold of magnetic island formation predicted by the conventional theory of spontaneous magnetic reconnection in a current sheet. In spite of the fact that the turbulence is driven by a short-wavelength instability due to a pressure gradient, the length of the magnetic island is the same order as the system size. The width of the island is several times the ion Larmor radius, and stronger turbulence causes wider magnetic islands. This suggests that the turbulence can trigger neoclassical tearing modes, which are the main nonlinear instability that limits the plasma pressure in magnetically confined plasmas. The long-wavelength magnetic island is formed by merging of small-scale magnetic islands. © 2010 American Institute of Physics. doi:10.1063/1.3463435


I. INTRODUCTION
2][3] These magnetic reconnection phenomena cause magnetic islands or plasmoids, and they can strongly affect the dynamics of a plasma. 1 Tearing instabilities are spontaneous magnetic reconnection that may occur in sheared magnetic configuration, such as a current sheet and a rational surface in a magnetic confinement. 4The conventional theory based on resistive magnetohydrodynamics ͑MHD͒ ͑Ref.4͒ shows that the instability is driven by the current density gradient of a current sheet and is caused by resistivity.When a current sheet is narrow, it is unstable for long-wavelength perturbation, and magnetic islands appear.The theory shows that the stability parameter of tearing instability ⌬Ј determines the instability threshold 4 and thus the formation of magnetic islands.
Typical time scale of tearing instability depends on the resistivity as 3/5 , 4 and its strong dependence on the resistivity is similar to the time scale of reconnection described by the Sweet-Parker model, 5,6 which predicts the scaling of 1/2 , where is the plasma resistivity.Fast magnetic reconnection, 1,3 which is independent of or weakly depends on the value of the resistivity, can be explained, for instance, by the appearance of Hall region around the reconnection point 1,3,[7][8][9] or by the appearance of turbulence.1][12][13][14][15][16][17][18][19][20] Theoretical studies show an anomalous magnetic flux diffusivity in Ohm's law 11 and the increase of reconnection rate in the presence of a stochastic component of the magnetic field. 12Numerical simulations of two-dimensional MHD suggest an enhancement of the magnetic reconnection by background turbulence 13 and a weaker dependence of reconnection on resistivity than the Sweet-Parker scaling. 20In magnetically confined plasmas, smallscale turbulence is driven by drift-wave instabilities, 21 and the turbulence can generate zonal flow and zonal magnetic field. 22The zonal magnetic field can affect the stability of tearing instabilities. 23Turbulence cannot only accelerate magnetic reconnection, but it also could change the threshold for magnetic island formation.In this paper, in order to understand the effect of turbulence on the threshold for magnetic island formation, we investigate turbulence in a current sheet, which is so thick that there is no spontaneous magnetic reconnection by numerical simulations based on an MHD model and on a twofluid model.We will demonstrate that long-wavelength magnetic islands are produced by turbulence in a current sheet even if the current sheet is so thick that the sheet is stable against tearing instabilities.The magnetic island is produced through energy transfer from small-scale turbulence.The threshold for magnetic island formation is usually given by the stability parameter of tearing instability; however, our results imply that the turbulence changes the threshold for magnetic island formation.
In our simulations turbulence is driven by shortwavelength instabilities.In order to include the shortwavelength instability we consider two models: one is a twodimensional MHD model including the interchange term; the other is a two-dimensional reduced two-fluid model.These two models describe typical short-wavelength instabilities.The former describes generalized Rayleigh-Taylor instabilities, and the latter describes drift-wave instabilities.Numerical simulations based on the two-fluid model show that turbulence driven by drift-wave instabilities causes longwavelength islands.The length of the magnetic island is the same order as the system size.The width of the island is several times the ion Larmor radius, and stronger turbulence causes wider magnetic islands.The long-wavelength magnetic island is formed by merging of small-scale magnetic islands.Numerical simulations based on the MHD model show that long-wavelength magnetic islands are produced through energy transfer from small-scale turbulence driven by short-wavelength instability due to a pressure gradient.
The paper is organized as follows.Two-fluid simulation results are shown in Sec.II.MHD simulation results are presented in Sec.III.Finally, summary is given in Sec.IV.

II. TWO-FLUID SIMULATION OF MAGNETIC ISLAND FORMATION
In this section, we examine turbulence driven by a typical drift-wave instability, which is the ion temperature gradient instability, and its effect on magnetic island formation by means of numerical simulations of a two-dimensional reduced set of two-fluid equations.

A. Simulation model
A reduced two-fluid plasma model is useful for the analysis of instabilities in strongly magnetized plasmas such as large aspect ratio magnetically confined plasmas with strong guide field. 24The well-known reduced two-fluid model 24 assumes that the inverse aspect ratio of a plasma column is small and that the guide field is much stronger than the magnetic field generated by a plasma current parallel to the guide field.Since the guide field is very strong, we can assume as follows: the plasma beta is small, plasma motion is slower than the fast magnetic wave, and variation along the guide field is smaller than variation normal to the field ͑see Ref. 24 or Sec.6.3 in Ref. 1͒.In addition, it is assumed that the spatial variations of density and temperature are small. 24In order to include ion temperature gradient instabilities, ion temperature fluctuation is included by using the Landau fluid closure, [25][26][27] while electron temperature fluctuation is omitted.A reduced two-fluid model in twodimensional slab geometry is obtained by assuming that the plasma is uniform along the strong guide field directed along the z-axis.The reduced model contains shear-Alfvén waves, sound waves, drift waves, and kinetic Alfvén waves.When the two-fluid parameter, ion skin depth d i , is zero, the electrostatic potential and the magnetic flux function obey the reduced MHD equations.The reduced two-fluid model is capable of describing nonlinear dynamics of tearing instability and ion temperature gradient instability.The model consists of the vorticity equation, the electron density equation, the parallel velocity equation, the generalized Ohm's law, and the ion temperature equation, where v In these equations, ⌽, , n, , p e = nT e , p = p i + p e are the electrostatic potential, the magnetic flux function, the electron density, the negative of current density, the parallel ion velocity, the parallel electron velocity, the ion temperature, the electron temperature, the ion pressure, the electron pressure, and the total pressure, respectively, where "parallel" means component along the guide magnetic field directed along the z-axis, and the electron temperature consists only of a given equilibrium part T e = T eq .The normalizations are and L are the ion skin depth, the plasma resistivity, the plasma beta, the ratio of specific heats, and the system size, respectively, where ⍀ i is the ion cyclotron frequency and i is the ion Larmor radius.The term including L = ͱ 8T eq / ͉B eq • ٌ͉ is the linear Landau damping term, which is obtained from the Hammett-Perkins closure. 28In order to clearly capture linear instabilities of the initial equilibrium, the field quantities are divided into the initial equilibrium part and the rest, as f͑x , y , t͒ = f eq ͑x͒ + f ˜͑x , y , t͒, and it is assumed that there are no initial equilibrium E ϫ B flow and no electron temperature fluctuation so that = eq + ˜, ⌽ = ⌽ ˜, n = n eq + n ˜, T i = T eq + T ˜i, T e = T eq , p i = n eq T eq + T eq n ˜+ n eq T ˜i, and p e = n eq T eq + T eq n ˜.Notice that f ˜can be the same order as f eq .The fluctuating component f ˜does not include k y = 0 component except ˜and ⌽ ˜so that the pressure gradient remains and zonal flow and zonal magnetic field are involved.In order to suppress high wavenumber perturbations the artificial dissipation term ٌ 2 f is added to each equation except in the Ohm's law, where the dissipation coefficient is =10 −6 .We numerically solve Eqs.͑1͒-͑5͒.Numerical settings are as follows.The simulation domain is rectangular in ͑x , y͒ plane and is periodic in the y direction.Time advancement is yielded by the fourth-order Runge-Kutta method.The Fourier decomposition is used in the y direction as f ˜= ͚ k y f ˆ͑x , k y , t͒exp͑ik y y͒, and the x derivative is approximated by a finite-difference method.When nonlinear terms are calculated, the pseudospectral method is used in the y direction, and a finite-difference method is used in the x direction.Notice that the nonlinear ion diamagnetic term in the vorticity equation is neglected in our simulations because of numerical difficulty caused by strong production Magnetic field configuration in ͑x , y͒ plane is the Harris current sheet, d eq ͑x͒ / dx = tanh͑x / L s ͒.The magnetic shear length L s of the Harris sheet is associated with the gradient of current density and controls the strength of spontaneous magnetic reconnection called tearing instability, and small L s , i.e., a thin current sheet, causes strong instability.The uniform current density model d eq ͑x͒ / dx = x / L s , which represents an infinitely thick current sheet, is also used.This profile has the neutral sheet at x = 0 and is stable against tearing instability.Following the MHD analysis, 4 the stability parameter of tearing instability ⌬Ј = ͑L / 2͓͒d ln ˜/ dx͔ 0− 0+ is obtained by solving the linearized form of B • ٌJ = B • ٌٌ 2 = 0.An equilibrium is stable ͑unstable͒ against tearing instability when ⌬Ј is negative ͑positive͒.In the following subsection, we consider the cases with ⌬Ј = −0.76 and Ϫ2.4 for the Harris sheet with L s = 0.3125 and 0.5, respectively, and the case with ⌬Ј = −3.4 for the uniform current density profile, where ⌬Ј is calculated for the longest wavelength in the domain.It is confirmed that these profiles are stable against tearing instabilities by numerically solving the linearized version of Eqs.͑1͒-͑5͒, in spite of the fact that two-fluid effects can affect tearing instability in the two-fluid model.

B. Magnetic islands caused by turbulence
We consider a plasma that is unstable against ion temperature gradient instability, which is a drift-wave instability and has short wavelengths characterized by the ion Larmor radius. 21The ion temperature gradient instabilities drive turbulence and cause anomalous heat transport in magnetically confined plasmas such as tokamaks. 21The instability is driven by an ion temperature gradient and is controlled by the parameter i = L n / L T , where L T =−͑d ln T eq / dx͒ −1 and L n =−͑d ln n eq / dx͒ −1 are temperature gradient and densitygradient lengths, respectively.Both density and temperature gradients are assumed to be uniform, and L n is the same as the system size, L n = L, in our simulations.The instability grows up exponentially at first; then the energy of the instability spreads in the Fourier space via nonlinear mode coupling.The instability also produces strong zonal flow, which is the k y = 0 mode, and the zonal flow then regulates the amplitude of the turbulence by means of its shear, 22 and the system reaches a quasi-steady state.
Figure 1 shows color map of electrostatic potential and equicontours of magnetic flux for the Harris sheet with L s / L = 0.3125, ⌬Ј = −0.76,and i = 3.5.The stability parameter ⌬Ј is calculated for the longest wavelength mode, k =1, where the wavenumber is given by k y = k / L. At t = 264, narrow magnetic islands appear at the neutral sheet, and electrostatic potential profile exhibits turbulent fluctuations.At t = 1188, magnetic islands become large, and zonal flow, which is uniform flow along the y-axis, appears.At the quasi-steady state, t = 1848, we clearly see magnetic islands appearing on the neutral sheet x = 0.The length of the mag-netic island is much longer than the ion Larmor radius and is the same order as the system size, and the width is several times the ion Larmor radius i / L =1/ 80.The shape of the magnetic island is deformed and is different from that of tearing instability because of turbulent fluctuations of the magnetic field.Note that the magnetic islands propagate toward the electron diamagnetic direction.Figure 2 shows energy spectrum of magnetic fluctuation.When the instability grows up exponentially at t = 66, the energy spectrum has a peak.After t Ϸ 130, the growth of instability is saturated and the energy of the short-wavelength instability spreads over the Fourier space via nonlinear mode coupling, and thus the energy spectrum becomes broad so that the long-wavelength modes are excited.Hence, the long-wavelength magnetic islands are formed.
Figure 3 shows the magnetic island width as a function of the stability parameter of tearing instability ⌬Ј for the case i = 3.5, where the width is time-averaged in the quasi-steady state.The stability parameter ⌬Ј is Ϫ2.4,Ϫ0.76, 0, 0.79, and 1.8 for the Harris sheet with L s = 0.5, 0.3125, 0.275, 0.25, and 0.225, respectively, and the parameter ⌬Ј is Ϫ3.4 for the uniform current density.The magnetic island width is several times the ion Larmor radius even if the current sheet is stable against tearing instability, ⌬Ј Ͻ 0. Hence, turbulence modifies the island formation threshold ⌬Ј = 0 given by the conventional theory. 4Figure 4 shows the dependence of the magnetic island width on the system size for the uniform current density profile that is stable against tearing instabilities, ⌬Ј = −3.4.The island width is about five times the ion Larmor radius and does not depend on the system size.
Figure 5 shows time trace of the reconnected flux that is the largest value of magnetic flux on the neutral sheet, x =0, for each temperature gradient parameter i that indicates the strength of the instability.When turbulence is strong, i = 4.5, the reconnected flux grows rapidly, and then it attains a large value and a quasi-steady state is formed.On the other hand, when the instability is weak, i = 2.5, the reconnected flux grows slowly, and it attains a small value after which it gets saturated.Thus, stronger turbulence causes faster reconnection and wider magnetic islands.
Figure 6 shows the dependence of the width on the turbulence energy E ‫ء‬ = ͚ k=1 N 1 2 ͉͐v Ek ͉ 2 dx in the quasi-steady state for the cases with ⌬Ј = −3.4 and Ϫ2.4,where the width is time-averaged in the quasi-steady state.Note that E ‫ء‬ does not include zonal flow component k = 0. Squares refer to the width for the uniform current density model with ⌬Ј = −3.4 and i = 2.5, 3.5, and 4.5.Circles refer to the width for the Harris current sheet with ⌬Ј = −2.4 and i = 2.5, 3.5, 4.5, and 5.5.The magnetic island width is large when turbulence is strong, W ϰ ͑E ‫ء‬ ͒ ␣ , where ␣ Ϸ 0.25.
We remark that the effect of zonal magnetic field, which changes the initial equilibrium magnetic field, on tearing instability is not responsible for the magnetic island formation in our simulations.This is because we have magnetic islands, which has the width W / i = 6.5 for the case ⌬Ј = −3.4 and i = 3.5 even if the zonal magnetic field is eliminated by setting ͗ ˜͘ = ͗J ˜͘ = 0, where ͗f͘ = ͐fdy / L is the average over y.

C. Formation process of magnetic islands
We examine the process of magnetic island formation.Figure 7 shows the time evolution of magnetic flux on the neutral sheet, x = 0, for the Harris sheet with ⌬Ј = −2.4,where the red region represents positive flux and blue represents negative flux, which represents the core parts of magnetic islands around O-points.Initially, the color patterns show  wave structure characterized by the wavelength of ion temperature gradient instability.After t Ϸ 200, large red and blue regions appear, and this means that long-wavelength magnetic islands are formed.Figure 8 shows the time evolution of the positions of X-points and O-points for the Harris sheet with ⌬Ј = −2.4.The red and black points represent X-points and O-points, respectively.Initially, there are many X-and O-points, and then some merge and disappear.In this process, small-scale magnetic reconnection may play a role, as described in Ref. 31.After the large red and blue regions appear in Fig. 7, sometimes X-points and O-points appear and disappear because the turbulent perturbation causes random magnetic reconnection.Figures 7 and 8 also show the propagation of magnetic islands.Initially, they propagate in the negative y direction, which is the ion diamagnetic direction, because the ion temperature gradient instability propagates in the ion diamagnetic direction.Then, after t = 200, the magnetic islands propagate in the positive y direction, which is the electron diamagnetic direction.
To understand the mechanism of long-wavelength island formation described above, we examine the role of each term in the generalized Ohm's law Eq.͑4͒.The long-wavelength ͑k =1͒ component of the generalized Ohm's law Eq.͑4͒ on the neutral sheet x = 0 is considered, • ٌ ˜1, and the propagation due to the equilibrium electron diamagnetic effect, i.e., ‫ء‬e effect, respectively.The last term, J ˜1, is the resistive diffusion term.Notice that the ͓⌽ ˜1 , eq ͔ and ͓p ˜e1 , eq ͔ terms vanish on the neutral sheet, x = 0. Figure 9 shows the time evolution of the absolute values of these terms for the Harris sheet with ⌬Ј = −2.4.The E ϫ B turbulent flow mixing and the electron diamagnetic turbulent flow mixing terms are comparable with the resistive diffusion term.Hence, these turbulent flow mixing terms play a role in producing long-wavelength magnetic islands.

III. MHD SIMULATION OF MAGNETIC ISLAND FORMATION
In this section, we demonstrate that turbulence caused by short-wavelength MHD instability also produces longwavelength magnetic islands and show universality of the island formation due to turbulence in a current sheet, which is stable against spontaneous magnetic reconnection.
We carry out MHD simulations of turbulence driven by short-wavelength instabilities due to a pressure gradient in a current sheet that is stable against tearing instability.Our two-dimensional MHD model consists of the vorticity equation, the Ohm's law, and the pressure equation, where g͓p͔ = ‫ץ‬ p / ‫ץ‬y is added in order to include the interchange instability ͑see Sec.4.5 in Ref. 1͒.In these equations, ⌽, , and p are the stream function, the magnetic flux function, and the pressure, respectively.Magnetic field configuration is assumed to be a uniform current density model, d eq ͑x͒ / dx = x / L s , which represents an infinitely thick current sheet.In the x direction, 1024 uniform grids are employed, and 343 Fourier modes are employed in the y direction.Parameters are set to be = 0.02, L s / L =3/ 4, L p / L =1/ 5, = 2.35ϫ 10 −5 , =1.ϫ 10 −6 , and p = 4.5ϫ 10 −6 , where L p and L s are the pressure-gradient length and the magnetic shear length, respectively.The stability parameter of the uniform current density profile is negative, ⌬Ј = −3.4.
It is confirmed that the profile is stable against tearing instability by numerically solving the linearized version of Eqs.͑7͒-͑9͒.
Here we show that long-wavelength magnetic islands are produced by turbulence due to short-wavelength MHD instabilities at the neutral sheet.A pressure gradient across the neutral sheet causes the interchange instability, which is a generalization of Rayleigh-Taylor instability.The instability is controlled by the parameter , i.e., the larger is, the more unstable equilibrium is.Initially, the instability grows up exponentially, and then a quasi-steady state is formed.Notice that the instability does not cause magnetic reconnection during the exponential growth.
Figure 10 shows color map of vorticity and equicontours of magnetic flux in the region −L / 4 Յ x Յ L / 4. Small-scale peaks of vorticity are merging, and small-scale magnetic islands appear at t = 262.Then, peaks of vorticity continue to merge, and magnetic islands grow at t = 438.The magnetic islands are twisted because of the interchange nature of the instability.We clearly see long-wavelength magnetic islands appearing on the neutral sheet x =0 at t = 700 and 875 even if there is no tearing instability.The length of the magnetic island is the same order as the system size.The shape of the magnetic island is different from that of tearing instability.Figure 11 shows energy spectrum of magnetic fluctuations.The energy spectrum has a peak during the exponential growth of the instability at t = 176.After the saturation of the instability t Ͼ 400, the energy of the instability spreads in the Fourier space via nonlinear mode coupling, i.e., the energy transfer from the typical scale of the instability excites both long and short wavelengths Fourier modes, and thus energy spectrum becomes broad.
f͔, and v di • ٌٌf =−ٌ͓f , p i ͔ are written in terms of the Poisson bracket ͓f , g͔ = ‫ץ‬ x f‫ץ‬ y g − ‫ץ‬ x g‫ץ‬ y f, where ٌ = ‫ץ͑‬ x , ‫ץ‬ y ͒ and B = B 0 e z − e z ϫ ٌ, v = e z ϫ ٌ⌽, and v di = e z ϫ ٌp i are the magnetic field, the E ϫ B flow velocity, and the ion diamagnetic velocity, respectively.

2 FIG. 1 .
FIG.1.͑Color online͒ Color map of electrostatic potential and equicontours of magnetic flux that indicates long-wavelength magnetic islands even if there is no spontaneous magnetic reconnection, ⌬Ј = −0.76.

FIG. 4 . 5 FIG. 5 .
FIG.4.Magnetic island width W as a function of the system size L for the case ⌬Ј = −3.4 and i = 3.5.
subscript 1 stands for k = 1 component.The left-hand side of the equation represents the production of longwavelength ͑k =1͒ magnetic flux on the neutral sheet.The first term in the right-hand side ͑v ˜• ٌ ˜͒1 = ͓⌽ ˜, ˜͔1 represents turbulent mixing by E ϫ B flow, where v ˜= e z ϫ ٌ⌽ ˜.The second term represents turbulent mixing by electron diamagnetic flow ͑v ˜de • ٌ ˜͒1 = ͓p ˜e , ˜͔1 , where v ˜de = e z ϫ ٌp ˜e.The third and fourth terms represent the propagation due to zonal flow, v ˜0

2 FIG. 10 .
FIG.10.͑Color online͒ Color map of vorticity and equicontours of magnetic flux that indicates long-wavelength magnetic islands even if there is no spontaneous magnetic reconnection, ⌬Ј = −3.4.
scale fluctuation, which is reported in Refs.29 and 30.The simulation box size is ͑x , y͒ = ͑L ,2L͒.In the x direction, 1024 uniform grids are employed, and 126 or 243 Fourier modes are employed in the y direction.Parameters are set to be ␤ = 0.01, L / i = 80, and = 2.35ϫ 10 −5 .