Auto merge of #142514 - LorrensP-2158466:miri-float-nondet-pow, r=<try>

Miri: handling of SNaN inputs in `f*::pow` operations

<!-- homu-ignore:start -->
<!--
If this PR is related to an unstable feature or an otherwise tracked effort,
please link to the relevant tracking issue here. If you don't know of a related
tracking issue or there are none, feel free to ignore this.

This PR will get automatically assigned to a reviewer. In case you would like
a specific user to review your work, you can assign it to them by using

    r? <reviewer name>
-->
<!-- homu-ignore:end -->

fixes [miri/#4286](rust-lang/miri#4286) and related to #138062 and [miri/#4208](rust-lang/miri#4208 (comment)).

For the following cases of the powf or powi operations, Miri returns either `1.0` or an arbitrary `NaN`:
- `powf(SNaN, 0.0)`
- `powf(1.0, SNaN)`
- `powi(SNaN, 0)`

Also added a macro in `miri/tests/pass/float.rs` which conveniently checks if both are indeed returned from such an operation.

Made these changes in the rust repo so I could test against stdlib, since these were impacted some time ago and were fixed in #138062. Tested with:
```fish
env MIRIFLAGS=-Zmiri-many-seeds ./x miri --no-fail-fast std core coretests -- f32 f64
```
This was successful. This does take a while, so I recommend using `--no-doc` and separate use of `f32` or `f64`

The pr is somewhat split up into 3 main commits, which implement the cases described above. The first commit also introduces the macro, and the last commit is just a global refactor of some things.

r? `@RalfJung`
try-job: x86_64-gnu-aux

Author: rust-bors[bot]

src/tools/miri/src/intrinsics/mod.rs

@@ -6,6 +6,7 @@ mod simd;
 use std::ops::Neg;
 
 use rand::Rng;
+use rand::rngs::StdRng;
 use rustc_abi::Size;
 use rustc_apfloat::ieee::{IeeeFloat, Semantics};
 use rustc_apfloat::{self, Float, Round};
@@ -191,7 +192,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
                 let [f] = check_intrinsic_arg_count(args)?;
                 let f = this.read_scalar(f)?.to_f32()?;
 
-                let res = fixed_float_value(intrinsic_name, &[f]).unwrap_or_else(||{
+                let res = fixed_float_value(this, intrinsic_name, &[f]).unwrap_or_else(||{
                     // Using host floats (but it's fine, these operations do not have
                     // guaranteed precision).
                     let host = f.to_host();
@@ -235,7 +236,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
                 let [f] = check_intrinsic_arg_count(args)?;
                 let f = this.read_scalar(f)?.to_f64()?;
 
-                let res = fixed_float_value(intrinsic_name, &[f]).unwrap_or_else(||{
+                let res = fixed_float_value(this, intrinsic_name, &[f]).unwrap_or_else(||{
                     // Using host floats (but it's fine, these operations do not have
                     // guaranteed precision).
                     let host = f.to_host();
@@ -312,7 +313,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
                 let f1 = this.read_scalar(f1)?.to_f32()?;
                 let f2 = this.read_scalar(f2)?.to_f32()?;
 
-                let res = fixed_float_value(intrinsic_name, &[f1, f2]).unwrap_or_else(|| {
+                let res = fixed_float_value(this, intrinsic_name, &[f1, f2]).unwrap_or_else(|| {
                     // Using host floats (but it's fine, this operation does not have guaranteed precision).
                     let res = f1.to_host().powf(f2.to_host()).to_soft();
 
@@ -330,7 +331,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
                 let f1 = this.read_scalar(f1)?.to_f64()?;
                 let f2 = this.read_scalar(f2)?.to_f64()?;
 
-                let res = fixed_float_value(intrinsic_name, &[f1, f2]).unwrap_or_else(|| {
+                let res = fixed_float_value(this, intrinsic_name, &[f1, f2]).unwrap_or_else(|| {
                     // Using host floats (but it's fine, this operation does not have guaranteed precision).
                     let res = f1.to_host().powf(f2.to_host()).to_soft();
 
@@ -349,7 +350,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
                 let f = this.read_scalar(f)?.to_f32()?;
                 let i = this.read_scalar(i)?.to_i32()?;
 
-                let res = fixed_powi_float_value(f, i).unwrap_or_else(|| {
+                let res = fixed_powi_float_value(this, f, i).unwrap_or_else(|| {
                     // Using host floats (but it's fine, this operation does not have guaranteed precision).
                     let res = f.to_host().powi(i).to_soft();
 
@@ -367,7 +368,7 @@ pub trait EvalContextExt<'tcx>: crate::MiriInterpCxExt<'tcx> {
                 let f = this.read_scalar(f)?.to_f64()?;
                 let i = this.read_scalar(i)?.to_i32()?;
 
-                let res = fixed_powi_float_value(f, i).unwrap_or_else(|| {
+                let res = fixed_powi_float_value(this, f, i).unwrap_or_else(|| {
                     // Using host floats (but it's fine, this operation does not have guaranteed precision).
                     let res = f.to_host().powi(i).to_soft();
 
@@ -489,56 +490,97 @@ fn apply_random_float_error_to_imm<'tcx>(
     interp_ok(ImmTy::from_scalar_int(res, val.layout))
 }
 
+/// Returns either a SNaN or a QNaN, with a randomly generated payload.
+fn random_nan<S: Semantics>(rng: &mut StdRng) -> IeeeFloat<S> {
+    if rng.random() {
+        IeeeFloat::<S>::snan(Some(rng.random()))
+    } else {
+        IeeeFloat::<S>::qnan(Some(rng.random()))
+    }
+}
+
 /// For the intrinsics:
 /// - sinf32, sinf64
 /// - cosf32, cosf64
 /// - expf32, expf64, exp2f32, exp2f64
 /// - logf32, logf64, log2f32, log2f64, log10f32, log10f64
 /// - powf32, powf64
 ///
+/// # Return
+///
 /// Returns `Some(output)` if the `intrinsic` results in a defined fixed `output` specified in the C standard
 /// (specifically, C23 annex F.10)  when given `args` as arguments. Outputs that are unaffected by a relative error
 /// (such as INF and zero) are not handled here, they are assumed to be handled by the underlying
 /// implementation. Returns `None` if no specific value is guaranteed.
+///
+/// # Note
+///
+/// For `powf*` operations of the form:
+///
+/// - `(SNaN)^(±0)`
+/// - `1^(SNaN)`
+///
+/// The result is implementation-defined:
+/// - musl returns for both `1.0`
+/// - glibc returns for both `NaN`
+///
+/// This discrepancy exists because SNaN handling is not consistently defined across platforms,
+/// and the C standard leaves behavior for SNaNs unspecified.
+///
+/// Miri chooses to adhere to both implementations and returns either one of them non-deterministically.
 fn fixed_float_value<S: Semantics>(
+    ecx: &mut MiriInterpCx<'_>,
     intrinsic_name: &str,
     args: &[IeeeFloat<S>],
 ) -> Option<IeeeFloat<S>> {
     let one = IeeeFloat::<S>::one();
-    match (intrinsic_name, args) {
+    Some(match (intrinsic_name, args) {
         // cos(+- 0) = 1
-        ("cosf32" | "cosf64", [input]) if input.is_zero() => Some(one),
+        ("cosf32" | "cosf64", [input]) if input.is_zero() => one,
 
         // e^0 = 1
-        ("expf32" | "expf64" | "exp2f32" | "exp2f64", [input]) if input.is_zero() => Some(one),
-
-        // 1^y = 1 for any y, even a NaN.
-        ("powf32" | "powf64", [base, _]) if *base == one => Some(one),
+        ("expf32" | "expf64" | "exp2f32" | "exp2f64", [input]) if input.is_zero() => one,
 
         // (-1)^(±INF) = 1
-        ("powf32" | "powf64", [base, exp]) if *base == -one && exp.is_infinite() => Some(one),
+        ("powf32" | "powf64", [base, exp]) if *base == -one && exp.is_infinite() => one,
 
-        // FIXME(#4286): The C ecosystem is inconsistent with handling sNaN's, some return 1 others propogate
-        // the NaN. We should return either 1 or the NaN non-deterministically here.
-        // But for now, just handle them all the same.
-        // x^(±0) = 1 for any x, even a NaN
-        ("powf32" | "powf64", [_, exp]) if exp.is_zero() => Some(one),
+        // 1^y = 1 for any y, even a NaN, *but* not a SNaN
+        ("powf32" | "powf64", [base, exp]) if *base == one => {
+            let rng = ecx.machine.rng.get_mut();
+            // Handle both the musl and glibc cases non-deterministically.
+            if !exp.is_signaling() || rng.random() { one } else { random_nan(rng) }
+        }
+
+        // x^(±0) = 1 for any x, even a NaN, *but* not a SNaN
+        ("powf32" | "powf64", [base, exp]) if exp.is_zero() => {
+            let rng = ecx.machine.rng.get_mut();
+            // Handle both the musl and glibc cases non-deterministically.
+            if !base.is_signaling() || rng.random() { one } else { random_nan(rng) }
+        }
 
         // There are a lot of cases for fixed outputs according to the C Standard, but these are mainly INF or zero
         // which are not affected by the applied error.
-        _ => None,
-    }
+        _ => return None,
+    })
 }
 
 /// Returns `Some(output)` if `powi` (called `pown` in C) results in a fixed value specified in the C standard
 /// (specifically, C23 annex F.10.4.6) when doing `base^exp`. Otherwise, returns `None`.
-fn fixed_powi_float_value<S: Semantics>(base: IeeeFloat<S>, exp: i32) -> Option<IeeeFloat<S>> {
-    match (base.category(), exp) {
-        // x^0 = 1, if x is not a Signaling NaN
-        // FIXME(#4286): The C ecosystem is inconsistent with handling sNaN's, some return 1 others propogate
-        // the NaN. We should return either 1 or the NaN non-deterministically here.
-        // But for now, just handle them all the same.
-        (_, 0) => Some(IeeeFloat::<S>::one()),
+// TODO: I'm not sure what I should document here about pown(1, SNaN) since musl and glibc do the same and the C standard is explicit here.
+fn fixed_powi_float_value<S: Semantics>(
+    ecx: &mut MiriInterpCx<'_>,
+    base: IeeeFloat<S>,
+    exp: i32,
+) -> Option<IeeeFloat<S>> {
+    match exp {
+        0 => {
+            let one = IeeeFloat::<S>::one();
+            let rng = ecx.machine.rng.get_mut();
+            Some(
+                // Handle both the musl and glibc powf cases non-deterministically.
+                if !base.is_signaling() || rng.random() { one } else { random_nan(rng) },
+            )
+        }
 
         _ => None,
     }

src/tools/miri/tests/pass/float.rs

@@ -1066,17 +1066,39 @@ pub fn libm() {
     assert_eq!((-1f32).powf(f32::NEG_INFINITY), 1.0);
     assert_eq!((-1f64).powf(f64::NEG_INFINITY), 1.0);
 
-    // For pow (powf in rust) the C standard says:
-    // x^0 = 1 for all x even a sNaN
-    // FIXME(#4286): this does not match the behavior of all implementations.
-    assert_eq!(SNAN_F32.powf(0.0), 1.0);
-    assert_eq!(SNAN_F64.powf(0.0), 1.0);
-
-    // For pown (powi in rust) the C standard says:
-    // x^0 = 1 for all x even a sNaN
-    // FIXME(#4286): this does not match the behavior of all implementations.
-    assert_eq!(SNAN_F32.powi(0), 1.0);
-    assert_eq!(SNAN_F64.powi(0), 1.0);
+    // Makes sure an operations returns both `1` and a `NaN` randomly.
+    macro_rules! test_snan_nondet {
+        ($pow_op:expr) => {{
+            let mut nan_seen = false;
+            let mut one_seen = false;
+
+            for _ in 0..64 {
+                let res = $pow_op;
+                nan_seen |= res.is_nan();
+                one_seen |= res == 1.0;
+
+                // little speedup
+                if nan_seen && one_seen { break; };
+            }
+
+            let op_as_str = stringify!($pow_op);
+
+            assert!(nan_seen && one_seen, "{} should return both `NaN` or `1.0` randomly", op_as_str);
+        }};
+    }
+
+    // x^(SNaN) = (1 | NaN)
+    test_snan_nondet!(f32::powf(SNAN_F32, 0.0));
+    test_snan_nondet!(f64::powf(SNAN_F64, 0.0));
+
+    // 1^(SNaN) = (1 | NaN)
+    test_snan_nondet!(f32::powf(1.0, SNAN_F32));
+    test_snan_nondet!(f64::powf(1.0, SNAN_F64));
+
+    // same as powf (keep it consistent):
+    // x^(SNaN) = (1 | NaN)
+    test_snan_nondet!(f32::powi(SNAN_F32, 0));
+    test_snan_nondet!(f64::powi(SNAN_F64, 0));
 
     assert_eq!(0f32.powi(10), 0.0);
     assert_eq!(0f64.powi(100), 0.0);