Is fallback `max`, `precise_max`, `min`, `precise_min` implemented correctly?

[ajakubowicz-canva]

This is a tracking issue that may require no action. I'm just logging it to investigate.

Currently fallback is implemented as (just looking at min):

    #[inline(always)]
    fn min_f32x4(self, a: f32x4<Self>, b: f32x4<Self>) -> f32x4<Self> {
        [
            f32::min(a[0usize], b[0usize]),
            f32::min(a[1usize], b[1usize]),
            f32::min(a[2usize], b[2usize]),
            f32::min(a[3usize], b[3usize]),
        ]
        .simd_into(self)
    }
    #[inline(always)]
    fn min_precise_f32x4(self, a: f32x4<Self>, b: f32x4<Self>) -> f32x4<Self> {
        [
            f32::min(a[0usize], b[0usize]),
            f32::min(a[1usize], b[1usize]),
            f32::min(a[2usize], b[2usize]),
            f32::min(a[3usize], b[3usize]),
        ]
        .simd_into(self)
    }

And wasm is similarly only using f32x4_min for both (although wasm should use pmin for precise_min).

I think we may need to address this so there's a difference between max and max_precise. Otherwise, why have both?

[LaurenzV]

The min method of scalar corresponds to precise_min IIRC, but since there is no instruction for fast_min, I just used the same. On NEON, two different instructions are used here. But if it turns out that they are both equally fast, maybe it's worth to just rename precise_min to min and use that everywhere.

[ajakubowicz-canva]

@LaurenzV Thank you

The min method of scalar corresponds to precise_min IIRC, but since there is no instruction for fast_min, I just used the same.

This made it click for me. I've been using the fallback as a source of truth for the WASM implementation, so I accidentally used the non-precise min/max variants in WASM for both operations. Thank you for clarifying, I can work on a fix! It may be what's causing the bug in vello_cpu + wasm + simd.

[raphlinus]

I don't have the full matrix in front of me, but the other relevant info is the behavior of _mm_min_ps and friends on x86. The documentation is a bit vague, just saying it isn't IEEE. My guess is that we should map min to this and have a less-optimized sequence for precise_min (very likely compares and blends). In other words, it is x86 that justifies having two different methods, otherwise I would be in favor of unification.

In implementing the robustness fixes for eval_err, I found that minimum behavior would be useful. I believe that's vmin_f32 on Neon, while vminnm_f32 is the IEEE one matching Rust min, but this would need to be double-checked. For fallback, this function is a bit tricky, as it's not been stabilized, though I think it should be reasonably easy to polyfill. It's not urgent though, at the moment that logic does an is_finite check which is probably no more expensive than a polyfilled version; the win is probably on Neon where there is an instruction of the desired semantics (and I think equivalent performance, though again don't know for sure).

See also DLTcollab/sse2neon#57 for discussion of this issue. They're dealing with somewhat of the reverse problem that we are, they consider the x86 version to be the source of truth and are concerned with faithfully implementing it on Neon. But it looks like good documentation of x86 semantics and very likely a source of test vectors.

[ajakubowicz-canva]

Thank you all! Adding one last reference that helped me get WASM passing the vello_sparse_tests test suite with fearless_simd enabled.

https://github.com/google/highway/blob/9d1fa11ad921b97357aed0cb277ded70ad9170f8/hwy/ops/wasm_128-inl.h#L899-L904

There is operator order significance with pmin and pmax on WASM. I'm working on a PR for this fix now and will add tests to prevent us from regressions. Will close this issue with that PR.

[ajakubowicz-canva]

Wow - the accurate simulation of _mm_max_ps using neon is pretty interesting.

These are the float implementations on Highway – which when I follow result in the passing of the vello_cpu test suite (using WASM SIMD):

• Wasm: wasm_f32x4_pmax(b.raw, a.raw);
• Neon: vminnmq_f32(a.raw, b.raw)
• x86_128: _mm_max_ps(a.raw, b.raw)

Following Highway's semantics gets us the following behavior between `neon` and `WASM` (click to expand code):

fn main() {
    // Run with `cargo run` (my machine is aarch64)
    #[cfg(target_arch = "aarch64")]
    unsafe {
        use std::arch::aarch64::*;

        let a = &[f32::NAN, -3.0, f32::INFINITY, 0.5];
        let b = &[1.0, f32::NAN, 7.0, f32::NEG_INFINITY];

        let a_vec = vld1q_f32(a.as_ptr());
        let b_vec = vld1q_f32(b.as_ptr());

        println!("Vector result: {:?}", vmaxnmq_f32(a_vec, b_vec));
        // Vector result: float32x4_t(1.0, -3.0, inf, 0.5)
    }

    // Built with: `cargo build --target=wasm32-wasip2`
    // Run with: `wasmtime ./target/wasm32-wasip2/debug/tmp_neon_tests.wasm`
    #[cfg(target_arch = "wasm32")]
    unsafe {
        use std::{arch::wasm32::*, mem};

        let a = f32x4(f32::NAN, -3.0, f32::INFINITY, 0.5);
        let b = f32x4(1.0, f32::NAN, 7.0, f32::NEG_INFINITY);

        let result: [f32; 4] = mem::transmute(f32x4_pmax(b, a));
        println!("Vector result: {:?}", result);
        // Vector result: [1.0, NaN, inf, 0.5]
    }
}

I don't have access to x86 but allegedly it should match the WASM behavior (without flipping arguments).

[valadaptive]

This can probably be closed because of #136, but even the fallback Rust-based min/max operations are not immune to architecture-specific weirdness. They say they always return the non-NaN operand if possible, but if either input is a signaling NaN, they return a quiet NaN on AArch64. See the discussion starting from this comment. We should document this.

[valadaptive]

The sNaN vs qNaN behavior has been added to the now-generated documentation.