summaryrefslogtreecommitdiffstats
path: root/libc-top-half/musl/src/math/hypotl.c
diff options
context:
space:
mode:
authorDaniel Baumann <daniel.baumann@progress-linux.org>2024-04-17 13:54:38 +0000
committerDaniel Baumann <daniel.baumann@progress-linux.org>2024-04-17 13:54:38 +0000
commit8c1ab65c0f548d20b7f177bdb736daaf603340e1 (patch)
treedf55b7e75bf43f2bf500845b105afe3ac3a5157e /libc-top-half/musl/src/math/hypotl.c
parentInitial commit. (diff)
downloadwasi-libc-8c1ab65c0f548d20b7f177bdb736daaf603340e1.tar.xz
wasi-libc-8c1ab65c0f548d20b7f177bdb736daaf603340e1.zip
Adding upstream version 0.0~git20221206.8b7148f.upstream/0.0_git20221206.8b7148f
Signed-off-by: Daniel Baumann <daniel.baumann@progress-linux.org>
Diffstat (limited to 'libc-top-half/musl/src/math/hypotl.c')
-rw-r--r--libc-top-half/musl/src/math/hypotl.c66
1 files changed, 66 insertions, 0 deletions
diff --git a/libc-top-half/musl/src/math/hypotl.c b/libc-top-half/musl/src/math/hypotl.c
new file mode 100644
index 0000000..479aa92
--- /dev/null
+++ b/libc-top-half/musl/src/math/hypotl.c
@@ -0,0 +1,66 @@
+#include "libm.h"
+
+#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
+long double hypotl(long double x, long double y)
+{
+ return hypot(x, y);
+}
+#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
+#if LDBL_MANT_DIG == 64
+#define SPLIT (0x1p32L+1)
+#elif LDBL_MANT_DIG == 113
+#define SPLIT (0x1p57L+1)
+#endif
+
+static void sq(long double *hi, long double *lo, long double x)
+{
+ long double xh, xl, xc;
+ xc = x*SPLIT;
+ xh = x - xc + xc;
+ xl = x - xh;
+ *hi = x*x;
+ *lo = xh*xh - *hi + 2*xh*xl + xl*xl;
+}
+
+long double hypotl(long double x, long double y)
+{
+ union ldshape ux = {x}, uy = {y};
+ int ex, ey;
+ long double hx, lx, hy, ly, z;
+
+ ux.i.se &= 0x7fff;
+ uy.i.se &= 0x7fff;
+ if (ux.i.se < uy.i.se) {
+ ex = uy.i.se;
+ ey = ux.i.se;
+ x = uy.f;
+ y = ux.f;
+ } else {
+ ex = ux.i.se;
+ ey = uy.i.se;
+ x = ux.f;
+ y = uy.f;
+ }
+
+ if (ex == 0x7fff && isinf(y))
+ return y;
+ if (ex == 0x7fff || y == 0)
+ return x;
+ if (ex - ey > LDBL_MANT_DIG)
+ return x + y;
+
+ z = 1;
+ if (ex > 0x3fff+8000) {
+ z = 0x1p10000L;
+ x *= 0x1p-10000L;
+ y *= 0x1p-10000L;
+ } else if (ey < 0x3fff-8000) {
+ z = 0x1p-10000L;
+ x *= 0x1p10000L;
+ y *= 0x1p10000L;
+ }
+ sq(&hx, &lx, x);
+ sq(&hy, &ly, y);
+ return z*sqrtl(ly+lx+hy+hx);
+}
+#endif