Zonal Flow as Pattern Formation
ZZonal Flow as Pattern Formation
Jeffrey B. Parker a) and John A. Krommes b) Princeton Plasma Physics Laboratory, Princeton University, Princeton, New Jersey 08543,USA (Dated: 31 October 2018)
Zonal flows are well known to arise spontaneously out of turbulence. We show that for statistically averagedequations of the stochastically forced generalized Hasegawa-Mima model, steady-state zonal flows and inho-mogeneous turbulence fit into the framework of pattern formation. There are many implications. First, thewavelength of the zonal flows is not unique. Indeed, in an idealized, infinite system, any wavelength within acertain continuous band corresponds to a solution. Second, of these wavelengths, only those within a smallersubband are linearly stable. Unstable wavelengths must evolve to reach a stable wavelength; this processmanifests as merging jets.Zonal flows (ZFs) — azimuthally symmetric, gener-ally banded, shear flows — are spontaneously generatedfrom turbulence and have been reported in atmospheric and laboratory plasma contexts. Recently, they havealso been observed in astrophysical simulations. In mag-netically confined plasmas, ZFs are thought to play acrucial role in regulation of turbulence and turbulenttransport.
A greater understanding of ZF behavior isvaluable for untangling a host of nonlinear processes inplasmas, including details of transitions between modesof low and high confinement.Zonal flows remain incompletely understood, even re-garding the basic question of the jet width (wavelength).In the plasma literature, one finds modulational or sec-ondary instability calculations of ZF generation, butthese cannot provide information on a saturated state.Other theories typically make an assumption of long-wavelength ZFs and leave the ZF scale as an undeter-mined parameter. Within geophysical contexts, variousauthors have attempted to relate the jet width or spacingto length scales that emerge from the vorticity equationby heuristically balancing the magnitudes of the Rossbywave term and the nonlinear advection. Those scales in-clude the Rhines scale and other, similar scales.
ARhines-like length scale is also obtained from argumentsbased on potential vorticity staircases.
However, nei-ther the heuristic Rhines estimates nor the paradigmof potential vorticity inversion and mixing generalize tomore complex situations involving realistic plasma mod-els. We are therefore motivated to seek a more system-atic approach to determining the ZF width that may offersuch a generalization.A related topic is the merging of jets. Coales-cence of two or more jets is ubiquitous in numericalsimulations.
The merging process occurs during theinitial transient period before a statistically steady stateis reached. It is clear that the merging is part of a dynam-ical process through which the ZF reaches its preferredlength scale, but the merging phenomeon has not been a) Electronic mail: [email protected] b) Electronic mail: [email protected] understood thus far.Our present work addresses these questions in the con-text of the stochastically forced generalized Hasegawa-Mima (GHM) equation, a model of magnetizedplasma turbulence in the presence of a background den-sity gradient. This model is mathematically similar tothe barotropic vorticity equation on a β plane. Ouranalysis is related to several recent works that focusedon that equation in the geophysical context.
Impor-tantly, numerical simulations of both models can displayemergence of steady ZFs. The GHM equation and the pa-rameterizations of forcing and dissipation that we use arenot realistic descriptions of plasma; however, the simplic-ity is an asset in understanding the qualitative behaviorof these systems.We study a statistical average of the flow. Statisticalapproaches enable one to gain physical insight by aver-aging away the details of the turbulent fluctuations andworking with smoothly varying quantities. Sometimesstatistical turbulence theories strive for quantitative ac-curacy, which requires rather complicated methods. Incontrast, our investigation is at a more basic level andconcerns the fundamental nature of ZFs interacting self-consistently with inhomogeneous turbulence.Within the statistical framework, we build uponrecent understanding of zonostrophic instability, inwhich homogeneous turbulence becomes unstable to ZFperturbations. Steady ZFs emerge from this bifurca-tion. We show that the bifurcation obeys a classic am-plitude equation, and therefore ZFs can be understoodas pattern formation.
Two important results followfrom the general properties of pattern-forming systems.First, the wavelength of the ZF is not unique. Indeed,in an idealized, infinite system, any wavelength within acertain continuous band corresponds to a steady-state so-lution. Second, of these wavelengths, only those withina smaller subband are linearly stable. Unstable wave-lengths must evolve to reach a stable wavelength. For un-stable jets of short (long) wavelength, this process man-ifests as merging (branching) jets.Our basic model is the 2D GHM equation, ∂ t w ( x, y ) + v · ∇ w − κ∂ y φ = ξ − µw − ν ( − h ∇ h w, (1)where φ = ( L n /ρ s ) eϕ/T e is the normalized electrostatic a r X i v : . [ phy s i c s . a o - ph ] O c t potential, L n is the density gradient scale length, ρ s is thesound radius, T e is the electron temperature, w = ∇ φ − ˆ αφ is the generalized vorticity and is related to ion gyro-center density fluctuations δn Gi by w = − ( L n /ρ s ) δn Gi /n where n is the background density, ˆ α is an operator suchthat in Fourier space ˆ α ( k ) = 0 if k y = 0 (ZF mode) andˆ α ( k ) = 1 if k y (cid:54) = 0 (drift wave mode), v = ˆz × ∇ φ is the E × B velocity, µ is a constant frictional drag, ν is the viscosity with hyperviscosity factor h , ξ is white-noise forcing, and κ is related to the density scale length.Lengths are normalized to ρ s and times are normalizedto the drift wave period ω − ∗ = ( L n /ρ s )Ω − i . These nor-malizations and scalings are convenient to make w , φ ,and the active length and time scales of order unity, andthey allow us to set κ = 1.The ZF behavior in numerical simulations of Eq. (1) isshown in Fig. 1(a). During the transient period, mergingjets are observed, while in the late time a statisticallysteady state is reached with stable unwavering jets.We restrict ourselves to the quasilinear (QL) approxi-mation of this system. To obtain the QL equations, weperform an eddy–mean decomposition, given by decom-posing all fields into a zonal mean and a deviation fromthe zonal mean, then neglect the eddy–eddy nonlineari-ties within the eddy equation. The QL approximationis not expected to be physically and quantitatively cor-rect in detail (though it may be in certain regimes );for example, material conservation of potential vorticity(in the undamped, undriven case) is lost. However, theQL model is useful because it exhibits the same basiczonal jet features as the full model, namely merging jetsand the formation of stable jets. Therefore, analysis ofthe QL model can provide a mathematical foundation forunderstanding and interpreting the physical behavior.We consider a statistical average of the QL system.In the presence of steady ZFs, a statistical homogeneityassumption is clearly invalid. Therefore, we allow theturbulence to be inhomogeneous in the direction ( x ) ofZF variation. The averaged equations, referred to as thesecond-order cumulant expansion (CE2), are ∂ t W + ( U + − U − ) ∂ y W − ( U (cid:48)(cid:48) + − U (cid:48)(cid:48)− ) (cid:18) ∇ + 14 ∂ x (cid:19) ∂ y C + [2 κ + ( U (cid:48)(cid:48) + + U (cid:48)(cid:48)− )] ∂ x ∂ x ∂ y C = F − µW − νD h W, (2a) ∂ t U + ∂ x ∂ x ∂ y C (0 , , x, t ) = − µU − ν ( − h ∂ hx U, (2b)where x and y represent two-point separations, x repre-sents the two-point average position (if the turbulencewere homogeneous, there would be no x dependence), W ( x, y | x, t ) and C ( x, y | x, t ) are the one-time, two-space-point correlation functions of vorticity and poten-tial, U ( x, t ) is the zonal flow velocity, U ± = U ( x ± x/ , t ), ∇ = ∂ x + ∂ y − F ( x, y ) is chosen to be isotropic, homo-geneous ring forcing, and D h is a hyperviscosity operator.There is a linear relation between W and C . Given the assumption that the stochastic forcing ξ iswhite (delta-correlated) noise, the only further assump- ω ∗ t x / ρ s t x (a)(b) FIG. 1. (a) Merging jets during the transient regime of equa-tion (1) (zonal-mean velocity is shown). (b) Merging behaviorin the amplitude equation (3) [Re A ( x, t ) is shown]. tions necessary for CE2 to be an exact description ofthe QL model are statistical homogeneity and ergodic-ity in the zonal ( y ) direction. This is because the QLmodel neglects the nonlinear eddy–eddy term that wouldgive rise to a closure problem. Alternatively, CE2 canbe regarded as a truncated statistical closure of the fullmodel. The CE2 equations exhibit important symmetries oftranslation and reflection, given by x → x + δx , ( x, x ) → ( − x, − x ), ( y, x ) → ( − y, − x ), and ( x, y ) → ( − x, − y ).Many studies of CE2 have been performedpreviously. Numerical simulations of CE2also exhibit merging jets. For Eq. (2) there always exists a homogeneous equilib-rium: W ( x, y ) = (2 µ + 2 νD h ) − F , U = 0. This equi-librium is stable in a certain regime of parameters. As acontrol parameter such as µ is varied, this homogeneousstate becomes zonostrophically unstable. Physically,zonostrophic instability occurs when dissipation is over-come by the mutually reinforcing processes of eddy tiltingby zonal flows and production of Reynolds stress forcesby tilted eddies. The eigenmode consists of perturbationsspatially periodic in x with zero real frequency, so thatzonostrophic instability arises as a Type I s instability ofhomogeneous turbulence. Zonostrophic instability withinCE2 may be thought of as a variant of modulational in-stability calculations of ZF generation.Just beyond the instability threshold, a bifurcationanalysis follows a standard procedure and involves a mul-tiscale perturbation expansion about the threshold. Let u be the state vector relative to the homogeneous equi-librium and let (cid:15) be a normalized control parameter. Theexpansion proceeds as u = (cid:15) / u + (cid:15)u + · · · . At firstorder one finds u = A ( x, t ) r + c.c., where c.c. denotescomplex conjugate and r ∼ e iq c x is the eigenmode thatis marginally stable at (cid:15) = 0. One determines a PDEfor the complex amplitude A as a solvability conditionat third order in the perturbation expansion. This am-plitude equation is constrained by the translation andreflection symmetries to take a universal form. Theamplitude equation, sometimes referred to as the realGinzburg-Landau equation, is ∂ t A ( x, t ) = A + ∂ x A − | A | A, (3)where all coefficients have been rescaled to unity. Thederivation of Eq. (3) from Eq. (2) will be reportedelsewhere. The amplitude equation (3) is well understood.
First, a steady-state solution exists for any wave num-ber within the continuous band − < k < A = αe ikx with | α | = 1 − k is asolution). Second, only solutions with k < / This is demonstrated in Fig. 1(b), wherean unstable solution that has been slightly perturbed un-dergoes merging behavior until a stable wave number isreached. The preceding qualitative behaviors are alsoexhibited by the CE2 system, as we now show.We proceed to find the steady-state solutions ofEq. (2). In the context of an infinite domain with noboundaries, these solutions are referred to as ideal states.Let q denote the basic ZF wave number of an ideal state.For a given q , we solve the time-independent form ofEq. (2) directly. This approach is distinct from time in-tegration of Eq. (2) to a steady state. Our procedure hastwo advantages for understanding the global structure ofthe system. First, we can specify precisely the q of thedesired solution. Second, we can solve directly for allsolutions, including unstable ones, rather than find onlythose which develop from time evolution.An ideal state is represented as a Fourier-Galerkin se-ries with coefficients to be determined : W ( x, y | x ) = M (cid:88) m = − M N (cid:88) n = − N P (cid:88) p = − P W mnp e imax e inby e ipqx , (4a) U ( x ) = P (cid:88) p = − P U p e ipqx . (4b)While the periodicity in x is desired, the correlation func-tion should decay in x and y ; periodicity in x and y arsiesfrom using the convenient Fourier basis. Thus, a and b ,unlike q , are numerical parameters. They represent thespectral resolution of the correlation function and shouldbe small enough to obtain an accurate solution.The CE2 symmetries allow us to seek a solution where U ( x ) = U ( − x ) and W ( x, y | x ) = W ( − x, − y | x ) = W ( x, − y | − x ) = W ( − x, y | − x ). These constraints,along with reality conditions, force U p to be real, U p = U − p , and W mnp = W ∗− m,n,p = W ∗ m, − n,p = W ∗ m,n, − p .We obtain a system of nonlinear algebraic equations forthe coefficients U p , W mnp by substituting the Galerkinseries into Eq. (2) and projecting onto the basis func-tions. To demonstrate the projection for Eq. (2a), let φ mnp = e imax e inby e ipqx . We project Eq. (2a) onto φ rst qρ s qρ s U U (b) U U E EN N E EN ND
FIG. 2. Zonal flow amplitude U , U as a function of idealstate wave number q at (a) µ = 0 .
21 ( R β = 1 .
48) and (b) µ = 0 .
19 ( R β = 1 . by operating with (cid:18) πa πb πq (cid:19) − (cid:90) π/a − π/a dx (cid:90) π/b − π/b dy (cid:90) π/q − π/q dx φ ∗ rst . (5)For instance, the term ( U + − U − ) ∂ y W projects to I rstp (cid:48) mnp U p (cid:48) W mnp , where repeated indices are summedover, I rstp (cid:48) mnp = inbδ n,s δ p (cid:48) + p − t, ( σ + − σ − ), σ ± =sinc( α ± π/a ), and α ± = ma − ra ± p (cid:48) q/
2. The other termsof Eq. (2a), as well as Eq. (2b), are handled similarly.The system of nonlinear algebraic equations is solvedwith a Newton’s method. Figure 2 shows the ZF am-plitude coefficients U p as functions of q at µ = 0 .
21 and µ = 0 .
19. Near the instability threshold, ideal statesexist at all q for which the homogeneous equilibrium iszonostrophically unstable [between the two lines labeled N in Fig. 2(a)]. Farther from threshold, there is a regionof q where the ideal-state solution seems to disappear[between the lines N and D in Fig. 2(b); see also Fig.3]. The values of the other parameters used are κ = 1, ν = 10 − , and h = 4. The forcing F ( k ) = 2 πεk f /δk for k f − δk < | k | < k f + δk , and is zero otherwise. We take k f = 1, δk = 1 /
8, and ε , which acts like a total energyinput rate, to be equal to 1.To investigate stability of the ideal states, we considerperturbations δW ( x, y | x, t ) and δU ( x, t ) about an equi-librium W, U and linearize Eq. (2). Since the underlyingequilibrium is periodic in x , the perturbations can beexpanded as a Bloch state : δW ( x, y | x, t ) = e σt e iQx (cid:88) mnp δW mnp e imax e inby e ipqx , (6a) δU ( x, t ) = e σt e iQx (cid:88) p δU p e ipqx , (6b)where Q is the Bloch wave number and can be taken tolie within the first Brillouin zone − q/ < Q ≤ q/
2. Wedo not use a Q x or Q y because as previously mentioned qρ s D i m e n s i o n l e ss P a r a m e t e r γ NSD
FIG. 3. Stability diagram for the CE2 equations. Abovethe neutral curve (N), the homogeneous turbulent state iszonostrophically unstable. Ideal states are stable within themarginal stability curve S . The stability curve is consistentwith the dominant ZF wavenumber from independent QL sim-ulations (crosses). The stationary ideal states vanish to theleft of D . Here, a = 0 . b = 0 . M = 20, N = 33, P = 5,and other parameters are given in the text. γ is varied bychanging µ while holding other parameters fixed. the periodicity in x and y is artificial. The perturbationequations are projected onto the basis functions in thesame way as in the ideal state calculation. This projec-tion results in a linear system at each Q for the coef-ficients δW mnp and δU p ; this determines an eigenvalueproblem for σ . The equilibrium is unstable if for any Q there are any eigenvalues with Re σ > γ = ε / κ / µ − / , an important di-mensionless parameter controlling the ZF dynamics. To vary γ , we change µ and hold other parameters fixedat their previous values. The stable ideal states existinside of the marginal stability curve marked S . Nearthe threshold, marginal stability is governed by the Eck-haus instability, a long-wavelength universal instability. Farther from threshold, the instability transitions intonew, nonuniversal instabilities; details will be reportedelsewhere. The ZFs are spontaneously generated for γ > .
53. For γ > .
53, the stability curve is consis-tent with the dominant ZF wavenumber observed in QLsimulations.Numerical simulations typically are done within a fi-nite domain. When periodic boundary conditions areused, our infinite-domain results are modified merely bythe discretization of wave numbers. This affects not onlythe possible equilibria, but also any perturbations andhence the stability boundaries too.For a time-evolving system, the exact q that is ulti-mately chosen within the stability balloon results from adynamical process and is not addressed in a systematicway by the present study.While the CE2 equations exhibit spontaneously gener-ated zonal flows, they neglect many physical effects. Animportant piece of physics missing from the CE2 equa-tions is the nonlinear eddy self-interaction, which clearlycannot be ignored in general. At least one particular in- stance of the qualitative failure of CE2 has been noted. Yet, the basic mathematical structure of the theorypresented here arises only from symmetry arguments andgeneral properties of the zonostrophic instability. If onewere to include the important physics neglected in CE2,those general symmetries and properties should remainintact. Therefore, we expect our qualitative conclusionsto likewise remain valid.In summary, by analyzing a second-order statisticalmodel of an ensemble of interacting zonal flows and tur-bulence, we have shown that zonal flows constitute pat-tern formation amid a turbulent bath. This continuesprevious work to provide a firm analytic understand-ing of zonal flow generation and equilibrium within CE2.We calculated the stability diagram of steady zonal jetsand explained the merging of jets as a means of attaininga stable wave number. In general, the use of statisticallyaveraged equations and the pattern formation methodol-ogy provide a path forward for further systematic inves-tigations of zonal flows and their interactions with turbu-lence. Further work should be done to understand howthis framework can shed light on practical problems in-volving realistic plasmas in the near-collisionless regime.We acknowledge useful discussions with Brian Farrell,Henry Greenside, Petros Ioannou, and Brad Marston.This material is based upon work supported by an NSFGraduate Research Fellowship and a US DOE Fusion En-ergy Sciences Fellowship. This work was also supportedby US DOE Contract DE-AC02-09CH11466. A. R. Vasavada and A. P. Showman, Rep. Prog. Phys. , 1935(2005). A. Fujisawa, Nucl. Fusion , 013001 (2009). A. Johansen, A. Youdin, and H. Klahr, Astrophys. J. , 1269(2009). Z. Lin, T. S. Hahm, W. W. Lee, W. M. Tang, and R. B. White,Science , 1835 (1998). P. H. Diamond, S.-I. Itoh, K. Itoh, and T. S. Hahm, PlasmaPhysics and Controlled Fusion , R35 (2005). B. N. Rogers, W. Dorland, and M. Kotschenreuther, Phys. Rev.Lett. , 5336 (2000). C. Connaughton, S. Nazarenko, and B. Quinn, EPL , 25001(2011). P. B. Rhines, J. Fluid Mech. , 417 (1975). G. K. Vallis and M. E. Maltrud, J. Phys. Oceanogr. , 1346(1993). S. Sukoriansky, N. Dikovskaya, and B. Galperin, J. Atmos. Sci. , 3312 (2007). D. G. Dritschel and M. E. McIntyre, J. Atmos. Sci. , 855(2008). R. K. Scott and D. G. Dritschel, Journal of Fluid Mechanics ,576 (2012). H.-P. Huang and W. A. Robinson, J. Atmos. Sci. , 611 (1998). R. K. Scott and L. M. Polvani, J. Atmos. Sci. , 3158 (2007). A. I. Smolyakov, P. H. Diamond, and M. Malkov, Phys. Rev.Lett. , 491 (2000). J. A. Krommes and C.-B. Kim, Phys. Rev. E , 8508 (2000). K. Srinivasan and W. R. Young, J. Atmos. Sci. , 1633 (2012). B. F. Farrell and P. J. Ioannou, J. Atmos. Sci. , 2101 (2003). B. F. Farrell and P. J. Ioannou, J. Atmos. Sci. , 3652 (2007). N. A. Bakas and P. J. Ioannou, Phys. Rev. Lett. , 224501(2013). N. C. Constantinou, P. J. Ioannou, and B. F. Farrell,
Emergenceand equilibration of jets in beta-plane turbulence: applications of Stochastic Structural Stability Theory , submitted to J. At-mos. Sci. J. B. Marston, E. Conover, and T. Schneider, J. Atmos. Sci. ,1955 (2008). S. M. Tobias, K. Dagon, and J. B. Marston, Astrophys. J. ,127 (2011). S. M. Tobias and J. B. Marston, Phys. Rev. Lett. , 104502(2013). J. A. Krommes, Phys. Rep. , 1 (2002). M. C. Cross and P. C. Hohenberg, Rev. Mod. Phys. , 851(1993). M. Cross and H. Greenside,
Pattern Formation and Dynamicsin Nonequilibrium Systems (Cambridge University Press, Cam-bridge, 2009). R. Hoyle,
Pattern Formation: An Introduction to Methods (Cambridge University Press, Cambridge, 2006). R. M. Clever and F. H. Busse, J. Fluid Mech. , 625 (1974). A. C. Newell, T. Passot, and M. Souli, J. Fluid Mech. , 187(1990). F. Bouchet, C. Nardini, and T. Tangarife, J. Stat. Phys. (2013). B. F. Farrell and P. J. Ioannou, J. Atmos. Sci. , 2444 (2009). J. B. Parker and J. A. Krommes,
Generation of zonal flowsthrough symmetry breaking of statistical homogeneity , submit-ted to New. J. Phys. C. T. Kelley,
Solving Nonlinear Equations with Newton’s Method (Society for Industrial and Applied Mathematics, Philadelphia,2003). B. Galperin, S. Sukoriansky, and N. Dikovskaya, Ocean Dynam-ics , 427 (2010). S. Danilov and D. Gurarie, Physics of Fluids16