PaperlayerView exact revision on arXiv

Source-linked AI summary

Time Without Timesteps: Simulating Coupled Dynamical Systems via Self-Consistency

Liyu Zerihun, Mark Shinyoung Lee

arXiv:2609.03358v1cs.LG

TL;DR

Coupled dynamical simulation is usually advanced sequentially through time, while the circular dependence among subsystem trajectories makes the whole coupled solution a fixed-point problem. The paper replaces subsystem integrations with learned trajectory operators and solves for self-consistency. It reports fewer sequential solver iterations, implicit gradients, and a Jacobian-based prediction of convergence, alongside limits from surrogate error and training coverage.

  • Problem

    In coupled systems, each subsystem’s driving signal depends on neighboring trajectories, making the trajectory-level consistency problem circular.

    Source: §1 · p. 1

  • Method

    The method replaces subsystem integration with learned trajectory operators and enforces self-consistency across their outputs through a fixed-point solve.

    Source: §1 · p. 2

  • Results

    The learned operator’s Jacobian spectral radius predicts the coupled solve’s convergence boundary, while operator accuracy and composability are distinct and do not improve together.

    Source: §5 · p. 10, §7 · p. 11

  • Takeaways & Limitations

    The approach makes trajectory-level simulation and differentiation through its fixed point practical in the demonstrated systems, but remains an early exploration rather than a classical-solver replacement.

    Source: §5 · p. 10, §6 · p. 10

    The study covers only small systems—20 oscillators and neuron pairs—and reports 1–4% surrogate error, with fixed-point conditioning amplifying that error.

    Source: §6 · p. 10

Abstract

from arXiv · show

Numerical simulation of dynamical systems is usually organized as a causal march through time: each state is computed from the previous one. We explore a different formulation for coupled systems. For each subsystem type we train a neural surrogate mapping a full driving trajectory and initial condition directly to a full output trajectory; following classical waveform relaxation, coupled systems are assembled by enforcing self-consistency among these trajectories: simulation becomes a fixed-point problem over complete trajectories rather than a stepwise rollout. On coupled van der Pol oscillators and Hodgkin-Huxley neuron networks, sequential depth becomes the number of solver iterations: 4-10 Newton iterations where the reference integrator takes 1500 steps. The gradient likewise loses its time recursion: it becomes a linear system solved by GMRES at memory independent of solver depth. A single scalar measured from the learned operator, the spectral radius of its Jacobian, predicts in advance where the coupled solve will converge; past that boundary, unrolled backpropagation diverges and a Neumann adjoint fails, while the implicit gradient remains correct to 0.04%. We report where the approach succeeds and where surrogate error degrades it.

1 Introduction

The paper reframes coupled simulation as a global consistency problem over whole trajectories: learned subsystem operators replace repeated integrations, and their outputs are coupled through a fixed point. This changes sequential computation from time steps to solver iterations and makes gradients a linear-system solve rather than a time-recursive calculation.

1 Introduction

  • The coupling is circular because each subsystem’s driving signal depends on the unknown trajectories of its neighbors.Waveform relaxation resolves this by integrating every subsystem over the full window on each sweep, making each iteration as expensive as a complete simulation.Source: §1 · pp. 1–2
  • TWT replaces each subsystem’s inner integration with a learned operator that maps a driving trajectory and initial condition to an output trajectory, then enforces mutual consistency across subsystems.The resulting physical trajectory is selected by dynamics and coupling constraints, rather than by the order in which a solver constructs it.Source: §1 · pp. 1–2
  • 4–10 Newton iterations solve trajectories represented over 1500 timesteps, replacing timestep-count sequential depth with iteration-count depth.Each fixed-point iteration evaluates the module operators independently.Source: §1 · p. 2
  • At a fixed point, the gradient is a linear system without time recursion, solvable by Krylov methods whose matrix-vector products apply the module operators.The paper contrasts this with time-stepping gradients, which either store intermediate states or sweep backward across the horizon.Source: §1 · p. 2
  • Time remains in the trajectory representation, but the solver no longer constructs trajectories through a local update at every timestep.For oscillators, trajectories use a Chebyshev coefficient basis; the paper presents TWT as a different paradigm, not a replacement for classical integration.Source: §1 · pp. 2–3

2 Background and Related Work

TWT combines waveform relaxation’s subsystem decomposition with learned trajectory operators and a self-consistency solve. It differs from other global-trajectory and neural-surrogate approaches by composing learned operators across coupled subsystems, while treating physical trajectories—not hidden network states—as the fixed point.

2.1 Waveform Relaxation

  • Waveform relaxation iterates whole subsystem waveforms, but each sweep requires full integrations; TWT replaces those integrations with learned operators, enabling Newton–Krylov forward solves and GMRES adjoints.Unlike Parareal, which decomposes the time axis, TWT decomposes by subsystem.Source: §2.1 · p. 3

2.2 Trajectories as Global Objects

  • Variational integrators and PINNs also use global trajectory views, but TWT amortizes subsystem solves into operators and composes them at inference through self-consistency.Variational integrators still use a discretized temporal scheme; PINNs optimize a solution for each problem instance.Source: §2.2 · p. 3

2.3 Neural Surrogates and Neural Operators

  • Unlike isolated neural operators, TWT places trajectory-to-trajectory surrogates inside a coupled graph and enforces mutual consistency, making operator Jacobian norms relevant alongside isolated accuracy.The two properties do not improve together.Source: §2.3 · p. 3

2.4 Differentiable Simulation and Equilibrium Models

  • Differentiable simulators and neural ODEs differentiate through time-stepping, whereas TWT borrows the fixed-point view of DEQs for coupled physical trajectories that can be checked against external ground truth.Time-stepping gradients inherit the forward solve’s sequential structure.Source: §2.4 · p. 4

3 Method

The method models coupled subsystems as trajectory operators whose outputs must satisfy coupling constraints simultaneously. It solves that fixed point with Newton–Krylov methods, differentiates it implicitly, and separates formulation or discretization error from learned-surrogate error.

3.1 Problem Formulation

  • Each subsystem evolves under a drive generated by neighboring states, so solving components independently is impossible until their mutually dependent input trajectories are determined.The coupling graph represents subsystem dependencies.Source: §3.1 · p. 4

3.2 Module-Level Trajectory Surrogates

  • A module-level surrogate approximates an isolated subsystem’s solution operator, returning a full output trajectory from a discretized input trajectory in one evaluation.Oscillator trajectories use Chebyshev coefficients, while neuron trajectories use raw samples; the representation is an implementation choice, not a change to the formulation.Source: §3.2 · p. 4
  • Training uses a rolling buffer of generated drives, but its coverage can miss inputs induced by coupling, causing errors initially misattributed to the method.The surrogate is a residual MLP trained with plain MSE; the paper flags training-distribution coverage as a practical concern.Source: §3.2 · p. 5

3.3 Self-Consistency Solve

  • At the coupled fixed point, every trajectory matches its module’s prediction under neighbor-induced inputs; Picard relaxation is contraction-limited, so the method uses Newton–Krylov with GMRES without forming the Jacobian.Newton converges where I − JF is nonsingular, without requiring contraction.Source: §3.3 · p. 5

Implicit differentiation.

  • Implicit differentiation solves a linear system without dependence on the forward solver’s path or time recursion; a Neumann-series adjoint restores the contraction requirement that Newton avoids.The paper instead solves the adjoint system by GMRES.Source: Implicit differentiation · p. 5

3.4 Exact-Surrogate Interpretation

  • With exact subsystem operators, a zero-residual fixed point matches the original coupled initial-value problem under standard existence and uniqueness conditions; learned operators make this correspondence approximate.The paper reports self-consistency residual, trajectory error, and an exact-integrator floor separately to distinguish formulation and discretization error from learning error.Source: §3.4 · p. 6
  • Two valid O(∆t) drive conventions produce different fixed points, with their discrepancy growing from 0.8% to 10.4% as coupling strength increases.An exact-operator check isolated this discretization mismatch from surrogate error.Source: §3.4 · p. 6

4 Experiments

Across oscillator and neuron systems, the experiments test exact and learned trajectory-level fixed points, their convergence, and gradients. Results show accurate coupled trajectories in tested regimes, a measured convergence boundary, and implicit gradients that remain valid beyond it.

4 Experiments

  • The experiments use coupled van der Pol oscillators and Hodgkin–Huxley neurons; held-out median relative surrogate errors are 3.05% and 2.61%, respectively.The neuron surrogate predicts only the coupling-relevant variable s(t), leaving voltage and gating variables latent.Source: §4 · p. 6

4.1 Validating the Fixed Point

  • With an exact RK4 module operator, the fixed point reproduces the coupled reference exactly; Picard iteration counts grow as the contraction condition binds.At kc = 2, exact-integrator Newton takes 1548 s, versus 13.7 s for Picard, because each GMRES iteration requires full subsystem integrations.Source: §4.1 · p. 6, Table 1 · p. 7

4.2 Coupled Trajectories

  • For kc ≤1.5, coupled oscillator errors are 1.0%–18.2% with all five seeds converged, versus feedforward errors above 100%.The solve takes 4–10 Newton iterations on trajectories spanning 1500 timesteps.Source: §4.2 · p. 7, Table 2 · p. 7

4.3 Contraction Predicts Convergence

  • The learned operator’s spectral radius predicts the coupled solve’s convergence boundary; iteration count, convergence rate, and trajectory error follow that radius across stiffness.The predicted boundary is k*c = 1/ρ, with ρ measured by power iteration on the learned operator.Source: §4.3 · p. 7, Figure 1 · p. 8

4.4 Implicit Differentiation

  • Past the contraction boundary, the implicit gradient remains correct to 0.04% against finite differences, while unrolled backpropagation diverges by a factor of 83 and the Neumann adjoint returns non-finite values.In random systems, parameter error falls from 83.2% to a median of 1.2% (worst case 2.1%), and output error falls from 105% to 1.3%.Source: §4.4 · pp. 8–9

4.5 Hodgkin–Huxley Neurons

  • True post-pulse neuron activity falls from 0.0527 to 0.0397 as recurrence strengthens, while feedforward activity stays at 0.0551; the fixed point tracks the changing activity but fails at wrec = 5.The exact-operator floor is near 1.9%, so most remaining error is attributed to the surrogate.Source: §4.5 · p. 9, Table 5 · p. 9

5 Discussion

The experiments support trajectory-level self-consistency as a valid formulation independently of learning, while showing that learned operators make Newton practical and that operator composability is distinct from isolated accuracy. They also expose a limit: a small self-consistency residual does not guarantee low physical error.

5 Discussion

  • The fixed point is the correct object independently of learning, and the exact-operator floor separates formulation error from surrogate error.This validation supports the formulation while clarifying what the surrogate contributes.Source: §5 · p. 10
  • Newton removes Picard's contraction condition for learned operators, but with an exact integrator it is roughly 100× slower than Picard.Learned operators make Newton affordable because their Jacobian-vector products are cheap.Source: §5 · p. 10
  • Operator accuracy and composability are distinct, need not improve together, and composability predicts end-to-end behavior before a coupled solve.The pattern appears across three independently trained models and across stiffness within one model; power iteration on the Jacobian measures composability.Source: §5 · p. 10
  • A self-consistency residual of 10^-12 can coexist with substantial physical error, so residuals and exact-operator floors are needed to interpret results.This occurs at kc = 2; the residual measures agreement under the learned surrogate, not satisfaction of the original differential equation.Source: §5 · p. 10

6 Limitations and Future Work

The demonstration is limited by small system scale, surrogate error, and training-distribution coverage, and it does not establish wall-clock speedup. Its experiments cover only ODE systems with linear or synaptic coupling; several other system classes remain untested.

6 Limitations and Future Work

  • The results cover 20 coupled oscillators and neuron pairs, with surrogate error of 1–4%; the exact-operator floors are 0.0% for oscillators and 1.9% for neurons.The authors have not demonstrated the approach beyond these small systems, and fixed-point conditioning can amplify surrogate error.Source: §6 · p. 10
  • No wall-clock speedup is reported because neither implementation is optimized, making that comparison unmeaningful.The reported structural comparison is 4–10 iterations against 1500 timesteps; practical speed claims require an optimized implementation.Source: §6 · p. 10
  • Training-distribution coverage was the dominant practical failure mode: experiments sampled drive magnitudes or recurrent weights absent from training data.Coupling determines the inputs operators encounter, so the induced input distribution is not known in advance.Source: §6 · pp. 10–11
  • The experiments cover ODEs with linear or synaptic coupling; PDEs, stiff multiscale systems, chaotic dynamics, discontinuities, and event-driven systems need separate investigation.Windowing the time axis, learned preconditioners, compressed trajectory representations, and spectral-radius-focused training are proposed as future directions.Source: §6 · p. 11

7 Conclusion

Time Without Timesteps recasts coupled simulation as self-consistency among learned trajectory operators, changing the sequential structure of both the forward solve and its gradient. The learned operator’s Jacobian predicts the convergence boundary, while the method remains a scoped computational paradigm rather than a universal replacement for numerical integration.

7 Conclusion

  • The coupled solve converges in 4–10 iterations on trajectories spanning 1500 timesteps, and its gradient becomes a linear system solvable at memory independent of solver depth.The operator Jacobian predicts the convergence boundary; past it, only the implicit gradient remains valid.Source: §7 · p. 11

Related papers

Loading 2609.03358v1…