From c05345f529b21d9e53fa29efa672be11c2f2af6c Mon Sep 17 00:00:00 2001 From: JJenkx Date: Sat, 24 Oct 2026 11:32:05 -0400 Subject: [PATCH 1/3] fix(libc): preserve acoshf log and reciprocal evaluation --- core/libs/prx/libc/src/Math.cpp | 4 ++-- docs/dev/debt/libc.md | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/core/libs/prx/libc/src/Math.cpp b/core/libs/prx/libc/src/Math.cpp index e635b0b1e1..4e53efaec8 100644 --- a/core/libs/prx/libc/src/Math.cpp +++ b/core/libs/prx/libc/src/Math.cpp @@ -1045,10 +1045,10 @@ float APS5_VABI acoshf_nid_postfix(float x) { if ((bits & 0x80000000u) != 0 || magnitude < 0x3f800000u) return (x - x) / (x - x); if (magnitude >= 0x7f800000u) return x + x; if (magnitude == 0x3f800000u) return 0.0f; - if (magnitude >= 0x4d800000u) return std::log(x) + ln2; + if (magnitude >= 0x4d800000u) return logf_nid_postfix(x) + ln2; if (magnitude > 0x40000000u) { const float square = x * x; - return std::log(2.0f * x - 1.0f / (x + std::sqrt(square - 1.0f))); + return logf_nid_postfix(2.0f * x + (-1.0f / (x + std::sqrt(square - 1.0f)))); } const float t = x - 1.0f; return AtanhfLog1p(t + std::sqrt(2.0f * t + t * t)); diff --git a/docs/dev/debt/libc.md b/docs/dev/debt/libc.md index 5171721bc8..586cbb6caa 100644 --- a/docs/dev/debt/libc.md +++ b/docs/dev/debt/libc.md @@ -111,7 +111,7 @@ - [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. - [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) uses FreeBSD msun's [`e_acoshf.c`](https://github.com/freebsd/freebsd-src/blob/5aa554c7e5b6fd9042986df07f12ac5687e41392/lib/msun/src/e_acoshf.c) range handling; the Windows UCRT differed from the `libc.prx` copies in PPSA07923 and PPSA13934 on 216,017 and 216,024 of 20 million inputs in the #9707 survey. The invalid-domain path raises `FE_INVALID` without setting `errno` as the source does. The exact system `libSceLibcInternal` behavior and FMA-contracted results remain unverified +- [acoshf](../../../core/libs/prx/libc/src/Math.cpp) (libc) does not reproduce bundled invalid-domain errno or x87 exception flags (#10610); 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. - [acoshl, asinhl, atanhl, atan2l, cbrtl, erfl, erfcl, expm1l, hypotl, log1pl, log2l, lgammal, tgammal, j0f, j1f, jnf, y0f, y1f and ynf](../../../core/libs/prx/libc/src/Math.cpp) (libc) return the host C library's results (the float Bessel functions through its double ones), not those of the console's own code, so the last bits can differ (`y0f`, `y1f`, `ynf` follow the console at `+-0`). As in the console's `libSceLibcInternal.sprx` (SDK 4.5x build of `J03950001`, 0x56410), `lgammal` does not store the sign in `signgam`; `hypotl` returns infinity for an infinite argument even when the other is a NaN, as C11 F.10.4.3 requires. [hypot3, hypot3f and hypot3l](../../../core/libs/prx/libc/src/Math.cpp) follow that module (0x51D30, 0x520D0, 0x52420): the first NaN of x, y and z is returned unchanged, any infinity gives positive infinity, a zero argument gives the two-argument hypot of the other two, and otherwise the largest magnitude is combined with the hypot of the other two - [asinl](../../../core/libs/prx/libc/src/Math.cpp) (libc) follows the `asinl` of the `libc.prx` bundled with PPSA28319 (disassembly), FreeBSD msun's [`e_asinl.c`](https://github.com/freebsd/freebsd-src/blob/9f0432041b8fb69832f10c32da50a67f34feeb61/lib/msun/src/e_asinl.c) with the ld80 [`invtrig.c`](https://github.com/freebsd/freebsd-src/blob/9f0432041b8fb69832f10c32da50a67f34feeb61/lib/msun/ld80/invtrig.c) coefficients, which the binary holds as the same 80-bit constants. For |x| > 1, the infinities and NaN, Sony's error helper raises `FE_INVALID` and sets `errno` to `EDOM` before returning `(x - x) / (x - x)`, so a NaN argument also sets `EDOM`; msun itself sets no `errno`. 20 million random arguments in (-1, 1) give the same bits as that code run on a Windows host. The result needs the x87 precision control at 64 bits, as FreeBSD amd64 starts threads (control word `0x37f`), but guest threads do not set it: on Windows a thread started by `_beginthreadex` runs at 53 bits (`0x27f`), where `asinl(0.75)` returns mantissa `0xd91a98ae3406e000` instead of `0xd91a98ae3406e041`. The system `libSceLibcInternal` that RetroArch (PPSA99169) imports is assumed to have the same code - [atanhf](../../../core/libs/prx/libc/src/Math.cpp) (libc) follows the `atanhf` of the `libc.prx` bundled with PPSA28319 (disassembly): FreeBSD msun's [`e_atanhf.c`](https://github.com/freebsd/freebsd-src/blob/59bc2f095289a10b8cfa55cbb2b45c3bd6c2d4ae/lib/msun/src/e_atanhf.c) calling msun's [`log1pf`](https://github.com/freebsd/freebsd-src/blob/59bc2f095289a10b8cfa55cbb2b45c3bd6c2d4ae/lib/msun/src/s_log1pf.c), whose paths for the positive arguments it passes are inlined here; the exported `log1pf` is still the host's. For |x| > 1, the infinities and NaN, Sony's error helper raises `FE_INVALID` and sets `errno` to `EDOM` before returning `(x - x) / (x - x)`, and for |x| = 1 it raises `FE_DIVBYZERO` and sets `errno` to `ERANGE` before returning `x / 0`; msun itself sets no `errno`. That code, run on a Windows host, gives the same bits, `errno` and floating-point exceptions as this implementation for all 2^32 arguments; it was not run on Linux. That the system `libSceLibcInternal`, which Dusklight (PPSA99305) imports while bundling a zero-return `libc.prx`, has the same code is assumed From 0e895b55be29357f6c4bbcd7a10c37965816941e Mon Sep 17 00:00:00 2001 From: JJenkx Date: Sat, 24 Oct 2026 11:29:18 -0400 Subject: [PATCH 2/3] test(libc): cover acoshf positive results across rounding modes --- core/libs/tests/GuestMath.cpp | 46 +++++++++++++++++++++++++++++++++++ 1 file changed, 46 insertions(+) diff --git a/core/libs/tests/GuestMath.cpp b/core/libs/tests/GuestMath.cpp index 7c6bd1e44e..4b2d8bcb55 100644 --- a/core/libs/tests/GuestMath.cpp +++ b/core/libs/tests/GuestMath.cpp @@ -937,6 +937,52 @@ static void CheckAcoshfMsunBoundaries() { Require(std::isinf(acoshf_nid_postfix(INFINITY)) && acoshf_nid_postfix(INFINITY) > 0.0f); Require(std::isnan(acoshf_nid_postfix(-INFINITY))); Require(std::isnan(acoshf_nid_postfix(std::numeric_limits::quiet_NaN()))); + const RestoreFloatingPointState saved; + const unsigned savedCsr = _mm_getcsr(); + const struct { + std::uint32_t input; + std::uint32_t result[4]; + unsigned raised; + } cases[] = { + {0x3f800000u, {0x00000000u, 0x00000000u, 0x00000000u, 0x00000000u}, 0u}, + {0x3f800001u, {0x39ffffffu, 0x39ffffffu, 0x3a000001u, 0x39ffffffu}, 32u}, + {0x3fffffffu, {0x3fa89213u, 0x3fa89212u, 0x3fa89214u, 0x3fa89212u}, 32u}, + {0x40000000u, {0x3fa89214u, 0x3fa89213u, 0x3fa89214u, 0x3fa89213u}, 32u}, + {0x40000001u, {0x3fa89215u, 0x3fa89214u, 0x3fa89215u, 0x3fa89214u}, 32u}, + {0x40000002u, {0x3fa89216u, 0x3fa89215u, 0x3fa89217u, 0x3fa89215u}, 32u}, + {0x40000003u, {0x3fa89217u, 0x3fa89217u, 0x3fa89218u, 0x3fa89217u}, 32u}, + {0x4d7fffffu, {0x41a0cf66u, 0x41a0cf65u, 0x41a0cf66u, 0x41a0cf65u}, 32u}, + {0x4d800000u, {0x41a0cf66u, 0x41a0cf64u, 0x41a0cf66u, 0x41a0cf64u}, 32u}, + {0x4d800001u, {0x41a0cf66u, 0x41a0cf65u, 0x41a0cf67u, 0x41a0cf65u}, 32u}, + {0x7f7fffffu, {0x42b2d4fcu, 0x42b2d4fbu, 0x42b2d4fdu, 0x42b2d4fbu}, 32u}, + }; + const int rounding[] = {FE_TONEAREST, FE_DOWNWARD, FE_UPWARD, FE_TOWARDZERO}; + float (APS5_VABI * volatile call)(float) = acoshf_nid_postfix; + for (unsigned mode = 0; mode != 4; ++mode) { + for (unsigned denorm = 0; denorm != 4; ++denorm) { + const unsigned control = 0x1f80u | (mode << 13u) | ((denorm & 1u) ? 0x40u : 0u) | ((denorm & 2u) ? 0x8000u : 0u); + for (std::uint16_t sticky : {0u, 5u, 32u}) { + for (const auto& item : cases) { + Require(std::fesetenv(FE_DFL_ENV) == 0 && std::fesetround(rounding[mode]) == 0); + Require(fesetexceptflag_nid_postfix(&sticky, 0x3f) == 0); + _mm_setcsr(control | sticky); + GuestFloatingPointEnvironment before; + Require(fegetenv_nid_postfix(&before) == 0 && (before.status & 0x3fu) == sticky); + errno = 1234; + const auto result = std::bit_cast(call(std::bit_cast(item.input))); + const int error = errno; + GuestFloatingPointEnvironment after; + Require(fegetenv_nid_postfix(&after) == 0); + if (result != item.result[mode] || after.mxcsr != (control | sticky | item.raised) || (after.status & 0x3fu) != sticky || after.control != before.control || error != 1234) { + std::fprintf(stderr, "acoshf positive: input=%08llx mode=%u denorm=%u sticky=%u bits=%08llx expected=%08llx csr=%x x87=%x errno=%d\n", static_cast(item.input), mode, denorm, sticky, static_cast(result), static_cast(item.result[mode]), after.mxcsr, after.status, error); + Require(false); + } + } + } + } + } + _mm_setcsr(savedCsr); + } static void CheckCoshfMsunBoundaries() { From 2fc4ca2fd0f951e29fbb2cbdfb538ccb81e907a1 Mon Sep 17 00:00:00 2001 From: JJenkx Date: Sun, 25 Oct 2026 07:53:43 -0400 Subject: [PATCH 3/3] test(libc): reject host log fallback in large acoshf inputs --- core/libs/tests/GuestMath.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/core/libs/tests/GuestMath.cpp b/core/libs/tests/GuestMath.cpp index 4b2d8bcb55..af198e192b 100644 --- a/core/libs/tests/GuestMath.cpp +++ b/core/libs/tests/GuestMath.cpp @@ -954,6 +954,7 @@ static void CheckAcoshfMsunBoundaries() { {0x4d7fffffu, {0x41a0cf66u, 0x41a0cf65u, 0x41a0cf66u, 0x41a0cf65u}, 32u}, {0x4d800000u, {0x41a0cf66u, 0x41a0cf64u, 0x41a0cf66u, 0x41a0cf64u}, 32u}, {0x4d800001u, {0x41a0cf66u, 0x41a0cf65u, 0x41a0cf67u, 0x41a0cf65u}, 32u}, + {0x4d8050f2u, {0x41a0d473u, 0x41a0d472u, 0x41a0d474u, 0x41a0d472u}, 32u}, {0x7f7fffffu, {0x42b2d4fcu, 0x42b2d4fbu, 0x42b2d4fdu, 0x42b2d4fbu}, 32u}, }; const int rounding[] = {FE_TONEAREST, FE_DOWNWARD, FE_UPWARD, FE_TOWARDZERO};