diff --git a/core/libs/prx/libc/src/Math.cpp b/core/libs/prx/libc/src/Math.cpp index d04aed2822..91e8efa251 100644 --- a/core/libs/prx/libc/src/Math.cpp +++ b/core/libs/prx/libc/src/Math.cpp @@ -1273,6 +1273,7 @@ float APS5_VABI cbrtf_nid_postfix(float x) { return static_cast(T); } float APS5_VABI remainderf_nid_postfix(float x, float y) { return std::remainder(x, y); } + float APS5_VABI nextafterf_nid_postfix(float x, float y) { const std::uint32_t ux = std::bit_cast(x); const std::uint32_t uy = std::bit_cast(y); @@ -1285,17 +1286,20 @@ float APS5_VABI nextafterf_nid_postfix(float x, float y) { return y; } if ((ux & 0x7fffffffu) == 0) { - return std::bit_cast((uy & 0x80000000u) | 1u); + const volatile float smallest = std::bit_cast((uy & 0x80000000u) | 1u); + volatile float product = smallest * smallest; + return product == smallest ? product : smallest; } - const bool towardZero = (x > 0.0f) == (x > y); + const bool towardZero = ((ux ^ uy) & 0x80000000u) != 0 || ux > uy; const std::uint32_t next = towardZero ? ux - 1 : ux + 1; const std::uint32_t exponent = next & 0x7f800000u; if (exponent == 0x7f800000u) { - errno = ERANGE; + _Feraise_nid_postfix(0x08); return x + x; } if (exponent == 0) { - errno = ERANGE; + volatile float product = x * x; + if (product != x) _Feraise_nid_postfix(0x10); } return std::bit_cast(next); } @@ -1311,17 +1315,20 @@ double APS5_VABI nextafter_nid_postfix(double x, double y) { return y; } if ((ux & 0x7fffffffffffffffull) == 0) { - return std::bit_cast((uy & 0x8000000000000000ull) | 1u); + const volatile double smallest = std::bit_cast((uy & 0x8000000000000000ull) | 1u); + volatile double product = smallest * smallest; + return product == smallest ? product : smallest; } - const bool towardZero = (x > 0.0) == (x > y); + const bool towardZero = ((ux ^ uy) & 0x8000000000000000ull) != 0 || ux > uy; const std::uint64_t next = towardZero ? ux - 1 : ux + 1; const std::uint64_t exponent = next & 0x7ff0000000000000ull; if (exponent == 0x7ff0000000000000ull) { - errno = ERANGE; + _Feraise_nid_postfix(0x08); return x + x; } if (exponent == 0) { - errno = ERANGE; + volatile double product = x * x; + if (product != x) _Feraise_nid_postfix(0x10); } return std::bit_cast(next); } diff --git a/core/libs/tests/GuestMath.cpp b/core/libs/tests/GuestMath.cpp index f9cffe9b6e..ea2219eb63 100644 --- a/core/libs/tests/GuestMath.cpp +++ b/core/libs/tests/GuestMath.cpp @@ -671,6 +671,69 @@ static void CheckRoundingNaNConsoleBits() { } } +template +static void CheckNextafterEnvironment(T (APS5_VABI *next)(T, T), UInt sign, UInt one, UInt minimum, UInt maximum, UInt infinity) { + const RestoreFloatingPointState saved; + const auto savedControl = _mm_getcsr(); + const struct { + UInt x, y, normalResult, dazResult; + unsigned normalFlags, dazFlags; + int normalError, dazError; + } cases[] = { + {0, one, 1, 0, 0x32, 0, 97, 97}, + {sign, sign | one, sign | 1, 0, 0x32, 0, 97, 97}, + {1, 0, 0, 0, 0x32, 0, ERANGE, 97}, + {1, one, 2, 2, 0x32, 0, ERANGE, 97}, + {sign | 1, sign | one, sign | 2, sign | 2, 0x32, 0, ERANGE, 97}, + {minimum, 0, minimum - 1, minimum - 1, 0x30, 0x30, ERANGE, ERANGE}, + {maximum, infinity, infinity, infinity, 0x28, 0x28, ERANGE, ERANGE}, + {1, 1, 1, 1, 0x02, 0, 97, 97}, + {one, one + 1, one + 1, one + 1, 0, 0, 97, 97}, + }; + for (unsigned denorm = 0; denorm != 4; ++denorm) { + const unsigned control = 0x1f80u | ((denorm & 1u) ? 0x40u : 0u) | ((denorm & 2u) ? 0x8000u : 0u); + const bool daz = (denorm & 1u) != 0; + for (const auto& item : cases) { + Require(std::fesetround(FE_TONEAREST) == 0); + Require(std::feclearexcept(FE_ALL_EXCEPT) == 0); + errno = 97; + _mm_setcsr(control); + const auto result = std::bit_cast(next(std::bit_cast(item.x), std::bit_cast(item.y))); + const auto flags = _mm_getcsr() & 0x3fu; + const int exceptions = std::fetestexcept(FE_ALL_EXCEPT); + const int error = errno; + const auto expectedFlags = daz ? item.dazFlags : item.normalFlags; + Require(result == (daz ? item.dazResult : item.normalResult)); + Require(flags == expectedFlags); + const int expectedExceptions = ((expectedFlags & 0x01u) ? FE_INVALID : 0) | + ((expectedFlags & 0x04u) ? FE_DIVBYZERO : 0) | + ((expectedFlags & 0x08u) ? FE_OVERFLOW : 0) | + ((expectedFlags & 0x10u) ? FE_UNDERFLOW : 0) | + ((expectedFlags & 0x20u) ? FE_INEXACT : 0); + Require(exceptions == expectedExceptions); + Require(error == (daz ? item.dazError : item.normalError)); + } + } + constexpr int rounding[] = {FE_DOWNWARD, FE_UPWARD, FE_TOWARDZERO}; + for (const int mode : rounding) { + Require(std::fesetround(mode) == 0); + const unsigned control = (_mm_getcsr() & 0x6000u) | 0x1f80u; + errno = 97; + Require(std::feclearexcept(FE_ALL_EXCEPT) == 0); + _mm_setcsr(control); + Require(std::bit_cast(next(std::bit_cast(maximum), std::bit_cast(infinity))) == (mode == FE_UPWARD ? infinity : maximum)); + Require((_mm_getcsr() & 0x3fu) == 0x28u); + Require(errno == ERANGE); + errno = 97; + Require(std::feclearexcept(FE_ALL_EXCEPT) == 0); + _mm_setcsr(control); + Require(std::bit_cast(next(std::bit_cast(UInt{1}), std::bit_cast(one))) == 2); + Require((_mm_getcsr() & 0x3fu) == 0x32u); + Require(errno == (mode == FE_UPWARD ? 97 : ERANGE)); + } + _mm_setcsr(savedControl); +} + static void CheckAsinl() { const struct { std::uint64_t inputMantissa; std::uint16_t inputSignExponent; std::uint64_t outputMantissa; std::uint16_t outputSignExponent; } cases[] = { {0, 0x0000, 0, 0x0000}, {0, 0x8000, 0, 0x8000}, @@ -1854,6 +1917,8 @@ int main() { CheckNextafterConsoleBits(); CheckTanhConsoleBits(); CheckRoundingNaNConsoleBits(); + CheckNextafterEnvironment(nextafterf_nid_postfix, 0x80000000u, 0x3f800000u, 0x00800000u, 0x7f7fffffu, 0x7f800000u); + CheckNextafterEnvironment(nextafter_nid_postfix, std::uint64_t{0x8000000000000000ull}, std::uint64_t{0x3ff0000000000000ull}, std::uint64_t{0x0010000000000000ull}, std::uint64_t{0x7fefffffffffffffull}, std::uint64_t{0x7ff0000000000000ull}); CheckAsinl(); Require(scalbn_nid_postfix(0.75, -2) == 0.1875); Require(scalbnf_nid_postfix(0.75f, 4) == 12.f); diff --git a/docs/dev/debt/libc.md b/docs/dev/debt/libc.md index f0c366cb14..ae8a6c9018 100644 --- a/docs/dev/debt/libc.md +++ b/docs/dev/debt/libc.md @@ -118,7 +118,7 @@ - [coshf](../../../core/libs/prx/libc/src/Math.cpp) (libc) follows FreeBSD msun's [`e_coshf.c`](https://github.com/freebsd/freebsd-src/blob/5aa554c7e5b6fd9042986df07f12ac5687e41392/lib/msun/src/e_coshf.c), using the `expm1f` and `expf` paths already inlined for `sinhf` below; the Windows UCRT differed from the `libc.prx` copies in PPSA07923 and PPSA13934 on 500,554 and 500,614 of 20 million inputs in the #9707 survey. The exact system `libSceLibcInternal` behavior and the FMA-contracted `libc.prx` results remain unverified - [expm1](../../../core/libs/prx/libc/src/Math.cpp) (libc) does not compare floating-point exceptions, and its double inputs are checked by sampling, not exhaustively. It is FreeBSD msun's [`s_expm1.c`](https://github.com/freebsd/freebsd-src/blob/release/11.4.0/lib/msun/src/s_expm1.c) with `huge` read at run time, and sets `ERANGE` only when a finite argument overflows to infinity, as the titles' `libc.prx` does. It gives the same bits and `errno` as the `expm1` bundled with PPSA07923 and PPSA05916 on 320 million random and 47,222 structured arguments on a Windows host; the FMA-contracted build bundled with PPSA13934 differs on 43,433 of the random ones. - [log1p](../../../core/libs/prx/libc/src/Math.cpp) (libc) does not compare floating-point exceptions, and its double inputs are checked by sampling, not exhaustively. It is FreeBSD msun's [`s_log1p.c`](https://github.com/freebsd/freebsd-src/blob/release/11.4.0/lib/msun/src/s_log1p.c). As the titles' `libc.prx` does, it sets `EDOM` at msun's NaN return (arguments below -1, `-inf` and negative NaNs) and leaves `errno` alone for `log1p(-1)`, where it returns `-inf`. It gives the same bits, NaN payloads included, and `errno` as the `log1p` bundled with PPSA07923 and PPSA05916 on 320 million random and 68,216 structured arguments on a Windows host. The FMA-contracted build bundled with PPSA13934 differs on 73,682 of the random ones. -- [nextafter and nextafterf](../../../core/libs/prx/libc/src/Math.cpp) (libc) do not compare floating-point exceptions, and their inputs are checked by sampling, not exhaustively. They are FreeBSD msun's [`s_nextafter.c`](https://github.com/freebsd/freebsd-src/blob/release/11.4.0/lib/msun/src/s_nextafter.c) and [`s_nextafterf.c`](https://github.com/freebsd/freebsd-src/blob/release/11.4.0/lib/msun/src/s_nextafterf.c), with `ERANGE` at msun's overflow and underflow returns as the titles' `libc.prx` sets it; a step away from zero returns early and sets none. The copies bundled with eleven titles give the same bits and `errno` on 4,000,784 sampled pairs per function on a Windows host. +- [nextafter and nextafterf](../../../core/libs/prx/libc/src/Math.cpp) (libc) have no physical-console or unmasked-exception validation. A Linux CPU probe compares result bits, `errno` and masked floating-point flags with the unchanged math code in FAR: Changing Tides PPSA03766 bundled `libc.prx` (`nextafter` 0x51b30, `nextafterf` 0x51c80, `_Feraise` 0x11da0): 16,864 calls across both precisions, four rounding modes and all DAZ/FTZ combinations, not exhaustive inputs. Earlier Windows sampling of eleven bundled libraries compares only bits and `errno`; exception behavior outside the FAR module remains unverified. - [tanh](../../../core/libs/prx/libc/src/Math.cpp) (libc) does not compare floating-point exceptions, and its double inputs are checked by sampling, not exhaustively. It is FreeBSD msun's [`s_tanh.c`](https://github.com/freebsd/freebsd-src/blob/release/11.4.0/lib/msun/src/s_tanh.c) over this file's msun `Expm1`, with `tiny` read at run time; it sets no `errno`, as the titles' `libc.prx` does not. It gives the same bits as the `tanh` bundled with PPSA07923 and PPSA05916 on 16 million random and 68,211 structured arguments on a Windows host; the FMA-contracted build bundled with PPSA13934 differs on 1,054 of them. - [sinhf](../../../core/libs/prx/libc/src/Math.cpp) (libc) follows the `sinhf` of the `libc.prx` bundled with PPSA28319 (disassembly): FreeBSD msun's [`e_sinhf.c`](https://github.com/freebsd/freebsd-src/blob/5aa554c7e5b6fd9042986df07f12ac5687e41392/lib/msun/src/e_sinhf.c) from before `__ldexp_expf`, which squares `expf(|x| / 2)` between `logf(FLT_MAX)` and the overflow threshold, calling msun's [`expm1f`](https://github.com/freebsd/freebsd-src/blob/59bc2f095289a10b8cfa55cbb2b45c3bd6c2d4ae/lib/msun/src/s_expm1f.c) and [`expf`](https://github.com/freebsd/freebsd-src/blob/59bc2f095289a10b8cfa55cbb2b45c3bd6c2d4ae/lib/msun/src/e_expf.c), whose paths for the arguments it passes are inlined here. Past the threshold Sony's error helper sets `errno` to `ERANGE` before returning `x * 1e37`; msun itself sets no `errno`. That code, run on a Windows host, gives the same bits as this implementation for every finite argument up to the threshold, of both signs; it was not run on Linux. That the system `libSceLibcInternal`, which Ruffle (PPSA68091) and RetroArch (PPSA99169) import while bundling a zero-return `libc.prx`, has the same code is assumed - [acoshf](../../../core/libs/prx/libc/src/Math.cpp) (libc) exact system `libSceLibcInternal` behavior and FMA-contracted copies remain unverified. Positive finite inputs use [msun](https://github.com/freebsd/freebsd-src/blob/5aa554c7e5b6fd9042986df07f12ac5687e41392/lib/msun/src/e_acoshf.c) range handling and this file's `logf`; the middle range divides -1 before addition under directed rounding. Bits, errno and SSE/x87 status and controls match the supplied PPSA03766/PPSA07923 copies on 236,592 positive finite calls (4,929 inputs, four rounding modes, DAZ/FTZ combinations and sticky flags); these samples are not exhaustive. Invalid-domain errno and SSE/x87 invalid flags are checked against two supplied bundled copies; those samples are not exhaustive either.