Constants of motion network
Authors: Muhammad Firmansyah Kasim
Yi Heng Lim
Venue: NeurIPS 2022 Accept (2022)
Keywords: neural ode
learning dynamics
constants of motion
energy conservation

Abstract: The beauty of physics is that there is usually a conserved quantity in an always-changing system, known as the constant of motion. Finding the constant of motion is important in understanding the dynamics of the system, but typically requires mathematical proficiency and manual analytical work. In this paper, we present a neural network that can simultaneously learn the dynamics of the system and the constants of motion from data. By exploiting the discovered constants of motion, it can produce better predictions on dynamics and can work on a wider range of systems than Hamiltonian-based neural networks. In addition, the training progresses of our method can be used as an indication of the number of constants of motion in a system which could be useful in studying a novel physical system.

Constants of motion network
M. F. Kasim & Y. H. Lim Machine Discovery Ltd. Oxford, United Kingdom
{muhammad, yi. hhw.g///acchhh--eddcceeerr...o
Abstract
The beauty of physics is that there is usually a conserved quantity in an always- changing system, known the constant of motion. Finding the constant of motion is important in understanding the dynamics of the system, but typically requires mathematical proficiency and manual analytical work. In this paper, present a neural network that can simultaneously learn the dynamics of the system and the constants of motion from data. By exploiting the discovered constants of motion, it can produce better predictions on dynamics and can work on a wider range of systems than Hamiltonian-bee based neural networks. In addition, the training progresses of our method can be used as an indication of the number of constants of motion system which could be useful in studying novel physical system.
Introduction
Noether's theorem [1] states that if system has symmetry, then there's constant of motion corresponding to it. As it turns out, constants of motion play a significant role in understanding the world around and are ubiquitous in almost every aspect of physics. Among of the most prominent examples are the conservation of energy, the conservation of momentum, as well as the conservation of angular momentum.
Historically, constants of motion were discovered by doing analytical works from observational data or from mathematical descriptions of the systems' equations of motion. For example, the law of conservation of energy was proposed by Châtelet based on Newton's 's work on classical mechanics [2). With the recent emergence of employing neural networks for scientific discovery and learning systems' behaviour from observational data ( 4 3 6) naturally raises a question, "can we find constants of motion dynamical systems from their data and exploit them to make better prediction?"
large body of literature has been moving into this direction by learning the Hamiltonian its variations [8] 9 of system. However, most of the previous works focus on Hamiltonian and its variations, they work well in conserving the Hamiltonian or energy. However, when the systems have other constants of motion, the Hamiltoniian-base works fail to discover and exploit those quantities. Here we present the "Constant Of Motion nETwork" (COMET) that can discover constants of motion of system and exploit them to make better prediction. In contrast to Hamiltonian-basee networks 66. 77 [8] 9], COMET not constrained to Hamiltonian systems and its coordinate choice, making it generally applicable to wider range of systems shown in Table In addition, also found that the training progress of COMET can used as an indication of how many independent constants of motion there are in the system (see section 6 which could be valuable hint in studying a novel physical system. Our implementation and experiments can be found in the public domair
or
36th Conference on Neural Information Processing Systems (NcurIPS 2022).

---

NODE
10
HNN
NSF
LNN
8
COMET
Conserves energy
Works general coordinates
Works on dissipative systems
Conserves other quantities
Table I:
Summary of the methods comparison. The compared methods are neural ODE, Hamiltonian neural network, neural symplectic form, Lagrangian neural network, and COMET (ours).
2
Related works
The simplest way to learn the dynamics of systems with neural networks is by using neural ordinary differential equation (NODE) [T]. NODE takes the full states of the system and produces the dynamics of the systems, i.e. the time derivative of the states. The simulation can then be run by solving the ODE from the output of NODE. As there is no inductive bias in the NODE, they typically struggle to conserve quantities that are important in some systems' dynamics, such as energy.
Hamiltonian neural network (HNN) 6 is attempt solve this conservation problem by learning the Hamiltonian and calculate the state dynamics from the Hamiltonian. has been shown that with HNN, one conserve the energy and produce better motion prediction in long time horizon. Due to its simplicity and the elegance of the idea, HNN has been applied on a wide range of tasks and neural network architectures 12 [3], and even on dissipative systems by adding dissipation term 14 [5]. HNN can also be combined with symplectic integrator [6. to produce a better result from the trajectory observations.
Despite its ability to conserve the energy, HNN is limited by the requirement of using canonical coordinates instead of arbitrary coordinates. Works on Lagrangian neural networks [18] 8 solve this limitation by learning the Lagrangian. Other attempts use neural symplectic form to learn the coordinate-free representation of Hamiltonian 7 or Poisson neural network learn the Poisson However, those works are still limited to Hamiltonian or Poisson systems.
system 9
3 COMET: constants of motion network
We start by denoting a set of states system as s R= where n,, the number of states. States are the internal parameters of a system that completely determines its dynamics without external influence. For example, in a classical particle motion, the particle's position and velocity constitute the states of the system. Without external influence, the change of the states typically depends on the states itself, i.e. ds/dt = $(s).
A constant of motion is a quantity that is conserved over the time in the system, like energy. In some systems, such integrable systems [19), there other quantities other than energy that conserved, for example, momentum or angular momentum. These constants of motion can typically be described as a function of the states of the system, so we denote it as c(s) E 18 with n. is the number of constants of motion. As their quantity is constant throughout the motion, their time derivative must be O, or de/dt- 0. By taking the dependency of c on s, the condition on can be written as
dc dc $ 0 0.
dt Os
where 0c/Os is n. X n, Jacobian matrix where each row of the matrix is the gradient of each constant of motion with respect to the states s. The equation above means that the state dynamics must be perpendicular to the gradient of each constant of motion.
To design a deep learning architecture that can simultaneously learn the constants of motion and learn the dynamics that conserve the constant of motion, we define two functions that depends on the states that can constructed with neural networks, $o and c(s). The function &n R + R the initial guess of the rate of change of the states. The function Th - Rd is the function that computes the constants of motion of the system. To ensure the constants of motion are conserved as in equation]]
compute the state dynamics by orthogonalizing the initial guess šo with respect to

---

the gradient of every constant of motion.
ortho (šo. Vcc. Ve2, .... 1),,)),
where e, is the i-th element of the constants of motion and ortho(a, V) is operation to orthogo- nalize the vector a to every vector in the set V.
3.1 Orthogonalization process
One way to produce an orthogonal vector against a set vectors by using QR decomposition, i.e.
A (Vei, Ve2, VCn. 1 $o)
Q, R = QR(A)
$ (c.ne) (. (ne,ma))
(3)
where Q(...)) is the last column of the matrix Q and ((.....) is the element at the last row and last column of the matrix R. The first row of the equations above shows construction of tall matrix A I x()+t)) where the first ne columns
the gradient of the constants of motion and the last column is the initial guess of the states rate of change, $p. QR decomposition is usually implemented using Householder transformation 20 which produces much smaller numerical error than the alternative Gram-Schmidt process 21 in practice.
The QR procedure above imposes constraint that the number of constants of motion must be less than the number of states, i.e. nhe f2.- This is agreement with the maximum number of independent constants of motion in an integrable system is ne n, 1 [19]. By using QR decomposition, COMET will try to find ne independent constants of motion. Reasoning behind this is given in Appendix ??.
3.2 Training loss function
We need train the trainable functions in COMET, (s) and c(s), that the state dynamics & from equation match the dynamics from the observation or training data, &. order to train COMET, the loss function in this case is constructed as
The
S0
W2
Vc, $oll?
(4)
where w. are the tunable regularization weights. The first term of the loss function is the standard Le error where the prediction must match the training data. The second term of the equation above is included to accelerate the training process by making the initial guess šo to be as close as possible to the actual value of the states' rate of change. The third term an additional regularization help the discovery of the constants of motion.
Learning constants of motion from data
To demonstrate the capability of COMET to simultancously learn both the dynamics and the constants of motion, we tested it in a variety of cases. For all the cases in this section, the training data were generated by simulating the dynamics of the system from f collected the states were added a Gaussian noise with standard deviation a = 0.05.
t 10. From the simulations, well the states rate of change, 3, which were calculated analytically and
There are 6 simple experiments performed to demonstrate the capability of COMET: (1) mass-spring, (2) 2D pendulum, (3) damped pendulum, (4) two body. (5) nonlinear spring, and (6) Lotka- Voltera The cases were selected to represent wide variety of cases. It includes cases with Hamiltonian in canonical coordinates (case 1, 4, 5), Hamiltonian with non-canonical coordinates (case 2, 6), a case with redundant states (case 2), dissipative system (case 3), and case with a moderate number of states (case 4). The details each case, including the number of constants of motion that we set for COMET training, well as the training setups described in appendix ??.
For each case, we compared the performance COMET with other methods: (1) simple neural ODE (NODE) [10] (2) Hamiltonian neural network (HNN) [6] with the coordinates given in each case below, (3) neural symplectic form (NSF) (7), and (4) Lagrangian neural network (LNN) [8. The
neural network architecture for each method is detailed in appendix ??.

---

Case
NODE [10
HNN 6
NSF 7
LNN
COMET
mass-spring
0.17
0.19
0.11
0.2210.13
0.17
0.12
0
0.09
0.10+0.16
0.09
2D pendulum
0.087
0.10
0.11 0.11 +.024
0.029 -0018 00011
0.18+0.37
-0.14
damped pendulum
-..4
110
fail
fail
0.007 +0.014
0,005
two body
460 460
0.49-6.0
+-400
fail
fail
0.42 +0.488
6.33
0..99
nonlinear spring
0.63 0..3
0.18+0.74
0.19
0.17 -0.14
0.23 0.23+0.44 +2300
0.33
0.11
0.15
0.18
Lotka- Volterra
0.12+0.36
0.65
0.080
+0.00
N/A
0.048+0.00 0.0411
0.59
0.071
Table 2:
Root mean squared error of 100 randomly initialized simulations for each case and each method. The main number is the median while the range represents the 95% percentile (i.e. lower and upper bounds are 2.5% and 97.5% percentiles, respectively). The bolded values the ones that give the best upper bound among other methods, while the italic values denote the second best. "fail" means that there integration failures with scipy's 's solve_ which makes it unable to integrate to = 100 in a reasonable time.
4.1 Results
For each case, we tested each method by running an- other 100 simulations from 0 to f 100 (10 times longer than the training) using 1000 sampled points with the initial condition randomly initialized above using different seed. The root mean squared errors of the state predictions are shown in table 2
Mass-spring (COMET)
Mass-spring (Thue)
0.4
velcity
From table we can see that COMET performs well across the test cases. In the mass-spring case, all meth- ods perform well. However, when it goes to the pen- dulum cases and Lotka-Volterra, HNN fails to predict the dynamics due to the chosen coordinates not being the canonical coordinates. Although NSF can perform reasonably well in 2D pendulum, it fails on the damped pendulum
0.4
0.2
3.0 Postion
0.2
0.4
0.4
022
0.0
0.2
04
Fostion
Lotka-Voiterra (COMET)
Lotka-Volterra (True)
1.75
Iredats
1.25
0.75
because it is not Hamiltonian sys- tem. COMET takes the advantage having constants of motion in the cases above that be exploited to guide the true trajectory, therefore, it can achieve better predictions of the dynamics regardless of the chosen coordinates and whether it conserves energy or not.
0.75 30
1.25 1.50 1.50 Prey
0.75 1.00 122 1.0 1.75
Prey
Figure 1: The contour plot of constant of motion discovered by COMET (left) com- pared to the true constant of motion (right) for mass-spring (top) and Lotka-Von Volterra (bottom) cases.
Figure D shows the discovered constant of motion in mass-spring and Lotka-Volterra cases. As can be seen from the figure, COMET can successfully discovered the constants of motion from the data. Figure 2 shows the evolution of the known constants of motion for every method in the mass-spring, 2D pendulum, and the two body cases. The periodic variation from the true constants of motion due to the added noise in the training data. In mass-spring case, figure shows that the HNN, NSF, and COMET conserves the energy while the NODE gets the energy decreasing
time.
different story can be found in figures bbb) where they show constants motion of pendulum in 2D coordinate. In this case. NODE and HNN fail to conserve the constants of motion, while COMET can conserve those constants of motion during long period time. The failure of HNN can be attributed to the state coordinates not being the canonical coordinates. This shows that COMET can discover the constants of motion with much less constraints the state coordinates than HNN.
For the two body case in Figure C)c,, we can see that NODE and NSF diverges quite quickly. The failure of NSF this case might be due to the added noise in the training data. Among the tested methods, only HNN and COMET conserves the energy. However, HNN is only designed for Hamiltonian or energy conservation, it fails to conserve other quantities, such as the momentum and angular momentum. COMET, the other hand, can successfully conserve those quantities.

---

Energy springo
Energy
Length
angle pends
1.0
0.10
1002
N00C
0.66
0.C
0.00
20
40
60
10c
20
40
80
100
20
60
80
100
0
60
103
Tlme
Time
Time
Tlmm
(a)
(b)
Energy tWN body!
x.momentum two body)
momeenuum ltwo body!
Angular morrettum tum rowp body?
1.0
0.4
02
0.
-0.2
0.6
-0.3
0.5
(c)
Figure 2: The constants of motion calculated for every method for (a) mass-spring, (b) 2D pendulum, and (c) two body. Please note that the integration of NSF and LNN for the two body case cannot be completed.
HAM
COMET
Figure 3: Motion trajectory of the simulated two body system using Neural ODE, HNN, NSF, and 20.
CÓMET (ours) from
Figure shows the trajectory of the two body system simulated using various methods tested here from
20. From the figure, we can see that only our method (COMET) that can produce closed trajectory. HNN produces almost closed trajectory, but it slightly deviates from the closed trajectory because it conserves only the energy, but does not necessarily conserve the other quantities. By exploiting many constants of motion possible, COMET reproduce the motion with the small error compared to the other methods.
Systems with external influences
One advantage of COMET is that it can easily work with systems with external influences, such as external forces. If the system has conserved quantities when the external influences are kept constant, then COMET with simple modification exploit them to get more accurate dynamics. The modification is just to make the initial guess of the dynamics and the constants of motion to depend the external influences well the states, iee. So(s, x) and c(s, x), where x E R is the external influence. The dynamics can still be calculated following the equation 2
be used learn the constants of motion and
We conducted an experiment using the 2D pendulum from the section 4 but with additional external force in the x-direction, Fz. The training data generated by having the external force with profile (t) = ag cos( a + a2) with uniformly-distributed random values of ag ~ 7 -0.5, 0.5), a1 ~ (0, 5). and a₂ ~ ((, 2x). The experiment was done similarly like in section 4 by adding the force as the input to the neural network for NODE, HNN, NSF, and LNN as well.
Figure&shord shows the constants of motion on the test system with constant external force. As seen on the figure, the values of the true constants motion produced by COMET is oscillating slightly around a constant offset, due to the added noise in the training. In contrast, NODE produces the shift of the energy values, NSF produces large oscillation even for the energy, LNN quickly diverges. Although HNN produce similar energy deviation with COMET, it has larger deviation other
5

---

Energy (forced 2d perdulum)
Length (forced 2d pendu lam)
Velocily angle iforced 23 pendulum)
0.50
0.25
Tue
MODE
1.4
0.4
12-
0.2
UNN
1.0
-0.25
CONET
0.8
-0.2
-0.50
0.6
-0.75
0.4
-0.6
40
60
100
40
60
100
40
100
Time
Time
Time
Figure 4: Constants of motion of the forced 2D pendulum case calculated using NODE, HNN, NSF, LNN, and COMET.
conserved quantities. This shows that even with external forces, COMET can still find and exploit the constants of motion for its dynamic predictions.
6 Finding the number of constants of motion
For most systems, the number of independent constants of motion is usually not known beforehand and not so obvious. Knowing the number of constants of motion can be useful in understanding the manifold dimension of the motion, however, this problem is not an easy problem to solve. During the research of this work, we observe that COMET's training progresses might provide a valuable indication to the number of constants of motion.
The parameter to look at is the first term in the loss function in equation 4 i.e
L
(5)
If we set the number of constants of motion greater than the true number, then that term could not get lower than certain value. is because & is constrained to be perpendicular to the constants of motion and there are excessive constants of motion, then it may not be able to match the value from experiments a certain accuracy.
We ran simple experiment to find the number of constants motion for known systems. Specifically, COMET was trained in the damped pendulum, two body, and 2D nonlinear spring cases from section 4 without added noise and ran for 3000 epochs. Those cases are known to have 2, 7, and 2 constants of motion out of 4, 8, and number of states, respectively. The numbers of constants of motion were scanned from 0 to the maximum value, n,, 1. Figure shows the value of C for the various cases with varying number of constants of motion.
From the figure 5 (top), we can see that once the number of constants of motion is set above a certain number, the value of CJ suddenly increases. This gives indication of the actual number of constants of motion. the system has the maximum number of constants motion, then the values of C1 will always be similar to the values with n. 0. Besides the final value of L1, the evolution value of L1 for various number of constants of motion can be an indication on the true number the constants of motion
in figure (5 (bottom). From figure 5 we see that the number of constants of motion for the damped pendulum case is 2, for the two body case is (the maximum number), and for the nonlinear spring case is 2.
6.1 Failure mode
This technique in determining the number of constants of motion depends on the ability of the neural network to find the constants of motion. Therefore, if the neural network is not expressive enough, it could fail to find the constants of motion and indicate lower number of constants of motion than should be.
Figure illlstrrasss this case where we scanned the number of constants of motion from 0 to the 2D nonlinear spring case where the neural network only has 50 hidden elements per layer instead of 250. It gives an indication that the number of constants of motion be instead of the true number
of 2.

---

Damped pendulum
Reletive mean values of L1
Two budy
Nonlinear spring (2D)
Reletive mean values of
10¹
Relative values
10
13a
10a
10a
10
10f
0.5
1.0
2.0
3.0
Number of CoMs
Number of CoMs
Nurbber of CoMs
Damped pendulum
Two body
Nonlinear spring (2D)
10-2
10-
--0
L1 value
10-1-
Fg
value
20
m=
10-3
0
0 =1
10"4.
10-1
10-1
200000 400000 600000
23000C 400000 600000
200200 400003 600000 Tra ninng step
Training
Tra Training step
Figure 5: (Top row) The relative mean values of L1 from equation for damped pendulum, two body, and nonlinear spring cases, with number of constants of motion The were scanned from 0 to n, 1. The relative values were calculated by dividing it by the value of and the error bars were respectively obtained by taking the mean and std from 5 COMETs trained with different random seeds. (Bottom row) The values of L1 during the training for various numbers of constants of motion for damped pendulum, two body. and nonlinear spring cases.
nhe O. The values
7 More complex cases
Simulating a system with infinite number of states - The previous examples only involve systems with finite countable number of states. To demonstrate the general applicability of COMET,
Nonlinear spring (2D)
10-1
valua
experiment simulation of systems with infinite (but discretized) number of states. Specifically, we trained the COMET learn the dynamics of shallow wave following Korteweg- De Vries (KdV) equation [22] 23] of u(z, t), øi
du
10000r20000n 300000 300000000000000000000000
52 8h The states in this case are the values of u along the a-axis which constitutes infinite number of states.
u Zv. Du
Training step
Figure 6: The failure of finding the number of constants of mo- tion using smaller network.
In our experiment, we simulate the behaviour of from x = 0 to L 5 with periodic boundary condition, sampled in 100 points with uniform spacing. We also set 8 = 0.00022 for numerical stability. The training dataset was generated by running 100 simulations with random initial condition from = 0 to 10 with 100 steps. The initial condition in the training dataset is u(x, 0) ag + dn cos (2 L + a2) where ao 01 and a2 randomly chosen within the range of (1.5, 2.5). (0, 1), and (0, 2Tr, respectively.
The neural network padding, followed by logsigmoid activation function. The pattern above was repeated times but without the activation function for the last one, using 250 channels in the hidden layers. The training was done as described section which takes about 5-7 hours on an NVIDIA T4 GPU, The number of channels in the input is 1 (only for u), and for the output it is 1 + ne where n. is the number of constants motion that we set. The first channel of the output is to represent the initial guess of the dynamics, (z).. The last n. channels are for what call as the constants of motion density, p; (z).
construeted with ID convolutional layers with kernel size 5 and circular
for = 1,..., Re. From pr(x), the constants of motion can be calculated as c; = Je2 (x(x) dx. Using the outputs from the network, the dynamics
be calculated following the equation 2
compared the performance NODE and COMET in solving the KdV equation for ( = 0 to 20. It is not obvious how to apply HNN, NSF, and LNN as the KdV equation has only as its states and not include velocity nor momentum. Figure shows the states u(x, 4) simulations using the true dynamics, NODE, and COMET. From the figure, can see that at f 20.
of the
7

---

4.0
4.0
4.0
COmS)
4.0
coms)
3.5
3.5-
35-
3.5
3.0
3.0
2.
2.5
2.0
1.5
1.5
1.0
1.0
1.0
10
1
2
3
Figure 7: Plot of u(x, t) at f 20 from simulations done by the true analytic expression, NODE, and COMET using 1 and 2 constants of motion. All simulations were initialized to the same initial condition. The simulation run by NODE already diverges & = 20.
the simulation done by NODE has diverged while COMET simulations are still intact. This shows that COMET can take advantage of the constants of motion to make its prediction more accurate.
Learning from pixels - Another experiment that we run to show COMET's capability is to learn the dynamics from pixels and see if it can exploit the constants of motion in the latent space of an auto-encoded In this case, we simulate the dynamics of boobbdy with gravitational interactions (like in section D and present the data as 30 30 pixels images that show the positions of the 2 objects.
Our model consists of an encoder, a COMET, and a decoder. The inputs to our model are 2 consecutive frames to capture the information about positions and velocities which we denote x. 'he encoder converts the pixel inputs to latent variables, s = fx().. The dynamics of the latent variables are then learned by COMET which will then be converted back to pixel images, X (s). by the decoder.
The loss function in this case is simply a sum of the COMET loss from equation ?? and the auto- encoder loss, i.e. Lae = ||x - x|. The observed state derivatives to compute the COMET loss in equation ?? is approximated by finite difference of the encoded pixel data from 2 consecutive time steps, i.e. & = (x + At)) f (x(t))) /tt. In contrast to HNN and LNN, there is no requirement in COMET to make half of the latent states to be the time derivative of the other half.
In this experiment, we compared the performance of COMET, NODE, HNN, and NSF to learn the dynamics of the latent variables of the auto-encoder. The auto-encoder architectures are kept to be the same for all methods. The number of latent states are picked arbitrarily to be 10. It is more than the true number of states which is 8. For COMET, the number of constants of motion for COMET was found by following the procedure in section 6 which is 9.
Figure shows a sample of the decoded images of the dynamics predicted by COMET, NODE, HNN, and NSF, compared to the ground truth. From the figure, we can see that the dynamics predicted by NODE diverges as soon as 20 and at some point the decoder produces non-sensible images. Although NSF and HNN can produce good images until the end of simulation at = 180, the dynamics predicted by NSF and HNN diverges from the ground truth. contrast, COMET can still match the dynamics of the ground truth simulation until the end of simulation. This shows COMET's capability in accurately learn the dynamics of latent states of an auto-encoder.
Discussions
Limitations - COMET works better if there is at least one constant of motion in the system. If there
constant of motion, then COMET works similarly like the neural ODE [10] Although presented a way to find out the number of constants of motion in section it still requires multiple
training processes and manual insight.
In the case of successful training, COMET sometimes produces dynamics that are stiffer than the true dynamics, although LNN and NSF more often produce stiffer dynamics. In a rare case, the dynamics from COMET are so stiff that the integration by scipy's 's solve_ivy cannot be completed in reasonable time. This only happens in the KdV case and did not happen in the other cases tested for this paper. We believe that the limitations above should be addressed to move forward.
Broader impact- - The impact of deep learning on physical sciences is expected to be similar to the impact of scientific computational method. Although it enables new fields of study, adds more point of failure. For example, if new unexpected result is discovered using deep learning methods, it
8

---

I- 0.0
E- 20.0
:- 45.0
I- 60.6
L- 80.0
I- 100.0
t- 1200
140.0
2 260.0
180
Croune frutts
CONET
N200
H2N
REF
Figure 8: Snapshots of the decoded images of the dynamics predicted by COMET, NODE, HNN,
and NSF, as well as the ground truth images.
could be true discovery or false discovery due to the failure/imperfection in training, unsuitable neural network architecture, bug in the code, among other things. Therefore, deep learning methods such COMET should be accompanied with other different and independent methods to confirm obtained results in scientific works.
Conclusions
We have shown that COMET can simultaneously learn the constants of motion and the dynamics of system from observational data. Because the assumption made by COMET (i.e. have constants of motion) is less strict than Hamiltonian-basee neural networks, it can be applied to a wider range of systems than the Hamiltonian-bast external influences. The training progresses of COMET can also give an indication on the number of constants motion in system. With all the advantages we presented, we believe that COMET can be valuable tool for scientific machine learning in the future.
neural networks, including dissipative systems and systems with
Acknowledgement
We would like to thank Sam Vinko and Brett Larder for their helpful comments improving the manuscript as well as Ayesha Chairannisa for proofreading the paper.
References
[1] Emmy Noether. Invariant variation problems. Transport theory and statistical physicx, 163):86-207,
1971.
[2] Ruth Hagengruber. Emilie du Châtelet between Leibniz and Newton, Springer, 2012.
John Jumper, Richard Evans, Alexander Pritzel, Tim Green, Michael Figumov, Olaf Ronneberger, Kathryn Tunyasuvunakool, Russ Bates, Augustin Zidek, Anna Potapenko, al. Highly accurate protein structure prediction with alphafold. Nature, 596(7873):583 2021.
[4] Li Li, Stephan Hoyer, Ryan Pederson, Ruoxi Sun, Ekin Cubuk, Patrick Riley, Kieron Burke, et Kohn-sham equations as regularizer: Building prior knowledge into machine-leamed physics. Physical
review letters, 126(3):036401. 2021.
Muhammad Kasim and Sam Vinko. Learning the exchange-correlation functional from nature with fully differentiable density functional theory. Physical Review Letters, 127(12):126403, 2021.
[6] Samuel Greydanus, Misko Dzamba, and Jason Yosinski. Hamiltonian neural networks. Advances in Neural Information Processing Systems, 32, 2019.

---

[7] Yuhan Chen, Takashi Matsubara, and Takaharu Yaguchi. Neural symplectic form: Learning hamiltonian equations on general coordinate systems. Advances in Neural Information Processing Systems, 34, 2021.
[8] Miles Cranmer, Sam Greydanus, Stephan Hoyer. Peter Battaglia, David Spergel, and Shirley Ho. La- grangian neural networks. arXiv preprint arXiv:2003.046330 2020.
Pengzhan Jin, Zhen Zhang. Joannis G Kevrekidis, and George Em Karniadakis. Learning poisson systems and trajectories of autonomous systems via poisson neural networks. IEEE Transactions on Neural Networks and Learning Systems, 2022.
[10] Ricky TQ Chen, Yulia Rubanova, Jesse Bettencourt, and David K Duvenaud. Neural ordinary differential equations. Advances in neural information processing systems, 31, 2018.
[11] Alvaro
Victor Bapst, Kyle Cranmer, and Peter Battaglia. Hamiltonian graph networks with ode integrators. arXiv preprint 1rX9v:1909. 127900 2019.
[12] Peter Toth, Danilo Jimenez Rezende, Andrew Jaegle, Sébastien Racanière, Alcksandar Botev, and Irina Higgins. Hamiltonian generative networks. arXiv preprint arXiv: 1909. 13789, 2019.
[13] Chen-Di Han, Bryan Glaz, Mulugeta Haile, and Ying-Cheng Lai. Adaptable hamiltonian neural networks. Physical Review Research, 3(2):22316,, 2021.
[14] Sam Greydanus and Andrew Sosanya. Dissipative hamiltonian neural networks: Learning dissipative and conservative dynamics separately. arXiv preprint arXiv:2201.100S3 2022.
[15] Yaofeng Desmond Zhong, Biswadip Dey, and Amit Chakraborty. Dissipative symoden: Encoding hamiltonian dynamics with dissipation and control into deep learning. arXiv preprint arXiv:2002. 08860, 2020.
[16] Zhengdao Chen, Jianyu 7hang, Martin Arjovsky, and Léon Bottou. Symplectic recurrent neural networks. In International Conference on Learning Representations, 2019.
[17] Yaofeng Desmond Zhong, Biswadip Dey, and Amit Chakraborty. Symplectic ode-net: Learning hamiltonian dynamics with control. In International Conference on Learning Representations, 2019.
[18] Michael Lutter, Christian Ritter, and Jan Peters. Deep lagrangian networks: Using physics as model prior for deep learning. arXiv preprint arXi::1907.04490 2019.
[19] Nigel Hitchin, Graeme B Segal, and Richard Samuel Ward. Integruble systems: Twistors. loop groups. and Riemann surfaces, volume 4. OUP Oxford, 2013.
[20] Alston S Houscholder. Unitary triangularization of nonsymmetric matrix. Journal of the ACM (JACM), 5(4):339-34) 342, 1958.
[21] Ward Cheney and David Kincaid. Linear algebra: Theory and applications. The Australian Mathematical Society, 110:544-550, 2009.
[22] Olivier Darrigol. Worlds flow: A history of hydrodynamics from the Bernoullis to Prandil. Oxford University Press, 2005.
[23] Norman J Zabusky and Martin D Kruskal. Interaction of" solitons" a collisionless plasma and the recurrence of initial states. Physical review letters, 15(6):240, 1965.
Checklist
1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope? [Yes]
(b) Did you deseribe the limitations of your work? [Yes] See sections 6.1 and [Yes] See section
(c) Did you discuss any potential negative societal impacts your work?
(d) Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yes]
2. If you are including theoretical results...
(a) Did you state the full set of assumptions of all theoretical results? [N/A]
(b) Did you include complete proofs of all theoretical results? [N/A]
3. If you experiments...
(a) Did you include the code, data, and instructions needed to reproduce the main experi- mental results (either in the supplemental material or as URL)? [Yes] See the footnote
at the end of section
10

---

(b) Did you specify all the training details (e.g., data splits, hyperparameters, how they were chosen)? [Yes] See section
(c) Did you report error bars (e.g., with respect to the random seed after running experi-
ments multiple times)? [Yes] See table and figure S
(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] See section 4
4. If you are using existing assets (c.g., code, data, models) or curating/relcusing new assets...
(a) your work uses existing assets, did you cite the creators? [N/A]
(b) Did you mention the license of the assets? [N/A]
(c) Did you include any new assets either in the supplemental material or as a URL? [N/A]
(d) Did you discuss whether and how consent was obtained from people whose data you' re using/curating? [N/A]
(e) Did you discuss whether the data you using/curating contains personally identifiable information or offensive content? [N/A]
5. If you used crowdsourcing or conducted research with human subjects...
(a) Did you include the full text of instructions given participants and screenshots, applicable? [N/A]
(b) Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals, applicable? [N/A]
(c) Did you include the estimated hourly wage paid to participants and the total amount spent participant compensation? [N/A]


Reviews:
[{"review_id": "koMUAhowAP", "reviewer": "Reviewer_yjkt", "summary": "In this paper, the authors proposed a novel method to solve dynamical systems through data. Instead of formulating the system Lagrangians or Hamiltonians, the method leverages the symmetry and conservation laws, thus does not assume any specific dynamical system forms. The proposed network predicts both the guess of the generalized velocity and the set of constants of the motions (COMs). Then, by exploring the orthogonality between the gradients of COMs and the generalized velocity, the final system velocity is projected through a QR decomposition.  The authors show that their method performs equally well or better than many baselines in many simulated dynamical systems. \n", "questions": "Besides the main comment, I have a few questions and suggestions for the authors:\n\n1, Explain a bit more on why QR decomposition is chosen as the way to orthogonalization. What other methods are available? Why this particular one (i.e. differentiability, computation simplicity, etc)? \n2, Explain a bit more on the “stiffness” of the dynamics. Why COMET, LNN, NST can produce stiffer dynamics?\n\n3, Figure 5 shows that COMET works even without all COMs included. Can the authors add an ablation study to show how well COMET performs with a partial set of COMs?", "limitations": "N/A", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "The strengths of the paper:\n\nThe method is quite novel yet simple, i.e. utilizes the system’s symmetry structure and conservation laws to better learn the system dynamics. \n\nThe proposed method can handle many different types of dynamical systems, including Hamiltonian systems, dissipative systems, and even infinite dimension PDEs. \n\nThe method outperforms many baselines in a wide range of tasks, and the generated trajectory has higher quality in some tasks (i.e. two body systems). \n\n\nThe method works well even with a partial list of constants of motions, as shown in Figure 5. \n\n\nThe weakness of the paper:\n\nMost of the experiments are conducted in relatively simple dynamical systems. It will be interesting to compare the method with Lagrangian and Hamiltonian based methods in higher dimensional mechanical systems (such as robots).  I suspect that LNN may outperform the proposed method in some of these cases. \n\nIn table 1, LNN is shown to not work with dissipative systems, which I am not convinced of.  It will be surprising that a Lagrangian based method cannot handle dissipations.\n\nThe proposed method, COMET, performs significantly worse than other baselines such as LNN in the 2D pendulum task, and even worse than NODE, which is surprising–According to the authors, if the system does not have any symmetry it will be the same as NODE. ", "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": "MopwIkWwmeS", "reviewer": "Reviewer_FmUK", "summary": "The paper is well written and easy to follow. It describes a methodology to learn a system of differential equations from data. The methodology can be used to capture the evolution of a given multi-dimensional time series s(t); possibly the time series can be influenced by an exogenous input (e.g. a force acting on the system). The working hypothesis is that the provided time series is the result of a constrained differential equation, i.e. there exists a function c(.) such that c(s) is constant on the trajectory s(t). \n\nThe proposed approach is tested on motion data (i.e. frictionless mass-spring, 2D pendulum, 2D damped pendulum, two body interactions, 2D nonlinear spring, Lotka-Volterra dynamics). Results show that the proposed method performs better than previous approaches (i.e. NODE, HNN, NSF, LNN) in terms of mean squared error on prediction. Some heuristics for choosing the number of constraints acting on the system are proposed; generalization to systems with infinite number of states is also proposed. ", "questions": "- ***How robust is the proposed methodology to the choice of the system state?*** Do the authors have evidence that the proposed approach could potentially work replacing the given system state with an embedding (such as the latent of a VAEs) learnt from observations which are connected to the state but possibly higher dimensional (e.g. an image of the pendulum)?\n\n- ***How robust is the proposed methodology to noise on the system derivative?*** The system state is always composed of the system generalized coordinates and their derivatives. As a result the system state derivative often contains the generalized coordinates accelerations (i.e. second order derivatives) which are extremely noisy when measured on real systems. Authors mentioned the states rate change was subject to a Gaussian noise with standard deviation of 0.05. Can authors elaborate on how this is representative of the typical noise observed in real systems when second order derivatives are approximated (e.g. finite differences)? \n\n- ***What is the role of $s_0(s)$ at test time?*** Looking at (4) and (2), it seems that the only role of $s_0$ is (at training time) to guide $s$ to $\\hat{s}$ and to be sure that at convergence $\\dot{s} \\in \\langle \\nabla c_1, \\dots, \\nabla c_n \\rangle^\\perp$. After training it is unclear whether the function $s_0$ is needed. Can author describe better the need to explicitly learn the function $s_0(s)$?", "limitations": "Not applicable. ", "rating": 3, "confidence": 4, "soundness": 2, "presentation": 3, "contribution": 1, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "The paper is well written and easy to read. The structure of the paper is sound and results have been presented with enough details to understand the main point of the submitted paper. The major but fundamental weakness of the paper is that it focuses on relatively simple and artificial examples (e.g. mass-spring, pendulum) which aren't of practical use but merely useful for assessing the soundness of the approach. Even though I can see the value of these examples when explaining Lagrangian or Hamiltonian physics (where we need to understand how the different analytical components interact), I see very limited value in using these artificial examples when explaining how to use neural-networks to approximate a given time series. By focusing on this artificial examples, authors fall short in addressing fundamental questions listed below.\n\n(1) ***Choosing the system state***. In classical mechanics  (e.g. Newton-Euler, Legrangian, Hamiltonian) the state of the system and its derived quantities (e.g. its derivative, the momentum) are fundamental quantities with specific properties. The choice of the state isn't trivial and it requires an analytical understanding of the system itself. In the proposed paper, authors assume that the state is given and often coincides with the one used in classical mechanics. However, in applications relevant to machine learning (as the one proposed), assuming that the state is given is quite a restrictive assumption. \n\n(2) ***Choosing the system state dimensionality and number of constraints***. Section 6 proposes an heuristic for finding the number of constraints of motion. The problem is quite fundamental and it deserves a more thorough investigation. Authors could try to get inspiration from solutions to a similar problem in the scope of system identification and model selection (e.g. MDL, minimum description language for model selection; AIC, Akaike information criteria). \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": "8uL7eZFYbra", "reviewer": "Reviewer_Jrbh", "summary": "In the context of learning models for dynamical systems, the authors propose to learn the constants of motion of the system by means of a neural mapping from states. The method is validated on several (albeit small) simulated dynamical systems, where it appears to improve upon previous approaches.", "questions": "The assumed number of constants of motion (n_c) appears to have a significant effect (Fig. 5). Maybe I missed it, but I couldn't find how these were selected in the benchmark examples?\n\nIt is a bit concerning that it degrades for 50 hidden neurons on the 2D non-linear spring system, with three layers there are still thousands of parameters in that network. Any idea why? Could this have something to do with biases induced by the QR decomposition? It is not obvious to me what kind of biases that will induce.\n", "limitations": "Mostly adequate.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "This paper was mostly clearly written (if occasionally abstract) and enjoyable to read. While mechanics is not my area of expertise, it proposes what seems like a novel idea for modelling of dynamical systems (learning constants of motion). It also empirically improves upon the results of the famous neural ODE paper and later developments. On the other hand, the evaluation is entirely on what seems like quite simple examples. \n\nConcerns:\n-----------------\nPresentation: \n- Some intuition of how the learning interacts with the integration scheme would have been helpful (c.f. NODE paper). It is quite abstract. \n\nExperiments: \n- The assumed number of constants of motion appears to have a significant effect (Fig. 5). How were these selected in the benchmark examples?\n\n- It is a bit concerning that it degrades for 50 hidden neurons on the 2D non-linear spring system, there are still thousands of parameters in that network. Any idea why? Could this have something to do with biases induced by the QR decomposition? It is not obvious to me what kind of biases that will induce.\n\n- As an engineer, I would have liked to see some example on how this approach can be used on some more complex real-world-like sysid problem.\n\nMinor presentation issues:\n\n\"grad\" -> writing out \"gradient\" in text would read better, or typeset it as a math operator (or nabla)\n\nl282: works similar/is comparable? to\n\nl290: which field, or just on science?\nl290: more point*s* \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: Accept