Auroral vortex street formed by the magnetosphere-ionosphere coupling instability

By performing three-dimensional nonlinear MHD simulations including Alfven eigenmode perturbations most unstable to the ionospheric feedback effects, we reproduced the auroral vortex street that often appears just before substorm onset. We found that an initially placed arc splits, intensifies, and rapidly deforms into a vortex street. We also found that there is a critical convection electric field for growth of the Alfven eigenmodes. The vortex street is shown to be a consequence of coupling between the magnetospheric Alfven waves carrying field-aligned currents and the ionospheric density waves driven by Pedersen/Hall currents.


Introduction
The problem of substorm onset has occupied the literature on solar-terrestrial physics for the past 50 years since Akasofu (1964), and the current understanding, as established by high-resolution ground and satellite optical observations (Donovan et al., 2006;Sakaguchi et al., 2009;Henderson, 2009), is that the auroral arc initially deforms into a vortex street on the scale of 30-70 km. It is clearly observed that the vortex street originates from a preexisting or newly produced arc that intensifies in ≈ 1 min and expands poleward over the course of 2 to 3 min (Lyons et al., 2002;Mende et al., 2009).
The vortex street has been interpreted in terms of instabilities in the plasma sheet, e.g., shear flow and ballooning instabilities (Voronkov et al., 1999). However, this interpretation is only supported by the faint expectation that factors affecting strong magnetic or pressure fluctuations come from the external domain. On the other hand, multiple satellite observations (Ohtani et al., 2002) have suggested that such a situation is not necessary and there is an alternative means of arc intensification. A simple scenario seems to be that an arc lying on a local field line becomes destabilized through changes in the global conditions, leading to connection to magnetotail plasma instabilities (cf. Haerendel, 2010;Henderson, 2009).
If the scenario is limited to explaining only deformation of the arc and not its poleward expansion, it can be viewed as a nonlinear evolution of shear Alfvén waves in the magnetosphere-ionosphere (MI) coupling with nonuniform active field lines. Two-dimensional magnetohydrodynamic (MHD) simulations of the feedback instability of Alfvén waves were performed in a local dipole flux tube by Lysak and Song (2002). The instability arises spontaneously due to a coupling between Alfvén waves and the electrostatic density waves in the ionosphere under a strong convection flow. More recently, there have been theoretical works that focused on the formation of field-aligned electric fields, ionospheric cavity resonances, accompanying density depletion, and filamentation of quiet arcs (e.g., Lysak and Song, 2008;Streltsov et al., 2012;Hasegawa et al., 2013). A three-dimensional simulation of feedback instability in the system indicated that the waves induce strong magnetic and flow shears to produce vortex structures around the magnetic equator (Watanabe, 2010). Various linear eigenmodes from low-frequency field line resonances to high-frequency ionospheric Alfvén resonances have been shown to become destabilized in a dipole magnetic field (Hiraki and Watanabe, 2011;hereafter, HW2011). These predictions are partially supported by evidence of the enhancement of the convection flow before substorm onset (Bristow and Jensen, 2007).

Y. Hiraki: Auroral vortex street formed by the MI coupling instability
In this study, we performed three-dimensional simulations of shear Alfvén waves in a full field line system with MI coupling, including an east-west-aligned arc. We report new results: (i) the initial arc splits and quickly deforms into a vortex street, (ii) there is a critical convection electric field for its growth, and (iii) the extent of arc intensification is controlled by the scale size of vortices. Unlike the previous studies on arc evolution starting from arbitrary setups (Lysak and Song, 2008;Streltsov et al., 2012), we solve equations describing the nonlinear evolution of the most unstable Alfvén eigenmode perturbations intrinsic in the full field line and show that arc deformation is a consequence of their growth. Starting from comprehensive studies of vortex street formation in the simplest system, we adopted a dipole magnetic field and dropped processes related to sharp v A cavities, two fluid effects, and field-aligned electric fields. Note that processes in the high-β plasma sheet are beyond the range of an approach based on magnetic perturbations.

Model description
In order to elucidate the physics involved in auroral structures, nonlinear evolution of shear Alfvén waves propagating along the dipole magnetic field B 0 can be modeled by using two-field reduced MHD equations (e.g., Chmyrev et al., 1988;Lysak and Song, 2002). The waves slightly slip ( /k ⊥ v A ) through the feedback coupling to density waves at the ionosphere. The system of interest is a field line of L ≈ 8.5 with a length of l ≈ 7 × 10 4 km and at a latitude of 70 • , where auroral arcs develop; note that it corresponds to the lower latitude in the tail magnetic field geometry. The field line position s is defined as s = 0 at the ionosphere and s = l at the magnetic equator. We set a local flux tube: a square of (l ⊥ × l ⊥ ) with l ⊥ = 10 −3 l ≈ 70 km at s = 0, a rectangle of (h ν l ⊥ × h ϕ l ⊥ ) at s, and (≈ 3300 km × ≈ 1700 km) at s = l using dipole metrics h ν (s) and h ϕ (s) with B 0 (s) = 1/h ν h ϕ (HW2011).
The electric field E is partitioned into a background convective part E 0 (⊥ B 0 ) and the Alfvénic perturbation The magnetic perturbation is expressed as B 1 = ∇ ⊥ ψ × B 0 . The equations at 0 < s ≤ l are written as The convective drift velocity v 0 = E 0 × B 0 /B 2 0 is set so that E 0 satisfies the equipotential condition, while v ⊥ = v 0 + v 1 (E 1 ), vorticity ω = ∇ 2 ⊥ φ, field-aligned current j = −∇ 2 ⊥ ψ, and ∇ = ∂ s + b 0 · ∇ ⊥ × ∇ ⊥ ψ. Suppose that the changes in the shape of the auroral arc are realized through changes in the variables (ω, j ) since E and its electron acceleration are dropped.
Ionospheric plasma motion including density waves is described by the two fluid equations. Considering the current dynamo layer (height of 100-150 km), we can assume that ions and electrons respectively yield the Pedersen drift v i = µ P E − D∇ ⊥ ln n i and the Hall drift v e = µ H E × B 0 /B 0 − j /en e , with µ P,H : mobilities and D: molecular diffusion coefficient. By integrating the continuity equations over the dynamo layer, the equations at s = 0 become (see HW2011 for details) Here, Rn e is a linearized recombination term, and the Hall mobility is normalized to be unity. We assume that j is carried by thermal electrons. Equation (4) includes the nonlinearity of the Pedersen and Hall current divergences. Also, diffusion and recombination in Eqs. (3) and (4) have an effective role in reducing the growth rates of Alfvén eigenmodes. We used the fourth-order central difference method in space and fourth-order Runge-Kutta-Gill method in time to solve Eqs. (1)-(4). The number of grids were (256,256,128) for the ν, ϕ, and s directions, respectively. The time resolution was changed in accord with the Courant condition: max(v 1 / x(s)) t < 0.25. The numerical viscosity ν v and resistivity η equaled 1 × 10 −7 /B 0 (s). Regarding the calculation domain x ⊥ (s = 0) ≡ [x, y], x and y pointed southward and eastward, respectively, in the Southern Hemisphere. We set a periodic boundary in the x ⊥ direction, e.g., at x, y = 0 and l ⊥ = 70 km (thus x ≈ 0.27 km) at the ionosphere s = 0; this is valid since we take a local flux tube approximation (L = const). It was confirmed that the boundary conditions do not affect the development of the initial arc shown in Sect. 3. An asymmetric boundary for the magnetic field ψ = 0 (or j = 0) was set at the magnetic equator s = l. At the ionospheric boundary of φ, Eq. (4) was solved using the multigrid-BiCGSTAB method.
For characteristic scales, the Alfvén velocity and transit time are set to be v A ≈ 1.5 × 10 3 km s −1 and τ A = l/v A ≈ 47 s, while l ⊥ ≈ 70 km, the magnetic field B 0 ≈ 5.7×10 −5 T, and the electron density n ≈ 3.8 × 10 3 cm −3 are values at the ionosphere s = 0. The drift velocity is v ⊥ = v A l ⊥ / l ≈ 1.5 km s −1 . We set v A = v A along s by using the dipole field B 0 (s) and a density profile n 0 (s). The ambient density at s = 0 is set to be n 0 = 10n ; note that the above v A was determined using the F-layer density (≈ 7 × 10 5 cm −3 ) and does not necessarily match this n 0 . The other values are the same as in Hiraki (2013): µ P /µ H = 0.5, P / A = n 0 µ P = 5, D = 4 × 10 5 m 2 s −1 , and R = 2 × 10 −3 s −1 .
We solved a linearized set of Eqs.
(1)-(4) to determine the eigenfunctions ( φ(s), ψ(s), n e (0)) and frequency of Alfvén waves as functions of the perpendicular wave number k ⊥ and the field-line harmonic number. For the above setting, we found that the fundamental mode with a frequency range of τ A ≈ π 2 − π has the maximum growth rate γ ≡ Im( )τ A /π . By fixing k x /2π = 1 (hereafter, k x = 1) that matches the arc form, the modes with γ switch from (1, 0), (2, 0), and (3, 0) as functions of the convection electric field E 0 ; γ = 0.1 corresponds to the timescale of τ ≈ 2.5 min.
(k x , k y ) = (1, 5) at E 0 = 20 mV m −1 to (1, 2) at 80 mV m −1 as shown in Fig. 1; the drift velocity of the feedback unstable mode, less than the convection speed v 0 , is proportional to the product of k ⊥ · E 0 , and preferable k ⊥ = (k x , k y ) is selected for finite E 0 (see HW2011). Here, the convection electric field E 0 is assumed to point poleward, so that the Pedersen current also points poleward (x), and the Hall current eastward (y). We assume a situation of the equatorward side of auroral oval where the background two-cell convection (east/westward E 0 ) deforms and the poleward E 0 dominates at the pre-midnight region. Note that, although there is no background current j 0 and E 0 is uniform in the perpendicular direction, the free energy stored as the convection drift v 0 causes the feedback instability (see also HW2011). We will study the E 0 dependence for the case of (k x , k y ) = (1, 5) in Sect. 3. We will discuss the case of (k x , k y ) = (1, 2) in Sect. 4 as well as (1, 0), (2, 0), and (3, 0) modes shown in Fig. 1. We performed a 3-D simulation to ascertain the growth of feedback eigenmodes ( φ, ψ, n e ), from an initially eastwest-aligned auroral arc, in the poleward convection field E 0 . Here, we make a note of the usage of "arc". Our MHD model did not treat the field-aligned electric field and electron energization that work for the luminosity of the real arc system. Application of our setting to the observed arc deformation, ignoring field-aligned electron acceleration and the related source process, is discussed in Sect. 4. The arc we treat in this paper still means an arc-like structure placed at the MI coupling system. The perturbed potentials and density are partitioned into the arc component and the feedback eigenmode with k ⊥ shown above, as (φ, ψ, n e ) = (φ a , ψ a , n ea ) + ( φ, ψ, n e ). The arc potential is yielded as the fundamental wave form of φ a (s) ∝ 1 B 0 (s) sin( π 2l s) while ψ a (s) = n ea = 0 for simplicity. The essence of our results was unchanged for the choice of ψ a and n ea since these quickly adjust to φ a . The perpendicular function of φ a is assumed to be Gaussian-like with l a ≈ 10 km and E a ≈ 20 mV m −1 at s = 0 (see Fig. 2); note that negative ω is quickly produced. The electric field points equatorward at the poleward edge and is reversed at the equatorward side. It is accompanied by a counterclockwise flow shear across the arc, though it is too weak to trigger some instability. The feedback eigenmode has an amplitude of | φ| = | ψ| = 10 −4 |φ a | at t = 0. Figure 2 shows the temporal variation in vorticity ω (t/τ A = 0.1, 4, 5.8, and 7) at the ionosphere s = 0 in the case of the feedback mode (k x , k y ) = (1, 5) and the poleward field E 0 = 60 mV m −1 . Note that since the results are shown in the frame of a convection drift v 0 , structures mainly move westward in the rest frame. As density waves related to the Pedersen current j P x propagate poleward, a new arc is produced by splitting through a current divergence between the wave and the initial arc. Once these arcs become dark, a vortex street forms in the poleward arc at t ≈ 3 min (panel b).

Results
Another vortex street develops from the equatorward arc (panel c) and expands into larger vortices (30-40 km) during the next ≈ 2 min (panel d). The k y mode related to the Hall current causes formation of the vortex streets propagating in the y direction. Since the density wave propagates to the direction of ionospheric currents (j P , j H ), the vortex street is distorted (asymmetric) in the x-y plane; it originates from the setting of µ P /µ H = 0.5 in this paper. The amplitude of the flow is max|v 1 | ≈ 0.32 km s −1 in panel b and grows to ≈ 0.75 km s −1 in panel d. In this case, the convection flow is v 0 ≈ 1.1 km s −1 . It is also clear that the flow in panel b is counterclockwise ( j H ∼ −v 0y ) at the poleward edge of the vortices and at the front of the j P -density waves. On the other hand, in panels c and d, a clockwise flow is added by negative ω. It is clear that the vortex street can be produced from an arc under a high E 0 without any background shear. Figure 3 shows the temporal variation in j (panels a-d) and n e (panels e-h) at s = 0 accompanied by ω in Fig. 2. The amplitude of the current is max(j > 0) ≈ 3.1 µA m −2 in panel b and grows to ≈ 30 µA m −2 in panel d (cf. Haerendel, 2010). We can easily find that j is almost out of phase with ω in the linear stage of feedback instability. Downward currents j < 0 are produced just on a vortex street at t/τ A = 4, connecting to the equatorward Pedersen current, and upward currents are induced at the equatorward arc. The vortex street with j < 0 (electron loss) quickly fades out. The equatorward arc deforms into the another vortex street with j > 0 (electron supply) at t/τ A = 5.8. It is noticeable that a new arc with j > 0 appears at the initial position (x ≈ 35 km) and the fragments are involved in the vortex street. This compound structure may correspond to the observed bead-like structure if an upward j represents the auroral luminosity.  60-70 km). The electron density is depleted by up to −30 % due to the net j < 0 at the poleward side, while a tongue-like structure of density enhancement forms at the equatorward side.
In order to help comprehend the underlying physics, in Fig. 4 we show the temporal variation in j P (panels a-d) and j H (panels e-h) along with j at s = 0. Plots are limited at (x, y) = (30-70 km, 35-70 km) for clearness, and current vectors have a relative intensity excluding the background components (j P0 , j H0 ) ∝ n 0 E 0 . There is a simple pattern at t/τ A = 0.1 where j P is converged at a placed arc (j > 0), while a clockwise j H flows on the both sides. When a vortex street forms at t/τ A = 4, the perturbed currents dominate over the arc-induced ones. The equatorward j P flows from the vortex street (j < 0) side to the arc (j > 0). A counterclockwise j H is also produced at the vortex streets (e.g., (x, y) = (45-50 km, 55-65 km)) due to their j < 0. After a new vortex street (j < 0) develops at the equatorward side of the arc, the polarized j H (−x, −y) is enhanced at t/τ A = 5.8 to produce strong j > 0 on the right side of each j < 0 through an increase in perturbed n e < 0 as shown in Fig. 3g.
The polarized j P (−x, y) is also enhanced in the middle of these pair currents j < 0 and > 0, e.g., (x, y) = (46 km, 61 km). Localized pair currents of j induce strong polarized j P and j H , e.g., at (x, y) = (55 km, 61 km) at t/τ A = 7, which is a characteristic of vortex street expansion. Figure 5 shows the field line distribution of the average vorticity ω (s) during certain periods along with its cross section at s = l; max|v 1 | ≈ 31 km s −1 at t/τ A = 7. The Alfvén wave propagates to the ionosphere at t/τ A = 1 but goes away at t/τ A = 2.5, which means the maximum ω at s = l. Some waves still remain at s = 0. Waves come back again to s = 0 by t/τ A = 4. We suppose that the apparent wave propagation time becomes longer (> 1) because the initial function is deformed, which means generation of a new wave on the way. The amplitude at s = 4-9 R E decreases while the vortex street forms during this period. Waves return to the magnetosphere during t/τ A = 4-5, and the amplitudes in the region of s > 3 R E increase. Although partial reflections continuously occur, the waves (or ω max ) are on the s = 0 side at t/τ A = 6 and on the s = l side at t/τ A = 7. The nonlinear coupling to the fast density waves (large E 0 ) Ann. Geophys., 33, 217-224  causes a rapid growth of Alfvén waves through slippage and partial reflection, resulting in a pileup of vortices. It is also clear that the vortex pattern changes between s = 0 and l; an equipotential mapping cannot be inferred at this scale of aurora. We further see that radially inward weak flows develop at x ≈ 1700 and ≈ 3000 km through the effect of arc splitting. Figure 6 shows the average electron density n e (t) for E 0 = 20-80 mV m −1 . When the vortex street forms (see Fig. 2), it increases by up to 10-15 % of n 0 = 10n . In the case of E 0 = 20 mV m −1 , arc splitting occurs just before every minimum of n e , and new arcs repeatedly vanish. Although shears of these arcs prevent the growth of the mode (k x , k y ) = (1, 5), eventually a north-south-elongated structure (not vortex street) grows at t/τ A ≈ 27. The knowledge we get here is that the change in the growth pattern within E 0 = 20-40 mV m −1 indicates the existence of a critical field E cr for deformation of the arc into the vortex street. Further sensitivity studies (not shown) confirmed that the results related to E cr are independent of the initial conditions of arc itself, i.e., width l a of 10-20 km, field E a of 20-80 mV m −1 , and polarity sgn(φ a ).

Discussion
Let us first mention an overarching problem with regard to the application of our modeled system to the real world where auroral mechanisms, i.e., field-aligned electron acceleration and the ionization by precipitating electrons, are at work from the first brightening of the onset arc. We suppose that these are generally important in the arc system but do not directly contribute to the arc deformation itself. The arc deformation could start from destabilization of the Alfvén waves, and field-aligned currents are amplified to be 4-20 µA m −2 in their developing stage as shown in Fig. 3. The high current density can yield a strong field-aligned voltage if some resistivity as electron inertia and the realistic deep cavities are included in the system. As our first step in a minimal model setup without these features, we addressed the pure nonlinear coupling between an arbitrarily placed arc structure and Alfvén eigenmodes. At the next step, we will investigate the generation of field-aligned electric fields, the related electron acceleration, and ionospheric source processes in our nonlinear MI feedback system. Also, we should address the relationship between j > 0 at s = 0 and j and E at the electron acceleration region. Since the assumption was made that field-aligned currents are carried only by thermal electrons at s = 0, regions of the upward current j > 0 at s = 0 do not necessarily represent the auroral luminosity in a rigorous sense; the associated electrons only flow up magnetic fluxes but cannot form the auroral luminosity.
It is noted that our initial setup of the auroral arc is not artificial. Although the imposed arc in the case of Fig. 2 is highly unstable, it is demonstrated to be stable under a weak convection electric field (e.g., E 0 = 20 mV m −1 ) as in Fig. 6 and discussions below. Our model captured the transition from stable to unstable regime of the arc system with feedback eigenmodes. The scenario is that the arc appears a long time (e.g., 30 min) before substorm onset under a weak (westward) E 0 , whose source is at the magnetospheric end. After the poleward component of E 0 is enhanced up to a critical level (see the next paragraph), the arc rapidly grows as shown in Figs. 2 and 3.
Let us describe the changes in the wave growth patterns in the convection electric field of E 0 = 20-40 mV m −1 (Fig. 6) in the context of arc splitting. The behavior can be understood by the linear growth rate of feedback eigenmodes in Fig. 1 because they certainly grow in spite of oscillatory motions of n e through Alfvén wave propagation (see also Fig. 5). Figure 1 shows that γ of the (2, 0) mode switches from negative to positive at 20-40 mV m −1 . On the other hand, γ of the (1, 5) mode increases by ≈ 1.5 times in this range, but decreases in the higher regime; as these tendencies are inconsistent with the changes in growth patterns, this mode is not considered to be the main trigger. γ of the (1, 2) mode monotonously increases, which affects its saturation level, but it is not directly related to the changes. Consequently, we can infer the underlying physics of the vortex street formation as follows. The initial arc splits to generate the (2, 0) mode. Alfvén wave amplification follows the growth of the (2, 0) mode, which cumulatively grows with the (1, 5) mode, forming into the vortex street. The critical field E cr ≈ 25 mV m −1 implies the transition of the poleward mode growth rate (i.e., γ of the (2, 0) mode in this case).
An increase in n e in the case of E 0 = 20 mV m −1 (Fig. 6) means that eastward modes (k y = 0) can grow, though it takes a long time (≈ 25 min), after a few arc splittings; remember that the mode is defined in the westward-moving frame of v 0 . However, these modes are decoupled from the southward modes during this growth, which is not consistent with observation of the arc changing into the vortex street (e.g., Sakaguchi et al., 2009). We conclude that this case is a product of certain ideal conditions, and the critical field for vortex formation is the value given in the previous paragraph.
Let us see to what extent the vortex street can brighten. Figure 7 shows the results for the feedback mode (k x , k y ) = (1, 2) and E 0 = 60 mV m −1 . As in the case of (k x , k y ) = (1, 5) (Figs. 2 and 3), the initial arc repeatedly splits until t/τ A ≈ 5, brightens, and deforms into a twin vortex after that. For the case of the (1, 5) mode, a vortex street appearing at the poleward arc fades out, but that from the equatorward arc develops. On the other hand, a vortex street ap-pearing at the poleward arc directly develops in this case. Strong upward j is produced at the poleward right edge of the vortices (not shown). It seems that the arc from which the vortex street forms depends on the growing eigenmodes; both cases were really observed (Lyons et al., 2002;Mende et al., 2009). It is also emphasized that, in every case of (1, 2) with E 0 ≥ 40 mV m −1 , the saturation level of n e shifts to a higher value, up to ≈ 38 % from n 0 , than the level of the (1, 5) cases. This behavior is not simply understood by the linear growth rates, i.e., γ (1, 2) < γ (1, 5) at E 0 = 40 mV m −1 , depicted in Fig. 1. We can guess that the (1, 2) mode is well suited for nonlinear coupling to the (2, 0) mode. What can be inferred from this result is that the extent of arc intensification related to vortices is controlled by the growing eastward feedback eigenmodes. The modes were mainly determined by conductances µ P /µ H and P / A along with E 0 .
Let us compare our results with the observed features of auroral arcs during substorm onset. Ground-based optical observations have indicated that the arcs often flap and split just before vortex street formation (e.g., Motoba et al., 2012;K. Hosokawa and K. Sakaguchi, personal communication, 2014). These behaviors imply that ionospheric currents and density waves are increased by their feedback coupling to Alfvénic j . The vortex sizes produced by a feedback growth of the (1, 5) and (1, 2) modes are ≈ 14 and 35 km, respectively. We should note that our simulation domain is a dipole field line (L ≈ 8.5) when the vortex sizes are compared to those seen in the onset arc as 40-80 km (Sakaguchi et al., 2009, and references therein). The field line length l can be extended if an actual magnetotail field line is considered. Then, l ⊥ ∝ l is also extended, and the above difference (factors of 2-3) between the vortex size in our model results and that in observations is not the main concern for explaining the behavior of the onset arc. Another point is on the dynamic range of optical (all-sky camera) observations. In case of the above paper, power spectra in scales of < 40 km are noisy, and, further, the scale estimation strongly depends on the zenith angle of the arc position. We urge that high-resolution measurements of the onset arc be made in order to estimate its real wavelength; multi-point stereo-observations are anticipated. There is another problem related to the growth timescale. The scale of the observations, τ = 1-2 min (e.g., Mende et al., 2009), is shorter than that of our calculation, τ = 3-5 min. The Alfvén transit time we use is as large as τ A ≈ 47 s, but the growth time was shown to be shorter in more realistic cases (τ A ≈ 25 s) with ionospheric v A cavities (Hiraki and Watanabe, 2012). Therefore, simulations in the cases can explain the small τ in the observations. However, there is a crucial difference between observations (Motoba et al., 2012;Hosokawa et al., 2013) and our calculations. The observed phase speed of bead structures is ≈ 3 km s −1 , but the calculated one is ≈ 0.75 km s −1 ; see the first paragraph of Sect. 3. We speculated that the difference by a factor of 4 may be reduced through inclusion of the enhancement of ionospheric conductivities. The ionization by precipitating electrons causes the enhancement by several factors and the amplitude of the triggered density waves may increase, and then the drift speed may also increase. This is really the next problem to be clarified in the series of our studies.
The next point is on the magnitude of E 0 . The critical field E cr ≈ 25 mV m −1 was determined on the basis of linear properties of Alfvén eigenmodes as mentioned above. From our previous studies (Hiraki andWatanabe, 2011, 2012), the drift frequency (∝ k ⊥ · E 0 ) with maximum γ does not so much depend on the profiles of B 0 (s) and v A (s) in cases of the Alfvén transit time τ A = 25 and 47 s. Thus, the value of E cr cannot be also largely changed when realistic B 0 and v A cavities are taken into account. The other factor that can affect E cr is ionization by nonthermal electrons. If we can consider that the term linearly depends on the field-aligned current, it does not largely change the wave linear properties and E cr . But, if it has a more steep function given by a product of E and j (cf. Lysak and Song, 2002), the value of E cr can be reduced; this should be examined in our next study like the problem mentioned just above. On conditions of some uncertainties in theoretical studies, the critical field E cr ≈ 25 mV m −1 (i.e., a flow of 0.43 km s −1 ) we present here agrees fairly well with the enhanced convection flows (≥ 0.5 km s −1 ) observed before onset (Bristow and Jensen, 2007). Further comparison provides a useful restriction for interpreting statistical data, as well as the direction of E 0 assumed to be poleward in this paper. We urge that detailed measurements of pre-conditioning of onset arcs should be undertaken in the future.
The results in this paper indicate that feedback modes are prevented from growing by the initial arc when E 0 < E cr but that they grow together into a vortex street when E cr < E 0 . The vortex street is a manifestation of the nonlinear evolution of Alfvén waves in the MI coupling system. The following problems remain to be solved: (i) how are waves trapped in the ionospheric cavity region to form a field-aligned electric field when the field-aligned current is amplified, and (ii) why does the active arc expand poleward in the next step? The second question is closely related to coupling with magnetotail high-β plasma dynamics (e.g., Henderson, 2009).

Conclusion
By performing three-dimensional MHD simulations including Alfvén eigenmode perturbations most unstable to the ionospheric feedback effects, we examined the auroral vortex street that often appears just before substorm onset. We found that (i) the initial arc splits, intensifies, and deforms into a vortex street through Alfvén wave amplification during their 2-3 bounce periods (3-5 min), (ii) the vortex street is characterized by an enhancement of polarized Pedersen/Hall currents j P,H due to localized pairs of field-aligned currents j > 0 and < 0, (iii) there is a critical convection electric field 224 Y. Hiraki: Auroral vortex street formed by the MI coupling instability of E cr (≈ 25 mV m −1 in the present setting) for growth of Alfvén eigenmodes, and (iv) the extent of arc intensification is controlled by the nonlinear behavior of finite k y modes with the arc-related k x modes. The results of our simulation indicate that the vortex street is a consequence of coupling between the shear Alfvén waves carrying j and the ionospheric density waves driven by j P,H .
Topical Editor L. Blomberg thanks the three anonymous referees for their help in evaluating this paper.