DEOS Floating-Point Determinism Specification

Version 2.0, December 2025

Summary

The kernel guarantees bit-identical floating-point computation across x86-64 CPUs. The way it does this is mostly by not using floating point. Neural network code uses software fixed-point in Q48.16. Residual FP operations (loading weights from disk, formatting numbers for humans) sit behind a pinned MXCSR and a compiler that does not emit the instructions whose semantics vary across CPUs.

This is not a clever technique. It is a discipline. Every part of the system that touches numbers is either pure integer arithmetic, or it is allowed only at well-defined I/O boundaries where the answer does not feed into the next computation.

Invariants the kernel enforces

MXCSR (the SSE/AVX control register)

The kernel sets MXCSR to 0x1F80 at boot and never changes it.

FTZ and DAZ off because they introduce hardware-specific denormal handling. Round-to-nearest because that is the IEEE 754 default and the only mode the rest of the system assumes.

x87 control word

Pinned to 0x037F: round to nearest, 80-bit extended precision, all exceptions masked.

The kernel avoids x87 instructions, but legacy code or libraries may use them. Setting the control word means even when they do, the behavior is the same as it would be on any other x86 CPU running with these defaults.

Compiler flags

# kernel/.cargo/config.toml
rustflags = [
    "-C", "target-feature=-fma",         # disable FMA instructions
    "-C", "target-feature=-avx512f,...", # disable AVX-512
]

FMA is disabled because a*b+c with FMA uses single rounding, without it uses double rounding, and not all x86 CPUs have FMA. Compiling with FMA disabled gives the same instructions everywhere.

AVX-512 is disabled because not all x86 CPUs have it. Staying on SSE2 plus AVX2 baseline keeps the instruction stream identical across machines.

Fixed-point arithmetic

Neural network code uses one type:

pub struct FixedPoint(i64);  // 48 integer bits, 16 fractional bits

impl FixedPoint {
    pub const SCALE: i64 = 65536;  // 2^16

    pub fn mul(self, other: Self) -> Self {
        let product = (self.0 as i128) * (other.0 as i128);
        FixedPoint((product >> 16) as i64)
    }
}

Integer arithmetic is bit-identical across all CPUs. No rounding mode, no denormal handling, no FMA variability. The cost is precision (1/65536 ≈ 0.000015 per step) and dynamic range (±2^47 before saturation). For transformer inference at reasonable hyperparameters this is fine, as the operational bounds section below shows.

Where f32 is allowed

Only at I/O boundaries.

// Loading pre-computed weights from external files
let weight = FixedPoint::from_f32(file.read_f32());

// Debug output for humans
serial_println!("Loss: {:.4}", loss.to_f32());

What it cannot do is feed back into computation:

// Wrong: turns a constant into a non-deterministic value
let half = FixedPoint::from_f32(0.5);

// Right: integer constant, bit-identical
let half = FixedPoint::HALF;  // raw value 32768

All common constants are pre-computed as raw integers:

pub const ZERO:       Self = FixedPoint(0);
pub const ONE:        Self = FixedPoint(65536);
pub const HALF:       Self = FixedPoint(32768);
pub const SQRT_2_PI:  Self = FixedPoint(52279);  // round(0.7978845608 * 65536)
pub const GELU_COEFF: Self = FixedPoint(2930);   // round(0.044715 * 65536)

Transcendentals

exp, log, tanh use compile-time lookup tables with linear interpolation.

static EXP_LUT: [i64; 257] = { /* computed at compile time */ };

pub fn exp_fixed(x: FixedPoint) -> FixedPoint {
    let idx = ((x.raw() + 524288) * 256) / (16 * 65536);
    let y0 = EXP_LUT[idx];
    let y1 = EXP_LUT[idx + 1];
    let delta = y1 - y0;
    let interp = (delta * frac) >> 16;
    FixedPoint::from_raw(y0 + interp)
}

LUT values are baked at compile time. Interpolation is integer only. Same inputs produce the same outputs everywhere.

Verification

A boot-time check that catches the case where something has clobbered the FPU state:

// kernel/src/fpu.rs
pub fn verify_fpu_state() -> bool {
    let mxcsr = get_mxcsr();
    let x87_cw = get_x87_cw();
    mxcsr == 0x1F80 && x87_cw == 0x037F
}

Callable at any point as a debug assertion. In release builds the kernel runs the check periodically as part of the housekeeping pass.

The claim

Given the same DEOS kernel binary, the same model weights stored as raw fixed-point, and the same input data as raw fixed-point, the kernel produces bit-identical output across all x86-64 CPUs.

Tested configurations

Known limits

f32 weight loading. Weights loaded from external f32 files may differ by one ULP across CPUs. The fix is to store weights as raw i64 fixed-point and never round through f32.

Denormal handling. If FTZ or DAZ were enabled, different CPUs would handle denormals differently. Both are disabled.

x87 instructions. If legacy code uses x87, 80-bit intermediates could cause divergence. The x87 control word is pinned, but avoiding x87 entirely is safer.

Numerics contract

Representation

Overflow behavior by operation

Operational bounds for matmul

The theoretical maximum (±2^47) is misleading. Real transformer values never approach it.

For d_model = 768:

Max single term:   262,144 × 655,360 ≈ 2^37
Max accumulation:  768 × 2^37        = 2^47
i128 range:                            2^127  (safe by factor of 2^80)

Enforcement comes from Xavier initialization, LayerNorm unit variance, and optional gradient clipping during training.

Softmax stability

pub fn softmax_fixed(x: &[FixedPoint]) -> Vec<FixedPoint> {
    // 1. Find max so exp does not overflow
    let max_val = x.iter().map(|v| v.raw()).max().unwrap_or(0);

    // 2. Subtract max before exp
    let shifted: Vec<i64> = x.iter()
        .map(|v| v.raw() - max_val)
        .collect();

    // 3. exp of shifted values (all ≤ 0, so exp ≤ 1)
    // 4. Sum and normalize in i128
}

After max-subtraction, all exp inputs are ≤ 0, so exp outputs are in [0, 1]. No overflow path.

LUT specification

exp domain and clamping

Error bounds

Monotonicity proof: LUT[i] < LUT[i+1] by construction; linear interpolation between adjacent points preserves order; therefore exp_fixed(a) < exp_fixed(b) whenever a < b.

tanh range

pub fn tanh_fixed(x: FixedPoint) -> FixedPoint {
    // tanh(x) = (exp(2x) - 1) / (exp(2x) + 1)
    // Range: [-1, 1], so fixed-point [-65536, 65536]
    // No overflow path
}

Revision history