Doctoral Research · Space Robotics Inspection with a Free-Flying Space Manipulator
A Doctoral Research Journal Aerospace Engineering

The omega_b “oscillation” is a forward-Euler instability of the velocity-damping term

Verdict. The persistent base angular-rate ripple flagged in the Jun-21 check-in is not a control or singularity pathology. It 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 too large for the loop’s stiffest damping mode. The continuous-time system is perfectly stable; the ripple is a pure discretization artifact. Pointing (z_b ~ 0.8 deg) and coverage are unharmed because the ripple is near-zero-mean and integrates out of the configuration. Confirmed four independent ways (below).

This corrects the earlier “damped-pseudo-inverse velocity reconstruction” framing: the gamma-regularised reconstruction is a passive, sign-stable map here and is exonerated. The driver is upstream, in the explicit integration of v_breve_dot = M_breve^{-1} RHS.

1. The loop

Per step the controller integrates the reduced (“breve”) velocity v̆ = [omega_b; nu_e] (9-dim) and then moves the configuration (breve_controller.py:524, com_controller.py:400-401):

v̆_dot = M̆^{-1} · RHS,    RHS = −C̆(v̆−v̆_des) − D̆(v̆−v̆_des) − Jxᵀ K̆ x̃ + coupling + posture + FF
v̆     ← v̆ + dt · v̆_dot                 # forward Euler, dt = 0.03 s
q      ← q ⊕ v̆-reconstructed · dt        # pin.integrate

M̆ is the reduced inertia (Γ⁻ᵀ M Γ⁻¹) restricted to the [omega_b; nu_e] block (robot.py:370); D̆, K̆ are the block gains from cfg.controller.gains (base.derate.gain = false, so the base block is full strength K=D=2).

The data (per-direction SVD analysis of the logged state.v, with Γ rebuilt from the logged q) showed the ripple is a clean period-2 in a well-conditioned direction of Γ (singular value ~0.77, not the smallest), with the near-null direction carrying almost none of it. A well-conditioned direction is not nulled by Γ, so the ripple is real task-space velocity motion, and a static SPD solve cannot self-oscillate — both point at the explicit integrator, not the reconstruction.

2. The derivation (period-2 condition)

Linearise about the desired velocity; let e = v̆ − v̆_des. The damping term dominates the unstable mode (the stiffness contributes omega·dt = 0.99 < 2, comfortably stable, so drop it):

\[\dot e = -\,\bar M^{-1}\bar D\, e .\]

Diagonalise M̆⁻¹D̆ = V diag(μ_i) V⁻¹ with μ_i > 0 (positive damping/inertia decay rates). In mode i the continuous solution e_i(t) = e_i(0)\,e^{-μ_i t} decays monotonically — unconditionally stable, time-constant τ_i = 1/μ_i.

Apply forward Euler, e^{k+1} = e^k + dt\,\dot e^k:

\[e_i^{k+1} = (1 - dt\,μ_i)\,e_i^k \equiv ρ_i\, e_i^k .\]

The amplification ρ_i = 1 - dt\,μ_i gives the full picture:

regime condition multiplier ρ_i behaviour
monotone decay 0 < dt μ_i < 1 0 < ρ < 1 stable, no overshoot
ringing 1 < dt μ_i < 2 −1 < ρ < 0 stable, decaying 2-cycle
marginal dt μ_i = 2 ρ = −1 undamped period-2
unstable dt μ_i > 2 ρ < −1 growing period-2 → limit cycle

So the period-2 instability condition is dt · μ_max > 2, i.e. τ_min < dt/2: the requested damping decays faster than two timesteps can resolve, so forward Euler overshoots equilibrium and reverses sign every step. The growth is bounded into a sustained limit cycle by the configuration-dependence of M̆(q) (the velocity clip v_max = 50 is far above the ~0.3 ripple and does not bind).

Critical timestep: dt* = 2 / μ_max. Stable iff dt < dt*.

Why μ_max is large near the singular config

M̆ = (Γ⁻ᵀ M Γ⁻¹)|_{[omega_b;nu_e]}. At full arm extension Γ⁻¹ has large entries, so M̆ develops a light mode (a small eigenvalue), making M̆⁻¹ large; with D̆ fixed at full strength, μ_max = max eig(M̆⁻¹D̆) is large. The singularity therefore enters through the reduced inertia, not the damped pseudo-inverse — which is why raising base damping (the obvious lever) made things worse, not better: it increases D̆, pushing dt μ further above 2.

3. Numerical confirmation (four independent ways)

Measured at a representative near-singular step (s_min_G ~ 0.032), M̆ exact from robot.py:370, gains from cfg.controller.gains; full script validation/omega_b_ripple_mechanism.py:

evidence dt = 0.03 (nominal) dt = 0.01 (causal test)
1. data — dominant-direction lag-1 autocorr −0.92 (period-2) +1.0 (smooth)
2. theory — dt · μ_max 3.84 > 2 1.28 < 2
2. theory — 18×18 discrete-map dominant eigenvalue −2.81 (real, period-2) +1.00
3. mode alignment, theory vs data direction |cos| = 0.92 —
4. per-step velocity jump ‖Δv‖ 0.409 0.0009
provenance — rebuilt σ_min vs logged s_min_G 3.8e-16 1.9e-16

dt* = 2/μ_max = 2/128 ≈ 0.0156 s: the nominal dt = 0.03 is ~2.5× over the stability limit; dt = 0.01 is below it. The causal test (re-running the mission at dt = 0.01, spec omega_ripple_dt010, on the cluster) removed the ripple exactly as predicted — every direction’s lag-1 autocorr flipped to ~+1.0. The lag-1 autocorr flip is the dt-independent proof (a smooth signal reads +1 at any timestep); the 450× amplitude drop is partly dt-scaled and only corroborates.

4. The principled fix (deferred)

Treat the stiff damping implicitly (backward Euler / IMEX — implicit on D̆, explicit on the rest):

\[(I + dt\,\bar M^{-1}\bar D)\,e^{k+1} = e^k \;\Rightarrow\; ρ_i = \frac{1}{1 + dt\,μ_i} \in (0,1)\ \forall\,μ_i>0,\ \forall\,dt .\]

Unconditionally stable, monotone, never oscillates — no timestep bound, no gain re-tuning. Alternatives: shrink dt below dt* (3× the steps, 3× the cost), or lower D̆ (fights tracking). Adopting any of these is a separate decision with its own baseline re-pin and is not done here; this note establishes the mechanism and the fix’s form, not its adoption.

References