Learning Deterministic and Stochastic Forced Hamiltonian Systems
Authors: Benedikt Brantner, Tomasz Tyranowski
Organizations: Max-Planck-Institut für Plasmaphysik Boltzmannstraße 2, 85748 Garching, Germany · Technische Universität München, Zentrum Mathematik Boltzmannstraße 3, 85748 Garching, Germany · University of Twente, Department of Applied Mathematics PO Box 217, 7500AE Enschede, The Netherlands
We develop a geometric framework for learning deterministic and stochastic forced Hamiltonian systems with neural networks. Motivated by the Lagrange-d'Alembert principle and the theory of variational integrators, we introduce the notion of a Lagrange-d'Alembert map and establish a Cr convergence theorem for first-order one-step methods. Building on these results, we propose Generalized Forced Hamiltonian Neural Networks (GFHNNs), a class of structure-preserving neural networks obtained by concatenating Lagrange-d'Alembert-Euler maps, and prove a universal approximation theorem for this architecture. We further extend the framework to parameter-dependent systems, leading to Parametric Generalized Forced Hamiltonian Neural Networks (PGFHNNs). By interpreting the multiple Stratonovich integrals appearing in the Stratonovich-Taylor expansion as parameters, the same framework can be applied to stochastic forced Hamiltonian systems whenever information about the underlying Wiener process is available. Our numerical experiments demonstrate that the proposed geometric architectures provide significantly improved long-time stability and accuracy compared to non-geometric residual neural networks, while requiring substantially less training data to achieve a comparable level of performance.
Figures & tables
Figure 6.1: Convergence of the training loss for a representative training run using data generated from the linearly damped harmonic oscillator with damping parameter ν=0.1 .
Figure 6.2: Comparison of trajectories generated by the GFHNN and ResNet flows for the linearly damped harmonic oscillator. The left column corresponds to the case ν=0.5 , while the right column shows the case ν=0.001 . In both cases, trajectories are generated from the same initial condition (qˉ0,pˉ0)=(−0.75,−0.75) . The first row displays the position q(t) , the second row the momentum p(t) , and the third row the Hamiltonian H(t) evaluated along the corresponding trajectories.
Figure 6.3: Comparison of trajectories generated by the GFHNN and ResNet flows for the linearly damped harmonic oscillator with ν=0.005 . The left column corresponds to the initial condition (qˉ0,pˉ0)=(−0.75,−0.75) , while the right column corresponds to (qˉ0,pˉ0)=(0.25,0.25) . In both cases, the first row displays the position q(t) , the second row the momentum p(t) , and the third row the Hamiltonian H(t) evaluated along the corresponding trajectories.
Figure 6.4: Averaged trajectory error ϵtraj(t) for trajectories generated by the GFHNN and ResNet flows for the linearly damped harmonic oscillator. For each value of the damping coefficient ν , the error is computed by averaging over 49 test initial conditions distributed in the square −1≤q,p≤1 in phase space, none of which were used during training. Each subplot corresponds to a distinct value of ν , as indicated in the legend.
Figure 6.5: Averaged Hamiltonian error ϵH(t) for trajectories generated by the GFHNN and ResNet flows for the linearly damped harmonic oscillator. The error is computed by averaging over the same set of 49 test initial conditions used in Figure 6.4 . Each subplot corresponds to a distinct value of the damping coefficient ν , as indicated in the legend.
Figure 6.6: Comparison of trajectories generated by the GFHNN and ResNet flows for the quadratically damped pendulum. The left column corresponds to the case ν=0.05 , while the right column shows the case ν=0.01 . In both cases, trajectories are generated from the same initial condition (qˉ0,pˉ0)=(−0.75,−0.75) . The first row displays the position q(t) , the second row the momentum p(t) , and the third row the Hamiltonian H(t) evaluated along the corresponding trajectories.
Figure 6.7: Averaged trajectory error ϵtraj(t) for trajectories generated by the GFHNN and ResNet flows for the quadratically damped pendulum. For each value of the damping coefficient ν , the error is computed by averaging over 49 test initial conditions distributed in the square −0.75≤q,p≤0.75 in phase space, none of which were used during training. Each subplot corresponds to a distinct value of ν , as indicated in the legend.
Figure 6.8: Averaged Hamiltonian error ϵH(t) for trajectories generated by the GFHNN and ResNet flows for the quadratically damped pendulum. The error is computed by averaging over the same set of 49 test initial conditions used in Figure 6.7 . Each subplot corresponds to a distinct value of the damping coefficient ν , as indicated in the legend. Also compare this to the case of the linearly damped harmonic oscillator in Figure 6.5 .
Figure 6.9: Comparison of trajectories generated by the PGFHNN and parametric ResNet flows for the time-dependent damped harmonic oscillator with V0=1 , Ω=3.5 , and ν=0.25 . The left column corresponds to the initial condition (qˉ0,pˉ0)=(−0.75,−0.75) , while the right column corresponds to (qˉ0,pˉ0)=(0.25,0.25) . In both cases, the first row displays the position q(t) , the second row the momentum p(t) , and the third row the Hamiltonian H(t) evaluated along the corresponding trajectories. Note that the plots corresponding to the ResNet with Ntrain=4761 and the PGFHNN overlap very closely and are therefore nearly indistinguishable in the figure.
Figure 6.10: Averaged trajectory ϵtraj(t) ( Top ) and Hamiltonian ϵH(t) ( Bottom ) errors for trajectories generated by the PGFHNN and parametric ResNet flows for the time-dependent damped harmonic oscillator with V0=1 , Ω=3.5 , and ν=0.25 . The errors are computed by averaging over 49 test initial conditions distributed in the square −1≤q,p≤1 in phase space, none of which were used during training. In order to reach a comparable accuracy, the non-structure-preserving ResNet needed a 7.6-times larger training data set than the geometric PGFHNN.
Figure 6.11: Comparison of trajectories generated by the PGFHNN and parametric ResNet flows for the forced Kubo oscillator with β=0.5 and ν=0.25 . The left column corresponds to the initial condition (qˉ0,pˉ0)=(−0.75,−0.75) , while the right column corresponds to (qˉ0,pˉ0)=(0,0.25) . In both cases, the first row displays the position q(t) , the second row the momentum p(t) , and the third row the Hamiltonian H(t) evaluated along the corresponding trajectories.
Figure 6.12: Comparison of trajectories generated by the PGFHNN and parametric ResNet flows for the forced Kubo oscillator with β=0.5 and ν=0.01 . The left column corresponds to the initial condition (qˉ0,pˉ0)=(−0.75,−0.75) , while the right column corresponds to (qˉ0,pˉ0)=(−0.25,−0.25) . In both cases, the first row displays the position q(t) , the second row the momentum p(t) , and the third row the Hamiltonian H(t) evaluated along the corresponding trajectories.
Figure 6.13: Averaged trajectory ϵtraj(t) ( Top ) and Hamiltonian ϵH(t) ( Bottom ) errors for trajectories generated by the PGFHNN and parametric ResNet flows for the forced Kubo oscillator with β=0.5 . The left column corresponds to the friction coefficient ν=0.25 , while the right column corresponds to ν=0.01 . In both cases, the errors are computed by averaging over 49 test initial conditions distributed in the square −1≤q,p≤1 in phase space, and 49 independent sample paths of the Wiener process, none of which were used during training.
Hamiltonian Neural Networks (HNNs) integrate physical priors into neural models by learning a system's Hamiltonian, improving generalization and sample efficiency. Identifying the system Hamiltonian from noisy observations of state variables is a challenging task. For simulations to faithfully reflect the long-term behavior of Hamiltonian systems, especially energy conservation, it is essential to use symplectic integrators, which preserve the system's geometric structure. This fidelity comes at a cost: implicit symplectic integrators are more computationally intensive and make backpropagation through the ODE solver non-trivial. However, by leveraging the fact that symplectic discretizations of the adjoint system yield the same sensitivities associated by backpropagation, we obtain an efficient method of training the Neural Network parameters. In our work, we explore this alternate method of HNN training under noisy observation of trajectories with our HNN model based on an implicit symplectic integrator. Computationally, a predictor-corrector based ODE solver and fixed point iteration help to mitigate the computational cost of the implicit timestepping, resulting in more efficient generation of gradient updates. We showcase the numerical advantage, in experiments, in system identification and energy preservation on a range of non-separable, chaotic systems and the efficient computation and memory complexity of our method. We also observe that the post-processing of the learned Hamiltonian using backward error analysis yields a modified Hamiltonian that is a more accurate approximation of the true Hamiltonian without the need to use more accurate discretizations of the flow map.
Computer Science Department, Indraprasth Institute of Information Technology, New Delhi, India 110020 · Athena Research and Innovation Center, Marousi, Greece 15125 · HSBC Holdings, Singapore
Hamiltonian Neural Networks (HNNs) parameterize conservative dynamics through a learned scalar Hamiltonian, providing an architectural prior that is absent from generic vector-field neural networks. We evaluate this prior under a controlled protocol in which an HNN and a parameter-matched feedforward baseline are trained on the same RK4-generated trajectories, use the same central-difference derivative targets and optimization settings, and are integrated at inference with the same RK4 scheme. Results are reported over five independent training seeds. On the nonlinear pendulum, the HNN reduces mean energy drift by 42-fold and mean trajectory MSE by 15.8-fold at T = 100, approximately 16 pendulum periods. Its energy drift also remains bounded and exhibits substantially lower seed-to-seed variability than the standard-network baseline. An energy-stratified analysis shows that the difference becomes more pronounced as trajectories explore more nonlinear regions of phase space. As an additional diagnostic, we examine an explicit Störmer--Verlet-style rollout of the learned HNN. Because the learned Hamiltonian is not constrained to the separable form H(q,p) = T(p) + V(q), the standard symplecticity guarantee of velocity Verlet does not directly apply. We further apply the same matched-integrator protocol to the three-dimensional Kepler two-body problem. The HNN again exhibits lower trajectory, energy, and angular-momentum drift than the parameter-matched baseline. These experiments provide a controlled study of how Hamiltonian parameterization affects long-horizon prediction and physical consistency across two conservative dynamical systems.
Machine learning of Hamiltonian dynamics has driven growing interest in Hamiltonian neural networks (HNNs), which encode Hamilton's equations of motion into the learning architecture. Despite this progress, it remains unknown whether such networks can predict dynamical regimes absent from their training data, in particular the broad chaotic sea that emerges beyond the observed parameter interval. We address this question using a parameter-aware random-feature Hamiltonian neural network (RF-HNN). Trained using data from only a small number of control-parameter values at which invariant tori dominate, the RF-HNN predicts autonomous long-time dynamics at unseen parameter values where mixed phase space develops and chaotic regions expand, with no data from that regime used in training or model selection. The method is demonstrated across four two-degree-of-freedom Hamiltonian families, including the Hénon-Heiles system. Using Poincaré-section geometry and finite-time Lyapunov exponents, we show that the RF-HNN reproduces the breakup of regular structures and the emergence and growth of chaotic regions, whereas conventionally trained HNNs with the same Hamiltonian structure remain too regular. These results show that what decides parameter extrapolation is not Hamiltonian structure alone but how the fitted Hamiltonian continues in the control parameter. To our knowledge, this is the first demonstration that a learned Hamiltonian can qualitatively extrapolate from predominantly regular dynamics into a broad chaotic sea absent from training.
Jaesung Choi
Center for Artificial Intelligence and Natural Sciences, Korea Institute for Advanced Study, Seoul 02455, Korea