summaryrefslogtreecommitdiff
path: root/sysdeps/ieee754/dbl-64/s_atan.c
diff options
context:
space:
mode:
Diffstat (limited to 'sysdeps/ieee754/dbl-64/s_atan.c')
-rw-r--r--sysdeps/ieee754/dbl-64/s_atan.c115
1 files changed, 59 insertions, 56 deletions
diff --git a/sysdeps/ieee754/dbl-64/s_atan.c b/sysdeps/ieee754/dbl-64/s_atan.c
index 556f5b216d..b948f503a6 100644
--- a/sysdeps/ieee754/dbl-64/s_atan.c
+++ b/sysdeps/ieee754/dbl-64/s_atan.c
@@ -1,7 +1,7 @@
/*
* IBM Accurate Mathematical Library
* written by International Business Machines Corp.
- * Copyright (C) 2001 Free Software Foundation
+ * Copyright (C) 2001, 2011 Free Software Foundation
*
* This program is free software; you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as published by
@@ -52,8 +52,11 @@ double __signArctan(double,double);
double atan(double x) {
- double cor,s1,ss1,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10,u,u2,u3,
- v,vv,w,ww,y,yy,z,zz;
+ double cor,s1,ss1,s2,ss2,t1,t2,t3,t7,t8,t9,t10,u,u2,u3,
+ v,vv,w,ww,y,yy,z,zz;
+#ifndef DLA_FMA
+ double t4,t5,t6;
+#endif
#if 0
double y1,y2;
#endif
@@ -78,44 +81,44 @@ double atan(double x) {
if (u<C) {
if (u<B) {
if (u<A) { /* u < A */
- return x; }
+ return x; }
else { /* A <= u < B */
- v=x*x; yy=x*v*(d3.d+v*(d5.d+v*(d7.d+v*(d9.d+v*(d11.d+v*d13.d)))));
- if ((y=x+(yy-U1*x)) == x+(yy+U1*x)) return y;
-
- EMULV(x,x,v,vv,t1,t2,t3,t4,t5) /* v+vv=x^2 */
- s1=v*(f11.d+v*(f13.d+v*(f15.d+v*(f17.d+v*f19.d))));
- ADD2(f9.d,ff9.d,s1,ZERO,s2,ss2,t1,t2)
- MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
- ADD2(f7.d,ff7.d,s1,ss1,s2,ss2,t1,t2)
- MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
- ADD2(f5.d,ff5.d,s1,ss1,s2,ss2,t1,t2)
- MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
- ADD2(f3.d,ff3.d,s1,ss1,s2,ss2,t1,t2)
- MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
- MUL2(x,ZERO,s1,ss1,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8)
- ADD2(x,ZERO,s2,ss2,s1,ss1,t1,t2)
- if ((y=s1+(ss1-U5*s1)) == s1+(ss1+U5*s1)) return y;
-
- return atanMp(x,pr);
+ v=x*x; yy=x*v*(d3.d+v*(d5.d+v*(d7.d+v*(d9.d+v*(d11.d+v*d13.d)))));
+ if ((y=x+(yy-U1*x)) == x+(yy+U1*x)) return y;
+
+ EMULV(x,x,v,vv,t1,t2,t3,t4,t5) /* v+vv=x^2 */
+ s1=v*(f11.d+v*(f13.d+v*(f15.d+v*(f17.d+v*f19.d))));
+ ADD2(f9.d,ff9.d,s1,ZERO,s2,ss2,t1,t2)
+ MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
+ ADD2(f7.d,ff7.d,s1,ss1,s2,ss2,t1,t2)
+ MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
+ ADD2(f5.d,ff5.d,s1,ss1,s2,ss2,t1,t2)
+ MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
+ ADD2(f3.d,ff3.d,s1,ss1,s2,ss2,t1,t2)
+ MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
+ MUL2(x,ZERO,s1,ss1,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8)
+ ADD2(x,ZERO,s2,ss2,s1,ss1,t1,t2)
+ if ((y=s1+(ss1-U5*s1)) == s1+(ss1+U5*s1)) return y;
+
+ return atanMp(x,pr);
} }
else { /* B <= u < C */
i=(TWO52+TWO8*u)-TWO52; i-=16;
z=u-cij[i][0].d;
yy=z*(cij[i][2].d+z*(cij[i][3].d+z*(cij[i][4].d+
- z*(cij[i][5].d+z* cij[i][6].d))));
+ z*(cij[i][5].d+z* cij[i][6].d))));
t1=cij[i][1].d;
if (i<112) {
- if (i<48) u2=U21; /* u < 1/4 */
- else u2=U22; } /* 1/4 <= u < 1/2 */
+ if (i<48) u2=U21; /* u < 1/4 */
+ else u2=U22; } /* 1/4 <= u < 1/2 */
else {
- if (i<176) u2=U23; /* 1/2 <= u < 3/4 */
- else u2=U24; } /* 3/4 <= u <= 1 */
+ if (i<176) u2=U23; /* 1/2 <= u < 3/4 */
+ else u2=U24; } /* 3/4 <= u <= 1 */
if ((y=t1+(yy-u2*t1)) == t1+(yy+u2*t1)) return __signArctan(x,y);
z=u-hij[i][0].d;
s1=z*(hij[i][11].d+z*(hij[i][12].d+z*(hij[i][13].d+
- z*(hij[i][14].d+z* hij[i][15].d))));
+ z*(hij[i][14].d+z* hij[i][15].d))));
ADD2(hij[i][9].d,hij[i][10].d,s1,ZERO,s2,ss2,t1,t2)
MUL2(z,ZERO,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
ADD2(hij[i][7].d,hij[i][8].d,s1,ss1,s2,ss2,t1,t2)
@@ -138,7 +141,7 @@ double atan(double x) {
i=(TWO52+TWO8*w)-TWO52; i-=16;
z=(w-cij[i][0].d)+ww;
yy=HPI1-z*(cij[i][2].d+z*(cij[i][3].d+z*(cij[i][4].d+
- z*(cij[i][5].d+z* cij[i][6].d))));
+ z*(cij[i][5].d+z* cij[i][6].d))));
t1=HPI-cij[i][1].d;
if (i<112) u3=U31; /* w < 1/2 */
else u3=U32; /* w >= 1/2 */
@@ -148,7 +151,7 @@ double atan(double x) {
t1=w-hij[i][0].d;
EADD(t1,ww,z,zz)
s1=z*(hij[i][11].d+z*(hij[i][12].d+z*(hij[i][13].d+
- z*(hij[i][14].d+z* hij[i][15].d))));
+ z*(hij[i][14].d+z* hij[i][15].d))));
ADD2(hij[i][9].d,hij[i][10].d,s1,ZERO,s2,ss2,t1,t2)
MUL2(z,zz,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
ADD2(hij[i][7].d,hij[i][8].d,s1,ss1,s2,ss2,t1,t2)
@@ -165,36 +168,36 @@ double atan(double x) {
}
else {
if (u<E) { /* D <= u < E */
- w=ONE/u; v=w*w;
- EMULV(w,u,t1,t2,t3,t4,t5,t6,t7)
- yy=w*v*(d3.d+v*(d5.d+v*(d7.d+v*(d9.d+v*(d11.d+v*d13.d)))));
- ww=w*((ONE-t1)-t2);
- ESUB(HPI,w,t3,cor)
- yy=((HPI1+cor)-ww)-yy;
- if ((y=t3+(yy-U4)) == t3+(yy+U4)) return __signArctan(x,y);
-
- DIV2(ONE,ZERO,u,ZERO,w,ww,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
- MUL2(w,ww,w,ww,v,vv,t1,t2,t3,t4,t5,t6,t7,t8)
- s1=v*(f11.d+v*(f13.d+v*(f15.d+v*(f17.d+v*f19.d))));
- ADD2(f9.d,ff9.d,s1,ZERO,s2,ss2,t1,t2)
- MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
- ADD2(f7.d,ff7.d,s1,ss1,s2,ss2,t1,t2)
- MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
- ADD2(f5.d,ff5.d,s1,ss1,s2,ss2,t1,t2)
- MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
- ADD2(f3.d,ff3.d,s1,ss1,s2,ss2,t1,t2)
- MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
- MUL2(w,ww,s1,ss1,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8)
- ADD2(w,ww,s2,ss2,s1,ss1,t1,t2)
- SUB2(HPI,HPI1,s1,ss1,s2,ss2,t1,t2)
- if ((y=s2+(ss2-U8)) == s2+(ss2+U8)) return __signArctan(x,y);
+ w=ONE/u; v=w*w;
+ EMULV(w,u,t1,t2,t3,t4,t5,t6,t7)
+ yy=w*v*(d3.d+v*(d5.d+v*(d7.d+v*(d9.d+v*(d11.d+v*d13.d)))));
+ ww=w*((ONE-t1)-t2);
+ ESUB(HPI,w,t3,cor)
+ yy=((HPI1+cor)-ww)-yy;
+ if ((y=t3+(yy-U4)) == t3+(yy+U4)) return __signArctan(x,y);
+
+ DIV2(ONE,ZERO,u,ZERO,w,ww,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
+ MUL2(w,ww,w,ww,v,vv,t1,t2,t3,t4,t5,t6,t7,t8)
+ s1=v*(f11.d+v*(f13.d+v*(f15.d+v*(f17.d+v*f19.d))));
+ ADD2(f9.d,ff9.d,s1,ZERO,s2,ss2,t1,t2)
+ MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
+ ADD2(f7.d,ff7.d,s1,ss1,s2,ss2,t1,t2)
+ MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
+ ADD2(f5.d,ff5.d,s1,ss1,s2,ss2,t1,t2)
+ MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
+ ADD2(f3.d,ff3.d,s1,ss1,s2,ss2,t1,t2)
+ MUL2(v,vv,s2,ss2,s1,ss1,t1,t2,t3,t4,t5,t6,t7,t8)
+ MUL2(w,ww,s1,ss1,s2,ss2,t1,t2,t3,t4,t5,t6,t7,t8)
+ ADD2(w,ww,s2,ss2,s1,ss1,t1,t2)
+ SUB2(HPI,HPI1,s1,ss1,s2,ss2,t1,t2)
+ if ((y=s2+(ss2-U8)) == s2+(ss2+U8)) return __signArctan(x,y);
return atanMp(x,pr);
}
else {
- /* u >= E */
- if (x>0) return HPI;
- else return MHPI; }
+ /* u >= E */
+ if (x>0) return HPI;
+ else return MHPI; }
}
}