Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere
EEPJ manuscript No. (will be inserted by the editor)
Mechanics and Thermodynamics of a New Minimal Model ofthe Atmosphere
Gabriele Vissio , , and Valerio Lucarini , CEN, Meteorological Institute, University of Hamburg, Hamburg, Germany Institute of Geosciences and Earth Resources (IGG) - National Research Council (CNR), Torino, Italy Department of Mathematics and Statistics, University of Reading, Reading, UK Centre for the Mathematics of Planet Earth, University of Reading, Reading, UKReceived: October 19, 2020
Abstract.
The understanding of the fundamental properties of the climate system has long benefitted from the use ofsimple numerical models able to parsimoniously represent the essential ingredients of its processes. Here we introducea new model for the atmosphere that is constructed by supplementing the now-classic Lorenz ’96 one-dimensionallattice model with temperature-like variables. The model features an energy cycle that allows for conversion betweenthe kinetic and potential forms and for introducing a notion of efficiency. The model’s evolution is controlled by twocontributions - a quasi-symplectic and a gradient one, which resemble (yet not conforming to) a metriplectic structure.After investigating the linear stability of the symmetric fixed point, we perform a systematic parametric investigationthat allows us to define regions in the parameters space where at steady state stationary, quasi-periodic, and chaoticmotions are realised, and study how the terms responsible for defining the energy budget of the system depend onthe external forcing injecting energy in the kinetic and in the potential energy reservoirs. Finally, we find preliminaryevidence that the model features extensive chaos. We also introduce a more complex version of the model that is able toaccommodate for multiscale dynamics and that features an energy cycle that more closely mimics the one of the Earth’satmosphere.
PACS. − a Nonlinear dynamics and chaos – 92.60.Bh General circulation – 05.70.Ln Nonequilibrium andirreversible thermodynamics – 02.60.Cb Numerical simulation; solution of equations The climate is a non-equilibrium system whose dynamics is primarily driven by the uneven absorption of solar radiation, whichis mainly absorbed near the surface and in the tropical latitudes, rather than aloft and in the mid-high latitudes, respectively. Thesystem reacts to such an inhomogeneity in the local energy input through a complex set of instabilities and feedbacks affectingits dynamical processes and thermodynamic and radiative fluxes. Such processes lead to an overall reduction in the temperaturegradients inside the the system and allow for the establishment of approximate steady state conditions [1,2].An example of the re-equilibration mechanism can be described as follows. The large scale energy transport is mainly per-formed by atmospheric disturbances in the form of synoptic and (to a lesser extent) planetary eddies which are, in turn, fuelledby the baroclinic conversion of available potential energy into kinetic energy, which is then dissipated by friction. In turn, thepresence of large low-high latitudes temperature gradients is responsible for the presence of a reservoir of available potentialenergy, which is continuously replenished thanks to the dishomogeneity of the radiative energy budget across the globe [1]. Thisis the core of the celebrated Lorenz energy cycle, which provides a powerful representation of the climate as an engine [3]. Thethermodynamic viewpoint on the climate allows to define its efficiency and irreversibility [4,5,6,7,8,9].The reconstruction, interpretation, and analysis of observative data (now facilitated by the recent advances in data science);analytical tools borrowed (mostly) from mathematics and physics; and numerical simulations all contribute to our understandingof the climate. This task is exceedingly demanding since the system features non-trivial variability on a vast range of temporaland spatial scales and, furthermore, our ability to observe it has changed enormously over time. Additionally, the presence ofperiodic as well as irregular fluctuations in the boundary conditions do not allow the climate to reach an exact steady state [10,11]. a Corresponding author.
Contact Email: [email protected] a r X i v : . [ n li n . C D ] O c t Gabriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere
One of the features of the numerical investigation of the climate system is the reliance on hierarchies of models. In otherterms, climate phenomena are investigated using a full range of models going from low dimensional ones to state-of-the-artEarth system models, able to represent with higher precision many aspects of the climate. See a discussion of the meaning andof the use of hierarchies of climate models in [12,13,14,11]. It is important to remark that, given the multiscale nature of theclimate system, the heterogeneity of its subdomains, and the number of the active physical, chemical, and biological processes,the endeavour of construcing a model able to directly simulate all of them appears as a Sisyphean task, whilst, instead, theparametrization of the effect of the unresolved scales on those that are explicitly simulated is an essential component of anyreasonable model of the Earth system[15,16,17].Low-order models have played and still play today a very important role for improving our understanding of the geophysicalflows. Apart from the landmark 3-dimensional model developed by Lorenz in 1963 [18] starting from the truncation of theequations describing the Rayleigh-Benard convection introduced by Saltzman [19], simple models have been key to scientificadvances in oceanography [20,21,22], dynamical meteorology [23,24,25], climate dynamics [26,27,28,29], turbulence [30,31,32,33] and convection [34,35], among others. Additionally, low-order models supplemented by stochastic forcings have alsoprovided the backbone of the stochastic theory of climate [36,37,38]. Aside from sheer mathematical-related curiosity, many ofthese models were created in order to shed some light on specific problems by using a physically-meaningful benchmark toolthat is easier to analyse mathematically and faster to simulate numerically with respect to more complex models.
Of special relevance for the present study is the now-celebrated Lorenz ’96 model [25,39] (hereafter L96), whose structure isbriefly recapitulated below. The model consists of a lattice of N gridpoints, whose state is described by a variable. The model hasperiodic boundary conditions, so that it can loosely be interpreted as describing the properties of the atmosphere along a latitudinalcircle. The model features in an extremely simplified - almost metaphorical - way the main processes of the atmosphere: forcing,dissipation, and advection. Two versions of the model have been proposed: a one-level version, where the dynamics takes place ona single length scale, and a two-level version, where the lattice is augmented in order to describe dynamical processes on smallerspatial and temporal scales. The two-level version of the L96 model has the especially attractive property that the time-scaleseparation between the fast and the slow variables can be controlled by modulating one parameter.The L96 model has rapidly gained relevance among geoscientists, physicists, and applied mathematicians, as it has become abenchmark testbed for parametrizations [40,41,42,43,44,45,46], for studying extreme events [47,48,49,50], for developing dataassimilation schemes [51,52,53,54], for developing ensemble forecasting techniques [55,56,57], for studying the properties ofLyapunov exponents and covariant Lyapunov vectors [58,59,60,61], for developing and testing ideas in nonequilibrium statisticalmechanics [62,63,64,65,66], and for investigating bifurcations [67,68,69,70,71,72,73]. By looking at these references, thereader can find a very thorough analysis of the properties of the L96 model.The one-level L96 model can be written as: dX k dt = X k − ( X k + − X k − ) − γ X k + F , (1)with boundary conditions: X k − N = X k + N = X k . (2) k = , . . . , N is the index of the gridpoints defining the lattice, the nonlinear term defines a nontrivial process of advection, F isan external forcing, and γ (usually taken with unitary value) modulates the dissipation. In the unforced and inviscid regime - i.e.setting F = γ =
0, the energy of the system, expressed as the sum of the squares of the variables, is conserved: dEdt = ddt N ∑ k = X k = N ∑ k = X k dX k dt = . (3)If γ = N (cid:29)
1, the model’s attractor is a fixed point for 0 ≤ F ≤ /
9. As F is increased, the fixed point X = X = . . . = X K = F loses stability as the system undergoes bifurcations leading to a quasi-periodic behavior for moderate values of F and chaoticbehaviour for F ≥ . F is changed: as discussed in detail in [70,71,72], the properties of the system depend on N in a very nontrivial way in theregime of moderate forcing. In the regime of strong forcing and developed turbulence, instead, some sort of universality emerges,as the L96 model is strongly chaotic when F ≥ N [65]. Byand large, the mechanism of instability of the L96 model boils down to exchanges of energy between the symmetric state andthe perturbations away from it, as particularly clear in the case of the linear instability analysis [39]. Nonetheless, the L96 modelclearly features only one form of energy, which we may refer to as kinetic . abriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere 3 In this paper we propose an extension of the L96 model whereby a second variable is attached to each gridpoint, representing,metaphorically, the local thermodynamical properties, which are advected by the dynamical variables of the L96 model, andundergo forcing and diffusion. The model proposed here features a meaningful definition of energy that includes the kineticpart already present in the L96 model plus a potential part associated to the fluctuations of the temperature in the domain. Thefundamental advantage of the new model proposed here and analysed in detail in Sect. 2 is that it features an energy cycle thatallows for conversion between the kinetic and potential forms of energy. Conceptually, this mirrors the change one has going fromone a one-layer quasi-geostrophic model, where only barotropic processes are possible, to a two-layer quasi-geostrophic model,which, instead, features a coupling between dynamical and thermodynamic processes via baroclinic conversion [74]. Note that theLorenz ’63 model [18] and, more completely so, its extensions to higher order modes [75] features a nontrivial energetics whereexchanges take place between potential and kinetic energy [76]. The energy of the system defines the symplectic component thatcontributes - together with the metric one - to defining the evolution of the model [47].The rest of the paper is structured as follows. Section 2 provides a thorough introduction to the one-level version of the modelpresented here. The evolution equations are presented together with the rationale on which the model is based. Additionally,a detailed analysis of its mechanical and thermodynamic properties is carried out. Finally, the linear stability analysis of thesymmetric fixed point is also described. In Sect. 3 we present the results of a large set of numerical integrations of the model,aimed at exploring its properties in a rather vast range of values of two main parameters, which control the input of energy in thekinetic and potential form. We discuss the thermodynamics of the model in terms of mean values and the fluctuations of the mainterms describing the energetics of the model. Such a physical characterisation of the model is complemented by the analysis ofhow the first Lyapunov exponent depends on the two considered parameters, in order to be able to separate the regions wherethe asymptotic dynamics of the system takes place in a regular vs in a strange attractor [77,67] We will discover a non-trivialinterplay between the two sources of forcings applied to the model. We then perform a preliminary analysis for assessing to whatextent the system obeys extensive chaos. In Sect. 4 we summarise the main features of the model and the results obtained so far,and propose future lines of investigations. As in the case of original L96 model, the model introduced here can be formulatedin a two-level fashion - see Appendix A -, with non trivial couplings among different levels and variables and with a fairlysophisticated energetics, which is conceptually rather similar to the one described by the Lorenz energy cycle in the atmosphere.The analysis of the properties of the two-level model will not be performed in this paper and will be the subject of future studies.
We want to extend the standard L96 model presented in the introduction by adding a second set of variables for all the gridpoints k = , . . . , N . The goal is to construct a toy model able to describe in a very simple yet conceptually correct way the interactionbetween dynamical and thermodynamical processes of the atmosphere. The evolution equations of the model we propose in thiscontribution are the following: dX k dt = X k − ( X k + − X k − ) − αθ k − γ X k + F , (4) d θ k dt = X k + θ k + − X k − θ k − + α X k − γθ k + G , (5)with k = , . . . , N . The variable θ k can be loosely interpreted as temperature at the grid-point k . The boundary conditions aredefined as X k − N = X k + N = X k , θ k − N = θ k + N = θ k , (6)The variable θ k undergoes a constant forcing G , a linear dissipation term, and a nonlinear term representing, loosely speaking,the advection performed by the X variables. Additionally, θ k and X k are linearly coupled through a term proportional to α . Thepurpose of this coupling is to represent, in a very simplified way, the effect of correlated thermal and dynamical fluctuationson the dynamics, which allow for an exchange between kinetic and potential energy associated with thermal fluctuations, asdiscussed below. The introduction of a term proportional to α is the only - yet important - modification in the dynamics of the X variables for this model as compared to the classical L96 model, see Eq. 1. In what follows, we consider F , G , α , γ ≥ X and the θ variables is constructed in such a way that in the unforced and inviscid limit ( F = G = γ =
0) the total energy of the system E = K + P = N ∑ k = ( X k + θ k ) Gabriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere given by the sum of its kinetic and potential components, is conserved: dEdt = dKdt + dPdt = N ∑ k = X k dX k dt + N ∑ k = θ k d θ k dt = , (7)whereas, in general, K and P are not separately conserved. The quadratic functional form of the potential energy is inspired by thefact that the available potential energy in the global circulation of the atmosphere is approximately proportional to the varianceof the temperature fluctuations [3,1,78]. The dynamical role of the function E is explored in the next section. By definition, the time derivative of any smooth observable Ψ ( X , . . . , X K , θ , . . . , θ k ) is obtained by applying the generator of theKoopman operator to Ψ , as follows: d Ψ d t = L [ Ψ ] . We will now show that linear operator L can be written as the sum of a contribution coming from a symplectic (indeed, quasi-symplectic, for the reasons detailed at the end of this section) term and a contribution coming from a gradient term. Indeed, wecan write: d Ψ d t = L [ Ψ ] = { Ψ , E } + (cid:104) Ψ , Γ (cid:105) (8)where { A , B } = − { B , A } is a suitably defined Poisson bracket for the functions A and B , while (cid:104) A , B (cid:105) = (cid:104) B , A (cid:105) gives the gradientcontribution. The evolution equations 4-5 are obtained by setting Ψ = X i and Ψ = θ i , i = , . . . , N , respectively.We have that (cid:104) A , B (cid:105) = ( ∂ X i A ) ( ∂ X i B ) + (cid:0) ∂ θ i A (cid:1) (cid:0) ∂ θ i B (cid:1) , where we use the Einstein convention for the indices. The function Γ defining the gradient contribution to the dynamics is: Γ = − γ E + F N ∑ k = X k + G N ∑ k = θ k = − γ E + F Ξ + G Θ , (9)where Ξ = ∑ Nk = X k and Θ = ∑ Nk = θ k . it is clear that such a component describes the irreversible dynamics as it vanishes in theunforced, inviscid limit F = G = γ = { A , B } = ( ∂ X i A ) J i j (cid:0) ∂ X j B (cid:1) + (cid:0) ∂ θ i A (cid:1) Y i j (cid:16) ∂ θ j B (cid:17) + ( ∂ X i A ) L i j (cid:16) ∂ θ j B (cid:17) + (cid:0) ∂ θ i A (cid:1) M i j (cid:0) ∂ X j B (cid:1) (10)where J i j = X i − δ i + , j − X j − δ i , j + Y i j = X i + δ i + , j − X i − δ i − , j L i j = − αδ i , j M i j = αδ i , j (11)It is clear that the energy E is the generator of time translations according to the symplectic contribution and the antisymmetryof the Poisson brackets enforces the corresponding conservation law already discussed in Eq. 7.We remark that, in the inviscid and unforced limit the system is not Hamiltonian because the Poisson brackets do not fulfillthe Jacobi identity { A , { B , C }} + { C , { A , B }} + { B , { C , A }} =
0. Because of this and of the fact that { Γ , E } (cid:54) =
0, the systemgiven in Eqs. 4-5 is not metriplectic, i.e. the standard generalisation of Hamiltonian system to the dissipative case. [79]. Notethat the dynamics of dissipative fluids is, instead, metriplectic [80], and so is the dynamics of the (extended) Lorenz ’63 model,which is in fact derived from the Rayleigh-Bénard equations through systematic modal truncation [47]. The lack of an underlyingHamiltonian skeleton confirms the well-known fact that the L96 model cannot be easily related to any model of fluid flows.
Using Eqs. 8-9 we obtain the time evolution of the energy of the system:d E d t = (cid:104) Ψ , Γ (cid:105) = − γ E + F Ξ + G Θ = Γ − γ E abriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere 5 KP I P I K D P D K C P , K
Fig. 1.
Energetics of the model presented in Eq. 4-5. We indicate the fluxes of energy in and out of the reservoirs of kinetic (K) and potential(P) energy. Dashed lines represent input of energy; dotted lines represent energy dissipation; and solid lines represent energy conversion terms.The direction of the arrows indicates a positive energy flux. See text for details. which implies that, at steady state 2 γ ¯ E = F ¯ Ξ + G ¯ Θ , where ¯ Ψ is the long term average of the quantity Ψ .We next analyse the separate budget of the kinetic and potential energy. By inserting K and P in Eq. 8 we obtain:d K d t = N ∑ k = X k d X k d t = I K − D K + C P , K I K = F Ξ , D K = γ K , C P , K = − N ∑ k = ( α X k θ k ) , (12)d P d t = N ∑ k = θ k d θ k d t = I P − D P − C P , K I P = G Θ , D P = γ P , (13)where I K ( I P ) is the rate of input of kinetic (potential) energy, D K ( D P ) is the dissipation rate of kinetic (potential) energy, and C P , K the conversion rate from potential to kinetic energy. We remark that the input and dissipation of energy in either kineticor potential form is due by the metric component of the dynamics. Instead, the conversion of energy between the potential andkinetic form is controlled by the Poisson brackets given in Eqs. 10-11. Nonetheless, the components J and Y of the Poissonbrackets, which describe advection, do not give any net contribution. Equations 12-13 describe the energetics of the modelpresented in this work, which is represented by the diagram shown in Fig. 1. One can draw a parallel between the energeticsof this model and the Lorenz energy cycle of the atmosphere, where, as well known, the input of energy comes almost entirelythrough the potential energy channel via baroclinic forcing associated with the differential heating of low versus high latituderegions [3,1]. The two-level version of the model introduced in this paper features an energetics that is conceptually closer to theone of the true atmosphere because it is able to describe energy cascades across scales on top of energy conversion processes,see Appendix A.At steady state conditions one has 2 γ ¯ K = F ¯ Ξ + ¯ C and 2 γ ¯ P = G ¯ Θ − ¯ C , which relate the size of the reservoirs of kinetic andpotential energy to intensity of the acting forcings and energy exchange. We can also introduce a notion of efficiency of thismodel η = ¯ C / ( F ¯ Ξ + G ¯ Θ ) = ¯ C / ( γ ¯ E ) , which relates the amount of energy exchanged between the two reservoirs of energy tothe total energy input. Since X k + θ k ≥ | X k θ k | ∀ k , we have that 2 E ≥ | C | / α . Therefore, | η | ≤ α / γ , which provides a constrainton the efficiency of the system. Note that η is positive if, on the average, energy is converted from potential to kinetic, andnegative otherwise.Note that K ≥ P ≥ G , γ > F =
0, one has ¯ C ≥ , whereas if F , γ > G =
0, ¯ C ≤ F = G =
0, instead, we have that ¯ K = ¯ P = ¯ E =
0, where the origin is a stable fixed point for the system. Note that this is, to a very good approximation, what applies to the climate system as a whole, because the geophysical fluids do not receiveany input of mechanical energy, apart from the very small lunar and solar tidal forcing. Gabriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere
We investigate the linear stability of the system analysed here around the fixed point corresponding to the symmetric solution X j = X k = const , 1 , . . . , j , k , . . . N and θ j = θ k = const , 1 , . . . , j , k , . . . N . By plugging this ansatz in Eqs. 4-5 one gets: X k = ˜ X = γ F − α G γ + α ∀ k , (14) θ k = ˜ θ = α F + γ G γ + α ∀ k . (15)Taking inspiration from [39] , we then investigate the linear stability of this solution by substituting X k = ˜ X + A exp ( σ t ) exp ( i κ k − i ω t ) and θ k = ˜ θ + B exp ( σ t ) exp ( i κ k − i ω t ) in Eqs. 4-5, where ˜ X and ˜ θ have been defined in Eqs. 14 and 15, respectively; A and B are complex constants, σ is a real number defining the growth rate (if positive) of the amplitude of the wave, whilst κ is thewavenumber and ω is the angular frequency of the wave. Neglecting terms that are quadratic in the wave amplitude, one obtains: ( σ − i ω ) A = ˜ XA ( exp ( i κ ) − exp ( − κ )) − α B − γ A (16) ( σ − i ω ) B =
2i ˜ XB sin ( κ ) +
2i ˜ θ A sin ( κ ) + α A − γ B , (17)We exclude the trivial solution A = B = A = b = B / A is indeed relevant). Weseparate real and imaginary part in the previous equations and obtain: σ = ˜ X ( cos ( κ ) − cos ( κ )) − α Re { b } − γ (18) ω = − ˜ X ( sin ( κ ) + sin ( κ )) + α Im { b } (19) σ Re { b } + ω Im { b } = − X Im { b } sin ( κ ) + α − γ Re { b } , (20) σ Im { b } − ω Re { b } = X Re { b } sin ( κ ) + θ sin ( κ ) − γ Im { b } , (21)where Re { x } and Im { x } are the real and imaginary part of the complex number x , respectively. The conditions leading to thebifurcation associated with the loss of stability of the fixed point given in Eqs. 14-15 can be derived by setting σ = ω , κ , Re { b } and Im { b } and a function of the parameters F , G , α , and γ . In the case ω (cid:54) =
0, the onset of theneutral wave corresponds to a Hopf bifurcation.Solving the previous Eqs. 18-21 and finding the expression of σ and ω as a function of κ and of the parameters F , G , α , and γ gives the dispersion relation of the waves. Additionally, obtaining the real and imaginary part of b allows for understanding therelative amplitude of the waves in the X and θ variables.Note that the linear stability analysis of the L96 model can be obtained by setting α = γ = θ variables. One then recovers the result first presented in [39] and discussed in greater detail in [70]. In whatfollows we consider F ≥
0; an analysis of the somewhat dynamics occurring for F < F such that the fixed point of the system loses stability and, correspondingly, to obtain thewavelength and frequency of the emerging neutral wave. One finds that, taking a continuum approximation ( N → ∞ ), the neutralwave is realised when F = F crit = /
9, where the critical wavenumber is κ = κ crit = arccos ( / ) , and the critical frequency is ω crit = − F crit ( sin ( κ crit )+ sin ( κ crit )) ≈ − .
29. If one assumes that the gridpoints are arranged like along a latitudinal circle wherethe longitude increases with the index of the gridpoints (note that the periodic boundary conditions of the system impose a toroidaltopology), we have that the crest of the neutral wave moves westward, because the phase velocity v p = ω crit / κ crit = ≈ − .
98 isnegative. Instead, the group velocity v g = ∂ ω crit / ∂ κ | κ = κ crit = − F crit ( cos ( κ crit )+ ( κ crit )) ≈ . > X and θ variables, it is hard to find for the model introduced in thispaper an explicit expression for the conditions supporting the presence of a neutral wave, also if one takes the special cases whereone between F and G vanishes. A simple solution is instead found if one takes γ = α = F = G , which implies ˜ X = θ = F . One then obtains the following results when imposing σ = Re { b } = − ω = Im { b } = √ κ = arcsin ( √ / F ) . This indicates that F crit = √
2, and κ crit = π / ω crit = √
2. Therefore, the phase velocity of the neutral wave v p = ω crit / κ crit = √ / π is positive, corresponding to aneastward motion of the wave crests. Since ω crit = F crit sin ( κ crit ) , we have v g = ∂ ω crit / ∂ κ | κ = κ crit = F crit ( cos ( κ crit )) =
0, implyingno net motion of the wave packets. See also the rather sophisticated analysis of the stability of the L96 model presented in [71,70,72].abriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere 7a) b)
Fig. 2.
First Lyapunov exponent of the one-level model with N =
36 gridpoints as a function of F and G . a) Detail for 0 ≤ F , G ≤
3. b) Fullrange for 0 ≤ F , G ≤
10. In all simulations γ = α = Many are the possible scientific questions one can address regarding the model introduced above. Building on the large literatureon the L96 model discussed in the introduction, and taking into account the extra features of the current model, we can mentionthe following lines of investigation: – Analysis of the bifurcations leading the system from fixed point to a periodic and quasi-periodic behaviour to a chaotic regimeas the forcing is increased; – Systematic investigation of the predictability of the system - e.g. analysis of the finite-time and asymptotic Lyapunov expo-nents and the corresponding covariant Lyapunov vectors as a function of the two forcing parameters F and G ; – Systematic investigation of the energetics of the system as a function of the two forcing parameters F and G ; – Analysis of the signal propagation through waves in the quasi-periodic and weakly chaotic regime; – Definition of asymptotic scaling laws for the properties if the system for large values of F and G ; – Detection and analysis of chaos extensivity as the number of gridpoints N → ∞ ; – Extension of the model to multiple scales and analysis of dynamics and of the energetics of scale-to-scale interaction.Obviously, it is impossible to address with a high level of detail all these aspects in the present paper. Rather than focusing on oneor few aspects among those above, since this is the first time this model is proposed to the scientific community, we will presentsome preliminary results that address partially each of the points above, in the hope of stimulating a reader into going in greaterdetail. Further studies by the authors that focus specifically in some of the aspects mentioned above will be reported elsewhere.All the numerical integrations are performed using a Dormand–Prince method with adaptive time step and a spin time of 100time units, with runs of 1000 time units. We make use of the Python module JiTCODE [81], an extension of SciPy’s ODE thatallows to numerically simulate ordinary differential equations, computing quantities of interest as Lyapunov exponents as well.All results have been double checked and confirmed using the MATLAB function ode45 where integrations are performed usingthe 4 th order Runge-Kutta integrator with adaptive time step [82]. A simple yet fundamentally correct way to characterise at qualitative level the dynamical properties of a system is to investigateto what extent its evolution is sensitive to its initial conditions. Roughly speaking, the first Lyapunov exponent of a systemmeasures the asymptotic rate of growth or decay of the distance between two orbits which are initialised in the attractor of thesystem at infinitesimal distance from each other [83]. Similarly, one can define the sum of the first p Lyapunov exponents asdefining the asymptotic average rate of growth or decay of the p -volume defined by p + Q dimensional continuous time dynamical system, it ispossible to compute Q Lyapunov exponents λ , . . . , λ Q , where the customary ordering is such that λ j ≤ λ k if k ≤ j of the [84].If λ < λ = Gabriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere quantitatively the rapidity with which two nearby trajectories diverge from each other. In this case, one has that there is at leastone λ j = j >
1, which corresponds to the direction of the flow [77,83].Figure 2 shows the estimate of λ for N =
36 as a function of F and G in the range 0 ≤ F , G ≤
10. We remark that the PythonJiTCODE module allows for the computation of the full spectrum of Lyapunov exponents using the algorithm proposed in [84].The system has a negative λ for small values of the forcings, as expected, see Fig. 2a). We remark that if G = λ ≤ F ≤ .
4, whereas for the L96 model F crit = /
9, indicating that presence of a mechanism of energy transfer between kinetic andpotential energy and the presence of a new channel of dissipation (for potential energy) leads to higher stability for the system.We observe that λ depends in a very nontrivial way on both F and G , as the system’s behaviour depends delicately on how theenergy is injected into it, because the dynamics of the X and θ variables is, in fact, quite distinct. It is extremely different toforce the system through kinetic or the potential energy channel. We also observe that the theoretical prediction of F crit = √ F = G agrees with what shown in Fig. 2a).Many other interesting features appear. Increasing F from zero to 3 while keeping G = .
7, the asymptotic dynamics ofthe system changes first from quasi-periodic to a fixed point, then again to quasi-periodic, which then alternates with chaoticbehaviour. Indeed, one can observe two complex tongue-like structures in Fig. 2a) for F ≥ . G ≥ .
5, which indicate thepresence of a very nontrivial set of bifurcations for that regions of the parameters’ space, defining the transition between thequasi-periodic behaviour - the light orange region - ad the chaotic regime - the dark orange and red region in Fig. 2a).Zooming out towards a larger range of values for F and G the intuitive argument that increasing either F or G makes thesystem less predictable becomes more robust, even though there are regions where a destructive interference is clear (in terms ofvalues of λ between the two forms of forcing, compare the two troughs near the diagonal in Fig. 2b).We remark that it is reasonable to expect that, as in the case of the L96 model [71,70,72], in the regime of moderate forcing theposition and nature of the bifurcations will depend delicately on the number of nodes N , so that one should expect modificationsespecially in Fig. 2a) when performing simulations for a value of N other than 36 considered here. Instead, as shown below inSect. 3.4, one finds some indication of universality associated with the continuum limit N → ∞ when sufficiently strong forcingis considered. It is useful to investigate the long-term average of the terms in Eqs.12-13 as a function of F and G , see Fig.3. The lack ofequivalence between applying forcing to the X vs to the θ variables is extremely clear by looking at the conversion term (panele). ¯ C is positive for large values of G and moderate values of F , and negative viceversa. The absolute value of ¯ C increases with F ( G ) if G ( F ) is kept constant. The zero isoline strongly deviates from the diagonal and indicate that if F = G there is a net transferof energy from kinetic to potential. The zero isoline of ¯ C coincides with the ridge in the value of λ shown in 2b), indicating thatthe condition of no net energy exchange between the two reservoirs of energy corresponds to a state where instabilities are ratherstrong. The zero-isoline of the efficiency η (see panel f), by definition, coincides with the one of ¯ C . The absolute value of theefficiency grows with the asymmetry of the forcing, and peaks for moderate intensity of either F or G , suggesting - see Sect. 3.5below - that the energy conversion becomes less efficient decreases when stronger forcings are considered.The behaviour of the other thermodynamical quantities is somehow unsurprising, as we have that both input and dissipationof kinetic (potential) energy increase with F ( G ). We remark that, once again, the response of the system to the two individualforcings is quantitatively different. It should be noted that, when one considers F ≤
2, for G ≥ G , the potential energy input is always positive - yet small. We highlight some qualitative features of the dynamics of the model that indicate the presence of wave-like structures amidstchaos in the regime of moderate forcing. Figure 4 shows some examples of evolution of the system of the system in the case of N =
36 sectors with F = G = F = G =
10 (panel b); and F = G =
10 (panel c). We are using a Hovmöller-type diagramme [85], where time is on the vertical axis and the variables X k and θ k , k = , . . . ,
36 are on the horizontal axis. Thisdiagramme is particularly well suited for appreciating wave-like structures, as it is easy to to visualise wave crests. If the forcingon the θ variables is switched off, the X variables behave similarly to the case of the L96 model, where, amidst chaos, the clearsignature of a westward propagating phase velocity can be found, as already observed in [39] and recently mentioned in [73].Ascan be guessed from the evolution equations, the θ variables feature weaker variability and similar pattern of the wave crests, asthey are advected by the X variables and receive energy from them. The situations is qualitatively similar when both the X and θ variables are forced, but, quite naturally, the fluctuations of the θ variables are stronger than in the previous case. Note that inthe case analysed here of F = G =
10, the wave crests travel in the opposite direction with respect to what we have found for theneutral wave emerging for F = G = √
2, see Sect. 2.3. Therefore, the presence of a turbulent background radically changes the abriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere 9a) b)c) d)e) f)
Fig. 3. a) F N ¯ Ξ = F N ∑ Nk = X k ; b) − γ N ¯ K = − γ N ∑ Nk = X k ; c) G N ¯ Θ = G N ∑ Nk = θ k ; d) − N γ ¯ P = − γ N ∑ Nk = θ k ; e) N ¯ C = − N N ∑ k = ( α X k θ k ) ;f) η = ¯ C / ( γ ¯ E ) as a function of F and G . In all simulations γ = α = kinematics of the waves. If, instead, the forcing acts on the θ variables only, the wave crests have a much less clear direction ofpropagation, both for the θ and for the X variables, where the latter feature a much lower variability, as expected. In other terms,the setup where F vanishes is characterised by absolute instability, with little or no advection of anomalies, whereas the othertwo cases are characterised by convective instability, where anomalies are spatially advected [86].We will further discuss in the following sections in more quantitative terms the differences emerging when forcing the X variables only, the θ variables only, or all variables. Fig. 4.
Slice of 10 time units of the evolution of the system with N =
36 for F = G = F = G =
10 (panel b); and F = G =
10 (panel c). Time is on the vertical axis. On the horizontal axis, the first 36 variables are the X k , k = ,...,
36, followed by the variables θ k , k = ,...,
36. Panel d): Extensivity of the system - spectrum of Lyapunov exponents as a function of the rescaled index x = ( j − / ) / ( N ) for F = G = F = G =
10 (black); F = G =
10 (blue). Solid lines: N =
72; Dots: N =
36; Crosses: N =
18. In all simulations, α = γ = Ruelle [87] proposed that systems with short-range interaction can feature extensive chaos, because large domains can be hi-erarchically partitioned into smaller, weakly interacting subdomains with similar properties. One way to test whether chaosextensivity is to analyse the finite-size scaling of the Lyapunov exponents. Specifically, one plots the obtained spectrum of Lya-punov exponents for different values of system size Q (in our case, Q = N ) against the rescaled index x = ( j − / ) / Q and testswhether a universal curve is obtained in the limit of large values of Q [88,65]. We remark that chaos extensivity implies that theratio between the Kaplan-Yorke dimension of the attractor [89], also referred to as Lyapunov dimension [67], and Q tends to aconstant as Q → ∞ .In order to prove convincingly the extensive nature of chaos in the system analyzed here, one should test such property for allvalues of F and G . This is beyond the scope of this paper. Yet, preliminary results do confirm extensivity for the three referencecases F = G = F = G = F = G =
10 shown in Panels a)-c) of Fig. 4. Indeed, we have here performed simulationswith N = N =
36, and N =
72 and, as shown in 4d), the Lyapunov exponents spectra seem to collapse to universal curvesas N grows. Indeed, the Lyapunov spectra plotted against their respective rescaled indices can hardly be visually distinguished.This is especially encouraging in view of the clear evidence for chaos extensivity in the L96 model [60,65]. As thoroughly analyzed in [65], in the one-layer L96 model the average energy per unit site scales to a very high degree ofapproximation as F . for large values of F . The origin of such a scaling law is still unknown. We report some preliminaryresults of scaling laws obtained for the current model in some special configurations of parameters. We have performed longintegrations (1000 time units) at steady state for three set of experiments: abriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere 11 F = j , j = , . . . , G = F = G = j , j = , . . . , G = j , j = , . . . , F = X variables only, on both the X and the θ variables, or onthe θ variables only, respectively. These are regimes of forcing where, see the case of the L96 model [65], one might expect thatchaos extensivity applies with a very good approximation, see Sect. 3.4. We obtain the following approximate asymptotic scalinglaws, which are rather accurate when F and/or G are larger than 256:1. ¯ K , ¯ E ∝ F . , ¯ Ξ ∝ F . , ¯ Θ ∝ F − . , ¯ P = − ¯ C / ∝ F . , λ ∝ F . ;2. ¯ K , ¯ E , ¯ P ∝ F . , ¯ P ≈ . K , ¯ C < , | ¯ C | ∝ F . , ¯ Ξ , ¯ Θ ∝ F . , ¯ Θ ≈ . Ξ , λ ∝ F . ;3. ¯ P , ¯ E ∝ G . , ¯ K = ¯ C / ∝ G . , ¯ Θ ∝ G . , ¯ Ξ ∝ G . , λ ∝ G . ;where the uncertainty is 0.01 for all the numbers above. As clear from these scaling relations, and in agreement with what onecould intuitively guess by looking at Fig. 4, it is rather different to force the system through the X or the θ variables, and theinterplay between the two reservoirs of energy is far from trivial. If the forcing is applied to only one set of variables, the energycycle is more enhanced, ceteris paribus , when the θ variables undergo the forcing. Indeed, the reservoir of total energy andthe conversion ¯ C of energy between the two forms are larger than for corresponding case of forcing applied uniquely to the X variables. The behaviour of the quantities ¯ Ξ and ¯ Θ is also extremely different in the two cases, implying a qualitatively differentway the forcing impacts the spatially-coherent fluctuations of the variables. If F = G , the ratio of the average size of the tworeservoirs of energy is a constant, with ¯ K being larger that ¯ P (and, correspondingly, ¯ Ξ being larger than ¯ Θ ), and the average fluxof energy goes from kinetic to potential. We remind that the dynamics of the F = G case is characterised by convective instability,similarly to the case where G =
0, compare 3a) and 3b). The reason why the forcing through the kinetic channel dominates inthe special case of F = G is still unclear and should be further investigated.In all cases, the amount of energy that is converted between the two forms becomes a negligible fraction of the total incomingenergy in the limit of large forcing. In other terms, the efficiency of the model η = ¯ C / ( γ ¯ E ) tends to zero if either F or G tendto infinity, even if ¯ C tends to infinity. It is then unsurprising that when considering the limit of large F , regardless of whether G is also increased, we obtain for the X variables results that are in agreement with what featured by the one-layer L96, comparewith [65]. At this regard, a useful piece of information is obtained by looking at the properties of λ in the large forcing limit.One obtains that in scenarios 1 and 2, λ ∝ F . , which is again in excellent agreement with what obtained for the L96 model(including the pre-exponential factor). Scenarios 1) and 2) seem like featuring a rather similar dynamics, the main differencebetween the two being the strength of the fluctuations of the θ variables; compare with panels a) and b) of Fig. 4.The growth of λ with G is slower in the case F is set to zero, where we have absolute instability. We can gain a qualitativeunderstanding of the different impact on λ of changes in the value of G vs F by comparing panels a), b) and c) of Fig. 4, whichnonetheless describe weaker regimes of forcing (what follows stands also in the case of stronger applied forcing). Simple and conceptual models have proved extremely useful for better understanding the dynamics of climate as a whole as wellas of its individual components. Indeed, their usefulness spans from being the testbeds for developing new methods in terms ofdata analysis, data assimulation, and model testing; to supporting the definition of new metrics for testing more complex models;to providing valuable insights in the basic active physical mechanisms and most prominent mathematical features.The model presented in this paper goes in this direction and has been constructed in order to provide a new layer of physicalcomplexity to the L96 model by adding a new variable to each gridpoint of the model. This variable can be loosely interpreted asa local temperature and allows for the establishment of a complex energetics for the system, encompassing energy input, output,and conversion. Two forms of energy are present in the system, a kinetic one and potential one. We are also able to introducea notion of efficiency, which is useful for studying the conversion of energy from one form to the other one. The energetics ofthe model is reminiscent of the one of the real atmosphere. Extending previous analyses, we have provided a fairly completeanalysis of the mechanics of the new model by separating a quasi-symplectic and a metric component to its dynamical structure.The energy of the system is used to construct the antisymmetric evolution operator, whose corresponding brackets are not truePoisson brackets because they do not obey the Jacobi identity, hence the symplectic structure is not complete.We have then performed a preliminary analysis of some of the key aspects of the new model by investigating how its propertieschange as a function of the two parameters that control the input of kinetic and potential energy.We have studied, in a special case,the Hopf bifurcation leading to the onset of the neutral wave from the fixed point solution. The interplay between the two forcingsis extremely non-trivial in the weak forcing regime, where much needs to be explored regarding the transition from fixed pointto quasi-periodic to chaotic asymptotic states, and one expects that the structure and position of the bifurcations might dependdelicately on the number of modes included in the system, similarly to the case of the L96 model. When considering regimesassociated with stronger forcing, the system exhibits extensive chaos, even if there is clear evidence of wave-like structuresemerging in the context of an overall strongly chaotic flow. Understanding the interplay between ordered wave-like structuresand turbulence seems of great interest.
The system reacts differently depending on how we force it. The nature of the flow is impacted because absolute vs convectiveinstability dominate if we force the system through the potential energy vs kinetic energy channel, respectively. The mechanismof energy conversion makes sure that also the variables that are not directly forced feature nontrivial variability. If the strength ofthe forcing is the same in the two channels, the kinetic energy channel ends up being more efficient: the dynamics is characterisedby convective instability, and, on the average, energy is transferred from the kinetic to the potential form. The reason for thisbehaviour is still unclear. Similarly to the case of the L96 model, it is possible to obtain accurate power laws describing howsome of the fundamental dynamical and thermodynamical properties of the system scale with the forcing parameters in the limitof very strong forcings.The analysis presented here is only a first step in the direction of better understanding the properties of this model, whichwe believe has the potential of being of great interest for investigations in areas like statistical physics, nonlinear dynamics, dataassimilation, mechanics, model reduction techniques, and extreme events.Finally, again along the lines of the L96 model, we have introduced a two-level version - see Appendix A - of the model,which allows for studying multi-scale dynamics and which features an energetics that resembles, conceptually, the one of theatmosphere, where the Lorenz energy cycle describes succinctly the input and output of energy in the kinetic and potential formas well as the conversion between the two forms and between energy compartments at small vs large scales. The study of theproperties of this model, which is a fortiori extremely promising in the fields above, will be carried out in a future work.
A Two-level new model
Here we propose an extension of the model able to represent multiscale dynamics and energy exchanges across scales. A de-tailed investigation of the properties of this model will be discussed elsewhere. We present the evolution equation and discussthe energetics of the model. Mimicking the structure of the classical two-layer L96 model, we define the following evolutionequations: dX k dt = X k − ( X k + − X k − ) − αθ k − γ X k + F − hcb J ∑ j = Y j , k , (22) dY j , k dt = cbY j + , k ( Y j − , k − Y j + , k ) − α c φ j , k − γ cY j , k + f + hcb X k , (23) d θ k dt = X k + θ k + − X k − θ k − + α X k − γθ k + G − hcb J ∑ j = φ j , k , (24) d φ j , k dt = cb ( Y j − , k φ j − , k − Y j + , k φ j + , k ) + α cY j , k − γ c φ j , k + g + hcb θ k , (25)with k = , ..., N ; j = , ..., J , where Y j , k ’s are j small-scale variables coupled with X k , similarly to the two-level L96 model, and φ j , k ’s are j small-scale variables similarly coupled with θ k . Additionally, the variables Y j , k and φ j , k are also mutually coupled.The constant hc / b determines the strength of the coupling between variables at different scale, while c defines the time scaleseparation between the two levels and b controls the relative amplitude of the fluctuations between large and small scales. Finally, f ( g ) defines the forcing on the variables Y j , k ( φ j , k ). The boundary conditions are the following: X k − N = X k + N = X k , Y j , k − N = Y j , k + N = Y j , k , Y j − J , k = Y j , k − , Y j + J , k = Y j , k + . θ k − N = θ k + N = θ k , φ j , k − N = φ j , k + N = φ j , k , φ j − J , k = φ j , k − , φ j + J , k = φ j , k + . (26)Note that, choosing α = θ and φ variables we obtain the classic two-layer L96 model.Below we show the equations for the fluxes of energy for X and θ alongside with those for Y and φ , as we have done inSection 2 for the one-level model.We define K L = N ∑ k = X k / K S = J ∑ j = N ∑ k = Y j , k / P L = N ∑ k = θ k /
2, and P S = J ∑ j = N ∑ k = φ j , k /
2, as the potential energy at large and small scales, respectively. One derives thefollowing equations for the various reservoirs of energy: ddt K L = I K L − C K L , K S + C P L , K L − D K L (27) ddt K S = I K S + C K L , K S + C P S , K S − D K S (28) ddt P L = I P L − C P L , P S − C P L , K L − D P L (29) ddt P S = I P S + C P L , P S − C P S , K S − D P S (30) abriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere 13 K L P L P S K S I P L D P L C P S , K S C P L , K L C P L , P S C K L , K S I K L I K S I P S D P S D K L D K S Fig. 5.
Diagram of the Lorenz energy cycle of the two-level model. We indicate the fluxes of energy between large (L) and small (S) scalesand between kinetic (K) and potential (P) energy. Dashed lines represent input of energy; dotted lines represent energy dissipation; and solidlines represent energy conversion terms. where I K L = F N ∑ k = X k I K S = f J ∑ j = N ∑ k = Y j , k D K L = γ K L D K S = γ cK S I P L = G N ∑ k = θ k I P S = g J ∑ j = N ∑ k = φ j , k D P L = γ P L D P S = γ cP S C K L , K S = hcb N ∑ k = J ∑ j = X k Y j , k C P L , P S = hcb N ∑ k = J ∑ j = θ k φ j , k C P L , K L = − α N ∑ k = X k θ k C P S , K S = − α c N ∑ k = J ∑ j = Y j , k φ j , k (31)The energy cycle of this model is depicted in Fig. 5 and is closely reminiscent of the Lorenz energy cycle of the atmosphere.Note that, if we add Eqs. 27 with 28, on the one hand, and Eqs. 29 with 30, on the other hand, we derive the energy budget forthe total kinetic and total potential energy, respectively. Instead, if we add Eqs. 27 with 29, on the one hand, and Eqs. 28 with 30,on the other hand, we derive the budgets for the total energy at large and small scale, respectively. Acknowledgements
GV wishes to thank C. Franzke for several useful discussions. VL wishes to thank R. Blender and G. Gallavotti for the manyinspiring conversations on the topics covered in this paper. VL acknowledges the support provided by the EU Horizon 2020project TiPES (grant No. 820970) and by the EPSRC project "Applied Nonautonomous Dynamical Systems: Theory, Methodsand Examples" (grant No. EP/T018178/1) The two authors have equally contributed to this paper. Scripts and data used for thepreparation of this paper can be found at: https://figshare.com/articles/journal_contribution/VissioLucarini2020.zip/12917984 . References
1. J.P. Peixoto and A.H. Oort.
Physics of Climate . American Institute of Physics, New York, NY, 1992.2. V. Lucarini, R. Blender, C. Herbert, F. Ragone, and S. Pascale. Mathematical and physical ideas for climate science.
Reviews of Geophysics ,52:809–859, 2014.3. Edward N Lorenz.
The nature and theory of the general circulation of the atmosphere , volume 218. World Meteorological OrganizationGeneva, 1967.4. O Pauluis. Sources and sinks of available potential energy in a moist atmosphere.
J. Atmos. Sci. , 64:2627–2641, 2007.5. V. Lucarini. Thermodynamic efficiency and entropy production in the climate system.
Phys. Rev. E , 80:021118, 2009.6. M. H. P. Ambaum.
Thermal Physics of the Atmosphere . Wiley, New York, 2010.7. V. Lucarini, K. Fraedrich, and F. Lunkeit. Thermodynamics of climate change: generalized sensitivities.
Atmospheric Chemistry andPhysics , 10(20):9729–9737, 2010.8. F. Laliberté, J. Zika, L. Mudryk, P. J. Kushner, J. Kjellsson, and K. Döös. Constrained work output of the moist atmospheric heat enginein a warming climate.
Science , 347(6221):540–543, 2015.9. V. Lembo, F. Lunkeit, and V. Lucarini. Thediato (v1.0) – a new diagnostic tool for water, energy and entropy budgets in climate models.
Geoscientific Model Development , 12(8):3805–3834, 2019.10. M. Ghil. A Mathematical Theory of Climate Sensitivity or, How to Deal With Both Anthropogenic Forcing and Natural Variability? InC.-P. Chang, M. Ghil, M. Latif, and J.M. Wallace, editors,
Climate Change: Multidecadal and Beyond , volume 6, pages 31–52. WorldScientific Publishing Co., Singapore, 2015.11. Michael Ghil and Valerio Lucarini. The physics of climate variability and climate change.
Rev. Mod. Phys. , 92:035002, Jul 2020.12. Stephen H. Schneider and Robert E. Dickinson. Climate modeling.
Reviews of Geophysics , 12(3):447–493, 1974.13. I. M. Held. The gap between simulation and understanding in climate modeling.
Bulletin of the American Meteorological Society ,86:1609–1614, 2005.14. V. Lucarini. Modeling complexity: the case of climate science. In U. Gohde, S. Hartmann, and J.H. Wolf, editors,
Models, Simulations,and the Reduction of Complexity , pages 229–254. De Gruyter, 2013.15. T.N. Palmer and P.D. Williams. Introduction. Stochastic physics and climate modelling.
Philosophical transactions. Series A, Mathemati-cal, physical, and engineering sciences , 366(1875):2421–7, 2008.16. C.L.E. Franzke, T.J. O’Kane, J. Berner, P.D. Williams, and V. Lucarini. Stochastic climate theory and modeling.
Wiley InterdisciplinaryReviews: Climate Change , 6(1):63–78, 2015.17. J. Berner, U. Achatz, L. Batté, L. Bengtsson, A. De La Cámara, H.M. Christensen, M. Colangeli, D.R.B. Coleman, D. Crommelin, S.I.Dolaptchiev, C.L.E. Franzke, P. Friederichs, P. Imkeller, H. Järvinen, S. Juricke, V. Kitsios, F. Lott, V. Lucarini, S. Mahajan, T.N. Palmer,C. Penland, M. Sakradzija, J.-S. Von Storch, A. Weisheimer, M. Weniger, P.D. Williams, and J.-I. Yano. Stochastic Parameterization:Towards a new view of Weather and Climate Models.
Bulletin of the American Meteorological Society , 98:3:565–588, 2017.18. E.N. Lorenz. Deterministic Nonperiodic Flow.
Journal of the Atmospheric Sciences , 20:130–141, 1963.19. B. Saltzman. Finite Amplitude Free Convection as an Initial Value Problem—I.
Journal of the Atmospheric Sciences , 19:329–341, 1962.20. H. Stommel. Thermohaline convection with two stable regimes of flow.
Tellus , 2:244–230, 1961.21. G. Veronis. An analysis of the wind-driven ocean circulation with a limited number of Fourier components.
Journal of the AtmosphericSciences , 20:577–593, 1963.22. C. Rooth. Hydrology and ocean circulation.
Progress in Oceanography , 11:131–149, 1982.23. J. G. Charney and J. G. DeVore. Multiple flow equilibria in the atmosphere and blocking.
J. Atmos. Sci. , 36:1205–1216, 1979.24. E.N. Lorenz. Irregularity: a Fundamental Property of the Atmosphere.
Tellus A: Dynamic Meteorology and Oceanography , 36:98–110,1984.25. Edward N. Lorenz. Predictability - a problem partly solved. In Tim Palmer and Renate Hagedorn, editors,
Predictability of Weather andClimate , pages 40–58. Cambridge University Press, 1996.26. M.I. Budyko. The effect of solar radiation variations on the climate of the earth.
Tellus , 21:611–619, 1969.27. W. D. Sellers. A global climatic model based on the energy balance of the earth atmosphere.
J. Appl. Meteorol. , 8:392–400, 1969.28. M. Ghil. Climate stability for a Sellers-type model.
J. Atmos. Sci. , 33:3–20, 1976.29. K. Fraedrich. Catastrophes and resilience of a zero-dimensional climate system with ice-albedo and greenhouse feedback.
QuarterlyJournal of the Royal Meteorological Society , 105(443):147–167, 1979.30. E. B. Gledzer. System of hydrodynamic type admitting two quadratic integrals of motion.
Dokl. Akad. Nauk SSSR , 209:1046–1048, 1973.31. L. Biferale, A. Lambert, R. Lima, and G. Paladin. Transition to Chaos in a Shell Model of Turbulence.
Physica D: Nonlinear Phenomena ,80:105–119, 1995.32. V.S. L’vov, E. Podivilov, A. Pomyalov, I. Procaccia, and D. Vandembroucq. Improved shell model of turbulence.
Physical Review E ,58:1811, 1998.33. L. Biferale. Shell Models of Energy Cascade in Turbulence.
Annual Review of Fluid Mechanics , 35:441–468, 2003.34. A. Brandenburg. Energy Spectra in a Model for Convective Turbulence.
Physical Review Letters , 69:605–608, 1992.35. J. Mingshun and L. Shida. Scaling behavior of velocity and temperature in a shell model for thermal convective turbulence.
PhysicalReview E , 56:441–446, 1997.36. K. Hasselmann. Stochastic climate models, part I. theory.
Tellus , 28(6):473–485, December 1976.37. B. Saltzman.
Dynamical Paleoclimatology: Generalized Theory of Global Climate Change . Academic Press New York, New York,November 2001.38. Imkeller, P. and von Storch, J.S.
Stochastic Climate Models . Birkhauser, Basel, 2001.39. Edward N. Lorenz. Designing Chaotic Models.
Journal of the Atmospheric Sciences , 62(5):1574–1587, 05 2005.abriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere 1540. D. Orrell. Model Error and Predictability over Different Timescales in the Lorenz ’96 Systems.
Journal of the Atmospheric Sciences ,60(17):2219–2228, 2003.41. D.S. Wilks. Effects of stochastic parametrizations in the Lorenz ’96 system.
Quarterly Journal of the Royal Meteorological Society ,131(606):389–407, 2005.42. Rafail Abramov. A simple stochastic parameterization for reduced models of multiscale dynamics.
Fluids , 1(1), 2016.43. H. M. Arnold, I. M. Moroz, and T. N. Palmer. Stochastic parametrizations and model uncertainty in the lorenz system.
PhilosophicalTransactions of the Royal Society A: Mathematical, Physical and Engineering Sciences , 371(1991):20110479, 2013.44. G. Vissio and V. Lucarini. A proof of concept for scale-adaptive parametrizations: the case of the Lorenz ’96 model.
Quarterly Journal ofthe Royal Meteorological Society , 144:63–75, 2018.45. G. Vissio. Statistical mechanical methods for parametrization in geophysical fluid dynamics.
Reports on Earth System Science , 212, 2018.46. A. Chattopadhyay, P. Hassanzadeh, and D. Subramanian. Data-driven predictions of a multiscale lorenz 96 chaotic system using machine-learning methods: reservoir computing, artificial neural network, and long short-term memory network.
Nonlinear Processes in Geo-physics , 27(3):373–389, 2020.47. R. Blender and V. Lucarini. Nambu representation of an extended lorenz model with viscous heating.
Physica D: Nonlinear Phenomena ,243(1):86 – 91, 2013.48. Tamás Bódai.
Extreme Value Analysis in Dynamical Systems: Two Case Studies , page 392–429. Cambridge University Press, 2017.49. A. E. Sterk and D. L. van Kekem. Predictability of extreme waves in the lorenz-96 model near intermittency and quasi-periodicity.
Complexity , 2017:9419024, 2017.50. Guannan Hu, Tamás Bódai, and Valerio Lucarini. Effects of stochastic parametrization on extreme value statistics.
Chaos: An Interdisci-plinary Journal of Nonlinear Science , 29(8):083102, 2019.51. A. Trevisan and F. Uboldi. Assimilation of Standard and Targeted Observations within the Unstable Subspace of the Observation – Analysis– Forecast Cycle System.
Journal of the Atmospheric Sciences , 61:103–113, 2004.52. A. Trevisan, D. Isidoro, and O. Talagrand. Four-dimensional variational assimilation in the unstable subspace and the optimal subspacedimension.
Quarterly Journal of the Royal Meteorological Society , 136:487–496, 2010.53. Guannan Hu and Christian Franzke. Data assimilation in a multi-scale model.
Mathematics of Climate and Weather Forecasting ,3:118–139, 12 2017.54. Julien Brajard, Alberto Carrassi, Marc Bocquet, and Laurent Bertino. Combining data assimilation and machine learning to emulatea dynamical model from sparse and noisy observations: A case study with the lorenz 96 model.
Journal of Computational Science ,44:101171, 2020.55. D. S. Wilks. Comparison of ensemble-mos methods in the lorenz ’96 setting.
Meteorological Applications , 13(3):243–256, 2006.56. Wansuo Duan and Zhenhua Huo. An Approach to Generating Mutually Independent Initial Perturbations for Ensemble Forecasts: Orthog-onal Conditional Nonlinear Optimal Perturbations.
Journal of the Atmospheric Sciences , 73(3):997–1014, 02 2016.57. N. Le Carrer and P. L. Green. A possibilistic interpretation of ensemble forecasts: experiments on the imperfect lorenz 96 system.
Advancesin Science and Research , 17:39–45, 2020.58. Diego Pazó, Ivan G. Szendro, Juan M. López, and Miguel A. Rodríguez. Structure of characteristic lyapunov vectors in spatiotemporalchaos.
Phys. Rev. E , 78:016209, Jul 2008.59. S. Hallerberg, D. Pazó, J.M. López, and M.A. Rodríguez. Logarithmic bred vectors in spatiotemporal chaos: Structure and growth.
PhysicalReview E - Statistical, Nonlinear, and Soft Matter Physics , 81(6):1–8, 2010.60. A. Karimi and M. R. Paul. Extensive chaos in the Lorenz-96 model.
Chaos: An Interdisciplinary Journal of Nonlinear Science , 20:043105,2010.61. M. Carlu, F. Ginelli, V. Lucarini, and A. Politi. Lyapunov analysis of multiscale dynamics: the slow bundle of the two-scale lorenz 96model.
Nonlinear Processes in Geophysics , 26(2):73–89, 2019.62. Rafail V. Abramov and Andrew J. Majda. New approximations and tests of linear fluctuation-response for chaotic nonlinear forced-dissipative dynamical systems.
Journal of Nonlinear Science , 18(3):303–341, 2008.63. V. Lucarini and S. Sarno. A statistical mechanical approach for the computation of the climatic response to general forcings.
NonlinearProcesses in Geophysics , 18(1):7–28, 2011.64. Valerio Lucarini. Stochastic perturbations to dynamical systems: A response theory approach.
Journal of Statistical Physics , 146(4):774–786, 2012.65. G. Gallavotti and V. Lucarini. Equivalence of Non-equilibrium Ensembles and Representation of Friction in Turbulent Flows : The Lorenz96 Model.
Journal of Statistical Physics , 156:1027–1065, 2014.66. Rafail V. Abramov. Leading order response of statistical averages of a dynamical system to small stochastic perturbations.
Journal ofStatistical Physics , 166(6):1483–1508, 2017.67. E. Ott.
Chaos in Dynamical Systems . Cambridge University Press, Cambridge, England, 1993.68. H. Broer, C. Simó, and R. Vitolo. Bifurcations and strange attractors in the Lorenz-84 climate model with seasonal forcing.
Nonlinearity ,15:1205–1267, 2002.69. D. Orrell and L.A. Smith. Visualising bifurcations in high dimensional systems: The spectral bifurcation diagram.
International Journalof Bifurcation and Chaos , 13(10):3015–3027, 2003.70. Dirk L. van Kekem and Alef E. Sterk. Travelling waves and their bifurcations in the lorenz-96 model.
Physica D: Nonlinear Phenomena ,367:38 – 60, 2018.71. D. L. van Kekem and A. E. Sterk. Wave propagation in the lorenz-96 model.
Nonlinear Processes in Geophysics , 25(2):301–314, 2018.72. Dirk L. van Kekem and Alef E. Sterk. Symmetries in the lorenz-96 model.
International Journal of Bifurcation and Chaos ,29(01):1950008, 2019.73. John Kerin and Hans Engler. On the lorenz ’96 model and some generalizations, 2020.6 Gabriele Vissio, Valerio Lucarini: Mechanics and Thermodynamics of a New Minimal Model of the Atmosphere74. James R. Holton.
An introduction to dynamic meteorology . International Geophysics Series. Elsevier Academic Press„ Burlington, MA, 4edition, 2004.75. Zhi-Min Chen and W.G. Price. On the relation between rayleigh-bénard convection and lorenz system.
Chaos, Solitons and Fractals ,28(2):571 – 578, 2006.76. Valerio Lucarini and Klaus Fraedrich. Symmetry breaking, mixing, instability, and low-frequency variability in a minimal lorenz-likesystem.
Phys. Rev. E , 80:026313, Aug 2009.77. J. Eckmann and D. Ruelle. Ergodic theory of chaos and strange attractors.
Reviews of Modern Physics , 57(4):617–656, 1985.78. R. Grotjahn.
Global Atmospheric Circulations: Observations and Theories . Oxford University Press, Oxford, 1993.79. Allan N. Kaufman. Dissipative hamiltonian systems: A unifying principle.
Physics Letters A , 100(8):419 – 422, 1984.80. Miroslav Grmela. Bracket formulation of dissipative fluid mechanics equations.
Physics Letters A , 102(8):355 – 358, 1984.81. G. Ansmann. Efficiently and easily integrating differential equations with JiTCODE, JiTCDDE, and JiTCSDE.
Chaos: An InterdisciplinaryJournal of Nonlinear Science , 28(4):043116, 2018.82. Lawrence F. Shampine and Mark W. Reichelt. The matlab ode suite.
SIAM Journal on Scientific Computing , 18(1):1–22, 1997.83. S. Strogatz.
Nonlinear Dynamics and Chaos: With Applications to Physics, Biology, Chemistry, and Engineering . Westview Press, 2 ed.,Boulder, United States of America, 2014.84. G. Benettin, L. Galgani, A. Giorgilli, and J.M. Strelcyn. Lyapunov Characteristic Exponents for smooth dynamical systems and forHamiltonian systems; A method for computing all of them.
Meccanica , 15:21–30, 1980.85. Ernest Hovmöller. The trough-and-ridge diagram.
Tellus , 1(2):62–66, 1949.86. P Huerre and P A Monkewitz. Local and global instabilities in spatially developing flows.
Annual Review of Fluid Mechanics , 22(1):473–537, 1990.87. David Ruelle. Large volume limit of the distribution of characteristic exponents in turbulence.
Communications in Mathematical Physics ,87(2):287–302, 1982.88. Kazumasa A. Takeuchi, Hugues Chaté, Francesco Ginelli, Antonio Politi, and Alessandro Torcini. Extensive and subextensive chaos inglobally coupled dynamical systems.
Phys. Rev. Lett. , 107:124101, Sep 2011.89. J.L. Kaplan and J.A. Yorke. Chaotic behavior of multidimensional difference equations.