crmvec

Why

A correctly rounded function has exactly one right answer for every input, so every correct implementation agrees, on every machine and every version. Vector math libraries aim for speed instead, within a few ulps, and their last bits move with the library, its version and the CPU's instruction set. Through PoCL, for example, glibc's libmvec returns something other than the correctly rounded result for 42% of random double exp inputs and 63% of double sin inputs: within OpenCL's bounds, and not reproducible across libraries.

crmvec gives those same entry points correctly rounded results, so a vectorized program computes the same bits on any machine, in any of the four rounding modes.

Its scalar preload does the same for the C library's libm. Against 14 LTS distributions, glibc's libm misrounds somewhere in 30 to 53 of the 54 functions measured, and its results depend on the CPU in up to 28 of them: what it changes on yours.

What it replaces

• glibc's libmvec.so.1 on x86-64 and aarch64 (checked natively on both);
• SLEEF's libsleefgnuabi.so.3 on aarch64;
• SLEEF's RVV libsleef.so.3 on riscv64 (checked under emulation and natively on a SpacemiT X60);
• since 0.9.0 on x86-64, and since 0.11.0 on aarch64, a program's scalar libm calls too: libcrpreload.so has CORE-MATH's 76 functions under the C library's names, preloaded by crmvec-run --libm, with glibc's errno and exception flags, on glibc 2.17 and later;
• since 0.11.0, ways to turn it on without changing the operating system: module files for clusters, a container recipe of one mounted file and one variable, and crmvec-stamp, which records in a run's outputs which math it used (docs/lab.md);
• since 0.10.0, on x86-64, a fast mode: fast/libmvec.so.1, run with crmvec-run --fast. It is not correctly rounded: each function stays within OpenCL's accuracy bound for it, at a median of 1.02 times glibc's time, with the same bits on every CPU with AVX2 and FMA. Since 0.12.0 those results are versioned: a release that changes any of them says so. Through PoCL, the OpenCL conformance tests' math checks pass with it for the functions PoCL vectorizes (with glibc's libmvec, exp and log fail).

26 functions, each in float and double:

sin cos tan asin acos atan atan2 sinh cosh tanh asinh acosh atanh exp exp2 exp10 expm1 log log2 log10 log1p pow cbrt hypot erf erfc

On x86 that is every symbol LLVM's vectorizer can call through libmvec, in SSE2, AVX, AVX2 and AVX-512 forms; a program needing a symbol the library lacks fails to link rather than falling back silently. OpenCL's other correctly rounded functions are there too, under their C23 names: sinpi, cospi, tanpi, asinpi, acospi, atanpi, atan2pi, lgamma, tgamma, rsqrt, powr, pown, and half-precision and bfloat16 arrays.

Using it

make lib                   # the libraries (gcc on x86)
make install PREFIX=/usr/local
crmvec-run your-program    # its vector math from crmvec
crmvec-run --libm your-program   # its scalar libm calls too (x86-64, aarch64)
crmvec-stamp OUTPUT_DIR          # record which math a run used
crmvec-run --fast your-program   # the fast mode: within OpenCL's bounds, not correctly rounded (x86-64)
pkg-config --cflags --libs crmvec

The libraries install to lib/crmvec/, a directory of their own, and replace the system's only for programs run with crmvec-run or linked with pkg-config's rpath, never system-wide. Recipes for Debian and Ubuntu, Fedora and Nix are in the repository, and it is on conda-forge: conda install -c conda-forge crmvec (crmvec-feedstock).

How it is checked

• Every float input: all 2³² inputs of each one-argument float function, through every entry point, against CORE-MATH.
• Doubles: 2³¹ random inputs per function (2³⁸ in round-to-nearest on a large machine), CORE-MATH's hard-to-round cases and edge values, against CORE-MATH; and a separate check against MPFR in all four rounding modes, through both x86 entry points.
• Under flush-to-zero, as -ffast-math sets it, and on emulated older CPUs (no AVX, no AVX2).
• Controls: each check has one that must fail, such as a deliberately wrong function or a rounding test switched off, so a check that tests nothing says so.
• make check runs a few minutes of this on every push, on x86-64 and natively on arm64.

Speed

Correct rounding costs speed. On one Zen 3 core (AVX2 entry points, nanoseconds per element), crmvec 0.7 took a median 3.3 times glibc's time, from 1.5 times (double asinh) to 9.9 times (double erfc); the table is that version's. glibc computes floats with single precision; correct rounding needs about nine bits more. Since 0.8.0, 18 float functions get them from float pairs on eight lanes, with a rounding test and CORE-MATH for the rare inputs it cannot decide, instead of double precision on four lanes; they take 22% to 90% of 0.7's time, with the same results (release notes). Every function but one is faster than calling scalar CORE-MATH in a loop.

The README has all 52 functions, the other CPUs (aarch64, riscv64, AVX-512) and how each was measured.

Limits

• The x86 library needs gcc; clang can build the checks but not the library.
• On riscv64 it is correct but slower than SLEEF.
• Natively checked on a few core types (Zen 3, Zen 4, Cascade Lake, Neoverse N1 and N2, SpacemiT X60), others under emulation.