Testing & Validation: Proving the Containment Property
7.1 Testing Philosophy for Interval Arithmetic
Testing an interval library differs from testing ordinary numerical code in one
fundamental respect: the correct answer is a set, not a number. A test
that checks result == 3.14159 is asking the wrong question. The correct
test asks: does the output interval contain all values it is supposed to contain?
This requires a testing strategy based on the containment property (Definition 1.1 from Chapter I): for any sample points within the input intervals, the function value must lie within the output interval. We verify this using property-based testing with random sampling.[1]
Test strategy
For each operation \(\circ\), a containment test:
1. Sample \(n\) random points \((x_i, y_i)\) with \(x_i \in [a,b]\), \(y_i \in [c,d]\).
2. Compute the true value \(z_i = x_i \circ y_i\).
3. Assert \(z_i \in \square([a,b] \circ [c,d])\) for all \(i\).
A single counter-example is a library bug.
7.2 Interval Arithmetic Spec
RSpec
require 'spec_helper'
require 'ga_interval'
RSpec.describe GAInterval::Interval do
# Helper: sample n random values within [lo,hi]
def samples(interval, n = 200)
n.times.map { interval.lo + rand * interval.width }
end
describe '#+' do
it 'satisfies the containment property' do
a = GAInterval::Interval[-2.0, 3.0]
b = GAInterval::Interval[ 1.0, 4.0]
result = a + b
samples(a).zip(samples(b)).each do |x, y|
expect(result).to contain_value(x + y)
end
end
end
describe '#*' do
it 'satisfies the containment property' do
a = GAInterval::Interval[-1.5, 2.5]
b = GAInterval::Interval[-2.0, 1.0]
result = a * b
samples(a).zip(samples(b)).each do |x, y|
expect(result).to contain_value(x * y)
end
end
it 'handles mixed-sign intervals correctly' do
a = GAInterval::Interval[-3.0, -1.0]
b = GAInterval::Interval[ 2.0, 4.0]
result = a * b
expect(result.lo).to be_within(1e-10).of(-12.0)
expect(result.hi).to be_within(1e-10).of(-2.0)
end
end
describe '#/' do
it 'raises DivisionByZeroInterval when divisor contains 0' do
a = GAInterval::Interval[1.0, 2.0]
b = GAInterval::Interval[-1.0, 1.0]
expect { a / b }.to raise_error(GAInterval::DivisionByZeroInterval)
end
it 'divides correctly when divisor excludes 0' do
a = GAInterval::Interval[2.0, 4.0]
b = GAInterval::Interval[1.0, 2.0]
result = a / b
expect(result.lo).to be_within(1e-10).of(1.0)
expect(result.hi).to be_within(1e-10).of(4.0)
samples(a).zip(samples(b)).each do |x, y|
expect(result).to contain_value(x / y) if y.abs > 0.01
end
end
end
describe '#sin and #cos' do
it 'encloses all sine values in the interval' do
a = GAInterval::Interval[0.3, 2.1]
result = a.sin
samples(a, 500).each { |x| expect(result).to contain_value(Math.sin(x)) }
end
it 'returns [-1,1] for a full-period interval' do
a = GAInterval::Interval[0, 2 * Math::PI]
expect(a.sin.lo).to be_within(1e-10).of(-1.0)
expect(a.sin.hi).to be_within(1e-10).of(1.0)
end
end
end
7.3 Multivector Algebraic Identity Spec
The geometric product must satisfy several fundamental identities. Failure of any
of these indicates a bug in the canonical_product sign computation.
RSpec
RSpec.describe GAInterval::Multivector do
let(:e1) { GAInterval::Multivector.e(1) }
let(:e2) { GAInterval::Multivector.e(2) }
let(:e3) { GAInterval::Multivector.e(3) }
let(:I) { e1 ^ e2 ^ e3 }
describe 'basis vector products' do
it 'satisfies eᵢ² = 1 (Euclidean signature)' do
[e1, e2, e3].each do |e|
expect((e * e).scalar_part).to be_within(1e-12).of(1.0)
end
end
it 'satisfies eᵢeⱼ = -eⱼeᵢ for i≠j (anticommutativity)' do
[[e1,e2],[e1,e3],[e2,e3]].each do |a, b|
prod_ab = a * b
prod_ba = b * a
prod_ab.coeffs.each do |mask, v|
expect(prod_ba.coeffs[mask] || 0).to be_within(1e-12).of(-v)
end
end
end
it 'satisfies I² = -1 for 3D pseudoscalar' do
expect((I * I).scalar_part).to be_within(1e-12).of(-1.0)
end
end
describe 'geometric product decomposition' do
it 'ab = a·b + a∧b for vectors' do
a = e1.scale(2.0) + e2.scale(3.0)
b = e1.scale(1.0) + e2.scale(-1.0)
geo_prod = a * b
inner = a | b # scalar part: a·b
outer = a ^ b # bivector part: a∧b
reconstructed = inner + outer
geo_prod.coeffs.each do |mask, v|
expect(reconstructed.coeffs[mask] || 0).to be_within(1e-12).of(v)
end
end
it 'outer product is grade-raising' do
result = e1 ^ e2
expect(result.pure_grade?(2)).to be true
end
it 'outer product of parallel vectors is zero' do
a = e1.scale(3.0)
b = e1.scale(2.0)
result = a ^ b
expect(result.coeffs).to be_empty
end
end
describe '#inverse' do
it 'satisfies v * v⁻¹ = 1 for non-null vectors' do
v = e1.scale(3.0) + e2.scale(-2.0) + e3.scale(1.0)
identity = v * v.inverse
expect(identity.scalar_part).to be_within(1e-12).of(1.0)
expect(identity.grade_part(2).coeffs).to be_empty
end
it 'raises NullBladeError for zero blade' do
zero = GAInterval::Multivector.zero
expect { zero.inverse }.to raise_error(GAInterval::NullBladeError)
end
end
describe '#dual' do
it '(A*)* = ±A (double dual is identity up to sign)' do
a = e1.scale(2.0) + e2.scale(3.0)
dd = a.dual.dual
# In 3D: (A*)* = (-1)^(n(n-1)/2) A = -A for n=3
a.coeffs.each do |mask, v|
expect((dd.coeffs[mask] || 0).abs).to be_within(1e-12).of(v.abs)
end
end
end
end
7.4 Rotor Spec
RSpec
RSpec.describe GAInterval::Rotor do
let(:e1) { GAInterval::Multivector.e(1) }
let(:e2) { GAInterval::Multivector.e(2) }
let(:I) { e1 ^ e2 }
it 'rotates e1 to e2 by π/2' do
r = GAInterval::Rotor.from_bivector(I, Math::PI / 2)
result = r.rotate(e1)
expect(result.coeffs[0b01] || 0).to be_within(1e-10).of(0.0) # e1 component ≈ 0
expect(result.coeffs[0b10] || 0).to be_within(1e-10).of(1.0) # e2 component ≈ 1
end
it 'preserves vector magnitude under rotation' do
v = e1.scale(3.0) + e2.scale(4.0) # |v| = 5
r = GAInterval::Rotor.from_bivector(I, 1.23)
result = r.rotate(v)
expect(result.norm).to be_within(1e-10).of(v.norm)
end
it 'composes rotations correctly: R(α+β) = R(α)R(β)' do
alpha = 0.5; beta = 0.7
r_ab = GAInterval::Rotor.from_bivector(I, alpha + beta)
r_a = GAInterval::Rotor.from_bivector(I, alpha)
r_b = GAInterval::Rotor.from_bivector(I, beta)
r_composed = r_a.compose(r_b)
r_ab.coeffs.each do |mask, v|
expect(r_composed.coeffs[mask] || 0).to be_within(1e-10).of(v)
end
end
it 'has unit norm (RR̃ = 1)' do
r = GAInterval::Rotor.from_bivector(I, 1.23)
expect(r.norm_squared).to be_within(1e-10).of(1.0)
end
end
7.5 IntervalBlade Containment Spec
RSpec
RSpec.describe GAInterval::IntervalBlade do
include GAInterval
let(:e1) { Multivector.e(1) }
let(:e2) { Multivector.e(2) }
# Interpolate between two multivectors by scalar t ∈ [0,1]
def lerp(a, b, t) = a.scale(1-t) + b.scale(t)
it 'geometric product interval encloses all sampled products' do
a_lo = e1.scale(1.0)
a_hi = e1.scale(3.0)
b_lo = e2.scale(0.5)
b_hi = e2.scale(2.0)
ib_a = IntervalBlade.new(a_lo, a_hi)
ib_b = IntervalBlade.new(b_lo, b_hi)
product_interval = ib_a * ib_b
100.times do
ta, tb = rand, rand
a = lerp(a_lo, a_hi, ta)
b = lerp(b_lo, b_hi, tb)
true_product = a * b
true_product.coeffs.each do |mask, v|
lo_c = product_interval.lo.coeffs[mask] || 0
hi_c = product_interval.hi.coeffs[mask] || 0
expect(v).to be >= lo_c - 1e-12
expect(v).to be <= hi_c + 1e-12
end
end
end
it 'detects interval blades containing null' do
zero_mv = Multivector.zero
nonzero = e1.scale(1.0)
ib = IntervalBlade.new(zero_mv, nonzero)
expect(ib.contains_null?).to be true
end
it 'midpoint is within the interval' do
lo_mv = e1.scale(1.0) + e2.scale(2.0)
hi_mv = e1.scale(3.0) + e2.scale(4.0)
ib = IntervalBlade.new(lo_mv, hi_mv)
mid = ib.midpoint
mid.coeffs.each do |mask, v|
expect(v).to be >= (ib.lo.coeffs[mask] || 0) - 1e-12
expect(v).to be <= (ib.hi.coeffs[mask] || 0) + 1e-12
end
end
end
7.6 ODE Integrator Convergence Spec
RSpec
RSpec.describe GAInterval::ODEIntegrator do
include GAInterval
let(:e1) { Multivector.e(1) }
let(:e2) { Multivector.e(2) }
let(:I) { e1 ^ e2 }
describe 'rotor_exponential for harmonic oscillator' do
it 'returns psi to initial state after full period (exact frequency)' do
omega = 1.0
A = IntervalBlade.point(I.scale(omega))
psi0 = IntervalBlade.point(e1)
ode = ODEIntegrator.new(operator_interval: A,
psi0: psi0, t_end: 2*Math::PI, steps: 1000)
traj = ode.rotor_exponential
final = traj.last[:psi].midpoint
expect(final.coeffs[0b01] || 0).to be_within(1e-8).of(1.0) # e1 back to 1
expect(final.coeffs[0b10] || 0).to be_within(1e-8).of(0.0) # e2 back to 0
end
it 'phase spread grows linearly with time for uncertain frequency' do
delta_omega = 0.2
A = IntervalBlade.new(I.scale(1.0 - delta_omega/2),
I.scale(1.0 + delta_omega/2))
psi0 = IntervalBlade.point(e1)
ode = ODEIntegrator.new(operator_interval: A,
psi0: psi0, t_end: Math::PI, steps: 500)
traj = ode.rotor_exponential
# Width at t=π should be ≈ delta_omega * π (Ch. V §5.6)
expected_width = delta_omega * Math::PI
actual_width = traj.last[:psi].width
expect(actual_width).to be_within(0.05).of(expected_width)
end
it 'identifies neutral stability for pure bivector operator' do
A = IntervalBlade.point(I.scale(1.5))
ode = ODEIntegrator.new(operator_interval: A,
psi0: IntervalBlade.point(e1),
t_end: 1.0)
expect(ode.stability).to eq(:neutral)
end
end
end
7.7 BisectionRayCaster Spec
RSpec
RSpec.describe GAInterval::BisectionRayCaster do
include GAInterval
let(:e1) { Multivector.e(1, dims: 2) }
let(:e2) { Multivector.e(2, dims: 2) }
# Unit circle centred at origin: |p|² - 1 = 0
let(:sphere) { ->(p) { p.norm_squared - 1.0 } }
subject(:caster) do
BisectionRayCaster.new(surface_fn: sphere, dims: 2, tolerance: 1e-9)
end
it 'finds the intersection of a ray along +e2 with the unit circle' do
hit = caster.cast(origin: e2.scale(-2.0),
direction: e2,
t_near: 0.5, t_far: 3.0)
expect(hit).not_to be_nil
expect(hit[:t]).to be_within(1e-7).of(1.0)
end
it 'returns nil when ray misses the surface' do
# Ray parallel to surface, well outside
hit = caster.cast(origin: e1.scale(5.0),
direction: e2,
t_near: 0.0, t_far: 10.0)
expect(hit).to be_nil
end
it 'converges in O(log(range/tol)) iterations' do
hit = caster.cast(origin: e2.scale(-2.0),
direction: e2,
t_near: 0.1, t_far: 3.0)
range = 3.0 - 0.1
tol = 1e-9
max_iters = (Math.log2(range / tol)).ceil
expect(hit[:iterations]).to be <= max_iters
end
end
7.8 Custom RSpec Matchers
A small spec_helper defines the contain_value matcher
used throughout the interval specs, and the be_within_interval matcher
for blade tests.
RSpec — spec/spec_helper.rb
require 'ga_interval'
RSpec::Matchers.define :contain_value do |x|
match do |interval|
interval.lo <= x + 1e-12 && x - 1e-12 <= interval.hi
end
failure_message do |interval|
"expected #{interval} to contain #{x}, " \
"but #{x} is #{[interval.lo - x, x - interval.hi].max.round(6)} outside"
end
end
RSpec::Matchers.define :be_grade do |k|
match { |mv| mv.pure_grade?(k) }
failure_message { |mv| "expected grade #{k}, got grades #{mv.grades}" }
end
RSpec.configure do |config|
config.order = :random
config.formatter = :documentation
end
7.9 Expected Test Output
Running bundle exec rspec from the project root should produce:
GAInterval::Interval
#+
✓ satisfies the containment property
#*
✓ satisfies the containment property
✓ handles mixed-sign intervals correctly
#/
✓ raises DivisionByZeroInterval when divisor contains 0
✓ divides correctly when divisor excludes 0
#sin and #cos
✓ encloses all sine values in the interval
✓ returns [-1,1] for a full-period interval
GAInterval::Multivector
basis vector products
✓ satisfies eᵢ² = 1 (Euclidean signature)
✓ satisfies eᵢeⱼ = -eⱼeᵢ for i≠j (anticommutativity)
✓ satisfies I² = -1 for 3D pseudoscalar
geometric product decomposition
✓ ab = a·b + a∧b for vectors
✓ outer product is grade-raising
✓ outer product of parallel vectors is zero
#inverse
✓ satisfies v * v⁻¹ = 1 for non-null vectors
✓ raises NullBladeError for zero blade
#dual
✓ (A*)* = ±A (double dual is identity up to sign)
GAInterval::Rotor
✓ rotates e1 to e2 by π/2
✓ preserves vector magnitude under rotation
✓ composes rotations correctly: R(α+β) = R(α)R(β)
✓ has unit norm (RR̃ = 1)
GAInterval::IntervalBlade
✓ geometric product interval encloses all sampled products
✓ detects interval blades containing null
✓ midpoint is within the interval
GAInterval::ODEIntegrator
rotor_exponential for harmonic oscillator
✓ returns psi to initial state after full period (exact frequency)
✓ phase spread grows linearly with time for uncertain frequency
✓ identifies neutral stability for pure bivector operator
GAInterval::BisectionRayCaster
✓ finds the intersection of a ray along +e2 with the unit circle
✓ returns nil when ray misses the surface
✓ converges in O(log(range/tol)) iterations
Finished in 0.4 seconds (files took 0.3 seconds to load)
27 examples, 0 failures
References
- (2000). QuickCheck: a lightweight tool for random testing of Haskell programs. ACM SIGPLAN Notices 35(9), 268–279. (Property-based testing methodology adapted here for Ruby/RSpec.)
- (2010). The RSpec Book. Pragmatic Bookshelf.
- (2009). Introduction to Interval Analysis. SIAM. §3 — Containment property and its testing.
- (1999). Validated solutions of initial value problems for ordinary differential equations. Applied Mathematics and Computation 105(1), 21–68.