user0/rust
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

fix mangitude of applied float error

Browse source
This commit is contained in:
Ralf Jung 2025-09-02 17:17:19 +02:00
parent d4f861ecb6
commit 9f0b2a2df4
6 changed files with 103 additions and 69 deletions
Download patch file Download diff file Expand all files Collapse all files

49
src/tools/miri/src/intrinsics/mod.rs
View file

@ -10,7 +10,7 @@ use rustc_abi::Size;
use rustc_apfloat::ieee::{IeeeFloat, Semantics};
use rustc_apfloat::{self, Float, Round};
use rustc_middle::mir;
use rustc_middle::ty::{self, FloatTy, ScalarInt};
use rustc_middle::ty::{self, FloatTy};
use rustc_span::{Symbol, sym};
use self::atomic::EvalContextExt as _;
@ -230,7 +230,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
let res = apply_random_float_error_ulp(
this,
res,
2, // log2(4)
4,
);
// Clamp the result to the guaranteed range of this function according to the C standard,
@ -274,7 +274,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
let res = apply_random_float_error_ulp(
this,
res,
2, // log2(4)
4,
);
// Clamp the result to the guaranteed range of this function according to the C standard,
@ -336,9 +336,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
// Apply a relative error of 4ULP to introduce some non-determinism
// simulating imprecise implementations and optimizations.
apply_random_float_error_ulp(
this, res, 2, // log2(4)
)
apply_random_float_error_ulp(this, res, 4)
});
let res = this.adjust_nan(res, &[f1, f2]);
this.write_scalar(res, dest)?;
@ -354,9 +352,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
// Apply a relative error of 4ULP to introduce some non-determinism
// simulating imprecise implementations and optimizations.
apply_random_float_error_ulp(
this, res, 2, // log2(4)
)
apply_random_float_error_ulp(this, res, 4)
});
let res = this.adjust_nan(res, &[f1, f2]);
this.write_scalar(res, dest)?;
@ -373,9 +369,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
// Apply a relative error of 4ULP to introduce some non-determinism
// simulating imprecise implementations and optimizations.
apply_random_float_error_ulp(
this, res, 2, // log2(4)
)
apply_random_float_error_ulp(this, res, 4)
});
let res = this.adjust_nan(res, &[f]);
this.write_scalar(res, dest)?;
@ -391,9 +385,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
// Apply a relative error of 4ULP to introduce some non-determinism
// simulating imprecise implementations and optimizations.
apply_random_float_error_ulp(
this, res, 2, // log2(4)
)
apply_random_float_error_ulp(this, res, 4)
});
let res = this.adjust_nan(res, &[f]);
this.write_scalar(res, dest)?;
@ -448,7 +440,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
}
// Apply a relative error of 4ULP to simulate non-deterministic precision loss
// due to optimizations.
let res = apply_random_float_error_to_imm(this, res, 2 /* log2(4) */)?;
let res = crate::math::apply_random_float_error_to_imm(this, res, 4)?;
this.write_immediate(*res, dest)?;
}
@ -486,31 +478,6 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
}
}
/// Applies a random ULP floating point error to `val` and returns the new value.
/// So if you want an X ULP error, `ulp_exponent` should be log2(X).
///
/// Will fail if `val` is not a floating point number.
fn apply_random_float_error_to_imm<'tcx>(
ecx: &mut MiriInterpCx<'tcx>,
val: ImmTy<'tcx>,
ulp_exponent: u32,
) -> InterpResult<'tcx, ImmTy<'tcx>> {
let scalar = val.to_scalar_int()?;
let res: ScalarInt = match val.layout.ty.kind() {
ty::Float(FloatTy::F16) =>
apply_random_float_error_ulp(ecx, scalar.to_f16(), ulp_exponent).into(),
ty::Float(FloatTy::F32) =>
apply_random_float_error_ulp(ecx, scalar.to_f32(), ulp_exponent).into(),
ty::Float(FloatTy::F64) =>
apply_random_float_error_ulp(ecx, scalar.to_f64(), ulp_exponent).into(),
ty::Float(FloatTy::F128) =>
apply_random_float_error_ulp(ecx, scalar.to_f128(), ulp_exponent).into(),
_ => bug!("intrinsic called with non-float input type"),
};
interp_ok(ImmTy::from_scalar_int(res, val.layout))
}
/// For the intrinsics:
/// - sinf32, sinf64
/// - cosf32, cosf64

2
src/tools/miri/src/machine.rs
View file

@ -1278,7 +1278,7 @@ impl<'tcx> Machine<'tcx> for MiriMachine<'tcx> {
ecx: &mut InterpCx<'tcx, Self>,
val: ImmTy<'tcx>,
) -> InterpResult<'tcx, ImmTy<'tcx>> {
crate::math::apply_random_float_error_to_imm(ecx, val, 2 /* log2(4) */)
crate::math::apply_random_float_error_to_imm(ecx, val, 4)
}
#[inline(always)]

67
src/tools/miri/src/math.rs
View file

@ -6,10 +6,6 @@ use rustc_middle::ty::{self, FloatTy, ScalarInt};
use crate::*;
/// Disturbes a floating-point result by a relative error in the range (-2^scale, 2^scale).
///
/// For a 2^N ULP error, you can use an `err_scale` of `-(F::PRECISION - 1 - N)`.
/// In other words, a 1 ULP (absolute) error is the same as a `2^-(F::PRECISION-1)` relative error.
/// (Subtracting 1 compensates for the integer bit.)
pub(crate) fn apply_random_float_error<F: rustc_apfloat::Float>(
ecx: &mut crate::MiriInterpCx<'_>,
val: F,
@ -17,11 +13,13 @@ pub(crate) fn apply_random_float_error<F: rustc_apfloat::Float>(
) -> F {
if !ecx.machine.float_nondet
|| matches!(ecx.machine.float_rounding_error, FloatRoundingErrorMode::None)
// relative errors don't do anything to zeros... avoid messing up the sign
|| val.is_zero()
{
return val;
}
let rng = ecx.machine.rng.get_mut();
// Generate a random integer in the range [0, 2^PREC).
// (When read as binary, the position of the first `1` determines the exponent,
// and the remaining bits fill the mantissa. `PREC` is one plus the size of the mantissa,
@ -37,43 +35,66 @@ pub(crate) fn apply_random_float_error<F: rustc_apfloat::Float>(
let err = r.scalbn(err_scale.strict_sub(F::PRECISION.try_into().unwrap()));
// give it a random sign
let err = if rng.random() { -err } else { err };
// multiple the value with (1+err)
(val * (F::from_u128(1).value + err).value).value
// Compute `val*(1+err)`, distributed out as `val + val*err` to avoid the imprecise addition
// error being amplified by multiplication.
(val + (val * err).value).value
}
/// [`apply_random_float_error`] gives instructions to apply a 2^N ULP error.
/// This function implements these instructions such that applying a 2^N ULP error is less error prone.
/// So for a 2^N ULP error, you would pass N as the `ulp_exponent` argument.
/// Applies an error of `[-N, +N]` ULP to the given value.
pub(crate) fn apply_random_float_error_ulp<F: rustc_apfloat::Float>(
ecx: &mut crate::MiriInterpCx<'_>,
val: F,
ulp_exponent: u32,
max_error: u32,
) -> F {
let n = i32::try_from(ulp_exponent)
.expect("`err_scale_for_ulp`: exponent is too large to create an error scale");
// we know this fits
let prec = i32::try_from(F::PRECISION).unwrap();
let err_scale = -(prec - n - 1);
apply_random_float_error(ecx, val, err_scale)
// We could try to be clever and reuse `apply_random_float_error`, but that is hard to get right
// (see <https://github.com/rust-lang/miri/pull/4558#discussion_r2316838085> for why) so we
// implement the logic directly instead.
if !ecx.machine.float_nondet
|| matches!(ecx.machine.float_rounding_error, FloatRoundingErrorMode::None)
// FIXME: also disturb zeros? That requires a lot more cases in `fixed_float_value`
// and might make the std test suite quite unhappy.
|| val.is_zero()
{
return val;
}
let rng = ecx.machine.rng.get_mut();
let max_error = i64::from(max_error);
let error = match ecx.machine.float_rounding_error {
FloatRoundingErrorMode::Random => rng.random_range(-max_error..=max_error),
FloatRoundingErrorMode::Max =>
if rng.random() {
max_error
} else {
-max_error
},
FloatRoundingErrorMode::None => unreachable!(),
};
// If upwards ULP and downwards ULP differ, we take the average.
let ulp = (((val.next_up().value - val).value + (val - val.next_down().value).value).value
/ F::from_u128(2).value)
.value;
// Shift the value by N times the ULP
(val + (ulp * F::from_i128(error.into()).value).value).value
}
/// Applies a random 16ULP floating point error to `val` and returns the new value.
/// Applies an error of `[-N, +N]` ULP to the given value.
/// Will fail if `val` is not a floating point number.
pub(crate) fn apply_random_float_error_to_imm<'tcx>(
ecx: &mut MiriInterpCx<'tcx>,
val: ImmTy<'tcx>,
ulp_exponent: u32,
max_error: u32,
) -> InterpResult<'tcx, ImmTy<'tcx>> {
let scalar = val.to_scalar_int()?;
let res: ScalarInt = match val.layout.ty.kind() {
ty::Float(FloatTy::F16) =>
apply_random_float_error_ulp(ecx, scalar.to_f16(), ulp_exponent).into(),
apply_random_float_error_ulp(ecx, scalar.to_f16(), max_error).into(),
ty::Float(FloatTy::F32) =>
apply_random_float_error_ulp(ecx, scalar.to_f32(), ulp_exponent).into(),
apply_random_float_error_ulp(ecx, scalar.to_f32(), max_error).into(),
ty::Float(FloatTy::F64) =>
apply_random_float_error_ulp(ecx, scalar.to_f64(), ulp_exponent).into(),
apply_random_float_error_ulp(ecx, scalar.to_f64(), max_error).into(),
ty::Float(FloatTy::F128) =>
apply_random_float_error_ulp(ecx, scalar.to_f128(), ulp_exponent).into(),
apply_random_float_error_ulp(ecx, scalar.to_f128(), max_error).into(),
_ => bug!("intrinsic called with non-float input type"),
};

7
src/tools/miri/tests/pass/float.rs
View file

@ -39,9 +39,8 @@ macro_rules! assert_approx_eq {
}};
($a:expr, $b: expr) => {
// accept up to 12ULP (4ULP for host floats and 4ULP for miri artificial error and 4 for any additional effects
// due to having multiple error sources.
assert_approx_eq!($a, $b, 12);
// accept up to 8ULP (4ULP for host floats and 4ULP for miri artificial error).
assert_approx_eq!($a, $b, 8);
};
}
@ -176,6 +175,7 @@ fn assert_eq_msg<T: PartialEq + Debug>(x: T, y: T, msg: impl Display) {
}
/// Check that floats have bitwise equality
#[track_caller]
fn assert_biteq<F: Float>(a: F, b: F, msg: impl Display) {
let ab = a.to_bits();
let bb = b.to_bits();
@ -189,6 +189,7 @@ fn assert_biteq<F: Float>(a: F, b: F, msg: impl Display) {
}
/// Check that two floats have equality
#[track_caller]
fn assert_feq<F: Float>(a: F, b: F, msg: impl Display) {
let ab = a.to_bits();
let bb = b.to_bits();

10
src/tools/miri/tests/pass/float_extra_rounding_error.rs
View file

@ -9,7 +9,7 @@ use std::hint::black_box;
fn main() {
let expected = cfg_select! {
random => 13, // FIXME: why is it 13?
random => 9, // -4 ..= +4 ULP error
max => 2,
none => 1,
};
@ -20,4 +20,12 @@ fn main() {
values.insert(val.to_bits());
}
assert_eq!(values.len(), expected);
if !cfg!(none) {
// Ensure the smallest and biggest value are 8 ULP apart.
// We can just subtract the raw bit representations for this.
let min = *values.iter().min().unwrap();
let max = *values.iter().max().unwrap();
assert_eq!(min.abs_diff(max), 8);
}
}

37
src/tools/miri/tests/pass/shims/x86/rounding-error.rs Normal file
View file

@ -0,0 +1,37 @@
// We're testing x86 target specific features
//@only-target: x86_64 i686
//! rsqrt and rcp SSE/AVX operations are approximate. We use that as license to treat them as
//! non-deterministic. Ensure that we do indeed see random results within the expected error bounds.
#[cfg(target_arch = "x86")]
use std::arch::x86::*;
#[cfg(target_arch = "x86_64")]
use std::arch::x86_64::*;
use std::collections::HashSet;
fn main() {
let mut vals = HashSet::new();
for _ in 0..50 {
unsafe {
// Compute the inverse square root of 4.0, four times.
let a = _mm_setr_ps(4.0, 4.0, 4.0, 4.0);
let exact = 0.5;
let r = _mm_rsqrt_ps(a);
let r: [f32; 4] = std::mem::transmute(r);
// Check the results.
for r in r {
vals.insert(r.to_bits());
// Ensure the relative error is less than 2^-12.
let rel_error = (r - exact) / exact;
let log_error = rel_error.abs().log2();
assert!(
rel_error == 0.0 || log_error < -12.0,
"got an error of {rel_error} = 2^{log_error}"
);
}
}
}
// Ensure we saw a bunch of different results.
assert!(vals.len() >= 50);
}