diff --git a/core/libs/prx/libc/src/Math.cpp b/core/libs/prx/libc/src/Math.cpp index d88fc6d00b..43c3546cce 100644 --- a/core/libs/prx/libc/src/Math.cpp +++ b/core/libs/prx/libc/src/Math.cpp @@ -1336,7 +1336,31 @@ double APS5_VABI nan_nid_postfix(const char* tagp) { double APS5_VABI modf_nid_postfix(double x, double* integral) { return std::modf(x, integral); } float APS5_VABI modff_nid_postfix(float x, float* integral) { return std::modf(x, integral); } -double APS5_VABI tanh_nid_postfix(double x) { return std::tanh(x); } +double APS5_VABI tanh_nid_postfix(double x) { + static volatile double tiny = 1.0e-300; + const double huge = 1.0e300; + const std::int32_t jx = static_cast(std::bit_cast(x) >> 32); + const std::int32_t ix = jx & 0x7fffffff; + if (ix >= 0x7ff00000) { + return jx >= 0 ? 1.0 / x + 1.0 : 1.0 / x - 1.0; + } + double z; + if (ix < 0x40360000) { + if (ix < 0x3e300000 && huge + x > 1.0) { + return x; + } + if (ix >= 0x3ff00000) { + const double t = Expm1(2.0 * std::fabs(x)); + z = 1.0 - 2.0 / (t + 2.0); + } else { + const double t = Expm1(-2.0 * std::fabs(x)); + z = -t / (t + 2.0); + } + } else { + z = 1.0 - tiny; + } + return jx >= 0 ? z : -z; +} float APS5_VABI tanhf_nid_postfix(float x) { return std::tanh(x); } long double APS5_VABI atanl_nid_postfix(long double x) { return std::atan(x); } long double APS5_VABI sinl_nid_postfix(long double x) { return std::sin(x); } diff --git a/core/libs/tests/GuestMath.cpp b/core/libs/tests/GuestMath.cpp index 8297d06458..9fd4dbd27a 100644 --- a/core/libs/tests/GuestMath.cpp +++ b/core/libs/tests/GuestMath.cpp @@ -630,6 +630,27 @@ static void CheckLog1pConsoleBits() { } } +static void CheckTanhConsoleBits() { + constexpr std::uint64_t cases[][2] = { + {0x0000000000000000ull, 0x0000000000000000ull}, {0x8000000000000000ull, 0x8000000000000000ull}, + {0x3e2fffffffffffffull, 0x3e2fffffffffffffull}, {0x3e30000000000000ull, 0x3e30000000000000ull}, + {0xbe30000000000000ull, 0xbe30000000000000ull}, {0x3fefffffffffffffull, 0x3fe85efab514f394ull}, + {0x3ff0000000000000ull, 0x3fe85efab514f394ull}, {0xbff0000000000000ull, 0xbfe85efab514f394ull}, + {0x3fe0000000000000ull, 0x3fdd9353d7568af3ull}, {0x4035ffffffffffffull, 0x3ff0000000000000ull}, + {0x4036000000000000ull, 0x3ff0000000000000ull}, {0xc036000000000000ull, 0xbff0000000000000ull}, + {0x4049000000000000ull, 0x3ff0000000000000ull}, {0x7fefffffffffffffull, 0x3ff0000000000000ull}, + {0x7ff0000000000000ull, 0x3ff0000000000000ull}, {0xfff0000000000000ull, 0xbff0000000000000ull}, + {0x7ff8000000000001ull, 0x7ff8000000000001ull}, {0xfff8000000000002ull, 0xfff8000000000002ull}, + {0x3f927ed0db8370b1ull, 0x3f927e4d1092a580ull}, {0x3f95d74e4cee78bfull, 0x3f95d67547d523a0ull}, + {0x3f18f1b062cb523cull, 0x3f18f1b06187faaeull}, {0x3f937129f25ae835ull, 0x3f937090dd8c9cb3ull}, + {0x3facb0d89f8aa91cull, 0x3faca92afe31a32bull}, {0x3f07f462f020acdeull, 0x3f07f462efd91530ull}, + }; + for (const auto& [input, expected] : cases) { + errno = 1234; + Require(std::bit_cast(tanh_nid_postfix(std::bit_cast(input))) == expected && errno == 1234); + } +} + 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}, @@ -1722,6 +1743,7 @@ int main() { CheckExpm1ConsoleBits(); CheckLog1pConsoleBits(); CheckNextafterConsoleBits(); + CheckTanhConsoleBits(); 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 c648774059..40e0f9edac 100644 --- a/docs/dev/debt/libc.md +++ b/docs/dev/debt/libc.md @@ -114,6 +114,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. - [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. +- [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) 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