1 files changed, 51 insertions, 0 deletions
diff --git a/pl/math/sv_hypot_1u5.c b/pl/math/sv_hypot_1u5.c
new file mode 100644
index 000000000000..cf1590e4b9ab
--- /dev/null
+++ b/pl/math/sv_hypot_1u5.c
@@ -0,0 +1,51 @@
+/*
+ * Double-precision SVE hypot(x) function.
+ *
+ * Copyright (c) 2023, Arm Limited.
+ * SPDX-License-Identifier: MIT OR Apache-2.0 WITH LLVM-exception
+ */
+
+#include "sv_math.h"
+#include "pl_sig.h"
+#include "pl_test.h"
+
+static const struct data
+{
+  uint64_t tiny_bound, thres;
+} data = {
+  .tiny_bound = 0x0c80000000000000, /* asuint (0x1p-102).  */
+  .thres = 0x7300000000000000,	    /* asuint (inf) - tiny_bound.  */
+};
+
+static svfloat64_t NOINLINE
+special_case (svfloat64_t sqsum, svfloat64_t x, svfloat64_t y, svbool_t pg,
+	      svbool_t special)
+{
+  return sv_call2_f64 (hypot, x, y, svsqrt_x (pg, sqsum), special);
+}
+
+/* SVE implementation of double-precision hypot.
+   Maximum error observed is 1.21 ULP:
+   _ZGVsMxvv_hypot (-0x1.6a22d0412cdd3p+352, 0x1.d3d89bd66fb1ap+330)
+    got 0x1.6a22d0412cfp+352
+   want 0x1.6a22d0412cf01p+352.  */
+svfloat64_t SV_NAME_D2 (hypot) (svfloat64_t x, svfloat64_t y, svbool_t pg)
+{
+  const struct data *d = ptr_barrier (&data);
+
+  svfloat64_t sqsum = svmla_x (pg, svmul_x (pg, x, x), y, y);
+
+  svbool_t special = svcmpge (
+      pg, svsub_x (pg, svreinterpret_u64 (sqsum), d->tiny_bound), d->thres);
+
+  if (unlikely (svptest_any (pg, special)))
+    return special_case (sqsum, x, y, pg, special);
+  return svsqrt_x (pg, sqsum);
+}
+
+PL_SIG (SV, D, 2, hypot, -10.0, 10.0)
+PL_TEST_ULP (SV_NAME_D2 (hypot), 0.71)
+PL_TEST_INTERVAL2 (SV_NAME_D2 (hypot), 0, inf, 0, inf, 10000)
+PL_TEST_INTERVAL2 (SV_NAME_D2 (hypot), 0, inf, -0, -inf, 10000)
+PL_TEST_INTERVAL2 (SV_NAME_D2 (hypot), -0, -inf, 0, inf, 10000)
+PL_TEST_INTERVAL2 (SV_NAME_D2 (hypot), -0, -inf, -0, -inf, 10000)