From c07686217324aecc33366d9f2a823d9990f39bdb Mon Sep 17 00:00:00 2001 From: Viviane Potocnik Date: Fri, 20 Oct 2023 21:39:28 +0200 Subject: [PATCH] math: add safe casts to log function --- ...-math-add-safe-casts-to-log-function.patch | 86 +++++++++++++++++++ sw/math/src/math/log2.c | 20 ++--- 2 files changed, 96 insertions(+), 10 deletions(-) create mode 100644 sw/deps/patches/musl/0009-math-add-safe-casts-to-log-function.patch diff --git a/sw/deps/patches/musl/0009-math-add-safe-casts-to-log-function.patch b/sw/deps/patches/musl/0009-math-add-safe-casts-to-log-function.patch new file mode 100644 index 000000000..2920149c5 --- /dev/null +++ b/sw/deps/patches/musl/0009-math-add-safe-casts-to-log-function.patch @@ -0,0 +1,86 @@ +From 246137e8768c1a7c9c74e76d494ae764c40fc85e Mon Sep 17 00:00:00 2001 +From: Viviane Potocnik +Date: Fri, 20 Oct 2023 21:38:16 +0200 +Subject: [PATCH] math: add safe casts to log function + +--- + src/math/log2.c | 20 ++++++++++---------- + 1 file changed, 10 insertions(+), 10 deletions(-) + +diff --git a/src/math/log2.c b/src/math/log2.c +index 1276ed4e..a0c669cf 100644 +--- a/src/math/log2.c ++++ b/src/math/log2.c +@@ -22,7 +22,7 @@ + /* Top 16 bits of a double. */ + static inline uint32_t top16(double x) + { +- return asuint64(x) >> 48; ++ return asuint64_safe(x) >> 48; + } + + double log2(double x) +@@ -32,14 +32,14 @@ double log2(double x) + uint32_t top; + int k, i; + +- ix = asuint64(x); ++ ix = asuint64_safe(x); + top = top16(x); +-#define LO asuint64(1.0 - 0x1.5b51p-5) +-#define HI asuint64(1.0 + 0x1.6ab2p-5) ++#define LO asuint64_safe(1.0 - 0x1.5b51p-5) ++#define HI asuint64_safe(1.0 + 0x1.6ab2p-5) + if (predict_false(ix - LO < HI - LO)) { + /* Handle close to 1.0 inputs separately. */ + /* Fix sign of zero with downward rounding when x==1. */ +- if (WANT_ROUNDING && predict_false(ix == asuint64(1.0))) ++ if (WANT_ROUNDING && predict_false(ix == asuint64_safe(1.0))) + return 0; + r = x - 1.0; + #if __FP_FAST_FMA +@@ -47,7 +47,7 @@ double log2(double x) + lo = r * InvLn2lo + __builtin_fma(r, InvLn2hi, -hi); + #else + double_t rhi, rlo; +- rhi = asdouble(asuint64(r) & -1ULL << 32); ++ rhi = asdouble_safe(asuint64_safe(r) & -1ULL << 32); + rlo = r - rhi; + hi = rhi * InvLn2hi; + lo = rlo * InvLn2hi + r * InvLn2lo; +@@ -67,12 +67,12 @@ double log2(double x) + /* x < 0x1p-1022 or inf or nan. */ + if (ix * 2 == 0) + return __math_divzero(1); +- if (ix == asuint64(INFINITY)) /* log(inf) == inf. */ ++ if (ix == asuint64_safe(INFINITY)) /* log(inf) == inf. */ + return x; + if ((top & 0x8000) || (top & 0x7ff0) == 0x7ff0) + return __math_invalid(x); + /* x is subnormal, normalize it. */ +- ix = asuint64(x * 0x1p52); ++ ix = asuint64_safe(x * 0x1p52); + ix -= 52ULL << 52; + } + +@@ -85,7 +85,7 @@ double log2(double x) + iz = ix - (tmp & 0xfffULL << 52); + invc = T[i].invc; + logc = T[i].logc; +- z = asdouble(iz); ++ z = asdouble_safe(iz); + kd = (double_t)k; + + /* log2(x) = log2(z/c) + log2(c) + k. */ +@@ -99,7 +99,7 @@ double log2(double x) + double_t rhi, rlo; + /* rounding error: 0x1p-55/N + 0x1p-65. */ + r = (z - T2[i].chi - T2[i].clo) * invc; +- rhi = asdouble(asuint64(r) & -1ULL << 32); ++ rhi = asdouble_safe(asuint64_safe(r) & -1ULL << 32); + rlo = r - rhi; + t1 = rhi * InvLn2hi; + t2 = rlo * InvLn2hi + r * InvLn2lo; +-- +2.31.1 + diff --git a/sw/math/src/math/log2.c b/sw/math/src/math/log2.c index 1276ed4e3..a0c669cf0 100644 --- a/sw/math/src/math/log2.c +++ b/sw/math/src/math/log2.c @@ -22,7 +22,7 @@ /* Top 16 bits of a double. */ static inline uint32_t top16(double x) { - return asuint64(x) >> 48; + return asuint64_safe(x) >> 48; } double log2(double x) @@ -32,14 +32,14 @@ double log2(double x) uint32_t top; int k, i; - ix = asuint64(x); + ix = asuint64_safe(x); top = top16(x); -#define LO asuint64(1.0 - 0x1.5b51p-5) -#define HI asuint64(1.0 + 0x1.6ab2p-5) +#define LO asuint64_safe(1.0 - 0x1.5b51p-5) +#define HI asuint64_safe(1.0 + 0x1.6ab2p-5) if (predict_false(ix - LO < HI - LO)) { /* Handle close to 1.0 inputs separately. */ /* Fix sign of zero with downward rounding when x==1. */ - if (WANT_ROUNDING && predict_false(ix == asuint64(1.0))) + if (WANT_ROUNDING && predict_false(ix == asuint64_safe(1.0))) return 0; r = x - 1.0; #if __FP_FAST_FMA @@ -47,7 +47,7 @@ double log2(double x) lo = r * InvLn2lo + __builtin_fma(r, InvLn2hi, -hi); #else double_t rhi, rlo; - rhi = asdouble(asuint64(r) & -1ULL << 32); + rhi = asdouble_safe(asuint64_safe(r) & -1ULL << 32); rlo = r - rhi; hi = rhi * InvLn2hi; lo = rlo * InvLn2hi + r * InvLn2lo; @@ -67,12 +67,12 @@ double log2(double x) /* x < 0x1p-1022 or inf or nan. */ if (ix * 2 == 0) return __math_divzero(1); - if (ix == asuint64(INFINITY)) /* log(inf) == inf. */ + if (ix == asuint64_safe(INFINITY)) /* log(inf) == inf. */ return x; if ((top & 0x8000) || (top & 0x7ff0) == 0x7ff0) return __math_invalid(x); /* x is subnormal, normalize it. */ - ix = asuint64(x * 0x1p52); + ix = asuint64_safe(x * 0x1p52); ix -= 52ULL << 52; } @@ -85,7 +85,7 @@ double log2(double x) iz = ix - (tmp & 0xfffULL << 52); invc = T[i].invc; logc = T[i].logc; - z = asdouble(iz); + z = asdouble_safe(iz); kd = (double_t)k; /* log2(x) = log2(z/c) + log2(c) + k. */ @@ -99,7 +99,7 @@ double log2(double x) double_t rhi, rlo; /* rounding error: 0x1p-55/N + 0x1p-65. */ r = (z - T2[i].chi - T2[i].clo) * invc; - rhi = asdouble(asuint64(r) & -1ULL << 32); + rhi = asdouble_safe(asuint64_safe(r) & -1ULL << 32); rlo = r - rhi; t1 = rhi * InvLn2hi; t2 = rlo * InvLn2hi + r * InvLn2lo;