fix: route the last three transcendentals through libm, and enforce it

CLAUDE.md requires transcendentals to go through libm rather than std,
because std delegates to the system math library and the two disagree by
one ULP often enough to change an iteration count in a fixed point.

Three production sites did not:

  factor/margin.rs:84  cavity.sigma().hypot(sigma)   12.136% of 1e6 inputs
  factor/margin.rs:92  f64::MIN_POSITIVE.ln()         3.437%
  factor/trunc.rs:98   f64::MIN_POSITIVE.ln()         same

`hypot` is the material one: it is on the path of every scored event, and
its divergence rate is HIGHER than the 9.7% the rule cites for `exp` as
its own justification. The two `ln` calls happen to agree bit-for-bit on
this host, which is exactly the platform dependence the rule exists to
remove.

The `hypot` choice itself was right and stays — the comment above it
explains why, and it is measured: naive sqrt(a^2 + b^2) overflows to inf
at 1e200 and flushes to zero at 1e-200 where hypot does neither. Only the
implementation moves.

tests/libm_rule.rs enforces it. The rule was stated plainly in CLAUDE.md
and still violated three times, so prose is evidently not sufficient. The
test strips `#[cfg(test)]` items by brace matching, plus comments and
string literals so prose is not mistaken for a call, then scans for std
method spellings. `sqrt` is exempt: IEEE 754 specifies it, so std and
libm cannot disagree.

Confirmed non-vacuous by reintroducing the `hypot` violation and watching
it fail with the offending line, then pass again on restore. Two further
tests pin the stripper itself, since a stripper that removed everything
would make the guard pass on anything.

Closes #63

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_011hcFjNDmHXZF8URGLku5zZ

src/factor/margin.rs

@@ -81,7 +81,7 @@ fn cavity_log_evidence(cavity: Gaussian, m_obs: f64, sigma: f64) -> f64 {
     // `hypot`, not `sqrt(a^2 + b^2)`: squaring overflows to infinity above a
     // sigma of ~1.3e154 and flushes to zero below ~1.5e-154, and `Gaussian`'s
     // constructors are public so a caller can reach both.
-    let combined_sigma = cavity.sigma().hypot(sigma);
+    let combined_sigma = libm::hypot(cavity.sigma(), sigma);
     let value = ln_pdf(m_obs, cavity.mu(), combined_sigma);
 
     // A degenerate cavity (infinite sigma) is the only way to reach a
@@ -89,7 +89,7 @@ fn cavity_log_evidence(cavity: Gaussian, m_obs: f64, sigma: f64) -> f64 {
     if value.is_finite() {
         value
     } else {
-        f64::MIN_POSITIVE.ln()
+        libm::log(f64::MIN_POSITIVE)
     }
 }
 

src/factor/trunc.rs

@@ -95,7 +95,7 @@ fn cavity_log_evidence(diff: Gaussian, margin: f64, tie: bool) -> f64 {
     if value.is_finite() {
         value
     } else {
-        f64::MIN_POSITIVE.ln()
+        libm::log(f64::MIN_POSITIVE)
     }
 }
 
@@ -203,7 +203,7 @@ mod tests {
             let approx = -0.5 * z * z - z.ln() - (2.0 * std::f64::consts::PI).sqrt().ln();
 
             assert!(
-                got < f64::MIN_POSITIVE.ln(),
+                got < libm::log(f64::MIN_POSITIVE),
                 "mu={mu}: {got} is still stuck on the old clamp floor"
             );
             assert!(

tests/libm_rule.rs

@@ -0,0 +1,175 @@
+//! The libm rule, enforced rather than asserted in prose.
+//!
+//! CLAUDE.md requires transcendentals to go through `libm`, not `std`:
+//!
+//! > IEEE 754 pins the basic operations and `sqrt` but says nothing about
+//! > `exp`/`log`/`erf`, and `std` delegates to the *system* math library —
+//! > measured, `f64::exp` and `libm::exp` disagree on 9.7% of inputs by one
+//! > ULP. Since inference is an iterative fixed point, one ULP can change an
+//! > iteration count.
+//!
+//! The rule was stated clearly and still violated in three production sites,
+//! one of them `hypot` on the path of every scored event — whose measured
+//! divergence, 12.1%, is *higher* than the `exp` figure the rule cites as its
+//! own justification. Prose is evidently not enough, so this is a test.
+//!
+//! Tests may use either, which the crate documents, so `#[cfg(test)]` blocks
+//! are excluded.
+
+use std::{fs, path::Path};
+
+/// Method-call spellings that reach the system math library.
+///
+/// `sqrt` is deliberately absent: IEEE 754 specifies it exactly, so `std` and
+/// `libm` cannot disagree. `abs`, `recip`, `powi` and `mul_add` are likewise
+/// exact or specified.
+const FORBIDDEN: &[&str] = &[
+    "exp", "exp2", "exp_m1", "ln", "ln_1p", "log", "log2", "log10", "powf", "sin", "cos", "tan",
+    "asin", "acos", "atan", "atan2", "sinh", "cosh", "tanh", "hypot", "cbrt", "erf", "erfc",
+];
+
+/// Strip `#[cfg(test)]` items by brace matching, plus comments and string
+/// literals, so a mention in prose is not mistaken for a call.
+fn production_code(source: &str) -> String {
+    let mut out = String::with_capacity(source.len());
+    let bytes: Vec<char> = source.chars().collect();
+    let mut i = 0;
+
+    while i < bytes.len() {
+        let rest: String = bytes[i..].iter().take(16).collect();
+
+        if rest.starts_with("#[cfg(test)]") {
+            // Skip to the opening brace of the guarded item, then past its
+            // matching close.
+            let mut j = i;
+            while j < bytes.len() && bytes[j] != '{' {
+                j += 1;
+            }
+            let mut depth = 0usize;
+            while j < bytes.len() {
+                match bytes[j] {
+                    '{' => depth += 1,
+                    '}' => {
+                        depth -= 1;
+                        if depth == 0 {
+                            j += 1;
+                            break;
+                        }
+                    }
+                    _ => {}
+                }
+                j += 1;
+            }
+            i = j;
+            continue;
+        }
+
+        if rest.starts_with("//") {
+            while i < bytes.len() && bytes[i] != '\n' {
+                i += 1;
+            }
+            continue;
+        }
+
+        if rest.starts_with("/*") {
+            i += 2;
+            while i + 1 < bytes.len() && !(bytes[i] == '*' && bytes[i + 1] == '/') {
+                i += 1;
+            }
+            i += 2;
+            continue;
+        }
+
+        if bytes[i] == '"' {
+            i += 1;
+            while i < bytes.len() && bytes[i] != '"' {
+                if bytes[i] == '\\' {
+                    i += 1;
+                }
+                i += 1;
+            }
+            i += 1;
+            continue;
+        }
+
+        out.push(bytes[i]);
+        i += 1;
+    }
+
+    out
+}
+
+fn rust_files(dir: &Path, out: &mut Vec<std::path::PathBuf>) {
+    for entry in fs::read_dir(dir).expect("read src") {
+        let path = entry.expect("dir entry").path();
+        if path.is_dir() {
+            rust_files(&path, out);
+        } else if path.extension().is_some_and(|e| e == "rs") {
+            out.push(path);
+        }
+    }
+}
+
+#[test]
+fn production_code_never_calls_a_std_transcendental() {
+    let mut files = Vec::new();
+    rust_files(Path::new("src"), &mut files);
+    assert!(files.len() > 10, "expected to find the crate's sources");
+
+    let mut offences = Vec::new();
+
+    for path in &files {
+        let source = fs::read_to_string(path).expect("read source");
+        let code = production_code(&source);
+
+        for (n, line) in code.lines().enumerate() {
+            for name in FORBIDDEN {
+                let needle = format!(".{name}(");
+                if line.contains(&needle) {
+                    offences.push(format!("{}:{}: {}", path.display(), n + 1, line.trim()));
+                }
+            }
+        }
+    }
+
+    assert!(
+        offences.is_empty(),
+        "production code must call libm, not std, for transcendentals \
+         (`sqrt` is exempt — IEEE 754 specifies it):\n{}",
+        offences.join("\n")
+    );
+}
+
+/// The stripper has to actually strip, or the test above passes vacuously.
+#[test]
+fn the_test_module_stripper_works() {
+    let source = r#"
+fn production() { let _ = libm::exp(1.0); }
+
+#[cfg(test)]
+mod tests {
+    fn allowed() { let x = 1.0f64.exp(); }
+}
+
+fn also_production() {}
+"#;
+    let code = production_code(source);
+    assert!(
+        code.contains("also_production"),
+        "stripped too much: {code}"
+    );
+    assert!(
+        !code.contains(".exp()"),
+        "failed to strip cfg(test): {code}"
+    );
+}
+
+/// And it must not strip a doc comment's worth of prose into oblivion, nor
+/// mistake prose for a call.
+#[test]
+fn prose_is_not_mistaken_for_a_call() {
+    let source = "/// Uses `x.exp()` in the docs.\nfn f() { let _ = libm::exp(1.0); }\n";
+    let code = production_code(source);
+    assert!(!code.contains(".exp()"), "doc comment leaked: {code}");
+    assert!(code.contains("libm::exp"), "stripped real code: {code}");
+}