Calculus of Clouds: Interval ODEs as Multivector Flows

5.1 From Scalar Systems to Multivector Equations

Consider the coupled linear system:

\[ \dot{x} = ax + by, \quad \dot{y} = cx + dy \]

In classical notation, this requires two scalar interval equations evolved in parallel, with coupling terms linking them at every step. In GA, we collapse the state into a single multivector:

\[ \Psi = x\,\mathbf{e}_1 + y\,\mathbf{e}_2 \]

and encode the dynamics in a GA operator \(A\) such that the system becomes:

\[ \dot{\Psi} = A\Psi \]
Definition 5.1 — Multivector ODE

The operator \(A = \frac{1}{2}(a+d) + \frac{1}{2}(a-d)\,\mathbf{e}_1\mathbf{e}_1 + \frac{1}{2}(b+c)\,\mathbf{e}_1\mathbf{e}_2 + \frac{1}{2}(b-c)\,\mathbf{e}_1\mathbf{e}_2^*\) encodes the matrix \(\bigl(\begin{smallmatrix}a&b\\c&d\end{smallmatrix}\bigr)\) as a GA element. The action \(A\Psi\) is the geometric product, unifying scaling, rotation, and shear in a single operation.[1,2]

The formal solution is the multivector exponential:

\[ \Psi(t) = e^{At}\Psi_0 = \left(\sum_{k=0}^{\infty} \frac{(At)^k}{k!}\right)\Psi_0 \]

When \(A\) is a pure bivector \(A = \omega\mathbf{I}\), the exponential is a rotor \(e^{\omega t \mathbf{I}}\), and the solution is pure rotation — a stable orbit. When \(A\) has a non-zero scalar part \(\sigma\), the solution grows or decays exponentially.[1]

5.2 The Interval Blade as Solution Cloud

With uncertain initial conditions \(\Psi_0 \in [\Psi_a, \Psi_b]\) and uncertain operator \(A \in [A_1, A_2]\), the solution at time \(t\) is an interval blade:

\[ \Psi(t) \in e^{[A_1,A_2]t}[\Psi_a, \Psi_b] \]

This is the "cloud" of all possible solution trajectories. Its geometric properties directly reveal the system's behaviour:

Stable orbit (rotating blade)

If \(\text{Re}(\lambda) = 0\) for all eigenvalues \(\lambda\) of \(A\), the blade rotates without growing. The pseudoscalar magnitude \(|\langle\Psi\tilde\Psi\rangle_0|\) remains bounded.

Divergent (growing blade)

If \(\text{Re}(\lambda) > 0\) for any eigenvalue, the blade's magnitude increases. The outer product \(\Psi \wedge \dot\Psi\) has positive scalar part, indicating expansion of the solution cloud.

Contracting (stable attractor)

If \(\text{Re}(\lambda) < 0\) for all eigenvalues, the blade shrinks toward the origin. \(|\Psi| \to 0\) as \(t \to \infty\), regardless of initial conditions. The system is globally asymptotically stable.

5.3 Stability via the Geometric Product

For the operator \(A = \sigma + \omega\mathbf{I}\) (scalar growth rate plus rotation), the stability criterion reduces to a simple geometric test:

\[ \frac{d}{dt}|\Psi|^2 = 2\langle\Psi\,A\Psi\rangle_0 = 2\sigma|\Psi|^2 \]

Hence:

\[ \text{Stable} \iff \sigma < 0 \iff \langle A + \tilde{A} \rangle_0 < 0 \]
Theorem 5.2 — Geometric stability criterion

The system \(\dot\Psi = A\Psi\) is stable if and only if the scalar part of \(A + \tilde{A}\) is negative, where \(\tilde{A}\) denotes the reverse of \(A\). For a pure bivector operator, \(A + \tilde{A} = 0\) and the system is neutrally stable. This characterisation holds for all grades of \(A\).[1, §6.3]

Interactive — phase plane: interval blade evolution
0
40
15
Fig. 5.1. The interval blade \(\Psi(t)\) in the phase plane. The ellipse represents the set of all possible solution states consistent with the initial uncertainty. Adjusting \(\sigma\) (growth) scales the ellipse radially; \(\omega\) (rotation) drives the orbital motion; \(\Delta\) sets the initial angular spread of the interval.

5.4 Propagating Uncertainty: The Wrapping Problem

A fundamental difficulty in interval ODE integration is the wrapping effect: axis-aligned interval boxes, when rotated, must be over-approximated by a new axis-aligned box, causing accumulated inflation of the enclosure.[3]

GA interval blades avoid this problem naturally: because a blade is an oriented subspace object, rotation by a rotor maps a blade interval to another blade interval of the same shape. There is no wrapping loss for rotations in GA. The wrapping effect is replaced by the precise geometric tracking of the oriented cloud.

// Classical interval box rotation (wrapping loss) R([x₁,x₂] × [y₁,y₂]) ⊆ [x₁',x₂'] × [y₁',y₂'] // over-approximation // GA interval blade rotation (exact for blades) R · [Ψ₁, Ψ₂] · R̃ = [RΨ₁R̃, RΨ₂R̃] // exact, no wrapping // Stability test via scalar part of reverse ⟨A + Ã⟩₀ < 0 ⟹ system is stable ⟨A + Ã⟩₀ = 0 ⟹ neutral orbit (pure rotation) ⟨A + Ã⟩₀ > 0 ⟹ divergent

5.5 Integration Methods for Interval Multivector ODEs

The geometric structure of GA enables integration methods that respect the blade structure throughout, avoiding the coordinate-based artefacts of classical interval Runge–Kutta schemes.

  • 1

    Rotor exponential method. When \(A\) is a bivector (antisymmetric operator), exact integration is possible via the rotor exponential \(e^{At}\). The interval is propagated by applying the interval rotor \([e^{A_1 t}, e^{A_2 t}]\).[4]

  • 2

    Geometric midpoint method. An analogue of the classical midpoint Runge–Kutta scheme formulated in GA. The midpoint step is a rotor sandwich, preserving the blade norm to second order.

  • 3

    Validated Lohner method in GA. An adaptation of Lohner's QR method[3] for interval ODEs, replacing QR decomposition with GA rotor factorisation. This eliminates wrapping for any rotation-dominated system.

CLASSICAL INTERVAL BOX t=0 R over-approximation Wrapping: box inflates at each rotation step GA INTERVAL BLADE t=0 RΨR̃ t=1 exact No wrapping: blade rotates without inflation
Fig. 5.2. Comparison of wrapping behaviour. Left: a classical interval box must be re-enclosed in an axis-aligned box after each rotation, causing progressive inflation. Right: a GA interval blade maps exactly under rotor action — the rotated blade has the same area as the original. After Nedialkov [3].

5.6 Example: Uncertain Harmonic Oscillator

Consider the harmonic oscillator \(\ddot{x} + \omega^2 x = 0\) with uncertain frequency \(\omega \in [\omega_1, \omega_2]\). Writing \(\Psi = x\,\mathbf{e}_1 + \dot{x}\,\mathbf{e}_2\), the system becomes:

\[ \dot\Psi = \underbrace{(\omega\,\mathbf{e}_1\mathbf{e}_2)}_{A}\Psi \]

Since \(A = \omega\mathbf{I}\) is a pure bivector, the solution is exactly:

\[ \Psi(t) = e^{\omega\mathbf{I}t}\Psi_0 \tilde{e}^{\omega\mathbf{I}t} \]

The interval blade evolves as:

\[ \Psi(t) \in [e^{\omega_1\mathbf{I}t}, e^{\omega_2\mathbf{I}t}]\Psi_0 \]

This is precisely the annular sector from Chapter II — the interval rotor applied to the initial state vector. The angular width of the sector grows as \((\omega_2 - \omega_1)t\), providing a precise geometric picture of how frequency uncertainty translates to phase uncertainty over time.

// Harmonic oscillator with ω ∈ [ω₁, ω₂] A = ω · I // pure bivector operator Ψ(t) = e^(ωIt) Ψ₀ // exact solution (no truncation) // Interval solution Ψ(t) ∈ [e^(ω₁·I·t), e^(ω₂·I·t)] · Ψ₀ // Angular spread at time t Δθ(t) = (ω₂ - ω₁) · t // grows linearly, no wrapping loss // Stability check ⟨A + Ã⟩₀ = 0 // neutral — consistent with conservative oscillator
!

Connection across chapters. The annular sector picture of Chapter II, the null blade tangency condition of Chapter III, and the phase-plane evolution of this chapter are all the same geometric object — the interval blade — viewed at different stages of the computation pipeline. Geometric algebra provides the unified language in which these apparently distinct phenomena share a common description.

5.7 Summary and Outlook

This guide has developed a single thesis across five chapters: that spatial uncertainty is fundamentally a geometric object, and that geometric algebra provides the natural language for its rigorous treatment.

The blade generalises the scalar interval to arbitrary dimension and orientation. The geometric product propagates uncertainty through rotations exactly, without the wrapping loss that afflicts coordinate-based methods. The Meet and Join implement intersection and span directly, enabling efficient CSG ray tracing. And the multivector ODE collapses systems of uncertain equations into a single geometric flow, whose stability is characterised by a simple test on the operator's grade structure.

Open directions include the development of validated GA integration methods for non-linear systems, the application of conformal GA interval blades to robot kinematics under joint uncertainty, and the integration of affine arithmetic[5] with the GA framework to reduce wrapping in non-rotational dynamics.

References

  1. Hestenes, D. (1999). New Foundations for Classical Mechanics (2nd ed.). Kluwer Academic. §6.3.
  2. Lasenby, A., Doran, C., & Gull, S. (1993). A multivector derivative approach to Lagrangian field theory. Foundations of Physics 23(10), 1295–1327.
  3. Nedialkov, N. S., Jackson, K. R., & Corliss, G. F. (1999). Validated solutions of initial value problems for ordinary differential equations. Applied Mathematics and Computation 105(1), 21–68.
  4. Doran, C., & Lasenby, A. (2003). Geometric Algebra for Physicists. Cambridge University Press. Chapter 6.
  5. de Figueiredo, L. H., & Stolfi, J. (2004). Affine arithmetic: concepts and applications. Numerical Algorithms 37(1–4), 147–158.
  6. Macdonald, A. (2017). A survey of geometric algebra and geometric calculus. Advances in Applied Clifford Algebras 27, 853–891.