summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorrander <rander.wang@intel.com>2017-03-31 14:20:32 +0800
committerYang Rong <rong.r.yang@intel.com>2017-04-17 16:08:49 +0800
commitd3f8a429505e018f4b4c7cf1fd7eba491657aaf9 (patch)
tree1ca0ccbb6ec7edfed6bd1a88664200e39f1c5d63
parente3b4bd1874376bd9cf58b3615b526f6841deb9ba (diff)
downloadbeignet-d3f8a429505e018f4b4c7cf1fd7eba491657aaf9.tar.gz
backend: add double version of hypot
cp from fdlibm and pass the cft after refined Signed-off-by: rander <rander.wang@intel.com> Tested-by: Yang Rong <rong.r.yang@intel.com>
-rw-r--r--backend/src/libocl/tmpl/ocl_math_common.tmpl.cl65
-rw-r--r--backend/src/libocl/tmpl/ocl_math_common.tmpl.h1
2 files changed, 66 insertions, 0 deletions
diff --git a/backend/src/libocl/tmpl/ocl_math_common.tmpl.cl b/backend/src/libocl/tmpl/ocl_math_common.tmpl.cl
index 0cf58005..4dd2caee 100644
--- a/backend/src/libocl/tmpl/ocl_math_common.tmpl.cl
+++ b/backend/src/libocl/tmpl/ocl_math_common.tmpl.cl
@@ -1444,6 +1444,71 @@ OVERLOADABLE double floor(double x)
}
}
+OVERLOADABLE double hypot(double x, double y)
+{
+ double a=x,b=y,t1,t2,y1,y2,w;
+ int j,k,ha,hb;
+
+ ha = __HI(x)&0x7fffffff; /* high word of x */
+ hb = __HI(y)&0x7fffffff; /* high word of y */
+ if(hb > ha) {a=y;b=x;j=ha; ha=hb;hb=j;} else {a=x;b=y;}
+ __setHigh(&a, ha); /* a <- |a| */
+ __setHigh(&b, hb); /* b <- |b| */
+ if((ha-hb)>0x3c00000) {return a+b;} /* x/y > 2**60 */
+ k=0;
+ if(ha > 0x5f300000) { /* a>2**500 */
+ if(ha >= 0x7ff00000) { /* Inf or NaN */
+ w = a+b; /* for sNaN */
+ if(((ha&0xfffff)|__LO(a))==0) w = a;
+ if(((hb^0x7ff00000)|__LO(b))==0) w = b;
+ return w;
+ }
+ /* scale a and b by 2**-600 */
+ ha -= 0x25800000; hb -= 0x25800000; k += 600;
+ __setHigh(&a, ha);
+ __setHigh(&b, hb);
+ }
+ if(hb < 0x20b00000) { /* b < 2**-500 */
+ if(hb <= 0x000fffff) { /* subnormal b or 0 */
+ if((hb|(__LO(b)))==0) return a;
+ t1=0;
+ __setHigh(&t1, 0x7fd00000); /* t1=2^1022 */
+ b *= t1;
+ a *= t1;
+ k -= 1022;
+ } else { /* scale a and b by 2^600 */
+ ha += 0x25800000; /* a *= 2^600 */
+ hb += 0x25800000; /* b *= 2^600 */
+ k -= 600;
+ __setHigh(&a, ha);
+ __setHigh(&b, hb);
+ }
+ }
+ /* medium size a and b */
+ w = a-b;
+ if (w>b) {
+ t1 = 0;
+ __setHigh(&t1, ha);
+ t2 = a-t1;
+ w = sqrt(t1*t1-(b*(-b)-t2*(a+t1)));
+ } else {
+ a = a+a;
+ y1 = 0;
+ __setHigh(&y1, hb);
+ y2 = b - y1;
+ t1 = 0;
+ __setHigh(&t1, ha+0x00100000);
+ t2 = a - t1;
+ w = sqrt(t1*y1-(w*(-w)-(t1*y2+t2*b)));
+ }
+ if(k!=0) {
+ t1 = 1.0;
+ __setHigh(&t1, __HI(t1) + (k<<20));
+ return t1*w;
+ } else return w;
+
+}
+
OVERLOADABLE double log(double x)
{
double ln2_hi = 6.93147180369123816490e-01, /* 3fe62e42 fee00000 */
diff --git a/backend/src/libocl/tmpl/ocl_math_common.tmpl.h b/backend/src/libocl/tmpl/ocl_math_common.tmpl.h
index cc94d4e4..a538b3c1 100644
--- a/backend/src/libocl/tmpl/ocl_math_common.tmpl.h
+++ b/backend/src/libocl/tmpl/ocl_math_common.tmpl.h
@@ -47,6 +47,7 @@ OVERLOADABLE double floor(double x);
OVERLOADABLE double fmax(double a, double b);
OVERLOADABLE double fmin(double a, double b);
OVERLOADABLE double fmod (double x, double y);
+OVERLOADABLE double hypot(double x, double y);
OVERLOADABLE double ldexp(double x, int n);
OVERLOADABLE double log(double x);
OVERLOADABLE double log2(double x);