Horizontal Contact Neural ODE

Integrate this module at the latent-dynamics component of a neural ODE, continuous-time SSM, or world model. Let x\in\mathbb{R}^d be the latent state, p\in\mathbb{R}^d an auxiliary cotangent state, and u\in\mathbb{R} a scalar energy or confidence state. Choose m\ll d vector fields X_i(x); they can be fixed from a known manifold or produced by a small neural network, with X(x)=[X_1(x),...,X_m(x)]. Parameterize A_ heta=R_ heta R_ heta^ op+\epsilon I so it is positive definite. At every ODE step compute q=X(x)^ op p, h=0.5*q^T*solve(A,q)-rho*softplus(-u), then compute dx=X(x)*solve(A,q). Obtain dp and du by automatic differentiation of h: dp=-grad_x(h)-(partial_u h)*p and du=p^T*dx-h. Use an adaptive RK solver or fixed-step RK4, and feed x_T and u_T to the task head. The paper’s mathematical guarantees used directly are horizontal motion, strict convexity from A\succ0, superlinear quadratic control cost, and bounded negative scalar derivative; derivatives and numerical trajectories are computed by autodiff. Initialize p with an encoder output and u with a learned scalar or zero. First run a cheap ablation on latent trajectory prediction for a synthetic Heisenberg-group or sphere dataset with known horizontal vector fields, comparing an MLP-ODE, a standard Hamiltonian neural ODE, and this module at matched parameters and function evaluations. Measure one-step error, 50-step rollout error, gradient norm variance, and held-out likelihood. A successful signal is lower long-horizon error and fewer exploding trajectories at equal function evaluations; a second test should check whether learned horizontal rank m<d improves sample efficiency rather than merely adding parameters.