Wave propagation in the Lorenz-96 model

In this paper we study the spatiotemporal properties of waves in the Lorenz-96 model and their dependence on the dimension parameter n and the forcing parameter F . For F > 0 the first bifurcation is either a supercritical Hopf or a double-Hopf bifurcation and the periodic attractor born at these bifurcations represents a traveling wave. Its spatial wave number increases linearly with n, but its period tends to a finite limit as n→∞. For F < 0 and odd n, the first bifurcation is again a supercritical Hopf bifurcation, but in this case the period of the traveling wave also grows linearly with n. For F < 0 and even n, however, a Hopf bifurcation is preceded by either one or two pitchfork bifurcations, where the number of the latter bifurcations depends on whether n has remainder 2 or 0 upon division by 4. This bifurcation sequence leads to stationary waves and their spatiotemporal properties also depend on the remainder after dividing n by 4. Finally, we explain how the double-Hopf bifurcation can generate two or more stable waves with different spatiotemporal properties that coexist for the same parameter values n and F .


Introduction
In this paper we study the Lorenz-96 model which is defined by the equations dx j dt = x j −1 (x j +1 − x j −2 ) − x j + F, j = 0, . . ., n − 1, (1) together with the periodic "boundary condition" implied by taking the indices j modulo n. The dimension n ∈ N and the forcing parameter F ∈ R are free parameters. Lorenz (2006) interpreted the variables x j as values of some atmospheric quantity in n equispaced sectors of a latitude circle, where the index j plays the role of "longitude". Hence, a larger value of n can be interpreted as a finer latitude grid. Lorenz also remarked that the vectors (x 0 , . . ., x n−1 ) can be interpreted as wave profiles, and he observed that for F > 0 sufficiently large these waves slowly propagate "westward", i.e., in the direction of decreasing j . Figure 1 shows a Hovmöller diagram illustrating two traveling waves with wave number 5 for dimension n = 24 and the parameter values F = 2.75 (in the periodic regime) and F = 3.85 (in the chaotic regime). The Lorenz-96 model was introduced as a tool for numerical experiments in predictability studies, rather than as a physically realistic model. Indeed, Lorenz (2006) wrote that "the physics of the atmosphere is present only to the extent that there are external forcing and internal dissipation, simulated by the constant and linear terms, while the quadratic terms, simulating advection, together conserve the total energy [. . . ]". The value of the Lorenz-96 model primarily lies in the fact that it has a very simple implementation in numerical codes while at the same time it can exhibit very complex dynamics for suitable choices of the parameters n and F . The famous Lorenz-63 model (Lorenz, 1963), which does have a clear physical interpretation (namely, as the Galerkin projection of a fluid dynamical model describing Rayleigh-Bénard convection), has two disadvantages. Firstly, it consists of only three ordinary differential equations. Secondly, for the classical parameter values the model has a Lyapunov spectrum (0.91, 0, −14.57), which makes the model very dissipative. Such properties are not typical for atmospheric models. In contrast, the dimension of the Lorenz-96 model can be chosen to be arbitrarily large, and for suitable values of the parameters the Lyapunov spectrum is similar to those observed in models obtained from discretizing partial differential equations. For those reasons the Lorenz- Table 1. Recent papers with applications of the Lorenz-96 model and the values of n that were used.
96 model has become a test model for a wide range of geophysical applications. Table 1 lists some recent papers with applications of the Lorenz-96 model. In most studies the dimension n is chosen ad hoc, but n = 36 and n = 40 appear to be popular choices. Many applications are related to geophysical problems, but the model has also attracted the attention of mathematicians working in the area of dynamical systems for phenomenological studies in high-dimensional chaos. Note that Eq. (1) is in fact a family of models parameterized by means of the discrete parameter n. An important question is to what extent both the qualitative and quantitative dynamical properties of Eq. (1) depend on n. For example, the dimension n has a strong effect on the predictability of large amplitudes of traveling waves in weakly chaotic regimes of the Lorenz-96 model . In general, the statistics of extreme events in dynamical systems strongly depend on topological properties and recurrence properties of the system (Holland et al., , 2016. Therefore, a coherent overview of the dependence of spatiotemporal properties on the parameters n and F is useful to assess the robustness of results when using the Lorenz-96 model in predictability studies. In this paper we address the question of how the spatiotemporal properties of waves, such as their period and wave number, in the Lorenz-96 model depend on the dimension n and whether these properties tend to a finite limit as n → ∞. We will approach this question by studying waves represented by periodic attractors that arise through a Hopf bifurcation of a stable equilibrium. Along various routes to chaos these periodic attractors can bifurcate into chaotic attractors representing irregular waves which "inherit" their spatiotemporal properties from the periodic attractor. For example, the wave shown in the left panel of Fig. 1 bifurcates into a three-torus attractor which breaks down and gives rise to the wave in the right panel. Note that both waves have the same wave number. Figure 2 shows power spectra of these waves, and clearly their dominant peaks are located at roughly the same period. Inheritance of spatiotemporal properties also manifests itself in a shallow water model studied by Sterk et al. (2010) in which a Hopf bifurcation (related to baroclinic instability) explains the observed timescales of atmospheric low-frequency variability. The Hopf bifurcation plays a key role in explaining the physics of low-frequency variability in many geophysical contexts. Examples are the Atlantic Multidecadal Oscillation (Te Raa and Dijkstra, 2002;Dijkstra et al., 2008;Frankcombe et al., 2009), the wind-driven ocean circulation (Simonnet et al., 2003a, b), and laboratory experiments (Read et al., 1992;Tian et al., 2001). In addition to traveling waves, such as illustrated in Fig. 1, we will also show the existence of stationary waves. In a recent paper by Frank et al. (2014) stationary waves have also been discovered in specific regions of the multi-scale Lorenz-96 model. Their paper uses dynamical indicators such as the Lyapunov dimension to identify the parameter regimes with stationary waves. Moreover, we will explain two bifurcation scenarios by which waves with different spatiotemporal properties coexist. This paper complements the results of our previous work (Van Kekem and Sterk, 2018), which considers only the case F > 0; these two papers together give a comprehensive picture of wave propagation in the Lorenz-96 model.
The remainder of this paper is organized as follows. In Sect. 2 we explain how to obtain an approximation of the periodic attractor born at a Hopf bifurcation, which enables us to derive spatiotemporal properties of waves in the Lorenz-96 model. In Sect. 3.1 we show that, for F > 0, periodic attractors indeed represent traveling waves as suggested by Lorenz. Also for F < 0 and odd values of n, periodic attractors represent traveling waves, as is demonstrated in Sect. 3.2. In Sect. 3.3, however, we show analytically that, for n = 6 and F < 0, stationary waves occur. By means of numerical experiments we show in Sect. 3.4 that stationary waves occur in general for even n and F < 0. In Sect. 4 we discuss the bifurcation scenarios by which stable waves with different spatiotemporal properties can coexist for the same values of the parameters n and F .

Hopf bifurcations
In this section we consider a general geophysical model in the form of a system of ordinary differential equations: In this equation, µ ∈ R is a parameter modeling external circumstances such as forcing. Assume that for the parameter value µ 0 the system has an equilibrium x 0 ; this means that f (x 0 , µ 0 ) = 0 and hence x 0 is a time-independent solution of Eq. (2). In the context of geophysics x 0 represents a steady flow, and its linear stability is determined by the eigenvalues of the Jacobian matrix Df (x 0 , µ 0 ). An equilibrium becomes unstable when eigenvalues of the Jacobian matrix cross the imaginary axis upon variation of the parameter µ. Dijkstra (2005) provides an extensive discussion of the physical interpretation of bifurcation behavior. Assume that Df (x 0 , µ 0 ) has two eigenvalues ±ωi on the imaginary axis. This indicates the occurrence of a Hopf bifurcation, i.e., the birth of a periodic solution from an equilibrium that changes stability. If the equilibrium x 0 is stable for µ < µ 0 and unstable for µ > µ 0 , then under suitable nondegeneracy conditions the Hopf bifurcation is supercritical, which means that a stable periodic orbit exists for µ > µ 0 (Kuznetsov, 2004). For small values of ε = √ µ − µ 0 the periodic orbit that is born at the Hopf bifurcation can be approximated by see Beyn et al. (2002). Without loss of generality we may assume that the corresponding complex eigenvectors u ± iv of the matrix Df (x 0 , µ 0 ) have unit length. In the context of geophysical applications this first-order approximation of the periodic orbit can be interpreted as a wave-like perturbation imposed on a steady mean flow. The spatiotemporal properties of this wave can now be determined by the vectors x 0 , u, v, and the frequency ω.

Waves in the Lorenz-96 model
In this section we study waves in the Lorenz-96 model and how their spatiotemporal characteristics depend on the parameters n and F .

Traveling waves for n ≥ 4 and F > 0
For all n ∈ N and F ∈ R the point x F = (F, . . ., F ) is an equilibrium solution of Eq. (1). This equilibrium represents a steady flow, and since all components are equal the flow is spatially uniform. The stability of x F is determined by the eigenvalues of the Jacobian matrix of Eq. (1). Note that the Lorenz-96 model is invariant under the symmetry x i → x i+1 while taking into account the periodic boundary condition. As a consequence the Jacobian matrix evaluated at x F is circulant, which means that each row is a right cyclic shift of the previous row, and so the matrix is completely determined by its first row. If we denote this row by (c 0 , c 1 , . . ., c n−1 ), then it follows from Gray (2006) that the eigenvalues of the circulant matrix can be expressed in terms of roots of unity ρ j = exp(−2π ij/n) as follows: An eigenvector corresponding to λ j is given by In particular, for the Lorenz-96 model in Eq. (1) the Jacobian matrix at x F has only three nonzero elements on its first row, viz. c 0 = −1, c 1 = F , c n−2 = −F . Hence, the eigenvalues λ j can be expressed in terms of n and F as follows: where the functions f and g are defined as For F = 0 the equilibrium x F is stable as Reλ j = −1 for all j = 0, . . ., n − 1. The real part of the eigenvalue λ j changes sign if the equation is satisfied. The graph of f in Fig. 3 shows that for F > 0 Eq. (6) can have at most four solutions. Since f is symmetric around x = π it follows that if j is a solution of Eq. (6) then so is n − j . This means that the equilibrium x F becomes unstable for F > 0 when either a pair or a double pair of eigenvalues becomes purely imaginary. The main result is summarized in the following theorem.
Theorem 1. Assume that n ≥ 4 and l ∈ N satisfies 0 < l < n 2 , l = n 3 . Then the lth eigenvalue pair (λ l , λ n−l ) of the trivial equilibrium x F crosses the imaginary axis at the parameter value F H (l, n) := 1/f (2π l/n) and thus x F bifurcates through either a Hopf or a double-Hopf bifurcation. A double-Hopf bifurcation, with two pairs of eigenvalues crossing the imaginary axis, occurs if and only if there exist l 1 , l 2 ∈ N such that cos 2π l 1 n + cos Otherwise, a Hopf bifurcation occurs. Moreover, the first Hopf bifurcation of x F is always supercritical.
The proof of Theorem 1 can be found in Van Kekem and Sterk (2018), in which also an expression for the first Lyapunov coefficient is derived which determines for which l the Hopf bifurcation is sub-or supercritical. Observe that Theorem 1 implies that a double-Hopf bifurcation occurs for n = 10m (with l 1 = m, l 2 = 3m, F = 2) and n = 12m (with l 1 = 2m, l 2 = 3m, F = 1). In Sect. 4 we will explain how double-Hopf bifurcations lead to the coexistence of two or more stable traveling waves with different wave numbers.
From the eigenvalues that cross the imaginary axis and the corresponding eigenvectors we can deduce the physical characteristics of the periodic orbit that arises after a Hopf bifurcation. When the lth eigenvalue pair (λ l , λ n−l ) crosses the imaginary axis we can write If we set ω = cot(π l/n), then according to Eq. (3) an approximation of the periodic orbit is given by This is indeed the expression for a traveling wave in which the spatial wave number and the period are given by respectively l and T = 2π/ω = 2π tan(π l/n). Thus the index of the eigenpair that crosses the imaginary axis determines the propagation characteristics of the wave. Note that Hopf bifurcations of an unstable equilibrium will result in an unstable periodic orbit. Therefore, not all waves that are guaranteed to exist by Theorem 1 will be visible in numerical experiments. Equation (6) implies that for F > 0 the first Hopf bifurcation occurs for the eigenpair (λ l , λ n−l ) with index In Appendix A1 it is shown that, except for n = 7, the integer l + 1 (n) satisfies the bounds which means that the wave number increases linearly with the dimension n. Since the function f has a maximum at which is consistent with Eq. (9). As a corollary we find that the period of this wave tends to a finite limit as n → ∞: Figure 4 shows a graph of the period and the wave number as a function of n. Note that the period settles down on the value T ∞ .

Traveling waves for odd n ≥ 4 and F < 0
Now assume that n is odd. For F < 0 Eq. (6) has precisely two solutions which implies that the first bifurcation of x F is a supercritical Hopf bifurcation. The index of the first bifurcating eigenpair (λ l , λ n−l ) follows by minimizing the value of the function f in Eq. (5): Again, the wave number increases linearly with n, but at a faster rate than in the case F > 0. Now the period of the wave is given by where the last equality follows from the computations in Appendix A2. This implies that contrary to the case F > 0 the period increases monotonically with n and does not tend to a limiting value as n → ∞.
Note that for even n and F < 0 the first bifurcation is not a Hopf bifurcation since λ n/2 = −1 − 2F is a real eigenvalue that changes sign at F = 1 2 . Surprisingly, the case n = 4 is not analytically tractable. The case n = 6 will be studied analytically in Sect. 3.3. In Sect. 3.4 we will numerically study the bifurcations for other values of n and F < 0.
3.3 Stationary waves for n = 6 and F < 0 We now consider the dimension n = 6. At F = − 1 2 the eigenvalue λ 3 changes sign. Note that the equilibrium x F cannot exhibit a saddle-node bifurcation since x F continues to exist for F < − 1 2 . Instead, at F = − 1 2 there must be a branching point which is either a pitchfork or a transcritical bifurcation. If we try for F < − 1 2 an equilibrium solution of the form x P = (a, b, a, b, a, b) then it follows that a and b must satisfy the equations Of course a = b = F is a solution to these equations, but this would lead to the already known equilibrium x F = (F, F, F, F, F, F ). There is an additional pair of solutions which is given by With these values of a and b we obtain two new equilibria x P ,1 = (a, b, a, b, a, b) and x P ,2 = (b, a, b, a, b, a) that exist for F < − 1 2 in addition to the equilibrium x F . This means that a pitchfork bifurcation occurs at F = − 1 2 .  Figure 4. As the equilibrium x F = (F, . . ., F ) loses stability through a (double-)Hopf bifurcation for F > 0 a periodic attractor is born which represents a traveling wave. The spatial wave number increases linearly with n, whereas the period tends to a finite limit.  Figure 5. As Fig. 1, but for two periodic attractors for n = 6 and F = −3.6. These attractors are born at Hopf bifurcations of the equilibria x P ,1 (a) and x P ,2 (b) at F = −7/2. Note that the waves do not travel "eastward" or "westward". The pitchfork bifurcation changed the mean flow which in turn changes the propagation of the wave.
As F decreases, each of the new equilibria x P ,1,2 may bifurcate again. We first consider the equilibrium x P ,1 for which the Jacobian matrix is given by Note that J is no longer circulant: in addition to shifting each row in a cyclic manner, the values of a and b also need to be interchanged. In particular, this means that the eigenvalues can no longer be determined by means of Eq. (4). Symbolic manipulations with the computer algebra package Mathematica (Wolfram Research, Inc., 2016) show that an eigenvalue crossing occurs for F = − 7 2 , in which case a = 1 2 (−1 + √ 6) and b = 1 2 (−1 − √ 6) so that the characteristic polynomial of J is given by This expression shows that J has two purely imaginary eigenvalues ±i √ 3 and the remaining four complex eigenvalues have a negative real part. Therefore the equilibrium x P ,1 undergoes a Hopf bifurcation at F = − 7 2 . Numerical experiments with Mathematica show that the matrix J − i √ 3I has a null vector of the form Hence, using Eq. (3) the periodic orbit can be approximated as Note that if ε = − 7 2 − F is sufficiently small, then x j (t) is always positive (resp. negative) for j = 0, 2, 4 (resp. j = 1, 3, 5). This implies that the periodic orbit represents a stationary wave rather than a traveling wave. The period of the wave is T = 2π/ √ 3 and the spatial wave number is 3. These spatiotemporal properties are clearly visible in the left panel of Fig. 5.
The computations for the equilibrium x P ,2 are similar and show that another Hopf bifurcation takes place at F = − 7 2 . This means that for F < − 7 2 there exists a second stable periodic orbit which coexists with the stable periodic orbit born at the Hopf bifurcation of x P ,1 . Its first-order approximation is almost identical to Eq. (11): only the numerators 1 − √ 6 and 1 + √ 6 need to be interchanged and therefore the complete expression will be omitted. Hence, the two coexisting stable waves that arise from the two Hopf bifurcations of the equilibria x P ,1 and x P ,2 have the same spatiotemporal properties, but they differ in the spatial phase which is indeed visible in the Hovmöller diagrams in Fig. 5. These results show how the pitchfork bifurcation changes the mean flow and hence also the propagation characteristics of the wave. In the next section we will explore spatiotemporal properties of waves for F < 0 and other even values of n.

Stationary waves for even n ≥ 4 and F < 0
The case n = 4 turns out to be more complicated than the case n = 6. If n = 4, then the equilibrium x F = (F, F, F, F ) undergoes a pitchfork bifurcation at F = − 1 2 since λ 2 = 0. Just as in the case n = 6 two new branches of equilibria appear which are given by where a, b are again given by Eq. (10). The Jacobian matrix at the equilibrium x P ,1 is given by For F = −3, in which case a = (−1+ √ 5)/2 and b = (−1− √ 5)/2, the characteristic polynomial of the matrix J is given by which implies that a real eigenvalue of J becomes zero at F = −3. For the equilibrium x P ,2 we obtain the same result. Since the equilibria x P ,1 and x P ,2 continue to exist for F < −3 a saddle-node bifurcation is ruled out. Numerical continuation using the software package AUTO-07p (Doedel and Oldeman, 2007) shows that again a pitchfork bifurcation takes place at F = −3. It is not feasible to derive analytic expressions for the new branches of equilibria as in Eq. (10). Continuation of the four branches while monitoring their stability indicates that at F ≈ −3.853 in total four Hopf bifurcations occur (one at each branch). Figure 6 shows the bifurcation diagrams for the cases n = 4 and n = 6. The question is whether the results described above persist for even dimensions n > 6. To that end we conducted the following numerical experiment. For all even dimensions 4 ≤ n ≤ 50 we used the software package AUTO-07p to numerically continue the equilibrium x F = (F, . . ., F ) for F < 0 while monitoring the eigenvalues to detect bifurcations. At each pitchfork bifurcation we performed a branch switch in order to follow the new branches of equilibria and detect their bifurcations. Once a Hopf bifurcation is detected we can compute the period of the wave as T = 2π/ω from the eigenvalue pair ±ωi. The results of this experiment reveal that the cases n = 4k and n = 4k + 2 are different both qualitatively and quantitatively.
If n = 4k + 2 for some k ∈ N, then one pitchfork bifurcation occurs at F = −0.5. This follows directly from Eq. (4) for the eigenvalues of the equilibrium x F : for even n we have λ n/2 = −1 − 2F which changes sign at F = − 1 2 . From the pitchfork bifurcation two new branches of stable equilibria emanate. Each of these equilibria is of the form with a > 0 and b < 0; the other equilibrium just follows by interchanging a and b. Each of the two equilibria undergoes a Hopf bifurcation, which leads to the coexistence of two stable waves. Figure 7a suggests that the value of F at which this bifurcation occurs is not constant, but tends to −3 as n → ∞. The period of the periodic attractor that is born at the Hopf bifurcation increases almost linearly with n: fitting the function T (n) = α+βn to the numerically computed periods gives α = 0.36 and β = 0.59 (see Fig. 7b). If n = 4k for some k ∈ N, then two Pitchfork bifurcations in a row occur at F = −0.5 and F = −3. After the second pitchfork bifurcation there are four branches of equilibria. Each of these equilibria is of the form (a, b, c, d, a, b, c, d, . . .), where a, b, c, and d alternate in sign; the other equilibria are obtained by applying a circulant shift. Each of the four stable equilibria undergoes a Hopf bifurcation at the same value of the parameter F , which leads to the coexistence of four stable waves. Figure 7a suggests that the value of F at which this bifurcation occurs is not constant, but tends to −3.64 as n → ∞. Contrary to the case n = 4k+2, the period of the periodic attractor that appears after the Hopf bifurcation settles down and tends to 1.92 as n → ∞.
In spite of the aforementioned quantitative differences between the cases n = 4k + 2 and n = 4k, the wave numbers depend in the same way on n in both cases. Equations (12) and (13) show that the n components of the equilibrium that undergoes the Hopf bifurcation alternate in sign. Therefore, sufficiently close to the Hopf bifurcation the components x 0 (t), . . ., x n−1 (t) of the periodic orbit will also alternate in Pitchfork Hopf Figure 6. Bifurcation diagrams obtained by continuation of the equilibrium x F = (F, . . ., F ) for F < 0 for the dimensions n = 4 and n = 6. Stable (unstable) branches are marked by solid (dashed) lines. For n = 4 two pitchforks in a row occur before the Hopf bifurcation, whereas for n = 6 only one pitchfork occurs before the Hopf bifurcation. The bifurcation diagram for n = 4k with k ∈ N (resp. n = 4k + 2) is qualitatively similar to the bifurcation diagram for n = 4 (resp. n = 6; see the main text). . Parameter values of the first Hopf bifurcation (a) and the periods of the periodic attractor (b) that appears after the Hopf bifurcation for F < 0 and even values of the dimension n. For clarity the cases n = 4k and n = 4k + 2 have been marked with different symbols in order to emphasize the differences between the two cases.
sign. Hence, the resulting stationary waves consists of n/2 "troughs" and "ridges", which means that their wave number equals n/2.

Multi-stability: coexistence of waves
The results of Sect. 3.4 show that for even n and F < 0 either two or four stable periodic orbits coexist for the same parameter values. This phenomenon is referred to as multistability in the dynamical systems literature. An overview of the wide range of applications of multi-stability in different disciplines of science is given by Feudel (2008).
Multi-stability also occurs when F > 0, but for a very different reason. For n = 12, Theorem 1 implies that the first bifurcation of the equilibrium x F = (F, . . ., F ) for F > 0 is not a Hopf bifurcation, but a double-Hopf bifurcation. Indeed, at F = 1 we have two pairs of purely imaginary eigenvalues, namely (λ 2 , λ 10 ) = (−i √ 3, i √ 3) and (λ 3 , λ 9 ) = (−i, i). Note that the double-Hopf bifurcation is a codimension-2 bifurcation which means that generically two parameters must be varied in order for the bifurcation to occur (Kuznetsov, 2004). However, symmetries such as those in the Lorenz-96 model can reduce the codimension of a bifurcation.
In previous work (Van Kekem and Sterk, 2018) we have introduced an embedding of the Lorenz-96 model in a two-  Figure 8. Bifurcation diagram of the two-parameter system (Eq. 14) in the (F, G) plane for n = 12. A double-Hopf bifurcation point is located at the point (F, G) = (1, 0) due to the intersection of two Hopf bifurcation lines. From this codimension-2 point two Neȋmark-Sacker bifurcation curves emanate which bound a lobeshaped region in which two periodic attractors coexist.
parameter family by adding a diffusion-like term multiplied by an additional parameter G: Note that by setting G = 0 we retrieve the original Lorenz-96 model in Eq. (1). Since the Jacobian matrix of Eq. (14) is again a circulant matrix we can use Eq. (4) to determine its eigenvalues: λ j = −1 − 2G(1 − cos(2πj/n)) + Ff (2πj/n) Also note that x F = (F, . . ., F ) remains an equilibrium solution of Eq. (14) for all (F, G). The Hopf bifurcations of x F described in Theorem 1 now occur along the lines and the intersection of two such lines leads to a double-Hopf bifurcation. Figure 8 shows a local bifurcation diagram of the twoparameter Lorenz-96 model in the (F, G) plane for n = 12 which was numerically computed using MATCONT (Dhooge et al., 2011). A double-Hopf point is located at (F, G) = (1, 0), which is indeed implied by Theorem 1. The normal form of a double-Hopf bifurcation depends on the values of two coefficients which determine the unfolding of the bifurcation. In total, there are 11 different bifurcation scenarios to consider. The normal form computation of the double-Hopf bifurcation for n = 12 and (F, G) = (1, 0) in Van Kekem and Sterk (2018) shows that the unfolding of this particular case is of "type I in the simple case" as described by Kuznetsov (2004). This means that from the double-Hopf point only two curves of Neȋmark-Sacker bifurcations emanate. In discrete-time dynamical systems, a Neȋmark-Sacker bifurcation is the birth of a closed invariant curve when a fixed point changes stability through a pair of complex eigenvalues crossing the unit circle in the complex plane. From a continuous-time system, such as Eq. (14), we can construct a discrete-time dynamical system by defining a Poincaré return map of a periodic orbit, in which case a Neȋmark-Sacker bifurcation refers to the birth of an invariant two-dimensional torus when the periodic orbit changes stability by a pair of Floquet multipliers crossing the unit circle in the complex plane.
In order to explain the dynamics in a neighborhood around the double-Hopf point, we now use Fig. 8 to describe the successive bifurcations that occur for G = 0.1 fixed and increasing F . At F = 1.1 the equilibrium x F becomes unstable through a supercritical Hopf bifurcation (the blue line given by G = F − 1) and a stable periodic orbit with wave number 2 is born. At F = 1.2 the now unstable equilibrium undergoes a second Hopf bifurcation (the red line given by G = 1 2 (F −1)) and an unstable periodic orbit with wave number 3 is born. The latter periodic orbit becomes stable at F ≈ 1.58 through a subcritical Neȋmark-Sacker bifurcation (orange curve) and an unstable two-dimensional invariant torus is born. Hence, for parameter values F > 1.58 two stable waves with wave numbers 2 and 3 coexist until one of these waves becomes unstable in a bifurcation. For all fixed values of 0 < G < 0.15 the same bifurcation scenario occurs, but the values of F are different. For −0.06 < G < 0 the roles of the two Hopf bifurcations and periodic orbits have to be interchanged.
The scenario described above shows how the presence of two subcritical Neȋmark-Sacker bifurcations emanating from a double-Hopf bifurcation determines a region of the (F, G) plane in which two stable periodic orbits coexist with an unstable two-dimensional invariant torus. We will refer to this region as the "multi-stability lobe". The scenario described above is not limited to the special case of the Lorenz-96 model, but occurs near a double-Hopf bifurcation of type I in any dynamical system (Kuznetsov, 2004).
Double-Hopf bifurcations are abundant in the twoparameter . The lines described in Eq. (16) have a different slope for all 0 < j < n/2 and j = n/3, and hence they mutually intersect each other. This implies that the number of double-Hopf points in the (F, G) plane grows quadratically with n, see Appendix A3. However, not all these points will have an influence on the dynamics: if x F is already unstable, then any dynamical object born through the double-Hopf bifurcation will also be unstable. In what follows, we only consider the double-Hopf bifurcations through which x F can change from stable to unstable. We can find such points as follows. Starting from the line in Eq. (16) with j = l + 1 (n) as defined by Eq. (8), we first compute double-Hopf points by computing the intersections with all other lines. From these intersections we select those that satisfy the condition max{Reλ j : j = 0, . . ., n − 1} = 0. Figure 9 shows the G coordinates of these double-Hopf points as a function of n. Clearly, for large n there exist double-Hopf points which are very close to the F axis, which suggests that the multi-stability lobe that emanates from such points can intersect the F axis and hence influence the dynamics of the original Lorenz-96 model for G = 0. Moreover, Fig. 9 shows that for n > 12 there are always two double-Hopf points by which x F can change from stable to unstable. It is then possible that two multi-stability lobes intersect each other, which leads to a region in the (F, G) plane in which at least three stable waves coexist. Figures 10-12 show bifurcation diagrams of three periodic orbits as a function of F for G = 0 for n = 40, 60, 80. For each periodic orbit the continuation is started from a Hopf bifurcation of the equilibrium x F . If x F is unstable, then so will be the periodic orbit. However, when the boundary of a multi-stability lobe is crossed, a Neȋmark-Sacker bifurcation occurs, by which a periodic orbit can gain stability. For specific intervals of the parameter F , three stable periodic orbits coexist. Since Fig. 9 shows that for large values of n the double-Hopf bifurcations are close to the F axis, we expect that the coexistence of three or more stable waves is typical for the Lorenz-96 model.
The double-Hopf bifurcation has been reported in many works on fluid dynamical models. A few examples are baroclinic flows (Moroz and Holmes, 1984), rotating cylinder flows (Marqués et al., 2002(Marqués et al., , 2003, Poiseuille flows (Avila et al., 2006), rotating annulus flows (Lewis and Nagata, 2003;Lewis, 2010), and quasi-geostrophic flows (Lewis and Na- The period of the orbit is plotted as a function of F . Stable (resp. unstable) orbits are indicated by solid (resp. dashed) lines. Circles denote Neȋmark-Sacker bifurcations and triangles denote period doubling bifurcations. The Hopf bifurcations generating the waves with wave numbers 8, 9, and 7 occur at respectively F = 0.894, F = 0.902, and F = 0.959. Clearly, for 1.15 < F < 2.79 three stable periodic orbits coexist. gata, 2005). In all of these examples, the coexistence of multiple waves is reported, where the nature of these waves depends on the specific model.

Conclusions
In this paper we have studied spatiotemporal properties of waves in the Lorenz-96 model and their dependence on the dimension n. For F > 0 the first bifurcation of the equilibrium x F = (F, . . ., F ) is either a supercritical Hopf or a double-Hopf bifurcation and the periodic attractor born at the Hopf bifurcation represents a traveling wave. The spatial wave number is determined by the index of the eigenpair that Period F Wave nr. 13 Wave nr. 12 Wave nr. 14 Figure 11. As Fig. 10, but for n = 60. For 1.01 < F < 2.03 three stable periodic orbits coexist. The Hopf bifurcations generating the waves with wave numbers 13, 12, and 14 occur at respectively F = 0.891, F = 0.894, and F = 0.923. crosses the imaginary axis and increases linearly with n, but the period tends to a finite limit as n → ∞. For F < 0 and n odd, the first bifurcation of x F is always a supercritical Hopf bifurcation and the periodic attractor that appears after the bifurcation is again a traveling wave. In this case the wave number equals (n − 1)/2 and the period is O(4n).
For n even and F < 0 the first bifurcation of x F is a pitchfork bifurcation which occurs at F = − 1 2 and leads to two stable equilibria. If n = 4k + 2 for some k ∈ N, then each of these equilibria undergoes a Hopf bifurcation which leads to the coexistence of two stationary waves. The role of the pitchfork bifurcation is to change the mean flow which in turn changes the propagation of the wave. If n = 4k for some k ∈ N, then two pitchfork bifurcations take place at F = − 1 2 and F = −3 before a Hopf bifurcation occurs, which leads to the coexistence of four stationary waves.
The occurrence of pitchfork bifurcations before the Hopf bifurcation leads to multi-stability, i.e., the coexistence of different waves for the same parameter settings. A second scenario that leads to multi-stability is via the double-Hopf bifurcation. For n = 12 the equilibrium x F loses stability through a double-Hopf bifurcation. By adding a second parameter G to the Lorenz-96 model we have studied the unfolding of this codimension-2 bifurcation. Two Neȋmark-Sacker bifurcation curves emanating from the double-Hopf point bound a lobe-shaped region in the (F, G) plane in which two stable traveling waves with different wave numbers coexist. For dimensions n > 12 we find double-Hopf bifurcations near the F axis, which can create two multistability lobes intersecting each other, and in turn this can lead to the coexistence of three stable waves coexisting for G = 0 and a range of F values. Hence, adding a parameter G to the Lorenz-96 model helps to explain the dynamics, which is observed in the original model for G = 0.
Our results provide a coherent overview of the spatiotemporal properties of the Lorenz-96 model for n ≥ 4 and F ∈ R. Since the Lorenz-96 model is often used as a model for testing purposes, our results can be used to select the most appropriate values of n and F for a particular application. The periodic attractors representing traveling or stationary waves can bifurcate into chaotic attractors representing irregular versions of these waves, and their spatiotemporal properties are inherited from the periodic attractor (see for example Figs. 1) and 2. This means that our results on the spatiotemporal properties of waves apply to broader parameter ranges of the parameter F than just in a small neighborhood of the Hopf bifurcation.
The results presented in this paper also illustrate another important point: both qualitative and quantitative aspects of the dynamics of the Lorenz-96 model depend on the parity of n. This phenomenon also manifests itself in discretized partial differential equations. For example, for discretizations of Burgers' equation, Basto et al. (2006) observed that for odd degrees of freedom the dynamics were confined to an invariant subspace, whereas for even degrees of freedom this was not the case. For the Lorenz-96 model the parity of n also determines the possible symmetries of the model. We will investigate these symmetries and their consequences on bifurcation sequences using techniques from equivariant bifurcation theory in forthcoming work .
Code availability. The scripts used for continuation with AUTO-07p are available upon request from Alef Sterk.