diff --git a/core/libs/prx/libc/src/Math.cpp b/core/libs/prx/libc/src/Math.cpp index 2aecbe9ab..9ecf79e14 100644 --- a/core/libs/prx/libc/src/Math.cpp +++ b/core/libs/prx/libc/src/Math.cpp @@ -7,6 +7,7 @@ #include #include #include +#include #include #include "prx/libc/include/General.hpp" @@ -24,6 +25,106 @@ std::int64_t RoundToInt64(TFloat value) { return static_cast(rounded); } +float SinhfExpm1(float x) { + constexpr float Ln2High = 6.9313812256e-01f, Ln2Low = 9.0580006145e-06f, InverseLn2 = 1.4426950216e+00f; + constexpr float Q1 = -3.3333212137e-2f, Q2 = 1.5807170421e-3f; + const std::uint32_t magnitude = std::bit_cast(x); + int k = 0; + float c = 0.0f; + if (magnitude > 0x3eb17218u) { + float high; + float low; + if (magnitude < 0x3f851592u) { + high = x - Ln2High; + low = Ln2Low; + k = 1; + } else { + k = static_cast(InverseLn2 * x + 0.5f); + const float t = static_cast(k); + high = x - t * Ln2High; + low = t * Ln2Low; + } + x = high - low; + c = (high - x) - low; + } + const float halfX = 0.5f * x; + const float halfSquare = x * halfX; + const float r1 = 1.0f + halfSquare * (Q1 + halfSquare * Q2); + const float t = 3.0f - r1 * halfX; + float e = halfSquare * ((r1 - t) / (6.0f - x * t)); + if (k == 0) return x - (x * e - halfSquare); + e = x * (e - c) - c; + e -= halfSquare; + if (k == 1) { + if (x < -0.25f) return -2.0f * (e - (x + 0.5f)); + return 1.0f + 2.0f * (x - e); + } + const float twoPowK = std::bit_cast(static_cast(0x7f + k) << 23); + const float oneMinusTwoPowMinusK = std::bit_cast(0x3f800000u - (0x1000000u >> k)); + return (oneMinusTwoPowMinusK - (e - x)) * twoPowK; +} + +float SinhfExp(float x) { + constexpr float Ln2High = 6.9314575195e-01f, Ln2Low = 1.4286067653e-06f, InverseLn2 = 1.4426950216e+00f; + constexpr float P1 = 1.6666625440e-1f, P2 = -2.7667332906e-3f; + const int k = static_cast(InverseLn2 * x + 0.5f); + const float kFloat = static_cast(k); + const float high = x - kFloat * Ln2High; + const float low = kFloat * Ln2Low; + x = high - low; + const float t = x * x; + const float c = x - t * (P1 + t * P2); + const float y = 1.0f - ((low - (x * c) / (2.0f - c)) - high); + if (k == 128) return y * 2.0f * 0x1p127f; + return y * std::bit_cast(static_cast(0x7f + k) << 23); +} + +static_assert(std::numeric_limits::digits == 64, "guest long double is the x87 80-bit format"); + +struct ExtendedParts { + std::uint64_t mantissa; + std::uint16_t signExponent; +}; + +ExtendedParts SplitExtended(long double value) { + unsigned char bytes[sizeof(long double)]; + std::memcpy(bytes, &value, sizeof(value)); + ExtendedParts parts; + std::memcpy(&parts.mantissa, bytes, sizeof(parts.mantissa)); + std::memcpy(&parts.signExponent, bytes + sizeof(parts.mantissa), sizeof(parts.signExponent)); + return parts; +} + +long double ClearLowMantissaWord(long double value) { + unsigned char bytes[sizeof(long double)]; + std::memcpy(bytes, &value, sizeof(value)); + std::memset(bytes, 0, sizeof(std::uint32_t)); + std::memcpy(&value, bytes, sizeof(value)); + return value; +} + +constexpr long double AsinNumeratorCoefficients[] = { + 1.66666666666666666631e-01L, -4.16313987993683104320e-01L, 3.69068046323246813704e-01L, -1.36213932016738603108e-01L, + 1.78324189708471965733e-02L, -2.19216428382605211588e-04L, -7.10526623669075243183e-06L, +}; +constexpr long double AsinDenominatorCoefficients[] = { + -2.94788392796209867269e+00L, 3.27309890266528636716e+00L, -1.68285799854822427013e+00L, 3.90699412641738801874e-01L, + -3.14365703596053263322e-02L, +}; +constexpr long double HalfPiHigh = 1.57079632679489661926e+00L; +constexpr long double HalfPiLow = -2.50827880633416601173e-20L; +constexpr long double QuarterPiHigh = 7.85398163397448309628e-01L; + +long double AsinNumerator(long double t) { + const auto& p = AsinNumeratorCoefficients; + return t * (p[0] + t * (p[1] + t * (p[2] + t * (p[3] + t * (p[4] + t * (p[5] + t * p[6])))))); +} + +long double AsinDenominator(long double t) { + const auto& q = AsinDenominatorCoefficients; + return 1.0L + t * (q[0] + t * (q[1] + t * (q[2] + t * (q[3] + t * q[4])))); +} + } extern "C" { @@ -98,6 +199,36 @@ double APS5_VABI trunc_nid_postfix(double x) { return std::trunc(x); } double APS5_VABI fma_nid_postfix(double x, double y, double z) { return std::fma(x, y, z); } long double APS5_VABI acosl_nid_postfix(long double x) { return std::acos(x); } +long double APS5_VABI asinl_nid_postfix(long double x) { + const ExtendedParts parts = SplitExtended(x); + const int exponent = parts.signExponent & 0x7fff; + if (exponent >= 0x3fff) { + if (exponent == 0x3fff && (parts.mantissa & 0x7fffffffffffffffull) == 0) return x * HalfPiHigh + x * HalfPiLow; + std::feraiseexcept(FE_INVALID); + errno = EDOM; + return (x - x) / (x - x); + } + if (exponent < 0x3ffe) { + if (exponent < 0x3fff - 32 && 1.0e300L + x > 1.0L) return x; + const long double t = x * x; + return x + x * (AsinNumerator(t) / AsinDenominator(t)); + } + const long double t = (1.0L - std::fabs(x)) * 0.5L; + const long double p = AsinNumerator(t); + const long double q = AsinDenominator(t); + const long double s = std::sqrt(t); + long double result; + if ((parts.mantissa >> 32) >= 0xf3333333u) { + result = HalfPiHigh - (2.0L * (s + s * (p / q)) - HalfPiLow); + } else { + const long double w = ClearLowMantissaWord(s); + const long double c = (t - w * w) / (s + w); + const long double r = p / q; + result = QuarterPiHigh - ((2.0L * s * r - (HalfPiLow - 2.0L * c)) - (QuarterPiHigh - 2.0L * w)); + } + return (parts.signExponent & 0x8000) != 0 ? -result : result; +} + long double APS5_VABI frexpl_nid_postfix(long double x, int* exponent) { if (!std::isfinite(x)) return x; return std::frexp(x, exponent); @@ -124,6 +255,27 @@ double APS5_VABI cos_nid_postfix(double x) { return std::cos(x); } double APS5_VABI sinh_nid_postfix(double x) { return std::sinh(x); } double APS5_VABI cosh_nid_postfix(double x) { return std::cosh(x); } float APS5_VABI coshf_nid_postfix(float x) { return std::cosh(x); } + +float APS5_VABI sinhf_nid_postfix(float x) { + const std::uint32_t magnitude = std::bit_cast(x) & 0x7fffffffu; + if (magnitude >= 0x7f800000u) return x + x; + const float half = std::signbit(x) ? -0.5f : 0.5f; + const float absolute = std::fabs(x); + if (magnitude < 0x41100000u) { + if (magnitude < 0x39800000u && 1.0e37f + x > 1.0f) return x; + const float t = SinhfExpm1(absolute); + if (magnitude < 0x3f800000u) return half * (2.0f * t - t * t / (t + 1.0f)); + return half * (t + t / (t + 1.0f)); + } + if (magnitude < 0x42b17217u) return half * SinhfExp(absolute); + if (magnitude <= 0x42b2d4fcu) { + const float w = SinhfExp(0.5f * absolute); + return half * w * w; + } + errno = ERANGE; + return x * 1.0e37f; +} + double APS5_VABI asinh_nid_postfix(double x) { return std::asinh(x); } float APS5_VABI asinhf_nid_postfix(float x) { return std::asinh(x); } double APS5_VABI acosh_nid_postfix(double x) { return std::acosh(x); } diff --git a/core/libs/tests/GuestMath.cpp b/core/libs/tests/GuestMath.cpp index f69166847..5a8e85f5f 100644 --- a/core/libs/tests/GuestMath.cpp +++ b/core/libs/tests/GuestMath.cpp @@ -35,6 +35,7 @@ int APS5_VABI abs_nid_postfix(int); double APS5_VABI trunc_nid_postfix(double); double APS5_VABI fma_nid_postfix(double, double, double); long double APS5_VABI acosl_nid_postfix(long double); +long double APS5_VABI asinl_nid_postfix(long double); long double APS5_VABI frexpl_nid_postfix(long double, int*); struct LibcFloatConstant { std::uint32_t bits[4]; }; extern LibcFloatConstant _FInf_nid_postfix; @@ -73,6 +74,7 @@ std::div_t APS5_VABI div_nid_postfix(int, int); double APS5_VABI sinh_nid_postfix(double); double APS5_VABI cosh_nid_postfix(double); float APS5_VABI coshf_nid_postfix(float); +float APS5_VABI sinhf_nid_postfix(float); double APS5_VABI asinh_nid_postfix(double); float APS5_VABI asinhf_nid_postfix(float); double APS5_VABI acosh_nid_postfix(double); @@ -197,6 +199,49 @@ static void CheckFmodErrno() { Require(fmod_nid_postfix(1e300, 3.0) == std::fmod(1e300, 3.0)); } +static long double MakeExtended(std::uint64_t mantissa, std::uint16_t signExponent) { + unsigned char bytes[sizeof(long double)]{}; + std::memcpy(bytes, &mantissa, 8); + std::memcpy(bytes + 8, &signExponent, 2); + long double value = 0; + std::memcpy(&value, bytes, sizeof(value)); + return value; +} + +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}, + {0x8000000000000000ull, 0x3fd7, 0x8000000000000000ull, 0x3fd7}, {0x8000000000000000ull, 0xbfd7, 0x8000000000000000ull, 0xbfd7}, + {0x8000000000000000ull, 0x3ffd, 0x815f4e705c3e68f3ull, 0x3ffd}, {0x8000000000000000ull, 0xbffd, 0x815f4e705c3e68f3ull, 0xbffd}, + {0xcccccccccccccccdull, 0x3ffb, 0xcd2493b01bc41e02ull, 0x3ffb}, + {0x8000000000000000ull, 0x3ffe, 0x860a91c16b9b2c23ull, 0x3ffe}, {0x8000000000000000ull, 0xbffe, 0x860a91c16b9b2c23ull, 0xbffe}, + {0xc000000000000000ull, 0xbffe, 0xd91a98ae3406e041ull, 0xbffe}, {0xe666666666666666ull, 0x3ffe, 0x8f549b82fd0b8fa3ull, 0x3fff}, + {0xf333333333333333ull, 0x3ffe, 0xa06a08acf9a04f78ull, 0x3fff}, {0xfd70a3d70a3d70a4ull, 0xbffe, 0xb6f1e3799ea8ecc3ull, 0xbfff}, + {0xffffffffffffffffull, 0x3ffe, 0xc90fdaa16c63cf01ull, 0x3fff}, + {0x8000000000000000ull, 0x3fff, 0xc90fdaa22168c235ull, 0x3fff}, {0x8000000000000000ull, 0xbfff, 0xc90fdaa22168c235ull, 0xbfff}, + }; + for (const auto& test : cases) { + *__error_nid_postfix() = 12345; + const long double result = asinl_nid_postfix(MakeExtended(test.inputMantissa, test.inputSignExponent)); + if (!ExtendedBits(result, test.outputMantissa, test.outputSignExponent) || *__error_nid_postfix() != 12345) { + std::fprintf(stderr, "Guest asinl failed for %016llx:%04x\n", static_cast(test.inputMantissa), test.inputSignExponent); + std::abort(); + } + } + Require(std::fabs(asinl_nid_postfix(0.5L) - 0.523598775598298873077107230547L) <= 0x1p-64L); + const long double domainErrors[] = { + 1.5L, -2.0L, MakeExtended(0x8000000000000001ull, 0x3fff), MakeExtended(0x8000000000000001ull, 0xbfff), + std::numeric_limits::max(), std::numeric_limits::infinity(), -std::numeric_limits::infinity(), + std::numeric_limits::quiet_NaN(), + }; + for (const long double input : domainErrors) { + std::feclearexcept(FE_ALL_EXCEPT); + *__error_nid_postfix() = 0; + Require(std::isnan(asinl_nid_postfix(input)) && *__error_nid_postfix() == EDOM && std::fetestexcept(FE_INVALID) != 0); + } + *__error_nid_postfix() = 0; +} + static void CheckIntegerConversions() { for (const long long numerator : {4294967301LL, -4294967301LL}) { for (const long long denominator : {3LL, -3LL}) { @@ -414,6 +459,38 @@ static void CheckHyperbolic() { Require(Near(acosh_nid_postfix(1e300), 691.4686750787736, 1e-15)); } +static void CheckSinhf() { + const float infinity = std::numeric_limits::infinity(); + const struct { std::uint32_t input; std::uint32_t output; } cases[] = { + {0x00000000u, 0x00000000u}, {0x80000000u, 0x80000000u}, {0x00000001u, 0x00000001u}, {0x39000000u, 0x39000000u}, + {0x3e800000u, 0x3e815667u}, {0xbf000000u, 0xbf056680u}, {0x3f800000u, 0x3f966cfeu}, {0xbf800000u, 0xbf966cfeu}, + {0x40000000u, 0x40681e7bu}, {0x41100000u, 0x457d38acu}, {0xc1a00000u, 0xcd675844u}, {0x42b17217u, 0x7effff84u}, + {0x42b20000u, 0x7f28e166u}, {0xc2b20000u, 0xff28e166u}, {0x42b2d4fcu, 0x7f7fffecu}, + {0x7f800000u, 0x7f800000u}, {0xff800000u, 0xff800000u}, + }; + for (const auto& test : cases) { + *__error_nid_postfix() = 12345; + const float result = sinhf_nid_postfix(std::bit_cast(test.input)); + if (std::bit_cast(result) != test.output || *__error_nid_postfix() != 12345) { + std::fprintf(stderr, "Guest sinhf failed for %08x\n", test.input); + std::abort(); + } + } + Require(std::isnan(sinhf_nid_postfix(std::numeric_limits::quiet_NaN())) && *__error_nid_postfix() == 12345); + Require(Near(sinhf_nid_postfix(1.f), 1.1752011936438014, 1e-6) && Near(sinhf_nid_postfix(-2.f), -3.626860407847019, 1e-6)); + Require(Near(sinhf_nid_postfix(0.25f), 0.25261231680816831, 1e-6) && Near(sinhf_nid_postfix(20.f), 242582597.70489514, 1e-6)); + for (const std::uint32_t input : {0x42b2d4fdu, 0x42c80000u, 0x7f7fffffu}) { + for (const std::uint32_t sign : {0u, 0x80000000u}) { + std::feclearexcept(FE_ALL_EXCEPT); + *__error_nid_postfix() = 0; + const float result = sinhf_nid_postfix(std::bit_cast(input | sign)); + Require(result == (sign != 0 ? -infinity : infinity) && *__error_nid_postfix() == ERANGE); + Require(std::fetestexcept(FE_OVERFLOW) != 0); + } + } + *__error_nid_postfix() = 0; +} + static void CheckFloatClassification() { Require(_FInf_nid_postfix.bits[0] == 0x7f800000u && _FNan_nid_postfix.bits[0] == 0x7fc00000u); for (int word = 1; word < 4; ++word) Require(_FInf_nid_postfix.bits[word] == 0 && _FNan_nid_postfix.bits[word] == 0); @@ -597,6 +674,7 @@ int main() { CheckExp2Errno(); CheckExactOperations(); CheckLongDouble(); + CheckAsinl(); Require(scalbn_nid_postfix(0.75, -2) == 0.1875); Require(scalbnf_nid_postfix(0.75f, 4) == 12.f); Require(nextafterf_nid_postfix(1.f, 2.f) == 0x1.000002p0f && nextafterf_nid_postfix(1.f, 0.f) == 0x1.fffffep-1f); @@ -679,6 +757,7 @@ int main() { Require(std::isinf(__powisf2_nid_postfix(-0.0f, -3)) && std::signbit(__powisf2_nid_postfix(-0.0f, -3))); Require(__powisf2_nid_postfix(2.0f, -127) == 0x1p-127f); CheckHyperbolic(); + CheckSinhf(); CheckHypotSpecialValue(); CheckFmodErrno(); } diff --git a/docs/dev/TechnicalDebt.md b/docs/dev/TechnicalDebt.md index d41233dfc..78d759cbd 100644 --- a/docs/dev/TechnicalDebt.md +++ b/docs/dev/TechnicalDebt.md @@ -445,6 +445,8 @@ Throughout the project, every function at every stage either **does exactly what - [nan and nanf](../../core/libs/prx/libc/src/Math.cpp) (libc) return fixed quiet NaN values (`0x7fc00000` and `0x7ff8000000000000`), ignoring the tag argument as observed in the console's `libc.prx` (both exports load a global constant via `mov rax,[rip+..]; vmovsd xmm0,[rax]; ret`); VLC (PPSA85300) stops on nanf at start-up after strsep. - [Draw input snapshots](../../core/libs/prx/libSceAgcDriver/Graphics/src/ShaderResources.cpp) (libSceAgcDriver): `PrepareDrawBindings` caps read-only input buffer snapshots at 64 MiB (overridable via `APS5_DRAW_SNAPSHOT_MAX_MIB`). Larger inputs (streaming heaps, bindless pools) or draws exceeding the cumulative budget skip the per-draw snapshot and bind in place, avoiding multi-gigabyte CPU copy stalls (such as Ratchet & Clank: Rift Apart PPSA01474 stalling at frame ~380 copying 36 GB every 5 s). If the title rewrites those large ranges on the CPU before GPU completion, the in-place binding observes the modification, as verified by `drawSnapshotCapTests`. - [pthread_attr_getinheritsched, pthread_attr_setinheritsched, pthread_attr_getscope, pthread_attr_setscope, pthread_attr_getstackaddr](../../core/libs/prx/libkernel/Pthread/Posix/Attr.cpp) (libkernel) follow FreeBSD 11 libthr ([thr_attr.c](https://github.com/freebsd/freebsd-src/blob/5ed7eb0d97ba4436218e810f61bd059acba984c2/lib/libthr/thread/thr_attr.c), [thr_create.c](https://github.com/freebsd/freebsd-src/blob/5ed7eb0d97ba4436218e810f61bd059acba984c2/lib/libthr/thread/thr_create.c)), not measured on the console: the scope defaults to `PTHREAD_SCOPE_SYSTEM` (2), only it and `PTHREAD_SCOPE_PROCESS` (0) are accepted and other values return `EINVAL`, and an inherit value other than `PTHREAD_INHERIT_SCHED` (4) or `PTHREAD_EXPLICIT_SCHED` (0) returns `ENOTSUP`. `pthread_attr_getstackaddr` returns the configured stack address pointer (null by default). A null output of `pthread_attr_getinheritsched` throws, since libthr writes through it unchecked. A thread keeps the scope it was created with, taken from the creating thread under `PTHREAD_INHERIT_SCHED`, and `pthread_attr_get_np` and `scePthreadAttrGet` report it, but every guest thread is a host thread at system scope, so process scope changes nothing +- [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 +- [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 - [SPIRV-Tools optimization](../../core/shader/recompiler/SpirvBackend/src/SpirvOptimizer.cpp) skips CCP and Simplification for modules with storage-only 8- or 16-bit types because the pinned folder produces illegal narrow constants. This also disables unrelated constant folding; the restricted full pipeline mirrors the pinned performance recipe. - [SRT signed comparison](../../core/shader/recompiler/Optimization/src/SrtWalker/SrtEvaluator.cpp) evaluates `SLessThan32` as signed 32-bit less-than, matching the existing SPIR-V backend and [SPIR-V `OpSLessThan`](https://registry.khronos.org/SPIR-V/specs/unified1/SPIRV.html#OpSLessThan). Signed boundary descriptor expressions are tested on the host; native RDNA2/PS5 execution and title behavior remain unverified. - On macOS the [write watch](../../core/libs/prx/libc/src/GuestWriteWatch.cpp) write-protects watched pages and records the first write to each from a SIGBUS/SIGSEGV handler, because macOS has no counterpart to `userfaultfd` write-protect or `MEM_WRITE_WATCH`. The protection is real, so: