RFR: JDK-8302028: Port fdlibm atan2 to Java
Joe Darcy
darcy at openjdk.org
Thu Feb 16 19:38:04 UTC 2023
On Thu, 16 Feb 2023 19:30:06 GMT, Joe Darcy <darcy at openjdk.org> wrote:
> Working down the porting list, next stop, atan2.
Diffs of the various ports, starting with the original C vs the transliteration port:
$ diff -w Atan2.c Atan2.translit.java
1c1,4
< /* __ieee754_atan2(y,x)
---
> /**
> * Returns the angle theta from the conversion of rectangular
> * coordinates (x, y) to polar coordinates (r, theta).
> *
27,34c30,31
<
< #include "fdlibm.h"
<
< #ifdef __STDC__
< static const double
< #else
< static double
< #endif
---
> static class Atan2 {
> private static final double
42,48c39
< #ifdef __STDC__
< double __ieee754_atan2(double y, double x)
< #else
< double __ieee754_atan2(y,x)
< double y,x;
< #endif
< {
---
> static double compute(double y, double x) {
51c42
< unsigned lx,ly;
---
> /*unsigned*/ int lx,ly;
100c91
< else z=atan(fabs(y/x)); /* safe to do y/x */
---
> else z=atan(Math.abs(y/x)); /* safe to do y/x */
103c94,96
< case 1: __HI(z) ^= 0x80000000;
---
> case 1:
> // original:__HI(z) ^= 0x80000000;
> z = __HI(z, __HI(x) ^ 0x80000000);
109a103
> }
And the transliteration vs the more idiomatic port:
$ diff -w Atan2.translit.java Atan2.fdlibm.java
30a31,32
> private Atan2() {throw new UnsupportedOperationException();}
>
33,37c35,37
< zero = 0.0,
< pi_o_4 = 7.8539816339744827900E-01, /* 0x3FE921FB, 0x54442D18 */
< pi_o_2 = 1.5707963267948965580E+00, /* 0x3FF921FB, 0x54442D18 */
< pi = 3.1415926535897931160E+00, /* 0x400921FB, 0x54442D18 */
< pi_lo = 1.2246467991473531772E-16; /* 0x3CA1A626, 0x33145C07 */
---
> pi_o_4 = 0x1.921fb54442d18p-1, // 7.8539816339744827900E-01
> pi_o_2 = 0x1.921fb54442d18p0, // 1.5707963267948965580E+00
> pi_lo = 0x1.1a62633145c07p-53; // 1.2246467991473531772E-16
44c44,45
< hx = __HI(x); ix = hx&0x7fffffff;
---
> hx = __HI(x);
> ix = hx & 0x7fff_ffff;
46c47,48
< hy = __HI(y); iy = hy&0x7fffffff;
---
> hy = __HI(y);
> iy = hy&0x7fff_ffff;
48,49c50,51
< if(((ix|((lx|-lx)>>31))>0x7ff00000)||
< ((iy|((ly|-ly)>>31))>0x7ff00000)) /* x or y is NaN */
---
> if (((ix | ((lx | -lx) >> 31)) > 0x7ff0_0000)||
> ((iy |((ly | - ly) >> 31)) > 0x7ff0_0000)) // x or y is NaN
51,52c53,55
< if(((hx-0x3ff00000)|lx)==0) return atan(y); /* x=1.0 */
< m = ((hy>>31)&1)|((hx>>30)&2); /* 2*sign(x)+sign(y) */
---
> if (((hx - 0x3ff0_0000) | lx) == 0) // x = 1.0
> return StrictMath.atan(y);
> m = ((hy >> 31) & 1)|((hx >> 30) & 2); // 2*sign(x) + sign(y)
54c57
< /* when y = 0 */
---
> // when y = 0
56,73c59,67
< switch(m) {
< case 0:
< case 1: return y; /* atan(+-0,+anything)=+-0 */
< case 2: return pi+tiny;/* atan(+0,-anything) = pi */
< case 3: return -pi-tiny;/* atan(-0,-anything) =-pi */
< }
< }
< /* when x = 0 */
< if((ix|lx)==0) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
<
< /* when x is INF */
< if(ix==0x7ff00000) {
< if(iy==0x7ff00000) {
< switch(m) {
< case 0: return pi_o_4+tiny;/* atan(+INF,+INF) */
< case 1: return -pi_o_4-tiny;/* atan(-INF,+INF) */
< case 2: return 3.0*pi_o_4+tiny;/*atan(+INF,-INF)*/
< case 3: return -3.0*pi_o_4-tiny;/*atan(-INF,-INF)*/
---
> return switch(m) {
> case 0, 1 -> y; // atan(+/-0, +anything) = +/-0
> case 2 -> Math.PI + tiny; // atan(+0, -anything) = pi
> default -> -Math.PI - tiny; // atan(-0, -anything) = -pi
> };
> }
> // when x = 0
> if ((ix | lx) == 0) {
> return (hy < 0)? -pi_o_2 - tiny : pi_o_2 + tiny;
74a69,78
>
> // when x is INF
> if (ix == 0x7ff0_0000) {
> if (iy == 0x7ff0_0000) {
> return switch(m) {
> case 0 -> pi_o_4 + tiny; // atan(+INF, +INF)
> case 1 -> -pi_o_4 - tiny; // atan(-INF, +INF)
> case 2 -> 3.0*pi_o_4 + tiny; // atan(+INF, -INF)
> default -> -3.0*pi_o_4 - tiny; // atan(-INF, -INF)
> };
76,80c80,85
< switch(m) {
< case 0: return zero ; /* atan(+...,+INF) */
< case 1: return -1.0*zero ; /* atan(-...,+INF) */
< case 2: return pi+tiny ; /* atan(+...,-INF) */
< case 3: return -pi-tiny ; /* atan(-...,-INF) */
---
> return switch(m) {
> case 0 -> 0.0; // atan(+..., +INF)
> case 1 -> -0.0; // atan(-..., +INF)
> case 2 -> Math.PI + tiny; // atan(+..., -INF)
> default -> -Math.PI - tiny; // atan(-..., -INF)
> };
82a88,90
> // when y is INF
> if (iy == 0x7ff0_0000) {
> return (hy < 0)? -pi_o_2 - tiny : pi_o_2 + tiny;
84,85d91
< /* when y is INF */
< if(iy==0x7ff00000) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
87c93
< /* compute y/x */
---
> // compute y/x
89,101c95,107
< if(k > 60) z=pi_o_2+0.5*pi_lo; /* |y/x| > 2**60 */
< else if(hx<0&&k<-60) z=0.0; /* |y|/x < -2**60 */
< else z=atan(Math.abs(y/x)); /* safe to do y/x */
< switch (m) {
< case 0: return z ; /* atan(+,+) */
< case 1:
< // original:__HI(z) ^= 0x80000000;
< z = __HI(z, __HI(x) ^ 0x80000000);
< return z ; /* atan(-,+) */
< case 2: return pi-(z-pi_lo);/* atan(+,-) */
< default: /* case 3 */
< return (z-pi_lo)-pi;/* atan(-,-) */
< }
---
> if (k > 60) { // |y/x| > 2**60
> z = pi_o_2+0.5*pi_lo;
> } else if (hx < 0 && k < -60) { // |y|/x < -2**60
> z = 0.0;
> } else { // safe to do y/x
> z = StrictMath.atan(Math.abs(y/x));
> }
> return switch (m) {
> case 0 -> z; // atan(+, +)
> case 1 -> -z; // atan(-, +)
> case 2 -> Math.PI - (z - pi_lo); // atan(+, -)
> default -> (z - pi_lo) - Math.PI; // atan(-, -), case 3
> };
A few notes: I decided to use expression switches, but that is not necessary. Also, in the final switch I used "-z" as a more direct idiom to negate z than the form in the original fdlibm.
-------------
PR: https://git.openjdk.org/jdk/pull/12608
More information about the core-libs-dev
mailing list