ulp_report.md

ULP Alignment Report

Auto-generated by pytest tests/test_all.py — aggregate summary of ulp_report.csv

Module	Dtype	Tests	0 ULP	≤tol ULP	Max ULP	Tol	Root Cause
norm.pdf	float64	21	6	15	2	3	std::exp vs scipy npy_exp (both libm, 1-2 ULP diff)
norm.pdf	float32	21	6	15	2	3	std::exp vs scipy npy_exp (both libm, 1-2 ULP diff)
norm.cdf	float64	14	14	0	0	0	Cephes erfc — std::exp same libm as scipy Cephes → 0 ULP
norm.cdf	float32	14	14	0	0	0	Cephes erfc — std::exp same libm as scipy Cephes → 0 ULP
norm.ppf	float64	11	11	0	0	0	Cephes ndtri — std::log/sqrt same libm as scipy Cephes → 0 ULP
norm.ppf	float32	11	11	0	0	0	Cephes ndtri — std::log/sqrt same libm as scipy Cephes → 0 ULP
integrate (trapezoid)	float64	2	2	0	0	0	Pure arithmetic → 0 ULP
integrate (trapezoid)	float32	2	2	0	0	0	Pure arithmetic → 0 ULP
integrate (simpson)	float64	1	1	0	0	0	Pure arithmetic → 0 ULP
integrate (simpson)	float32	1	1	0	0	0	Pure arithmetic → 0 ULP
linalg.solve	float64	107	3	104	31715	atol=1e-10	Eigen3 partialPivLu vs LAPACK gesv: different pivot/decomp paths
linalg.solve	float32	107	4	103	88991	atol=1e-10	Eigen3 partialPivLu vs LAPACK gesv: different pivot/decomp paths
cdist	float64	4	4	0	0	0	Pure arithmetic, float32 promoted to float64 → 0 ULP
cdist	float32	4	4	0	0	0	Pure arithmetic, float32 promoted to float64 → 0 ULP
KDTree	float64	2	2	0	0	0	scipy's own ckdtree, internal float64 → 0 ULP
KDTree	float32	2	2	0	0	0	scipy's own ckdtree, internal float64 → 0 ULP
gaussian_filter1d	float64	9	9	0	0	0	Kernel built in Python numpy, C++ convolution only → 0 ULP
gaussian_filter1d	float32	9	9	0	0	0	Kernel built in Python numpy, C++ convolution only → 0 ULP
medfilt	float64	4	4	0	0	0	Sort-based median → 0 ULP
medfilt	float32	4	4	0	0	0	Sort-based median → 0 ULP
Rotation	float64	886	886	0	0	0	Delegates to scipy.spatial.transform.Rotation → 0 ULP
Rotation	float32	886	886	0	0	0	Delegates to scipy.spatial.transform.Rotation → 0 ULP

Key Takeaways

1. norm.pdf is the only API with non-zero ULP differences (max 2 ULP). numpy::exp() → std::exp() (libm) vs scipy's npy_exp (also libm). Same library, different compiler flags (-ffp-contract=off -ffloat-store). pdf(0), pdf(1), pdf(-1), extreme values, and tiny-scale are exact (0 ULP); non-trivial inputs show 1-2 ULP.

2. norm.cdf / norm.ppf are 0 ULP. Cephes erfc/ndtri use std::exp/log/sqrt directly. scipy's own Cephes is compiled with #define exp npy_exp — both resolve to the same libm library.

3. linalg.solve tolerates up to atol=1e-10 absolute difference. Eigen3 partialPivLu and LAPACK gesv are fundamentally different LU decomposition algorithms. Well-conditioned small matrices (2×2, identity) are 0 ULP.

4. All other APIs (integrate, cdist, KDTree, gaussian_filter1d, medfilt, Rotation) are 0 ULP — either pure arithmetic, float64-internal, or delegated to scipy's own implementation.