Fix incorrect rounding for subnormals in widening FMA

[Shnatsel]

Adopts the formally proven algorithm from the 2008 paper "Emulation of FMA and correctly-rounded sums:
proved algorithms using rounding to odd" by Sylvie Boldo and Guillaume Melquiond, https://guillaume.melquiond.fr/doc/08-tc.pdf

Fixes #1262

The common path is unchanged, only the handling of subnormals is altered.

I've added regression tests in-tree, and additionally validated the implementation via random testing on subnormals.

[Shnatsel]

The added unit test involving f128 on 32-bit ARM failed due to ABI bugs: my_f64 as f128 as well as .narrow() produce garbage values by reading from the wrong reigster. I've worked around that by disabling the f128 unit test on ARM, but 32-bit ARM probably shouldn't declare f128_enabled at all.

The instruction count on benchmarks also regresses, but that's expected given that we're doing more work for to ensure correctness. The bulk of it is under an if that only fires for subnormals, which are slow even in hardware so it shouldn't be a big deal.

[tgross35]

Thank you for the issue and the fix, this looks good at first pass but I need to do some number crunching.

To confirm: #1262 mentions x86-64 specifically, but this is present on all arches without hardware fmaf right?

The added unit test involving f128 on 32-bit ARM failed due to ABI bugs: my_f64 as f128 as well as .narrow() produce garbage values by reading from the wrong reigster. I've worked around that by disabling the f128 unit test on ARM, but 32-bit ARM probably shouldn't declare f128_enabled at all.

Happen to have a repro? I'm not aware of any ABI bugs on this platform, everything we know of is tracked at https://github.com/rust-lang/rust/blob/ef20314466010b8b9259ec5f86230c530ca08661/compiler/rustc_codegen_llvm/src/llvm_util.rs#L392-L412.

[Shnatsel]

Thank you! I'd appreciate a thorough review. I fit the algorithm from the paper into the existing structure with raising floating-point flags as best I could, but I'm not 100% confident I got it right.

To confirm: #1262 mentions x86-64 specifically, but this is present on all arches without hardware fmaf right?

Yes. I've double-checked and 32-bit ARM on qemu fails the f32 rounding test too. I mentioned x86_64 to avoid confusion around non-conformant i586 (x87) floats.

Happen to have a repro?

Yes, check out commit 8b351a6 in this branch. The logs from CI runs show the failures clearly. It's not just rounding failures, the result is just garbage.

[tgross35]

I meant would you be able to reproduce the ABI issue with a smaller snip of code? I wrote a quick test and annotated what I think is going on at https://rust.godbolt.org/z/9PjYETqf6. It's a little weird that Rust and C have different return ABIs, but everything seems otherwise correct there.

[Shnatsel]

Oh! Yes, this fails in qemu for me on 32-bit ARM, both targets that failed in CI, while passing on x86_64:

#![feature(f128)]

fn main() {
    let x = std::hint::black_box(1.0_f64);
    let y = std::hint::black_box(x as f128);
    let bits = y.to_bits();

    assert_eq!(
        bits,
        0x3fff_0000_0000_0000_0000_0000_0000_0000,
        "got {bits:#034x}",
    );
}

[tgross35]

Whoops, I missed that it's the hf targets that fail, but the codegen still looks fine https://rust.godbolt.org/z/hPd196xbE. Digging...

[Shnatsel]

I don't know enough ARM assembly to help here, but here's what my friendly neighbourhood LLM thinks of this:

LLM analysis

I ran it with:

rustc +nightly reproducers/arm-f128-hardfloat-abi.rs \
  --target armv7-unknown-linux-gnueabihf \
  -C linker=arm-linux-gnueabihf-gcc \
  -o /tmp/arm-f128-abi

qemu-arm -L /usr/arm-linux-gnueabihf /tmp/arm-f128-abi

It passes on x86_64 but fails on ARM hard-float. The incorrect bits can vary because the callee reads the wrong registers.

To expose the problematic call:

rustc +nightly reproducers/arm-f128-hardfloat-abi.rs \
  --target armv7-unknown-linux-gnueabihf \
  -C linker=arm-linux-gnueabihf-gcc \
  -O --emit=asm \
  -o /tmp/arm-f128-abi.s

rg -C 8 __extenddftf2 /tmp/arm-f128-abi.s

The assembly shows rustc placing the f64 in d0:

vldr d0, [...]
bl   __extenddftf2

But compiler-builtins defines __extenddftf2 as extern "aapcs", which reads the argument from r0/r1.

Verified with:

rustc 1.99.0-nightly (84b36a78a 2026-08-06)
LLVM 23.1.0

Please take it with a grain of salt as always.

[tgross35]

I think that's pretty accurate, wrote up #1271 to follow up more

[Shnatsel]

I've spent some time minifying the diff to the algorithm because my previous changes messed up the exception flag handling, e.g. on these values:

a = f32::from_bits(0x5f78_0000);
b = f32::from_bits(0x5f84_2108);
c = f32::from_bits(0x8000_0001);

In the first iteration of this PR that raised the overflow flag because of the early narrowing.

I believe the state as of e254558 deviates from the original algorithm as little as possible: it hoists one variable out of a branch (which gets const-folded anyway), adds a check for subnormals and adds one additional check required for correctness.

It took me a while to fully grasp the algorithm and realize this simple fix is possible 😅