FINDE: Neural Differential Equations for Finding and Preserving Invariant Quantities
Authors: Takashi Matsubara
Takaharu Yaguchi
Venue: NeurIPS 2022 Submitted (2022)
Keywords: neural ordinary differential equations
first integral
conservetaion law

Abstract: Neural networks have shown promise for modeling dynamical systems from data. Recent models, such as Hamiltonian neural networks, have been designed to ensure known geometric structures of target systems and have shown excellent modeling accuracy. However, in most situations where neural networks learn unknown systems, their underlying structures are also unknown. Even in such cases, one can expect that target systems are associated with first integrals (a.k.a. invariant quantities), which are quantities remaining unchanged over time. First integrals come from the conservation laws of system energy, momentum, and mass, from constraints on states, and from other features of governing equations. By leveraging projection methods and discrete gradient methods, we propose first integral-preserving neural differential equations (FINDE). The proposed FINDE finds and preserves first integrals from data, even in the absence of prior knowledge about the underlying structures. Experimental results demonstrate that the proposed FINDE is able to predict future states of given systems much longer and find various quantities consistent with well-known first integrals of the systems in a unified manner.

Neural Differential Equations for Finding and Preserving Invariant Quantities
FINDE:
Anonymous Author(s)
Affiliation
Address
email
Abstract
Neural networks have shown promise for modeling dynamical systems from data. Recent models, such as Hamiltonian neural networks, have been designed to ensure known geometric structures of target systems and have shown excellent modeling accuracy. However, in most situations where neural networks learn unknown systems, their underlying structures also unknown. Even in such cases, one can expect that target systems are associated with first integrals (akka. in- variant quantities), which are quantities remaining unchanged integrals come from the conservation laws of system energy, momentum, and mass, from constraints states, and from other features of governing equations. By leveraging projection methods and discrete gradient methods, we propose first integral-presern neural differential equations (FINDE). The proposed FINDE finds and preserves first integrals from data, even in the absence of prior knowl- edge about the underlying structures. Experimental results demonstrate that the proposed FINDE is able to predict future states of given systems much longer and find various quantities consistent with well-known first integrals of the systems unified manner.
2
4
6
time. First
8
10
11
12
13
14
15
18
Introduction
18
Although neural networks have achieved remarkable results in image and natural language pro-
cessing [1 28 they have also been actively investigated for modeling dynamical systems [41).
20
Target systems include the chemical dynamics to accelerate computer simulations [46), the climate dynamics for climate change prediction and weather forecasting 47
52), and the physical dynamics 22 of vehicles and robots for optimal control [41] Their history dates back to at least the 1990s, and
23 many approaches have been proposed far (see
55 for example). Recently, 24 neural ordinary differential equation (NODE) has redefined neural networks for continuous-time
12
35
49
25 dynamics S A target system is described by an ordinary differential equation (ODE) d
dr u = f(t, u). where denotes the system state. Then, NODE replaces the vector field with neural network and employs a numerical integrator to obtain a solution uut)..
27
28
Most world systems are associated with first integrals (a.k.a. invariant quantities), which are quantities remaining unchanged over time [27]. a system has first integral V(u). the solution ul()) for the initial condition (0)) remains at contour line V(u(()) (u(()) over time. Many previous studies have attempted to learn a target system accurately by incorporating prior knowledge about first integrals. Greydanus et al. (26] proposed Hamiltonian neural network (HNN), which
29
31
33 employs a neural network to approximate Hamilton's 's equation, thereby conserving the system energy called the Hamiltonian. Finzi al. [19] proposed neural network architectures that conserve lincar 35 and angular momenta by utilizing the graph structure. Finzi et al. (00 also extended HNN to a system with holonomic constraints, which lead to first integrals such pendulum length. Matsubara et al.
Submitted to 36th Conference on Neural Information Processing Systems (NeurIPS 2022). Do not distribute.

---

Table 1: Comparison between Related Studies on Preservation of First Integrals.
learning invariants
exact conservation
monentum
constraint
energy
musS
NODE
HNN
LicCony
19
DGNet
38
CHNN
20
continuous FINDE (proposed)
discrete FINDE
(proposed)
/
37
proposed model that preserves the total mass of a discretized partial differential equation PDE.. These studies have demonstrated that neural network with more prior knowledge about first integrals predicts the dynamics of the target system more accurately. See Table for comparison.
38
39
40
Previous studies have mainly attempted to preserve known first integrals. However, in situations where neural network learns unknown target system, is naturally expected that first integrals associated with the target system are also unknown, and is not clear which of the above methods are available. Given the above, this study proposes First (FINDE) to find and preserve first integrals from data. FINDE has the following ad advantages.
41
42
43
Neural Differential Equation
44
Learning First Integrals For modeling continuous-tine dynamics with known first integrals, many studies have designed architectures or operations of neural networks [13 type of first integral, one dedicated method was proposed. However, the properties of a target system are generally unknown in practice. In contrast, the proposed FINDE finds various kinds of first integrals from data in unified manner and preserves them in predictions. A symbolic regression confirms that the learned first integrals are consistent with well-known first integrals of target systems.
46
19 20
26
38 For each
48
49
50
51
The proposed FINDE be combined with previously proposed neural networks designed to preserve known first integrals, such as HNN. Therefore, FINDE available in various situations.
Combination with Known First Integrals
52
53
54 Exact Preservation of First Integrals Even if first integral is associated with a continuous-time system, it is destroyed after the system discretized time for computer simulations. This true even when using a symplectic integrator, which preserves the system energy only approximately (27]
55
56
57 By leveraging discrete gradients (38), the discrete-time version of FINDE preserves first integrals exactly (up to rounding errors) in discrete time and further improves the prediction performance.
58
2
Background and Related Work
59
80
First Integrals Let us consider time-invariant differential system d dimensional manifold M. where denotes the system state and f M - T.M represents vector field on the manifold M. The manifold M can be M = S x IR for a pendulum. In this 63 paper, we suppose the manifold M be a Eucleadian space N for simplicity.
di f(u) on an N-
81
82
Definition 1 (first integral). A quantity V M R is referred d d de = f(u) if it remains constant along with any solution u(t), iee.
first integral of system dt V(u)== =
85
If differential system मे du- f(u) is associated with K functionally independent first integrals 87 Vi..., Vx. the solution (1) given an initial value ug stays at the (N - K )dddmenniona submanifold
88
= (== Vi(n)= = (uo)). VK(u)) = Vx(uo)).
(1)
89 The tangent space T,.M' C T.M of the submanifold M' M at a point u is the orthogonal complement to the space spanned by the gradients VVi (u) of the first integrals Vn, for k = K,
71 that is,
TuM' {w T.M: VV)(u) 0 for k 1,...,Ki K}
(2)

---

72 If a quantity V& is first integral of the system di u = f(u), the time-derivative f at point u is on the
d
73 tangent space Tu M', being orthogonal to the gradient VV; of the first integral VE. Then, it holds
74 that di (u)) = VV (w) मe (u) f(u)= 0.
One of the most well-known first integrals is the Hamiltonian H, which represents the system energy 76 of a Hamiltonian system. Noether's 's theorem states that a continuous symmetry of system leads to a
77 conservation law (and hence first integral) [27]: Hamiltonian system is symmetric translation in time and conserves the Hamiltonian. Symmetries to translation and rotation in space lead to the conservation of linear and angular momenta. Not all first integrals are related to symmetries. A pendulum can be expressed in Cartesian coordinates, and then the rod length constrains the mass position. This kind of constraint is called a holonomic constraint and leads a first integral. model for disease spreading called an cirepptble-infected-rcooo
78
79
80
81
82
(SIR) model and the dynamics of 83 chemical reactions have the total mass (population) as first integral. Also for system described by a PDE, the total mass is sometimes first integral 23 See Appendix A for theoretical classification
84
a5 of dynamics.
86
First Integrals in Numerical Analysis For computer simulations, a differential system is dis- cretized in time and solved by numerical integration. Then, the geometric structures of the system often destroyed, and most first integrals are no longer preserved. A common remedy is sym- plectic integrator, which preserves the symplectic structure and integrates a Hamiltonian system accurately [27]. However, Ge-Marsden theorem states that symplectic integrator conserves the Hamiltonian only approximately [56]. Hence, many numerical schemes have also been investigated for preserving first integrals exactly, while they cannot preserve the symplectic structure.
87
89
90
91
92
93 Let superscript s denote the state u² or time t" at s-th time step, and At" = = - - denote a time step size. A projection method predicts next state ū u- from the current state u² using numerical integrator and projects it onto the submanifold M obtaining the projected state u++1 that preserves the first integrals V2 [24] (see also [27] Section IV.4]). In particular, the projected state is obtained by solving the optimization problem
95
97 u ++1
u ++1 arg min I /0+1 - - 1 subject to Vi((u = V (w)) for k = I,.... K.
(3)
A local coordinate method defines coordinate system to the neighborhood of the current state u² and integrates a differential equation on it [43 (see also (27 Section IV:5)). A discrete gradient method defines discrete analogue to a given differential system and integrates it in discrete time (6
39
100
101
23
(45.. This method eliminates numerical errors caused by temporal discretization and is used preserve the Hamiltonian exactly (up to rounding errors) in discrete time.
29
102
103
Except for DGNet, which used discrete gradients to preserve the Hamiltonian (38. all the above methods have never been applied to neural networks due to difficulties that we will introduce later. To our best knowledge, the discrete- time version of FINDE is the first projection method for dynamical systems modeled using neural networks.
105
106
107
Preservation of First Integrals by Neural Networks NODE defines an ODE using neural net- work in the most general way with no associated first integrals X NODE is a universal approximator to ODEs [51], and it can approximate any ODE with arbitrary accuracy if there is an infinite amount
108
109
110
of training data. In practice, the amount of training data limited, and prior knowledge about the
target system is helpful for learning (see (48) for the case with convolutional neural networks). HNN assumes the target system to be a Hamiltonian system in the canonical form (26). HNN guarantees various properties of Hamiltonian systems by definition, including the conservation of the energy and the preservation of the symplectic structure in continuous time [27]. Some studies employed a symplectic integrator for HNN to preserve the energy and symplectic structure with smaller numerical errors 10]. LicConv and EMLP-HNN employed neural network architectures with translational and rotational symmetries to preserve momenta [19 21. CHNN incorporates a known holonomic constraint in the dynamics 20 Deep conservation extracts latent dynamics of a PDE system and
112
113
114
116
118
preserves quantity of interest by forcing its flux to be zero [4]. HNN++ also guarantees the conservation of the mass in PDE systems by using a coefficient matrix derived from differential
120
121
operators (88)
122
Several studies proposed neural networks to learn Lyapunov functions, which are expected to be non-increasing over time, in contrast to first integrals [37 50). If the state moves in the direction of
3

---

124
increasing the function, it is projected onto or moved inside the counter line of the gradient of the Lyapunov function. Their idea is similar to the continuous time version of FINDE but limited to a single non-increasing quantity in continuous time. On the other hand, our proposed FINDE preserves multiple quantities both continuous and discrete time.
125
126
127
128
Previous studies aimed to preserve known first integrals. Moreover, except for DGNet [38], all the above methods suffer from numerical errors caused by temporal discretization. In contrast, our proposed FINDE learns first integrals from data and can eliminate discretization errors.
129
130
3
First Integral-Preserving Neural Differential Equation
131
132
The main purpose is to find and preserve first integrals from data by neural networks. We suppose that a target system has at least K unknown functionally independent first integrals. Even when a NODE learns the target system, it is not guaranteed to learn these first integrals. Hence, we introduce a neural network with K outputs, each of which is expected to learn one of first integrals expressed as Vv R + R for k = ....., K. We denote the set of first integrals by vector Vx(u))) Then, the submanifold M' is defined using the neural
133
134
135
V (u)) = (V) V/2 (
137
network V as in Eq. ①..
139
Because there is no way to define local coordinates on such submanifolds, a local coordinate method is not applicable. When using a projection method, the optimization problem in Eq. (3) should be solved at every training iteration as well as the prediction phase. Optimization problems are computationally expensive, and common libraries for neural networks do not provide backpropagation algorithms for optimization problems 42 Until a recent study has proposed an algorithm [88.
140
141
142
143
144
there was no way to obtain discrete gradients of neural networks. Because of these difficulties, no methods for preserving first integrals have been applied to neural networks. By leveraging a projection method and a discrete gradient method, we propose FINDE as follows.
145
147
3.1
Continuous FINDE: Time-Derivative Projection Method
First, we propose dime-derivative projection method called continuous FINDE (cFINDE) for neural networks, which projects the time-derivative onto the tangent space M.. While it still suffers from numerical errors, it is sufficient to find first integrals from data.
149
150
We suppose that neural network called base model defines the time-derivative f: RN - R.V. Then, we define the iim-deervative f of the cFINDE d dr u = f(u) as
152
f(u)= f(u) Ek=1 d (v) f(()) - M(u)) T((x), and X(u) (Ai(u) 22(u) ... Ak(u)) If V
(4)
where A is a Lagrange multiplier, M- remains constant,
eV Du
154
0
(u(t)) M()) # M()f((u) M (())(f(u)
(()
(5)
where 0- = (0 ... T)T. By transforming Fq. (5). we obtain the Lagrange multiplier X(u)=
155
(M) T)
M (w) (u). By eliminating it, the cFINDE d uu fw) is given by
156
s(u) = - (I(")))( where Y(u)= = (MM)) (M Remark 1 (continuous- first integral preservation). The cFINDE de u f() preserves all first integrals Vi k 1,.... K in continuous time, ie., ae Vı = 0.
(6)
158
The base model f can be NODE, an HNN, or other models depending on available prior knowledge. Also, if a first integral is already known, one can use it directly as one of first integrals Ve instead of learning it using neural network. Note that even though the base model f is an HNN, due to projection, the cFINDE ʃ is no longer a Hamiltonian system the strict sense.
159
160
162
Compared the base model f. the cFINDE requires the additional computation of the neural network V. several matrix multiplications, and an inverse operation. The inverse operation needs computa- tional cost of ( which is not costly if the number K of first integrals small. For satisfying the constraints and geometric structures, many previous models also need the inverse operation, such
163
164
185
166
167
as Lagrangian neural network (LNN) 13, neural symplectic form 99, and CHNN 20
The algorithm proposed in 2 might work, but it is outside the scope of this paper.

---

168
3.2
Discrete FINDE: Discrete-Time Projection Method
189
To eliminate numerical errors caused by temporal discretization, we employ discrete gradients and propose a projection method called discrete FINDE (dFINDE).
170
171
45]. Recall that a gradient VV of a function V R + R can be regarded as a function I + IR that satisfies the T d de Analogously, a discrete gradient V is defined as follows.
discrete gradient VV discrete analogue to a gradient VV 63
29
chain rule d de V(u)= = (v)
173
Definition 2 (discrete gradient). A discrete gradient VV of a function V RN R is function RN x RN + RN that satisfies
175
V(v) - V(u) = VV T)T(u u) and VV(u, u) VV(u).
(7)
The first condition is a discrete analogue to the chain rule when replacing the time-derivatives d V ät and d u with finite differences (V(u) V(u)) and (v - u), respectively, and the second condition ensures the consistency with the ordinary gradient VV. A discrete gradient V is not uniquely determined and has been obtained manually. Recently, the automatic discrete differentiation algorithm (ADDA) has been proposed (38), which obtains a discrete gradient of a neural network in a similar way to the automatic differentiation algorithm [ 42] The discrete gradient is defined in discrete time, and hence a numerical integration using the discrete gradient is free from numerical errors caused by temporal discretization. See Appendix B and the references [6] 23 (88 for more details.
176
177
178
180
182
Following 11 15), we introduce discrete analogue to the tangent space T,,M called the discrete tangent space T.. In particular, for pair (v,, u) E M' of points, it is defined as
184
185
Tw = fw E RR,, (v, u) w = 0 for A: 1,...,K).
(8)
If the finite difference (u++1 u²) between the predicted and current states is on the discrete tangent space ( the first integrals V&2 are preserved because V (x=*+) ) - Vu((") =
187
VV Vu* u²) (u*+1 u²) 0. Note that similar concepts defined in different ways are also
189
referred to as discrete tangent spaces 16).
Let & denote a discrete-time base model that satisfies #-++_/" (((w"(A)))) where ū++1 denotes Ar"
190
the predicted state. We assume that the base model & is composed of continuous base model and numerical integrator. Then, the dFINDE u=+1 -"' v'uu+1, u²; At*) given by AP
192
b(u*t1 u²: At") = (uu;; At*) - M Muttt u²) X(u
199
where Mu )²) = ( the preservation of the first integrals Va leads
As is the case continuous time,
V(u++))-V V M(u+1 u²) Ar
t Muu+1, u (u A1))
(10)
195
Substituting Eq. and eliminating the Lagrange multiplier X, obtain
v(u=+1 u²; AP))
Y (u u²) (u""; Atr) where = -M M T
- (11)
Remark (discrete-time first integral preservation). The dFINDE ===l_u" vb(u-! uu: At preserves all first integrals Vi k I,..., K in discrete time, i.e., Vi(u+t1) Va(u²) 0.
198
Due to projection, dFINDE can be regarded as projection method using discrete gradients. For the base model vi the continuous- time base model ff can be NODE, an HNN, or other models, and the numerical integrator can be Runge-Kutta method, the leapfrog integrator, others.
199
Because dFINDE is an implicit method, it is computationally expensive for prediction. However, the next state u 0+ is given for training, and the ADDA can explicitly obtain the discrete gradient w.r.t. the pair (u++1 u²) as well as its computational graph. Thus, dFINDE can be computed explicitly and optimized by standard backpropagation algorithms. Moreover, we suppose that dFINDE projects
202
203
204
the finite difference v only at every time step, whereas cFINDE projects the time-derivative ff at every substep inside numerical integrator. Therefore, dFINDE is less computationally expensive
205
206
than cFINDE for training. In contrast, a typical projection method requires much computational
208
cost to solve an optimization problem for training, and standard backpropagation algorithms are not
applicable to it.
5

---

Table 2: Datasets, Dynamies, and First Integrals.
First Integrals
Dataset
Dynamics
N
Energy Momentum
Mass
Constraint
Twa-body problem
Canonical Hamiltonian
8
Discretized KdV equation
non-canonical Hamiltonian
50
Double pendulum
Poisson
8
FitzHugh -Naguoo model
Dirac
4
Remark 3 (trainability). The dFINDE can be trained using the standard backpropagation algorithm, whereas a straightforward application of a projection method cannot.
211
4
Experiments
212
213
4.1
Experimental Settings
214
Target Systems We evaluated FINDE and base models using datasets associated with first integrals, summarized in Table A gravitational two-body problem (2-body) on 2-dimensional configuration space is a typical Hamiltonian system in the canonical form. In addition to the total energy, it has first integrals related to symmetries in space, namely, the linear and angular momenta. The Korteweg-De Vries (KdV) equation is a PDE model of shallow water waves. This is a Hamiltonian system in a non-canonical form and has the Hamiltonian, total mass, and many other quantities as first integrals. double pendulum (2-pend) is Hamiltonian system in polar coordinates. However, we transformed it to Cartesian coordinates; it of two rods work as holonomic constraints and lead to four first integrals. The model is biological neuron model an electric circuit, which exhibits rapid and transient change of voltage called a spike. As an electric circuit, the currents through and voltages applied to the inductor and capacitor can be regarded system states, and the states are constrained by the circuit topology and Kirchhoff's current and voltage laws. Then, this system has a state of four elements and two first integrals. Due energy dissipation the resistor, the model is not a Poisson system, but one can find a Dirac structure [53]. See Appendix C for more details.
215
216
218
220
longer Hamiltonian system but Poisson system. The lengths
222
224
226
227
228
229
Implementation We implemented the proposed FINDE and evaluated under the following settings. We implemented all codes by modifying the officially released codes of HNN (26] and DGNet [38 We used Python v3.8. 12 with packages scipy v1.7.3, pytorch v1.10.2, torchdiffeq v0. I.I.I.I functorch v0 preview, and gplearn v0.4.2. We used the Dormand-Pringt method (dopri5) [[]]
231
233
as the numerical integrator, unless otherwise stated. All experiments were performed on a single NVIDIA A100 provided by (ANONYMOUS PROVIDER).
Following HNN [26] and DGNet 388. we represented the first integrals V. NODE, and HNN H using fully-connected neural networks with two hidden layers. Each hidden layer had 200 units and preceded a hyperbolic tangent activation function. Each weight matrix was initialized as an orthogonal matrix. The input was the state u, and the output represented the first integrals V for FINDE, time-derivative f for NODE, and the Hamiltonian H for HNN. For the KdV dataset, we used a I I-dimensional convolutional neural network (CNN), each of whose layers had a kernel size of 3. The double pendulum is second-order system, implying that the time-derivative de of the position 241
235
236
237
240
242
is known as the velocity v. Hence, we treated only the acceleration d: v as the output to learn. This assumption slightly improved the absolute performances but did not change the relative trends.
243
244
used the /-step error as the loss function to be minimized. particular, enrrr (MSE) between the ground truth state uer and the state "fred predicted from the previous step
the mean squared
246
The base model and FINDE were jointly trained using the Adam optimizer (33 with the parameters (B1, (2)) = (0.9, 0.999) and a batch size of 200. The learning rate was initialized 10-3 and decayed to zero with a cosine annealing [36
247
Apache-22.0 License)
httpps///ggaab...c////eettttttte http:///thtaiioooeaaeettttttttta (MIT License)

---

249
Evaluation Metric As an evaluation metric, we used the 1-step error, which is identical to the loss function. We displayed it at the scale of x x10-9 The lower this indicator, the better, as emphasized by 4- While several studies used the MSEs of the state or system energy over the whole time series [26 38), we consider these indicators are misleading, as pointed in several studies 4 For example, in the case of a periodic orbit, an orbit that is correctly learned except for slight difference
250
252
253
254 in angular velocity will have the same MSE as an orbit that never moves from its initial position. Instead, we used the valid prediction time (VPT) 4 32 [54). VPT denotes the time point divided by
296
the length S of time series at which the MSE of the predicted state u ffred exceeds a given threshold for the first time in an initial value problem, that is,
257
V VT U p UGT) = Y/2 arg max., { MMS((U pred. uGr) < 0 for all s S szy.
(12)
258
To obtain VPTs, we normalized each element of state to have zero mean and unit variance in the training data and set 0 to 0.01. The higher this indicator, the better, emphasized by 1. Because of the "spiking" behavior of the FitzHugh-Nagumo model, a small error in phase is regarded as a large error in state. To measure the qualitative performance, we calculated VPTs by allowing for a delay and advance of up to 5 steps.
259
280
261
262
analytical cFINDE
leapfrog dFINDE
263
4.2
First Integral Preservation for Hamiltonian System
1.0
Before learning first integrals from data, we first evaluated FINDE as a numerical integrator using a known mass-spring
state
0.0
285
system. The system has the state u (q v) T, the dynamics d
-1.0
and 4e -q. and the system energy E(q.v) (q2 ( Using the initial value (1.0 0.0) T and the time step size A = 0.2, we solved the initial value problem of the true
0.500
q
energy
288
269
0.495
270
ODE using the leapfrog integrator. We applied FINDE with the
step
50
271
true system energy E' as the first integral V. Note that no neural
Figure 1: Integration of a known mass-spring system.
272
networks nor training were involved.
273
The results with the analytical solution are shown in Fig. | The upper panel shows that the time series predicted by comparison methods overlap each other and are apparently almost identical. However, the lower panel shows that the energy obtained from the states predicted by the leapfrog integrator is fluetuating. The same is true for the case with the cFINDE. This is because the symplectic integrator and the cFINDE suffer from numerical errors caused by temporal discretization. In contrast, the dFINDE preserves the energy accurately. This is because, at every step, the dFINDE projects the state
274
275
276
277
278
(g v) onto the discrete tangent space T( MM' Although a smaller step size reduces numerical errors, this result demonstrates the advantage of dFINDE.
279
280
4.3
Learning First Integrals from Data of Hamiltonian System
We evaluated FINDE learning from the 2-body dataset. We used HNN as the base model f. We found that the FINDE got better performances if it did not treat the Hamiltonian H of the HNN as of first integrals Vi. The medians and standard deviations of 5 trials summarized in the leftmost column of Table 3 The cFINDE achieved better VPT' than the vanilla HNN with K = 1
283
284
285
286
to 2, and the performance was suddenly degraded for K 3. The dFINDE showed similar trend with slightly better performances. The HNN with FINDE found two first integrals in addition to the Hamiltonian H of the HNN. Even though two-body problem is a Hamiltonian system that HNN can learn, the prior knowledge that there exist first integrals other than the Hamiltonian H can be a clue better learning. The HNN with FINDE got worse 1-step errors, suggesting that without FINDE, HNN overfitted short-term change and had difficulty predicting long-tern term dynamics.
287
288
289
290
291
292
We performed symbolic regression of first integrals V learned by the neural network. For K = 2. the learned first integrals V were identical to the linear momenta in the x- and y-directions up to affine transformation in most cases. See Appendix D for more details.
293
294
295
We depict example results in Fig. In the absence of FINDE, the mass positions (21, (I2, V2) became inaccurate a short time and the center-of-gravity position (Iee Be) = ( 21 21 deviated rapidly. The HNN with cFINDE accurately predicted the state for longer period. Even after errors in the mass positions became non-negligible, errors the center-of-gravitt position were still small. We show the absolute errors averaged over all trials in Fig. B In each of r- and y-directions.
296
2
298
7

---

Table 3: Results of FINDE.
2-body HNN
KdV
2-pend
Model
K
1-step]
VPT1
1-step|
VPT
1-step.
VPT1
1-step].
VPT1
base model
5.17
2057
0.362
#0026
5.59 40.00
0.339
20.8s
0.82 #0:2
=00
0.1
=0035
73.66-12.19
0.236 20053
1
7.10 #1.25
0.374
10.066
6.24
10.44
0.371
10.08
0.75
10.04
0.156
20,042
54.18
s.12
0.127
40.148
2
7.78 11:39
0,450
#0052
2.59 20:11
0.608
20.ss
0.73 20:05
0.198
20:088
37.03
=181
0.437
20.04
cFINDE
3
> 10a
0.147
146"
3.19 20:37
0.730
0.69
0.411
20.093
> 109
0.007 100077
10a
0.101
10
3.65 40.30
0.641
10.77
0.77 10:07
0.395
20.083
5
> 10a
0.080
40.014
4.68 +0.43
0.601
+0,009
0,80
+0.07
0,585
+0.097
6
>> 10
0.070
40.019
7.79 40.51
0.425
H0,007
12.53
+0.00
0.005
+0.000
1
7.01 11:06
0.379
10.040
11.61
16.60
0.288
10.083
0.75
10.10
0.152
20.017
47.07
=6.03
0.117 10.122
2
7.03
11.00
0,475
40.022
2.70 :0:26
0.598
10.059
0.74
10,05
0.271
10.111
33.24
=3.40
0,455 =0.022
dFINDE
3
54.78
136.9
0.309
10.024
3.78 10:27
0.636
10.024
0.69
10.05
0,447
10.081
319.70
=9.11
0.049 10.00
10a
0.102
10015
3.48 +0.32
0.780
+0.059
0.71 +0.03
0.454
10:060
5
> 10
0.086 2001
5.26
+0.15
0.718
+0.08
0.86 +0:09
0.591
20.087
6 >> 103
0.059 #0017
9.60 #3,61
0.573 #0.121 58.88 #22.98 0.037
=0:039
A standard deviation follows + symbol. Underlines indicate results better than the base models" results, and bolded fonts indicate the best results. denotes that some trials failed in training because of the underflow of the step size. A dash denotes a case we did not try.
true (Fier (Ne,, Be)
(11. vn
3 272
(122, 1/2)
V/. 3/2
V
0.5
HNN
10a
+cFINDE 10a
ground truth
HNN
+cFINDE
Figure 3: Mean absolute errors of states for the 2-body dataset with or without cFINDE.
Figure 2: Example results of the 2-body dataset.
the HNN without FINDE produced errors in the center-of-gres of-ravty position Te (or ye) and those in the mass positions ø1, I2 (or y1, 3/2) at almost the same level. In contrast, when the cFINDE is present, errors in the center-of-gravie gravity position were much smaller than those in the mass positions, implying that errors in one mass position canceled out errors in the other mass position.
301
303
304
Therefore, we conclude that FINDE not only had better prediction accuracy but also found and
305
preserved linear momenta (which are related to symmetries in space) more accurately despite not having prior knowledge about symmetries.
307
4.4
Learning First Integrals from Data of Unknown Systems
It is often unclear whether target system is a Hamiltonian system or not, but one can expect that the target system has several first integrals. We evaluated FINDE using NODE as the base model. We summarized the results in Table 3
309
311
For the KdV dataset, the NODE with FINDE got much better 1-step errors and VPTs for wide range of K. Figure shows an example result. The top panels show that the prediction results were apparently similar. The bottom panels summarize mean absolute errors in states u. total mass Σg ug. and energy. In the absence of FINDE, the NODE increased all of its errors in proportion to time. With the <FINDE, the error in total mass increased at the point where the two solitons collided but then returned to the original level. Although the calculation is slightly inaccurate, the cFINDE learned
312
313
314
315
316
317
to preserve the total mass. The rightmost panel shows that the error energy continued to increase for K = 2, but it stayed within small range for K = 3. These results suggest that the first or second quantity learned by the cFINDE was total mass, the third quantity was system energy, and the
319
remaining quantity may correspond to one of the many first integrals of the KdV equation.
8

---

50
04
Ground Truth
NODE
+cFINDE K=2
+cFINDE K=
10
Truth
10-
NODE
K = 3
10
in state 10a error in total mass 10
in energy 10
Figure 4: Example results of the KdV dataset. (top) Predicted states. Red belts denote moving solitons. (bottom) Mcan absolute emrors.
Truth
NODE
K=1
2
0-
Ground Truth NODE
+<FINDE
Figure 5: Example results of the 2-pend dataset for 2,000 steps.
2000
2000
Figure 6: Example results of the
dataset.
321
For the 2-pend dataset, the NODE with FINDE got better 1-step errors and VPTs for K = 1 to 5 except for the 1-step error of the dFINDE with K 5. In addition to the system energy, the double pendulum has two holonomic constraints on the position, which lead to two additional constraints involving the velocity (see Appendix [ooor details). Thus, is reasonable that the NODE with FINDE
322
323
324
325
got the best VPTs for K 5 first integrals and totally failed when assuming K >> 5 first integrals. As exemplified in Fig. 5 the NODE without FINDE did not preserve the lengths of rods, making
326
327
the states deviate gradually. See Appendix D.2 for the case when actual constraints are known. For
328
the SeyzHugh-Nagumme dataset, the NODE with FINDE got much better 1-step errors and VPTs for
periodic orbit, and only the reproduced such dynamics. On the other hand, the state did not stay at limited region without FINDE and converged to wrong equilibrium with the cFINDE for
K- 2. As exemplified in Fig. 6 the ground truth state converged
330
NODE with cFINDE for
332
K = 1. For K = 1, the sole quantity V12 may have tried to learn both of the two first integrals and
remained under-trained. In these two cases, FINDE found all first integrals; K 5 for the 2-pend dataset and K = for the SitzHuh-Naaumm dataset.
334
Conclusion
335
336
This study proposed first
neural differential equation (FINDE). FINDE projects the time evolution onto the submanifold defined using the (discrete) gradients of first integrals represented by neural network. With an appropriate number of assumed first integrals, FINDE predicted future states more accurately than base models. Not only that, FINDE found and preserved the system energy and the total mass as first integrals, first integrals related to symmetrics in space, and first integrals led by constraints in unified manner. Therefore, FINDE has the potential to make a scientific discovery by revealing unknown propertics of target dynamical systems.
338
340
342
343
The step errors were on the order of 10-5 to 10-4 in absolute error, being much larger than the numerical error tolerance of 10 used in the experiments; numerical errors were negligible compared to modeling errors. However, the dFINDE tended to get VPTs better than the cFINDE despite the fact that its advantage is to eliminate numerical errors caused by temporal discretization. This result suggests that method leading to smaller numerical errors results in model with smaller modeling emror.. Similar tendencies have been observed previous works 10 38 and these results may form
344
345
346
348
new frontier for integrating numerical and modeling errors.

---

References
350
[1] Abadi, M.. Agarwal, A., Barham, P., Brevdo, E., Chen, Z.. Citro, C.,. Corrado, G. S.,, Davis, A., Dean, J., Devin, M., Ghemawat, S., Goodfellow, I., Harp, A., Irving, G., Isard, M. Jozefowicz. R.. Jia, Y., Kaiser. L., Kudlur, M., Levenberg, J., Mané, D., Schuster, M., Monga, R., Moore, S., Murray, D., Olah, C., Shlens, J. Steiner, B., Sutskever, I., Talwar, K., Tucker, P., Vanhoucke, V., Vasudevan, V., Viégas, F,, Vinyals, O.,
352
355
354
Warden, P., Wattenberg, M., Wicke, M., Yu, Y., and Zheng, X. (2016). TensorFlow: Large-scale machine learning on heterogencous systems. USENIX Symposium on Operating Systems Design and Implementation (OSDI).
356
357
358
[2] Bai, S., Kolter, J. 7., and Koltun, V. (2019). Deep Equilibrium Models. Advances in Neural Information Processing Systems (NeurIPS).
359
360
[3] Barrett, D. G. T. and Dherin, B. (2021). Implicit Gradient Regularization. In International Conference on Learning Representations (ICLR).
381
382
[4] Botev, A., Jaegle, A., Wirnsberger, P.. Hennes, D., and Higgins, L. (2021). Which priors matter? Benchmark- ing models for learning latent dynamics. In Advances in Neural Information Processing Systems (NeurlPS) Track Datasets and Benchmarks.
364
[S] Cao. Y., Fang, Z., Wu, Y., Zhou, D. X., and Gu, Q. (2021). Towards Understanding the Spectral Bias of Deep Learning. International Joint Conference on Artificial Intelligence (IJCAI), pages 2205-2211.
366
387
[6] Celledoni, E., Grimm, V., McLachlan, R.. McLaren, D.. O' O'Neale, D., Owren, B., and Quispel. G. (2012). Preserving energy resp. dissipation in numerical PDEs using the "Average Vector Field" method. Journal of Computational Physicx, (31)0)66777-6778
388
389
370
[7] Chen, S., Billings, S. A., and Grant, P. M. (1990). Non-linear system identification using neural networks. International Journal of Control, 51(6):1191-122 1214.
372
[8] Chen, T. Q., Rubanowa, Y., Bettencourt, J., Duvenaud, D., Chen, R. T. Q., Rubanova, Y., Bettencourt, J., and Duvenaud, D. (2018). Neural Ordinary Differential Equations. In Advances Neural Information Processing
373
374
Systems (NeurIPS), pages 1-19.
375
[9]
Chen, Y., Matsubara, T., and Yaguchi, T. (2021). Neural Symplectic Form: : Learning Hamiltonian Equations on General Coordinate Systems. In Advances in Neural Information Processing Systems (NeurlPS).
376
377
[10] Chen, Z., Zhang, J., Arjovsky, M., and Bottou, L. (2020). Symplectic Recurrent Neural Networks. In International Conference on Learning Representations (ICLR), pages 1-23.
378
379
[1I] Christiansen, S. H., Munthe- Kaas, H. Z., and Owren, B. (2011). Topics in structure praeturr- discretiza- tion. Acta Numerica, 20:1-119.
[12] Clouse, D. S., Giles, C. L., Home, B. G., and Cottrell, G. W. (1997). Time-delay neural networks: representation and induction of finite-state machines. IEEE Transactions Neural Networks, 865):1065-70.
382
[13] Cranmer, M., Greydanus, S., Hoyer, S., Battaglia, P., Spergel, D., and Ho, S. (2020). Lagrangian Neural Networks. In ICLR Deep Differential Equations Workshop, pages 1-9.
384
385
[14] Cuell, C. and Patrick, G. W. (2009). Geometric discrete analogues of tangent bundles and constrained Lagrangian systems. Journal of Geometry and Physics, 59(7)99769997.
386
387
[15] Dahlby, M., Owren, B., and Yaguchi, T. (2011). Preserving multiple first integrals by discrete gradients. Journal of Physics A: Mathematical and Theoretical, 44(30).
388
389
[16] Dehmamy, N., Walters, R., Liu, Y., Wang, D., and Yu, R. (2021). Automatic Symmetry Discovery with Lie Algebra Convolutional Network. In Advances in Neural Information Processing Systems (NeurIPS), number 2018, pages 1-30.
391
[17] Devlin, J., Chang, M.-W., Lee, K., and Toutanova, K. (2018). BERT: Pre-training of Deep Bidirectional Transformers for Language Understanding. arXiv, pages 1-15.
393
394
[18] Dormand, R. and Prince, P. J. (1986). A reconsideration of some embedded Runge-Kutta formulae. Journal of Computational and Applied Mathematics, 15(2)2203-211.
395
396
[19] Finzi, M., Stanton, S., Izmailov, P., and Wilson, A. G. (2020a). Generalizing Convolutional Neural Networks for Equivariance to Lie Groups on Arbitrary Continuous Data. In International Conference on Machine Learning (ICML), pages 3146-3157.
387
10

---

389
[20] Finzi, M., Wang, K. A., and Wilson, A. G. (2020b). Simplifying Hamiltonian and Lagrangian Neural Networks via Explicit Constraints. In Advances in Neural Information Processing Systems (NeurIPS).
[21] Finzi, M., Welling, M., and Wilson, A. G. (2021). A Practical Method for Constructing Equivariant Multilayer Perceptrons for Arbitrary Matrix Groups. In International Conference on Machine Learning (ICML).
403
404
[22] Furihata, D. (2001). A stable and conservative finite difference scheme for the Numerische Mathematik, 87(4):675-66 699.
equation.
405
408
[23] Furihata, D. and Matsuo, T. (2010). Discrete Variational Derivative Method: A Senncture-Preeeringg Numerical Method Partial Disferential Equations. Chapman and Hall/CRC.
407
408
[24] Gear, C. W. (1986). Maintaining Solution Invariants the Numerical Solution ODE SIAM Journal on Scientific and Statistical Computing, 7(3):744-743.
409
410
[25] Gonzalez, O. (1996). Time integration and discrete Hamiltonian systems, Journal of Nonlinear Science, 6(5):4494467
411
412
[26] Greydanus. S., Dzamba, M., and Yosinski, J. (2019). Hamiltonian Neural Networks. In Advances in Neural Information Processing Systems (NeurIPS), pages 1-16.
413
414
[27] Hairer, E., Lubich, C., and Wanner, G. (2006). Geometric Numerical Integration: Structure- Prestrrvnn Algorithms for Ordinary Differential Equations, volume 31 of Springer Series in Computational Mathematics. Springer-Verlie Berlin/lleidelberg.
415
417
[28] He, K., Zhang, X., Ren, S., and Sun, J. (2016). Deep Residual Learning for Image Recognition. In IEFF Conference on Computer Vision and Pattern Recognition (CVPR), pages 1-9.
418
[29] Hong, J., Zhai, S., and Zhang, J. (2011). Disercte Gradient Approach Stochastic Differential Equations with Conserved Quantity. SIAM Journal on Numerical Analysis, 49(5):2017-2
420
[30]
Izhikevich, E. M. and FitzHugh. R. (2006). FitzHughhNaguno model.
422 [31] Jin, P., Zhang. Z., Kevrekidis, I. G., and Karniadakis, G. E. (2020a). Learning Poisson systems and trajectories of autonomous systems via Poisson neural networks. pages 1-12.
[32] Jin, P., Zhu, A., Karniadakis, G. E., and Tang, Y. (2020b). Symplectic networks: Intrinsic structure- preserving networks for identifying Hamiltonian systems. Neural Nenworks, 132:166-179.
425
[33] Kingma, D. P. and Ba. J. (2015). Adam: A Method for Stochastic Optimization. In International Conference on Learning Representations (ICLR), pages 1-15.
427
[34] Lee, K. and Carlberg, K. (2021). Deep Conservation: latent-dynamics model for exact satisfaction of physical conservation laws. In AAAI Conference on Artificial Intelligence (AAAI).
429
430
[35] Levin, A. and Narendra, K. (1995). Recursive identification using feedforward neural networks. International Journal of Control, 61(3):333-547.
431
432
[36] Loshchilov, and Hutter, F. (2017). SGDR: Conference on Learning Representations (ICLR), pages 1-16.
Stochastic gradient descent with warm restarts. International
433
434
[37] Manek, G. and Kolter, (2019). Learning Stable Deep Dynamics Models. In Advances in Neural Information Processing Systems (NeurIPS), pages 1-9.
435
436
[38] Matsubara, T., Ishikawa, A., and Yaguchi, T. (2020), Deep Energy- Based Modeling of Discrete-Time Physics. In Advances in Neural Information Processing Systems (NeurIPS).
437
438
[39] Miura, R. M., Gardner, C. S., and Kruskal, M. D. (1968). Korteweg-de Vries equation and generalizations. II. Existence of conservation laws and constants of motion. Journal of Mathematical Physics, 98):
439
440
[40] Narendra, K. S. and Parthasarathy, K. (1990). Identification and Control of Dynamical Systems Using Neural Networks. IEEE Transactions on Neural Networks, 1(1):4-27 27.
441
442
[41] Nelles, O. (2001). Nonlinear System Identification. Springer Berlin Heidelberg, Berlin, Heidelberg.
443
[42] Paszke, A., Chanan, G., Lin, Z., Gross, S., Yang. E., Antiga, L., and Devito, Z. (2017). Automatic differentiation in PyTorch. In Autodiff Workshop on Advances in Neural Information Processing Systems. pages 1-4.
444

---

448
[43] Potra, F. A. and Yen, J. (1991). Implicit numerical integration for euler-lagrange equations via tangent space parametrization. Mechanics Structures and Machines, 19(1):77-98.
448
[44] Quispel, G. R. and Capel, H. W. (1996). Solving ODEs numerically while preserving a first integral. Physics Letters, Section A: General, Atomic and Solid State Physics, 218(3-6):225 22-22222-222.
450
[45] Quispel, G. R. and Turner, G. S. (1996). Discrete gradient methods for solving ODEs numerically while preserving first integral. Journal of Physics A: Mathematical and General, 29(13).
[46] Raff, I.., Komanduri, R., Hagan, M., and Bukkapamamm S. (2012). Neural Networks in Chemical Reaction Dynamics.
454
[47] Rasp, S., Dueben, P. D., Scher, S., Weyn, J. A., Mouatadid, S., and Thuerey, N. (2020). WeatherBench: A
Benchmark Data Set for Data-Driven Weather Forecasting. Journal of Advances in Modeling Earth Systems, 12(11).
456
[48] Sannai, A., Imaizumi, M., and Kawano, M. (2021). Improved Generalization Bounds of Group Invariant Equivariant Deep Networks via Quotient Feature Spaces. Conference on Uncertuinty Artificial Intelligence (UAI).
458
459
460
[49] Sjöberg. J., Hjalmarsson, H., Ljung, L. (1994). Neural Networks System Identification. IFAC Proceedings Volumes, 2793):35.-38.
481
462
[50] Takeishi, N. and Kawahara, Y. (2020). Learning dynamics models with stable invariant sets. In AAAI Conference on Artificial Intelligence (AAAI).
463
464
[SI] Teshima, T., Tojo, K., Ikeda, M., Ishikawa, I., and Oono, K. (2020). Universal Approximation Property of Neural Ordinary Differential Equations. In NeurlPS Workshop on Disferential Geometry meets Deep Learning (DifGGeo4DL..
485
486
467
[52] Trigo, R. M. and Palutikof, J. P. (1999). Simulation of daily temperatures for climate change scenarios over Portugal: A neural network model approach. Climate Research, 13(1):45-5.
488
469
[53] van der Schaft, A. and Jeltsema, D. (2014). Port-Hamiltonian Systems Theory; An Introductory Overview. Foundations and Trends in Systems and Control, 1()):1333778.
470
471
[54] Vlachas, P. R., Pathak, J., Hunt, B. R., Sapsis, T. P., Girvan, M., Ott, E.,, and Koumoutsakos, P. (2020). Backpropagation algorithms and Reservoir Computing in Recurrent Neural Networks for the forecasting of complex spatiotemporal dynamics. Neural Networks, 126:191-217.
472
474 [55] Wang. Y. J. and Lin, C. T. (1998). Kutte-Kutta neural network for identification of dynamical systems in high accuracy. IEEE Transactions on Neural Networks, 9(2):294 307.
476
[56] Zhong, G. and Marsden, J. E. (1988). Lie-Poisson Hamilton-Jacobi Theory and Lie-Poisson Integrators. Physics Letters A, 133(3):3-8.
[57] Zhong, Y. D., Dey, B., and Chakraborty, A. (2020). Dissipative SymODEN: Encoding Hamiltonian Dynamics with Dissipation and Control into Deep Learning. arXiv, pages 1-6.
12

---

Checklist
480
1. For all authors...
482
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope? [Yes] The theoretical contributions are summarized in Remarks D and 2 The performance improvements were validated numerically in Table and visually in Figs. o
483
484
485
486
(b) Did you describe the limitations of your work? [Yes] We have discussed an increase in computational complexity at the bottoms of Sections and (3.2] We also presented
487
488
the limitations in Appendix D.2 while their situations were originally outside the scope of the proposed method.
(c) Did you discuss any potential negative societal impacts of your work? [N/A] No societal impact is supposed.
491
492
(d) Have you read the ethics review guidelines and ensured that your paper conforms to them? [N/A]
498
We have read the guidelines carefully, but no ethical impact is supposed.
484
2.
If you are including theoretical results...
495
(a) Did you state the full set of assumptions of all theoretical results? [Yes] We have introduced the background of the proposed method and provided the full set of as- sumptions in Section [ Even though a slight modification may make the proposed method available general manifold, we have clearly stated that our theoretical and experimental results were limited to the finiteddimensional Eucleadian spaces.
496
497
498
489
500
(b) Did you include complete proofs of all theoretical results? [Yes] We have provided a proof just before each of Remarks D and 2
501
3. If you ran experiments..
503
(a) Did you include the code, data, and instructions needed to reproduce the main exper- imental results (either in the supplemental material or as a URL)? [Yes] We have enclosed the source code to reproduce all experiments in supplemental materials.
504
505
508
(b) Did you specify all the training details (e.g., data splits, hyperparameters, how they were chosen)? [Yes] We have provided the software and hardware environment, network architectures, and hyperparameters in Section 4.1 We have also provided detailed hyperparameters to generate datasets in Appendix
507
508
509
510
(c) Did you report error bars (e-g., with respect to the random seed after running experi- ments multiple times)? [Yes] We have summarized the standard deviations over five trials in T
512
513
(d) Did you include the total amount of compute and the type of resources used (e.g., type of GPUs, internal cluster, or cloud provider)? [Yes] We have provided the hardware environment in Section 1 but we have anonymized the cloud service providers to avoid a potential violation of the double-blind policy.
515
516
517
4. If you are using existing assets (e.g., code, data, models) or curating/relcasing new assets...
518
(a) If your work uses existing assets, did you cite the creators? [Yes] We have made the source code for all experiments by modifying the source codes of HNN (26 and DGNet [38]. We have cited these references and added links to respective repositories in the footnotes.
520
522
(b) Did you mention the license of the assets? [Yes] We have verified that the source codes of HNN [26) and DGNet [38] are provided in Apache-2.0 License and MIT License, respectively. We have clearly stated these facts in the footnotes.
523
524
525
(c) Did you include any new assets either in the supplemental material or as a URL? [Yes]
526
We have enclosed the source code to generate datasets in supplemental materials.
(d) Did you discuss whether and how consent was obtained from people whose data you're asing/curating? [N/A]
$28
$29
(c) Did you discuss whether the data you are using/curating contains personally identifiable information or offensive content? [N/A]
531
5. If you used crowdsourcing or conducted research with human subjects...
532
(a) Did you include the full text of instructions given to participants and screenshots, applicable? [N/A]
13

---

534
(b) Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals, applicable? [N/A]
535
536
(c) Did you include the estimated hourly wage paid to participants and the total amount
spent participant compensation? [N/A]
14


Reviews:
[{"review_id": "qRc9Er1NHg-", "reviewer": "Reviewer_x8Jk", "summary": "This paper proposes a novel neural model for predicting the dynamics of a system by finding multiple unknown conserved quantities (first integral) from data. The discrete derivative-based formulation eliminates numerical errors, and the model can be trained to strictly satisfy the conservation laws. Experiments on physical systems with multiple conserved quantities show that the proposed model can capture multiple conserved quantities and predict the dynamics of the state more accurately than the base model.", "questions": "Q1. When the number K of conserved quantities is unknown, do you have a way to determine K?\n\nQ2. It would be better to specify the objective function and/or the learning algorithm. I was wondering if it is really appropriate to use 1-step error for learning. In [R1,R2], learning is performed using time-window, which might be useful for this study.\n\nQ3. The evaluation experiments are carefully conducted, but I was concerned that the only comparison method is the base model. If possible, it would be better to have a comparison with at least one of the latest models (e.g., DGNet).\n\nQ4. How robust is the proposed model when varying the number of trajectories and time resolution in the training data? It might work well for sparse data because multiple conserved quantities are considered.\n\n[R1] Kevin Course, Trefor Evans, and Prasanth Nair. Weak form generalized Hamiltonian learning. In Advances in Neural Information Processing Systems, volume 33, pages 18716–18726, 2020.\n\n[R2] Yaofeng Desmond Zhong, Biswadip Dey, and Amit Chakraborty. Symplectic ODE-Net: Learn- ing Hamiltonian dynamics with control. In International Conference on Learning Representa- tions, 2020.", "limitations": "We believe that learning dynamics via a data-driven approach is effective in the following ways: 1) it allows us to construct simulators even when the mathematical model is unknown, and 2) it can lead to scientific discoveries (e.g., conservation laws). In order to clarify the significance of this research, we would like to ask for clarification of possible applications. In particular, please share your thoughts on how this could lead to scientific discoveries. If there is a gap, it would be good to mention it as a limitation.\n", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "- Originality\n  - Pros 1: The problem setting of inferring multiple unknown conserved quantities at the same time is new and interesting.\n  - Pros 2: Inforporating the prior knowledge of the first integral is reasonable.\n  - Pros 3: The formulation based on discrete derivatives is useful.\n\n- Quality\n  - Pros 1: This submission is technically sound.\n  - Pros 2: Experiments with physical systems with multiple conserved quantities have shown the effectiveness of the proposed model.\n  - Cons 1: The description of the learning algorism is somewhat lacking.\n  - Cons 2: Only a comparison with the base model is presented.\n\n- Clarity\n  - The presentation is excellent.\n\n- Significance\n  - Pros 1: Data-driven Approaches to modeling physical systems have attracted considerable attention. It is useful because it has the potential to construct simulators from data for physical phenomena for which mathematical models are unknown.\n  - Cons 1: The application is not stated in detail.", "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": "TCLS6rnLhAQ", "reviewer": "Reviewer_kyoV", "summary": "The authors propose modify the Neural ODE model to conserve the first integrals of the physical system, where the  first integrals are either known or learned by a neural network.", "questions": "- Can the authors clarify the difference between the cFINDE and dFINDE? cFINDE is derived for continuous functions f(u) and V(u), but in practice, it surely uses an integrator to compute the ODE solution. dFINDE also uses a numerical integrator for Psi: “We assume that the base model ψ is composed of a continuous-time base model and a numerical integrator.” Is the only difference between the models in the computation of the grad_V (u)?\n\n- In case first integrals are learned, how do you prevent it from learning degenerate solutions, e.g. predicting zero function for all inputs u? My understanding is that function f(u) is trained together with the learned first integrals V(u). In the best case, the function f(u) can already preserve all or some of the first integrals, and V(u) can only learn a degenerate function.\n\n- In Figure 1 authors demonstrate that leapfrog integrator does not preserve the total energy, while dFINDE does, while dFINDE uses Dormand–Prince method (dopri5)  (mentioned in section 4.1). Presumably, the ability to conserve the integral of dFINDE heavily depends on the integrator used. It is unfair to compare dFINDE equipped with dopri5 versus leapfrog integrator. Can authors provide the results or hypothesise about the following comparisons: 1) leapfrog versus dFINDE with leapfrog integrator 2) dopri5 versus on dFINDE with dopri5 integrator?\n\n- Have authors tested the model on the systems with different masses and total energy than on the training set to verify that the learned constraints indeed learn to compute the  total energy/mass that need to be conserved?\n", "limitations": "I strongly suggest the authors to acknowledge the limitations of their work.\n\nSome limitations worth mentioning in the paper:\n- The ability of the cFINDE and dFINDE to preserve the first integrals depends on the numerical integrators. dFINDE also depends on the time discretisation delta_t. \n- The results in Table 3 suggest that the model is sensitive to K — the number of the first integrals chosen. Moreover, with K>3 the predictions “explode” in case of 2-body in \n- FitzHugh–Nagumo systems.\n- By construction, enforces that certain quantities of the system have to be conserved exactly. The model is not applicable for the systems that are conserving, e.g. having dissipative energy.\n", "rating": 6, "confidence": 5, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "The authors propose a neat way to introduce the conservation constraint by adding the Lagrange multipliers to the derivative function of the Neural ODE. The first integrals are conserved exactly, which, to my knowledge, previous models do not provide. The first integrals can be learned from data, and the authors show that they correspond to the known conservation, such as conservation of mass and energy. It is notoriously hard to preserve the first integrals of the system in the learned physical simulations, and this work provides an elegant way to do so.", "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": "IUo9KR1PrFt", "reviewer": "Reviewer_f3su", "summary": "This paper proposes to learn the first integrals of an ODE system by a neural network and projects\n the time evolution of the dynamical system onto the submanifold defined using the (discrete) gradients of the learned first integrals. With an appropriate number of assumed first integrals, FINDE predicted future states more accurately than base models.", "questions": "1. Can the authors explain more clearly how the implicit dFINDE can be implemented in an explicit way?\n\n2. I don't understand this sentence: \"Not only that, FINDE found and preserved the system energy and the total mass as first integrals\".  Does it mean that FINDE can't guarantee to preserve other types of first integrals?", "limitations": "The limitations have be addressed in the Weaknesses summary.", "rating": 5, "confidence": 2, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths:\n\n1. This paper addresses an important topic in predicting an unknown physical system, that is symmetry, or the associated conservation law;\n\n2. The Discrete-Time Projection Method for preserving the first integral of the discretized dynamical system seems to be novel;\n\n3. The empirical results are through and sound.\n\nWeaknesses:\n\n1. The methodology section is hard to follow, no intuition is given after the formal introduction of the method;\n\n2. The contribution of this article is entangled, the authors should separate the question of finding the right first integrals and making more accurate future predictions of the dynamical system. It's hard to say that making a better prediction comes from learning the right first integrals.\n", "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": "FI-kuxazeKQ", "reviewer": "Reviewer_fysu", "summary": "The objective of the paper is to propose a method that enables to find (and preserve) invariant quantities in the observations of a dynamical system. \n\nThe authors proposed in their work  the joint learning of a derivative ($f$) approximating the dynamics of the observations and parameterized first integrals $V=(V_1,..., V_k)$  that are used between each prediction step. Several technical developments are conducted to enforce the respect of the first integrals.\n\n\n### POST REBUTTAL ###\n\nI would like to thanks the authors for their rebuttal statements. The proposed explanations details and submission modification for the unclear parts are satisfying (notably on the projection and the discretization errors).  Therefore, I maintain my score.\n", "questions": "1) Section 3.1 provides the dynamics that respects the invariants given $\\hat{f}$ and lagrange multipliers. I suggest for clarity to present the optimization problem associated with these $\\lambda$. \n\n2) Can the authors also detail why the continuous version is subject to numerical errors compared to the discrete ones ?\n\n3) What is the method used to solve eq.3 (at each iteration) ? Or is the integration performed using the “modified” (i.e c/d Finde) version of the learned dynamics ?\n\n4) L.144: it is unclear to me why gradients are not computable. Indeed solving eq.3 can be thought of as a meta-learning problem (yet the optimization would not be exact). \n\n5) I get that the whole point of the method is to provide a way to “respect” first integrals, i.e. to “alter” an initial integrator so that $(V_1,..., V_k)$  are constant. However, since it is done pointwise it allows potential spikes in the learned first integral as for instance illustrated in fig.4. Can the authors comment on this and potential implications (for instance it severely arms the interpretability of the learned first integral )? Can the authors think of any additional constraints or inductive biases to palliate such drawback ?\n\n6) What is the influence of the hyper-parameters on the learning mainly on the relative influence of the LR of the networks associated with the first integral ? \n", "limitations": "see questions.", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "S:  The paper and its objective are clearly presented and very relevant to the community and the development of further methods.\nMoreover,  the paper is well-written and the technical contributions of the paper are clear.\n\n\nW: I identify two major weaknesses. \nThe first one lies in the presentation:  the loss objective is not clearly written in the main text which slightly hampers the understanding of the text. Indeed, the authors remain quite far from the application to focus on their technical contribution.Far from being essential, I nonetheless highly recommend the inclusion of an “algorithm” section to detail how the first integrals are learned practically using for instance NN. \n\nThe experiments are conducted using low dimensional data (and it presents a clear case of the potential of the proposition in this case) but the application of the proposition to  higher dimension data is not answered.\n\n", "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: Reject