[Libclc-dev] [PATCH 2/9] Implement acosh builtin
Tom Stellard
thomas.stellard at amd.com
Thu Mar 12 13:11:39 PDT 2015
This implementation was ported from the AMD builtin library
and has been tested with piglit, OpenCV, and the ocl conformance tests.
---
generic/include/clc/clc.h | 1 +
generic/include/clc/math/acosh.h | 24 +++++++
generic/include/clc/math/acosh.inc | 23 +++++++
generic/lib/SOURCES | 2 +
generic/lib/math/acosh.cl | 127 +++++++++++++++++++++++++++++++++++++
generic/lib/math/ep_log.cl | 94 +++++++++++++++++++++++++++
generic/lib/math/ep_log.h | 29 +++++++++
7 files changed, 300 insertions(+)
create mode 100644 generic/include/clc/math/acosh.h
create mode 100644 generic/include/clc/math/acosh.inc
create mode 100644 generic/lib/math/acosh.cl
create mode 100644 generic/lib/math/ep_log.cl
create mode 100644 generic/lib/math/ep_log.h
diff --git a/generic/include/clc/clc.h b/generic/include/clc/clc.h
index 9c09e67..269d30d 100644
--- a/generic/include/clc/clc.h
+++ b/generic/include/clc/clc.h
@@ -33,6 +33,7 @@
/* 6.11.2 Math Functions */
#include <clc/math/acos.h>
+#include <clc/math/acosh.h>
#include <clc/math/asin.h>
#include <clc/math/atan.h>
#include <clc/math/atan2.h>
diff --git a/generic/include/clc/math/acosh.h b/generic/include/clc/math/acosh.h
new file mode 100644
index 0000000..f1bb0a5
--- /dev/null
+++ b/generic/include/clc/math/acosh.h
@@ -0,0 +1,24 @@
+/*
+ * Copyright (c) 2014, 2015 Advanced Micro Devices, Inc.
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to deal
+ * in the Software without restriction, including without limitation the rights
+ * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included in
+ * all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+ * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+ * THE SOFTWARE.
+ */
+
+#define __CLC_BODY <clc/math/acosh.inc>
+#include <clc/math/gentype.inc>
diff --git a/generic/include/clc/math/acosh.inc b/generic/include/clc/math/acosh.inc
new file mode 100644
index 0000000..41ec429
--- /dev/null
+++ b/generic/include/clc/math/acosh.inc
@@ -0,0 +1,23 @@
+/*
+ * Copyright (c) 2014, 2015 Advanced Micro Devices, Inc.
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to deal
+ * in the Software without restriction, including without limitation the rights
+ * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included in
+ * all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+ * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+ * THE SOFTWARE.
+ */
+
+_CLC_OVERLOAD _CLC_DECL __CLC_GENTYPE acosh(__CLC_GENTYPE x);
diff --git a/generic/lib/SOURCES b/generic/lib/SOURCES
index c4df016..29c40d8 100644
--- a/generic/lib/SOURCES
+++ b/generic/lib/SOURCES
@@ -62,11 +62,13 @@ integer/sub_sat_if.ll
integer/sub_sat_impl.ll
integer/upsample.cl
math/acos.cl
+math/acosh.cl
math/asin.cl
math/atan.cl
math/atan2.cl
math/copysign.cl
math/cos.cl
+math/ep_log.cl
math/exp.cl
math/exp10.cl
math/half_rsqrt.cl
diff --git a/generic/lib/math/acosh.cl b/generic/lib/math/acosh.cl
new file mode 100644
index 0000000..cc10dd4
--- /dev/null
+++ b/generic/lib/math/acosh.cl
@@ -0,0 +1,127 @@
+/*
+ * Copyright (c) 2014 Advanced Micro Devices, Inc.
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to deal
+ * in the Software without restriction, including without limitation the rights
+ * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included in
+ * all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+ * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+ * THE SOFTWARE.
+ */
+
+#include <clc/clc.h>
+
+#include "ep_log.h"
+#include "math.h"
+#include "../clcmacro.h"
+
+_CLC_OVERLOAD _CLC_DEF float acosh(float x) {
+ uint ux = as_uint(x);
+
+ // Arguments greater than 1/sqrt(epsilon) in magnitude are
+ // approximated by acosh(x) = ln(2) + ln(x)
+ // For 2.0 <= x <= 1/sqrt(epsilon) the approximation is
+ // acosh(x) = ln(x + sqrt(x*x-1)) */
+ int high = ux > 0x46000000U;
+ int med = ux > 0x40000000U;
+
+ float w = x - 1.0f;
+ float s = w*w + 2.0f*w;
+ float t = x*x - 1.0f;
+ float r = sqrt(med ? t : s) + (med ? x : w);
+ float v = (high ? x : r) - (med ? 1.0f : 0.0f);
+ float z = log1p(v) + (high ? 0x1.62e430p-1f : 0.0f);
+
+ z = ux >= PINFBITPATT_SP32 ? x : z;
+ z = x < 1.0f ? as_float(QNANBITPATT_SP32) : z;
+
+ return z;
+}
+
+_CLC_UNARY_VECTORIZE(_CLC_OVERLOAD _CLC_DEF, float, acosh, float)
+
+#ifdef cl_khr_fp64
+#pragma OPENCL EXTENSION cl_khr_fp64 : enable
+
+_CLC_OVERLOAD _CLC_DEF double acosh(double x) {
+ const double recrteps = 0x1.6a09e667f3bcdp+26; // 1/sqrt(eps) = 9.49062656242515593767e+07
+ //log2_lead and log2_tail sum to an extra-precise version of log(2)
+ const double log2_lead = 0x1.62e42ep-1;
+ const double log2_tail = 0x1.efa39ef35793cp-25;
+
+ // Handle x >= 128 here
+ int xlarge = x > recrteps;
+ double r = x + sqrt(fma(x, x, -1.0));
+ r = xlarge ? x : r;
+
+ int xexp;
+ double r1, r2;
+ __clc_ep_log(r, &xexp, &r1, &r2);
+
+ double dxexp = xexp + xlarge;
+ r1 = fma(dxexp, log2_lead, r1);
+ r2 = fma(dxexp, log2_tail, r2);
+
+ double ret1 = r1 + r2;
+
+ // Handle 1 < x < 128 here
+ // We compute the value
+ // t = x - 1.0 + sqrt(2.0*(x - 1.0) + (x - 1.0)*(x - 1.0))
+ // using simulated quad precision.
+ double t = x - 1.0;
+ double u1 = t * 2.0;
+
+ // (t,0) * (t,0) -> (v1, v2)
+ double v1 = t * t;
+ double v2 = fma(t, t, -v1);
+
+ // (u1,0) + (v1,v2) -> (w1,w2)
+ r = u1 + v1;
+ double s = (((u1 - r) + v1) + v2);
+ double w1 = r + s;
+ double w2 = (r - w1) + s;
+
+ // sqrt(w1,w2) -> (u1,u2)
+ double p1 = sqrt(w1);
+ double a1 = p1*p1;
+ double a2 = fma(p1, p1, -a1);
+ double temp = (((w1 - a1) - a2) + w2);
+ double p2 = MATH_DIVIDE(temp * 0.5, p1);
+ u1 = p1 + p2;
+ double u2 = (p1 - u1) + p2;
+
+ // (u1,u2) + (t,0) -> (r1,r2)
+ r = u1 + t;
+ s = ((u1 - r) + t) + u2;
+ // r1 = r + s;
+ // r2 = (r - r1) + s;
+ // t = r1 + r2;
+ t = r + s;
+
+ // For arguments 1.13 <= x <= 1.5 the log1p function is good enough
+ double ret2 = log1p(t);
+
+ ulong ux = as_ulong(x);
+ double ret = x >= 128.0 ? ret1 : ret2;
+
+ ret = ux >= 0x7FF0000000000000 ? x : ret;
+ ret = x == 1.0 ? 0.0 : ret;
+ ret = (ux & SIGNBIT_DP64) != 0UL | x < 1.0 ? as_double(QNANBITPATT_DP64) : ret;
+
+ return ret;
+}
+
+_CLC_UNARY_VECTORIZE(_CLC_OVERLOAD _CLC_DEF, double, acosh, double)
+
+#endif
diff --git a/generic/lib/math/ep_log.cl b/generic/lib/math/ep_log.cl
new file mode 100644
index 0000000..3c2c62c
--- /dev/null
+++ b/generic/lib/math/ep_log.cl
@@ -0,0 +1,94 @@
+/*
+ * Copyright (c) 2014 Advanced Micro Devices, Inc.
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to deal
+ * in the Software without restriction, including without limitation the rights
+ * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included in
+ * all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+ * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+ * THE SOFTWARE.
+ */
+
+#ifdef cl_khr_fp64
+
+#include <clc/clc.h>
+#include "ep_log.h"
+#include "math.h"
+#include "tables.h"
+
+#pragma OPENCL EXTENSION cl_khr_fp64 : enable
+
+#define LN0 8.33333333333317923934e-02
+#define LN1 1.25000000037717509602e-02
+#define LN2 2.23213998791944806202e-03
+#define LN3 4.34887777707614552256e-04
+
+#define LF0 8.33333333333333593622e-02
+#define LF1 1.24999999978138668903e-02
+#define LF2 2.23219810758559851206e-03
+
+_CLC_DEF void __clc_ep_log(double x, int *xexp, double *r1, double *r2)
+{
+ // Computes natural log(x). Algorithm based on:
+ // Ping-Tak Peter Tang
+ // "Table-driven implementation of the logarithm function in IEEE
+ // floating-point arithmetic"
+ // ACM Transactions on Mathematical Software (TOMS)
+ // Volume 16, Issue 4 (December 1990)
+ int near_one = x >= 0x1.e0faap-1 & x <= 0x1.1082cp+0;
+
+ ulong ux = as_ulong(x);
+ ulong uxs = as_ulong(as_double(0x03d0000000000000UL | ux) - 0x1.0p-962);
+ int c = ux < IMPBIT_DP64;
+ ux = c ? uxs : ux;
+ int expadjust = c ? 60 : 0;
+
+ // Store the exponent of x in xexp and put f into the range [0.5,1)
+ int xexp1 = ((as_int2(ux).hi >> 20) & 0x7ff) - EXPBIAS_DP64 - expadjust;
+ double f = as_double(HALFEXPBITS_DP64 | (ux & MANTBITS_DP64));
+ *xexp = near_one ? 0 : xexp1;
+
+ double r = x - 1.0;
+ double u1 = MATH_DIVIDE(r, 2.0 + r);
+ double ru1 = -r * u1;
+ u1 = u1 + u1;
+
+ int index = as_int2(ux).hi >> 13;
+ index = ((0x80 | (index & 0x7e)) >> 1) + (index & 0x1);
+
+ double f1 = index * 0x1.0p-7;
+ double f2 = f - f1;
+ double u2 = MATH_DIVIDE(f2, fma(0.5, f2, f1));
+
+ double2 tv = USE_TABLE(ln_tbl, (index - 64));
+ double z1 = tv.s0;
+ double q = tv.s1;
+
+ z1 = near_one ? r : z1;
+ q = near_one ? 0.0 : q;
+ double u = near_one ? u1 : u2;
+ double v = u*u;
+
+ double cc = near_one ? ru1 : u2;
+
+ double z21 = fma(v, fma(v, fma(v, LN3, LN2), LN1), LN0);
+ double z22 = fma(v, fma(v, LF2, LF1), LF0);
+ double z2 = near_one ? z21 : z22;
+ z2 = fma(u*v, z2, cc) + q;
+
+ *r1 = z1;
+ *r2 = z2;
+}
+
+#endif
diff --git a/generic/lib/math/ep_log.h b/generic/lib/math/ep_log.h
new file mode 100644
index 0000000..7f99ac6
--- /dev/null
+++ b/generic/lib/math/ep_log.h
@@ -0,0 +1,29 @@
+/*
+ * Copyright (c) 2014 Advanced Micro Devices, Inc.
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a copy
+ * of this software and associated documentation files (the "Software"), to deal
+ * in the Software without restriction, including without limitation the rights
+ * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the Software is
+ * furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included in
+ * all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+ * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+ * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+ * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+ * THE SOFTWARE.
+ */
+
+#ifdef cl_khr_fp64
+
+#pragma OPENCL EXTENSION cl_khr_fp64 : enable
+
+_CLC_DECL void __clc_ep_log(double x, int *xexp, double *r1, double *r2);
+
+#endif
--
2.0.4
More information about the Libclc-dev
mailing list