1*412f47f9SXin Li /*
2*412f47f9SXin Li * Double-precision SVE hypot(x) function.
3*412f47f9SXin Li *
4*412f47f9SXin Li * Copyright (c) 2023, Arm Limited.
5*412f47f9SXin Li * SPDX-License-Identifier: MIT OR Apache-2.0 WITH LLVM-exception
6*412f47f9SXin Li */
7*412f47f9SXin Li
8*412f47f9SXin Li #include "sv_math.h"
9*412f47f9SXin Li #include "pl_sig.h"
10*412f47f9SXin Li #include "pl_test.h"
11*412f47f9SXin Li
12*412f47f9SXin Li static const struct data
13*412f47f9SXin Li {
14*412f47f9SXin Li uint64_t tiny_bound, thres;
15*412f47f9SXin Li } data = {
16*412f47f9SXin Li .tiny_bound = 0x0c80000000000000, /* asuint (0x1p-102). */
17*412f47f9SXin Li .thres = 0x7300000000000000, /* asuint (inf) - tiny_bound. */
18*412f47f9SXin Li };
19*412f47f9SXin Li
20*412f47f9SXin Li static svfloat64_t NOINLINE
special_case(svfloat64_t sqsum,svfloat64_t x,svfloat64_t y,svbool_t pg,svbool_t special)21*412f47f9SXin Li special_case (svfloat64_t sqsum, svfloat64_t x, svfloat64_t y, svbool_t pg,
22*412f47f9SXin Li svbool_t special)
23*412f47f9SXin Li {
24*412f47f9SXin Li return sv_call2_f64 (hypot, x, y, svsqrt_x (pg, sqsum), special);
25*412f47f9SXin Li }
26*412f47f9SXin Li
27*412f47f9SXin Li /* SVE implementation of double-precision hypot.
28*412f47f9SXin Li Maximum error observed is 1.21 ULP:
29*412f47f9SXin Li _ZGVsMxvv_hypot (-0x1.6a22d0412cdd3p+352, 0x1.d3d89bd66fb1ap+330)
30*412f47f9SXin Li got 0x1.6a22d0412cfp+352
31*412f47f9SXin Li want 0x1.6a22d0412cf01p+352. */
SV_NAME_D2(hypot)32*412f47f9SXin Li svfloat64_t SV_NAME_D2 (hypot) (svfloat64_t x, svfloat64_t y, svbool_t pg)
33*412f47f9SXin Li {
34*412f47f9SXin Li const struct data *d = ptr_barrier (&data);
35*412f47f9SXin Li
36*412f47f9SXin Li svfloat64_t sqsum = svmla_x (pg, svmul_x (pg, x, x), y, y);
37*412f47f9SXin Li
38*412f47f9SXin Li svbool_t special = svcmpge (
39*412f47f9SXin Li pg, svsub_x (pg, svreinterpret_u64 (sqsum), d->tiny_bound), d->thres);
40*412f47f9SXin Li
41*412f47f9SXin Li if (unlikely (svptest_any (pg, special)))
42*412f47f9SXin Li return special_case (sqsum, x, y, pg, special);
43*412f47f9SXin Li return svsqrt_x (pg, sqsum);
44*412f47f9SXin Li }
45*412f47f9SXin Li
46*412f47f9SXin Li PL_SIG (SV, D, 2, hypot, -10.0, 10.0)
47*412f47f9SXin Li PL_TEST_ULP (SV_NAME_D2 (hypot), 0.71)
48*412f47f9SXin Li PL_TEST_INTERVAL2 (SV_NAME_D2 (hypot), 0, inf, 0, inf, 10000)
49*412f47f9SXin Li PL_TEST_INTERVAL2 (SV_NAME_D2 (hypot), 0, inf, -0, -inf, 10000)
50*412f47f9SXin Li PL_TEST_INTERVAL2 (SV_NAME_D2 (hypot), -0, -inf, 0, inf, 10000)
51*412f47f9SXin Li PL_TEST_INTERVAL2 (SV_NAME_D2 (hypot), -0, -inf, -0, -inf, 10000)
52