test/java/lang/StrictMath/FdlibmTranslit.java

*** 71,80 ****
--- 71,141 ----
      public static double hypot(double x, double y) {
          return Hypot.compute(x, y);
      }
  
      /**
+      * cbrt(x)
+      * Return cube root of x
+      */
+     public static class Cbrt {
+         // unsigned
+         private static final int B1 = 715094163; /* B1 = (682-0.03306235651)*2**20 */
+         private static final int B2 = 696219795; /* B2 = (664-0.03306235651)*2**20 */
+ 
+         private static final double C =  5.42857142857142815906e-01; /* 19/35     = 0x3FE15F15, 0xF15F15F1 */
+         private static final double D = -7.05306122448979611050e-01; /* -864/1225 = 0xBFE691DE, 0x2532C834 */
+         private static final double E =  1.41428571428571436819e+00; /* 99/70     = 0x3FF6A0EA, 0x0EA0EA0F */
+         private static final double F =  1.60714285714285720630e+00; /* 45/28     = 0x3FF9B6DB, 0x6DB6DB6E */
+         private static final double G =  3.57142857142857150787e-01; /* 5/14      = 0x3FD6DB6D, 0xB6DB6DB7 */
+ 
+         public static strictfp double compute(double x) {
+             int     hx;
+             double  r, s, t=0.0, w;
+             int sign; // unsigned
+ 
+             hx = __HI(x);           // high word of x
+             sign = hx & 0x80000000;             // sign= sign(x)
+             hx  ^= sign;
+             if (hx >= 0x7ff00000)
+                 return (x+x); // cbrt(NaN,INF) is itself
+             if ((hx | __LO(x)) == 0)
+                 return(x);          // cbrt(0) is itself
+ 
+             x = __HI(x, hx);   // x <- |x|
+             // rough cbrt to 5 bits
+             if (hx < 0x00100000) {               // subnormal number
+                 t = __HI(t, 0x43500000);          // set t= 2**54
+                 t *= x;
+                 t = __HI(t, __HI(t)/3+B2);
+             } else {
+                 t = __HI(t, hx/3+B1);
+             }
+ 
+             // new cbrt to 23 bits, may be implemented in single precision
+             r = t * t/x;
+             s = C + r*t;
+             t *= G + F/(s + E + D/s);
+ 
+             // chopped to 20 bits and make it larger than cbrt(x)
+             t = __LO(t, 0);
+             t = __HI(t, __HI(t)+0x00000001);
+ 
+ 
+             // one step newton iteration to 53 bits with error less than 0.667 ulps
+             s = t * t;          // t*t is exact
+             r = x / s;
+             w = t + t;
+             r= (r - t)/(w + r);  // r-s is exact
+             t= t + t*r;
+ 
+             // retore the sign bit
+             t = __HI(t, __HI(t) | sign);
+             return(t);
+         }
+     }
+ 
+     /**
       * hypot(x,y)
       *
       * Method :
       *      If (assume round-to-nearest) z = x*x + y*y
       *      has error less than sqrt(2)/2 ulp, than