Learning Chaotic Dynamics in Dissipative Systems
Authors: Zongyi Li
Miguel Liu-Schiaffini
Nikola Borislavov Kovachki
Kamyar Azizzadenesheli
Burigede Liu
Kaushik Bhattacharya
Andrew Stuart
Anima Anandkumar
Venue: NeurIPS 2022 Accept (2022)
Keywords: Dissipative Chaotic systems
operator learning
invariant statistics
attractor learning

Abstract: Chaotic systems are notoriously challenging to predict because of their sensitivity to perturbations and errors due to time stepping. Despite this unpredictable behavior, for many dissipative systems the statistics of the long term trajectories are governed by an invariant measure supported on a set, known as the global attractor; for many problems this set is finite dimensional, even if the state space is infinite dimensional. For Markovian systems, the statistical properties of long-term trajectories are uniquely determined by the solution operator that maps the evolution of the system over arbitrary positive time increments. In this work, we propose a machine learning framework to learn the underlying solution operator for dissipative chaotic systems, showing that the resulting learned operator accurately captures short-time trajectories and long-time statistical behavior. Using this framework, we are able to predict various statistics of the invariant measure for the turbulent Kolmogorov Flow dynamics with Reynolds numbers up to $5000$.

Learning Dissipative Dynamics in Chaotic Systems
Zongyi Li Caltech zongyili@caltech.edd
Miguel Liu-Schiaffiii Caltech m mliuchh@calltchh edu
Nikola Kovachki
NVIDIA nkovachki @nvidia. com
Burigede Liu University of Cambridge b1377@eng cam. uk
Kamyar Azizzadenesheli
Kaushik Bhattacharya
NVIDIA
Caltech bhatta@caltech. edu
kamyara@nvidia.com com
Andrew Stuart Caltech astuart@caltech.com edu
Anima Anandkumar
Caltech
anima@caltech.edi edu
Abstract
Chaotic systems notoriously challenging to predict because of their sensitivity perturbations and errors due to time stepping. Despite this unpredictable behavior, for many dissipative systems the statisties of the long term trajectories by an invariant measure supported on a set, known as the global attractor; for many problems this set is finite dimensional, even if the state space is infinite dimensional. For Markovian systems, the statistical properties of long-term trajectories are uniquely determined by the solution operator that maps the evolution of the system over arbitrary positive time increments. In this work, we propose a machine learning framework learn the underlying solution operator for dissipative chaotic systems, showing that the resulting learned operator accurately captures short-time trajectories and long long-time time statistical behavior. Using this framework, we are able to predict various statistics of the invariant measure for the turbulent Kolmogorov Flow dynamics with Reynolds numbers up to 5000.
governed
1 Introduction
Machine learning methods for chaotic systems. Chaotic systems are characterized by strong instabilities. Small perturbations to the system, such the initialization, lead to errors which accumulate during time-evolution, due to positive Lyapunov exponents. Such instability makes chaotic systems challenging, both for mathematical analysis and numerical simulation. Because of the intrinsic instability, is infeasible for any method to capture the exact trajectory of a chaotic system for long periods. Therefore, prior works mainly focus fitting recurrent neural networks (RNN) on extremely short trajectories or learning step-wise projection from a randomly generated evolution using reservoir computing (RC) (1.) Another approach this problem is to work in the space of probability densities propagated by, or invariant under, the dynamics [ 77. the Frobenius- Perron approach; a second, related, approach is consider the adjoint of this picture and study the propagation of functionals on the state space of the system, the Koopman operator approach and Koopman
g However, even when the state space finite dimensional, the Frobenius- approaches require approximation of a linear operator on a function space; when the state space deffiite-dimensionaa approximation of this operator is particularly challenging. These previous
"Equal contribution
The code is available ottpat/tttt........occcc.o.oootr
36th Conference on Neural Information Processing Systems (NcurIPS 2022).

---

attempts are able to push the limits of faithful prediction to moderate periods on low dimensional ordinary differential equations (ODEs), e.g. the Lorenz-63 system, or on one-dimensional partial differential equations (PDEs), e.g. the Muramoto-Sivashinsee (KS) equation. However, they are less effective at modeling more complicated turbulent systems such as the Navier-Stokes equation (NS), especially over long time periods. Indeed, because of positive Lyapunov exponents, we cannot expect predictions of long-time trajectories of such chaotic systems to be successful. Instead, we take a new perspective: we capture statistical properties of long trajectories by accurately learning the solution operator, mapping initial condition to solution at a short time later.
Invariants in chaos. Despite their instability, many chaotic systems exhibit certain reproducible statistical properties, such as the auto-correlation and, for PDEs, the energy spectrum. Such properties remain the same for different realizations of the initial condition [0]. This is provably the case for the Lorenz-63 model equation and the two-dimensional Navier-Stokes equation (Kolmogorov flows) [T3.. Dissipativity is physical property of many natural systems. Dissipativity may be encoded mathematically as the existence of a compact set into which all bounded sets of initial conditions evolve in a finite time, and therafter remain inside; the global attractor captures all possible long time behaviour of the system by mapping this compact set forward in time [3 [4]. There is strong empirical evidence that many dissipative systems are ergodic, i.c., there exists invariant measure which is supported on the global attractor. For autonomous systems, our focus here, the solution operator for the dynamical system maps the initial condition to the solution at any positive time. By fixing an arbitrary positive time and composing this map with itself many times, predictions can be made far into the future. This compositional property follows from the Markov semigroup property of the solution operator It is natural to ask that approximations of the solution operator inherit dissipativity of the
[ 12 and empirically holds for many dissipative PDEs, such as the KS
13 14. true solution operator, and that invariant sets and the invariant measure, which characterize the global attractor, are well-approximmted; this problem is well-studied in classical numerical analysis [15] [6] and here we initiate the study of inheriting dissipativity for machine-learner solution operators. By learning solution operator, we are able to quickly and accurately generate an approximate attractor and estimate its invariant measure for variety of chaotic systems that are of interest to the physics and applied mathematics communities 1 23 Prior work in learning complex dynamics from data [24 2
learning the attractors of dynamical systems), but all lack principled data-driven modeling paradigm to enforce it. In this paper we work within a precise mathematical setting encompassing a wide variety of practical problems [13), and we demonstrate two principled approaches to enforce dissipativity. One of these leads to provably dissipative models, something not previously achieved in the literature for our definition of dissipativity. B
has labeled long-term stability a desirable feature (e.g.. 28 29
Neural operators. To learn the solution operators for PDEs, one must model the time-evolution of functions in dimenteoadi function spaces. This is especially challenging when we need to generate long trajectories since even a small error accumulates over multiple compositions of the learned operator, potentially causing an exponential blow-up or a collapse due to the space's s high dimensionality. Because we study the evolution of functions in time, we propose to use neural operators 30 31 a recently developed operator learning method. Neural operators are deep learning models that are maps between function spaces. Neural operators generalize conventional neural networks which are maps between finiteddimensional spaces and subsume neural networks when limited to fixed grids. The input functions to neural operators can be represented in any discretization, and the output functions can be evaluated at any point in the domain. Neural operator remedies the mesh-dependent nature of mnite-dimensionall neural network models such as RNNs. CNNs, and RC. Neural operators are guaranteed to universally approximate any operator in a mesh independent manner BT, and hence, can capture the solution operator of chaotic systems. This approximation guarantee and the absorption of trajectories by the global attractor makes it possible to accurately follow it over long time horizons, allowing access the invariant measure of chaotic systems.
Our contributions. In this work, we formulate a machine learning framework for chaotic systems exploiting their dissipativity and Markovian properties. We propose the Markov neural operator (MNO) and train it given only one-step evolution data from a chaotic system. By composing the
'There are other, more restrictive, definitions of dissipativity; we work with the widely adopted definition from R 4
2

---

Dascline model [||u]|> 10")
Discipative model (levle rv 20)
log distaince from altrrccton meenn
dissipative attractor outer bound
Figure I: Dynamic evolution of the Markov neural operator for Kolmogorov flow systems Starting from initial conditions near the attractor, with and without dissipativity regularization. The baseline model has dissipativity regularization while the dissipative model has the regularization enforced during training. Baseline model blows up, whereas the dissipative model returns to the attractor, The dissipative model trained using the Fourier neural operator architecture in the manner shown in Figure 2
learned operator over long horizon, we accurately approximate the global attractor of the system its performance, we study the statistics of
Our architecture is outlined in Figure In order the associated invariant measure such as the Fourier spectrum, the spectrum of the proper orthogonal decomposition (POD), the point-wise distribution, the auto-correlation, and other domain -specifccc statistics such as the turbulence kinetic energy and the dissipation rate. Furthermore we study the behavior of our leaned operator over long horizons and ensure that it does not blow up or collapse but rather accurately follows the global attractor. In this work:
We theoretically prove that, under suitable conditions. the MNO can approximate the underlying solution operator of chaotic PDEs, while conventional neural networks lack such strong guarantees.
We impose dissipativity by augmenting the data on an outer shell to enforce that the dynamic evolution stays close to the attractor. As an additional safeguard, we impose hard dissipativity constraint on MNO predictions far from the attractor and the data augmentation shell. We show this is crucial for learning in chaotic regime as demonstrated by Figure I The resulting system remains stable against large perturbations.
We study the choice of time steps for training the MNO, demonstrating that the error follows valley-shaped phenomenon (cf.
53 35
This gives rise to recipe for choosing the optimal time
step for accurate learning.
We show that standard mean square error (MSE) type losses for training are not adequate, and the models often fail capture the higher frequency information induced from the derivatives. We investigate various Sobolev losses in operator learning. We show that using Sobolev norms for training captures higher-order derivatives and moments, as well high frequency details [36 37) This is similar in spirit pre-multiplying the spectrum by the wavenumber, approach commonly used in the study of fluid systems 38).
We investigate multiple existing deep learning architectures, including U-Net 39 memory convolution neural networks (LSTM-CNN) 40. and gated recurrent units (GRU) 411, in place of the neural operator learn the solution operator. We show MNO provides order of magnitude lower error on all loss functions studied. Furthermore, we show that MNO outperforms the aforementioned neural network architectures on all statisties studied.
long short-ters term
In summary, we propose a principled approach to the learning dissipative chaotic systems. Dissipativ- ity is enforced through data augmentation and through post-processing of the learned model. The methodology allows for treatment of infinite dimensional dynamical systems (PDEs) and leads to methods with short term trajectory accuracy and desirable long-term statistical properties.
2 Problem setting
We consider potentially infinite dimensional dynamical systems where the phase space is Banach space and. in particular, function space Lipschitz domain D Ed (for finite dimensional

---

Prediction
Global attractor
Markov neural operator
Sobolev loss
Dissipativity regulurization shell
uuft 3
Loss function
Ground Truth
Dissipativity constraint
Dissipativity regularization
off
Markov neural operator
Prediction
Figure 2: Markov neural operator (MNO): learn global dynamics from local data
Learning the MNO from the local time-evolution data with the Sobolev loss and dissipativity regularization. u(t) is t'th time step of the chaotic system. Sobolev norms of various orders are used compute the step-wise loss. Dissipativity regularization is computed by drawing random sample u(t) from the dissipativity shell to make sure that expectation next time step prediction ū(t 1) dissipates norm.
systems, 2/ will be a Euclidean space). We are interested in the initial-value problem
du
di (t) = F(u(t),
n(0) = uo. E 00)
(I)
for initial conditions uo E u where F is usually a non-linear operator. We will assume, given some appropriate boundary conditions on 0D when applicable, the solution u(() E U exists and is unique for all times E (0, oc). When making the spatial dependence explicit, if it is present, we will write u(z, t) to indicate the evaluation l(t) for any E D. We define the family of operators S; : M + as mapping ug + u(f) for any t 2 0, and note that, since ( is autonomous, S, satisfies the Markov property i.c. (s,))
u(s + t) for any s, O. We adopt the viewpoint of casting time-dependent PDEs into function space ODEs I.. as this leads to the semigroup approach to evolutionary PDEs which underlies our learning methodology.
One approach to this problem to learn F itself (e.g., 44)) However, infinite dimensions, F is unbounded operator and learning unbounded operators is considerably more challenging than learning bounded operators; this is proven 46 for linear operators. We avoid these issues by instead learning the solution operator S, for some fixed & O; this map is bounded as an operator on Zt. Furthermore, large time statistical behaviour can be obtained more efficiently by using the solution operator and not resolving dynamics at timescales smaller than time t.
Dissipativity. Systems for which there exists some fixed, bounded, positively positively-invariint invariant set such that for any bounded B C U4 there is some time t" = t- (B) beyond which the dynamics of any trajectory starting B enters and remains in are known as dissipative systems [13] 14]. The set E is known as the absorbing set of the system. For such systems, the global attractor A, defined subsequently, characterized as the w-limit set of F.. particular, for any initial condition uo E U the trajectory u(t) approaches A as f -> xo. In this work, we consider dissipative dynamical systems where there exist some a 2 0 and 3 > O such that d
|u|2= (u, F(w)) < <a aa-Bll B1|ul
dt
for all U. It can be shown that systems which satisfy this inequality dissipative 14 with the absorbing set E, an open ball of radius /a/B + for any 0. There several well-known this paper we consider the
examples of dynamical systems that satisfy the above inequality. fiite-dimennional Lorenz-63 system and the aefnite-dimensione cases of the Kuramoto-Sivashinsky and 2D incompressible Navier-Stokes equations, in the form of Kolmogorov flows 1B
Global Attractors. A compact, invariant set A is called global attractor if, for any bounded set B C 24 and any e > 0 there exists time 1" L" (e( B) such that SBB contained within an 6- neighborhood of A for all t 2 t". Formally it is the w-limit set of any absorbing set; this characterizes

---

it uniquely II |4
under assumption 2 Although unique, the w-limit set of a single initial condition may differ across initial conditions [4] In this paper we concentrate on the ergodic setting where the same w- w-limit set is seen for almost all initial conditions with respect to the invariant measure. Many PDEs arising in physics such as reaction-diffusion equations (dynamics of biochemical systems) or and Navier-Stokes equation (fluid mechanics) are dissipative and possess
the
a global attractor which is often fiite-dimennional [3.. Therefore, numerically characterizing the attractor is important problem in scientific computing with many potential applications.
Data distribution. For many applications, exact form for the possible initial conditions
is not available; is therefore convenient to use a stochastic model to describe the initial states. To that end, let (Ao be probability measure l/ and assume that all possible initial conditions to ( come as samples from μo i.c. ug ~ Po- Then any possible state of the dynamic i after some time
t > O can be thought of as being distributed according to the pushforward measure f4 == - i.e. u(t) ~ ft. Therefore as the dynamic evolves, so does the type of likely functions that result. This further complicates the problem of long time predictions since training data may only be obtained up to finite time horizons hence the model will need the ability to predict not only on data that out-of-sample but also
Ergodic systems. To alleviate of the previously presented challenges, consider ergodic systems. Roughly speaking, system is ergodic if there exists an invariant measure which is unchanged by pushforward under S, and time averages of functionals on state space U converge to averages against μ.. For dissipative systems the invariant measure is supported on the global attractor A and together the pair (A, μt) capture the large time dynamics and their statistics. Proving ergodicity is hard for deterministic systems. Rigorous proofs have been developed for the Lorenz-63 model [22: empirical evidence of ergodicity may be seen much more widely, including for the Kuramoto-Sivasee Sivashinsky and Navier-Stokes equations studied in this paper. For stochastic differential equations ergodicity easier to establish; see for details.
Ergodicity mitigates learning a model that is able to predict out-of-distributiin since both the input and the output of Sr. an approximation to Sa. will approximately be distributed according to Furthermore, we may use S), to learn about μ since sampling it simply corresponds to running the dynamic forward. Indeed, we need only generate data on a finite time horizon in order to learn $n. and, once learned, we may use it to sample yu indefinitely by repeatedly composing Sn, with itself. Having samples of μ then allows us to compute statistics which characterize the long term behavior of the system and therefore the global attractor A. This strategy avoids the issue of accumulating errors in long term trajectory predictions since are only interested in the property that Sy (u())) p.
3 Learning the Markov neural operator in chaotic dynamics
We propose the Markov neural operator (MNO), method for learning the underlying solution operators of autonomous, dissipative, chaotic dynamical systems. In particular, we approximate the operator mapping the solution from the current to the next step S₁ : uut) + h(t+ We approximate the solution operator Sμ, an element of the underlying continuous time semigroup {Se [0, oo)) using a neural operator as in Figure 2 in the finite dimensional case we use feedforward neural network. See Appendix B. for background on the solution operator and semigroup.
Sobolev norm as loss We incorporate Sobolev norms in the training process to better capture invariant statistics of the learned operator. Given the ground truth operator S₁ and the learned operator Sn and f f= = Sa(u) Sn(u), compute the step-wise loss in the Sobolev norm,
Il = Σ (1 I
In particular, we use p = 2, where the Sobolev norm can easily computed in Fourier space:
llfl = Σ (1+n²+ ..................... 2 f(n)
........
(4)
n=-
where f is the Fourier series of f. In practice, use k = 0, 1,2 for the training.
5

---

Long-term predictions. Having access to the map Sh, its semigroup properties allow for approxi- mating long time trajectories of ① by repeatedly composing Sr, with its own output. Therefore, for any n E N, we compute u(nh) as follows,
w(nh) ~ Sh (wo) = (ŠA o .... ( (uo).
(5)
times
The above semigroup formulation can be applied with various choices the backbone model for Sj.- In general, prefer models that be evaluated quickly and have approximation guarantees the per-step error can be controlled. Therefore, we choose the standard feed-forward neural network [48] for ODE systems, and the Fourier neural operator (49] for infinite dimensional PDE systems.
For the the neural operator parametric class, we prove the following theorem regarding the MNO. The result states that our construction can approximate trajectories of dinite-dimeesional dynamical systems arbitrary well. The proof is given in Appendix B.2
Theorem 1. Let K C be a compact set and assume that, for some h > O, the solution operator SA U -u associated to the dynamic ① is locally Lipschitz. Then, for any N and >> there exists a neural operator Sn U + u such that
[|u(kh) - (uo) <.
sup sup NOEK kE(1......
Theorem [ indicates that if the backbone model is rich enough, it can approximate many chaotic systems for arbitrarily long periods. For minee-dimensional systems, the same theorem holds with feed- feed-forward neural networks instead neural operators. We note that standard neural networks such RNNs and CNNs do not possess such approximation theorems in the infnite-dimensionall setting.
Invariant statistics. A useful application of the solution operators is to estimate statistics of the invariant measure of chaotic system. Assume the target system is ergodic and there exists invariant measure j as discussed in Section 2
An invariant statistic is defined as
To = -[[ G(u) du(u)= lim
T- T G(u(t)) dt f"
for any functional G U + Rd. Examples include the L norm, any spectral coefficients, and the spatial correlation, as well as problem-specinee statistics such as the turbulence kinetic energy and dissipation rates in fluid flow problems. Given the property (5) and using the ergodicity from 16 the approximate model S, be used estimate any invariant statistic simply by computing E G(S) (u))) where nh and > 0 large enough.
Encouraging dissipativity via regularization. In practice, training data for dissipative systems is typically drawn from trajectories close to the global attractor, so a priori there is no guarantee of learned model's behavior far from the attractor. Thus, if we seek to learn the global attractor and invariant statistics of a dynamical system, it is crucial that we encourage learning dissipativity. To learn S₁ u(t) + u(t + ft) we propose the following loss function, found by supplementing the Sobolev loss with dissipativity-inducii regularization term:
/ [1()(u)) - Xulli dv(u).
(w) S (w)li dp(())p(n))))) +
step-wise Subuley loss
dissiputivity regalarization
Here 4 is the underlying distribution of the training data, is a loss weighting hyperparameter, and O < A < 1 is some constant factor for scaling down (i.e., enforcing dissipativity) inputs drawn from probability measure v. We choose v to be a uniform probability distribution supported on some shell with fixed inner and outer radii from the origin in Lf; these are chosen so that the support of v disjoint from the attractor. Our dissipative regularization term scales down by constant but in principle alternative dissipative cost functionals can be used.

---

(a) No dissipativity
(b) Diss. regularization
(c) Diss. post post-processing
Figure 3: Dissipativity regularization on the Lorenz 63 system flow maps.
Red points are training data on the attractor. The dissipativity regularization in (b) is imposed by augmenting the data the blue shell. Post-processing with 100 (the orange circle) and 0.1 in eq.
Enforcing dissipativity via post-processing. While the previous dissipativity regularization forces a dissipative map in practice, it does not readily yield to theoretical guarantees of dissipativity. For additional safeguard against model instability, enforce hard dissipativity constraint far from the attractor and from the shell where is supported, resulting in provably dissipative dynamics.
In detail, post-process the model: whenever the dynamic moves out of an priori defined stable region, we switch to the second model that pushes the dynamic back. The new model combines the learned model S,, and the safety model v, via a threshold function p
(u, ((],= (1 p(Il
where V/ is some dissipative map and is a partition of unity. For simplicity we define
V (u) Au
p(]|u]) + ((ll--a)
(9)
where a is the effective transition radius between S, and v and 3 controls the transition rate. Note that this choice of V/ is consistent with the regularization term in the loss
By its construction, the post-processed model dissipative. However, its dynamics (in particular the unlearnable transition between Šn, and I/) are not smooth as the dissipative regularization model, which learn the whole dynamic together, as shown in Figure Combining regularization and post-processing results dissipative model with both smooth dynamics and theoretical
quarrainteggaaaeerr guarantees.
Experiments
We evaluate our approach on the finte-diimensional, chaotic Lorenz-63 system well as the chaotic 1D Buramoto-Sivashinst and 2D favier-Stokes equations. all cases we show that encouraging dissipativity is crucial for capturing the global attractor and evaluating statistics of the invariant measure. To the best our knowledge, we showcase the first machine learning method able to predict the statistical behavior of the incompressible Navier-Stokes equation in strongly chaotic regime.
4.1 Lorenz-63 system
To motivate and justify our framework for learning chaotic systems in the anfintte-dimensiona setting, we first apply our framework on the simple yet still highly chaotic Lorenz-63 ODE system, widely-studied (50] simplified model for atmospheric dynamics given by
0g = a 1 Uz). úy = - My Ug U2.
i = U2 - bu. b(r + a).
(10)
We use the canonical parameters (a. b.) = (10, 8/3, 28)
51I Since the solution operator of
Lorenz-63 system is finite-dimensional, we learn it by training a feedforward neural network on single trajectory with h = 0.05s the Lorenz attractor. Figure shows that dissipativity regularization produces predictions that isotropically point towards the attractor. Observe that the our network is also dissipative outside the regularization shell.
7

---

LSTM
GRU
2000
toom
comeettio
ene
2.01
1000
109
-100
62
10
20
wazerumber
40
60
100
500
1000
1500
2000
(a)
(b)
(c)
Figure 4: Choice of time step (a), Fourier spectrum (b), and spatial correlation (c) for KS equations. Both learning solution operator and learning identity residual induce smaller per step error for smaller h. When learned models are composed to generate longer trajectories, we observe that learning residual is advantageous. (b) Fourier spectrum of predicted attractor. All models capture Fourier modes with magnitude larger than (1), while MNO is more accurate on the tail. (c) Spatial correlation of the attractor, averaged in time. MNO more accurate on near-range correlation, but all models miss long-range correlation.
We empirically find that dissipativity can be encouraged without significantly affecting the model's step-wing wise relative L error and the learned statistical properties of the attractor. We also find that varying the dissipativity hyperparameters (loss weight a, scaling factor X, and shell radius) does not significantly affect the step-wise or dissipative error. Details in Appendix A.
4.2 Kuramoto-Sivashinkky equation
We consider the following one- iheediiensional KS equation,
du du
0tu
on [0, L] x (0, oo), on [0, L].
ōt dr d==
Fx
u(:, 0) uo,
where the spatial domain [0, L] is equipped with periodic boundary conditions. We study the impact the time step / has on learning. Our study shows that when the time steps are too large, the correlation chaotic and hard to capture. But counter-intuitively, when the time steps too small, the evolution also hard to capture. In this case, the model's 's input and output are very close, and the identity map will be local minimum. We thus use the MNO also learn the time-derivative or residual. Figure 4a shows the results for varying h and when MNO is used to learn either the identity residual or the solution operator. We observe that the residual model has better per-step error and accumulated error at smaller h. When the time step is large, there is no difference in modeling the residual.
As shown in Figures 4b and 4C
compare the performance of the MNO model against LSTM and GRU that we use to model the evolution operator of the KS equation with h = 1s. We observe that MNO model accurately recovers the Fourier spectrum of KS equation. For all other statistics (see Appendix A.2) all models perform similarly, except on the velocity distribution where MNO performs best. We emphasize that of these statistics are very challenging capture and most machine learning approaches in the literature thus far fail to do so (see Appendix A.2 for details).
4.3 Kolmogorov Flow
We consider wvo-dimensional Kolmogorov flow form of the Navier-Stokes equations) for a viscous, incompressible fluid,
Ou Ōt
1 Re
= u Vu Vp
Au si
V u = 0, on [0, 2m x (0, oo)
(11)
with initial condition ul-, 0) up where denotes the velocity, the pressure, and Re > 0 is the Reynolds number. We encourage dissipativity during training with the criterion described eq. 7
with À = 0.5 and being uniform probability distribution supported on shell around the origin.

---

100
love disqpation
10000
167
0.2
8000
0.2
42
2000
-3
108
0.02 21dr
53 100
(b)
(a)
(d)
Figure 5: The learned attractor of Kolmogorov flow (Re 500) and choice of Sobolev loss for KF. (a-c) The 10000 time steps trajectory generated by MNO projected onto the first two principal components. Each point corresponds snapshot the attractor. The vorticity field for two points is shown. (d) Velocity error for models trained on stream function, velocity, and vorticity using Sobolev loss of different orders.
vartiity
workilyy
157
12a0
1071
30"
predicice 6-050 acoling
prdiicoon
arodddes
60 acalng
5-000 acaing
6- scalrg
Figure 6: Simulation of Kolmogorov flow (Re = 5000) with the dissipative MNO model. The first column corresponds to the initial condition and is sampled from Gaussian random field. Around = 10 to f 20, we see the energy injected from the source term sin (4y) 12 and also the energy transfers from higher frequencies to the lower frequencies. The dotted line is the Kolmogorov energy cascade rate.
We test the effect of encouraging dissipativity in the turbulent (and blow-up prone) Re 500 setting, where we observe that a non-dissipative model blows up when composed with itself multiple times the initial condition is perturbed slightly from the attractor (Figure D]. even though the model achieves relatively low 1.2 error. In contrast, we empirically observe that the dissipative MNO does not blow up and its composed predictions returns to the attractor even when the initial condition is perturbed.
Simulation on turbulent Re 5000 case. The proposed dissipative MNO model can stably learn complex chaotic dynamics at high Reynolds numbers. As shown in Figure 6 the model captures the energy spectrum that converges the Kolmogorov energy cascade rate of k 5/3 The learned model is dissipative with the dissipative regularization, without compromising the accuracy as shown in Table
5 Even given a O(10))-scaled perturbation, the dynamic quickly returns to the attractor the system. In contrast, the model without dissipative regularization will stay outside the attractor given even a (110-scelle perturbation.
Accuracy with respect to various norms. We study performance of MNO along with other predictive architectures. U-Net [39 and LSTM-CNN A0, on modeling the vorticity w in the relatively non-tur cnnturbulent Re 40 case The MNO achieves one order of magnitude better accuracy compared to this architectures. As shown in Table 3] we train each model using the balanced 1²(== H H3, and H²2 losses, defined as the sum of the relative 1 loss grouped by each order of derivative. We measure the error with respect to the standard (unbalanced) norms. The MNO with H² loss consistently achieves the smallest error on vorticity for all of the 1.2, H1, and H² norms. However, the L loss achieves the smallest error the turbulence kinetic energy (TKE), while the

---

H' loss achieves the smallest error on the dissipation Vorticity fields for NS equation (Re = 40) are shown Figure II The figure indicates that the predicted trajectories (b) (c) (d) share the same behaviors as in the ground truth (a), indicating that the MNO model is stable.
Estimating the attractor and invariant statisties. We compose MNO 10000 times obtain the global attractor, and we compute the PCA (POD) basis of these 10000 snapshots and project them onto the first two components, As shown Figure Sa we obtain a cycle-shaped attractor. The true attractor has dimension on the order of (1000) 13.. If the attractor is a high-dimensional sphere, then most of its mass concentrates around the equator. Therefore, when projected to low-dimension, the attractor will have ring shape. We see that most points are located the ring, while few points are located in the center. The points in the center have high dissipation, implying they are intermittent states. In Figure 5b add the time axis. While the trajectory jumps around the cycle, observe there is a rough period of 2000s. We perform the same PCA on the training data, showing the same behavior. Furthermore, in Figure 12 we present invariant statistics for the NS equation (Re 40), computed from a large number of samples from our approximation the invariant measure. We find that the MNO consistently outperforms all other models in accurately capturing these statistics.
Derivative orders. Roughly speaking, vorticity is the derivative of velocity; velocity is the derivative of the stream function. Therefore stream function as 2, 1, and () respectively. Combining vorticity, velocity, and stream function, with I.2, H and H² loss, we have in total the order of derivatives ranging from 0 to 4. We observe, in general, it is best practice to keep the order of derivatives in the model at a number slightly higher than that of the target quantity. For example, shown in Figure 5d when querying the velocity (first-order quantity), it is best to use second-order (modeling velocity plus H loss or modeling vorticity plus L loss). This is further illustrated in Table for the Re 500 case. In general, using a higher order of derivatives as the loss will increase the power of the model and capture the invariant statisties more accurately. However, a higher-order of derivative means higher irregularity. It in turn requires a higher resolution for the model to resolve and for computing the discrete Fourier transform. This trade-off again suggests it is best to pick Sobolev norm not too low or too high.
denote the order of derivative of vorticity, velocity, and
Conclusion
In this work, propose machine learning framework that trains from local data and enforces dissipative dynamics. Experiments also show MNO predicts the attractor which shares the same distribution and invariant statistics as the true function space trajectories. The simulations achieved by the MNO have the potential to further our understanding of many physical phenomena and the mathematical models that underlie them.
Furthermore, the MNO shows great potential for modeling partially observed systems or studying systems that exhibit bifurcations, both of which are of great interest for engineering applications. The study of non- non-eroodic dissipative systems and their basins of attraction is also an important direction for future study.
10

---

Acknowledgements
gratefully acknowledges the financial support from the Kortschak Scholars, PIMCO Fellows, and Amazon A14Science Fellows programs. M. Liu-Schiaffini is supported by the Stephen Adelman Memorial Endowment. A. Anandkumar is supported in part by Bren endowed chair. K. Bhattacharya, N. B. Kovachki, B. Liu, A. M. Stuart gratefully acknowledge the financial support of the Army Research Laboratory through the Cooperative Agreement Number W91 12-0022. A. M. Stuart is also grateful to the US Department of Defense for support as Vannevar Bush Faculty Fellow. Research was sponsored by the Army Research Laboratory and was accomplished under Cooperative Agreement Number W91 12-2-12--00022 A part of this work took place when K. Azizzadeneshelii was at Purdue University. The views and conclusions contained this document are those of the authors and should not be interpreted representing the official policies, either expressed implied, of the Army Research Laboratory or the U.S. Government. The U.S. Government is authorized to reproduce and distribute reprints for Government purposes notwithstanding any copyright notation herein.
References
[1] Pantelis R Vlachas, Wonmin Byeon, Zhong Y Wan, Themistoklis P Sapsis, and Petros Koumout- sakos. Data-driven forecasting of high-dimensional chaotic systems with long short-term memory networks. Proceedings of the Royal Society A: Mathematical,, Physical and Engineer- 2018.
ing Sciences,
[2] Jaideep Pathak, Brian Hunt, Michelle Girvan, Zhixin Lu, and Edward Ott. Model-free prediction of large spatiotemporally chaotic systems from data: A reservoir computing approach. Physical review letters, 120(2)x02402,, 2018.
[3] Pantelis R Vlachas, Jaideep Pathak, Brian R Hunt, Themistoklis Sapsis, Michelle Girvan, Edward Ott, and Petros Koumoutsakos. Backpropagation algorithms and reservoir computing recurrent neural networks for the forecasting of complex spatiotemporal dynamics. Neural Networks, 2020.
[4] Huawei Fan, Junjie Jiang, Chun Zhang, Xingang Wang. and prediction of chaotic systems with machine learning. Physical Review Research, 2(1):012080, 2020.
Lai. Long-term
[5] Jaideep Pathak. Zhixin Lu, Brian R Hunt, Michelle Girvan, and Edward Ott. Using machine learning to replicate chaotic attractors and calculate лyaрunov exponents from data. Chaos: An Interdisciplinary Journal of Nonlinear Science, 12100, 2017.
[6] Michael Dellnitz and Oliver Junge. On the approximation of complicated dynamical behavior. SIAM Journal on Numerical Analysis, 36(2):491-515 1999.
Peter Deufthard, Michael Dellnitz, Oliver Junge, and Christof Schütte. Computation of essen- tial molecular dynamics by subdivision techniques. In Computational molecular dynamics: challenges, methods, ideas, pages 98-115, Springer, 1999.
[8] Romeo Alexander and Dimitrios Giannakis. Operator- framework for forecasting non- linear time series with kernel analog techniques. Physica D: Nonlinear Phenomena, 409:132520, 2020.
[9] Eurika Kaiser, J Nathan Kutz, and Steven L Brunton. Data-driven discovery of koopman eigenfunctions for control. Machine Learning: Science and Technology, 2(3):335033, 2021.
[10] James Gleick. Chaos: Making science. Open Road Media, 2011.
[11] Warwick Tucker. The lorenz attractor exists. Comptes Rendus de l'Académie des Sciences Series I Mathematics, 1977122 1999.
[12] Mark Holland and Ian Melbourne. Central limit theorems and invariance principles for Lorenz attractors. Journal of the London Mathematical Society, 76(2):345-364 10 2007.

---

[13] Roger Temam. Infinite-dinn dinniit-diiesioonall dynamical systems in mechanics and physics, volume 68. Springer Science Business Media, 2012.
[14] A. Stuart and A.R. Humphries. Dynamical Systems and Numerical Analysis. Cambridge Monographs on Applie. Cambridge University Press, 1998.
[15] AR Humphries and AM Stuart. Runge-kutte methods for dissipative and gradient dynamical systems. SIAM journal on numerical analysis, (1)5:1445222 -185,, 1994.
[16] Andrew M Stuart. Numerical analysis of dynamical systems. Acta numerica, 3:467-572, 1994.
[17] Kai Fukami, Koji Fukagata, and Kunihiko Taira. Machine-leaminn--uusee spatio-temporal super resolution reconstruction of turbulent flows. Journal of Fluid Mechanics, 909, 2021.
[18] José Cardesa, Alberto Vela-Martin, and Javier Jiménez. The turbulent cascade in five dimen- sions. Science, 357(6353),782-784, 2017.
[19] Gary Chandler and Rich Kerswell. Invariant recurrent solutions embedded in a turbulent two-dimensional kolmogorov flow. Journal of Fluid Mechanics, 722:554 595, 2013.
[20] Péter Koltai and Stephan Weiss. Diffusion maps embedding and transition matrix analysis of the large-scale flow structure in turbulent rayleigh- bénard convection, Nonlinearity, 33(4):1723, 2020.
[21] Jason J Bramburger, J Nathan Kutz, and Steven L Brunton. Data-driven stabilization of periodic orbits. IEEE Access, 9:43504-43521, 2021.
[22] Jacob Page, Michael P Brenner, and Rich R Kerswell. Revealing the state space of turbulence using machine learning. Physical Review Fluids, 6(3):034402, 2021.
[23] Victor Churchill and Dongbin Xiu. Deep learning of chaotic systems from partially-observed data. arXiv preprint arXiv:2205.0833: 2022.
[24] Ashesh Chattopadhyay, Jaideep Pathak, Ebrahim Nabizadeh, Wahid Bhimji, and Pedram Hassanzadeh. Long-term stability and generalization of dnservationalllyconstti data-driven models for geophysical turbulence. arXiv preprint arXiv:2205.0460) 2022.
stochastic
[25] Jaideep Pathak, Shashank Subramanian, Peter Harrington, Sanjeev Raja, Ashesh Chatropadhyay. Morteza Mardani, Thorsten Kurth, David Hall, Zongyi Li, Kamyar Azizzadenesheli, et al. Fourcastnet: A global data-driven high-resolution weather model using adaptive fourier neural operators. arXiv preprint arXiv: 141 2022.
[26] Ashesh Chattopadhyay, Mustafa Mustafa, Pedram Hassanzadeh, and Karthik Kashinath. Deep spatial transformers for autoregressive data-driven forecasting of geophysical turbulence. In Proceedings of the 10th International Conference on Climate Informatics, pages 106-112, 2020.
[27] Jonathan A Weyn, Dale R Durran, and Rich Caruana. Improving data-driven global weather prediction using deep convolutional neural networks on a cubed sphere. Journal of Advances in Modeling Earth Systems, 2000 09, 2020.
[28] Dominik Schmidt, Georgia Koppe, Zahra Monfared, Max Beutelspacher, and Daniel Durstewitz. Identifying nonlinear dynamical systems with multiple time scales and long- rangee dependencies. arXiv preprini arXiv: 190.0334,, 2019.
[29] Manuel Brenner, Florian Hess, Jonas M Mikhaeil, Leonard F Bereska, Zahra Monfared, Po-Chen Kuo, and Daniel Durstewitz. Tractable dendritic mns for reconstructing nonlinear dynamical systems. In International Conference on Machine Learning, pages 2292-2320. PMLR, 2022.
[30] Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar, Neural operator: Graph kernel network for partial differential equations. arXiv preprint arXiv:2003.03485 2020.
[31] Nikola B. Kovachki, Zongyi Li, Liu Burigede, Kamyar Azizzadenesheli, Kaushik Bhattacharya, Andrew M. Stuart, and Anima Anandkumar. Neural operator: Learning maps between function spaces. In Preparation, 2021.
12

---

[32] Georgy Anatolevich Sviridyuk. On the general theory of operator semigroups. Russian Mathematical Surveys, 49(4):45-74, 1994.
[33] Yuying Liu, Nathan Kutz, and Steven L Brunton. Hierarchical deep learning of multi- scale differential equation time-steppers. Philosophical Transactions of the Royal Society A, (00 2022.
[34] Ashesh Chattopadhyay, Mustafa Mustafa, Pedram Hassanzadeh, Eviatar Bach, and Karthik Kashinath. Towards physics-inspired data-driven weather forecasting: integrating data assimila- tion with deep trantfoaltrmmsormmerr based u-net in study with era5. Geoscientific Model Development, (155)2221-22377 2022.
[35] Matthew E Levine and Andrew M Stuart. A framework for machine learning of model error in dynamical systems. arXiv preprint arXiv:2107.06655 2021. To appear in Communications of the AMS.
[36] Wojeiech Marian Czamecki, Simon Osindero, Max Jaderberg, Grzegorz Šwirszez, and Razvan Pascanu. Sobolev training for neural networks. arXiv preprint
2017.
[37] Alex Beatson, Jordan Ash, Geoffrey Roeder, Tianju Xue, and Ryan P Adams. Learning com- posable energy surrogates for pde order reduction. Advances in Neural Information Processing Systems, 33, 2020.
[38] Alexander J Smits, Beverley J McKeon, and Ivan Marusic. High-reynolds number wall turbulence. Annual Review of Fluid Mechanics, 43:353-378 2011.
[39] Olaf Ronneberger, Philipp Fischer, and Thomas Brox. U-net: Convolutional networks for biomedical image segmentation. In International Conference on Medical image computing and computer-assion asssssted intervention, pages 234 241. Springer, 2015.
[40] Xingjian Shi, Zhourong Chen, Hao Wang, Dit-Yan Yeung, Wai-Kin Wong, and Wang-chun Woo. Convolutional Istm network: A machine learning approach for precipitation nowcasting. arXiv preprint 106.:05006002144 2015.
[41] Junyoung Chung, Caglar Gulcehre, KyungHyun Cho, and Yoshua Bengio. Empirical evaluation of gated recurrent neural networks sequence modeling. arXiv preprint arXiv:1412.3555,
2014.
[42] Joseph Bakarji, Kathleen Champion, J Nathan Kutz, and Steven L Brunton. Discovering governing equations from partial measurements with deep delay autoencoders. arXiv preprint arXiv:2201.051366 2022.
[43] Kathleen Champion, Bethany Lusch, J Nathan Kutz, and Steven L Brunton. Data-driven discovery of coordinates and governing equations. Proceedings of the National Academy of Sciences, 116(45):22445-22455 2019.
[44] Daniel E Shea, Steven I. Brunton, and Nathan Kutz. Sindy-bvp: Sparse identification of nonlinear dynamics for boundary value problems. Physical Review Research, 3(2):023255, 2021.
[45] Tian Qi Chen, Yulia Rubanova, Jesse Bettencourt, and David K Duvenaud. Neural ordinary differential equations. Advances in neural information processing systems, pages 6571-6582 2018.
[46] Maarten V de Hoop, Nikola B Kovachki, Nicholas Nelsen, and Andrew M Stuart. Conver- gence rates for learning linear operators from noisy data. arXiv preprint arXiv:2 12515, 2021.
[47] Grigorios Pavliotis and Andrew Stuart. Multiscale Methods: Averaging and Homogenization, volume 53. 01 2008
[48] Jooyoung Park and Irwin W Sandberg. Universal approximation using citiial-bais-fuutttiii networks. Neural computation, 3(2):246-277, 1991.
13

---

[49] Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar. Fourier neural operator for parametric partial differen- tial equations. arXiv preprint arXiv:2010.08895, 2020.
[50] Colin Sparrow. The Lorenz equations: bifurcations, chaos, and strange attractors, volume 41. Springer Science & Business Media, 2012.
[51] Edward N Lorenz. Deterministic nonperiodic flow. Journal of atmospheric sciences, 20(2):130- 141, 1963.
[52] Aly-Khan Kassam and Lloyd N. Trefethen. Fourth-order time-stepping for stiff pdes. SIAM Journal on Scientific Computing, 26(4):1214-1233, 2005.
[53] Gabriel J Lord, Catherine E Powell, and Tony Shardlow. An introduction to computational stochastic PDEs, volume 50. Cambridge University Press, 2014.
[54] Felix A Gers, Jürgen Schmidhuber, Fred Cummins. Learning to forget: Continual prediction with Istm. 1999.
[55] Nikola Kovachki, Samuel Lanthaler, and Siddhartha Mishra. On universal approximation and error bounds for fourier neural operators. Journal of Machine Learning Research, 22:1-76, 2021.
[56] Andrew Stuart, Perturbation theory for infinite dimensional dynamical systems. 1995.
14


Reviews:
[{"review_id": "gmdLSnUEqSB", "reviewer": "Reviewer_3spS", "summary": "This paper proposes a machine learning framework, named Markov neural operator (MNO), to learn the Markov operator of dissipative chaotic systems. Different from the conventional neural networks, they use Sobolev norms in operator learning and add dissipativity losses to ensure dissipativity. They experimentally show that the MNO can accurately approximate the global attractor and estimate various statistics of the invariant measure for dissipative chaotic systems. Besides, they provide a theoretical guarantee that this model is rich enough to approximate many chaotic dynamical systems for an arbitrarily long period.", "questions": "1. How the choice of the inner and outer radii of the shell of the dissipativity losses influences the performance of the MNO? Whether the choice significantly affects the MNO’s step-wise error and the estimation of the statistical properties of the attractor? How to choose the radii based on the training data?\n2. How to balance the influences of the Sobolev and the dissipativity losses? Is there a hyperparameter to adjust their proportion in the loss function?\n3. Does the choice of lambda significantly influence the performance of the MNO? How to choose the optimal lambda according to the training data?\n4. Is it possible to provide a method to judge whether the MNO can be directly applied to a certain target system?\n", "limitations": "The MNO can be used only when the dynamical system is Markovian, but actually we do not know whether that assumption holds in a model-free task. This paper lacks a method to judge whether the MNO can be directly applied to the certain target system.", "rating": 7, "confidence": 4, "soundness": 2, "presentation": 3, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "This is a novel method to approximate the global attractor and estimate the invariant statistics of dissipative chaotic systems. They pioneered the combination of the Sobolev losses and the dissipativity losses, which improves the performance for preserving the invariant statistics. The dissipativity losses they proposed can enforce the dynamics close to the attractor, which improves the robustness against large perturbations. They also investigate the effect of time steps and different Sobolev losses in operator learning, which makes this work more systematic. However, there are still some confusing details. The choice of the hypermeter lambda and the inner and outer radii of the shell in the dissipativity losses remains to be discussed. There is no description of the proportion of the Sobolev and the dissipativity losses in the loss function. In addition, sometimes we do not know whether the target system is Markovian in a model-free task. As a result, this paper lacks a method to judge whether the MNO can be directly applied to a certain target system.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "SUPNpaxV3Zs", "reviewer": "Reviewer_NdmN", "summary": "This article discusses robust estimation of chaotic dynamical systems, in particular their invariant statistics, using Markov neural operators (MNOs). In case of simple ODE systems, feedforward NN  instantiate the MNO, while Fourier operators are used for spatially extended (PDE) systems. The major innovation in the present work is the addition of a ‘dissipativity loss’ that enforces attracting behavior within a shell around the data, and thereby tends to enforce globally attracting properties. Essentially, the dissipativity loss tends to penalize larger flow vectors away from the data (why this works was not fully clear to me, however). This enables to learn highly challenging systems like the strongly turbulent Kolmogorov flow.", "questions": "1) If I understood correctly, a crucial premise of this approach is that there is indeed only one globally attracting limit set. This will be violated in most biological systems at least, for instance in molecular biology or neuroscience, as well as in many other complex dynamical systems of current interest, like climate models. For many real-world data from complex systems we would not know this in advance. Isn’t this a serious limitation of the approach?\nAlso, how does the vector field far from the attractor actually compare to that of the true system under this approach?\n\n2) I felt that a lot of relevant literature was left aside (all the work by the Brunton group or Edward Ott’s group for instance, https://arxiv.org/abs/2201.05136, https://arxiv.org/abs/1904.02107, https://arxiv.org/abs/1710.07313, https://arxiv.org/abs/2110.07238, https://arxiv.org/abs/1712.09707, https://arxiv.org/abs/1910.03471, https://arxiv.org/abs/2110.05266, https://arxiv.org/abs/2005.10756, https://arxiv.org/abs/1707.01146, https://arxiv.org/abs/2207.02542, or the recent developments surrounding neural ODEs and PDEs, e.g. https://iopscience.iop.org/article/10.1088/1367-2630/abeb90). The authors are not the first to put the focus on invariant statistics; power spectra, Lyapunov spectra, or geometrical statistics have been discussed before in the context of reconstructing dynamical systems. Likewise, although the novelty about Theorem 1 may be its extension to spatio-temporal operators (function spaces), similar ideas I believe appeared before in the works of Funahashi and Nakamura or, more recently, by Raginsky and colleagues (https://proceedings.mlr.press/v120/hanson20a.html and references therein). Also, much stronger competitors are likely to be found among these previous methods than the plain LSTMs and GRUs used in Fig. 4 or the models used in Table 1. Many groups (like Brunton or Vlachas) have tested their approaches on similar systems like the Lorenz-63, Kuramoto-Sivashinsky or Navier-Stokes systems. These, in my mind, would have been the more relevant comparisons. This is especially relevant when the claim is made that other ML approaches fail to capture invariant statistics (lines 255-256).\n\n3) If I understood correctly, the specific regularization solution is partly also motivated by the fact that iterating the MNO across many time steps (eq. 3) may easily lead to divergence. But this it seems to me is really a standard problem in RNN training, especially if the loss is only based on consecutive time steps. Specific algorithms to deal with this issue had been suggested before in the context of dynamical systems, e.g. https://arxiv.org/abs/2110.07238. Would be nice to know how these approaches differ in performance \\& properties.\n\n4) More methodological and algorithmic details in the main text would be appreciated. What is the loss function, apart from the novel regularization term? What precisely is the setup and parameterization of the forward NNs and the Fourier NOs, how were hyper-parameters of the architectures determined (especially $\\lambda$ in the  dissipativity loss)?\n\nMinor points:\n\n- In my understanding RNNs also learn a Markov operator that is used as in eq. 3, what is the difference?\n\n- Fig. 3 provides just single examples, likewise it appears Fig. 4. It would be good to see statistically more comprehensive evaluations with error bands on many different training runs. Also, it would be  good to see reconstructions of the actual attractors as well (the Lorenz is barely visible in Fig. 3). \n\n- I found it somewhat difficult to follow parts of sect. 4.2 and Fig. 4, simply cos in my mind too few details were provided. What, for instance, is the ‘residual model’?\n", "limitations": "The major limitation I pointed out above was not discussed.\n", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths:\nIn general I find this a sophisticated, concisely written and well presented contribution to the field from which I learned various things. This specific idea of enforcing dissipativity seems novel and useful to me.\n\nWeaknesses:\n... at the same time, for many systems this global form of dissipativity seems to be a quite strong and often too restrictive assumption, see my detailed Q below.\nMuch related lit. is not discussed, some methodological details are missing, for some results only single examples but no stats are presented.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "Q1-_MnTeh6_", "reviewer": "Reviewer_MCtx", "summary": "This paper proposes a machine learning framework (MNO) by using neural networks or Fourier neural operators to learn the underlying Markov operator of dissipative chaotic systems. Numerical results show that MNO has lower loss and outperforms the other neural network models, including U-Net, LSTM-CNN, and GRU.", "questions": "1. Could you provide the formula of Sobolev Loss and Dissipativity Loss in Figure 2. I guess that the Dissipativity Loss is given by Equation 5. Then, is the ground truth exactly match the Equation 5? If not, does \"lower loss\" make sense? I think it might be more appropriate to treat it as a regular term.\n\n2. How about using continuous models to learn continuous dynamics？For example, Neural ODE or its structure-preserving extensions [1]. In [1], they also consider learning Lorenz-63 system. For infinite dimensional PDE systems, in my view, we can directly replace the neural network in a continuous model by a neural operator to obtain a infinite continuous model. \n\n[1] Haijun. Yu, Xinyuan. Tian, Weinan. E, Qianxiao. Li, Onsagernet: Learning stable and interpretable dynamics using a generalized onsager principle, PHYSICAL REVIEW FLUIDS 6, 114402 (2021)\n", "limitations": "Yes, the authors addressed the limitations.", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "This paper proposes to learn the underlying Markov operator of dissipative chaotic systems using deep learning. Several architectures are investigated in detail. The results is interesting, however, I'm not sure whether learning Markov operator of dissipative chaotic systems is relevant for deep learning. In my view, this paper is more suitable for some dynamical journals, such as SIAM Journal on Applied Dynamical Systems, Chaos, ... .", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "FWQNWdTEl-", "reviewer": "Reviewer_Hc1s", "summary": "This is a very interesting paper and addresses an important challenge in terms of long-term stability of high-dimensional turbulent flow. The authors show that adding a dissipative term in the loss and training it with sobolev loss would improve the stability of the autoregressive data-driven models (here they show that a Markov neural operator) is the best skeleton model. \n", "questions": "The paper discusses and investigates realistic metrics pertinent to the problem and I've thoroughly enjoyed reading it. However, there are still some questions/clarifications that remain for the authors to justify their claims. Herein, I summarize some of the main issues that I find. Beyond that, there is a large class of literature in weather/climate modeling and fluid dynamics that talks about these issues at large; I suggest that the authors look at some of the data-driven weather prediction papers that talk about long-term instability. Adding that context to this paper would make it even more meaningful for readers who would actually care about why this is an important problem for the physical sciences. \n\n\n1. Line 2. -- The issue is not just with time-stepping, but with the fact that because of a positive Lyapunov exponent the initial perturbation grows exponentially. \n\n2. Line 19-20. There seems to be other work (one at least) that looks at long-term statistics of chaotic systems, e.g. PDFs including the tails that capture extreme events. Even more so, it turns out RCs can do both short- and long-term seamlessly. https://npg.copernicus.org/articles/27/373/2020/. Surely, there could be other papers that do the same. The authors' statement is not accurate here. They are however correct when they say that RCs definitely do not scale for more complex problems. Even so, the authors are still looking at toy systems that RCs excel at (except for the NS case). \n\n3. Line 25-26 What's the rationale behind calling it \"ill-posed\", citations to relevant papers are essential. \n\n4. Line 40. Typo. charges->changes\n\n5. For Lorenz 63, simpler methods are long-term stable, even for KS, and for other more complex Lorenz-type systems as shown above. The first two examples are not interesting or relevant when we want to talk about long-term statistics. 2D Kolomogrov is a great example, and there are many other relevant problems that have been studied for long-term stability along those lines, e.g., QG turbulence, other climate models, and real weather data. There is a big community working on this specific problem. The authors need to look at that literature:\n(i) https://arxiv.org/abs/2205.04601\n(ii) https://arxiv.org/abs/2202.11214\n(ii) https://arxiv.org/abs/2202.07575\n(iv) https://dl.acm.org/doi/10.1145/3429309.3429325\n(v) https://agupubs.onlinelibrary.wiley.com/doi/full/10.1029/2020MS002109\n\n\nIn the purely ML community, there are a few other papers as well (who tried and improved stability but eventually couldn't make it infinitely stable) for Kolomogrov turbulence (more complicated cases than the current paper). The manuscript would benefit from putting the problem in the context of this bigger goal. \n https://arxiv.org/abs/2112.15275\n\n6. The effect of \"h\" is very interesting in the KS problem. There are a few other papers that have found results that might be contradicting this specific paper, or maybe not if I have misunderstood. But either way, that's good. It generates discussion and the importance of exploring this problem more carefully. However, I feel that the authors should acknowledge those papers as well. While not done in the context of operator modeling, they have been done on more complex, realistic, and practical problems of chaos and turbulence. Here are a few of those papers :\n\n(i) https://royalsocietypublishing.org/doi/10.1098/rsta.2021.0200 (Section 2)\n(ii) https://agupubs.onlinelibrary.wiley.com/doi/epdf/10.1029/2020MS002203 (Figure 2)\n(ii) https://gmd.copernicus.org/articles/15/2221/2022/ (Section 4.3)\n\n7. I appreciate how the authors interpret stability as something beyond just running the model and not getting \"NaN\". Stability means that the flow is both physical and non-drifting and has the right variability. However, I am still not 100% convinced that this model can serve as a real long-term stable model for a real complex flow. Here are a few things the authors can show to remove my doubts:\n    --- While the model captures the spectrum quite well, it deviates in the high-wavenumbers, an essential reason for instability (Figure 4b). Can the authors show how the spectrum deviates as we move forward in time? Basically, take the best model and plot the output's spectrum at various time steps during prediction and compare it with the truth. Furthermore, I don't think that KS is that complicated a system, wherein \"spectral bias\" (the fact that you lose high wavenumbers as you integrate in time) is that big of a problem. However, in the N-S case with high Re, this would become a definite problem especially because you have inverse cascade. Also, for KS, I think Pathak et al, PRL, had long-term statistics, without doing anything? Does this KS have the same dimension or is it higher?\n   ----  It's great that the authors look at the PC1 autocorrelation curve. It must be pointed out that the authors' model does not capture it, or am I mistaken? I also understand that it's difficult to capture PC1 autocorrelation in KS because, in reality, I don't think PC1 in this system has any physical meaning (e.g., in a real system like Earth's climate PC1 has connections to ENSO and so on..). Still the fact that the model doesn't go absolutely haywire is promising. Yet, the authors should clarify this\n  ---  I don't think the PDF is also captured correctly. Please plot PDF in semilog, so that the tails (extreme events) are more prominent, that's what matters. But here, I think even the bulk is wrong in some portions. It's still impressive that you have the right range of values. \n\n8.  For the N-S case, I have some major confusion. In the main paper (not the appendix), the spectrum is not plotted. In the appendix, I see some images with Re=40. Re=40 is not nearly turbulent enough. That's hardly interesting in the context of turbulence and it's long-term predictability. Even Re=500 is not a realistic set-up. I must say here, that I understand the difficulty to scale up the turbulent flow to high Re, it's both difficult computationally, and probably very hard to predict. But based on the context of the problem, if this model really is going to serve practical problems in terms of long-term stability, the 2D N-S case should have Re=10000 or higher. That's the regime when we have all the complexities of this flow including when the effect of inverse cascade would be prominent. Here're a few papers that look at NS cases in realistic regimes.  \n    (i) https://www.pnas.org/doi/10.1073/pnas.2101784118\n    (ii) https://www.sciencedirect.com/science/article/pii/S0021999122001528\n    (iii) https://www.cambridge.org/core/journals/journal-of-fluid-mechanics/article/abs/subgrid-modelling-for-twodimensional-turbulence-using-neural-networks/10EDED1AEAA52C35F3E3A3BB6DC218C1\n\n9. For the N-S case, the energy spectrum at different time steps needs to be shown for the best model and the baseline. That would clearly show how well the spectrum is captured. It would also only be interesting for a case where the flow is turbulent (Re>10000)\n\n10. Again, the PDFs should be shown in semilog plots to highlight the tails. Also, what's going on with the PDF of KE for the N-S case? It's not really getting captured right? In this case, as well, there's no dominant mode of variability, so PC1 does not have any meaning right? But still, plot the true PC1 autocorrelation plot to show how well/poorly it's captured. \n\n11. The authors might benefit from looking at another system. A system that has a dominant mode and coherent structures while being turbulent enough. I think the discussion around stability would have far more importance there with PC1, autocorr, etc even if they are not perfectly captured. \n\nFinally, despite my long review and many questions, I enjoyed reading the paper. The authors have done a great job in addressing the question of stability \"correctly\"  by looking at the right set of metrics. While there are these lingering questions, once properly addressed this would be an important contribution to the ML and physics community. ", "limitations": "N/A", "rating": 5, "confidence": 5, "soundness": 3, "presentation": 3, "contribution": 4, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "The major strength of this paper is that it addresses a difficult problem in the physical sciences that have a lot of stakeholders and look at pertinent and physically justified metrics to support their claim. \n\nHowever, the major weakness of this work is that 2 of the 3 systems considered are not realistic enough to show issues with instability (the problem at hand) and the 3rd system is also not the best realistic case. Other papers have looked at more complex systems. The weakness still does not take away a very good set of analysis and transparent results that the authors have shown. Still it needs some work before getting accepted. ", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}]

Author rebuttal:

Decision: Accept