Crate dreamwell_math 

dreamwell-math — Classical Mathematics for Simulation

Pure mathematical primitives for the Dreamwell engine and research ecosystem. No engine dependencies. No allocation on hot path. Standalone usable.

Modules

• Complex: Complex number arithmetic (add, mul, conj, exp_i, norm)

• Spatial: Distance metrics, Gaussian kernels, AABB extraction

• Wave: Traveling wave, standing wave, wave packets, dispersion

• Fibonacci: Golden ratio, phyllotaxis disc, Fibonacci sphere

• Interpolation: Lerp, exponential decay, smoothstep

• Random: Deterministic splitmix64 PRNG

• Linear Algebra: Complex matrix multiply (CGEMM), tiled blocking

• Eigenvalue: Jacobi rotation solver for Hermitian matrices

• Matrix Exponential: Scaling-and-squaring with Padé approximant

Clean Compute Compliance

• Zero heap allocation in caller-provides-output functions
• All functions are #[inline] where beneficial
• No external dependencies
• Deterministic: same inputs always produce same outputs

Modules

eigen

Eigenvalue Decomposition — Jacobi rotation solver for Hermitian matrices.

linalg

Linear Algebra — Complex matrix operations for quantum simulation.

matrix_exp

Matrix Exponential — Scaling-and-squaring with Padé approximant.

Structs

Complex

Complex number z = re + i·im.

Constants

HALF_NEIGHBORHOOD

Forward half-neighborhood offsets for broadphase pair generation. 14 cells: self + 13 lexicographically forward neighbors. Produces zero-duplicate pairs by construction (Causal Spatial Mask).

Functions

adiabatic_fraction

Adiabatic force distribution: returns the fractional force to apply on the current step within an adiabatic window.

angular_frequency

Angular frequency from frequency: ω = 2πf.

atomic_f32_add

Atomic f32 addition via compare-exchange on AtomicU32.

baumgarte_correction

Baumgarte positional correction for penetration resolution. correction = max(0, depth - slop) * pct / inv_mass_sum

clamp01

Saturating clamp to [0, 1].

contact_impulse

Contact impulse magnitude (sequential impulse formula). j = -(1 + e) * v_along_normal / inv_mass_sum

coulomb_friction_clamp

Coulomb friction clamp: limits tangential impulse to mu * |normal impulse|. jt_clamped = clamp(jt, -|j_normal| * mu, |j_normal| * mu)

dispersion_linear

Linear dispersion relation: ω = c·k (non-dispersive medium).

distance_3d

Euclidean distance between two 3D points.

distance_sq_3d

Squared Euclidean distance between two 3D points. Avoids sqrt for comparison (cheaper than full distance).

ease_in_out_cubic

Cubic ease-in-out: smooth acceleration and deceleration. t in [0, 1]. Returns smooth curve from 0 to 1.

exp_decay

Exponential decay toward target: current + (target - current) · (1 - e^{-rate·dt}).

fibonacci_sphere

Fibonacci sphere: N points uniformly distributed on a unit sphere using the golden ratio spiral.

fnv1a_3d

FNV-1a hash of a 3D integer cell coordinate. Used for spatial hash grid broadphase. Deterministic, fast, low collision rate.

gaussian_kernel

Gaussian radial basis function: exp(-d²/σ²). Returns 1.0 at center, decays to ~0 at distance 3σ. Used for spatial attention and relevance scoring.

golden_angle

The golden angle θ = 2π/φ² ≈ 2.3999632297 radians ≈ 137.508°. Used in phyllotaxis for optimal packing on discs and spheres.

golden_ratio

The golden ratio φ = (1 + √5) / 2 ≈ 1.6180339887.

lerp

Linear interpolation: a + (b - a) · t.

lerp3

3D linear interpolation.

phyllotaxis_disc

Phyllotaxis disc: N points on a flat disc with golden angle spacing. Radius grows as √(index) for uniform area density.

remap

Range remapping: maps value from [in_min, in_max] to [out_min, out_max].

slerp

Quaternion spherical linear interpolation (slerp). Interpolates between two unit quaternions [x, y, z, w] along the great arc.

smoothstep

Smoothstep: Hermite interpolation between edge0 and edge1. Returns 0 below edge0, 1 above edge1, smooth curve between.

splitmix64

Splitmix64: deterministic PRNG. Returns (value, next_state).

splitmix64_f32

Splitmix64 -> f32 in [0, 1). Deterministic.

standing_wave

Standing wave: y = A·sin(kx)·cos(ωt). Created by superposition of two counter-propagating traveling waves.

traveling_wave

Traveling wave: y = A·sin(kx - ωt). k = wave number (2π/λ), ω = angular frequency (2πf). Propagates in +x direction.

traveling_wave_2d

2D traveling wave with perpendicular tapering: y = A·sin(k·along - ωt)·cos(k·perp·taper)

union_find_merge

Union-find merge by rank. O(alpha(n)) amortized.

union_find_root

Union-find root with path compression. O(alpha(n)) amortized. Used for island splitting in the Parallel Jacobi solver.

vec3_cross

3D cross product.

vec3_dot

3D dot product.

vec3_len

3D vector length.

vec3_normalize

3D vector normalize (returns zero vector if length < epsilon).

vec3_scale

3D vector scale.

wave_number

Wave number from wavelength: k = 2π/λ.

wave_packet_gaussian

Gaussian wave packet: ψ(x,t) = exp(-(x-x0-v_g·t)²/(4σ²)) · exp(i(k0·x - ω0·t)) Returns the real part (amplitude envelope × carrier).