contrib/math: preserve negative zero in Expm1

[EylonKrause]

Summary

Expm1(-0.0) returns +0.0 instead of -0.0, losing the sign of zero. C11 Annex F (F.10.3.3) requires expm1(±0) = ±0, and scalar std::expm1 complies — so this is a signed-zero correctness bug (observable via std::signbit, or 1.0 / expm1(-0.0) giving +inf instead of -inf).

Cause

On the |x| < ln2/2 polynomial path, ExpPoly computes MulAdd(poly, x*x, x). For x = -0.0 this is 0.5 * (+0.0) + (-0.0), and IEEE addition of oppositely-signed zeros yields +0.0, so the sign bit is dropped. Unlike the odd functions in this file (Sin, Tan, Sinh, …), Expm1 never reapplies the sign of x, so the negative zero is not restored.

Fix

Reapply x's sign on the polynomial path using the sign-bit mask already defined in the function:

const V z = IfThenElse(Lt(Abs(x), kLn2Over2), Or(y, And(x, kNegZero)),
                       Sub(impl.LoadExpShortRange(d, Add(y, kOne), q), kOne));

expm1(x) shares the sign of x for all x, so OR-ing in x's sign bit is a no-op for every input except x == -0.0, which it corrects to -0.0. +0.0, small nonzero x, the large-|x| branch, and the -1 clamp are unchanged, on every SIMD target.

Testing

• Added TestExpm1SignedZero, a bit-exact check that Expm1(±0) returns ±0 — the existing ULP comparator treats +0 == -0 and cannot detect this (which is why it went unnoticed).
• Before the fix the new test fails (Expm1(-0) yields +0); after, it passes.
• Full math_test passes (104 tests) across every enabled SIMD target.

Disclosure: this contribution was authored with an AI coding assistant (Claude) and reviewed before submission.

[jan-wassenberg]

Thanks again for fixing! Rather than Or+And, please use CopySign(y, x) - this is more efficient on some platforms.

[jan-wassenberg]

Once again thanks for fixing! Rather than Or+And, please use CopySign(y, x) - this is more efficient on some platforms.