diff --git a/ci/docker/aarch64-unknown-linux-gnu/Dockerfile b/ci/docker/aarch64-unknown-linux-gnu/Dockerfile index a7b23cb9..7fa06b28 100644 --- a/ci/docker/aarch64-unknown-linux-gnu/Dockerfile +++ b/ci/docker/aarch64-unknown-linux-gnu/Dockerfile @@ -3,7 +3,7 @@ FROM ubuntu:24.04 RUN apt-get update && \ apt-get install -y --no-install-recommends \ gcc libc6-dev ca-certificates \ - gcc-aarch64-linux-gnu libc6-dev-arm64-cross \ + gcc-aarch64-linux-gnu m4 make libc6-dev-arm64-cross \ qemu-user-static ENV TOOLCHAIN_PREFIX=aarch64-linux-gnu- diff --git a/ci/docker/i686-unknown-linux-gnu/Dockerfile b/ci/docker/i686-unknown-linux-gnu/Dockerfile index 3b0bfc0d..37e206a8 100644 --- a/ci/docker/i686-unknown-linux-gnu/Dockerfile +++ b/ci/docker/i686-unknown-linux-gnu/Dockerfile @@ -2,4 +2,4 @@ FROM ubuntu:24.04 RUN apt-get update && \ apt-get install -y --no-install-recommends \ - gcc-multilib libc6-dev ca-certificates + gcc-multilib m4 make libc6-dev ca-certificates diff --git a/ci/docker/x86_64-unknown-linux-gnu/Dockerfile b/ci/docker/x86_64-unknown-linux-gnu/Dockerfile index 15723ab5..c84a31c5 100644 --- a/ci/docker/x86_64-unknown-linux-gnu/Dockerfile +++ b/ci/docker/x86_64-unknown-linux-gnu/Dockerfile @@ -2,4 +2,4 @@ FROM ubuntu:24.04 RUN apt-get update && \ apt-get install -y --no-install-recommends \ - gcc libc6-dev ca-certificates + gcc m4 make libc6-dev ca-certificates diff --git a/ci/run.sh b/ci/run.sh index 30265e51..94612adc 100755 --- a/ci/run.sh +++ b/ci/run.sh @@ -35,6 +35,22 @@ case "$target" in *) extra_flags="$extra_flags --features libm-test/build-musl" ;; esac +# Configure which targets test against MPFR +case "$target" in + # MSVC cannot link MPFR + *windows-msvc*) ;; + # FIXME: MinGW should be able to build MPFR, but setup in CI is nontrivial. + *windows-gnu*) ;; + # Targets that aren't cross compiled work fine + # FIXME(ci): we should be able to enable aarch64 Linux here once GHA + # support rolls out. + x86_64*) extra_flags="$extra_flags --features libm-test/test-multiprecision" ;; + # i686 works fine, i586 does not + i686*) extra_flags="$extra_flags --features libm-test/test-multiprecision" ;; + # Apple aarch64 is native + aarch64*apple*) extra_flags="$extra_flags --features libm-test/test-multiprecision" ;; +esac + # FIXME: `STATUS_DLL_NOT_FOUND` testing macros on CI. # case "$target" in diff --git a/crates/libm-test/Cargo.toml b/crates/libm-test/Cargo.toml index 703524bc..72ac5723 100644 --- a/crates/libm-test/Cargo.toml +++ b/crates/libm-test/Cargo.toml @@ -10,18 +10,21 @@ default = [] # Generate tests which are random inputs and the outputs are calculated with # musl libc. test-musl-serialized = ["rand"] +test-multiprecision = ["dep:az", "dep:rug"] # Build our own musl for testing and benchmarks build-musl = ["dep:musl-math-sys"] [dependencies] anyhow = "1.0.90" +az = { version = "1.2.1", optional = true } libm = { path = "../.." } libm-macros = { path = "../libm-macros" } musl-math-sys = { path = "../musl-math-sys", optional = true } paste = "1.0.15" rand = "0.8.5" rand_chacha = "0.3.1" +rug = { version = "1.26.1", optional = true, default-features = false, features = ["float", "std"] } [target.'cfg(target_family = "wasm")'.dependencies] # Enable randomness on WASM diff --git a/crates/libm-test/src/gen/random.rs b/crates/libm-test/src/gen/random.rs index 601ef4f1..c73937aa 100644 --- a/crates/libm-test/src/gen/random.rs +++ b/crates/libm-test/src/gen/random.rs @@ -7,7 +7,7 @@ use rand::{Rng, SeedableRng}; use rand_chacha::ChaCha8Rng; use super::CachedInput; -use crate::GenerateInput; +use crate::{CheckCtx, GenerateInput}; const SEED: [u8; 32] = *b"3.141592653589793238462643383279"; @@ -40,9 +40,10 @@ static TEST_CASES_JN: LazyLock = LazyLock::new(|| { let mut cases = (&*TEST_CASES).clone(); // These functions are extremely slow, limit them - cases.inputs_i32.truncate((NTESTS / 1000).max(80)); - cases.inputs_f32.truncate((NTESTS / 1000).max(80)); - cases.inputs_f64.truncate((NTESTS / 1000).max(80)); + let ntests_jn = (NTESTS / 1000).max(80); + cases.inputs_i32.truncate(ntests_jn); + cases.inputs_f32.truncate(ntests_jn); + cases.inputs_f64.truncate(ntests_jn); // It is easy to overflow the stack with these in debug mode let max_iterations = if cfg!(optimizations_enabled) && cfg!(target_pointer_width = "64") { @@ -105,11 +106,10 @@ fn make_test_cases(ntests: usize) -> CachedInput { } /// Create a test case iterator. -pub fn get_test_cases(fname: &str) -> impl Iterator +pub fn get_test_cases(ctx: &CheckCtx) -> impl Iterator where CachedInput: GenerateInput, { - let inputs = if fname == "jn" || fname == "jnf" { &TEST_CASES_JN } else { &TEST_CASES }; - - CachedInput::get_cases(inputs) + let inputs = if ctx.fname == "jn" || ctx.fname == "jnf" { &TEST_CASES_JN } else { &TEST_CASES }; + inputs.get_cases() } diff --git a/crates/libm-test/src/lib.rs b/crates/libm-test/src/lib.rs index 2abe7f60..13b76d6c 100644 --- a/crates/libm-test/src/lib.rs +++ b/crates/libm-test/src/lib.rs @@ -1,4 +1,6 @@ pub mod gen; +#[cfg(feature = "test-multiprecision")] +pub mod mpfloat; mod num_traits; mod special_case; mod test_traits; @@ -14,14 +16,18 @@ pub type TestResult = Result; // List of all files present in libm's source include!(concat!(env!("OUT_DIR"), "/all_files.rs")); -/// ULP allowed to differ from musl (note that musl itself may not be accurate). +/// Default ULP allowed to differ from musl (note that musl itself may not be accurate). const MUSL_DEFAULT_ULP: u32 = 2; -/// Certain functions have different allowed ULP (consider these xfail). +/// Default ULP allowed to differ from multiprecision (i.e. infinite) results. +const MULTIPREC_DEFAULT_ULP: u32 = 1; + +/// ULP allowed to differ from muls results. /// /// Note that these results were obtained using 400,000,000 rounds of random inputs, which /// is not a value used by default. pub fn musl_allowed_ulp(name: &str) -> u32 { + // Consider overrides xfail match name { #[cfg(x86_no_sse)] "asinh" | "asinhf" => 6, @@ -42,6 +48,27 @@ pub fn musl_allowed_ulp(name: &str) -> u32 { } } +/// ULP allowed to differ from multiprecision results. +pub fn multiprec_allowed_ulp(name: &str) -> u32 { + // Consider overrides xfail + match name { + "asinh" | "asinhf" => 2, + "acoshf" => 4, + "atanh" | "atanhf" => 2, + "exp10" | "exp10f" => 3, + "j0" | "j0f" | "j1" | "j1f" => { + // Results seem very target-dependent + if cfg!(target_arch = "x86_64") { 4000 } else { 800_000 } + } + "jn" | "jnf" => 1000, + "lgamma" | "lgammaf" | "lgamma_r" | "lgammaf_r" => 16, + "sinh" | "sinhf" => 2, + "tanh" | "tanhf" => 2, + "tgamma" => 20, + _ => MULTIPREC_DEFAULT_ULP, + } +} + /// Return the unsuffixed version of a function name; e.g. `abs` and `absf` both return `abs`, /// `lgamma_r` and `lgammaf_r` both return `lgamma_r`. pub fn canonical_name(name: &str) -> &str { diff --git a/crates/libm-test/src/mpfloat.rs b/crates/libm-test/src/mpfloat.rs new file mode 100644 index 00000000..44962d11 --- /dev/null +++ b/crates/libm-test/src/mpfloat.rs @@ -0,0 +1,368 @@ +//! Interfaces needed to support testing with multi-precision floating point numbers. +//! +//! Within this module, the macros create a submodule for each `libm` function. These contain +//! a struct named `Operation` that implements [`MpOp`]. + +use std::cmp::Ordering; + +use az::Az; +use rug::Assign; +pub use rug::Float as MpFloat; +use rug::float::Round::Nearest; +use rug::ops::{PowAssignRound, RemAssignRound}; + +use crate::Float; + +/// Create a multiple-precision float with the correct number of bits for a concrete float type. +fn new_mpfloat() -> MpFloat { + MpFloat::new(F::SIGNIFICAND_BITS + 1) +} + +/// Set subnormal emulation and convert to a concrete float type. +fn prep_retval(mp: &mut MpFloat, ord: Ordering) -> F +where + for<'a> &'a MpFloat: az::Cast, +{ + mp.subnormalize_ieee_round(ord, Nearest); + (&*mp).az::() +} + +/// Structures that represent a float operation. +/// +/// The struct itself should hold any context that can be reused among calls to `run` (allocated +/// `MpFloat`s). +pub trait MpOp { + /// Inputs to the operation (concrete float types). + type Input; + + /// Outputs from the operation (concrete float types). + type Output; + + /// Create a new instance. + fn new() -> Self; + + /// Perform the operation. + /// + /// Usually this means assigning inputs to cached floats, performing the operation, applying + /// subnormal approximation, and converting the result back to concrete values. + fn run(&mut self, input: Self::Input) -> Self::Output; +} + +/// Implement `MpOp` for functions with a single return value. +macro_rules! impl_mp_op { + // Matcher for unary functions + ( + fn_name: $fn_name:ident, + CFn: $CFn:ty, + CArgs: $CArgs:ty, + CRet: $CRet:ty, + RustFn: fn($fty:ty,) -> $_ret:ty, + RustArgs: $RustArgs:ty, + RustRet: $RustRet:ty, + fn_extra: $fn_name_normalized:expr, + ) => { + paste::paste! { + pub mod $fn_name { + use super::*; + pub struct Operation(MpFloat); + + impl MpOp for Operation { + type Input = $RustArgs; + type Output = $RustRet; + + fn new() -> Self { + Self(new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0.assign(input.0); + let ord = self.0.[< $fn_name_normalized _round >](Nearest); + prep_retval::(&mut self.0, ord) + } + } + } + } + }; + // Matcher for binary functions + ( + fn_name: $fn_name:ident, + CFn: $CFn:ty, + CArgs: $CArgs:ty, + CRet: $CRet:ty, + RustFn: fn($fty:ty, $_fty2:ty,) -> $_ret:ty, + RustArgs: $RustArgs:ty, + RustRet: $RustRet:ty, + fn_extra: $fn_name_normalized:expr, + ) => { + paste::paste! { + pub mod $fn_name { + use super::*; + pub struct Operation(MpFloat, MpFloat); + + impl MpOp for Operation { + type Input = $RustArgs; + type Output = $RustRet; + + fn new() -> Self { + Self(new_mpfloat::<$fty>(), new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0.assign(input.0); + self.1.assign(input.1); + let ord = self.0.[< $fn_name_normalized _round >](&self.1, Nearest); + prep_retval::(&mut self.0, ord) + } + } + } + } + }; + // Matcher for ternary functions + ( + fn_name: $fn_name:ident, + CFn: $CFn:ty, + CArgs: $CArgs:ty, + CRet: $CRet:ty, + RustFn: fn($fty:ty, $_fty2:ty, $_fty3:ty,) -> $_ret:ty, + RustArgs: $RustArgs:ty, + RustRet: $RustRet:ty, + fn_extra: $fn_name_normalized:expr, + ) => { + paste::paste! { + pub mod $fn_name { + use super::*; + pub struct Operation(MpFloat, MpFloat, MpFloat); + + impl MpOp for Operation { + type Input = $RustArgs; + type Output = $RustRet; + + fn new() -> Self { + Self(new_mpfloat::<$fty>(), new_mpfloat::<$fty>(), new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0.assign(input.0); + self.1.assign(input.1); + self.2.assign(input.2); + let ord = self.0.[< $fn_name_normalized _round >](&self.1, &self.2, Nearest); + prep_retval::(&mut self.0, ord) + } + } + } + } + }; +} + +libm_macros::for_each_function! { + callback: impl_mp_op, + skip: [ + // Most of these need a manual implementation + fabs, ceil, copysign, floor, rint, round, trunc, + fabsf, ceilf, copysignf, floorf, rintf, roundf, truncf, + fmod, fmodf, frexp, frexpf, ilogb, ilogbf, jn, jnf, ldexp, ldexpf, + lgamma_r, lgammaf_r, modf, modff, nextafter, nextafterf, pow,powf, + remquo, remquof, scalbn, scalbnf, sincos, sincosf, + ], + fn_extra: match MACRO_FN_NAME { + // Remap function names that are different between mpfr and libm + expm1 | expm1f => exp_m1, + fabs | fabsf => abs, + fdim | fdimf => positive_diff, + fma | fmaf => mul_add, + fmax | fmaxf => max, + fmin | fminf => min, + lgamma | lgammaf => ln_gamma, + log | logf => ln, + log1p | log1pf => ln_1p, + tgamma | tgammaf => gamma, + _ => MACRO_FN_NAME_NORMALIZED + } +} + +/// Implement unary functions that don't have a `_round` version +macro_rules! impl_no_round { + // Unary matcher + ($($fn_name:ident, $rug_name:ident;)*) => { + paste::paste! { + // Implement for both f32 and f64 + $( impl_no_round!{ @inner_unary [< $fn_name f >], (f32,), $rug_name } )* + $( impl_no_round!{ @inner_unary $fn_name, (f64,), $rug_name } )* + } + }; + + (@inner_unary $fn_name:ident, ($fty:ty,), $rug_name:ident) => { + pub mod $fn_name { + use super::*; + pub struct Operation(MpFloat); + + impl MpOp for Operation { + type Input = ($fty,); + type Output = $fty; + + fn new() -> Self { + Self(new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0.assign(input.0); + self.0.$rug_name(); + prep_retval::(&mut self.0, Ordering::Equal) + } + } + } + }; +} + +impl_no_round! { + fabs, abs_mut; + ceil, ceil_mut; + floor, floor_mut; + rint, round_even_mut; // FIXME: respect rounding mode + round, round_mut; + trunc, trunc_mut; +} + +/// Some functions are difficult to do in a generic way. Implement them here. +macro_rules! impl_op_for_ty { + ($fty:ty, $suffix:literal) => { + paste::paste! { + pub mod [] { + use super::*; + pub struct Operation(MpFloat, MpFloat); + + impl MpOp for Operation { + type Input = ($fty, $fty); + type Output = $fty; + + fn new() -> Self { + Self(new_mpfloat::<$fty>(), new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0.assign(input.0); + self.1.assign(input.1); + self.0.copysign_mut(&self.1); + prep_retval::(&mut self.0, Ordering::Equal) + } + } + } + + pub mod [] { + use super::*; + pub struct Operation(MpFloat, MpFloat); + + impl MpOp for Operation { + type Input = ($fty, $fty); + type Output = $fty; + + fn new() -> Self { + Self(new_mpfloat::<$fty>(), new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0.assign(input.0); + self.1.assign(input.1); + let ord = self.0.pow_assign_round(&self.1, Nearest); + prep_retval::(&mut self.0, ord) + } + } + } + + pub mod [] { + use super::*; + pub struct Operation(MpFloat, MpFloat); + + impl MpOp for Operation { + type Input = ($fty, $fty); + type Output = $fty; + + fn new() -> Self { + Self(new_mpfloat::<$fty>(), new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0.assign(input.0); + self.1.assign(input.1); + let ord = self.0.rem_assign_round(&self.1, Nearest); + prep_retval::(&mut self.0, ord) + } + } + } + + pub mod [] { + use super::*; + pub struct Operation(MpFloat); + + impl MpOp for Operation { + type Input = ($fty,); + type Output = ($fty, i32); + + fn new() -> Self { + Self(new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0.assign(input.0); + let (sign, ord) = self.0.ln_abs_gamma_round(Nearest); + let ret = prep_retval::<$fty>(&mut self.0, ord); + (ret, sign as i32) + } + } + } + + pub mod [] { + use super::*; + pub struct Operation(i32, MpFloat); + + impl MpOp for Operation { + type Input = (i32, $fty); + type Output = $fty; + + fn new() -> Self { + Self(0, new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0 = input.0; + self.1.assign(input.1); + let ord = self.1.jn_round(self.0, Nearest); + prep_retval::<$fty>(&mut self.1, ord) + } + } + } + + pub mod [] { + use super::*; + pub struct Operation(MpFloat, MpFloat); + + impl MpOp for Operation { + type Input = ($fty,); + type Output = ($fty, $fty); + + fn new() -> Self { + Self(new_mpfloat::<$fty>(), new_mpfloat::<$fty>()) + } + + fn run(&mut self, input: Self::Input) -> Self::Output { + self.0.assign(input.0); + self.1.assign(0.0); + let (sord, cord) = self.0.sin_cos_round(&mut self.1, Nearest); + ( + prep_retval::<$fty>(&mut self.0, sord), + prep_retval::<$fty>(&mut self.1, cord) + ) + } + } + } + } + }; +} + +impl_op_for_ty!(f32, "f"); +impl_op_for_ty!(f64, ""); + +// Account for `lgamma_r` not having a simple `f` suffix +pub mod lgammaf_r { + pub use super::lgamma_rf::*; +} diff --git a/crates/libm-test/src/special_case.rs b/crates/libm-test/src/special_case.rs index df263d74..dac7a349 100644 --- a/crates/libm-test/src/special_case.rs +++ b/crates/libm-test/src/special_case.rs @@ -58,20 +58,6 @@ impl MaybeOverride<(f32,)> for SpecialCase { ctx: &CheckCtx, ) -> Option { if ctx.basis == CheckBasis::Musl { - if ctx.fname == "acoshf" && input.0 < -1.0 { - // acoshf is undefined for x <= 1.0, but we return a random result at lower - // values. - return XFAIL; - } - - if ctx.fname == "sincosf" { - let factor_frac_pi_2 = input.0.abs() / f32::consts::FRAC_PI_2; - if (factor_frac_pi_2 - factor_frac_pi_2.round()).abs() < 1e-2 { - // we have a bad approximation near multiples of pi/2 - return XFAIL; - } - } - if ctx.fname == "expm1f" && input.0 > 80.0 && actual.is_infinite() { // we return infinity but the number is representable return XFAIL; @@ -82,15 +68,40 @@ impl MaybeOverride<(f32,)> for SpecialCase { // doesn't seem to happen on x86 return XFAIL; } + } - if ctx.fname == "lgammaf" || ctx.fname == "lgammaf_r" && input.0 < 0.0 { - // loggamma should not be defined for x < 0, yet we both return results - return XFAIL; - } + if ctx.fname == "acoshf" && input.0 < -1.0 { + // acoshf is undefined for x <= 1.0, but we return a random result at lower + // values. + return XFAIL; + } + + if ctx.fname == "lgammaf" || ctx.fname == "lgammaf_r" && input.0 < 0.0 { + // loggamma should not be defined for x < 0, yet we both return results + return XFAIL; } maybe_check_nan_bits(actual, expected, ctx) } + + fn check_int( + input: (f32,), + actual: I, + expected: I, + ctx: &CheckCtx, + ) -> Option> { + // On MPFR for lgammaf_r, we set -1 as the integer result for negative infinity but MPFR + // sets +1 + if ctx.basis == CheckBasis::Mpfr + && ctx.fname == "lgammaf_r" + && input.0 == f32::NEG_INFINITY + && actual.abs() == expected.abs() + { + XFAIL + } else { + None + } + } } impl MaybeOverride<(f64,)> for SpecialCase { @@ -117,15 +128,40 @@ impl MaybeOverride<(f64,)> for SpecialCase { // musl returns -0.0, we return +0.0 return XFAIL; } + } - if ctx.fname == "lgamma" || ctx.fname == "lgamma_r" && input.0 < 0.0 { - // loggamma should not be defined for x < 0, yet we both return results - return XFAIL; - } + if ctx.fname == "acosh" && input.0 < 1.0 { + // The function is undefined for the inputs, musl and our libm both return + // random results. + return XFAIL; + } + + if ctx.fname == "lgamma" || ctx.fname == "lgamma_r" && input.0 < 0.0 { + // loggamma should not be defined for x < 0, yet we both return results + return XFAIL; } maybe_check_nan_bits(actual, expected, ctx) } + + fn check_int( + input: (f64,), + actual: I, + expected: I, + ctx: &CheckCtx, + ) -> Option> { + // On MPFR for lgamma_r, we set -1 as the integer result for negative infinity but MPFR + // sets +1 + if ctx.basis == CheckBasis::Mpfr + && ctx.fname == "lgamma_r" + && input.0 == f64::NEG_INFINITY + && actual.abs() == expected.abs() + { + XFAIL + } else { + None + } + } } /// Check NaN bits if the function requires it @@ -142,6 +178,11 @@ fn maybe_check_nan_bits(actual: F, expected: F, ctx: &CheckCtx) -> Opt return SKIP; } + // MPFR only has one NaN bitpattern; allow the default `.is_nan()` checks to validate. + if ctx.basis == CheckBasis::Mpfr { + return SKIP; + } + // abs and copysign require signaling NaNs to be propagated, so verify bit equality. if actual.to_bits() == expected.to_bits() { return SKIP; @@ -158,9 +199,10 @@ impl MaybeOverride<(f32, f32)> for SpecialCase { _ulp: &mut u32, ctx: &CheckCtx, ) -> Option { - maybe_skip_min_max_nan(input, expected, ctx) + maybe_skip_binop_nan(input, expected, ctx) } } + impl MaybeOverride<(f64, f64)> for SpecialCase { fn check_float( input: (f64, f64), @@ -169,47 +211,86 @@ impl MaybeOverride<(f64, f64)> for SpecialCase { _ulp: &mut u32, ctx: &CheckCtx, ) -> Option { - maybe_skip_min_max_nan(input, expected, ctx) + maybe_skip_binop_nan(input, expected, ctx) } } /// Musl propagates NaNs if one is provided as the input, but we return the other input. // F1 and F2 are always the same type, this is just to please generics -fn maybe_skip_min_max_nan( +fn maybe_skip_binop_nan( input: (F1, F1), expected: F2, ctx: &CheckCtx, ) -> Option { - if (ctx.canonical_name == "fmax" || ctx.canonical_name == "fmin") - && (input.0.is_nan() || input.1.is_nan()) - && expected.is_nan() - { - return XFAIL; - } else { - None + match ctx.basis { + CheckBasis::Musl => { + if (ctx.canonical_name == "fmax" || ctx.canonical_name == "fmin") + && (input.0.is_nan() || input.1.is_nan()) + && expected.is_nan() + { + XFAIL + } else { + None + } + } + CheckBasis::Mpfr => { + if ctx.canonical_name == "copysign" && input.1.is_nan() { + SKIP + } else { + None + } + } } } impl MaybeOverride<(i32, f32)> for SpecialCase { fn check_float( input: (i32, f32), - _actual: F, - _expected: F, + actual: F, + expected: F, ulp: &mut u32, ctx: &CheckCtx, ) -> Option { - bessel_prec_dropoff(input, ulp, ctx) + match ctx.basis { + CheckBasis::Musl => bessel_prec_dropoff(input, ulp, ctx), + CheckBasis::Mpfr => { + // We return +0.0, MPFR returns -0.0 + if ctx.fname == "jnf" + && input.1 == f32::NEG_INFINITY + && actual == F::ZERO + && expected == F::ZERO + { + XFAIL + } else { + None + } + } + } } } impl MaybeOverride<(i32, f64)> for SpecialCase { fn check_float( input: (i32, f64), - _actual: F, - _expected: F, + actual: F, + expected: F, ulp: &mut u32, ctx: &CheckCtx, ) -> Option { - bessel_prec_dropoff(input, ulp, ctx) + match ctx.basis { + CheckBasis::Musl => bessel_prec_dropoff(input, ulp, ctx), + CheckBasis::Mpfr => { + // We return +0.0, MPFR returns -0.0 + if ctx.fname == "jn" + && input.1 == f64::NEG_INFINITY + && actual == F::ZERO + && expected == F::ZERO + { + XFAIL + } else { + bessel_prec_dropoff(input, ulp, ctx) + } + } + } } } diff --git a/crates/libm-test/src/test_traits.rs b/crates/libm-test/src/test_traits.rs index c24ac6e4..deb83788 100644 --- a/crates/libm-test/src/test_traits.rs +++ b/crates/libm-test/src/test_traits.rs @@ -52,6 +52,8 @@ impl CheckCtx { pub enum CheckBasis { /// Check against Musl's math sources. Musl, + /// Check against infinite precision (MPFR). + Mpfr, } /// A trait to implement on any output type so we can verify it in a generic way. diff --git a/crates/libm-test/tests/compare_built_musl.rs b/crates/libm-test/tests/compare_built_musl.rs index 208b8e28..5a118f7c 100644 --- a/crates/libm-test/tests/compare_built_musl.rs +++ b/crates/libm-test/tests/compare_built_musl.rs @@ -29,8 +29,8 @@ macro_rules! musl_rand_tests { fn [< musl_random_ $fn_name >]() { let fname = stringify!($fn_name); let ulp = musl_allowed_ulp(fname); - let cases = random::get_test_cases::<$RustArgs>(fname); let ctx = CheckCtx::new(ulp, fname, CheckBasis::Musl); + let cases = random::get_test_cases::<$RustArgs>(&ctx); for input in cases { let musl_res = input.call(musl::$fn_name as $CFn); diff --git a/crates/libm-test/tests/multiprecision.rs b/crates/libm-test/tests/multiprecision.rs new file mode 100644 index 00000000..f8d94a16 --- /dev/null +++ b/crates/libm-test/tests/multiprecision.rs @@ -0,0 +1,71 @@ +//! Test with "infinite precision" + +#![cfg(feature = "test-multiprecision")] + +use libm_test::gen::random; +use libm_test::mpfloat::{self, MpOp}; +use libm_test::{CheckBasis, CheckCtx, CheckOutput, TupleCall, multiprec_allowed_ulp}; + +/// Implement a test against MPFR with random inputs. +macro_rules! multiprec_rand_tests { + ( + fn_name: $fn_name:ident, + CFn: $CFn:ty, + CArgs: $CArgs:ty, + CRet: $CRet:ty, + RustFn: $RustFn:ty, + RustArgs: $RustArgs:ty, + RustRet: $RustRet:ty, + attrs: [$($meta:meta)*] + ) => { + paste::paste! { + #[test] + $(#[$meta])* + fn [< multiprec_random_ $fn_name >]() { + type MpOpTy = mpfloat::$fn_name::Operation; + + let fname = stringify!($fn_name); + let ulp = multiprec_allowed_ulp(fname); + let mut mp_vals = MpOpTy::new(); + let ctx = CheckCtx::new(ulp, fname, CheckBasis::Mpfr); + let cases = random::get_test_cases::<$RustArgs>(&ctx); + + for input in cases { + let mp_res = mp_vals.run(input); + let crate_res = input.call(libm::$fn_name as $RustFn); + + crate_res.validate(mp_res, input, &ctx).unwrap(); + } + } + } + }; +} + +libm_macros::for_each_function! { + callback: multiprec_rand_tests, + attributes: [ + // Also an assertion failure on i686: at `MPFR_ASSERTN (! mpfr_erangeflag_p ())` + #[ignore = "large values are infeasible in MPFR"] + [jn, jnf], + ], + skip: [ + // FIXME: MPFR tests needed + frexp, + frexpf, + ilogb, + ilogbf, + ldexp, + ldexpf, + modf, + modff, + remquo, + remquof, + scalbn, + scalbnf, + + // FIXME: test needed, see + // https://github.com/rust-lang/libm/pull/311#discussion_r1818273392 + nextafter, + nextafterf, + ], +}