discretization_stability_omega_b

Base-Angular-Velocity Ripple Is a Forward-Euler Instability, Not a Control or Singularity Pathology

Statement

The persistent ripple observed in the base angular velocity ${\mathbf{\omega }}_b$ of the
coordinated free-flying controller is a period-2 (Nyquist) limit cycle produced by integrating
the reduced-impedance loop’s velocity-damping term with explicit forward Euler at a timestep
larger than the loop’s stiffest damping mode can resolve. It is a pure discretization artifact:
the underlying continuous-time loop is unconditionally stable, and the ripple is near-zero-mean so
it integrates out of the configuration (pointing ${\tilde{\mathbf{x}}}_b\sim 0.8^{\circ }$ and
coverage are unharmed).

Linearizing the reduced velocity error $\mathbf{e}=\breve{\mathbf{v}}-{\breve{\mathbf{v}}}_d$
about the desired velocity, the damping term dominates the unstable mode (stiffness contributes
$\omega dt=0.99<2$, comfortably stable):

$$\dot{\mathbf{e}}=-{\breve{\mathbf{M}}}^{-1}\breve{\mathbf{D}}\mathbf{e}.$$

Diagonalizing ${\breve{\mathbf{M}}}^{-1}\breve{\mathbf{D}}=\mathbf{V}\mathrm{diag}({\mu }_i){\mathbf{V}}^{-1}$
with decay rates ${\mu }_i>0$, the continuous solution $e_i(t)=e_i(0)e^{-{\mu }_{i}t}$ decays
monotonically. Forward Euler with step $dt$ replaces this with the scalar recurrence

$$e_i^{k+1}=(1-dt{\mu }_i)e_i^k\equiv {\rho }_i{e}_i^k,{\rho }_i=1-dt{\mu }_i,$$

so the discrete mode is unstable (a growing period-2 cycle, ${\rho }_i<-1$) exactly when

$$\boxed{dt{\mu }_{\max }>2}⟺dt>d{t}^{\star }\equiv \frac{2}{{\mu }_{\max }},{\mu }_{\max }=\max \mathrm{eig}\left({\breve{\mathbf{M}}}^{-1}\breve{\mathbf{D}}\right).$$

At a representative near-singular step (${\sigma }_G\sim 0.032$) the largest decay rate is
${\mu }_{\max }\approx 128{\mathrm{s}}^{-1}$, giving $d{t}^{\star }\approx 0.0156\mathrm{s}$. The nominal
$dt=0.03\mathrm{s}$ sits $\sim 2.5\times$ above the stability limit ($dt{\mu }_{\max }\approx 3.84$);
at $dt=0.01\mathrm{s}$ ($dt{\mu }_{\max }\approx 1.28<2$) the ripple vanishes.

Why ${\mu }_{\max }$ blows up near the singular configuration. Since
$\breve{\mathbf{M}}=({\mathbf{\Gamma }}^{-\top }\mathbf{M}{\mathbf{\Gamma }}^{-1}){|}_{[{\mathbf{\omega }}_b;{\mathbf{\nu }}_e^{\oplus }]}$,
at full arm extension ${\mathbf{\Gamma }}^{-1}$ has large entries, $\breve{\mathbf{M}}$ develops a
light mode (small eigenvalue), and ${\breve{\mathbf{M}}}^{-1}$ grows. With $\breve{\mathbf{D}}$
held at full strength, ${\mu }_{\max }=\max \mathrm{eig}({\breve{\mathbf{M}}}^{-1}\breve{\mathbf{D}})$
grows. The singularity therefore enters through the reduced inertia $\breve{\mathbf{M}}$, not
the damped pseudo-inverse — which is why the obvious lever (raising base damping) makes things
worse: it increases $\breve{\mathbf{D}}$, pushing $dt{\mu }_{\max }$ further past 2.

This exonerates two suspects. The $\gamma$-regularized velocity reconstruction is a passive,
sign-stable SPD map that cannot self-oscillate, and the singularity / conditioning layer
(singularity_threshold_cascade) is downstream of the artifact.
The driver is upstream, in the explicit integration of $\dot{\breve{\mathbf{v}}}={\breve{\mathbf{M}}}^{-1}\mathrm{RHS}$.

Assumptions

• Free-flying regime (fully-actuated 6-DOF base). The result is about the time-discretization of
  the reduced attitude+EE impedance loop, not about momentum conservation or the free-floating GJM
  machinery.
• The reduced state is $\breve{\mathbf{v}}=[{\mathbf{\omega }}_b;{\mathbf{\nu }}_e^{\oplus }]$ and the
  loop integrates $\dot{\breve{\mathbf{v}}}={\breve{\mathbf{M}}}^{-1}\mathrm{RHS}$ with explicit
  forward Euler at fixed $dt$, followed by a configuration step.
• Linearization about ${\breve{\mathbf{v}}}_d$ with the damping term dominating the stiffest mode;
  the stiffness contribution ($\omega dt=0.99<2$) is stable and dropped from the instability
  analysis.
• $\breve{\mathbf{M}}$ is SPD and $\breve{\mathbf{D}}$ is SPD, so ${\breve{\mathbf{M}}}^{-1}\breve{\mathbf{D}}$
  has real positive eigenvalues ${\mu }_i$ (real, diagonalizable spectrum — no complex damping modes).
• Numerical values (${\mu }_{\max }\approx 128$, $d{t}^{\star }\approx 0.0156$) are configuration-specific,
  measured at ${\sigma }_G\sim 0.032$ with $\breve{\mathbf{M}}$ exact and gains from the controller
  config (base.derate.gain = false, so the base block runs at full strength).
• The velocity clip ($v_{\max }=50$) does not bind — it is far above the $\sim 0.3$ ripple
  amplitude — so the limit cycle is bounded by the configuration-dependence of $\breve{\mathbf{M}}(\mathbf{q})$,
  not by saturation.

Proof sketch

Full derivation: omega_b_forward_euler_instability.md.

1. Reduce the loop. Write the per-step update $\breve{\mathbf{v}}\leftarrow \breve{\mathbf{v}}+dt{\breve{\mathbf{M}}}^{-1}\mathrm{RHS}$;
   linearize about ${\breve{\mathbf{v}}}_d$ in the error $\mathbf{e}=\breve{\mathbf{v}}-{\breve{\mathbf{v}}}_d$;
   keep the dominant damping term $\dot{\mathbf{e}}=-{\breve{\mathbf{M}}}^{-1}\breve{\mathbf{D}}\mathbf{e}$.
2. Diagonalize and discretize. Spectral decomposition ${\breve{\mathbf{M}}}^{-1}\breve{\mathbf{D}}=\mathbf{V}\mathrm{diag}({\mu }_i){\mathbf{V}}^{-1}$;
   per-mode forward Euler gives $e_i^{k+1}=(1-dt{\mu }_i)e_i^k$.
3. Read off the stability map. ${\rho }_i=1-dt{\mu }_i$: monotone decay ($0<dt{\mu }_i<1$), decaying
   ring ($1<dt{\mu }_i<2$), marginal period-2 ($dt{\mu }_i=2$), growing period-2 limit cycle
   ($dt{\mu }_i>2$). Hence the critical step $d{t}^{\star }=2/{\mu }_{\max }$.
4. Locate ${\mu }_{\max }$. Trace its growth to ${\breve{\mathbf{M}}}^{-1}$ blowing up near the
   singular configuration via ${\mathbf{\Gamma }}^{-1}$, with $\breve{\mathbf{D}}$ fixed.
5. Numerical confirmation — four independent ways (full table in the derivation):
   (i) the dt-independent lag-1 autocorrelation of the dominant data direction flips
   $-0.92\to +1.0$ between $dt=0.03$ and $dt=0.01$ (period-2 $\to$ smooth);
   (ii) the theory ratio $dt{\mu }_{\max }=3.84>2$ vs $1.28<2$, and the $18\times 18$ discrete-map
   dominant eigenvalue $-2.81$ (real, period-2) vs $+1.00$;
   (iii) mode alignment theory-vs-data $|\cos |=0.92$;
   (iv) a causal re-run of the mission at $dt=0.01$ (spec omega_ripple_dt010) removed the
   ripple exactly as predicted ($\sim 450\times$ amplitude drop, every direction’s autocorr
   flipping to $\sim +1.0$).

The principled fix (deferred — a separate decision with its own baseline re-pin) is to treat the
stiff damping implicitly (backward Euler / IMEX): $(\mathbf{I}+dt{\breve{\mathbf{M}}}^{-1}\breve{\mathbf{D}}){\mathbf{e}}^{k+1}={\mathbf{e}}^k$
gives ${\rho }_i=1/(1+dt{\mu }_i)\in (0,1)$ for all ${\mu }_i>0$ and all $dt$ — unconditionally stable,
monotone, no timestep bound, no gain re-tuning.

Source / provenance

• Ours — the project’s own numerical finding. Full derivation, code references, and the
  four-way confirmation table:
  omega_b_forward_euler_instability.md.
• Verification script: validation/omega_b_ripple_mechanism.py; causal re-run spec
  omega_ripple_dt010; provenance and verdicts in logs/logs_Jun21_26/CLAIMS.md (omega_b entries)
  and validation/INSIGHTS.md.
• Code: the forward-Euler velocity step and reduced-inertia $\breve{\mathbf{M}}$ live in the GNC
  controller and utils/robot.py (cited in the derivation).
• Notation: notation.md.

Caveats

• ${\mu }_{\max }\approx 128$ and $d{t}^{\star }\approx 0.0156\mathrm{s}$ are configuration-specific
  (measured at ${\sigma }_G\sim 0.032$); $d{t}^{\star }=2/{\mu }_{\max }$ shrinks as the system approaches
  singularity, so a single global timestep that is safe in open configurations can violate the bound
  near full extension.
• The analysis drops the stiffness term ($\omega dt=0.99$). This is justified at the nominal
  operating point but is not a proof for arbitrarily large stiffness gains; with much stiffer
  $\breve{\mathbf{K}}$ the stiffness mode could itself approach the $dt\omega =2$ bound.
• The limit-cycle amplitude is bounded by $\breve{\mathbf{M}}(\mathbf{q})$ nonlinearity, not by
  the velocity clip; if operating conditions push the ripple toward $v_{\max }$ the bounding mechanism —
  and the near-zero-mean property that protects pointing — would have to be re-checked.
• This note establishes the mechanism and the form of the fix; it does not adopt the
  implicit/IMEX integrator, shrink $dt$, or lower $\breve{\mathbf{D}}$. Any of those is a separate
  change requiring its own baseline re-pin.
• The result is about discretization of the coordinated impedance loop and is orthogonal to the
  Lyapunov stability of the continuous cascade
  (coordinated_control_lyapunov_stability): continuous
  stability is necessary but not sufficient for stability under explicit integration.

Related

• Sibling result: coordinated_control_lyapunov_stability —
  the continuous-time stability the discrete loop must preserve; this page shows where explicit
  integration breaks it.
• Sibling result: singularity_threshold_cascade — the conditioning
  layer this finding exonerates; the ripple is upstream of ${\sigma }_G$ scheduling, but ${\mu }_{\max }$
  (hence $d{t}^{\star }$) is driven by the same near-singular ${\mathbf{\Gamma }}^{-1}$ blow-up.
• Sibling result: circumcentroidal_decoupling — supplies
  $\mathbf{\Gamma }$, the hat/breve distinction, and the reduced inertia $\breve{\mathbf{M}}$.
• Topics: coordinated_control ·
  cascaded_systems.
• Notation: notation.md.