Input-to-State Stable Coupled Oscillator Networks for Closed-form Model-based Control in Latent Space

Authors: maximilian  stölzle
della  cosimo
  santina
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 19:28:33

Abstract: Even though a variety of methods have been proposed in the literature, efficient and effective latent-space control (i.e., control in a learned low-dimensional space) of physical systems remains an open challenge.
We argue that a promising avenue is to leverage powerful and well-understood closed-form strategies from control theory literature in combination with learned dynamics, such as potential-energy shaping.
We identify three fundamental shortcomings in existing latent-space models that have so far prevented this powerful combination: (i) they lack the mathematical structure of a physical system, (ii) they do not inherently conserve the stability properties of the real systems, (iii) these methods do not have an invertible mapping between input and latent-space forcing.
This work proposes a novel Coupled Oscillator Network (CON) model that simultaneously tackles all these issues. 
More specifically, (i) we show analytically that CON is a Lagrangian system - i.e., it possesses well-defined potential and kinetic energy terms. Then, (ii) we provide formal proof of global Input-to-State stability using Lyapunov arguments.
Moving to the experimental side, we demonstrate that CON reaches SoA performance when learning complex nonlinear dynamics of mechanical systems directly from images.
An additional methodological innovation contributing to achieving this third goal is an approximated closed-form solution for efficient integration of network dynamics, which eases efficient training.
We tackle (iii) by approximating the forcing-to-input mapping with a decoder that is trained to reconstruct the input based on the encoded latent space force.
Finally, we leverage these three properties and show that they enable latent-space control. We use an integral-saturated PID with potential force compensation and demonstrate high-quality performance on a soft robot using raw pixels as the only feedback information.

Main paper: 
 introduction: 
 Learning how the environment evolves around us from high-dimensional observations (i.e., world models [1] ) is essential for achieving both artificial and physical intelligence [2] . For example, world models are required for effectively planning an artificial/robotic agent's actions in complex and unstructured environments [3] . However, learning such dynamics directly in high-dimensional observation space is usually intractable. Seminal works have shown that we can leverage autoencoders to compress the state information into a low-dimensional latent space [4, 5] in which it is much more feasible to learn the dynamics [6, 7, 8, 9, 10] . However, strong limitations still persist when it comes to using these learned models to generate low-level intelligence. < l a t e x i t s h a 1 _ b a s e 6 4 = " y B f 9 s g i f h l a X Z u X w j 7 e K c D q 8 m Q Q = " > A A A B 6 3 i c b V B N S 8 N A E J 3 U r 1 q / q h 6 9 L B b B U 0 l E q s e C F 4 8 V 7 A e 0 o W y 2 m 3 b p b h J 2 J 0 I J / Q t e P C j i 1 T / k z X / j p s i N O E + x E d K R E K R t F K D z h w B + W K W 3 U X I O v E y 0 k F c j Q H 5 a / + M G Z p x B U y S Y 3 p e W 6 C f k Y 1 C i b 5 r N R P D U 8 o m 9 A R 7 1 m q a M S N n y 1 O n Z E L q w x J G G t b C s l C / T 2 R 0 c i Y a R T Y z o j i 2 K x 6 c / E / r 5 d i e O N n Q i U p c s W W i 8 J U E o z J / G 8 y F J o z l F N L K N P C 3 k r Y m G r K 0 K Z T s i F 4 q y + v k / Z V 1 a t V a / f X l U Y t j 6 M I Z 3 A O l + B B H R p w B 0 1 o A Y M R P M M r v D n S eC o f H z S N n G q G W + x W M a 6 G 1 D D p V C 8 h Q I l 7 y a a 0 y i Q v B N M b u d + 5 4 l r I 2 L 1 i N O E + x E d K R E K R t F K D z j w B u W K W 3 U X I O v E y 0 k F c j Q H 5 a / + M G Z p x B U y S Y 3 p e W 6 C f k Y 1 C i b 5 r N R P D U 8 o m 9 A R 7 1 m q a M S N n y 1 O n Z E L q w x J G G t b C s l C / T 2 R 0 c i Y a R T Y z o j i 2 K x 6 c / E / r 5 d i e O N n Q i U p c s W W i 8 J U E o z J / G 8 y F J o z l F N L K N P C 3 k r Y m G r K 0 K Z T s i F 4 q y + v k / Z V 1 a t V a / f X l U Y t j 6 M I Z 3 A O l + B B H R p w B 0 1 o A Y M R P M M r v D n S eC o f H z S N n G q G W + x W M a 6 G 1 D D p V C 8 h Q I l 7 y a a 0 y i Q v B N M b u d + 5 4 l r I 2 L 1 i N O E + x E d K R E K R t F K D z i Y D M o V t + o u Q N a J l 5 M K 5 G g O y l / 9 Y c z S i C t k k h r T 8 9 w E / Y x q F E z y W a m f G p 5 Q N q E j 3 r N U 0 Y g b P 1 u c O i M X V h m S M N a 2 F J K F + n s i o 5 E x 0 y i w n R H F s V n 1 5 u J / X i / F 8 M b P h E p S 5 I o t F 4 W p J B i T + d 9 k K D R n K K e W U K a F v Z W w M d W U o Ug O B w z n 3 c G 9 O m E h h 0 H W / n b X 1 j c 2 t 7 c J O c X d v / + C w d H T c M n G q G W + y W M a 6 E 1 L D p V C 8 i Q I l 7 y S a 0 y i U v B 2 O b 2 d + + 5 F r I 2 L 1 g J O E B x E d K j E Q j K K V f L 8 f Y + Y / i W k F L 3 q l s l t 1 5 y C r x M t J G X I 0 e q U vs i r N i z L 2 q R v X e w M y 4 = " > A A A B 7 X i c b V D L S g N B E O y N r x h f U Y 9 e F o M Q L 2 F X J H o M e P E Y w T w g W c L s Z D Y Z M z u z z P Q K M e Q f v H h Q x K v / 4 8 2 / c Z L s Q R M L G o q q b r q 7 w k R w g 5 7 3 7 e T W 1 j c 2 t / L b h Z 3 d v f 2 D 4 u F R 0 6 h U U 9 a g S i j d D o l h g k v W Q I 6 C t R P N S B w K 1 g p H N z O / 9 c i 0 4 U r e 4 z h h Q U w G k k e c E r R S 8 6 m M P f + 8 V y x 5 F W 8 O d 5 X 4 G S l B h n q v + N X t K 5 r G T C I V x J i O 7 y U Y T I h G T g W b F r q p Y Q m h I z J g H U s l i Z k J J v N r p + 6 Z V f p u p L Q t i e 5 c / T 0 x I b E x 4 z i 0 n T H B o V n 2 Z u J / X i f F 6 D q Y c J m k y C R d L I p S 4 a J y Z 6 + 7 f a 4 Z R T G 2 h F D N 7 a 0 u H R J N K N q A C j Y E f / n l V d K 8 q P j V S v X u s l S r Z n H k 4 Q R O o Q w + X E E N b q E O DC E + C 9 R P g i w o L B z f B s = " > A A A B 7 X i c b V D L S g N B E O y N r x h f U Y 9 e F o M Q L 2 F X J H o M e P E Y w T w g W c L s Z D Y Z M z u z z P Q K M e Q f v H h Q x K v / 4 8 2 / c Z L s Q R M L G o q q b r q 7 w k R w g 5 7 3 7 e T W 1 j c 2 t / L b h Z 3 d v f 2 D 4 u F R 0 6 h U U 9 a g S i j d D o l h g k v W Q I 6 C t R P N S B w K 1 g p H N z O / 9 c i 0 4 U r e 4 z h h Q U w G k k e c E r R S 8 6 m M v d F 5 r 1 j y K t 4 c 7 i r x M 1 K C D P V e 8 a v b V z S N m U Q q i D E d 3 0 s w m B C N n A o 2 L X R T w x J C R 2 T A O p Z K E j M T T O b X T t 0 z q / T d S G l b E t 2 5 + n t i Q m J j x n F o O 2 O C Q 7 P s z c T / v E 6 K 0 X U w 4 T J J k U m 6 W B S l w k X l z l 5 3 + 1 w z i m J s C a G a 2 1 t d O i S a U L Q B F W w I / v L L q 6 R 5 U f GM x N k C + W n Z E L q w x J l G j 7 F J K F + j u R 0 9 i Y a R z a y Z j i 2 K x 6 c / E / r 5 d h d B v k Q q U Z c s W W H 0 W Z J J i Q + e V k K D R n K K e W U K a F 3 Z W w M d W U o e 2 n Y k v w V k 9 e J + 2M x N k C + W n Z E L q w x J l G j 7 F J K F + j u R 0 9 i Y a R z a y Z j i 2 K x 6 c / E / r 5 d h d B v k Q q U Z c s W W H 0 W Z J J i Q + e V k K D R n K K e W U K a F 3 Z W w M d W U o e 2 n Y k v w V k 9 e J + 2= " > A A A C H X i c b V D L S s N A F J 3 4 r P E V d e k m W A R X J S l S X R b c u K x g H 9 C E M p n c t E M n k z A z K Z T Q H 3 H j r 7 h x o Y g L N + L f O E 0 j a O u B g c M 5 5 3 L n n i B l V C r H + T L W 1 j c 2 t 7 Y r O + b u 3 v 7 B o X V 0 3 J F J J g i 0 S c I S 0 Q u w B E Y 5 t B V V D H q p A B w H D L r B + G b u d y c g J E 3 4 v Z q m 4 M d 4 y G l E C V Z a G l i X X g B D y n M C X I G Y m S 0 B I S U K Q s 8 z k 0 C C m B R B a X r A w 5 / U w K o 6 N a e A v U r c k l R R i d b A + v D C h G S x H i c M S 9 l 3 n V T 5 O R a K E g Y z 0 8 s k p J i M 8 R D 6 m n I c g / T z 4 r q Z f a 6 V 0 I 4 S o R 9 X d q H + n s h x L O U 0 D n Q y x m o k l 7 2 5 + J / X z 1 R 0 7 e e U p 5 k C T h a L o o z Z K r H n V d k h F U A U m 2 q C i a D 6 r z Y Z Y Y F 1 P U K a u g R 3 + e R V 0 q n X 3 E a t c V e v N h t l H R V 0 i s 7 Q B X L R F W q i W 9 R C b U T Q A 3 p C L + j V e D S e j T f j f R F d M 8 q Z E / Q H x u c 3 D n 2 j H w = = < / l a t e x i t >Predicted observations < l a t e x i t s h a 1 _ b a s e 6 4 = " + 4 / w K f T r P1 S 5 O p N b G 0 o o K Q T K K k w = " > A A A B 8 3 i c b V B N S 8 N A E J 3 U r 1 q / q h 6 9 L B a h X k o i U j 0 W v H i s Y D + g C W W z 3 b Z L N 5 u w O x F K 6 N / w 4 k E R r / 4 Z b / 4 b t 2 0 O 2 v p g 4 P H e D D P z w k Q K g 6 7 7 7 R Q 2 N r e 2 d 4 q 7 p b 3 9 g 8 O j 8 v F J 2 8 S p Z r z F Y h n r b k g N l 0 L x F g q U v J t o T q N Q 8 k 4 4 u Z v 7 n S e u j Y j V I 0 4 T H k R 0 p M R Q M I p W 8 v 0 x x S y e V b H v X f b L F b f m L k D W i Z e T C u R o 9 s t f / i B m a c Q V M k m N 6 X l u g k F G N Q o m + a z k p 4 Y n l E 3 o i P c s V T T i J s g W N 8 / I h V U G Z B h r W w r J Q v 0 9 k d H I m G k U 2 s 6 I 4 t i s e n P x P 6 + X 4 v A 2 y I R K U u S K L R c N U 0 k w J v M A y E B o z l B O L a F M C 3 s r Y W O q K U M b U 8 m G 4 K 2 + v E 7 a V z W v X q s / X F c a 9 T y O I p z B O V T B g x t o w D 0 0 o Q U M E n i G V 3 h z U u f F e X c+ l q 0 F J 5 8 5 h T 9 w P n 8 A c j W R R g = = < / l a t e x i t > ô(t1) < l a t e x i t s h a 1 _ b a s e 6 4 = " j Q f We encode the initial observation o(t0) and the input u(t) into latent space where we leverage the Coupled Oscillator Network (CON) to predict future latent states. Finally, we decode both the latent-space torques τ (t) and the predicted latent states z(t).p I V W J k K J N c q E w Q B p S e Y j g K m 4 = " > A A A B 8 3 i c b V B N S 8 N A E J 3 U r 1 q / q h 6 9 L B a h X k o i U j 0 W v H i s Y D + g C W W z 3 b Z L N 5 u w O x F K 6 N / w 4 k E RA i P J 1 q 4 / 4 H i M 9 K U = " > A A A C G H i c b V D L S g M x F M 3 U V x 1 f o y 7 d D L Z C 3 d S Z L q r L g g g u K 9 g H d E r J Z G 7 b 0 E w y J B m h D P 0 M N / 6 K G x e K u O 3 O v z F 9 C F o 9 E D i c c y 4 3 9 4 Q J o 0 p 7 3 q e V W 1 v f 2 N z K b 9 s 7 u 3 v 7 B 8 7 h U V O J V B J o E M G E b I d Y A a M c G p p q B u 1 E A o 5 D B q 1 w d D 3 z W w 8 g F R X 8 X o 8 T 6 M Z 4 w G m f E q y N 1 H M u g h A G l G c E u A Y 5 s W 8 4 E R H I I L C L Q X 1 I S + K 8 a A f A o + 9 A z y l 4 Z W 8 O 9 y / x l 6 S A l q j 3 n G k Q C Z L G Z p w w r F T H 9 x L d z b D U l D C Y 2 E G q I M F k h A f Q M Z T j G F Q 3 m x 8 2 c c + M E r l 9 I c 3 j 2 p 2 r P y c y H C s 1 j k O T j L E e q l V v J v 7 n d V L d v + p m l C e p B k 4 W i / o p c 7 V w Z y 2 5 E Z V A N B s b g o m k 5 q 8 u G W K J i e l A 2 a Y E f / X k v 6 R Z K f v V c v W u U q h V l 3 X k 0 Q k 6 R S X k o 0 t U Q 7 e o j h q I o E f 0 j F 7 R m / V k v V j v 1 s c i m r O W M 8 f o F 6 z p F 6 f 2 n 3 8 = < / l a t e x i t > Encoder (o) < l a t e x i t s h a 1 _ b a s e 6 4 = " R s 4 v V J G f / Y 0 h N K S c v f A 2 l 1 z A V Z 4 = " > A A A C G H i c b V D L S s N A F J 3 U V 4 2 v q E s 3 w V a o m 5 p 0 U V 0 W d O G y g n 1 A E 8 p k c t s O n U z C z E S o o Z / h x l 9 x 4 0 I R t 9 3 5 N 0 4 f g r Y e G D i c c y 5 3 7 g k S R q V y n C 8 j t 7 a + s b m V 3 z Z 3 d v f 2 D 6 z D o 6 a M U 0 G g Q W I W i 3 a A J T D K o a G o Y t B O B O A o Y N A K h t d T v / U A Q t K Y 3 6 t R A n 6 E + 5 z 2 K M F K S 1 3 r w g u g T 3 l G g C s Q Y / M G S BQ i r g V z 2 p u J / X i d V v S s / o z x J F X A y X 9 R L m a 1 i e 9 q S H V I B R L G R J p g I q v 9 q k w E W m O g O p K l L c J d P X i X N S t m t l q t 3 l U K t u q g j j 0 7 Q K S o h F 1 2 i G r p F d d R A B D 2 h F / S Gshortcomings, such as a limited planning horizon and slow control rates (MPC and gradient-based approaches), sample inefficiency (RL), or they pose a requirement for learning linear dynamics [6] (LQR), which is not even possible for systems that are inherently non-linearizable [22] . One interesting avenue is to leverage model-based control approaches, such as potential shaping [23, 24, 10] , for effective and computationally efficient control in latent space [25] . For these techniques to be feasible, the dynamical model needs to fulfill four characteristics: (i) the dynamics need to have the mathematical structure of physical systems, (ii) conserve the stability properties of real systems, (iii) the latent state needs to be relatively low-dimensional, and (iv) there needs to exist a well-defined, invertible mapping between the input and the forcing in latent space. However, existing model structures that are used for learning latent dynamics [26] do not meet all of these criteria. Relevant examples are Multilayer Perceptrons (MLPs), Neural ODEs (NODEs) [27, 28] , many variants of Recurrent Neural Networks (RNNs) (e.g., LSTMs [29] , Gated Recurrent Units (GRUs) [30] , etc.), and physics-informed neural networks (e.g., Lagrangian Neural Networks (LNNs) [31, 32, 10] , Hamiltonian Neural Networks (HNNs)) [33] . For example, MLPs do not have a physical interpretation and do not provide an invertible mapping of the forcing generated by the input, NODEs are usually not easily stabilizable [34] , most RNNs require a relatively high-dimensional latent space (i.e., many hidden states), and energy-shaping control approaches based on LNNs [10] do not come with any formal stability guarantees.In recent years, oscillatory networks [35, 36, 37, 38, 39] have been shown to exhibit state-of-the-art performance on time sequence modeling tasks while being parameter-efficient, thus fulfilling our requirement (iii). Consequently, we believe that they are a promising option for control-oriented dynamics learning in latent space. Still, these models do not fulfill the remaining requirements that we have listed above. Despite being an interpretable combination of harmonic oscillators, they do not have the structure of a physical system -i.e., they do not possess a well-defined energy function. Moreover, only local stability [35, 37] has been shown, with sufficient conditions that appear to be very stringent. Finally, in addition to training an encoder that maps inputs to latent-space forcing, we propose also training a decoder that learns to reconstruct inputs based on latent-space forcing. This enables us to easily switch between inputs and forcing, which is essential when implementing control strategies.We resolve all the above-mentioned challenges by proposing Coupled Oscillator Networks (CONs), a new formulation of a coupled oscillator network that is inherently Input-to-State Stability (ISS) stable, for learning the dynamics of physical systems and subsequently exploiting its structure for model-based control in latent space. The network consists of damped, harmonic oscillators connected through elastic springs, damping elements, and a neuron-like coupling force and can be excited by a nonlinear actuation term. We identify a transformation into a set of coordinates from which we can derive the networks' kinetic and potential energy. This allows us to leverage Lyapunov arguments [40] for proving the global asymptotic stability of the unforced system and ISS stability for the forced system under relatively mild assumptions on the network parameters. Even though we constrain the dynamics to a very specific structure, we demonstrate (a) the CON network achieves similar performance as NODEs when learning the dynamics of unactuated, mechanical systems with two orders of magnitude fewer parameters and (b) that the proposed model achieves, for the complex task of learning the actuated, highly nonlinear dynamics of continuum soft robots [41, 42] directly from pixels, a 60 % lower prediction error than Coupled Oscillatory Recurrent Neural Network (coRNN) [35] and reaches the SoA performance across all techniques that we tested. Finally, we show some initial results that the proposed CON model is also able to learn the latent dynamics of Partial Differential Equations (PDEs), in this case containing reaction-diffusion [9, 43] dynamics.Subsequently, we derive an approximate closed-form solution, that is, in parameter regimes in which the linear, decoupled dynamics dominate transient, more accurate than numerical integrators with comparable computational requirements and which increases training speed by 2x with a small decrease in prediction accuracy. Finally, as we can derive the system's potential energy, we can leverage potential shaping [23, 24] to derive a controller that combines an integral-saturated PID controller with a feedforward term compensating potential forces. As the feedback acts on a wellshaped potential field, tuning the feedback gains becomes very simple and out-of-the-box, and the controller exhibits a faster response time and a 26 % lower trajectory tracking Root Mean Squared Error (RMSE) than a pure feedback controller based on a latent NODE [27] model.The proposed methodology is particularly well-suited for learning the latent dynamics of mechanical systems with continuous dynamics, dissipation, and a single, attractive equilibrium point. Examples of such systems include many soft robots, deformable objects with dominant elastic behavior, Lagrangian systems immersed in a dominant potential field, or locally other mechanical systems such as robotic manipulators, legged robots, etc. For these systems, we can fully leverage the structural prior of the proposed latent dynamics, including the integrated stability guarantees. If the system is actuated, the learned dynamics can be subsequently exploited for model-based control, as demonstrated in Sec. 5.The code associated with this paper is available on GitHub1 .2 Input-to-State Stable (ISS) Coupled Oscillator Networks (CONs)Formulation. The integral component to (coupled) oscillatory RNNs [35, 36, 37, 38] are onedimensional, potentially damped, harmonic oscillators, which are described by their statey i = [x i ẋi ] T ∈ R 2, where x i and ẋi are the position and velocity of the oscillator, respectively. Then, the oscillator's dynamics are defined by the following Equation of Motion (EOM)EQUATIONHere, m i is the mass, κ i is the stiffness, and d i is the damping coefficient of the damped harmonic oscillator. F i (t) ∈ R is a (possibly time-dependent) external forcing term acting on the mass.Even though the state is extremely low dimensional and the number of parameters is small, this single, damped harmonic oscillator can already exhibit a variety of (designable) behaviors: The expressions ω n,i = κi mi andζ i = di 2√ κi mi let us determine the natural frequency and the damping factor, respectively and allow us to design the transient behavior. For example, ω n,i lets us isolate a spectrum of the input signal F i (t) [37] and ζ i determines the damping regime: underdamped (ω n,i < 1), critically damped (ω n,i = 1), overdamped (ω n,i > 1). Furthermore, as (damped) harmonic oscillators are omnipresent in nature (and especially in physical systems), they have been intensively studied and are well understood (e.g., characteristics, closed-form solutions, etc.). In this work, we will exploit some of these properties and knowledge to learn stable (latent) dynamics efficiently. By intercoupling damped harmonic oscillators, we can drastically increase the expressiveness of the dynamical system [35, 37, 38] while preserving some of the intuition and understanding we have for these systems. In this work, we propose a ISS-stable CON consisting of n damped harmonic oscillators that are coupled through both linear and nonlinear terms. The networks' state is defined as y = x T ẋT T ∈ R 2n and its dynamics can be formulated as a 2 nd -order Ordinary Differential Equation (ODE)EQUATIONwhere K, D ∈ R n×n are the linear stiffness and damping matrices, respectively. The neuron-inspired term tanh(W x(t) + b) with W ∈ R n×n , b ∈ R n provides nonlinear coupling between the harmonic oscillators. The network is excited by the time-dependent input u(t) ∈ R m through the possibly nonlinear mapping g : R m → R n . Specifically, we consider in this work a formulation where an input-dependent matrix B(u) ∈ R n×m projects the input u(t) to a time-dependent forcing on the oscillators: τ = g(u) = B(u) u.Here B(u) could, for example, be parametrized by a MLP.We specifically designed the network architecture such that (i) the system exhibits a unique and isolated equilibrium and (ii) we can derive expressions for the kinetic and potential energies. These two features allow us to (a) prove Global Asymptotic Stability (GAS) and ISS stability using an established procedure based on strict Lyapunov arguments [44, 45] , and (b) implement model-based controller based on potential shaping.One key insight of this work is that in the coordinates x(t), ẋ(t), we cannot derive a potential as the hyperbolic force tanh(W x(t)+b) is not symmetric. Therefore, we propose a coordinate transformation into W-coordinates:yw(t) = xw(t) ẋw(t) = W x(t) W ẋ(t) ∈ R 2n. The coordinate transformation is valid if its Jacobian is full-rank, which is the case if rank(W ) = n. In W-coordinates, the dynamics can be rewritten asEQUATIONwithKw = K W , Dw = D W and Mw = W -1 .A difference of this formulation compared to prior work [35, 36, 37, 38] is that (i) the forcing produced by the input term τ = g(u) is fully separated from the forcing produced by the elastic coupling terms Kw, and (ii) the generalized force is symmetric, which we prove in Appendix A.1, allowing us to define a potential energy expression, which we can later on leverage for stability analysis and control.The equilibria ȳw = xT Proof. The proof is straightforward and provided in Appendix A.2.Next, we introduce a mapping into the tilde coordinates ỹw = yw -ȳw. The residual dynamics (w.r.t. the equilibrium ȳw) can now be stated asEQUATIONIn the following, we will write ∥A∥ to denote the induced norm of matrix A and λm(A), λM(A) to refer to its minimum and maximum Eigenvalue respectively.Global Asymptotic Stability (GAS) for the unforced system. We first consider the unforced system with τ = g(u) = 0, ∀t ∈ [t0, t∞) and strive to prove global asymptotic stability [40] for the attractor x. We propose a strict Lyapunov candidate with skewed level sets [45] Vµ(ỹw) = 1 2 ỹT w PV ỹw + n i=1 xw,i 0 tanh(xw,i + σ + bi) dσ - n i=1 xw,i 0 tanh(xw,i + bi) dσ, = 1 2 ỹT w PV ỹw + n i=1 (lcosh(xw,i + xw,i + bi) -lcosh(xw,i + bi) -tanh(xw,i + bi) xw,i) , with PV = Kw µ Mw µ M T w Mw ∈ R 2n×2n , lcosh(•) = log(cosh(•)),and µ > 0.(EQUATIONLemma 2. The scalar function Vµ(ỹw) defined in (5) is continuously differentiable and verifies the conditionVµ(0) = 0. Furthermore, let Mw, Kw ≻ 0. Now, if we choose 0 < µ < √ λm(Mw) λm(Kw) ∥Mw∥ := µV, then Vµ(ỹw) > 0 ∀ ỹw ∈ R 2n \ {0}. Additionally, then Vµ(ỹw) is radially unbounded as ∥ỹw∥ → ∞ ⇒ Vµ(ỹw) → ∞.Proof. We provide the proof in Appendix A.3 and demonstrate that the bounds on µ are required for the Lyapunov candidate to be positive-definite.Theorem 1. Let Mw, Kw and Dw be positive definite and suppose the system be unforced: g(u(t)) = 0. Then, ỹw = 0 is globally asymptotically stable for the system dynamics defined (4) such that Vµ(ỹw) < 0, ∀ ỹw ∈ R 2n \ {0}.Proof. First, we show that ỹw = 0 is an equilibrium of (4):EQUATIONLemma 2 states that we can always choose µ such that (5) is strict Lyapunov function (e.g., Lipschitz continuous, zero-valued at ỹw = 0, positive definite, and radially unbounded) [40, 45] . We now evaluate the time-derivative of Vµ(ỹw) in the case of an unforced system (i.e., g(u(t)) = 0):EQUATIONwhere we exploited the force balance at equilibrium Kw xw = -tanh(xw + b) and Lemma 7 for defining the upper bound on Vµ(ỹw). Lemma 8 states that P V ≻ 0 for 0 < µ < µ V. Similarly, Lemma 2 requests that 0 < µ < µV. Indeed, both conditions can always be fulfilled by choosing µ ∈ (0, min{µV, µ V}). With P V ≻ 0 ⇔ λm (P V) > 0 [46] , we can state that Vµ(ỹw) < 0 ∀ ỹw ∈ R 2n \ {0} and conclude that the unforced system is globally asymptotically stable around ỹw = 0.Global Input-to-State Stability (ISS) for the forced system. We now take the forcing g(u) into account again and demonstrate that the system states remain proportionally bounded to the initial conditions and as a function of the supremum of the input forcing.Theorem 2. Suppose Mw, Kw, Dw ≻ 0, 0 < θ < 1, and that we choose 0 < µ < min{µV, µ V}}. Then, (4) is globally Input-to-State Stable (ISS) such that the solution ỹw(t) verifies ∥ỹw∥2 ≤ β (∥ỹw(t0)∥2, t -t0) + γ supEQUATIONwhere β(r, t) ∈ KL, γ(r) =(1+µ 2 ) λ M (P V ) r 2 +4 θ √ n √ 1+µ 2 λm(P V ) r θ 2 λm(P V ) λ 2 m (P V ) ∈ K.Proof. The proof is provided in Appendix A.5. We first demonstrate that the ISS-Lyapunov function is bounded from both sides by K∞ functions. Subsequently, we derive the energy dissipation of the forced system and establish attracting regions as a function of the norm of the forcing. Outside these regions, it is ensured that the system has a minimal rate of decay and, therefore, converges exponentially fast into the attracting region. 3 An approximate closed-form solution for the rollout of CON To predict future system states, we need to integrate the ODE in Eq. ( 2), with the solution given by y(t k+1 ) = yt k + t k+1 t k f (y(t ′ ), u(t ′ )) dt ′ .Unfortunately, a closed-form solution for the nonlinear dynamics f (y, u) does not (yet) exist. Therefore, we traditionally need to revert to (high-order) numerical ODE solvers that are computationally very expensive and introduce additional memory overhead [47] . This considerably increases the training time of models involving such continuous-time dynamics. While the computational time can be reduced by increasing the (minimum) time step of the integrator, this comes at the expense of an integration error, and we lose (part of) the theoretical guarantees and practical characteristics that the nominal ODE provides. In this work, we take an alternative approach by splitting the problem into (i) decoupled linear dynamics that can be cheaply and precisely integrated using a closed-form solution and (ii) the residual, coupled nonlinear dynamics, which we integrate numerically at a slower time scale:EQUATIONwhere κ = diag(K11 . . . Knn), d = (D11 . . . Dnn) are the diagonal components of the stiffness and damping matrices, respectively, and F ∈ R n is a constant, external forcing term on the oscillators.For a short-time-interval δt, we now approximate (9) asEQUATION)For a scalar 2 nd -order, linear ODE of form ẏi = f ld,i (x(t ′ ), F (t k )), a well-known, closed-form solution [48] exists. We exploit this characteristic by formulating the approximate solution asEQUATIONand denote fCFA-CON : R n × R m → R n as the Closed-Form Approximation of the Coupled Oscillator Network (CFA-CON) model. The implicit assumption behind (10) is that f ẍ,ld (y) ≻ f ẍ,nld (y, u) (i.e., the linear, decoupled dynamics dominate the nonlinear, coupled, time-varying dynamics). We refer the interested reader to Appendix B for derivation and implementation details, where we summarize the integration procedure in Algorithm 1. We also provide qualitative results for the integration accuracy in Fig. 2 and quantitative results for the integration accuracy and computational speed-up w.r.t. numerical integrators in Appendix B. 
learning control-oriented latent dynamics from pixels: 
 We now move towards learning latent dynamical models based on CON and CFA-CON. CONs are an ideal fit for learning latent dynamics as they guarantee that the latent states stay bounded.We assume to have access to observations in the form of images o ∈ R ho×wo×co , where co denotes the number of channels. Please note that this could also be other high-dimensional observations such as LiDAR scans, point clouds, etc. We now leverage an encoder-decoder architecture to map these high-dimensional observations into a compressed latent space: The encoder Φ : R ho×wo×co → R nz with nz ≪ ho wo identifies a lowdimensional latent representation z ∈ R nz of the images. The decoder Ψ : R nz → R ho×wo×co approximates the inverse operation by reconstructing an image ô ∈ R ho×wo×co based on the latent representation. To promote the learning of a smooth and monotonic mapping into latent space, we specifically choose to implement the autoencoder here as a β-Variational Autoencoder (VAE) [5, 49] . Instead of just statically reconstructing the image ô(t k ), we are interested in predicting future observations ô(t k+l ), where l ∈ 1 . . . N . For this, we train a 2 nd -order dynamical model that is, when integrated, able to predict future latent representations z(t k+l ). This requires us to define a latent state ξ(t) = z T (t) żT (t) T ∈ R 2nz consisting of the latent representation and latent velocity ż(t) ∈ R nz .We now rely on CON with n = nz oscillators to provide us with the latent state derivative ξ = fw(yw(t), u(t)), where we defined ξ = yw, and z = xw. To ensure stability, we make use of the Cholvesky decomposition to ensure that Mw, Kw and Dw always remain positive definite (see Theorem 2). It is important to note that we train the encoder, decoder, and dynamical model all jointly. Please refer to Appendix C for more implementation details.Training. It is important to remember that because we are using a β-VAE [5, 49] , the image encoding becomes stochastic, and the encoder neural network actually outputs µz(o), 2 log(σz)(o) ∈ R nz . After executing the reparametrization trick as z(t k ) ∼ N (µz(t k ), σ 2 z (t k )), we formulate the loss function, evaluated on each trajectory consisting of N time-steps, asEQUATIONModel RMSE M-SP+F [26] ↓ RMSE S-P+F [26] ↓ RMSE D-P+F [26] ↓ RMSE CS ↓ RMSE PCC-NS-2 ↓ RMSE PCC-NS-3 ↓ RNN 0.2739 ± 0.0057 0.2378 ± 0.0352 0.1694 ± 0.0004 0.1011 ± 0.0009 0.1373 ± 0.0185 0.2232 ± 0.0075 GRU [30] 0.0267 ± 0.0033 0.1457 ± 0.0078 0.1329 ± 0.0005 0.1125 ± 0.0100 0.0951 ± 0.0021 0.2148 ± 0.0196 coRNN [35] 0.0265 ± 0.0002 0.1333 ± 0.0044 0.1324 ± 0.0016 0.2537 ± 0.0018 0.2504 ± 0.0899 0.2474 ± 0.0018 NODE [27] 0.0264 ± 0.0010 0.1260 ± 0.0013 0.1324 ± 0.0024 0.2415 ± 0.0021 0.1867 ± 0.0561 0.3373 ± 0.0565 MECH-NODE 0.0328 ± 0.0034 0.1650 ± 0.0205 0.1710 ± 0.0111 0.2494 ± 0.0028 0.1035 ± 0.0012 0.1900 ± 0.0024 CON-S (our) 0.0303 ± 0.0053 0.1303 ± 0.0064 0.1323 ± 0.0018 0.1993 ± 0.0646 0.0996 ± 0.0012 0.1792 ± 0.0038 CON-M (our) 0.0303 ± 0.0053 0.1303 ± 0.0064 0.1323 ± 0.0018 0.1063 ± 0.0027 0.1008 ± 0.0006 0.1785 ± 0.0023 CFA-CON (our) 0.0313 ± 0.0026 0.1352 ± 0.0073 0.1307 ± 0.0012 0.1462 ± 0.0211 0.1124 ± 0.0025 0.1803 ± 0.0003 Table 1 : Benchmarking of CON and CFA-CON at learning latent dynamics against baseline methods. The first three datasets, based on [26] , contain samples of a mass-spring with friction (M-SP + F), a single pendulum with friction (S-P + F), and a double pendulum with friction (D-P + F) (all without system inputs). The CS dataset considers a continuum soft robot consisting of one segment with three constant planar strains. The PCC-NS-2 and PCC-NS-3 datasets contain trajectories of a continuum soft robot made of two and three piecewise constant curvature segments, respectively. We choose the latent dimensions of the models as nz = 4, nz = 4, and nz = 12 for the M-SP + F, S-P + F, and D-P + F datasets, and nz = 8, nz = 12, and nz = 12 for the PCC-NS-2, PCC-NS-3, and CS soft robotic datasets. We report the mean and standard deviation over three different random seeds. whereẑ(t k ) is predicted by ξ(t k ) = t k t=t 0 f ξ (ξ(t ′ ), u(t ′ )) dt ′ , and ξ(t0) = z T (t0) żT (t0) T . Here, z(t0)is given by the encoder, and ż(t0) is approximated using finite differences in image-space (see Appendix C.5 for more details). β, λ ⃗ o , λz ∈ R are loss weights.Models. We train the CON with the input-to-forcing mapping g(u) = B(u) u, where B(u) is parametrized by a MLP with a hyperbolic tangent activation function applied in between layers. We report results for two variants of the CON model: for the medium-sized CON-M and small-sized CON-S, the MLP consists of five and two layers with a hidden dimension of 30 and 12, respectively. The model CFA-CON uses the same architecture as CON-M. We compare against several popular latent space model architectures: The NODE model uses a MLP with an hyperbolic activation functions and predicts ξ(t) = fNODE(ξ(t), u(t)).To make the comparison fair, we parametrize the NODE's MLP in the same fashion as for CON-M. The MECH-NODE integrates prior knowledge towards learning 2 nd -order mechanical ODEs and, therefore, predicts z(t) = fMECH-NODE(ξ(t), u(t)). Furthermore, we consider multiple autoregressive models: RNN, GRU, and coRNN and let them parameterize the following transition function: ξ(t k+1 ) = far(ξ(t k ), u(t k )) As common in the relevant literature [26] , we allow the autoregressive models to perform multiple time step transitions before predicting the next sample. For the autoencoder, we use a vanilla Convolutional Neural Network (CNN). More details can be found in Appendix C.Datasets. We consider in total six datasets that are based on simulations of unactuated mechanical systems, and actuated continuum soft robots. The first three, mechanical dataset are based on the work of Botev et al. [26] and contain video sequences of a mass-spring system with friction (M-SP+F), a single pendulum with friction (S-P+F), and a double pendulum with friction (D-P+F). Continuum soft robots have theoretically infinite Degree of Freedom (DOF), evolve with highly nonlinear and often time-dependent dynamical behaviors, and are notoriously difficult to model from first principles [50] . For that reason, it is a very interesting proposition if we could learn latent-space dynamical models directly from video [51] and later leverage them for control [52] . Therefore, we generate three datasets based on the Piecewise Constant Strain (PCS) soft robot model. CS considers one segment with constant strain and is modeled using three configuration variables. PCC-NS-2 and PCC-NS-3 only consider bending deformations and contain soft robots with two and three segments, respectively. For all datasets, we render images with a resolution of 32x32px of the system's state. More information on the datasets can be found in Appendix C.1.We tune all hyperparameters for each model and dataset separately using Optuna [53] .The results in Tab. 1 show that the NODE model slightly outperforms the CON network on the M-SP+F and S-P+F datasets. However, as the datasets do not consider system inputs, we can remove the input mapping from all models (e.g., RNN, GRU, coRNN, CON, and CFA-CON). With that adjustment, the CON network has the fewest parameters among all models, particularly two orders of magnitude less than the NODE model. Therefore, we find it very impressive that the CON network is roughly on par with the NODE model. For the D-P+F dataset, we can conclude that the CFA-CON model offers the best performance across all methods. Finally, most of the time, the CON & CFA-CON networks outperform the other baseline methods that have more trainable parameters. Actuated continuum soft robot datasets: The results in Tab. 1 show that CON-M matches the performance of the state-of-the-art methods across all experiments. In the case of PCC-NS-3, CON-M even decreases the RMSE error by 6 % w.r.t. the closest baseline method (MECH-NODE). Impressively, the performance is not reduced (but instead often even improved) compared to other models that offer a much larger design space for learning the dynamics (e.g. NODE). Furthermore, CON-S and CFA-CON often only exhibit slightly lower performance than CON-M, even though they have significantly fewer parameters and consider an approximated solution, respectively. Supplementary results (e.g., more evaluation metrics) can be found in Appendix D. We also conduct on the PCC-NS-2 dataset an analysis concerning the effect of the latent dimension on the performance (see Fig. 3 ). For this experiment, all hyperparameters were tuned for nz = 8, and we observe that the CON models have a much-improved consistency and smaller variance w.r.t. the baseline methods when the latent dimensionality is increased.The results for a dataset based on the 1 st -order reaction-diffusion PDE [9] are presented in Apx. D.1.Question: Does the paper fully disclose all the information needed to reproduce the main experimental results of the paper to the extent that it affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?Answer: [Yes] Justification: All the code for the experiments has been open-sourced on GitHub. Furthermore, we provide a detailed description of the implementation details in Appendix C. 
exploiting the dynamic structure for latent-space control: 
 We consider the problem setting of guiding the system towards a desired observation o d by providing a sequence of inputs u(t k ) such that at time tN , the actual observation o(tN ) matches o d . A relatively simple way would be to encode the desired observation into latent space z d = Φ(o d ) and then to design a feedback controller (e.g., PID) in latent space: g(u) = PID(z d -z(t), -ż). Unfortunately, several challenges appear: first, it is not clear how the latent-space force τ = g(u) can be mapped back to an actual input u(t) as the inverse input-to-forcing mapping g -1 is generally not known. Furthermore, relying entirely on a PID controller has several well-known drawbacks, such as poor and slow transient behavior, steady-state errors (in case the integral gain is chosen to be zero), and instability for high proportional and integral gains. We take inspiration from potential shaping strategies [23, 24] , which are widely used for effectively controlling (elastic) robots, and, therefore, combine a feedforward term compensating the latent-space potential forces with an integral-saturated, PID-like feedback term. For mapping the desired forcing τ back to an input u(t), we train an forcing decoder η : R n → R m that approximates g -1 . Specifically, we consider here the structure u = η(τ ) = E(τ ) τ , where E ∈ R m×n is parameterized by an MLP.The latent-space control law is given byτ (t) = g(u) = Kw z d + tanh(z d + b)Feedforward term: compensation of potential forces+ Kp (z d -z) -K d ż + Ki t 0 tanh(υ(z d (t ′ ) -z(t ′ )))dt ′Feedback term: P-satI-D (13) where Kp, Ki, K d ∈ R n×n are the proportional, integral, and derivative control gains, respectively. As integral terms can often lead to instability when applied to nonlinear systems [54] , we adopt an integral term saturation [55] with the associated dimensionless gain υ ∈ R, which ensures that the integral error added at each time step is bounded to the interval (-1, 1). Subsequently, τ is decoded to the input as u(t) = η(τ ) = E(τ ) τ . For training this decoder, we add a reconstruction loss to (12) :Lu(t k ) = λu MSE(u(t k ), û(t k )) = λu MSE (u(t k ), η(g(u(t k )))).Experimental setup. We train a CON model with two latent variables (nz = 2) on the PCC-NS-2 dataset.Analog to the input encoder mapping B(u), the forcing decoder mapping E(τ ) is parametrized by an MLP consisting of five layers with hidden dimension 30 and a hyperbolic tangent activation function. The CON model achieves an RMSE of 0.1628 on the test set. We benchmark two controllers on the simulated continuum soft robot consisting of two segments: (i) a pure P-satI-D controller (i.e., the feedback term in (13) ) that leverages the smooth mapping into the latent representation enabled by the CON dynamic model and the β-VAE, and (ii) a P-satI-D+FF (i.e., (13) ) that exploits the structure of the CON dynamics by compensating for the potential forces. The stable closed-loop system dynamics made the control gain tuning very easy, and we selected < l a t e x i t s h a 1 _ b a s e 6 4 = " j f q H c M 0 B S P +J P f E D h e J H a x E 0 R m w = " > A A A C K H i c b V B d S w J B F J 3 t 0 7 Y v q 8 d e h i Q w R N n 1 w X p L C K Q X w T A z c E V m x 6 s N z s 4 u M 7 O B L P 6 c X v o r v U Q U 4 W u / p P E j K O 3 A w L n n 3 s O d e / y I M 6 U d Z 2 y t r K 6 t b 2 y m t u z t n d 2 9 / f T B 4 Z 0 K Y 0 m h Q U M e y n u f K O B M Q E M z z e E + k k A C n 0 P T H 1 x N + s 1 H k I q F 4 l Y P I 2 g H p C 9 Y j 1 G i j d R J X 3 o + 9 J l I K A g N c m R X i V L 5 e i S Z 6 H u e n c M V y e h k 0 h T Z a r 5 e y 1 X O b A 9 E 9 8 f Q S W e c g j M F X i b u n G T Q H L V O + s 3 r h j Q O j J 1 y s 6 z l O p F u J 0 R q R j m M b C 9 W E B E 6 I H 1 o G S p I A K q d T A 8 d 4 V O j d H E v l O Y J j a f q b 0 d C A q W G g W 8 m A 6 I f 1 G J v I v 7 X a 8 W 6 d 9 F O m I h i D Y L O F v V i j n W I J 6 n h L p N A N R 8 a Q q h k 5 q + Y P h B J q M l A 2 S Y E d / H k Z X J X L L i l Q u m m m C k X 5 3 G k 0 D E 6 Q V n k o n N U R t e o h h q I o i f 0 g t 7 R h / V s v V q f 1 n g 2 u m L N P U f o D 6 y v b + l P p T Q = < / l a t e x i t > Mass-Spring + Friction (M-SP+F)< l a t e x i t s h a 1 _ b a s e 6 4 = " 0 u c u Z 3 5 y g j t 1 m l M d t H l S i 8 X 1 1 m Q = " > A A A C K H i c b V D L S g M x F M 3 4 d n x V X b o J F k E X l p k K 6 k 6 h G 1 d S H 1 W h U 0 o m v V O D m W R I 7 g h l 6 O e 4 8 V f c i C j i 1 i 8 x r R V 8 H Q i c n H M P y T 1 x J o X F I H j z x s Y n J q e m Z 2 b 9 u f m F x a X S 8 s q F 1 b n h 0 O B a a n M V M w t S K G i g Q A l X m Q G W x h I u 4 5 v a w L + 8 B W O F V u f Y y 6 C V s q 4 S i e A M n d Q u H U Q x d I U q O C g E 0 / d r W q F Q e Z 5 G k X + m E 6 S n O t b o L p v 1 W m 3 7 + G x 7 Z 8 u P Q H W + A u 1 S O a g E Q 9 C / J B y R M h m h 3 i 4 9 R R 3 N 8 9 T F u W T W N s M g w 1 b B D A o u o e 9 H u Y W M 8 R v W h a a j i q V g W 8 V w 0 T 7 d c E q H J t q 4 o 5 A O 1 e + J g q X W 9 t L Y T a Y M r + 1 v b y D + 5 z V z T P Z b h V B Z j q D 4 5 0 N J L i l q O m i N d o Q B j r L n C O N G u L 9 S fV R A d w 0 1 k a Q F X J K O T e X v J l 8 + r h c q Z G 1 j / J 9 X O 5 r y i N w V e J v 6 c 5 N A c 1 X b 2 N e j E 1 E Q 2 T j l R q u l 7 i W 6 l R G p G O Y z c w C h I C B 2 Q H j Q t F S Q C 1 U q n 2 4 7 w q V U 6 u B t L e 4 T G U / V 3 I i W R U s M o t J M R 0 X 2 1 6 E 3 E / 7 y m 0 d 2 r V s p E Y j Q I O n u o a z j W M Z 5 U h z t M A t V 8 a A m h k t m / Y t o n k l D b g X J t C f 7 i y s v k o V T 0 L 4 o X d 6 X c d W l e R w Y d o x O U R z 6 6 R N f o F l V R D V H 0 h F 7 Q G LA A A C H X i c b V D L S s N A F J 3 U V 4 2 v q k s 3 w S K 4 s S R F q s u C G 9 1 V s Q 9 o Q p l M b t q h k 0 m Y m S g l 9 E f c + C t u X C j i w o 3 4 N 0 7 T C t p 6 Y O B w z r n c u c d P G J X K t r + M w t L y y u p a c d 3 c 2 N z a 3 i n t 7 r V k n A o C T R K z W H R 8 L I F R D k 1 F F Y N O I g B H P o O 2 P 7 y Y + O 0 7 E J L G / F a N E v A i 3 O c 0 p A Q r L f V K p 6 4 P f c o z A l y B G J s 3 g N n J f S x Y 4 L r m F Z c J F X n S d I E H P 6 l e q W x X 7 B z W I n F m p I x m a P R K H 2 4 Q k z T S 4 4 R h K b u O n S g v w 0 J R w m B s u q m E B J M h 7 k N X U 4 4 j k F 6 W X z e 2 j r Q S W G E s 9 O P K y t X f E x m O p B x F v k 5 G W A 3 k v D c R / / O 6 q Q r P v Y z y J F X A y X R R m D J L x d a k K i u g A o h i I 0 0 w E V T / 1 S I D L D D R H U h T l + D M n 7 x I W t W K U 6 v U r q v l e n V W R x E d o E N 0 j B x 0 h u r o E j V Q E x H 0 g J 7 Q C 3 o 1 H o 1 n 4 8 1 4 n 0 Y L x m x m H / 2 B 8 f k N l Z C i 0 Q = = < / l a t e x i t >< l a t e x i t s h a 1 _ b a s e 6 4 = " q N L g n 8C t v x 2 R S 5 k a u X 0 A 8 Q b S 8 6 w = " > A A A C H X i c b V D L S s N A F J 3 U V 4 2 v q E s 3 g 0 V w V Z I i 1 W V B F 7 q r Y m u h C W U y u W 2 H T i Z h Z i K U 0 B 9 x 4 6 + 4 c a G I C z f i 3 z h 9 C N p 6 Y e B w z j 3 3 z j 1 h y p n S r v t l F Z a W V 1 b X i u v 2 x u b W 9 o 6 z u 9 d U S S Y p N G j C E 9 k K i Q L O B D Q 0 0 x x a q Q Q S h x z u w s H 5 W L + 7 B 6 l Y I m 7 1 M I U g J j 3 B u o w S b a i O c + K H 0 G M i p y A 0 y J F 9 A y I C C Z H v 2 x d E m 8 E a X x k L 2 L 4 R f r o 6 T s k t u 5 P C i 8 C b g R K a V b 3 j f P h R Q r P Y 2 C k n S r U 9 N 9 V B T q R m l M P I 9 j M F K a E D s 6 h t o C A x q C C f X D f C R 4 a J c D e R 5 g m N J + x v R 0 5 i p Y Z x a D p j o v t q X h u T / 2 n t T H f P g p y J N N M g 6 H R R N + N Y J 3 g c F Y 6 Y B K r 5 0 A B C J T N / x b R P J K E m A 2 W b E L z 5 k x d B s 1 L 2 q u X q d a V U q 8 z i K K I D d I i O k Y d O U Q 1 d o j p q I I o e 0 B N 6 Q a / W o / V s v V n v 0 9 a C N f P s o z 9 l f X 4 D 2 e m i W g = = < / l a t e x i t > Rendered Dataset ImageWhereas evolution has integrated self-assembly, self-organization and periodic RD processes at the nanoscale to produce beautiful nanomachinery, there is clearly a lot to do to combine these three processes in nanoscale abiological systems. We believe that the advances outlined in this Review have laid the foundation for investigating self-assembly, self-organization and RD at the nano-and microscales in a holistic manner.Recent studies of the prototypical BZ oscillating chemical reaction in nano-and microdroplets have provided considerable insight into how emergent properties arise from RD processes. Moreover, future challenging task would be to render such DNA objects with desired functionalities as well 47 . From an instrumentation point of view, the rapid progress of tools such as cryo-electron microscopy 103 and X-ray crystallography 104 has greatly enhanced our capability to study self-assembly, self-organization and RD processes.Finally, observing that a key feature of life is its rhythmic nature 105 , we take chemical oscillation as an example to direct the reader towards a set of speci c problems that merit attention from researchers. First, new chemical oscillators: the BZ reaction has revealed many important insights about oscillating reactions and may still serve as a good toy model; how can we use these principles to develop new types of chemical oscillator for controlling RD processes 106,107 or regulating RD processes in living systems? Second, communication, computation and memory: we know that BZ oscillators can communicate and exhibit certain features of computation 108 and memory systems, so how can we develop synthetic systems, with intermediate complexity between chemical oscillators and the brain, but not based on silicon electronics, that can perform such tasks?ird, motion: nature has evolved ATP motors by integrating selfassembly, self-organization and RD processes, so what do we need to construct an analogous e cient chemomechanical engine? Fourth, periodicity: how can nanoscale oscillations that result in micro-and macroscale phenomena help us to develop RD processes for controlling the fate of cells and treating diseases? (Although some progress has been made towards this latter goal, we are still at a very early stage before we can treat diseases 109-111 .) Fi h, multiscale aspects: although there is an observable, phenomenological link between a variety of nanoscale processes and microscale phenomena, how can we rationalize how large-scale behaviour emerges from an understanding of nanoscale processes and interactions?We hope that this short Review will stimulate scientists to develop new, nanoscale, synthetic systems for a wide range of functions and inspire them to shi from studying single molecules (or objects) to multiple entities, from homotypic to heterotypic molecular interactions, from single to multiple events, from static to dynamic patterns, and from structure to function. Of course, it is unlikely that the integration of self-assembly, self-organization and RD processes at the nanoscale will produce systems that capture all the features of living systems. But exploration along these lines will undoubtedly not only lead to a new paradigm of doing nanoscience and to the development of new nanotechnologies, but may also help address fundamental questions, such as how life evolves from the coupling of chemical reactions 1,112-114 and other types of nonlinear behaviour 115,116 . Z Z e x J Q r W a A h t 1 O M 7 D d j C x w / H 0 = " > A A A C J 3 i c b V D L S g M x F M 3 U V x 1 f V Z d u g k X Q h W W m i + p K C n b h s h b 7 g E 4 p m f R O G 5 r J D E l G K E P / x o 2 / 4 k Z Q E V 3 6 J 6 Y v U O u F w L n n 3 M P N P X 7 M m d K O 8 2 l l V l b X 1 j e y m / b W 9 s 7 u X m 7 / o K G i R F K o 0 4 h H s u U T B Z w J q G u m O b R i C S T 0 O T T 9 4 f V E b 9 6 D V C w S d 3 o U Q y c k f c E C R o k 2 V D d 3 5 f n Q Z y K l I D T I s Y 1 x D Q i d a J 5 n m g o L g k Q t u t P a e e X M 9 k D 0 F v P d X N 4 p O N P C y 8 C d g z y a V 7 W b e / F 6 E U 1 C Y 6 e c K N V 2 n V h 3 U i I 1 o x z G t p c o i A k d k j 6 0 D R Q k B N V J p 3 e O 8 Y l h e j i I p H l C 4 y n 7 0 5 G S U K l R 6 J v J k O i B + q t N y P + 0 d q K D y 0 7 K R J x o E H S 2 K E g 4 1 h G e h I Z 7 T A L V f G Q A o Z K Z v 2 I 6 I N I E Z a K 1 T Q j u 3 5 O X Q a N Y c E u F 0 m 0 x X y 7 O 4 8 i i I 3 S M T p G L L l A Z 3 a A q q i O K H t A T e k V v 1 q P 1 b L 1 b H 7 P R j D X 3 H K J f Z XW K W h h h p J q h S X c d O t J d R q T k T O C 3 1 U o U J Z W M 6 x K 6 h E Q 1 R e d n 8 0 C k 5 M c q A B L E 0 F W k y V 5 c n M h o q N Q l 9 0 x l S P V K / v Z n 4 l 9 d N d X D t Z T x K U o 0 R W y w K U k F 0 T G Z f k w G X y L S Y G E K Z 5 O Z W w k Z U U q Z N N i U T w v e n 5 H / S O q 8 6 b t V t X F Z q b h 5 H E Y 7 g G E 7 B g S u o w S 3 U o Q k M E B 7 h G V 6 s e + v JW i r N J Y U m T X k q O y F R w J m A p m a a Q y e T Q J K Q Q z s c X U / 9 9 j 1 I x V J x p 8 c Z 9 B I S C x Y x S r S R + o 4 b h B A z U V A Q G u T E x r i R a s M Z 4 V g N S c Z E b N s B i M F y o u 9 U 3 K o 7 A 1 4 l n n / h + t h b K B W 0 Q K P v f A S D l O a J W a e c K N X 1 3 E z 3 C i I 1 o x w m d p A r y A g d k R i 6 h g q S g O o V s 2 Q T f G a U A Y 5 S a Z 7 Q e K a u b h Q k U W q c h G Y y I X q o f n t T 8 S + v m + v o q l c w k e U m L Z 0 f i n K O d Y q n N e E B k 0 A 1 H x t C q G T m r 5 g O i S T U d K B s U 8 I y K f 6 f t G p V z 6 / 6 t 7 V K 3 V / U U U Y n 6 B S d I w 9 d o j q 6 Q Q 3 U R B Q 9 o Cf 0 g l 6 t R + v Z e r P e 5 6 M l a 7 F z j H 7 A + v w G e 9 q g l g = = < / l a t e x i t >< l a t e x i t s h a _ b a s e = " u f dV O n c B d s D T J m o m s i h M G k = " > A A A C J i c d V D L S s N A F J r P E V d e l m s A h u L E m F E o K b l y q G A f J Q y m d y Q y e T M D M R S u j f u P F X A g q o k v / x O k L q u i F C d z z u X e e K U M V d N a W l Z X V s v b N i b W s u e f l l m a R Q o w l P Z D M k C j g T U N N M c i m E k g c c m i E g u x r g H q V g i v Q w h X Z M e o J G S X a U B n M g i h x R O Q W i Q I x v j G I N P l U p o R A E h u g m k j L R w x H Q J A J p B y C i u b / j F N S O y m C D z / z P W x N O K a F b V j v M S R A n N Y j N O O V G q b m p b u d E a k Y j O w g U A D g P W g Y K E o N q M / R / j Y M N H H t N C w m O J G T W K l h H B p n T H R f / d b G F a K P d i b O R J q Z + l U T f j W C d H B q O m A S q + d A A Q i U z t L a J J Q k G y T Q j z T / H / o F u e X J v y X K / s j g I R E f o B H n o H F X Q N a q i G q L o A T h V / R m P V r P r v M b U u W b O Z A / S j r K v v t e l w = = < / l a t e x i t >Latent-space forcing decoder < l a t e x i t s h a 1 _ b a s e 6 4 = " D q B Y l r 2 5 3 h P 9 U s 3 g a s 2 v O w a U J 7 A = " > A A A C G H i c d V D L S g M x F M 3 4 r O O r 6 t J N s A h u W m c q j C 4 L u t B d B f u A t p R M 5 k 4 b m s k M S U Y o Q z / D j b / i x o U i b r v z b 0 x f U E U P B A 7 n n M v N P X 7 C m d K O 8 2 W t r K 6 t b 2 z m t u z t n d 2 9 / f z B Y V 3 F q a R Q o z G P Z d M n C j g T U N N M c 2 g m E k j k c 2 j 4 g + u J 3 3 g E q V g s H v Q w g U 5 E e o K F j B J t p G 7 + v O 1 D j 4 m M g t A g R z b G 1 a I i + q 5 4 g 0 O A w C d 0 Y L d B B I t A N 1 9 w S s 4 U e J m 4 3 o X j Y X e u F N A c 1 W 5 + 3 A 5 i m k Z m n H K i V M t 1 E t 3 J i N S M c h j Z 7 V R B Y r a Q H r Q M F S Q C 1 c m m h 4 3 w q V E C H M b S P K H x V F 2 e y E i k 1 D D y T T I i u q 9 + e x P x L 6 + V 6 v C q k z G R p B o E n S 0 K U 4 5 1 j C c t 4 Y B J o J o P D S F U M v N X T P t E E m o 6 U L Y p Y X E p / p / U y y X X K 3 n 3 5 U L F m 9 e R Q 8 f o B J 0 h F 1 2 i C r p F V V R D F D 2 h F / S G 3 q 1 n 6 9 X 6 s D 5 n 0 R V r P n O E f s A a f w O u X 5 + M < / l a t e x i t > P-satI-D feedback < l a t e x i t s h a 1 _ b a s e 6 4 = " j W B i s V 3 a Y + 4 Z N O q V w k i 9 M L c R R J E = " > A A A C i H i c d V F b T 9 s w F H b C G K z s U u B x L 0 e r G K 3 Y q o R L G Q 9 I i L 1 M 2 g t I K y A 1 J X I c h 1 o 4 T m S f I L V R f s v + E 2 / 8 G 5 x S RX i m R C E b R U m H z b 4 C 0 C I O U 4 k i n Z R J V c P g V f v 8 j 8 g q C b 9 C e X D 4 S c Q X f Y d K x 1 5 N T P H U K 4 g z L S Q V b c 4 q w i l A Y e p c I t p A a t Y M i N 0 J m 6 n n O N m 7 W K S f 1 2 + k 8 8 b g Z N l t e 1 5 s a z A O / t + P 1 w J 8 x L T K z k 7 B 5 a x t h R c o V M k m N G f h e j s O S a h R M 8 q o R F I b n l F 3 T K z 6 w U N G U m 2 E 5 H W Q F G 5 a J I c m 0 P Q p h y s 5 H l D Q 1 Z p z a K W 3 U P Z r / t Z p 8 S R s U m P w Y l k L l B X L F H g o l h Q T M o N 4 K x E J z h n J s A W V a 2 F 6 B j a i m D O 3 u G n Y I j z + F 1 8 H Z d t f v d X uh i S j 1 G G + y Y 1 f I r j Q h S J u / Q g z f U g t e i B w 7 r n 3 k H u P E z A q l W m + G K m Z 2 b n 5 h f R i Z m l 5 Z X U t u 7 5 R l X 4 o M K l g n / m i 7 i B J G O W k o q h i p B 4 I g j y H k Z r T O 0 v 6 t V s i J P X 5 t e o H p O m h D q c u x U h p q Z W 9 s B U K W 7 a H V F d 4 k e v G 8 B R e / t R 3 M b Q P 4 P 3 N u G 7 H c B 9 q B + / u T o n O X i u b M / P m A H C S W M V D s w i t k Z I D I 5 R b 2 S e 7 7 e P Q I 1 x h h q R s W G a g m h E S i m J G 4 o w d S h I g 3 E M d 0 t C U I 4 / I Z j Q 4 O Y Y 7 W m l D 1 x f 6 c Q U H 6 q Q j Q p 6 U f c / R k 8 m a c r q X i H / 1 G q F y T 5 o R 5 U G o C M f D j 9 y Q Q e X D J D / Y p o J g x f q a I C y o 3 h X i L h I I K 5 1 y R o c w v h T + T 6 q F v F X M F 6 + O c q X C K I 4 0 2 A L b Y B d Y 4 B ig Z m R c F L e F 1 z J V L H 0 D 2 B O R m Y V a t u z Q A X i e 2 c W g 6 0 C 6 U K C r Q G 5 o c 7 j L E K S S Q x Q 0 L 0 b S u R X o a 4 p J i R a d l V g i Q I T 9 C I 9 D W N U E i E l 8 2 O n 8 I j r Q x h E H N d k Y Q z d X E i Q 6 E Q a e j r z h D J s f j t 5 e J f X l / J 4 M L L a J Q o S S I 8 X x Q o B m U M 8 y T g k H K C J U s 1 Q Z h T f S v E Y 8 Q R l j q v s g 7 h + 1 P 4 P + k 0 6 r Z T d 2 7 P q s 1 G E U c J H I B D U A M 2 O A d N c A N a o A 0 w S M E j e A Y v x o P x Z L w a b / P W J a O Y q Y A f M N 6 / A A d K k x A = < / l a t e x i t > u = E(⌧ ) ⌧ < l a t e x i t s h a 1 _ b a s e 6 4 = " E + / U + h U h T z j x O L T X C 9 C r R S 0 z s K Q = " > A A A B 7 n i c d V B N S 8 N A E J 3 U r 1 q / q h 6 9 L B b B U 0 k q R I 8 F E T x W s B / Q h r L Z b N q l mI r R N J J e q F 2 J N O R O 0 b Z j h t J c q i p O Q 0 2 4 4 u Z 7 7 3 Q e q N J P i 3 k x T G i R 4 J F j M C D Z W 6 t 4 I I i O q h t W a W 3 c X Q K v E 8 y 9 c H 3 m F U o M C r W H 1 Y x B J k i V U G M K x 1 n 3 P T U 2 Q Y 2 U Y 4 X R W G W S a p p h M 8 I j 2 L R U 4 o T r I F + f O 0 J l V I h R L Z U s Y t F B X J 3 K c a D 1 N Q t u Z Y D P W v 7 2 5 + J f X z 0 x 8 F e R M p J m h g i w X x R l H R q L 5 7 y h i i h L D p 5 Z g o p i 9 F Z E x V p g Y m 1 D F h v D 9 K f q f d B= " > A A A C H n i c d V D L S g M x F M 3 4 r O O r 6 t J N s A i u y k z F 0 W W h G 5 c V 7 A M 6 Q 8 m k t 2 1 o J j M k m U I Z + i V u / B U 3 L h Q R X O n f m L Z T q K I H A i f n n k N y T 5 h w p r T j f F l r 6 x u b W 9 u F H X t 3 b / / g s H h 0 3 F R x K i k 0 a M x j 2 Q 6 J A s 4 E N D T T H N q J B B K F H F r h q D a b t 8 Y g F Y v F v Z 4 k E E R k I F i f U a K N 1 C 1 e + S Ec V m h L 2 h C m U x u 2 6 G T S Z i Z F E r o j 7 j x V 9 y 4 U M S F G / F v n L 6 g i h 4 Y O J x z L n f u C V P O l H a c L 2 t t f W N z a 7 u w Y + / u 7 R 8 c F o + O m y r J J I U G T X g i 2 y F R w J m A h m a a Q z u V Q O K Q Q y s c 3 k 7 9 1 g i k Y o m o 6 3 E K Q U z 6 g v U Y J d p I 3 e K V H 0 K f i Z y C 0 C A n N s Z 1 I v u g f d / Q J F Q g R 7 O k 7 Y O I l q l u s e S U n R n w K n G 9 S 8 f D 7 k I p o Q V q 3 e KI r R N J J e q F 2 J N O R O 0 b Z j h t J c q i p O Q 0 2 4 4 u Z 7 7 3 Q e q N J P i 3 k x T G i R 4 J F j M C D Z W 6 t 4 I I i O q h t W a W 3 c X Q K v E 8 y 9 c H 3 m F U o M C r W H 1 Y x B J k i V U G M K x 1 n 3 P T U 2 Q Y 2 U Y 4 X R W G W S a p p h M 8 I j 2 L R U 4 o T r I F + f O 0 J l V I h R L Z U s Y t F B X J 3 K c a D 1 N Q t u Z Y D P W v 7 2 5 + J f X z 0 x 8 F e R M p J m h g i w X x R l H R q L 5 7 y h i i h L D p 5 Z g o p i 9 F Z E x V p g Y m 1 D F h v D 9 K f q f d BW K W h h h p J q h S X c d O t J d R q T k T O C 3 1 U o U J Z W M 6 x K 6 h E Q 1 R e d n 8 0 C k 5 M c q A B L E 0 F W k y V 5 c n M h o q N Q l 9 0 x l S P V K / v Z n 4 l 9 d N d X D t Z T x K U o 0 R W y w K U k F 0 T G Z f k w G X y L S Y G E K Z 5 O Z W w k Z U U q Z N N i U T w v e n 5 H / S O q 8 6 b t V t X F Z q b h 5 H E Y 7 g G E 7 B g S u o w S 3 U o Q k M E B 7 h G V 6 s e + v JV o m h B Q 1 H V T X e X n w i u t G 1 / W o W V 1 b X 1 j e J m a W t 7 Z 3 e v v H / Q U n E q G T Z Z L G L Z 8 a l C w S N s a q 4 F d h K J N P Q F t v 3 x z c x v P 6 B U P I 7 u 9 C R B L 6 T D i A e c U W 2 k x l m / X L G r 9 h x k m T j u h e 0 S J 1 c q k K P e L 3 / 0 B j F L Q 4 w 0 E 1 S p r m M n 2 s u o 1 J w J n J Z 6 q c K E s j E d Y t f Q i I a o v G x + 6 J S c G G V A g l i a i j S Z q 8 s T G Q 2 V m o S + 6 Q y p H q n f 3 k z 8 y + u m O r j 2 M h 4 l q c a I L R Y F q S A 6 J r O v y Y B L Z F p M D K F M c n M r Y S M q K d M m m 5 I J 4 f t T 8 j 9 p n V c d t + o 2 L i s 1 N 4 + j C E d w D K f g w B X U 4 B b q 0 A Q G C I / w D C / W v f V k v V p v iC u V B i j R H t Y / g l j i j B N h M E N a D z w 3 N W G O l K G Y k V k t y D R J E Z 6 g E R l Y K h A n O s z n m W f w x C o x T K S y T x g 4 V 5 c 3 c s S 1 n v L I T h Y J 9 W + v E P / y B p l J r s K c i j Q z R O D F o S R j 0 E h Y F A B j q g g 2 b G o J w o0 Q c D j / d m m J k X p o J r 4 7 q f z s r q 2 v r G Z m G r u L 2 z u 7 d f O j h s 6 i R T D B s s E Y l q h 1 S j 4 B I b h h u B 7 V Q h j U O B r X B 0 M / N b D 6 g 0 T + S d G a c Y x H Q g e c Q Z N V a q J 7 1 S 2 a 2 4 c 5 B l 4 v k X r k + 8 X C l D j l q v 9 N H t J y y L U R o m q N Y d z 0 1 N M K H K c C Z w W u x m G l P K R n S A H U s l j V E H k / m h U 3 J q l T 6 J E m V L G j J X l y c m N N Z 6 H I e 2 M 6 Z m q H 9 7 M / E v r 5 O Z 6 D q Y c J l m B i V b L I o y Q U x C Z l + T P l f I j B h b Q p n i 9 l b C h l R R Z m w 2 R R v C 9 6 f k f 9 I 8 r 3 h + x a 9 f l q t + H k c B j u E E z s C D K 6 j C L d S g A Q w Q H u E Z X p0 Q c D j / d m m J k X p o J r 4 7 q f z s r q 2 v r G Z m G r u L 2 z u 7 d f O j h s 6 i R T D B s s E Y l q h 1 S j 4 B I b h h u B 7 V Q h j U O B r X B 0 M / N b D 6 g 0 T + S d G a c Y x H Q g e c Q Z N V a q Z 7 1 S 2 a 2 4 c 5 B l 4 v k X r k + 8 X C l D j l q v 9 N H t J y y L U R o m q N Y d z 0 1 N M K H K c C Z w W u x m G l P K R n S A H U s l j V E H k / m h U 3 J q l T 6 J E m V L G j J X l y c m N N Z 6 H I e 2 M 6 Z m q H 9 7 M / E v r 5 O Z 6 D q Y c J l m B i V b L I o y Q U x C Z l + T P l f I j B h b Q p n i 9 l b C h l R R Z m w 2 R R v C 9 6 f k f 9 I 8 r 3 h + x a 9 f l q t + H k c B j u E E z s C D K 6 j C L d S g A Q w Q H u E Z X pG F 2 0 p u M m X 0 w M y v E J V / g x Y M i X v 0 k b / 6 N k 2 S F K F r Q U F R 1 0 9 3 l J 4 I r b d u f V m F l d W 1 9 o 7 h Z 2 t r e 2 d 0 r 7 x + 0 V J x K h k 0 W i 1 h 2 f K p Q 8 A i b m m u B n U Q i D X 2 Bb X 9 8 P f P b 9 y g V j 6 N b P U n Q C + k w 4 g F n V B u p 8 d A v V + y q P Q d Z J o 5 7 b r v E y Z U K 5 K j 3 y x + 9 Q c z S E C P N B F W q 6 9 i J 9 j I q N W c C p 6 V e q j C h b E y H 2 D U 0 o i E q L 5 s f O i U n R h m Q I J a m I k 3 m 6 v J E R k O l J q F v O k O q R + q 3 N x P / 8 r q p D q 6 8 j E d J q j F i i 0 Y / r c r S 8 s r q W n W 9 t r G 5 t b 1 T 3 9 3 r q i S T h H Z I w h N 5 G 2 B F O Y t p R z P N 6 W 0 q K R Y B p 7 1 g f F n 4 v X s q F U v i G z 1 J q S / w M G Y R I 1 g b y X u 4 8 w T W I y n y c D q o N + y m P Q N a J I 5 7 a r v I K Z U G l G g P 6 h 9 e m J B M 0 F g T j p X q O 3 a q / R x L z Q i n 0 5 q X K Z p i M s Z D 2 j c 0 x o I q P 5 9 l n q I j o4 Q o S q R 5 s U Y z d X E j x 0 K p i Q j M Z J F Q / f Y K 8 S + v n + n o w s 9 Z n G a a x m R + K M o 4 0 g k q C k A h k 5 R o P j E E E 8 l M V k R G W G K i T U 0 1 U 8 L 3 T 9 H / p H v S d N y m e 3 3 W a L l l H V U 4 g E M 4 B g f O o Q V X 0 I Y O E E j h E Z 7 h x cq s J + v V e p u P V q x y Z x 9 + w H r / A r L e k h o = < / l a t e x i t > z d Furthermore, we compare the control performance of our model-based controllers with a baseline control strategy based on the MECH-NODE (nz = 2) that achieves an error of 0.1104 on the test set. First, we utilize the same P-satI-D feedback controller as for the CON model to generate the control action τ (t) in latent space. As the MECH-NODE uses an MLP to parameterize the function ξ = f ξ (ξ, u), we cannot easily map τ (t) into an input u(t). Therefore, we linearize the latent space dynamics w.r.t to the input as f ξ,ac (ξ, u) = f ξ (ξ, 0) + A(ξ) u, where A(ξ) = ∂f ξ ∂u (ξ, 0) is computed using autodiff. Then, u(t) = A T (ξ) τ (t). After tuning the control gains, we choose Kp = 0.001, Ki = 0.02, K d = 1e-5, υ = 1.A A A B + H i c d V D L S g M x F L 3 j s 9 Z H R 1 2 6 C R b B V Z m p M L o s u H F Z w T 6 g H U o m T d v Q J D M k G W E c + i V u X C j i 1 k 9 x 5 9 + Y t i N U 0 Q M X TWe train all models (e.g., MECH-NODE, CON) on three different random seed and choose the best model instance. Please refer to Appendix E for more details on the model selection. We generate a trajectory of 7 setpoints, where q d (tj) ∼ U(-5π, 5π) rad/m ∈ R 2 is a sampled configuration of the soft robot. Then, we render an image o d (tj) that represents the target observation for the controller and encode it into latent space to retrieve z d ∈ R 2 . At time step k, we render an image o(t k ) of the robot's current configuration q(t k ) and encode the image. Subsequently, we evaluate the control law and apply the decoder u(t k ) = η(τ (t k ))), which is finally passed to the simulator that integrates the ground-truth dynamics to the next time-step t k+1 considering the actuation u(t k ). The controller runs at 100 Hz, and we simulate the ground-truth dynamics with a Dopri5 ODE integrator at a time-step of 1e-5 s.Results. As an evaluation metric, we consider the RMSE between the actual and the reference trajectory.The P-satI-D applied to the MECH-NODE model (baseline) achieves an RMSE of 2.88 rad/m w.r.t. to the desired configuration q d (but unknown to the algorithm). The P-satI-D CON controller, which does not exploit the learned latent dynamics for control, exhibits an RMSE of 4.08 rad/m w.r.t. to the desired configuration q d . The P-satI-D+FF controller exhibits an RMSE of 2.12 rad/m w.r.t. to the desired configuration q d . We also visualize the closed-loop trajectories in Fig. 5 and as sequences of stills in Apx. E. We conclude that the nicely structured latent space generated by the β-VAE allows the P-satI-D controller to effectively regulate the system towards the setpoint, although the response time is rather slow. The P-satI-D+FF controller is able to exploit the structure of the CON model through its potential shaping feedforward term. With that, CON P-satI-D+FF exhibits a faster response time and a 26 % lower RMSE than the MECH-NODE P-satI-D baseline. We provide results for the control of an actuated damped harmonic oscillator in Apx. D.Question: Does the paper provide open access to the data and code, with sufficient instructions to faithfully reproduce the main experimental results, as described in supplemental material? Answer: [Yes] Justification: The code associated with this paper is available on GitHub2 . It allows the user to generate the datasets, run the hyperparameters selection, train the models, and generats result plots based on training checkpoints.Guidelines:• The answer NA means that paper does not include experiments requiring code.• Please see the NeurIPS code and data submission guidelines (https://nips.cc/public/ guides/CodeSubmissionPolicy) for more details. • While we encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for a new open-source benchmark).• The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https://nips.cc/ public/guides/CodeSubmissionPolicy) for more details. • The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc. • The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only a subset of experiments are reproducible, they should state which ones are omitted from the script and why. • At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable). • Providing as much information as possible in supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted. 
conclusion and limitations: 
 Conclusion. In this work, we propose a new formulation for a coupled oscillator that is inherently input-tostate stable. Additionally, we identify a closed-form approximation, that is able to simulate the network dynamics more accurately compared to numerical ODE integrators with similar computational costs. When learning latent dynamics with CON, we observe that the performance is on par or slightly better compared to SoA methods such as RNNs, NODEs, etc., even though we constrained the solution space to a ISS-stable coupled oscillator structure. Furthermore, we point out that the performance of the CON models is more consistent across latent dimensions compared to the baselines and improved when not specifically tuned for a given dimension. Furthermore, as seen in Tab. 8, the closed-form approximation achieves, with the same number of model parameters, similar accuracies and double the training speed w.r.t. to the continuous-time model. Finally, we demonstrate that even a simple PID-like latent-space controller can effectively regulate the system to a setpoint. By exploiting the Figure 5 : Latent-space control of a continuum soft robot (simulated using two piecewise constant curvature segments) at following a sequence of setpoints: The upper two rows show the performance of a pure P-satI-D feedback controller operating in latent space z learned with the MECH-NODE and CON models, respectively. The lower row displays the results for a latent space controller based on the CON model that additionally also compensates for the learned potential forces.network structure and compensating for potential forces, regulation performance can be greatly improved, and response time decreased by more than 55 %.Limitations. While we think our proposed method shows great potential and opens interesting avenues for future research, there exist certain limitations. For example, the proposed method of learning (latent) dynamics implicitly assumes that the underlying system adheres to the Markov property (e.g., the full state of the system is observable), that a system with mechanical structure can approximate it, and that it has an isolated, globally asymptotically stable equilibrium. This is, for example, the case for many mechanical systems (e.g., some continuum soft robots, deformable objects, and elastic structures) with continuous dynamics, convex elastic behavior, dissipation, and whose time-dependent effects (e.g., viscoelasticity, hysteresis) are negligible. Even if these conditions are not met globally, the method can be applied to model the local behavior around an asymptotic equilibrium point of the system (e.g., robotic manipulators, legged robots) with added stability benefits for out-of-distribution samples. Alternatively, the method could be extended to relax some of these assumptions, e.g., by allowing for multiple equilibria, zero damping, or by incorporating additional terms to capture discontinuous dynamics (e.g., stick-slip models) or period motions (e.g., limit cycles such as the Van der Pol oscillator). The proposed method might not be suitable for some physical systems, such as nonholonomic systems, partially observable systems, or systems with non-Markovian properties. Examples of such systems include mobile robots and systems with hidden states or delayed observations. Furthermore, the approximated closed-form solution shows the best integration for situations where linear, decoupled dynamics dominate the transient. For dominant nonlinear, coupled forces, the performance of CFA-CON degrades, and it might be better to revert to numerical integration of the CON ODE. Finally, the control works exceptionally well in the setting where the latent dimension equals the input dimension. We hypothesize that this enables the method to identify a diffeomorphism between the input and the latent-space forcing. Still not investigated is how the performance could degrade if nz > m (or nz < m for that matter).Question: Does the paper specify all the training and test details (e.g., data splits, hyperparameters, how they were chosen, type of optimizer, etc.) necessary to understand the results? Answer: [Yes] Justification: In Appendix C, we include implementation details such as the used libraries, algorithms, optimizers, and evaluation procedures. All hyperparameters (e.g., learning rate, loss weights, weight decay, etc.) can be found in the code on GitHub (specifically, in the sweep_generic_dynamics_autoencoder.py Python script that is placed in the examples/sweeping folder). Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in the appendix, or as supplemental material. 
experiment statistical significance: 
 Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments? Answer: [Yes] Justification: We run all experiments with various initializations (i.e., different random seeds), and each result table/plot is accompanied by a description of how the variability of results is captured. Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors). • It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). Justification: This paper primarily involves fundamental research, and the presented application of predicting and controlling the future evolutions of dynamical systems does not directly have any broader societal impact. Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)? Answer: [NA] Justification: This work by itself does not pose any risk for misuse. Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. 
 
 Appendix: In the following, we will write ∥A∥ to denote the induced norm of matrix A and λm(A), λM(A) to refer to its minimum and maximum Eigenvalue respectively.< l a t e x i t s h a 1 _ b a s e 6 4 = " j x C w L a u G K a v r 6 Y e 9 V Q z g i I 7 M x b s = " > A A A C M H i c d V B N S y N B E O 3 x M 8 Z d N + r R S 2 M Q s q C h J y 4 a b 4 I H B T 0 o m C h k w t D T q W h j T 8 / Q X a O E Y X 6 S F 3 + K X h Q U 8 e q v c C a J o M v u g 4 b X r 1 5 R V S + I l b T I 2 J M z M T k 1 P T N b m i v P / / i 5 8 K u y u N S 2 U W I E t E S k I n M e c A t K a m i h R A X n s Q E e B g r O g q u 9 o n 5 2 D c b K S J / i I I Z u y C + 0 7 E v B M Z f 8 y r 4 X A P K a p 9 p g k H o o V Q / S Q e Z 7 I c d L E 6 Y 3 W Q 1 9 9 p t 6 p j C s U 9 w Y f a W m Q 4 v g K j 0 8 y s p + p c r q O 8 1 N 5 j L K 6 m y I g m y 6 r M G o O 1 a q Z I x j v 3 L v 9 S K R h K B R K G 5 t x 2 U x d l N u U A o F W d l L L M R c X P E L 6 O R U 8 x B s N x 0 e n N G 1 X O n R f m T y p 5 E O 1 a 8 d K Q + t H Y R B 7 i y W t H / X C v F f t U 6 C / W Y 3 l T p O E L Q Y D e o n i m J E i / R o T x o Q q A Y 5 4 c L I f F c q L r n h A v O M i x A + L 6 X / J + 1 G 3 d 2 q b 5 3 8 q e 4 2 x n G U y A p Z J T X i k m 2 y S w 7 I M W k R Q W 7 J A 3 k m L 8 6 d 8 + i 8 O m 8 j 6 4 Q z 7 l k m 3 + C 8 f w A 2 / q k A < / l a t e x i t > (kỹ w (t 0 )k, t t 0 ) 2 KL < l a t e x i t s h a 1 _ b a s e 6 4 = " A s S l D 7 J g i w n l t D b v B 3 d Y 2 O l B U x c = " > A A A B / X i c d V D L S s N A F J 3 U V 6 2 v + N i 5 G S y C q 5 J U i S 4 L b l x W s A 9 o Q p h M J u 3 Q y Y O Z i R J D 8 F f c u F D E r f / h z r 9 x 0 k a o o g c G D u f c y z 1 z v I R R I Q 3 j U 6 s t L a + s r t X X G x u b W 9 s 7 + u 5 e X 8 Q p x 6 S H Y x b z o Y c E Y T Q i P U k l I 8 O E E x R 6 j A y 8 6 W X p D 2 4 J F z S O b m S W E C d E 4 4 g G F C O p J F c / s C V l P s m z w r V D J C c 8 z O 8 K V 2 8 a L W M G u E h M 6 9 S w o F k p T V C h 6 + o f t h / j N C S R x A w J M T K N R D o 5 4 p J i R o q G n Q q S I D x F Y z J S N E I h E U 4 + S 1 / A Y 6 X 4 M I i 5 e p G E M 3 V x I 0 e h E F n o q c k y o f j t l e J f 3 i i V w Y W T 0 y h J J Y n w / F C Q M i h j W F Y B f c o J l i x T B G F O V V a I J 4 g j L F V h D V X C 9 0 / h / 6 T f b p l W y 7 o + a 3 b a V R 1 1 c A i O w A k w w T n o g C v Q B T 2 A w T 1 4 B M / g R X v Q n r R X 7 W 0 + W t O q n X 3 w A 9 r 7 F 8 o P l h A = < / l a t e x i t > ỹw (a) GAS: g(u) = 0 < l a t e x i t s h a 1 _ b a s e 6 4 = " r T Q B n 8 9 H P H z i f Y s o G 8 u v 5 Y s y Q w g = " > A A A C a H i c d V F d b 9 M w F H X C 1 + j 4 C F 9 C i B d r F V q q Q Z W U r Y O 3 S b z w O C T a T a q r y H F v U m t 2 E u w b U B V F / E f e + A GN g J k o V W n O U 2 5 B y Q J m K F H B e W W A 6 1 T B W X r x r t P P P o O x s i w + 4 q a C p e Z 5 I T M p O D o q C b 6 y F J C H T M 3 B I G U o 1 Q q a T Z s w z X F t d P O l D T G J R p S Z L u E l x V f b 4 w F l O d e a U 6 Y g Q x p S Z u s q a Z z W M Z 8 o 7 v d 7 2 3 f O w z r E / d F V I 7 f J f I 2 j Q R I M o / E 0 P j w + i m g 0 j r b R g d e T t / E R j X t m S P o 4 T Y J v b F W K W k O B Q n F r F 3 F U 4 b L h B q V Q 0 A 5 Y b a H i 4 o L n s H C w 4 B r s s t k a 1 d I X j l n R r D R u F U i 3 7 O 8 V D d f W b n T q M r v 5 7 d 9 a R / 5 L W 9 S Y v V k 2 s q h q h E J c X p T V i m J J O 9 f p S h o Q q D Y O c G G k e y s V a 2 6 4 Q P c 3 n Q l X k 9 L / g / l k H E / H 0 w + H w 5 N J b 8 c O e U 7 2 S E h i c k x O y H t y S m Z E k B / e r v f Y e + L 9 9 A P / q f / s M t X 3 + p p H 5 I / w 9 3 4 B R o C 2 N Q = = < / l a t e x i t > (kỹ w (t 0 )k, t t 0 ) + ✓ sup t0t 0 t kg(u(t 0 ))k ◆ < l a t e x i t s h a 1 _ b a s e 6 4 = " j x C w L a u G K a v r 6 Y e 9 V Q z g i I 7 M x b s = " > A A A C M H i c d V B N S y N B E O 3 x M 8 Z d N + r R S 2 M Q s q C h J y 4 a b 4 I H B T 0 o m C h k w t D T q W h j T 8 / Q X a O E Y X 6 S F 3 + K X h Q U 8 e q v c C a J o M v u g 4 b X r 1 5 R V S + I l b T I 2 J M z M T k 1 P T N b m i v P / / i 5 8 K u y u N S 2 U W I E t E S k I n M e c A t K a m i h R A X n s Q E e B gZ I x j v 3 L v 9 S K R h K B R K G 5 t x 2 U x d l N u U A o F W d l L L M R c X P E L 6 O R U 8 x B s N x 0 e n N G 1 X O n R f m T y p 5 E O 1 a 8 d K Q + t H Y R B 7 i y W t H / X C v F f t U 6 C / W Y 3 l T p O E L Q Y D e o n i m J E i / R o T x o Q q A Y 5 4 c L I f F c q L r n h A v O M i x A + L 6 X / J + 1 G 3 d 2 q b 5 3 8 q e 4 2 x n G U y A p Z J T X i k m 2 y S w 7 I M W k R Q W 7 J A 3 k m L 8 6 d 8 + i 8 O m 8 j 6 4 Q z 7 l k m 3 + C 8 f w A 2 / q k A < / l a t e x i t > (kỹ w (t 0 )k, t t 0 ) 2 KL < l a t e x i t s h a 1 _ b a s e 6 4 = " A s S l D 7 J g i w n l t D b v B 3 d Y 2 O l B U x c = " > A A A B / X i c d V D L S s N A F J 3 U V 6 2 v + N i 5 G S y C q 5 J U i S 4 L b l x W s A 9 o Q p h M J u 3 Q y Y O Z i R J D 8 F f c u F D E r f / h z r 9 x 0 k a o o g c G D u f c y z 1 z v I R R I Q 3 j U 6 s t L a + s r t X X G x u b W 9 s 7 + u 5 e X 8 Q p x 6 S H Y x b z o Y c E Y T Q i P U k l I 8 O E E x R 6 j A y 8 6 W X p D 2 4 J F z S O b m S W E C d E 4 4 g G F C O p J F c / s C V l P s m z w r V D J C c 8 z O 8 K V 2 8 a L W M G u E h M 6 9 S w o F k p T V C h 6 + o f t h / j N C S R x A w J M T K N R D o 5 4 p J i R o q G n Q q S I D x F Y z J S N E I h E U 4 + S 1 / A Y 6 X 4 M I i 5 e p G E M 3 V x I 0 e h E F n o q c k y o f j t l e J f 3 i i V w Y W T 0 y h J J Y n w / F C Q M i h j W F Y B f c o J l i x T B G F O V V a I J 4 g j L F V h D V X C 9 0 / h / 6 T f b p l W y 7 o + a 3 b a V R 1 1 c A i O w A k w w T n o g C v Q B T 2 A w T 1 4 B M / g R X v Q n r R X 7 W 0 + W t O q n X 3 w A 9 r 7 F 8 o P l h A = < / l a t e x i t > ỹw (b) ISSFigure 6 : Illustration of global asymptotic stability for the unforced system with g(u) = 0 and input-to-state stability for the forced system, where the black dashed line denotes the input u(t).We introduce some expressions often used throughout this section: The gradient of Vµ(ỹw) w.r.t. the residual coordinate ỹw is given byEQUATION)Next, the Hessian of the Lyapunov candidate can be derived asEQUATIONwhereEQUATION) Furthermore, the Schur complement of PV is given byEQUATIONA.1 Potential force and energy EQUATION)Proof. First, we take the derivative of U(xw):EQUATION)The Hessian of the potential is given byEQUATION)As Kw ≻ 0 ⇒ Kw = K T w , we can easily show that the potential force is symmetric:EQUATIONA. EQUATIONBased on the assumption Mw ≻ 0, Kw ≻ 0, we can state ∥Mw∥ λw(Kw) > 0. Therefore, the critical case for the lower bound on λm(SPEQUATIONConsequently, the Eigenvalue sensitivity theorem [46] demands that SP V ≻ 0. Proof. By inspecting the expressions for PV ≻ 0 and HV(xw) ≻ 0 ∀xw ∈ R n in Equations ( 5) and ( 15), respectively, it can be easily seen that HV(xw) ⪰ PV ∀xw ∈ R. As Lemma 4 states that the Schur complement of PV is positive definite, it follows that HV(xw) ⪰ PV ≻ 0. We strive to find the critical points (i.e., minimas and maximas) r of hV,th (r) and, for this, analyze where the first derivative of hV,th (r) is zeroEQUATIONwhich is the case only for r = 0. Next, we compute the second derivative at r asEQUATIONThus, hV,th (r) is convex and its global minimum at r = 0 takes the value h V,th (0) = 0. As a result, h V,th (xw) is also positive semi-definite. Proof.Step 1: It can be easily seen that Vµ(ỹw) in ( 5) is smooth and continuously differentiable.Step 2: Proof that Vµ(0) = 0.EQUATION)Step 3: Proof that the Lyapunov candidate is positive definite; i.e., Vµ(ỹw) > 0 ∀ỹw ∈ R n \ {0}.As the gradient of the Lyapunov candidate, as defined in (14) , is zero for ỹw = 0:EQUATIONỹw = 0 is a critical point of Vµ(ỹw). According to Lemma 5, the Hessian in ( 15) is positive-definite [56] :HV(ỹw) ≻ 0 ∀ỹw ∈ R 2n .With that, ( 5) is convex and its global minimum is at ỹw = 0, where Vµ(0) = 0. In summary, we state Vµ(ỹw) > 0 ∀ỹw ∈ R n \ {0}.Step 4: Proof that the Lyapunov candidate is radially unbounded: i.e., ∥ỹw∥ → ∞ ⇒ Vµ(ỹw) → ∞. Lemma 6 is exploited for identifying a lower bound on Vµ(ỹw):EQUATION)Lemma 5 tells us that PV ≻ 0 and with that λm(PV) > 0. Therefore, if ∥ỹw∥ → ∞, it also follows that Vµ(ỹw) → ∞.The Lemmas introduced below are used to prove Theorem 1.Lemma 7. Suppose xw, xw, b ∈ R n and n ∈ N + . Then, the function h V,th (xw) defined asEQUATIONis positive semi-definite.Proof. Proving h V,th (xw) ≥ 0 is equivalent to proving that the scalar function h V,th (r) = (tanh(r + a) -tanh(a)) r ≥ 0 ∀ r, a ∈ R, where we set r = xw,i and a = xw,i + bi. Now, we expand the hyperbolic tangent:EQUATIONas the denominator e 2r+4a + e 2r+2a + e 2a + 1 > 0 ∀ r ∈ R and as sign 2 e 2a e 2r -1 = sign(r). For example, e 2r -1 ≥ 0 ∀ r ≥ 0. Analog, e 2r -1 < 0 ∀ r < 0.Lemma 8. Let Mw ≻ 0, Kw ≻ 0, and Dw ≻ 0. Also, let µ ∈ R + be chosen such that 0 < µ < Proof. The Schur complement of P V is given byEQUATIONThe lower bound on the smallest Eigenvalue of SP V can be identified asEQUATIONAs Kw, Dw ≻ 0, we know that ∥Dw∥ 2 λm(Kw) > 0. Therefore, the case µ = λm(Dw)λm(Mw)+ ∥Dw ∥ 2 4λm (Kw ) := µ V determinesthe lower bound on λm(SP V ):EQUATIONWe conclude, based on the Eigenvalue sensitivity theorem of symmetric matrices [46] , that SP V ≻ 0 and with that P V ≻ 0 [56] .A.5 Proof of Theorem 2: Proof of Input-to-State Stability (ISS)EQUATIONProof. Proving h V,th (xw) ≤ 2 |xw| is equivalent to proving that the scalar function hV,th (r) = Next, we demonstrate that the slope of 2 |r| is always larger than the magnitude of the slope of hV,th (r):EQUATIONAdditionally, hV,th (0) = 2 |0| = 0. We conclude that hV,th (r) ≤ 2|r| ∀ r ∈ R and with that,h V,th (xw) ≤ 2 |xw| ∀ xw ∈ R n .Lemma 10. Let Mw ≻ 0 and Kw ≻ 0. Then, (5) is bounded by the two scalar, class K∞ functions α1(r) = 1 2 λm(PV) r 2 and α2(r) = 1 2 λM(PV) r 2 + 2√ n r: α1(∥ỹw∥ 2 2 ) ≤ Vµ(ỹw) ≤ α2(∥ỹw∥) 2 2 .Proof. With Lemma 2, we already showed that Vµ(ỹw) is a Lyapunov candidate. Now, we additionally also verify the conditions for ISS-Lyapunov candidates [40] .Step 1: Establishing bounds on Vµ(ỹw).We first identify the lower bound of Vµ(ỹw) by leveraging Lemma 6:EQUATIONSimilarly, we derive an upper bound for Vµ(ỹw) exploiting Lemma 9.EQUATIONStep 2: Proof that α1(r), α2(r) belong to class K∞. According to Lemma 5, PV ≻ 0 and with that λm(PV) > 0. First, we analyze the behavior of α1(r): as it is strictly increasing and α1(0) = 0, it belongs to class K. Furthermore, we can evaluate limr→∞ α1(r) = ∞. Therefore, α1(r) ∈ K∞ [40] . α2(r) is also strictly increasing for r ∈ [0, ∞), α2(0) = 0, and it is radially unbounded as limr→∞ α2(r) = ∞. For that reason, α2(r) ∈ K∞ as well.B Appendix on an approximate closed-form solution for coupled oscillator networks As introduced in (1), we consider the linear dynamics of a 1D forced harmonic oscillator with state yi =EQUATIONwhere Fi(t) ∈ R is the externally applied force acting on the oscillator.The characteristic equation for the unforced dynamics (i.e., Fi(t) = 0) can be stated as [48]EQUATIONwhere ωn,i = √ κi and ζi = d i 2 √ κ i are the natural frequency and the damping factor of the ith homogeneous oscillator, respectively. This harmonic oscillator exhibits three regimes: underdamped (ζi < 1), critically damped (ζi = 1), and overdamped regime (ζi > 1).We approximate the forcing using the Heavyside function H(t): Fi(t) = Fi(t k ) H(t), where Fi(t k ) is the constant external forcing as computed by (10) . The solution for ζi ̸ = 1 is given by [48] EQUATIONwhere δt = t k+1 -t k , αi = ζi ωn,i, and βi = ωn,i 1 -ζ 2 i . After enforcing the initial conditions xi(t k ), xi(t k ), the integration constantsEQUATIONcan be identified with ∆λi = λi,2 -λi,1 = -2 βi j, where j is the imaginary value. While we could derive the solution for the critically damped case di < 2 √ κi separately, we instead approximate it in our network dynamics with (51) by setting ∆λi ∼ = sign(-2 βi j) ϵ |2 βi j| < ϵ -2 βi j |2 βi j| ≥ ϵ , where ϵ ∈ R + ≪ 1 is a small, positive value.We can now leverage the closed form solution to the evolution of a single, decoupled damped harmonic oscillator of (51) to solve the integral in ( 11)EQUATIONwithκ = diag(K11 . . . Knn), d = (D11 . . . Dnn), ωn = √ κ, ζ = d 2 √ κ , α = ζ ⊙ ωn, β = ωn 1 -ζ 2 , F (t k ) = g(u(t k )) -(K -κ) x(t k ) -(D -d) ẋ(t k ) -tanh (W x(t k ) + b) , c1 = x(t k ) - F (t k ) κ , c2 = 1 β ẋ(t k ) + α ⊙ x(t k ) - F (t k ) κ .(54) We summarize the approach of integrating/rolling out the CFA-CON dynamics in Algorithm 1.Algorithm 1 Rollout of CFA-CON.Inputs: initial state y(t 0 ), input sequence {u(t 0 ), . . . u(t k ), . . . u(t N )} Outputs: state sequence {y(t 0 ), . . .y(t k ), . . . y(t N )} 1: k ← 0 2: while k ≤ N do 3: (x(t k ), ẋ(t k )) ← y(t k ) 4: F (t k ) ← g(u(t k )) -(K -κ) x(t k ) -(D -d) ẋ(t k ) -tanh (W x(t k ) + b) 5: ω n , ζ ← √ κ, d 2 √ κ▷ Compute the characteristics of the decoupled harmonic oscillators.6:α, β ← ζ ⊙ ω n , ω n 1 -ζ 2 7: c 1 ← x(t k ) -F (t k )κ ▷ Compute integration constants using initial conditions.8:c 2 ← 1 β ẋ(t k ) + α ⊙ x(t k ) -F (t k ) κ 9: δt = t k+1 -t k ▷ Set time step. 10:▷ Update state with approximated closed-form solution.11:x(t k+1 ) ← (c 1 ⊙ cos(β δt) + c 2 ⊙ sin(β δt)) ⊙ e -α δt + F κ 12: ẋ(t k+1 ) ← -((c 1 ⊙ α -c 2 ⊙ β) cos(β δt) + (c 1 ⊙ β + c 2 ⊙ α) sin(β δt)) ⊙ e -α δtk ← k + 1 ▷ Update time index. 14: end while B.3 Approximation bounds for CFA-CON Lemma 11 demonstrates how, for the particular case of no external input and linearly decoupled oscillators (which we are always free to choose), we can establish bounds on the approximation error when using the closed-form solution instead of the ground-truth coupled oscillator dynamics. Lemma 11. Suppose that the network is unforced with g(u(t)) = 0 and that K = diag(κ1, . . . , κn), D = diag(d1, . . . , dn) such that the oscillators are not linearly coupled. Then, given any t ≥ 0, and the initial state y(0) ∈ R 2n , the error between the continuous dynamics ẍ(t) of (2) and the approximated dynamics x(t) in (10) , is bounded by ∥ẍ(t) -x(t)∥ ≤ 2.Proof. (2), (10) and F = -f ẍ,nld (y(0), 0) give usEQUATIONIn Table 3 , we present a comparison of CFA-CON with several other strategies for integrating nonlinear dynamics, such as CON. Following the implicit assumption made in Section B, we consider the case of g(u) = 0, K = diag(κ1, . . . , κn) and D = diag(d1, . . . , dn) but with the hyperbolic coupling between the oscillators active (i.e., a full W matrix). Integrating the dynamics at a very small time step (i.e., δt = 5e-5 s) with a high-order ODE solver would give us a very accurate solution, but this is computationally infeasible in practice. We, therefore, regard this as the upper bound on the accuracy of the solution. A feasible solution would be to implement either a high-order solver such as Tsit5 at a larger integration time-step, e.g., δt = 1e-1 s) or a low-order solver with a slightly smaller integration time step, e.g., δt = 5e-2 s). Therefore, we also benchmark these options. We also benchmark an implementation specialized on the underdamped case (i.e., ζi < 1): Closed-Form Approximation of the Underdamped Coupled Oscillator Network (CFA-UDCON). This specialized implementation allows us to avoid using complex numbers in the algorithm and reduces the number of computations necessary for calculating the approximated solution. As a result, we see a considerable increase in the sim-time to real-time factor.We perform the integration error benchmark over 100 different network configurations, all consisting of 50 oscillators (n = 50): First, we sample the natural frequency of the ith oscillator from a uniform distribution as ωn,i ∼ U(0.05 Hz, 0.5 Hz), then we sample κi ∼ U(0.2 N/m, 2 N/m) such that K = diag(κ1, . . . , κn) ≻ 0, which lets us determine each mass mi= κ i ω 2 n,i. Next, the damping ratio is determined as ζi ∼ U(0.1, 0.9) and ζi ∼ U (0.1, 2.0) for the underdamped and general case, respectively. As a result, D = diag(d1, . . . , dn) ≻ 0 with di = 2 ζi √ mi κi given. Finally, by leveraging the Cholesky decomposition, we sample a W ≻ 0 and bi ∼ U (-1, 1). We compute the estimation error of all integrated trajectories with respect to the high-precision solution (i.e., Tsitouras' 5/4 method (Tsit5) at δt = 5e-4 s). For this, we compute the RMSE for each 60 s trajectory and then take the mean and standard deviation across the 100 different network configurations.The simulation vs. real-time factor is computed as the simulated rollout duration per second of computational time. For this, we let each method simulate a 60 s trajectory for 100 times and record the minimum run time on an Intel Core i7-10870H CPU (single core) over 10 trials. Because of computational constraints, we simulated with the high-precision Tsit5 solver the trajectory only 5 times.The results in Table 3 show that CFA-CON is 30 % more accurate than the Euler integrator at half of the speed. Compared against the Tsit5 integrator, CFA-CON exhibits a 1.56x speed increase while being significantly less accurate. For the underdamped case with ζ < 1, the specialized implementation CFA-UDCON is 14.8 % faster and at the same time 32 % more accurate than the Euler integrator. Furthermore, CFA-UDCON is 3.7x faster and significantly less accurate than the Tsit5 integrator. We can conclude that in the pure rollout setting (i.e., no backpropagation involved) for a generic CON, the CFA-CON does not show clear advantages to an appropriately tuned Euler or Tsit5 solver. However, the specialized version CFA-UDCON demonstrates a 2.4x speed-up at no reduction of accuracy vs. CFA-CON for underdamped oscillator networks.Method RMSE [m] ↓ RMSE ζ < 1 [m] ↓ Complexity ↓ Sim. time Real time ↑ CON with Tsit5 at δt = 5e-5 s n/a n/a O n log 2 7 ph δt = O(3.5e11) 5.68xCON with Tsit5 at δt = 1e-1 s 5e-5 ± 1e-5 8e-6 ± 1e-5O n n log 2 7 ph δt = O(1.8e8) 11310xCON with Euler at δt = 5e-2 s 0.010 ± 0.003 0.022 ± 0.005EQUATIONCFA-CON (our) with δt = 1e-1 s 0.007 ± 0.002 0.015 ± 0.003O n log 2 7 h δt = O(3.5e7) 17680x CFA-UDCON (our) with δt = 1e-1 s n/a 0.015 ± 0.003 O n log 2 7 h δt = O(3.5e7) 41900xTable 3 : Benchmarking of various methods for integrating the CON dynamics. The RMSE is computed with respect to the Tsitouras' 5/4 method (Tsit5) (i.e., extremely high accuracy but also extremely high computational complexity). We denote with n the number of oscillators in the network (in this case n = 50), with p the order of the numerical ODE solver, and with δt the time-step. When stating the complexity, we refer to h = tN -t0 as the rollout horizon in seconds. In this case, we report the results for a horizon of h = 60 s. The RMSE column states the RMSE of the various integration strategies with respect to the CON with Tsit5 at δt = 5e-5s solution, which we consider to be the ground-truth. The RMSE ζ < 1 computes the same metrics, but this time for a dataset that contains only underdamped oscillators. The Sim. time Real time column states the ratio between the duration of the simulation achieved (in seconds) per second of real-time (i.e., computational time). We report the mean and standard deviation of the RMSE over 100 different network configurations.coRNN. In the 1 st -order version, we define the transition function asEQUATIONCFA-CON. We adapt a 1 st -order version of Alg. 1 for predicting the time evolution in latent-spaceEQUATIONwhere the closed-form solution for the integral is given byEQUATIONFor 2 nd -order systems and when integrating the evolution of the latent state ξ(t) = z T (t) żT (t)T in time,we need to have access to an initial latent velocity ż(t0) such that we can roll out the latent state ξ(t) in time. A naive approach to estimating such an initial latent velocity would be to encode multiple (at least two) images of the system at the start of the trajectory into latent space and then perform numerical differentiation (e.g., finite differences) in latent space. However, we found the resulting ż(t k ) to be relatively noisy and susceptible to small encoding errors. Instead, we propose to perform numerical differentiation in image space and then map this velocity into latent space using the encoder's Jacobian. First, we estimate the image-space velocity at t k using finite differences:ȯ(t k ) ≈ o(t k+1 )-o(t k-1 ) t k+1 -t k-1. The latent velocity is then estimated as ż(t k ) = ∂Φ ∂o (o(t k )) ȯ(t k ), where ∂Φ ∂o is obtained with forward-mode automatic differentiation.We implement the network dynamics and the neural networks (e.g., encoder, decoder, and MLPs) in JAX [60] and Flax [68] , respectively. When training or inferring time-continuous dynamical models (e.g., NeuralODE, CON), we rely on Diffrax [47] for numerical integration of the ODE using the Dormand-Prince's 5/4 method [69] (Dopri5). For the numerical integration of both the time-continuous and the time-discrete models (e.g., RNN, coRNN, CFA-CON), we use an integration time-step δt of 0.025 s and 0.01 s for the Toy Physics [26] and soft robotic datasets, respectively.Because of the GPU memory constraints, we limit ourselves to a batch size of 30 and 80 trajectories for the Toy Physics [26] and soft robotic datasets, respectively. We implement a learning rate schedule consisting of a warm-up (5 epochs) and a cosine annealing [70] period (remaining epochs). We employ an AdamW optimizer [71, 72] with β1 = 0.9, β2 = 0.999 for updating both neural network weights (e.g., encoder, decoder) and parameters of the dynamical model (e.g., K, Dw, Mw, etc.).Before training, we conduct a hyperparameter selection study using Optuna [53] . For this, we leverage a Tree-Structured Parzen Estimator [73] for identifying hyperparameters such as the base learning rate, the weight decay, the loss function weights, and model-specific hyperparameters such as the number of MLP layers, the hidden dimension of the MLP layers, the γ and ϵ values for the coRNN model etc. that minimize the RMSE of the predicted images. To reduce computational requirements, we employ the Asynchronous Successive Halving Algorithm [74] to stop unpromising trials early.Similar to other publications in the field [75, 76, 77] , we state the RMSE, the Peak Signal-to-Noise Ratio (PSNR) and the Structural Similarity Index Measure (SSIM) [78] between the ground-truth image o ∈ R ho×wo and the predicted image image ô ∈ R ho×wo . We use the separated test set for all evaluation results.The RMSE between the two images is given by RMSE(o, ô) = We report the full set of quantitative results, including the additional evaluation metrics PSNR and SSIM, in Tables 7 to 9 . For PCC-NS-2 in Tab. 8, we additionally also recorded the training steps per second on an Nvidia RTX 3090 GPU with a batch size of 80 (which leads to 8080 images per batch). We plot the results of a sweep across the latent dimensions for the additional evaluation metrics PSNR and SSIM in Fig. 8 . Correspondingly, we visualize the number of trainable parameters of each model vs. the latent dimension in Fig. 9 . In Figs. 10 11 12 13 14 15 16 , we present sequences of stills for the rollout of the trained CON(-M) models on the various datasets. Table 8 : Benchmarking of CON and CFA-CON at learning latent dynamics on the PCC-NS-2 (soft robot with two constant curvature segments) dataset. For all models, a latent dimension of nz = 8 is chosen. CON-S and CON-M are small and medium-sized versions of the CON model, respectively. MECH-NODE is a NODE with prior knowledge about the mechanical structure of the system (i.e., dx dt = ẋ). We report the mean and standard deviation over three different random seeds. Furthermore, we state the number of parameters of each latent dynamics model and the training steps per second on a Nvidia RTX 3090 GPU. Each batch contains 80 trajectories and 8080 images of resolution 32x32px in total. Finally, we report the inference time averaged over 5000 runs for performing a rollout of 2.02 s (while encoding and decoding all images along the trajectory) on an Nvidia RTX 3090 GPU with a batch size of 1.As all previous examples exampled ODEs, we strive to test the proposed approach also on a system that is governed by PDEs. Specifically, we consider the Reaction-diffusion (R-D) dataset as introduced in Apx. C.1.3.To address the unactuated nature of the dataset, we remove, analog to the M-SP+F, S-P+F, and D-P+F datasets, the input-to-state mapping parameters of the dynamical models (e.g., the B(u) and E(τ ) MLPs for the CON models). Furthermore, the PDE describing the system dynamics is of 1 st -order. Therefore, we leverage the 1 st -order versions of the latent dynamics as specified in Apx. C.4.We report the metrics of the test set evaluations in Tab. 10 . Furthermore, we also present a sequence of stills of the rollout of a trained latent dynamics CON model in Fig. 17 . We find it impressive that CON with its strong stability guarantees can accurately model the dynamics of a high-dimensional PDE system. dataset. For all models, a latent dimension of nz = 4 is chosen. As this dataset does not consider inputs, we remove all parameters in the RNN, GRU, coRNN, CON, and CFA-CON models related to the input mapping. Also, as the reaction-diffusion system is governed by 1 st -order PDE dynamics, we use specialized, 1 st -order version of the CON, CFA-CON, and coRNN dynamics. We report the mean and standard deviation over three different random seeds and the number of parameters of each latent dynamics model. Target q d P-satI-D q D+FF q P-satI-D+FF q (a) Configuration q(t) ∈ R 2 The orange line visualizes the behavior of a feedforward controller with only a damping term applied in feedback (D+FF). The green line shows the performance of our proposed combination of feedback and feedforward terms (P-satI-D+FF). The dotted and solid lines show the reference and actual values, respectively. For each setpoint, we randomly sample a desired shape q d and render the corresponding image o d . This image is then encoded to a target latent z d . The controller then computes a latent-space torque F d , which is decoded to an input u. Finally, we provide this input to the simulator, which performs a roll-out of the closed-loop dynamics. Important: The robot's configuration (i.e., the first-principle, minimal-order state) is solely used for generating a target image and simulating the closed-loop system.For the control experiments, we train instances of the MECH-NODE and CON-M models with latent dimension nz = 2 and with neural network weights initialized with three different random seeds. For MECH-NODE, we choose the model with the lowest validation loss (seed 0).For the CON network, we found that model-based control does not perform as well when the latent stiffness Γw (as visualized in Fig. 21(a) ) is significantly larger along one of the Eigenvectors than along the other one. Therefore, we evaluate the Eigenvalues of the learned stiffness matrix in W-coordinates after training: λ1,2(Γw). Particularly, we choose the seed that minimizes the normalized standard deviation of the EigenvaluesEQUATIONAdditional results for the P-satI-D feedback controller based on the MECH-NODE and CON models are provided in Fig. 23 , Fig. 24 , respectively and for the P-satI-D+FF controller based on the CON model in Fig. 25 . Sequences of stills for the CON P-satI-D+FF controller are provided in Fig. 22 . 23 : Latent-space control of a continuum soft robot (simulated using two PCC segments) following a sequence of setpoints with a pure P-satI-D feedback controller operating in a 2D latent space learned with the MECH-NODE model. The CON model weights are initialized using a random seed of 0. The dotted and solid lines show the reference and actual values, respectively. For each setpoint, we randomly sample a desired shape q d and render the corresponding image o d . This image is then encoded to a target latent z d . The controller then computes a latent-space torque F d , which is decoded to an input u. Finally, we provide this input to the simulator, which performs a roll-out of the closed-loop dynamics. Important: The robot's configuration (i.e., the first-principle, minimal-order state) is solely used for generating a target image and simulating the closed-loop system. Latent-space control of a continuum soft robot (simulated using two PCC segments) following a sequence of setpoints with a pure P-satI-D feedback controller operating in a 2D latent space learned with the CON model. The CON model weights are initialized using a random seed of 0. The dotted and solid lines show the reference and actual values, respectively. For each setpoint, we randomly sample a desired shape q d and render the corresponding image o d . This image is then encoded to a target latent z d . The controller then computes a latent-space torque F d , which is decoded to an input u. Finally, we provide this input to the simulator, which performs a roll-out of the closed-loop dynamics. Important: The robot's configuration (i.e., the first-principle, minimal-order state) is solely used for generating a target image and simulating the closed-loop system. : Latent-space control of a continuum soft robot (simulated using two PCC segments) following a sequence of setpoints with a pure P-satI-D+FF feedback & feedforward controller operating in a 2D latent space learned with the CON model. The CON model weights are initialized using a random seed of 0. The dotted and solid lines show the reference and actual values, respectively. For each setpoint, we randomly sample a desired shape q d and render the corresponding image o d . This image is then encoded to a target latent z d . The controller then computes a latent-space torque F d , which is decoded to an input u. Finally, we provide this input to the simulator, which performs a roll-out of the closed-loop dynamics. Important: The robot's configuration (i.e., the first-principle, minimal-order state) is solely used for generating a target image and simulating the closed-loop system.F Extended discussion on future applications and limitations F.1 Systems for which we would expect the proposed method to work Mechanical systems with continuous dynamics, dissipation, and a single, attractive equilibrium point. The proposed method is a very good fit for mechanical systems with continuous dynamics, dissipation, and a single, attractive equilibrium point. In this case, the real system and the latent dynamics share the energetic structure and stability guarantees. Examples of such systems include many soft robots, deformable objects with dominant elastic behavior, and other mechanical structures with elasticity.Local modeling of (mechanical) systems that do not meet the global assumptions. Even if the global assumptions of the proposed method are not met, the method can still be applied to model the local behavior around a local asymptotic equilibrium point of the system (i.e., in the case of multi-stability). For example, the method could be used to model the behavior of a robotic leg locally in contact with the ground, a cobot's interaction with its environment, etc.Mechanical systems without dissipation. The proposed method would currently not work well for mechanical systems without any dissipation, as (a) the original system will likely not have a globally asymptotically stable equilibrium point, and more importantly, (b) we currently force the damping learned in latent space to be positive definite. However, these systems are not common in practice as friction and other dissipation mechanisms are omnipresent, and the proposed method can learn very small damping values (e.g., the mass-spring+friction system). A possible remedy could be to relax the positive definiteness of the damping matrix in the latent space, allowing for zero damping. This would allow the method to work for systems without dissipation, such as conservative systems. Examples of such systems include a mass-spring system without damping, the n-body problem, etc.Systems with discontinuous dynamics. The proposed method might underperform for systems with highly discontinuous dynamics, such as systems with impacts, friction, or other discontinuities. In these cases, the latent dynamics might not capture the real system's behavior accurately, and the control performance of feedforward + feedback will very likely be worse than pure feedback. Again, the method should be able to capture local behavior well. A possible remedy for learning global dynamics could be to augment the latent dynamics with additional terms that capture the discontinuities, such as contact and friction models (e.g., stick-slip friction).Systems with multiple equilibrium points. The original system having multiple equilibria conflicts with the stability assumptions underlying the proposed CON latent dynamics. In this case, as, for example, seen on the pendulum+friction and double pendulum + friction results, the method might work locally but will not be able to capture the global behavior of the system. A possible remedy could be to relax the global stability assumptions of the CON network. For example, the latent dynamics could be learned in the original coordinates of CON while allowing W also to be negative definite. This would allow the system to have multiple equilibria & attractors. Examples of such systems include a robotic arm under gravity, pendula under gravity, etc.Systems with periodic behavior. The proposed method will likely not work well for systems with periodic behavior, as they do not have a single, attractive equilibrium point. Examples of such systems include a massspring system with a periodic external force, a pendulum with a periodic external force, some chemical reactions, etc. Again, it is likely possible to apply the presented method to learning a local behavior (i.e., not completing the full orbit). A possible remedy could be to augment the latent dynamics with additional terms that capture the periodic behavior, such as substituting the harmonic oscillators with Van der Pol oscillators to establish a limit cycle or a supercritical Hopf bifurcation.F.3 Systems for which we would not expect the proposed method to work Nonholonomic systems. The proposed method likely would not work well for nonholonomic systems, as both structure (e.g., physical constraints) and stability characteristics would not be shared between the real system and the latent dynamics. Examples of such systems include vehicles, a ball rolling on a surface, and many mobile robots.Partially observable and non-markovian systems. As the CON dynamics are evaluated based on the latent position and velocity encoded by the observation of the current time step and the observation-space velocity, we implicitly assume that the system is (a) fully observable and (b) fulfills the Markov property. This assumption might not hold for partially observable systems, such as systems with hidden states or systems with delayed observations. Examples of such cases include settings where the system is partially occluded or in situations without sufficient (camera) perspectives covering the system. Furthermore, time-dependent material properties, such as viscoelasticity or hysteresis, that are present and significant in some soft robots and deformable objects are not captured by the method in its current formulation.• All the theorems, formulas, and proofs in the paper should be numbered and cross-referenced.• All assumptions should be clearly stated or referenced in the statement of any theorems.• The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch to provide intuition. • Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material. • Theorems and Lemmas that the proof relies upon should be properly referenced.• The answer NA means that the paper does not include experiments.• If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not. • If the contribution is a dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable. • Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed. • While NeurIPS does not require releasing code, the conference does require all submissions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example (a) If the contribution is primarily a new algorithm, the paper should make it clear how to reproduce that algorithm. (b) If the contribution is primarily a new model architecture, the paper should describe the architecture clearly and fully. (c) If the contribution is a new model (e.g., a large language model), then there should either be a way to access this model for reproducing the results or a way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset). (d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results.https://github.com/tud-phi/uncovering-iss-coupled-oscillator-networks-from-pixelshttps://github.com/tud-phi/uncovering-iss-coupled-oscillator-networks-from-pixels 


Decision: The paper presents a very interesting combination of traditional control methods and machine learning, applied to robotic control (non-soft robot results were added at rebuttal time). Unfortunately there were only 2 reviews but both are very enthusiastic and I concur with their assessment. The authors should make sure that the final version reflects the clarifications provided in the rebuttal.

Reviews: review: paper_summary: paper_summary: - This paper proposes Coupled Oscillator Networks (CONs); Networks consisting of coupled one-dimensional dampened harmonic oscillators, coupled through a neuron-like connection.
- It is shown that under some constraints, the unforced coupled oscillators have a single, globally stable equilibrium and the forced system is globally ISS stable.
- As the CON has no closed-form solution, the authors propose to split the solution in one decoupled linear, analytically solvable term and a residual non-linear term which is. For this system, a closed-form approximation is derived (CON-CFA).
- High dimensional observations are encoded via a VAE into a compressed latent space that is trained to reconstruct future observations that are predicted via a CON latent dynamics model. This system is trained jointly.
- “CONs are an ideal fit for learning latent dynamics as they guarantee that the latent states stay bounded.”
- The approaches (CON with two different latent dimensions,S and M, and CON-CFA) are evaluated on a simulated soft-robots dataset and benchmarked against other latent-space dynamics methods based on neural ODEs ore autoregressive models.
- When predicting system state, CON-M performs on par with sota methods, CON-S and CON-CFA perform almost as good, while cheaper computationally.
- For control, potential-shaping strategies are combined with PID, resulting in a latent-space control law. A forcing decoder is then trained to to predict control inputs. This approach strongly outperforms a simple latent space PID controller.
strength_weakness: strength_weakness: - Limited experimental evaluation on a single simulated soft robot.
- It is unclear how the obtained results generalize to other robots, in particular with less smooth contact dynamics.
questions: questions: - Have you tried the method on other (non-soft) robots? Do you have any intuition on how your method would perform on environments with many contact forces?
- Have you considered non-physical systems?
limitations: limitations: Many limitations of the method are clearly stated.
paper_summary: paper_summary: In this paper, the authors introduce Coupled Oscillator Networks (CONs), extending the idea of coRNNs but with the difference that the generalized force is symmetric, enabling the definition of a potential energy expression which can be exploited for energy shaping based control methods. It is proven, that this model is asymptotically stable and, given an input, input-to-state-stable. It is shown that CONs provide similar performance compared to other state of the art methods. Then, the proposed network is used together with a variational autoencoder to learn a low-dimensional representation of the dynamics of a soft robot. Finally, a PID controller with an energy shaping part is designed to steer the soft robot to a given desired state. The main contribution is the new network structure which allows us to exploit the potential energy expression.
strength_weakness: strength_weakness: - Theorem 1 seems to be trivial as it elaborates on the (nonlinear) coupling of passive systems (damped mass-spring systems) which is always passive since it does not “produce” any energy and thus, is inherently stable. (Same for Theorem 2). In this sense, the theoretical contribution seems limited. Is there anything I missed?
- The statement in Line 266 “PID controller has several well-known drawbacks, such […] steady-state errors (in case the integral gain is chosen to be zero)” is misleading as in this case we would say it’s a PD controller.
- The main contribution seems to be the new network structure that allows to exploit the potential energy for energy shaping methods. The VAE and controller are existing methods. Thus, a more detailed elaboration on the performance and limitation of CONs would be beneficial as it is evaluated for soft robotics data sets only.
questions: questions: - Can you add the inference time for the different methods in Table 1? If that’s not possible, could you give some general comments on the inference time of CONs?
- Figure 3 visualizes that the controller with feed-forward part leads to heavy oscillations in the systems. Is that due to a poorly tuned controller or do you see the reason in the CON model?

UPDATE: New score after discussion
limitations: limitations: Yes

Comments: title: Response to Reviewer W9L3 Part 3
Response to Reviewer W9L3 Part 2
Response to Reviewer W9L3 Part 1
Technical response to Reviewer Ys99: Example 2
Technical response to Reviewer Ys99: Example 1
General response to Reviewer Ys99
comment: ### Systems for which we would not expect the proposed method to work

- **Nonholonomic systems:** The proposed method likely would not work well for nonholonomic systems, as both structure (e.g., physical constraints) and stability characteristics would not be shared between the real system and the latent dynamics. Examples of such systems include vehicles, a ball rolling on a surface, and many mobile robots.

- **Partially observable and non-markovian systems:** As the CON dynamics are evaluated based on the latent position and velocity encoded by the observation of the current time step and the observation-space velocity, we implicitly assume that the system is (a) fully observable, and (b) fulfills the Markov property. This assumption might not hold for partially observable systems, such as systems with hidden states or systems with delayed observations. Examples of such cases include settings where the system is partially occluded or in situations without sufficient (camera) perspectives covering the system. Furthermore, time-dependent material properties, such as viscoelasticity or hysteresis, that are present and significant in some soft robots and deformable objects are not captured by the method in its current formulation.
## Planned additions to the Appendix

Furthermore, we could envision adding a more detailed discussion (similar to the one below), or alternatively a table, to the Appendix:

### Systems for which we would expect the proposed method to work

- **Mechanical systems with continuous dynamics, dissipation, and a single, attractive equilibrium point:** The proposed method is a very good fit for mechanical systems with continuous dynamics, dissipation, and a single, attractive equilibrium point. In this case, the real system and the latent dynamics share both the energetic structure and stability guarantees. Examples of such systems include many soft robots, deformable objects with dominant elastic behavior, other mechanical structures with elasticity, etc.

- **Local modeling of (mechanical) systems that do not meet the global assumptions:** Even if the global assumptions of the proposed method are not met, the method can still be applied to model the local behavior around a local asymptotic equilibrium point of the system (i.e., in the case of multi-stability). For example, the method could be used to model the behavior of a robotic leg locally in contact with the ground, a cobot's interaction with its environment, etc.

### Systems for which we could envision the proposed method to work under (minor) modifications

- **Mechanical systems without dissipation:** The proposed method would currently not work well for mechanical systems without any dissipation, as (a) the original system will likely not have a globally asymptotically stable equilibrium point, and more importantly, (b) we currently force the damping learned in latent space to be positive definite. However, these systems are not common in practice as friction and other dissipation mechanisms are omnipresent, and the proposed method can learn very small damping values (e.g., the mass-spring+friction system).
A possible remedy could be to relax the positive definiteness of the damping matrix in the latent space, allowing for zero damping. This would allow the method to work for systems without dissipation, such as conservative systems. Examples of such systems include a mass-spring system without damping, the n-body problem, etc.

- **(Mechanical) systems with discontinuous dynamics:** The proposed method might underperform for systems with highly discontinuous dynamics, such as systems with impacts, friction, or other discontinuities. In these cases, the latent dynamics might not capture the real system's behavior accurately, and the control performance of feedforward + feedback will very likely be worse than pure feedback. Again, the method should be able to capture local behavior well. A possible remedy for learning global dynamics could be to augment the latent dynamics with additional terms that capture the discontinuities, such as contact and friction models (e.g., stick-slip friction).

- **(Mechanical) systems with multiple equilibrium points:** The original system having multiple equilibria conflicts with the stability assumptions underlying the proposed CON latent dynamics. In this case, as, for example, seen on the pendulum+friction and double pendulum + friction results, the method might work locally but will not be able to capture the global behavior of the system. A possible remedy could be to relax the global stability assumptions of the CON network. For example, the latent dynamics could be learned in the original coordinates of CON while allowing $W$ also to be negative definite. This would allow the system to have multiple equilibria & attractors. Examples of such systems include a robotic arm under gravity, pendula under gravity, etc.

- **(Mechanical) systems with periodic behavior:** The proposed method will likely not work well for systems with periodic behavior, as they do not have a single, attractive equilibrium point. Examples of such systems include a mass-spring system with a periodic external force, a pendulum with a periodic external force, some chemical reactions, etc. Again, it is likely possible to apply the presented method to learning a local behavior (i.e., not completing the full orbit). A possible remedy could be to augment the latent dynamics with additional terms that capture the periodic behavior, such as substituting the harmonic oscillators with Van der Pol oscillators to establish a limit cycle or a supercritical Hopf bifurcation.
> _RW9L3:_ I agree with the authors' view that their method might underperform for systems with highly discontinuous dynamics. Since this paper is framed as a method for prediction and control of systems, with soft robots being only one application area (rather than a method for the control of soft robots), I would appreciate a lucid analysis of classes of systems/application areas in which the method might perform strongly/weakly (including a clear reasoning why; and potentially a remedy left to future work). Please consider including this, both from a systems theory and an application standpoint into the paper.

We appreciate the Reviewer's feedback and suggestion. We fully agree with the Reviewer that specific categories or examples of systems for which the proposed method might be suitable or unsuitable should be discussed in the paper. In the following, we will detail the planned changes.

## Planned changes to the Introduction

First, we will add the following sentence to the end of the introduction (i.e., Section 1):

> The proposed methodology is particularly well-suited for learning the latent dynamics of mechanical systems with continuous dynamics, dissipation, and a single, attractive equilibrium point. Examples of such systems include many soft robots, deformable objects with dominant elastic behavior, other mechanical structures with elasticity, or locally for other mechanical systems such as robotic manipulators, legged robots etc. For these systems, we can fully leverage the structural prior of the proposed latent dynamics including the integrated stability guarantees. If the system is actuated, the learned dynamics can be subsequently exploited for model-based control, as demonstrated in Section 5.

## Planned changes to the _Limitations_ section

We will also extend the _Limitations_ section (i.e., Section 6):

> **Limitations.** While we think our proposed method shows great potential and opens interesting avenues for future research, there are currently certain limitations. For example, the proposed method of learning (latent) dynamics implicitly assumes that the underlying system adheres to the Markov property (e.g., the full state of the system is observable), that it can be approximated by a system with mechanical structure, and that it has an isolated, globally asymptotically stable equilibrium. This is, for example, the case for many mechanical systems (e.g., some continuum soft robots, deformable objects, and elastic structures) with continuous dynamics, convex elastic behaviour, dissipation and whose time-dependent effects (e.g., viscoelastcity, hysteresis) are negligible. Even if these conditions are not met globally, the method can be applied to model the local behavior around an asymptotic equilibrium point of the system (e.g., robotic manipulators, legged robots) with added stability benefits for out-of-distribution samples. Alternatively, the method could be extended to relax some of these assumptions, e.g., by allowing for multiple equilibria, zero damping, or by incorporating additional terms to capture discontinuous dynamics (e.g., stick-slip models) or period motions (e.g., limit cycles such as the Van der Pol oscillator). For some physical systems, such as nonholonomic systems, partially observable systems, or systems with non-Markovian properties, the proposed method might not be suitable. Examples of such systems include mobile robots and systems with hidden states or delayed observations. Finally, the application of this method to non-physical systems, such as financial systems, social networks, or other complex systems, is out of the scope of this work. _THE DISCUSSION OF OTHER LIMITATIONS OF THIS PAPER WHICH ARE ALREADY MENTIONED IN THE INITIAL PAPER SUBMISSION CONTINUES HERE..._
Thank you for the response.

**Application of CON to non-soft robots**

I agree with the authors' view that their method might underperform for systems with highly discontinuous dynamics. Since this paper is framed as a method for prediction and control of systems, with soft robots being only one application area (rather than a method for the control of soft robots), I would appreciate a lucid analysis of classes of systems/application areas in which the method might perform strongly/weakly (including a clear reasoning why; and potentially a remedy left to future work). Please consider including this, both from a systems theory and an application standpoint into the paper.

**Application of CON to non-physical systems**

Thank you for explicitly limiting the application area to physical systems.
Thank you for your response! I hope that you don't see my direct nature as rude because that was not at all my intention.

1. Yes, I agree
2. and 3. No need to rename it but I definitely like the perspective that you intentionally designed the network so that you can find a "simple" proof. That would be great to emphasize in the paper.

Especially thank you for the new example. I appreciate your work here and now understand better the details. I will raise my score by 3 points.
## Example 2: A system where the W-coordinate transformation does not help

The proof underlying Theorem 1 is only (relatively) simple as the coordinate transformation $x_\mathrm{w} = W x$ helps us to identify the potential energy function in the $\mathcal{W}$ coordinates, where the hyperbolic nonlinearity operates elementwise. However, this is not always possible. For example, consider the system $\ddot{x} + K x + D \dot{x} + W^{\mathrm{T}} \tanh(Wx + b) = 0$, where $K$, $D$, $W \succ 0$ are positive definite matrices. 

Unfortunately, the unactuated dynamics in the \mathcal{W} coordinates are now given by 

$M_w \ddot{x}_w + K_w x_w + D_w \dot{x}_w+ W^T \tanh(x_w + b) = 0$, 

for which we cannot easily derive a potential energy function (i.e., integrate), as the hyperbolic term $W^\mathrm{T} \tanh(x_\mathrm{w} + b)$ is not (easily) separable.

Remark: In contrast to Example 1, this system actually has a single equilibrium point. 

This example motivates why we chose the proposed CON network architecture such that we can (easily) derive the kinetic and potential energy terms, subsequently prove GAS & ISS, and perform model-based control.
> Finally, it's the same story for deriving stability from the single oscillators to the network. Of course you can make the stiffness matrix negative definite to achieve a non-stable network but again, if all parameters are positive (definite), meaning defined so that you would expect the system to be stable, I think that the system is stable. I'd change me opinion if you can provide a counter example where the parameters are positive (definite) but would lead to an unstable system / network.

We are now ready to go back to our technical discussion, with the discussed counterexamples. These, hopefully, also better illustrate our point #3. Essentially, our argument here is that even positive-definite stiffness & damping matrices are not sufficient to claim global asymptotic stability. This should answer the Reviewer's direct technical question.

In the following, we will give two examples of modifications to the CON network (in the original coordinates) that would have, according to the Reviewer's argument, appeared to be stable but for which global stability proof would be hard to prove.

## Example 1: A system with multiple equilibria and without a valid potential energy

We consider the slightly modified network dynamics $\ddot{x} + K x + D \dot{x} + W^{-1} \tanh(Wx + b) = 0$, where $K, D, W \succ 0$ are positive definite matrices.

The equilibria of this system are given by the characteristic equation $K \bar{x} + W^{-1} \tanh(W \bar{x} + b) = 0$.
For the system to be globally asymptotically stable, it would need to have a single equilibrium point. However, if we simulate the system with $K = [[5.0, -2.2], [-2.2, 1.0] \succ 0$ (positive-definite as symmetric and positive eigenvalues $0.027$ and $5.973$), $D = \mathrm{diag}(0.2, 0.2) \succ 0$, $W = [[1.0, -2.2], [-2.2, 5.0]] \succ 0$ (positive-definite as symmetric and positive eigenvalues $0.027$ and $5.973$), we notice that the system is actually bistable with two attractors at $\bar{x}_1 = [-212.6, -475.2]$ and $\bar{x}_2 = [212.6, 475.2]$. This is illustrated in the time series plot and phase portrait attached at https://anonymous.4open.science/r/neurips24-20062-rebuttal-7770. Therefore, the system is **not** globally asymptotically stable.

For $\tau_\mathrm{pot} = -K x -W^{-1} \tanh(Wx + b)$ to be a valid potential force, it would need to satisfy the property $\frac{\partial \tau_\mathrm{pot}}{\partial x} = \left ( \frac{\partial \tau_\mathrm{pot}}{\partial x} \right )^\mathrm{T}$. Therefore, we derive $\frac{\partial \tau_\mathrm{pot}}{\partial x} = -K - W^{-1} \mathrm{diag}(\mathrm{sech}^2(Wx + b)) W$. Its transpose is given by $\left ( \frac{\partial \tau_\mathrm{pot}}{\partial x} \right )^\mathrm{T} = -K - W^\mathrm{T} \mathrm{diag}(\mathrm{sech}^2(Wx + b)) W^{-\mathrm{T}}$. Therefore, $\tau_\mathrm{pot}$ only stems from a potential iff $W^{-\mathrm{T}} = W \rightarrow W^\mathrm{T} W = \mathbb{I}$ (i.e., $W$ is orthogonal). This is not the case for a general positive-definite $W$.

This example motivates that even if the matrices are positive-definite, the system can still have multiple equilibria (i.e., lose global asymptotic stability) and not have a valid potential energy function.
> Finally, it's the same story for deriving stability from the single oscillators to the network. Of course you can make the stiffness matrix negative definite to achieve a non-stable network but again, if all parameters are positive (definite), meaning defined so that you would expect the system to be stable, I think that the system is stable. I'd change me opinion if you can provide a counter example where the parameters are positive (definite) but would lead to an unstable system / network.

We appreciate the direct nature of the Reviewer's comment, which rises indeed a compelling point with which however we still disagree. In our next response, we will try again to convince the Reviewer of our argument, although proving that something is _not trivial_ is quite a challenging task, as triviality is very subjective.

So, before delving into it, we would like to take the chance of zooming out. Indeed, at the moment, we see the risk that we are debating on a relatively narrow point while it may be that we could essentially agree on several key aspects. Or at least, we believe we can work to find a common ground on these while we continue in parallel our discussion.

1. We believe that Theorem 2 is the main _stability proof_ contribution of the paper. To the best of our knowledge, we are the first in this community/problem setting to provide explicit input-to-state stability guarantees, including convergence rates into the region of attraction (Theorem 2). We believe that the proof of Theorem 2 is far from trivial.

Does the Reviewer disagree on this point?

2. In this sense, we agree with the Reviewer that the proof of Theorem 1 is comparably simpler (although not trivial!) and less impactful than the one of Theorem 2.

Would renaming "_Theorem 1_" as a "_Lemma_" or even as a "_Proposition_" better reflect the opinion of the Reviewer on the importance of this contribution?

3. Finally, we realize now that we did not make a great job, with the manuscript and with our previous answers, in conveying the message that: the relatively simplicity of the proof of Theorem 1 is no accident.
The proof of Theorem 1 could be derived via _standard_ arguments from the control of mechanical systems, _because_ we designed the CON network in a particular way. So the Theorem 1 statement with its relative simplicity implicitly stresses this point.

If the Reviewer agrees on this point, we are happy to include in the revision any remarks or changes that they see fit to better bring this point home.

Specifically, the proof in its current form is possible because (a) a coordinate transformation exists that allows us to identify a Lyapunov candidate, and (b) the system has a potential energy. We see these two design choices a contribution by themselves as they allow the proof to be simpler as it could have been otherwise.

We want to stress that already small modifications to the network formulation and the underlying assumptions in Theorem 1 would have made the proof much more difficult (or even impossible).
Thank your for your response! You are right that passivity does not imply *asymptotic* stability but most often stability. That said, your counterexamples do not convince me because of the following:

1. Example: It is pretty obvious that for a negative W the system's equilibrium can be unstable as you introduce a negative stiffness coefficient in the system. Thus, I'd still still say that the result seems trivial as the systems seems always stable (asymptotic stable) for positive damping and stiffness coefficients.

2. Example: Your analysis is wrong. Following Khalil, the storage function must be a "positive semidefinite function" that is obviously not the case in your example for a negative K. I don't want to look into the details, but I think that for the linear system as shown in example 2, passivity leads to (asymptotic) stability. Please provide a valid counter example if my believe is wrong.

Finally, it's the same story for deriving stability from the single oscillators to the network. Of course you can make the stiffness matrix negative definite to achieve a non-stable network but again, if all parameters are positive (definite), meaning defined so that you would expect the system to be stable, I think that the system is stable. I'd change me opinion if you can provide a counter example where the parameters are positive (definite) but would lead to an unstable system / network.