Bifurcation in electrostatic resistive drift wave turbulence
Abstract
The Hasegawa–Wakatani equations, coupling plasma density and electrostatic potential through an approximation to the physics of parallel electron motions, are a simple model that describes resistive drift wave turbulence. We present numerical analyses of bifurcation phenomena in the model that provide new insights into the interactions between turbulence and zonal flows in the tokamak plasma edge region. The simulation results show a regime where, after an initial transient, drift wave turbulence is suppressed through zonal flow generation. As a parameter controlling the strength of the turbulence is tuned, this zonal flow dominated state is rapidly destroyed and a turbulence-dominated state re-emerges. The transition is explained in terms of the Kelvin-Helmholtz stability of zonal flows. This is the first observation of an upshift of turbulence onset in the resistive drift wave system, which is analogous to the well-known Dimits shift in turbulence driven by ion temperature gradients.
I Introduction
Fusion plasmas and other turbulent flows in quasi-two-dimensional (2D) geometry can undergo spontaneous transitions to a turbulence-suppressed regime. In plasmas they are known as L–H (low-to-high confinement) transitions and are studied intensively because they effectively enhance the confinement, through suppression of anomalous or turbulent particle and heat fluxes. It is now widely accepted that emergent zonal flows are crucial to achieving confinement improvement Diamond et al. (2005). The L–H transition is associated with nonlinearly self-generated poloidal shear or zonal flows tra () in the tokamak edge region, which comprises the transition zone from inner hot core plasma to the outer cold scrape-off layer. Zonal flows reduce anomalous transport by absorbing energy from drift waves and by shearing apart eddies which mediate turbulent transport, and thus play a key role in its regulation.
In this paper we present the results of analytic and numerical investigations of transitions between turbulence-dominated and zonal-flow-dominated regimes, using the Hasegawa–Wakatani (HW) model Hasegawa and Wakatani (1983); Wakatani and Hasegawa (1984) for electrostatic resistive drift wave turbulence in 2D slab geometry. We find that bifurcations in the model correspond to the onset of drift wave turbulence, the generation of zonal flows, and the re-emergence of turbulence as the zonal flows become unstable, and observe that this is drift wave turbulence analog of the Dimits shift Dimits et al. (2000) in ion temperature gradient turbulence.
Three energetic subsystems interact to produce the complexity observed in L–H transition dynamics: the kinetic energy of turbulence, the kinetic energy of shear flows, and the potential energy contained in density or pressure gradients. The three major governing processes are generation of turbulence by drift waves, self-organization of zonal flows, and destabilization of the zonal flows. The instabilities that lead to these changes correspond to bifurcations of equilibrium solutions of model equations. If a tunable parameter crosses a stability threshold the qualitative nature of the solution changes. We say that a primary instability occurs at a linear stability threshold of the equilibrium with zero background flow, which physically corresponds to the onset and growth of drift waves. Theoretical Biskamp and Kaifen (1985) and experimental Klinger et al. (1997) studies have indicated that the generation of drift wave turbulence in plasmas may occur by the Ruelle-Takens mechanism Ruelle and Takens (1971), in which a limit cycle generated by a Hopf bifurcation undergoes a Niemark-Sacker bifurcation to a torus, which may undergo one or more bifurcations to higher-dimensional tori before the motion becomes chaotic.
However, to complicate this generic turbulence onset scenario, in plasmas zonal flows will be generated beyond the primary threshold due to an instability of the drift waves, effectively suppressing drift wave activity. This instability causing the zonal flow onset is termed a secondary instability. We can consider the turbulence to be well-developed at the secondary instability; i.e., for heuristic purposes we assume the Ruelle-Takens sequence to have already occurred.
A strong candidate for this secondary instability mechanism is modulational instability Guzdar et al. (2001); Dewar and Abdullatif (2007), a special case of nonlinear mode coupling whereby modulation of a small scale monochromatic wave can transfer energy non-locally to a longer wavelength structure due to the ponderomotive force effect leading to excitation of zonal flows. One might also expect an inverse energy cascade, endemic to quasi two-dimensional flows in general, whereby local mode coupling channels energy into large scale structures.
A different mechanism for this secondary instability that generates zonal flows is Kelvin–Helmholtz (KH) instability Rogers et al. (2000); Jenko et al. (2000). In this scenario the KH instability may be driven by radially elongated drift wave eigenmodes. The KH mode of the drift waves necessarily possess a zonal flow component, and provide a natural mechanism for the zonal flow growth.
As the zonal flows become more energetic they are subject to tertiary instability which breaks up the coherent zonal structuring of the flow into turbulent small scale eddies via KH instabilities of the zonal flows. The small scale turbulence may again cohere via secondary instabilities. These interactions are schematized in Fig. 1.
Nonlinear interactions between zonal flows and drift waves results in an upshift of the boundary in parameter space for the tertiary onset of turbulence. This is known as the Dimits shift in ion-temperature-gradient (ITG) driven turbulence, and the turbulence suppressed regime was mapped by gyrokinetic and gyrofluid simulations Dimits et al. (2000).
The simplest approach that captures the essential physics underlying the problem is low-dimensional dynamical modeling and analysis Sugama and Horton (1995); Ball et al. (2002); Ball (2005); Kolesnikov and Krommes (2005), which can provide a very economical tool to predict the transition. However, the tradeoff with such highly coarse-grained modeling is that it necessarily whites out information, and may therefore miss important physics and predict unphysical singular behavior Ball (2005). Thus we require validation of the low-dimensional modeling results by computational simulations of finer models.
The HW model Hasegawa and Wakatani (1983); Wakatani and Hasegawa (1984) was developed to investigate anomalous edge transport due to collisional drift waves, and has been widely studied Hasegawa and Wakatani (1987); Horton (1999); Pedersen et al. (1996); Camargo et al. (1995); Gang et al. (1989). It includes the effects of inhomogeneous background density and parallel electron dynamics described by Ohm’s law. The density gradient drives the drift waves, which are destabilized by the parallel electron resistivity. Convective nonlinearity regulates the linear growth of the resistive drift wave instability, and a quasi-stationary state is achieved where the resistive coupling balances the input. The HW model is particularly simple yet includes the essential physics for studying the self-consistent generation of turbulence and growth and decay of coherent macroscopic structures such as zonal flows Hasegawa and Wakatani (1987), even though it does not describe physics that can be important in specific situations, such as magnetic curvature, magnetic shear, and electromagnetic effects.
We emphasize that the parallel electron motion is important for generation, stabilization, and destabilization of zonal flows. The parallel electron response given by the generalized Ohm’s law leads to resistive coupling between the electrostatic potential and the density fluctuations. In toroidal geometry this coupling does not act on the flux-averaged parts Dorland and Hammett (1993), and in the original or unmodified HW model we do not observe zonal flows. Modification of the resistive coupling term, described in Sec. II, enables the generation of zonal flows. This corresponds to the difference between the ITG and the ETG (electron-temperature-gradient) cases discussed by Jenko et al. Jenko et al. (2000), who found that suppression of the secondary KH instability in the ETG case, due to the adiabatic electron response, is removed in the ITG limit.
In Sec. II, we describe the HW model and discuss the treatment of parallel electron motions. Linear stability analysis of the zero-flow background is also given to calculate transition points in parameter space. Numerical simulation results are given in Sec. III. We carry out a systematic parameter survey to locate the transition from a zonal-flow-dominated state to a turbulent state. To examine the hypothesis that this transition may be ascribed to the tertiary KH instability of the zonal flow, we study the KH stability of the generated zonal flows in the HW model in Sec. IV and compare the KH stability threshold with the transition boundary determined by simulation. Discussions and conclusions are presented in Sec. V.
Ii Modified Hasegawa–Wakatani Model
The physical setting of the HW model may be considered as the edge region of a tokamak plasma of nonuniform density and in a constant equilibrium magnetic field . Following the drift wave ordering Hasegawa and Mima (1977), the ion vorticity ( is the electrostatic potential, is the 2D Laplacian) and the density fluctuations are governed by the equations
(1) | ||||
(2) |
where is the Poisson bracket, is the dissipation coefficient. The background density is assumed to have an unchanging exponential profile: . Electron parallel motion is determined by Ohm’s law with electron pressure ,
(3) |
assuming electron temperature to be constant (isothermal electron fluid). This relation gives the coupling between and through the adiabaticity operator appearing in Eqs. (1) and (2). In our 2D setting becomes a constant coefficient when acting on the drift wave components of and by the replacement , where is a length characteristic of the drift waves’ phase variation along the field lines. However, for the zonal flow components, this resistive coupling term must be treated carefully because zonal components of fluctuations ( modes) do not contribute to the parallel current Dorland and Hammett (1993). Recalling that turbulence in the tokamak edge region, where there is strong magnetic shear, is considered here, should always coincide with because any potential fluctuation on the flux surface is neutralized by parallel electron motion. Let us define zonal and non-zonal components of a variable as
where is the periodic length in , and remove the contribution by the zonal components in the resistive coupling term in Eqs. (1) and (2). Subtraction of the zonal components from the resistive coupling term yields the modified HW (MHW) equations,
(4) | ||||
(5) |
Evolutions of the zonal components can be extracted from Eqs. (4) and (5) by averaging in the direction:
where stands for and .
Wakatani and Hasegawa found Wakatani and Hasegawa (1984) that excitations of waves having that maximizes the linear growth rate (for given and ) are most likely to occur, since the plasma can choose any parallel wavenumber (). Using the parallel wave number of the maximum growth rate, is given by . This also gives for the zonal mode.
The MHW model spans two limits with respect to the adiabaticity parameter . In the adiabatic limit (collisionless plasma), the non-zonal component of electron density obeys the Boltzmann relation , and the equations are reduced to the Hasegawa–Mima equation Hasegawa and Mima (1977). In the hydrodynamic limit , the equations are decoupled. The vorticity is determined by the 2D Navier-Stokes equation, and the density becomes a passive scalar. The advantage of our choice of as a free parameter is the capability for treating the limits in a unified manner.
The variables in Eqs. (4) and (5) have been normalized by
where is the ion sound Larmor radius ( is the ion sound velocity in the cold ion limit), is the fluctuating part of the density.
In the adiabatic, ideal limit (, ) the MHW system has two dynamical invariants, the energy and the potential enstrophy ,
(6) |
where , which constrain the fluid motion. Conservation laws are given by
(7) |
where fluxes and dissipations are given by
Unlike the Hasegawa–Mima model which is an energy-conserving system, the MHW model has an energy source . Due to the parallel resistivity, and can fluctuate out of phase which produces non-zero . The system can absorb free energy contained in the background density profile through the resistive drift wave instability.
Note that the same conservation laws hold for the original, unmodified original HW (OHW) model, Eqs. (1) and (2), except that is defined by both zonal and non-zonal components; . In the OHW model, the zonal modes as well as the non-zonal modes suffer resistive dissipation.
We present the linear stability analysis for the zero background (the primary instability). Beyond this stability threshold we expect excitation of drift waves. Since the zonal modes have linearly decaying solutions, we only consider the form (). Linearization of Eqs. (4) and (5) around the zero equilibrium () yields the dispersion relation,
(8) |
where we defined , , and the drift frequency . Solutions to the dispersion relation (8) are given by
, . In the limit where , it is readily proved that one of the growth rates is positive if is finite, thus unstable. However, there exists a range of where the drift wave instability is suppressed. The stability threshold is given by
(9) |
and is depicted in Fig. 2. The first unstable mode shown in the figure is the mode. Below this threshold, an initial perturbation damps out and nothing happens. If we choose the parameters in the region beyond the threshold, more than one mode starts to grow linearly until the nonlinear terms set in. The left panel shows how many modes are excited for given parameters. Most unstable modes are on axis.
Iii Simulation Results
The HW equations are solved in a doubly periodic square slab domain with box size where the lowest wavenumber (). The equations are discretized on grid points by the finite difference method. Arakawa’s method is used for evaluation of the Poisson bracket Arakawa (1966). The time stepping algorithm is the third order explicit linear multistep method Karniadakis et al. (1991). We examine the effects of the parameters and on the nonlinearly saturated state, and fix throughout this paper.
We start simulations by imposing small amplitude random perturbations. The perturbations grow linearly in the initial phase and generate drift waves, then the drift waves undergo secondary instabilities which excite zonal flows until nonlinear saturation occurs. In the saturated state, we observe that . We compare the MHW and the OHW models by showing the spatial behavior of the saturated electrostatic potential in Fig. 3, and the time evolution of the total kinetic energy, the zonal component of the kinetic energy, and the cross-field transport in Fig. 4. From Fig. 3 we see that zonally elongated structures of the electrostatic potential are generated in the MHW model, while rather isotropic vortices are generated in the OHW model. From Fig. 4 we see that growth of the drift waves is not changed by the modification, but that in the MHW model the zonal flows saturate at a higher amplitude (because the modification removes the unphysical resistive dissipation of the zonal modes). In fact, in the MHW model, the zonal flows carry nearly all the kinetic energy in the final state — they have absorbed nearly all the energy from the drift waves. In both models, the cross-field transport initially increases as the turbulent kinetic energy level increases, but in the MHW model it begins to fall as zonal flows absorb the drift wave energy. The build-up of the zonal flow in the MHW model and the resulting transport suppression highlight the importance of the difference between the MHW and the original HW model in the nonlinear regime Numata et al. (2007).
Let us show how the parameters and affect the saturated state in the MHW model. In Fig. 5, we plot the ratio of the kinetic energy of the zonal flow () to the total kinetic energy () against and . It is clearly seen that there are two types of saturated states. One is a zonal-flow-dominated state where turbulence is almost completely suppressed, and the other is an isotropic turbulence-dominated state. The zonal-flow-dominated state suddenly jumps to the turbulent state in a narrow range of the parameter space. If we strongly drive the drift wave instability by increasing , the system is likely to reach the turbulent state. From the dependence on , we can see that zonal flows are generated in the adiabatic regime () while isotropic flows are generated in the hydrodynamic regime (). These results are compatible with the properties of the Hasegawa–Mima model and of hydrodynamic flows as discussed in the next section.
Let us assume that the generated zonal flows in the direction can be expressed by a sinusoidal profile,
(10) |
The amplitude and wavenumber are determined from the simulation results. To estimate we plot the average wavenumber of the generated zonal flow
(11) |
in Fig. 6, and amplitude of the zonal flow in Fig. 7. The average wavenumbers are small and rather insensitive to the parameters. This illustrates a feature of 2D flows, which tend to generate large scale structures. The wavenumber of a stable zonal flow is typically (corresponding to ). The amplitudes of zonal flows are roughly proportional to and are independent of .
Iv Stability of Zonal Flow
We examine the stability of the zonal flows obtained from the numerical simulations, and compare the stability threshold and the transition point in this section. We consider the perturbation around the zonal flow background. The electrostatic potential and the density are decomposed as , and where gives the background flow in the direction. By linearizing the MHW equations, we obtain an eigenvalue equation containing the effect of and ,
(12) |
We neglect the viscosity. The density fluctuation is determined by
(13) |
We solve the eigenvalue equation by the standard shooting method in the domain . The boundary is assumed to be rigid for simplicity.
iv.1 Hydrodynamic and adiabatic limit
Before going to the analysis of the HW case, we briefly review the results in two limits: the hydrodynamic limit () and the adiabatic limit ().
In the limit, we recover the Rayleigh eigenvalue equation for neutral fluids,
(14) |
The well-known Rayleigh’s inflection point theorem demands existence of an inflection point for the instability Lord Rayleigh (1879). The necessary and sufficient condition is also known for this case. Tollmien Tollmien (1935) showed the existence of a marginally stable eigenfunction satisfying where is the inflection point. satisfies,
(15) |
The solution is given by
(16) |
and the critical wave number is
(17) |
If , the marginally stable wave number exists. It should be noted that Tollmien does not exclude the possibility that the marginally stable mode is isolated. However, perturbation analysis around the marginally mode shows the existence of solutions smoothly connected to the marginal solution Lin (1945); Drazin and Reid (1981).
A similar analysis can be applied to the adiabatic limit,
(18) |
if . The marginally stable eigenfunction satisfies,
(19) |
The solution is identical with the previous case, but the critical wave number is slightly modified to
(20) |
The necessary and sufficient condition of the flow shear for instability is .
We can judge the stability by finding the critical wavenumber. We consider the flow given by Eq. (10) with . The critical wavenumber exists only in the hydrodynamic limit for the given profile. On the other hand, the given flow is stable in the adiabatic limit. The difference of the two stability conditions (17) and (20) comes not from but from the strong coupling between and , and reflects the stabilizing effect of adiabatic parallel electron motion.
Figure 8 shows the imaginary parts of the eigenvalues for case in the hydrodynamic limit. The eigenvalues are pure imaginary in this limit because of antisymmetry of the flow []. Another property in this limit is the scale invariance. The eigenvalues do not depend on and , but are determined by .
We set . Or, in other words, is normalized out by considering . The eigenvalue problem of the given flow profile with has the same eigenvalues as that of the flow with in the half domain (solid line). The critical wavenumber for this curve is given by . In addition, we find another branch of solutions (broken line) which continue to exist until .
Next, let us consider the effect of . Since the critical wavenumber does not exist for the profile with , we examine a profile having stronger flow shear by setting , and take for simplicity. In this setting the marginal wavenumber exists ().
Figure. 9 shows the eigenvalues obtained in the adiabatic limit for and . is also normalized by . seems independent of . Thus the same stability condition still holds for finite, but not too large, . As we see from the figure, the growth rate decreases with increasing and disappears for large even though exists. We need another condition for . Multiplying (18) by complex conjugate of and integrating over the domain, we obtain
(21) |
If , must be satisfied somewhere in the domain Kuo (1949). Applying this condition to our assumed flow profile, we obtain the condition for the instability. This gives only a necessary condition for the instability, but provides a good estimate [Fig. 9 (b)].
If we find the eigenvalue and the corresponding eigenfunction , the complex conjugate of is also an eigenvalue and the corresponding eigenfunction is given by the complex conjugate of . Thus, we can always restrict our quest for eigenvalues in the upper half plane of the complex plane without loss of generality. This greatly simplifies the situation because we can neglect the continuous spectrum on the real axis.
iv.2 Hasegawa–Wakatani case (intermediate value of )
Unlike the previous cases, the complex conjugate of an eigenvalue is not a eigenvalue if we include finite . In this case we must solve for negative as well. Moreover, there exist two continuous spectra in this case:
(22) |
Both represent convective transport due to the background flow. One of them is damped by the resistivity. These continua may interact with the point spectrum. Thus the situation is much more complicated in the intermediate case compared with the adiabatic and hydrodynamic limits.
First, we show the effect of and neglect effect of . We consider for simplicity. Figure 10 shows the imaginary parts of the eigenvalues. Three different cases, and the dependence of the positive branches, are shown. The continuous spectra are shown by thick solid lines. In the case, two branches from the case (dotted line) are also shown for reference, so that it is seen that is slightly shifted downwards for finite . As decreases, the upper (unstable) branch intersects the continuous spectrum at marginal stability, and there exists a gap (interval in ) occupied by the two continuous spectra before this branch continues as a stable mode. The eigenfunctions belonging to the eigenvalues in the point spectrum close to this gap become singular.
For increasing , we observe the positive eigenvalues disappear at . In addition to the two stable branches seen at , at another stable branch has appeared in the small region. By further increase of we find that the lower two branches merge. Beyond this merging point, finite real parts appear, and the eigenmode starts to travel in the direction.
Our concern is to determine the stability threshold in the - plane. Next, we consider the effect of in addition to . Since always appears in the form of and is small in the vicinity of the threshold, the effect of is rather minor. does not significantly affect the behavior of the eigenvalues except that controls the amplitude of flow. As we stated earlier, the parameters are normalized by , , , in the shooting calculation, where is proportional to .
Finally, we summarize the shooting calculation by showing the bifurcation diagram in - plane together with the numerically obtained results. The only excitable mode that can be resolved in the numerical simulation is the mode, which is the first unstable mode of the primary instability (see Sec. II). In Fig. 11, we show the stability threshold of mode for the primary instability (resistive drift wave instability) and the tertiary instability (KH instability). Each mark in the figure denotes a numerically obtained saturated state: , , represent respectively the zonal-flow-dominated, transitional, and turbulence-dominated states. In these states zonal flows contain more than 90%, 20-90%, and less than 20% of the total kinetic energy, respectively. The qualitative tendency of the thresholds in the bifurcation diagram shows agreement between the numerical simulations and the KH analysis, i.e. increasing () is stabilizing (destabilizing). Zonal-flow-dominated states are observed in between the primary and the tertiary instability thresholds. The emergence of a turbulent state is shifted from the primary threshold to the tertiary threshold due to the turbulence suppression effect of the zonal flow, which is analogous to the Dimits shift observed in ITG turbulence.
The reasons for the quantitative discrepancy between the boundary of the zonal and the turbulent states may be because of the simplification made in the KH analysis; the simplified flow profile, the boundary condition and viscosity may also affect the results.
V Conclusion
In summary, we have analyzed bifurcation phenomena in two-dimensional resistive drift wave turbulence. First, we have performed numerical simulations of the modified HW model to study bifurcation structures in a two-parameter (-) space. We have shown that, in the MHW model, zonal flows are self-organized and suppress turbulence and turbulent transport over a range of parameters beyond the linear stability threshold for resistive drift waves. By performing a systematic parameter survey, we have found that such zonal-flow-dominated states suddenly disappear as a threshold is crossed, being replaced by a turbulence-dominated state.
The threshold of the onset of turbulence has been compared with the linear stability threshold of an assumed laminar zonal flow profile. Simple theoretical predictions in limiting cases explain the qualitative tendency of the stability of the zonal flow. determines the amplitude of the zonal flows, thus, large destabilizes the zonal flows. On the other hand, the adiabatic response of parallel electrons given by stabilizes them. Numerical analysis of the eigenvalue problem determining the stability of the assumed zonal flow profile in the HW model also confirms this trend. The constructed bifurcation diagram in the - plane for the HW model confirms the scenario of the onset of turbulence in the drift wave/zonal flow system being due to the disruption of zonal flows by KH instability.
The HW model considered here is a particularly simple model, but includes the essential physics of interactions between turbulence and coherent structures. This system exhibits many other interesting phenomena, but in this paper we have focused on the effect of the linear driving term and the parallel electron response (including the resistivity). To do so, we set the viscosity very small. In this case the zonal flow survives for a very long time. However, when the viscosity comes into play, the zonal flows are damped rapidly, and the turbulence grows again until zonal flows can be nonlinearly excited and the cycle repeats. Thus the system exhibits predator-prey oscillatory behavior.
Vi Acknowledgments
The authors would like to acknowledge B. D. Scott for providing the simulation code used in this work. We thank P. N. Guzdar, W. Dorland, C. Tebaldi and J. A. Krommes for useful discussions. This work is supported by the Australian Research Council.
References
- P. H. Diamond, S.-I. Itoh, K. Itoh, and T. S. Hahm, Plasma Phys. Control. Fusion 47, R35 (2005).
- The L–H transition may also be induced by externally generated flows, for example by edge biasing [see R. J. Taylor, M. L. Brown, B. D. Fried et al., Phys. Rev. Lett. 63, 2365 (1989)], but our interest here is focused on internally or self-generated shear flows.
- A. Hasegawa and M. Wakatani, Phys. Rev. Lett. 50, 682 (1983).
- M. Wakatani and A. Hasegawa, Phys. Fluids 27, 611 (1984).
- A. M. Dimits, G. Bateman, M. A. Beer, B. I. Cohen, W. Dorland, G. W. Hammett, C. Kim, J. E. Kinsey, M. Kotschenreuther, A. H. Kritz, et al., Phys. Plasmas 7, 969 (2000).
- D. Biskamp and H. Kaifen, Phys. Fluids 28, 2172 (1985).
- T. Klinger, A. Latten, A. Piel, G. Bonhomme, T. Pierre, and T. Dudok de Wit, Phys. Rev. Lett. 79, 3913 (1997).
- D. Ruelle and F. Takens, Commun. Math. Phys. 20, 167 (1971).
- P. N. Guzdar, R. G. Kleva, and L. Chen, Phys. Plasmas 8, 459 (2001).
- R. L. Dewar and R. F. Abdullatif, in Proceedings of the CSIRO/COSNet Workshop on Turbulence and Coherent Structures, Canberra, Australia, 10-13 January 2006, edited by J. P. Denier and J. S. Frederiksen (World Scientific, Singapore, 2007), vol. 6 of World Scientific Lecture Notes in Complex Systems, pp. 415–430.
- B. N. Rogers, W. Dorland, and M. Kotschenreuther, Phys. Rev. Lett. 85, 5336 (2000).
- F. Jenko, W. Dorland, M. Kotschenreuther, and B. N. Rogers, Phys. Plasmas 7, 1904 (2000).
- H. Sugama and W. Horton, Plasma Phys. Control. Fusion 37, 345 (1995).
- R. Ball, R. L.Dewar, and H. Sugama, Phys. Rev. E 66, 066408 (2002).
- R. Ball, Phys. Plasmas 12, 090904 (2005).
- R. A. Kolesnikov and J. A. Krommes, Phys. Plasmas 12, 122302 (2005).
- A. Hasegawa and M. Wakatani, Phys. Rev. Lett. 59, 1581 (1987).
- W. Horton, Rev. Mod. Phys. 71, 735 (1999).
- T. S. Pedersen, P. K. Michelsen, and J. J. Rasmussen, Plasma Phys. Control. Fusion 38, 2143 (1996).
- S. J. Camargo, D. Biskamp, and B. D. Scott, Phys. Plasmas 2, 48 (1995).
- F. Y. Gang, B. D. Scott, and P. H. Diamond, Phys. Fluids B 1, 1331 (1989).
- W. Dorland and G. W. Hammett, Phys. Fluids B 5, 812 (1993).
- A. Hasegawa and K. Mima, Phys. Rev. Lett. 39, 205 (1977).
- A. Arakawa, J. Comput. Phys. 1, 119 (1966).
- G. E. Karniadakis, M. Israeli, and S. A. Orszag, J. Comput. Phys. 97, 414 (1991).
- R. Numata, R. Ball, and R. L. Dewar, in Frontiers in Turbulence and Coherent Structures: Proceedings of the CSIRO/COSNet Workshop on Turbulence and Coherent Structures, Canberra, Australia, 10-13 January 2006, edited by J. P. Denier and J. S. Frederiksen (World Scientific, Singapore, 2007), vol. 6 of World Scientific Lecture Notes in Complex Systems, pp. 431–442.
- Lord Rayleigh, Proc. London Math. Soc. 11, 57 (1879).
- W. Tollmien, Nachr. Ges. Wiss. Göttingen, Math.-Phys. Kl. 50, 79 (1935).
- C. C. Lin, Quart. Appl. Math. 3, 218 (1945).
- P. G. Drazin and W. H. Reid, Hydrodynamic Stability (Cambridge University Press, Cambridge, 1981).
- H. L. Kuo, J. Met. 6, 105 (1949).