一、前言
关于apache的commons-math3开源包中的org.apache.commons.math3.util.FastMath数字处理类，对三角函数sinh、tanh、acosh、atanh等处理、对随机数random生成、指数函数exp、x的y次幂（次方）pow等。
二、源码说明
package org.apache.commons.math3.util;@b@@b@import java.io.PrintStream;@b@@b@public class FastMath@b@{@b@ public static final double PI = 3.141592653589793D;@b@ public static final double E = 2.718281828459045D;@b@ static final int EXP_INT_TABLE_MAX_INDEX = 750;@b@ static final int EXP_INT_TABLE_LEN = 1500;@b@ static final int LN_MANT_LEN = 1024;@b@ static final int EXP_FRAC_TABLE_LEN = 1025;@b@ private static final double LOG_MAX_VALUE = StrictMath.log(1.7976931348623157E+308D);@b@ private static final boolean RECOMPUTE_TABLES_AT_RUNTIME = 0;@b@ private static final double LN_2_A = 0.6931470632553101D;@b@ private static final double LN_2_B = 1.173046352508235E-007D;@b@ private static final double[][] LN_QUICK_COEF = { { 1.0D, 5.669184079525E-024D }, { -0.25D, -0.25D }, { 0.3333333134651184D, 1.986821492305628E-008D }, { -0.25D, -6.663542893624021E-014D }, { 0.199999988079071D, 1.192105680146323E-008D }, { -0.166666656732559D, -7.800414592973399E-009D }, { 0.142857134342194D, 5.650007086920087E-009D }, { -0.1250253021717072D, -7.44321345601866E-011D }, { 0.1111380755901337D, 9.219544613762692E-009D } };@b@ private static final double[][] LN_HI_PREC_COEF = { { 1.0D, -6.032174644509064E-023D }, { -0.25D, -0.25D }, { 0.3333333134651184D, 1.986816177772435E-008D }, { -0.2499999701976776D, -2.957007209750105E-008D }, { 0.1999995410442352D, 1.583099333206127E-010D }, { -0.1662487983703613D, -2.603382435519167E-008D } };@b@ private static final int SINE_TABLE_LEN = 14;@b@ private static final double[] SINE_TABLE_A = { 0.0D, 0.1246747374534607D, 0.2474039494991303D, 0.366272509098053D, 0.4794255495071411D, 0.5850973129272461D, 0.6816387176513672D, 0.7675435543060303D, 0.8414709568023682D, 0.9022675752639771D, 0.9489846229553223D, 0.9808930158615112D, 0.9974949359893799D, 0.9985313415527344D };@b@ private static final double[] SINE_TABLE_B = { 0.0D, -4.068233003401932E-009D, 9.755392680573412E-009D, 1.998799458285729E-008D, -1.090293811300796E-008D, -3.99867839389446E-008D, 4.23719669792332E-008D, -5.207000323380292E-008D, 2.800552834259E-008D, 1.883511811213715E-008D, -3.599736051276557E-009D, 4.116164446561962E-008D, 5.061467454812738E-008D, -1.012902791249686E-009D };@b@ private static final double[] COSINE_TABLE_A = { 1.0D, 0.9921976327896118D, 0.9689123630523682D, 0.9305076599121094D, 0.8775825500488281D, 0.8109631538391113D, 0.7316888570785523D, 0.6409968137741089D, 0.5403022766113281D, 0.4311765432357788D, 0.3153223395347595D, 0.194547712802887D, 0.07073719799518585D, -0.05417713522911072D };@b@ private static final double[] COSINE_TABLE_B = { 0.0D, 3.443971723674285E-008D, 5.865827662008209E-008D, -3.799979508385053E-008D, 1.184154459111628E-008D, -3.43338934259355E-008D, 1.179526864021679E-008D, 4.438921624363781E-008D, 2.925681159240093E-008D, -2.643711263204181E-008D, 2.286050914396312E-008D, -4.813899778443457E-009D, 3.672517058035558E-009D, 2.021743975633808E-010D };@b@ private static final double[] TANGENT_TABLE_A = { 0.0D, 0.1256551444530487D, 0.2553419470787048D, 0.3936265707015991D, 0.546302437782288D, 0.7214844226837158D, 0.9315965175628662D, 1.197421550750732D, 1.55740761756897D, 2.092571258544922D, 3.009569644927979D, 5.041914939880371D, 14.101419448852539D, -18.430862426757813D };@b@ private static final double[] TANGENT_TABLE_B = { 0.0D, -7.877917738262007E-009D, -2.585766856747989E-008D, 5.224033637135667E-009D, 5.206150291559893E-008D, 1.830718859967703E-008D, -5.761879374977071E-008D, 7.848361555046424E-008D, 1.070859325039445E-007D, 1.782725712942381E-008D, 2.893485277253286E-008D, 3.166009922273796E-007D, 4.983191803254889E-007D, -3.356118100840571E-007D };@b@ private static final long[] RECIP_2PI = { 2935890503282001226L, 9154082963658192752L, 3952090531849364496L, 9193070505571053912L, 7910884519577875640L, 113236205062349959L, 4577762542105553359L, -5034868814120038111L, 4208363204685324176L, 5648769086999809661L, 2819561105158720014L, -4035746434778044925L, -302932621132653753L, -2644281811660520851L, -3183605296591799669L, 6722166367014452318L, -3512299194304650054L, -7278142539171889152L };@b@ private static final long[] PI_O_4_BITS = { -3958705157555305932L, -4267615245585081135L };@b@ private static final double[] EIGHTHS = { 0.0D, 0.125D, 0.25D, 0.375D, 0.5D, 0.625D, 0.75D, 0.875D, 1.0D, 1.125D, 1.25D, 1.375D, 1.5D, 1.625D };@b@ private static final double[] CBRTTWO = { 0.629960524947437D, 0.7937005259840998D, 1.0D, 1.259921049894873D, 1.587401051968199D };@b@ private static final long HEX_40000000 = 1073741824L;@b@ private static final long MASK_30BITS = -1073741824L;@b@ private static final int MASK_NON_SIGN_INT = 2147483647;@b@ private static final long MASK_NON_SIGN_LONG = 9223372036854775807L;@b@ private static final double TWO_POWER_52 = 4503599627370496.0D;@b@ private static final double TWO_POWER_53 = 9007199254740992.0D;@b@ private static final double F_1_3 = 0.3333333333333333D;@b@ private static final double F_1_5 = 0.2D;@b@ private static final double F_1_7 = 0.1428571428571429D;@b@ private static final double F_1_9 = 0.111111111111111D;@b@ private static final double F_1_11 = 0.09090909090909091D;@b@ private static final double F_1_13 = 0.07692307692307693D;@b@ private static final double F_1_15 = 0.06666666666666667D;@b@ private static final double F_1_17 = 0.05882352941176471D;@b@ private static final double F_3_4 = 0.75D;@b@ private static final double F_15_16 = 0.9375D;@b@ private static final double F_13_14 = 0.928571428571429D;@b@ private static final double F_11_12 = 0.9166666666666666D;@b@ private static final double F_9_10 = 0.9D;@b@ private static final double F_7_8 = 0.875D;@b@ private static final double F_5_6 = 0.8333333333333334D;@b@ private static final double F_1_2 = 0.5D;@b@ private static final double F_1_4 = 0.25D;@b@@b@ private static double doubleHighPart(double d)@b@ {@b@ if ((d > -Precision.SAFE_MIN) && (d < Precision.SAFE_MIN))@b@ return d;@b@@b@ long xl = Double.doubleToRawLongBits(d);@b@ xl &= -1073741824L;@b@ return Double.longBitsToDouble(xl);@b@ }@b@@b@ public static double sqrt(double a)@b@ {@b@ return Math.sqrt(a);@b@ }@b@@b@ public static double cosh(double x)@b@ {@b@ double t;@b@ if (x != x) {@b@ return x;@b@ }@b@@b@ if (x > 20.0D) {@b@ if (x >= LOG_MAX_VALUE)@b@ {@b@ t = exp(0.5D * x);@b@ return (0.5D * t * t);@b@ }@b@ return (0.5D * exp(x));@b@ }@b@ if (x < -20.0D) {@b@ if (x <= -LOG_MAX_VALUE)@b@ {@b@ t = exp(-0.5D * x);@b@ return (0.5D * t * t);@b@ }@b@ return (0.5D * exp(-x));@b@ }@b@@b@ double[] hiPrec = new double[2];@b@ if (x < 0.0D)@b@ x = -x;@b@@b@ exp(x, 0.0D, hiPrec);@b@@b@ double ya = hiPrec[0] + hiPrec[1];@b@ double yb = -(ya - hiPrec[0] - hiPrec[1]);@b@@b@ double temp = ya * 1073741824.0D;@b@ double yaa = ya + temp - temp;@b@ double yab = ya - yaa;@b@@b@ double recip = 1.0D / ya;@b@ temp = recip * 1073741824.0D;@b@ double recipa = recip + temp - temp;@b@ double recipb = recip - recipa;@b@@b@ recipb += (1.0D - yaa * recipa - yaa * recipb - yab * recipa - yab * recipb) * recip;@b@@b@ recipb += -yb * recip * recip;@b@@b@ temp = ya + recipa;@b@ yb += -(temp - ya - recipa);@b@ ya = temp;@b@ temp = ya + recipb;@b@ yb += -(temp - ya - recipb);@b@ ya = temp;@b@@b@ double result = ya + yb;@b@ result *= 0.5D;@b@ return result;@b@ }@b@@b@ public static double sinh(double x)@b@ {@b@ double t;@b@ double result;@b@ double[] hiPrec;@b@ double ya;@b@ double yb;@b@ boolean negate = false;@b@ if (x != x) {@b@ return x;@b@ }@b@@b@ if (x > 20.0D) {@b@ if (x >= LOG_MAX_VALUE)@b@ {@b@ t = exp(0.5D * x);@b@ return (0.5D * t * t);@b@ }@b@ return (0.5D * exp(x));@b@ }@b@ if (x < -20.0D) {@b@ if (x <= -LOG_MAX_VALUE)@b@ {@b@ t = exp(-0.5D * x);@b@ return (-0.5D * t * t);@b@ }@b@ return (-0.5D * exp(-x));@b@ }@b@@b@ if (x == 0.0D) {@b@ return x;@b@ }@b@@b@ if (x < 0.0D) {@b@ x = -x;@b@ negate = true;@b@ }@b@@b@ if (x > 0.25D) {@b@ hiPrec = new double[2];@b@ exp(x, 0.0D, hiPrec);@b@@b@ ya = hiPrec[0] + hiPrec[1];@b@ yb = -(ya - hiPrec[0] - hiPrec[1]);@b@@b@ double temp = ya * 1073741824.0D;@b@ double yaa = ya + temp - temp;@b@ double yab = ya - yaa;@b@@b@ double recip = 1.0D / ya;@b@ temp = recip * 1073741824.0D;@b@ double recipa = recip + temp - temp;@b@ double recipb = recip - recipa;@b@@b@ recipb += (1.0D - yaa * recipa - yaa * recipb - yab * recipa - yab * recipb) * recip;@b@@b@ recipb += -yb * recip * recip;@b@@b@ recipa = -recipa;@b@ recipb = -recipb;@b@@b@ temp = ya + recipa;@b@ yb += -(temp - ya - recipa);@b@ ya = temp;@b@ temp = ya + recipb;@b@ yb += -(temp - ya - recipb);@b@ ya = temp;@b@@b@ result = ya + yb;@b@ result *= 0.5D;@b@ }@b@ else {@b@ hiPrec = new double[2];@b@ expm1(x, hiPrec);@b@@b@ ya = hiPrec[0] + hiPrec[1];@b@ yb = -(ya - hiPrec[0] - hiPrec[1]);@b@@b@ double denom = 1.0D + ya;@b@ double denomr = 1.0D / denom;@b@ double denomb = -(denom - 1.0D - ya) + yb;@b@ double ratio = ya * denomr;@b@ double temp = ratio * 1073741824.0D;@b@ double ra = ratio + temp - temp;@b@ double rb = ratio - ra;@b@@b@ temp = denom * 1073741824.0D;@b@ double za = denom + temp - temp;@b@ double zb = denom - za;@b@@b@ rb += (ya - za * ra - za * rb - zb * ra - zb * rb) * denomr;@b@@b@ rb += yb * denomr;@b@ rb += -ya * denomb * denomr * denomr;@b@@b@ temp = ya + ra;@b@ yb += -(temp - ya - ra);@b@ ya = temp;@b@ temp = ya + rb;@b@ yb += -(temp - ya - rb);@b@ ya = temp;@b@@b@ result = ya + yb;@b@ result *= 0.5D;@b@ }@b@@b@ if (negate) {@b@ result = -result;@b@ }@b@@b@ return result;@b@ }@b@@b@ public static double tanh(double x)@b@ {@b@ double result;@b@ double[] hiPrec;@b@ double ya;@b@ double yb;@b@ double na;@b@ double nb;@b@ double daa;@b@ double dab;@b@ double ratio;@b@ double ratioa;@b@ double ratiob;@b@ boolean negate = false;@b@@b@ if (x != x) {@b@ return x;@b@ }@b@@b@ if (x > 20.0D) {@b@ return 1.0D;@b@ }@b@@b@ if (x < -20.0D) {@b@ return -1.0D;@b@ }@b@@b@ if (x == 0.0D) {@b@ return x;@b@ }@b@@b@ if (x < 0.0D) {@b@ x = -x;@b@ negate = true;@b@ }@b@@b@ if (x >= 0.5D) {@b@ hiPrec = new double[2];@b@@b@ exp(x * 2.0D, 0.0D, hiPrec);@b@@b@ ya = hiPrec[0] + hiPrec[1];@b@ yb = -(ya - hiPrec[0] - hiPrec[1]);@b@@b@ na = -1.0D + ya;@b@ nb = -(na + 1.0D - ya);@b@ double temp = na + yb;@b@ nb += -(temp - na - yb);@b@ na = temp;@b@@b@ double da = 1.0D + ya;@b@ double db = -(da - 1.0D - ya);@b@ temp = da + yb;@b@ db += -(temp - da - yb);@b@ da = temp;@b@@b@ temp = da * 1073741824.0D;@b@ daa = da + temp - temp;@b@ dab = da - daa;@b@@b@ ratio = na / da;@b@ temp = ratio * 1073741824.0D;@b@ ratioa = ratio + temp - temp;@b@ ratiob = ratio - ratioa;@b@@b@ ratiob += (na - daa * ratioa - daa * ratiob - dab * ratioa - dab * ratiob) / da;@b@@b@ ratiob += nb / da;@b@@b@ ratiob += -db * na / da / da;@b@@b@ result = ratioa + ratiob;@b@ }@b@ else {@b@ hiPrec = new double[2];@b@@b@ expm1(x * 2.0D, hiPrec);@b@@b@ ya = hiPrec[0] + hiPrec[1];@b@ yb = -(ya - hiPrec[0] - hiPrec[1]);@b@@b@ na = ya;@b@ nb = yb;@b@@b@ double da = 2.0D + ya;@b@ double db = -(da - 2.0D - ya);@b@ double temp = da + yb;@b@ db += -(temp - da - yb);@b@ da = temp;@b@@b@ temp = da * 1073741824.0D;@b@ daa = da + temp - temp;@b@ dab = da - daa;@b@@b@ ratio = na / da;@b@ temp = ratio * 1073741824.0D;@b@ ratioa = ratio + temp - temp;@b@ ratiob = ratio - ratioa;@b@@b@ ratiob += (na - daa * ratioa - daa * ratiob - dab * ratioa - dab * ratiob) / da;@b@@b@ ratiob += nb / da;@b@@b@ ratiob += -db * na / da / da;@b@@b@ result = ratioa + ratiob;@b@ }@b@@b@ if (negate) {@b@ result = -result;@b@ }@b@@b@ return result;@b@ }@b@@b@ public static double acosh(double a)@b@ {@b@ return log(a + sqrt(a * a - 1.0D));@b@ }@b@@b@ public static double asinh(double a)@b@ {@b@ double absAsinh;@b@ boolean negative = false;@b@ if (a < 0.0D) {@b@ negative = true;@b@ a = -a;@b@ }@b@@b@ if (a > 0.167D) {@b@ absAsinh = log(sqrt(a * a + 1.0D) + a);@b@ } else {@b@ double a2 = a * a;@b@ if (a > 0.097D)@b@ absAsinh = a * (1.0D - a2 * (0.3333333333333333D - a2 * (0.2D - a2 * (0.1428571428571429D - a2 * (0.111111111111111D - a2 * (0.09090909090909091D - a2 * (0.07692307692307693D - a2 * (0.06666666666666667D - a2 * 0.05882352941176471D * 0.9375D) * 0.928571428571429D) * 0.9166666666666666D) * 0.9D) * 0.875D) * 0.8333333333333334D) * 0.75D) * 0.5D);@b@ else if (a > 0.036D)@b@ absAsinh = a * (1.0D - a2 * (0.3333333333333333D - a2 * (0.2D - a2 * (0.1428571428571429D - a2 * (0.111111111111111D - a2 * (0.09090909090909091D - a2 * 0.07692307692307693D * 0.9166666666666666D) * 0.9D) * 0.875D) * 0.8333333333333334D) * 0.75D) * 0.5D);@b@ else if (a > 0.0036D)@b@ absAsinh = a * (1.0D - a2 * (0.3333333333333333D - a2 * (0.2D - a2 * (0.1428571428571429D - a2 * 0.111111111111111D * 0.875D) * 0.8333333333333334D) * 0.75D) * 0.5D);@b@ else@b@ absAsinh = a * (1.0D - a2 * (0.3333333333333333D - a2 * 0.2D * 0.75D) * 0.5D);@b@@b@ }@b@@b@ return ((negative) ? -absAsinh : absAsinh);@b@ }@b@@b@ public static double atanh(double a)@b@ {@b@ double absAtanh;@b@ boolean negative = false;@b@ if (a < 0.0D) {@b@ negative = true;@b@ a = -a;@b@ }@b@@b@ if (a > 0.15D) {@b@ absAtanh = 0.5D * log((1.0D + a) / (1.0D - a));@b@ } else {@b@ double a2 = a * a;@b@ if (a > 0.08699999999999999D)@b@ absAtanh = a * (1.0D + a2 * (0.3333333333333333D + a2 * (0.2D + a2 * (0.1428571428571429D + a2 * (0.111111111111111D + a2 * (0.09090909090909091D + a2 * (0.07692307692307693D + a2 * (0.06666666666666667D + a2 * 0.05882352941176471D))))))));@b@ else if (a > 0.031D)@b@ absAtanh = a * (1.0D + a2 * (0.3333333333333333D + a2 * (0.2D + a2 * (0.1428571428571429D + a2 * (0.111111111111111D + a2 * (0.09090909090909091D + a2 * 0.07692307692307693D))))));@b@ else if (a > 0.003D)@b@ absAtanh = a * (1.0D + a2 * (0.3333333333333333D + a2 * (0.2D + a2 * (0.1428571428571429D + a2 * 0.111111111111111D))));@b@ else@b@ absAtanh = a * (1.0D + a2 * (0.3333333333333333D + a2 * 0.2D));@b@@b@ }@b@@b@ return ((negative) ? -absAtanh : absAtanh);@b@ }@b@@b@ public static double signum(double a)@b@ {@b@ return ((a > 0.0D) ? 1.0D : (a < 0.0D) ? -1.0D : a);@b@ }@b@@b@ public static float signum(float a)@b@ {@b@ return ((a > 0.0F) ? 1.0F : (a < 0.0F) ? -1.0F : a);@b@ }@b@@b@ public static double nextUp(double a)@b@ {@b@ return nextAfter(a, (1.0D / 0.0D));@b@ }@b@@b@ public static float nextUp(float a)@b@ {@b@ return nextAfter(a, (1.0D / 0.0D));@b@ }@b@@b@ public static double random()@b@ {@b@ return Math.random();@b@ }@b@@b@ public static double exp(double x)@b@ {@b@ return exp(x, 0.0D, null);@b@ }@b@@b@ private static double exp(double x, double extra, double[] hiPrec)@b@ {@b@ double intPartA;@b@ double intPartB;@b@ int intVal;@b@ double result;@b@ if (x < 0.0D) {@b@ double result;@b@ intVal = (int)(-x);@b@@b@ if (intVal > 746) {@b@ if (hiPrec != null) {@b@ hiPrec[0] = 0.0D;@b@ hiPrec[1] = 0.0D;@b@ }@b@ return 0.0D;@b@ }@b@@b@ if (intVal > 709)@b@ {@b@ result = exp(x + 40.19140625D, extra, hiPrec) / 2.850400951440118E+017D;@b@ if (hiPrec != null) {@b@ hiPrec[0] /= 2.850400951440118E+017D;@b@ hiPrec[1] /= 2.850400951440118E+017D;@b@ }@b@ return result;@b@ }@b@@b@ if (intVal == 709)@b@ {@b@ result = exp(x + 1.494140625D, extra, hiPrec) / 4.455505956692757D;@b@ if (hiPrec != null) {@b@ hiPrec[0] /= 4.455505956692757D;@b@ hiPrec[1] /= 4.455505956692757D;@b@ }@b@ return result;@b@ }@b@@b@ ++intVal;@b@@b@ intPartA = ExpIntTable.access$000()[(750 - intVal)];@b@ intPartB = ExpIntTable.access$100()[(750 - intVal)];@b@@b@ intVal = -intVal;@b@ } else {@b@ intVal = (int)x;@b@@b@ if (intVal > 709) {@b@ if (hiPrec != null) {@b@ hiPrec[0] = (1.0D / 0.0D);@b@ hiPrec[1] = 0.0D;@b@ }@b@ return (1.0D / 0.0D);@b@ }@b@@b@ intPartA = ExpIntTable.access$000()[(750 + intVal)];@b@ intPartB = ExpIntTable.access$100()[(750 + intVal)];@b@ }@b@@b@ int intFrac = (int)((x - intVal) * 1024.0D);@b@ double fracPartA = ExpFracTable.access$200()[intFrac];@b@ double fracPartB = ExpFracTable.access$300()[intFrac];@b@@b@ double epsilon = x - intVal + intFrac / 1024.0D;@b@@b@ double z = 0.04168701738764507D;@b@ z = z * epsilon + 0.1666666505023083D;@b@ z = z * epsilon + 0.500000000004269D;@b@ z = z * epsilon + 1.0D;@b@ z = z * epsilon + -3.940510424527919E-020D;@b@@b@ double tempA = intPartA * fracPartA;@b@ double tempB = intPartA * fracPartB + intPartB * fracPartA + intPartB * fracPartB;@b@@b@ double tempC = tempB + tempA;@b@@b@ if (extra != 0.0D)@b@ result = tempC * extra * z + tempC * extra + tempC * z + tempB + tempA;@b@ else {@b@ result = tempC * z + tempB + tempA;@b@ }@b@@b@ if (hiPrec != null)@b@ {@b@ hiPrec[0] = tempA;@b@ hiPrec[1] = (tempC * extra * z + tempC * extra + tempC * z + tempB);@b@ }@b@@b@ return result;@b@ }@b@@b@ public static double expm1(double x)@b@ {@b@ return expm1(x, null);@b@ }@b@@b@ private static double expm1(double x, double[] hiPrecOut)@b@ {@b@ if ((x != x) || (x == 0.0D)) {@b@ return x;@b@ }@b@@b@ if ((x <= -1.0D) || (x >= 1.0D))@b@ {@b@ double[] hiPrec = new double[2];@b@ exp(x, 0.0D, hiPrec);@b@ if (x > 0.0D)@b@ return (-1.0D + hiPrec[0] + hiPrec[1]);@b@@b@ double ra = -1.0D + hiPrec[0];@b@ double rb = -(ra + 1.0D - hiPrec[0]);@b@ rb += hiPrec[1];@b@ return (ra + rb);@b@ }@b@@b@ boolean negative = false;@b@@b@ if (x < 0.0D) {@b@ x = -x;@b@ negative = true;@b@ }@b@@b@ int intFrac = (int)(x * 1024.0D);@b@ double tempA = ExpFracTable.access$200()[intFrac] - 1.0D;@b@ double tempB = ExpFracTable.access$300()[intFrac];@b@@b@ double temp = tempA + tempB;@b@ tempB = -(temp - tempA - tempB);@b@ tempA = temp;@b@@b@ temp = tempA * 1073741824.0D;@b@ double baseA = tempA + temp - temp;@b@ double baseB = tempB + tempA - baseA;@b@@b@ double epsilon = x - intFrac / 1024.0D;@b@@b@ double zb = 0.008336750013465571D;@b@ zb = zb * epsilon + 0.04166666387918665D;@b@ zb = zb * epsilon + 0.1666666666674539D;@b@ zb = zb * epsilon + 0.4999999999999999D;@b@ zb *= epsilon;@b@ zb *= epsilon;@b@@b@ double za = epsilon;@b@ double temp = za + zb;@b@ zb = -(temp - za - zb);@b@ za = temp;@b@@b@ temp = za * 1073741824.0D;@b@ temp = za + temp - temp;@b@ zb += za - temp;@b@ za = temp;@b@@b@ double ya = za * baseA;@b@@b@ temp = ya + za * baseB;@b@ double yb = -(temp - ya - za * baseB);@b@ ya = temp;@b@@b@ temp = ya + zb * baseA;@b@ yb += -(temp - ya - zb * baseA);@b@ ya = temp;@b@@b@ temp = ya + zb * baseB;@b@ yb += -(temp - ya - zb * baseB);@b@ ya = temp;@b@@b@ temp = ya + baseA;@b@ yb += -(temp - baseA - ya);@b@ ya = temp;@b@@b@ temp = ya + za;@b@@b@ yb += -(temp - ya - za);@b@ ya = temp;@b@@b@ temp = ya + baseB;@b@@b@ yb += -(temp - ya - baseB);@b@ ya = temp;@b@@b@ temp = ya + zb;@b@@b@ yb += -(temp - ya - zb);@b@ ya = temp;@b@@b@ if (negative)@b@ {@b@ double denom = 1.0D + ya;@b@ double denomr = 1.0D / denom;@b@ double denomb = -(denom - 1.0D - ya) + yb;@b@ double ratio = ya * denomr;@b@ temp = ratio * 1073741824.0D;@b@ double ra = ratio + temp - temp;@b@ double rb = ratio - ra;@b@@b@ temp = denom * 1073741824.0D;@b@ za = denom + temp - temp;@b@ zb = denom - za;@b@@b@ rb += (ya - za * ra - za * rb - zb * ra - zb * rb) * denomr;@b@@b@ rb += yb * denomr;@b@ rb += -ya * denomb * denomr * denomr;@b@@b@ ya = -ra;@b@ yb = -rb;@b@ }@b@@b@ if (hiPrecOut != null) {@b@ hiPrecOut[0] = ya;@b@ hiPrecOut[1] = yb;@b@ }@b@@b@ return (ya + yb);@b@ }@b@@b@ public static double log(double x)@b@ {@b@ return log(x, null);@b@ }@b@@b@ private static double log(double x, double[] hiPrec)@b@ {@b@ if (x == 0.0D)@b@ return (-1.0D / 0.0D);@b@@b@ long bits = Double.doubleToRawLongBits(x);@b@@b@ if (((((bits & 0x0) != 0L) || (x != x))) && (x != 0.0D)) {@b@ if (hiPrec != null) {@b@ hiPrec[0] = (0.0D / 0.0D);@b@ }@b@@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if (x == (1.0D / 0.0D)) {@b@ if (hiPrec != null) {@b@ hiPrec[0] = (1.0D / 0.0D);@b@ }@b@@b@ return (1.0D / 0.0D);@b@ }@b@@b@ int exp = (int)(bits >> 52) - 1023;@b@@b@ if ((bits & 0x0) == 0L)@b@ {@b@ if (x == 0.0D)@b@ {@b@ if (hiPrec != null) {@b@ hiPrec[0] = (-1.0D / 0.0D);@b@ }@b@@b@ return (-1.0D / 0.0D);@b@ }@b@@b@ bits <<= 1;@b@ while ((bits & 0x0) == 0L) {@b@ --exp;@b@ bits <<= 1;@b@ }@b@@b@ }@b@@b@ if ((((exp == -1) || (exp == 0))) && (x < 1.01D) && (x > 0.99D) && (hiPrec == null))@b@ {@b@ double xa = x - 1.0D;@b@ double xb = xa - x + 1.0D;@b@ double tmp = xa * 1073741824.0D;@b@ double aa = xa + tmp - tmp;@b@ double ab = xa - aa;@b@ xa = aa;@b@ xb = ab;@b@@b@ double[] lnCoef_last = LN_QUICK_COEF[(LN_QUICK_COEF.length - 1)];@b@ double ya = lnCoef_last[0];@b@ double yb = lnCoef_last[1];@b@@b@ for (int i = LN_QUICK_COEF.length - 2; i >= 0; --i)@b@ {@b@ aa = ya * xa;@b@ ab = ya * xb + yb * xa + yb * xb;@b@@b@ tmp = aa * 1073741824.0D;@b@ ya = aa + tmp - tmp;@b@ yb = aa - ya + ab;@b@@b@ double[] lnCoef_i = LN_QUICK_COEF[i];@b@ aa = ya + lnCoef_i[0];@b@ ab = yb + lnCoef_i[1];@b@@b@ tmp = aa * 1073741824.0D;@b@ ya = aa + tmp - tmp;@b@ yb = aa - ya + ab;@b@ }@b@@b@ aa = ya * xa;@b@ ab = ya * xb + yb * xa + yb * xb;@b@@b@ tmp = aa * 1073741824.0D;@b@ ya = aa + tmp - tmp;@b@ yb = aa - ya + ab;@b@@b@ return (ya + yb);@b@ }@b@@b@ double[] lnm = lnMant.access$400()[(int)((bits & 0x0) >> 42)];@b@@b@ double epsilon = (bits & 0xFFFFFFFF) / (4503599627370496.0D + (bits & 0x0));@b@@b@ double lnza = 0.0D;@b@ double lnzb = 0.0D;@b@@b@ if (hiPrec != null)@b@ {@b@ double tmp = epsilon * 1073741824.0D;@b@ double aa = epsilon + tmp - tmp;@b@ double ab = epsilon - aa;@b@ double xa = aa;@b@ double xb = ab;@b@@b@ double numer = bits & 0xFFFFFFFF;@b@ double denom = 4503599627370496.0D + (bits & 0x0);@b@ aa = numer - xa * denom - xb * denom;@b@ xb += aa / denom;@b@@b@ double[] lnCoef_last = LN_HI_PREC_COEF[(LN_HI_PREC_COEF.length - 1)];@b@ double ya = lnCoef_last[0];@b@ double yb = lnCoef_last[1];@b@@b@ for (int i = LN_HI_PREC_COEF.length - 2; i >= 0; --i)@b@ {@b@ aa = ya * xa;@b@ ab = ya * xb + yb * xa + yb * xb;@b@@b@ tmp = aa * 1073741824.0D;@b@ ya = aa + tmp - tmp;@b@ yb = aa - ya + ab;@b@@b@ double[] lnCoef_i = LN_HI_PREC_COEF[i];@b@ aa = ya + lnCoef_i[0];@b@ ab = yb + lnCoef_i[1];@b@@b@ tmp = aa * 1073741824.0D;@b@ ya = aa + tmp - tmp;@b@ yb = aa - ya + ab;@b@ }@b@@b@ aa = ya * xa;@b@ ab = ya * xb + yb * xa + yb * xb;@b@@b@ lnza = aa + ab;@b@ lnzb = -(lnza - aa - ab);@b@ }@b@ else@b@ {@b@ lnza = -0.1662488244041857D;@b@ lnza = lnza * epsilon + 0.1999995412025452D;@b@ lnza = lnza * epsilon + -0.24999999976775D;@b@ lnza = lnza * epsilon + 0.3333333333332802D;@b@ lnza = lnza * epsilon + -0.5D;@b@ lnza = lnza * epsilon + 1.0D;@b@ lnza *= epsilon;@b@ }@b@@b@ double a = 0.6931470632553101D * exp;@b@ double b = 0.0D;@b@ double c = a + lnm[0];@b@ double d = -(c - a - lnm[0]);@b@ a = c;@b@ b += d;@b@@b@ c = a + lnza;@b@ d = -(c - a - lnza);@b@ a = c;@b@ b += d;@b@@b@ c = a + 1.173046352508235E-007D * exp;@b@ d = -(c - a - 1.173046352508235E-007D * exp);@b@ a = c;@b@ b += d;@b@@b@ c = a + lnm[1];@b@ d = -(c - a - lnm[1]);@b@ a = c;@b@ b += d;@b@@b@ c = a + lnzb;@b@ d = -(c - a - lnzb);@b@ a = c;@b@ b += d;@b@@b@ if (hiPrec != null) {@b@ hiPrec[0] = a;@b@ hiPrec[1] = b;@b@ }@b@@b@ return (a + b);@b@ }@b@@b@ public static double log1p(double x)@b@ {@b@ if (x == -1.0D) {@b@ return (-1.0D / 0.0D);@b@ }@b@@b@ if (x == (1.0D / 0.0D)) {@b@ return (1.0D / 0.0D);@b@ }@b@@b@ if ((x > 9/E-007D) || (x < -9/E-007D))@b@ {@b@ double xpa = 1.0D + x;@b@ double xpb = -(xpa - 1.0D - x);@b@@b@ double[] hiPrec = new double[2];@b@ double lores = log(xpa, hiPrec);@b@ if (Double.isInfinite(lores)) {@b@ return lores;@b@ }@b@@b@ double fx1 = xpb / xpa;@b@ double epsilon = 0.5D * fx1 + 1.0D;@b@ return (epsilon * fx1 + hiPrec[1] + hiPrec[0]);@b@ }@b@@b@ double y = (x * 0.3333333333333333D - 0.5D) * x + 1.0D;@b@ return (y * x);@b@ }@b@@b@ public static double log10(double x)@b@ {@b@ double[] hiPrec = new double[2];@b@@b@ double lores = log(x, hiPrec);@b@ if (Double.isInfinite(lores)) {@b@ return lores;@b@ }@b@@b@ double tmp = hiPrec[0] * 1073741824.0D;@b@ double lna = hiPrec[0] + tmp - tmp;@b@ double lnb = hiPrec[0] - lna + hiPrec[1];@b@@b@ double rln10a = 0.4342944622039795D;@b@ double rln10b = 1.969927233546363E-008D;@b@@b@ return (1.969927233546363E-008D * lnb + 1.969927233546363E-008D * lna + 0.4342944622039795D * lnb + 0.4342944622039795D * lna);@b@ }@b@@b@ public static double log(double base, double x)@b@ {@b@ return (log(x) / log(base));@b@ }@b@@b@ public static double pow(double x, double y)@b@ {@b@ double ya;@b@ double yb;@b@ double tmp1;@b@ double[] lns = new double[2];@b@@b@ if (y == 0.0D) {@b@ return 1.0D;@b@ }@b@@b@ if (x != x) {@b@ return x;@b@ }@b@@b@ if (x == 0.0D) {@b@ long bits = Double.doubleToRawLongBits(x);@b@ if ((bits & 0x0) != 0L)@b@ {@b@ long yi = ()y;@b@@b@ if ((y < 0.0D) && (y == yi) && ((yi & 1L) == 1L)) {@b@ return (-1.0D / 0.0D);@b@ }@b@@b@ if ((y > 0.0D) && (y == yi) && ((yi & 1L) == 1L))@b@ return -0.0D;@b@@b@ }@b@@b@ if (y < 0.0D)@b@ return (1.0D / 0.0D);@b@@b@ if (y > 0.0D) {@b@ return 0.0D;@b@ }@b@@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if (x == (1.0D / 0.0D)) {@b@ if (y != y)@b@ return y;@b@@b@ if (y < 0.0D)@b@ return 0.0D;@b@@b@ return (1.0D / 0.0D);@b@ }@b@@b@ if (y == (1.0D / 0.0D)) {@b@ if (x * x == 1.0D) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if (x * x > 1.0D)@b@ return (1.0D / 0.0D);@b@@b@ return 0.0D;@b@ }@b@@b@ if (x == (-1.0D / 0.0D)) {@b@ long yi;@b@ if (y != y) {@b@ return y;@b@ }@b@@b@ if (y < 0.0D) {@b@ yi = ()y;@b@ if ((y == yi) && ((yi & 1L) == 1L)) {@b@ return -0.0D;@b@ }@b@@b@ return 0.0D;@b@ }@b@@b@ if (y > 0.0D) {@b@ yi = ()y;@b@ if ((y == yi) && ((yi & 1L) == 1L)) {@b@ return (-1.0D / 0.0D);@b@ }@b@@b@ return (1.0D / 0.0D);@b@ }@b@ }@b@@b@ if (y == (-1.0D / 0.0D))@b@ {@b@ if (x * x == 1.0D) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if (x * x < 1.0D)@b@ return (1.0D / 0.0D);@b@@b@ return 0.0D;@b@ }@b@@b@ if (x < 0.0D)@b@ {@b@ if ((y >= 9007199254740992.0D) || (y <= -9007199254740992.0D)) {@b@ return pow(-x, y);@b@ }@b@@b@ if (y == ()y)@b@ {@b@ return (((()y & 1L) == 0L) ? pow(-x, y) : -pow(-x, y));@b@ }@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if ((y < 7/E+298D) && (y > -7/E+298D)) {@b@ tmp1 = y * 1073741824.0D;@b@ ya = y + tmp1 - tmp1;@b@ yb = y - ya;@b@ } else {@b@ tmp1 = y * 9.313225746154785E-010D;@b@ double tmp2 = tmp1 * 9.313225746154785E-010D;@b@ ya = (tmp1 + tmp2 - tmp1) * 1073741824.0D * 1073741824.0D;@b@ yb = y - ya;@b@ }@b@@b@ double lores = log(x, lns);@b@ if (Double.isInfinite(lores)) {@b@ return lores;@b@ }@b@@b@ double lna = lns[0];@b@ double lnb = lns[1];@b@@b@ double tmp1 = lna * 1073741824.0D;@b@ double tmp2 = lna + tmp1 - tmp1;@b@ lnb += lna - tmp2;@b@ lna = tmp2;@b@@b@ double aa = lna * ya;@b@ double ab = lna * yb + lnb * ya + lnb * yb;@b@@b@ lna = aa + ab;@b@ lnb = -(lna - aa - ab);@b@@b@ double z = 0.008333333333333333D;@b@ z = z * lnb + 0.04166666666666666D;@b@ z = z * lnb + 0.1666666666666667D;@b@ z = z * lnb + 0.5D;@b@ z = z * lnb + 1.0D;@b@ z *= lnb;@b@@b@ double result = exp(lna, z, null);@b@@b@ return result;@b@ }@b@@b@ public static double pow(double d, int e)@b@ {@b@ if (e == 0)@b@ return 1.0D;@b@ if (e < 0) {@b@ e = -e;@b@ d = 1.0D / d;@b@ }@b@@b@ int splitFactor = 134217729;@b@ double cd = 134217729.0D * d;@b@ double d1High = cd - cd - d;@b@ double d1Low = d - d1High;@b@@b@ double resultHigh = 1.0D;@b@ double resultLow = 0.0D;@b@@b@ double d2p = d;@b@ double d2pHigh = d1High;@b@ double d2pLow = d1Low;@b@@b@ while (e != 0)@b@ {@b@ if ((e & 0x1) != 0)@b@ {@b@ tmpHigh = resultHigh * d2p;@b@ double cRH = 134217729.0D * resultHigh;@b@ double rHH = cRH - cRH - resultHigh;@b@ double rHL = resultHigh - rHH;@b@ tmpLow = rHL * d2pLow - tmpHigh - rHH * d2pHigh - rHL * d2pHigh - rHH * d2pLow;@b@ resultHigh = tmpHigh;@b@ resultLow = resultLow * d2p + tmpLow;@b@ }@b@@b@ double tmpHigh = d2pHigh * d2p;@b@ double cD2pH = 134217729.0D * d2pHigh;@b@ double d2pHH = cD2pH - cD2pH - d2pHigh;@b@ double d2pHL = d2pHigh - d2pHH;@b@ double tmpLow = d2pHL * d2pLow - tmpHigh - d2pHH * d2pHigh - d2pHL * d2pHigh - d2pHH * d2pLow;@b@ double cTmpH = 134217729.0D * tmpHigh;@b@ d2pHigh = cTmpH - cTmpH - tmpHigh;@b@ d2pLow = d2pLow * d2p + tmpLow + tmpHigh - d2pHigh;@b@ d2p = d2pHigh + d2pLow;@b@@b@ e >>= 1;@b@ }@b@@b@ return (resultHigh + resultLow);@b@ }@b@@b@ private static double polySine(double x)@b@ {@b@ double x2 = x * x;@b@@b@ double p = 2.755381745227222E-006D;@b@ p = p * x2 + -0.0001984126965958651D;@b@ p = p * x2 + 0.008333333333329196D;@b@ p = p * x2 + -0.1666666666666667D;@b@@b@ p = p * x2 * x;@b@@b@ return p;@b@ }@b@@b@ private static double polyCosine(double x)@b@ {@b@ double x2 = x * x;@b@@b@ double p = 2.479773539153719E-005D;@b@ p = p * x2 + -0.001388888868903988D;@b@ p = p * x2 + 0.04166666666662117D;@b@ p = p * x2 + -0.4999999999999999D;@b@ p *= x2;@b@@b@ return p;@b@ }@b@@b@ private static double sinQ(double xa, double xb)@b@ {@b@ int idx = (int)(xa * 8.0D + 0.5D);@b@ double epsilon = xa - EIGHTHS[idx];@b@@b@ double sintA = SINE_TABLE_A[idx];@b@ double sintB = SINE_TABLE_B[idx];@b@ double costA = COSINE_TABLE_A[idx];@b@ double costB = COSINE_TABLE_B[idx];@b@@b@ double sinEpsA = epsilon;@b@ double sinEpsB = polySine(epsilon);@b@ double cosEpsA = 1.0D;@b@ double cosEpsB = polyCosine(epsilon);@b@@b@ double temp = sinEpsA * 1073741824.0D;@b@ double temp2 = sinEpsA + temp - temp;@b@ sinEpsB += sinEpsA - temp2;@b@ sinEpsA = temp2;@b@@b@ double a = 0.0D;@b@ double b = 0.0D;@b@@b@ double t = sintA;@b@ double c = a + t;@b@ double d = -(c - a - t);@b@ a = c;@b@ b += d;@b@@b@ t = costA * sinEpsA;@b@ c = a + t;@b@ d = -(c - a - t);@b@ a = c;@b@ b += d;@b@@b@ b = b + sintA * cosEpsB + costA * sinEpsB;@b@@b@ b = b + sintB + costB * sinEpsA + sintB * cosEpsB + costB * sinEpsB;@b@@b@ if (xb != 0.0D) {@b@ t = ((costA + costB) * (1.0D + cosEpsB) - (sintA + sintB) * (sinEpsA + sinEpsB)) * xb;@b@@b@ c = a + t;@b@ d = -(c - a - t);@b@ a = c;@b@ b += d;@b@ }@b@@b@ double result = a + b;@b@@b@ return result;@b@ }@b@@b@ private static double cosQ(double xa, double xb)@b@ {@b@ double pi2a = 1.570796326794897D;@b@ double pi2b = 6.123233995736766E-017D;@b@@b@ double a = 1.570796326794897D - xa;@b@ double b = -(a - 1.570796326794897D + xa);@b@ b += 6.123233995736766E-017D - xb;@b@@b@ return sinQ(a, b);@b@ }@b@@b@ private static double tanQ(double xa, double xb, boolean cotanFlag)@b@ {@b@ int idx = (int)(xa * 8.0D + 0.5D);@b@ double epsilon = xa - EIGHTHS[idx];@b@@b@ double sintA = SINE_TABLE_A[idx];@b@ double sintB = SINE_TABLE_B[idx];@b@ double costA = COSINE_TABLE_A[idx];@b@ double costB = COSINE_TABLE_B[idx];@b@@b@ double sinEpsA = epsilon;@b@ double sinEpsB = polySine(epsilon);@b@ double cosEpsA = 1.0D;@b@ double cosEpsB = polyCosine(epsilon);@b@@b@ double temp = sinEpsA * 1073741824.0D;@b@ double temp2 = sinEpsA + temp - temp;@b@ sinEpsB += sinEpsA - temp2;@b@ sinEpsA = temp2;@b@@b@ double a = 0.0D;@b@ double b = 0.0D;@b@@b@ double t = sintA;@b@ double c = a + t;@b@ double d = -(c - a - t);@b@ a = c;@b@ b += d;@b@@b@ t = costA * sinEpsA;@b@ c = a + t;@b@ d = -(c - a - t);@b@ a = c;@b@ b += d;@b@@b@ b = b + sintA * cosEpsB + costA * sinEpsB;@b@ b = b + sintB + costB * sinEpsA + sintB * cosEpsB + costB * sinEpsB;@b@@b@ double sina = a + b;@b@ double sinb = -(sina - a - b);@b@@b@ a = b = c = d = 0.0D;@b@@b@ t = costA * 1.0D;@b@ c = a + t;@b@ d = -(c - a - t);@b@ a = c;@b@ b += d;@b@@b@ t = -sintA * sinEpsA;@b@ c = a + t;@b@ d = -(c - a - t);@b@ a = c;@b@ b += d;@b@@b@ b = b + costB * 1.0D + costA * cosEpsB + costB * cosEpsB;@b@ b -= sintB * sinEpsA + sintA * sinEpsB + sintB * sinEpsB;@b@@b@ double cosa = a + b;@b@ double cosb = -(cosa - a - b);@b@@b@ if (cotanFlag)@b@ {@b@ double tmp = cosa; cosa = sina; sina = tmp;@b@ tmp = cosb; cosb = sinb; sinb = tmp;@b@ }@b@@b@ double est = sina / cosa;@b@@b@ temp = est * 1073741824.0D;@b@ double esta = est + temp - temp;@b@ double estb = est - esta;@b@@b@ temp = cosa * 1073741824.0D;@b@ double cosaa = cosa + temp - temp;@b@ double cosab = cosa - cosaa;@b@@b@ double err = (sina - esta * cosaa - esta * cosab - estb * cosaa - estb * cosab) / cosa;@b@ err += sinb / cosa;@b@ err += -sina * cosb / cosa / cosa;@b@@b@ if (xb != 0.0D)@b@ {@b@ double xbadj = xb + est * est * xb;@b@ if (cotanFlag) {@b@ xbadj = -xbadj;@b@ }@b@@b@ err += xbadj;@b@ }@b@@b@ return (est + err);@b@ }@b@@b@ private static void reducePayneHanek(double x, double[] result)@b@ {@b@ long shpi0;@b@ long shpiA;@b@ long shpiB;@b@ long inbits = Double.doubleToRawLongBits(x);@b@ int exponent = (int)(inbits >> 52 & 0x7FF) - 1023;@b@@b@ inbits &= 4503599627370495L;@b@ inbits |= 4503599627370496L;@b@@b@ ++exponent;@b@ inbits <<= 11;@b@@b@ int idx = exponent >> 6;@b@ int shift = exponent - (idx << 6);@b@@b@ if (shift != 0) {@b@ shpi0 = (idx == 0) ? 0L : RECIP_2PI[(idx - 1)] << shift;@b@ shpi0 |= RECIP_2PI[idx] >>> 64 - shift;@b@ shpiA = RECIP_2PI[idx] << shift | RECIP_2PI[(idx + 1)] >>> 64 - shift;@b@ shpiB = RECIP_2PI[(idx + 1)] << shift | RECIP_2PI[(idx + 2)] >>> 64 - shift;@b@ } else {@b@ shpi0 = (idx == 0) ? 0L : RECIP_2PI[(idx - 1)];@b@ shpiA = RECIP_2PI[idx];@b@ shpiB = RECIP_2PI[(idx + 1)];@b@ }@b@@b@ long a = inbits >>> 32;@b@ long b = inbits & 0xFFFFFFFF;@b@@b@ long c = shpiA >>> 32;@b@ long d = shpiA & 0xFFFFFFFF;@b@@b@ long ac = a * c;@b@ long bd = b * d;@b@ long bc = b * c;@b@ long ad = a * d;@b@@b@ long prodB = bd + (ad << 32);@b@ long prodA = ac + (ad >>> 32);@b@@b@ boolean bita = (bd & 0x0) != 0L;@b@ boolean bitb = (ad & 0x80000000) != 0L;@b@ boolean bitsum = (prodB & 0x0) != 0L;@b@@b@ if (((bita) && (bitb)) || ((((bita) || (bitb))) && (!(bitsum))))@b@ {@b@ prodA += 1L;@b@ }@b@@b@ bita = (prodB & 0x0) != 0L;@b@ bitb = (bc & 0x80000000) != 0L;@b@@b@ prodB += (bc << 32);@b@ prodA += (bc >>> 32);@b@@b@ bitsum = (prodB & 0x0) != 0L;@b@@b@ if (((bita) && (bitb)) || ((((bita) || (bitb))) && (!(bitsum))))@b@ {@b@ prodA += 1L;@b@ }@b@@b@ c = shpiB >>> 32;@b@ d = shpiB & 0xFFFFFFFF;@b@ ac = a * c;@b@ bc = b * c;@b@ ad = a * d;@b@@b@ ac += (bc + ad >>> 32);@b@@b@ bita = (prodB & 0x0) != 0L;@b@ bitb = (ac & 0x0) != 0L;@b@ prodB += ac;@b@ bitsum = (prodB & 0x0) != 0L;@b@@b@ if (((bita) && (bitb)) || ((((bita) || (bitb))) && (!(bitsum))))@b@ {@b@ prodA += 1L;@b@ }@b@@b@ c = shpi0 >>> 32;@b@ d = shpi0 & 0xFFFFFFFF;@b@@b@ bd = b * d;@b@ bc = b * c;@b@ ad = a * d;@b@@b@ prodA += bd + (bc + ad << 32);@b@@b@ int intPart = (int)(prodA >>> 62);@b@@b@ prodA <<= 2;@b@ prodA |= prodB >>> 62;@b@ prodB <<= 2;@b@@b@ a = prodA >>> 32;@b@ b = prodA & 0xFFFFFFFF;@b@@b@ c = PI_O_4_BITS[0] >>> 32;@b@ d = PI_O_4_BITS[0] & 0xFFFFFFFF;@b@@b@ ac = a * c;@b@ bd = b * d;@b@ bc = b * c;@b@ ad = a * d;@b@@b@ long prod2B = bd + (ad << 32);@b@ long prod2A = ac + (ad >>> 32);@b@@b@ bita = (bd & 0x0) != 0L;@b@ bitb = (ad & 0x80000000) != 0L;@b@ bitsum = (prod2B & 0x0) != 0L;@b@@b@ if (((bita) && (bitb)) || ((((bita) || (bitb))) && (!(bitsum))))@b@ {@b@ prod2A += 1L;@b@ }@b@@b@ bita = (prod2B & 0x0) != 0L;@b@ bitb = (bc & 0x80000000) != 0L;@b@@b@ prod2B += (bc << 32);@b@ prod2A += (bc >>> 32);@b@@b@ bitsum = (prod2B & 0x0) != 0L;@b@@b@ if (((bita) && (bitb)) || ((((bita) || (bitb))) && (!(bitsum))))@b@ {@b@ prod2A += 1L;@b@ }@b@@b@ c = PI_O_4_BITS[1] >>> 32;@b@ d = PI_O_4_BITS[1] & 0xFFFFFFFF;@b@ ac = a * c;@b@ bc = b * c;@b@ ad = a * d;@b@@b@ ac += (bc + ad >>> 32);@b@@b@ bita = (prod2B & 0x0) != 0L;@b@ bitb = (ac & 0x0) != 0L;@b@ prod2B += ac;@b@ bitsum = (prod2B & 0x0) != 0L;@b@@b@ if (((bita) && (bitb)) || ((((bita) || (bitb))) && (!(bitsum))))@b@ {@b@ prod2A += 1L;@b@ }@b@@b@ a = prodB >>> 32;@b@ b = prodB & 0xFFFFFFFF;@b@ c = PI_O_4_BITS[0] >>> 32;@b@ d = PI_O_4_BITS[0] & 0xFFFFFFFF;@b@ ac = a * c;@b@ bc = b * c;@b@ ad = a * d;@b@@b@ ac += (bc + ad >>> 32);@b@@b@ bita = (prod2B & 0x0) != 0L;@b@ bitb = (ac & 0x0) != 0L;@b@ prod2B += ac;@b@ bitsum = (prod2B & 0x0) != 0L;@b@@b@ if (((bita) && (bitb)) || ((((bita) || (bitb))) && (!(bitsum))))@b@ {@b@ prod2A += 1L;@b@ }@b@@b@ double tmpA = (prod2A >>> 12) / 4503599627370496.0D;@b@ double tmpB = (((prod2A & 0xFFF) << 40) + (prod2B >>> 24)) / 4503599627370496.0D / 4503599627370496.0D;@b@@b@ double sumA = tmpA + tmpB;@b@ double sumB = -(sumA - tmpA - tmpB);@b@@b@ result[0] = intPart;@b@ result[1] = (sumA * 2.0D);@b@ result[2] = (sumB * 2.0D);@b@ }@b@@b@ public static double sin(double x)@b@ {@b@ boolean negative = false;@b@ int quadrant = 0;@b@@b@ double xb = 0.0D;@b@@b@ double xa = x;@b@ if (x < 0.0D) {@b@ negative = true;@b@ xa = -xa;@b@ }@b@@b@ if (xa == 0.0D) {@b@ long bits = Double.doubleToRawLongBits(x);@b@ if (bits < 0L)@b@ return -0.0D;@b@@b@ return 0.0D;@b@ }@b@@b@ if ((xa != xa) || (xa == (1.0D / 0.0D))) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if (xa > 3294198.0D)@b@ {@b@ double[] reduceResults = new double[3];@b@ reducePayneHanek(xa, reduceResults);@b@ quadrant = (int)reduceResults[0] & 0x3;@b@ xa = reduceResults[1];@b@ xb = reduceResults[2];@b@ } else if (xa > 1.570796326794897D) {@b@ CodyWaite cw = new CodyWaite(xa);@b@ quadrant = cw.getK() & 0x3;@b@ xa = cw.getRemA();@b@ xb = cw.getRemB();@b@ }@b@@b@ if (negative) {@b@ quadrant ^= 2;@b@ }@b@@b@ switch (quadrant)@b@ {@b@ case 0:@b@ return sinQ(xa, xb);@b@ case 1:@b@ return cosQ(xa, xb);@b@ case 2:@b@ return (-sinQ(xa, xb));@b@ case 3:@b@ return (-cosQ(xa, xb));@b@ }@b@ return (0.0D / 0.0D);@b@ }@b@@b@ public static double cos(double x)@b@ {@b@ int quadrant = 0;@b@@b@ double xa = x;@b@ if (x < 0.0D) {@b@ xa = -xa;@b@ }@b@@b@ if ((xa != xa) || (xa == (1.0D / 0.0D))) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ double xb = 0.0D;@b@ if (xa > 3294198.0D)@b@ {@b@ double[] reduceResults = new double[3];@b@ reducePayneHanek(xa, reduceResults);@b@ quadrant = (int)reduceResults[0] & 0x3;@b@ xa = reduceResults[1];@b@ xb = reduceResults[2];@b@ } else if (xa > 1.570796326794897D) {@b@ CodyWaite cw = new CodyWaite(xa);@b@ quadrant = cw.getK() & 0x3;@b@ xa = cw.getRemA();@b@ xb = cw.getRemB();@b@ }@b@@b@ switch (quadrant)@b@ {@b@ case 0:@b@ return cosQ(xa, xb);@b@ case 1:@b@ return (-sinQ(xa, xb));@b@ case 2:@b@ return (-cosQ(xa, xb));@b@ case 3:@b@ return sinQ(xa, xb);@b@ }@b@ return (0.0D / 0.0D);@b@ }@b@@b@ public static double tan(double x)@b@ {@b@ double result;@b@ boolean negative = false;@b@ int quadrant = 0;@b@@b@ double xa = x;@b@ if (x < 0.0D) {@b@ negative = true;@b@ xa = -xa;@b@ }@b@@b@ if (xa == 0.0D) {@b@ long bits = Double.doubleToRawLongBits(x);@b@ if (bits < 0L)@b@ return -0.0D;@b@@b@ return 0.0D;@b@ }@b@@b@ if ((xa != xa) || (xa == (1.0D / 0.0D))) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ double xb = 0.0D;@b@ if (xa > 3294198.0D)@b@ {@b@ double[] reduceResults = new double[3];@b@ reducePayneHanek(xa, reduceResults);@b@ quadrant = (int)reduceResults[0] & 0x3;@b@ xa = reduceResults[1];@b@ xb = reduceResults[2];@b@ } else if (xa > 1.570796326794897D) {@b@ CodyWaite cw = new CodyWaite(xa);@b@ quadrant = cw.getK() & 0x3;@b@ xa = cw.getRemA();@b@ xb = cw.getRemB();@b@ }@b@@b@ if (xa > 1.5D)@b@ {@b@ double pi2a = 1.570796326794897D;@b@ double pi2b = 6.123233995736766E-017D;@b@@b@ double a = 1.570796326794897D - xa;@b@ double b = -(a - 1.570796326794897D + xa);@b@ b += 6.123233995736766E-017D - xb;@b@@b@ xa = a + b;@b@ xb = -(xa - a - b);@b@ quadrant ^= 1;@b@ negative ^= true;@b@ }@b@@b@ if ((quadrant & 0x1) == 0)@b@ result = tanQ(xa, xb, false);@b@ else {@b@ result = -tanQ(xa, xb, true);@b@ }@b@@b@ if (negative) {@b@ result = -result;@b@ }@b@@b@ return result;@b@ }@b@@b@ public static double atan(double x)@b@ {@b@ return atan(x, 0.0D, false);@b@ }@b@@b@ private static double atan(double xa, double xb, boolean leftPlane)@b@ {@b@ int idx;@b@ boolean negate = false;@b@@b@ if (xa == 0.0D) {@b@ return ((leftPlane) ? copySign(3.141592653589793D, xa) : xa);@b@ }@b@@b@ if (xa < 0.0D)@b@ {@b@ xa = -xa;@b@ xb = -xb;@b@ negate = true;@b@ }@b@@b@ if (xa > 16331239353195370.0D) {@b@ return (((negate ^ leftPlane)) ? -1.570796326794897D : 1.570796326794897D);@b@ }@b@@b@ if (xa < 1.0D) {@b@ idx = (int)((-1.716814692820414D * xa * xa + 8.0D) * xa + 0.5D);@b@ } else {@b@ double oneOverXa = 1.0D / xa;@b@ idx = (int)(-((-1.716814692820414D * oneOverXa * oneOverXa + 8.0D) * oneOverXa) + 13.07D);@b@ }@b@ double epsA = xa - TANGENT_TABLE_A[idx];@b@ double epsB = -(epsA - xa + TANGENT_TABLE_A[idx]);@b@ epsB += xb - TANGENT_TABLE_B[idx];@b@@b@ double temp = epsA + epsB;@b@ epsB = -(temp - epsA - epsB);@b@ epsA = temp;@b@@b@ temp = xa * 1073741824.0D;@b@ double ya = xa + temp - temp;@b@ double yb = xb + xa - ya;@b@ xa = ya;@b@ xb += yb;@b@@b@ if (idx == 0)@b@ {@b@ double denom = 1.0D / (1.0D + (xa + xb) * (TANGENT_TABLE_A[idx] + TANGENT_TABLE_B[idx]));@b@@b@ ya = epsA * denom;@b@ yb = epsB * denom;@b@ } else {@b@ double temp2 = xa * TANGENT_TABLE_A[idx];@b@ za = 1.0D + temp2;@b@ zb = -(za - 1.0D - temp2);@b@ temp2 = xb * TANGENT_TABLE_A[idx] + xa * TANGENT_TABLE_B[idx];@b@ temp = za + temp2;@b@ zb += -(temp - za - temp2);@b@ za = temp;@b@@b@ zb += xb * TANGENT_TABLE_B[idx];@b@ ya = epsA / za;@b@@b@ temp = ya * 1073741824.0D;@b@ double yaa = ya + temp - temp;@b@ double yab = ya - yaa;@b@@b@ temp = za * 1073741824.0D;@b@ double zaa = za + temp - temp;@b@ double zab = za - zaa;@b@@b@ yb = (epsA - yaa * zaa - yaa * zab - yab * zaa - yab * zab) / za;@b@@b@ yb += -epsA * zb / za / za;@b@ yb += epsB / za;@b@ }@b@@b@ epsA = ya;@b@ epsB = yb;@b@@b@ double epsA2 = epsA * epsA;@b@@b@ yb = 0.07490822288864472D;@b@ yb = yb * epsA2 + -0.09088450866185192D;@b@ yb = yb * epsA2 + 0.1111109594231331D;@b@ yb = yb * epsA2 + -0.1428571423679182D;@b@ yb = yb * epsA2 + 0.1999999999992358D;@b@ yb = yb * epsA2 + -0.3333333333333329D;@b@ yb = yb * epsA2 * epsA;@b@@b@ ya = epsA;@b@@b@ temp = ya + yb;@b@ yb = -(temp - ya - yb);@b@ ya = temp;@b@@b@ yb += epsB / (1.0D + epsA * epsA);@b@@b@ double za = EIGHTHS[idx] + ya;@b@ double zb = -(za - EIGHTHS[idx] - ya);@b@ temp = za + yb;@b@ zb += -(temp - za - yb);@b@ za = temp;@b@@b@ double result = za + zb;@b@ double resultb = -(result - za - zb);@b@@b@ if (leftPlane)@b@ {@b@ double pia = 3.141592653589793D;@b@ double pib = 1.224646799147353E-016D;@b@@b@ za = 3.141592653589793D - result;@b@ zb = -(za - 3.141592653589793D + result);@b@ zb += 1.224646799147353E-016D - resultb;@b@@b@ result = za + zb;@b@ resultb = -(result - za - zb);@b@ }@b@@b@ if ((negate ^ leftPlane)) {@b@ result = -result;@b@ }@b@@b@ return result;@b@ }@b@@b@ public static double atan2(double y, double x)@b@ {@b@ if ((x != x) || (y != y)) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if (y == 0.0D) {@b@ double result = x * y;@b@ double invx = 1.0D / x;@b@ double invy = 1.0D / y;@b@@b@ if (invx == 0.0D) {@b@ if (x > 0.0D)@b@ return y;@b@@b@ return copySign(3.141592653589793D, y);@b@ }@b@@b@ if ((x < 0.0D) || (invx < 0.0D)) {@b@ if ((y < 0.0D) || (invy < 0.0D))@b@ return -3.141592653589793D;@b@@b@ return 3.141592653589793D;@b@ }@b@@b@ return result;@b@ }@b@@b@ if (y == (1.0D / 0.0D)) {@b@ if (x == (1.0D / 0.0D)) {@b@ return 0.7853981633974483D;@b@ }@b@@b@ if (x == (-1.0D / 0.0D)) {@b@ return 2.356194490192345D;@b@ }@b@@b@ return 1.570796326794897D;@b@ }@b@@b@ if (y == (-1.0D / 0.0D)) {@b@ if (x == (1.0D / 0.0D)) {@b@ return -0.7853981633974483D;@b@ }@b@@b@ if (x == (-1.0D / 0.0D)) {@b@ return -2.356194490192345D;@b@ }@b@@b@ return -1.570796326794897D;@b@ }@b@@b@ if (x == (1.0D / 0.0D)) {@b@ if ((y > 0.0D) || (1.0D / y > 0.0D)) {@b@ return 0.0D;@b@ }@b@@b@ if ((y < 0.0D) || (1.0D / y < 0.0D))@b@ return -0.0D;@b@@b@ }@b@@b@ if (x == (-1.0D / 0.0D))@b@ {@b@ if ((y > 0.0D) || (1.0D / y > 0.0D)) {@b@ return 3.141592653589793D;@b@ }@b@@b@ if ((y < 0.0D) || (1.0D / y < 0.0D)) {@b@ return -3.141592653589793D;@b@ }@b@@b@ }@b@@b@ if (x == 0.0D) {@b@ if ((y > 0.0D) || (1.0D / y > 0.0D)) {@b@ return 1.570796326794897D;@b@ }@b@@b@ if ((y < 0.0D) || (1.0D / y < 0.0D)) {@b@ return -1.570796326794897D;@b@ }@b@@b@ }@b@@b@ double r = y / x;@b@ if (Double.isInfinite(r)) {@b@ return atan(r, 0.0D, x < 0.0D);@b@ }@b@@b@ double ra = doubleHighPart(r);@b@ double rb = r - ra;@b@@b@ double xa = doubleHighPart(x);@b@ double xb = x - xa;@b@@b@ rb += (y - ra * xa - ra * xb - rb * xa - rb * xb) / x;@b@@b@ double temp = ra + rb;@b@ rb = -(temp - ra - rb);@b@ ra = temp;@b@@b@ if (ra == 0.0D) {@b@ ra = copySign(0.0D, y);@b@ }@b@@b@ double result = atan(ra, rb, x < 0.0D);@b@@b@ return result;@b@ }@b@@b@ public static double asin(double x)@b@ {@b@ if (x != x) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if ((x > 1.0D) || (x < -1.0D)) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if (x == 1.0D) {@b@ return 1.570796326794897D;@b@ }@b@@b@ if (x == -1.0D) {@b@ return -1.570796326794897D;@b@ }@b@@b@ if (x == 0.0D) {@b@ return x;@b@ }@b@@b@ double temp = x * 1073741824.0D;@b@ double xa = x + temp - temp;@b@ double xb = x - xa;@b@@b@ double ya = xa * xa;@b@ double yb = xa * xb * 2.0D + xb * xb;@b@@b@ ya = -ya;@b@ yb = -yb;@b@@b@ double za = 1.0D + ya;@b@ double zb = -(za - 1.0D - ya);@b@@b@ temp = za + yb;@b@ zb += -(temp - za - yb);@b@ za = temp;@b@@b@ double y = sqrt(za);@b@ temp = y * 1073741824.0D;@b@ ya = y + temp - temp;@b@ yb = y - ya;@b@@b@ yb += (za - ya * ya - 2.0D * ya * yb - yb * yb) / 2.0D * y;@b@@b@ double dx = zb / 2.0D * y;@b@@b@ double r = x / y;@b@ temp = r * 1073741824.0D;@b@ double ra = r + temp - temp;@b@ double rb = r - ra;@b@@b@ rb += (x - ra * ya - ra * yb - rb * ya - rb * yb) / y;@b@ rb += -x * dx / y / y;@b@@b@ temp = ra + rb;@b@ rb = -(temp - ra - rb);@b@ ra = temp;@b@@b@ return atan(ra, rb, false);@b@ }@b@@b@ public static double acos(double x)@b@ {@b@ if (x != x) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if ((x > 1.0D) || (x < -1.0D)) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ if (x == -1.0D) {@b@ return 3.141592653589793D;@b@ }@b@@b@ if (x == 1.0D) {@b@ return 0.0D;@b@ }@b@@b@ if (x == 0.0D) {@b@ return 1.570796326794897D;@b@ }@b@@b@ double temp = x * 1073741824.0D;@b@ double xa = x + temp - temp;@b@ double xb = x - xa;@b@@b@ double ya = xa * xa;@b@ double yb = xa * xb * 2.0D + xb * xb;@b@@b@ ya = -ya;@b@ yb = -yb;@b@@b@ double za = 1.0D + ya;@b@ double zb = -(za - 1.0D - ya);@b@@b@ temp = za + yb;@b@ zb += -(temp - za - yb);@b@ za = temp;@b@@b@ double y = sqrt(za);@b@ temp = y * 1073741824.0D;@b@ ya = y + temp - temp;@b@ yb = y - ya;@b@@b@ yb += (za - ya * ya - 2.0D * ya * yb - yb * yb) / 2.0D * y;@b@@b@ yb += zb / 2.0D * y;@b@ y = ya + yb;@b@ yb = -(y - ya - yb);@b@@b@ double r = y / x;@b@@b@ if (Double.isInfinite(r)) {@b@ return 1.570796326794897D;@b@ }@b@@b@ double ra = doubleHighPart(r);@b@ double rb = r - ra;@b@@b@ rb += (y - ra * xa - ra * xb - rb * xa - rb * xb) / x;@b@ rb += yb / x;@b@@b@ temp = ra + rb;@b@ rb = -(temp - ra - rb);@b@ ra = temp;@b@@b@ return atan(ra, rb, x < 0.0D);@b@ }@b@@b@ public static double cbrt(double x)@b@ {@b@ long inbits = Double.doubleToRawLongBits(x);@b@ int exponent = (int)(inbits >> 52 & 0x7FF) - 1023;@b@ boolean subnormal = false;@b@@b@ if (exponent == -1023) {@b@ if (x == 0.0D) {@b@ return x;@b@ }@b@@b@ subnormal = true;@b@ x *= 18014398509481984.0D;@b@ inbits = Double.doubleToRawLongBits(x);@b@ exponent = (int)(inbits >> 52 & 0x7FF) - 1023;@b@ }@b@@b@ if (exponent == 1024)@b@ {@b@ return x;@b@ }@b@@b@ int exp3 = exponent / 3;@b@@b@ double p2 = Double.longBitsToDouble(inbits & 0x0 | (exp3 + 1023 & 0x7FF) << 52);@b@@b@ double mant = Double.longBitsToDouble(inbits & 0xFFFFFFFF | 0x0);@b@@b@ double est = -0.01071469073319593D;@b@ est = est * mant + 0.0875862700108075D;@b@ est = est * mant + -0.3058015757857271D;@b@ est = est * mant + 0.7249995199969751D;@b@ est = est * mant + 0.5039018405998234D;@b@@b@ est *= CBRTTWO[(exponent % 3 + 2)];@b@@b@ double xs = x / p2 * p2 * p2;@b@ est += (xs - est * est * est) / 3.0D * est * est;@b@ est += (xs - est * est * est) / 3.0D * est * est;@b@@b@ double temp = est * 1073741824.0D;@b@ double ya = est + temp - temp;@b@ double yb = est - ya;@b@@b@ double za = ya * ya;@b@ double zb = ya * yb * 2.0D + yb * yb;@b@ temp = za * 1073741824.0D;@b@ double temp2 = za + temp - temp;@b@ zb += za - temp2;@b@ za = temp2;@b@@b@ zb = za * yb + ya * zb + zb * yb;@b@ za *= ya;@b@@b@ double na = xs - za;@b@ double nb = -(na - xs + za);@b@ nb -= zb;@b@@b@ est += (na + nb) / 3.0D * est * est;@b@@b@ est *= p2;@b@@b@ if (subnormal) {@b@ est *= 3.814697265625E-006D;@b@ }@b@@b@ return est;@b@ }@b@@b@ public static double toRadians(double x)@b@ {@b@ if ((Double.isInfinite(x)) || (x == 0.0D)) {@b@ return x;@b@ }@b@@b@ double facta = 0.01745329052209854D;@b@ double factb = 1.997844754509471E-009D;@b@@b@ double xa = doubleHighPart(x);@b@ double xb = x - xa;@b@@b@ double result = xb * 1.997844754509471E-009D + xb * 0.01745329052209854D + xa * 1.997844754509471E-009D + xa * 0.01745329052209854D;@b@ if (result == 0.0D)@b@ result *= x;@b@@b@ return result;@b@ }@b@@b@ public static double toDegrees(double x)@b@ {@b@ if ((Double.isInfinite(x)) || (x == 0.0D)) {@b@ return x;@b@ }@b@@b@ double facta = 57.2957763671875D;@b@ double factb = 3.145894820876798E-006D;@b@@b@ double xa = doubleHighPart(x);@b@ double xb = x - xa;@b@@b@ return (xb * 3.145894820876798E-006D + xb * 57.2957763671875D + xa * 3.145894820876798E-006D + xa * 57.2957763671875D);@b@ }@b@@b@ public static int abs(int x)@b@ {@b@ int i = x >>> 31;@b@ return ((x ^ (i ^ 0xFFFFFFFF) + 1) + i);@b@ }@b@@b@ public static long abs(long x)@b@ {@b@ long l = x >>> 63;@b@@b@ return ((x ^ (l ^ 0xFFFFFFFF) + 1L) + l);@b@ }@b@@b@ public static float abs(float x)@b@ {@b@ return Float.intBitsToFloat(0x7FFFFFFF & Float.floatToRawIntBits(x));@b@ }@b@@b@ public static double abs(double x)@b@ {@b@ return Double.longBitsToDouble(0xFFFFFFFF & Double.doubleToRawLongBits(x));@b@ }@b@@b@ public static double ulp(double x)@b@ {@b@ if (Double.isInfinite(x))@b@ return (1.0D / 0.0D);@b@@b@ return abs(x - Double.longBitsToDouble(Double.doubleToRawLongBits(x) ^ 1L));@b@ }@b@@b@ public static float ulp(float x)@b@ {@b@ if (Float.isInfinite(x))@b@ return (1.0F / 1.0F);@b@@b@ return abs(x - Float.intBitsToFloat(Float.floatToIntBits(x) ^ 0x1));@b@ }@b@@b@ public static double scalb(double d, int n)@b@ {@b@ if ((n > -1023) && (n < 1024)) {@b@ return (d * Double.longBitsToDouble(n + 1023 << 52));@b@ }@b@@b@ if ((Double.isNaN(d)) || (Double.isInfinite(d)) || (d == 0.0D))@b@ return d;@b@@b@ if (n < -2098)@b@ return ((d > 0.0D) ? 0.0D : -0.0D);@b@@b@ if (n > 2097) {@b@ return ((d > 0.0D) ? (1.0D / 0.0D) : (-1.0D / 0.0D));@b@ }@b@@b@ long bits = Double.doubleToRawLongBits(d);@b@ long sign = bits & 0x0;@b@ int exponent = (int)(bits >>> 52) & 0x7FF;@b@ long mantissa = bits & 0xFFFFFFFF;@b@@b@ int scaledExponent = exponent + n;@b@@b@ if (n < 0)@b@ {@b@ if (scaledExponent > 0)@b@ {@b@ return Double.longBitsToDouble(sign | scaledExponent << 52 | mantissa); }@b@ if (scaledExponent > -53)@b@ {@b@ mantissa |= 4503599627370496L;@b@@b@ long mostSignificantLostBit = mantissa & 1L << -scaledExponent;@b@ mantissa >>>= 1 - scaledExponent;@b@ if (mostSignificantLostBit != 0L)@b@ {@b@ mantissa += 1L;@b@ }@b@ return Double.longBitsToDouble(sign | mantissa);@b@ }@b@@b@ return ((sign == 0L) ? 0.0D : -0.0D);@b@ }@b@@b@ if (exponent == 0)@b@ {@b@ while (mantissa >>> 52 != 1L) {@b@ mantissa <<= 1;@b@ --scaledExponent;@b@ }@b@ ++scaledExponent;@b@ mantissa &= 4503599627370495L;@b@@b@ if (scaledExponent < 2047)@b@ return Double.longBitsToDouble(sign | scaledExponent << 52 | mantissa);@b@@b@ return ((sign == 0L) ? (1.0D / 0.0D) : (-1.0D / 0.0D));@b@ }@b@@b@ if (scaledExponent < 2047)@b@ return Double.longBitsToDouble(sign | scaledExponent << 52 | mantissa);@b@@b@ return ((sign == 0L) ? (1.0D / 0.0D) : (-1.0D / 0.0D));@b@ }@b@@b@ public static float scalb(float f, int n)@b@ {@b@ if ((n > -127) && (n < 128)) {@b@ return (f * Float.intBitsToFloat(n + 127 << 23));@b@ }@b@@b@ if ((Float.isNaN(f)) || (Float.isInfinite(f)) || (f == 0.0F))@b@ return f;@b@@b@ if (n < -277)@b@ return ((f > 0.0F) ? 0.0F : -0.0F);@b@@b@ if (n > 276) {@b@ return ((f > 0.0F) ? (1.0F / 1.0F) : (1.0F / -1.0F));@b@ }@b@@b@ int bits = Float.floatToIntBits(f);@b@ int sign = bits & 0x80000000;@b@ int exponent = bits >>> 23 & 0xFF;@b@ int mantissa = bits & 0x7FFFFF;@b@@b@ int scaledExponent = exponent + n;@b@@b@ if (n < 0)@b@ {@b@ if (scaledExponent > 0)@b@ {@b@ return Float.intBitsToFloat(sign | scaledExponent << 23 | mantissa); }@b@ if (scaledExponent > -24)@b@ {@b@ mantissa |= 8388608;@b@@b@ int mostSignificantLostBit = mantissa & 1 << -scaledExponent;@b@ mantissa >>>= 1 - scaledExponent;@b@ if (mostSignificantLostBit != 0)@b@ {@b@ ++mantissa;@b@ }@b@ return Float.intBitsToFloat(sign | mantissa);@b@ }@b@@b@ return ((sign == 0) ? 0.0F : -0.0F);@b@ }@b@@b@ if (exponent == 0)@b@ {@b@ while (mantissa >>> 23 != 1) {@b@ mantissa <<= 1;@b@ --scaledExponent;@b@ }@b@ ++scaledExponent;@b@ mantissa &= 8388607;@b@@b@ if (scaledExponent < 255)@b@ return Float.intBitsToFloat(sign | scaledExponent << 23 | mantissa);@b@@b@ return ((sign == 0) ? (1.0F / 1.0F) : (1.0F / -1.0F));@b@ }@b@@b@ if (scaledExponent < 255)@b@ return Float.intBitsToFloat(sign | scaledExponent << 23 | mantissa);@b@@b@ return ((sign == 0) ? (1.0F / 1.0F) : (1.0F / -1.0F));@b@ }@b@@b@ public static double nextAfter(double d, double direction)@b@ {@b@ if ((Double.isNaN(d)) || (Double.isNaN(direction)))@b@ return (0.0D / 0.0D);@b@ if (d == direction)@b@ return direction;@b@ if (Double.isInfinite(d))@b@ return ((d < 0.0D) ? -1.797693134862316E+308D : 1.7976931348623157E+308D);@b@ if (d == 0.0D) {@b@ return ((direction < 0.0D) ? -4.940656458412465E-324D : 4.9E-324D);@b@ }@b@@b@ long bits = Double.doubleToRawLongBits(d);@b@ long sign = bits & 0x0;@b@ if ((((direction < d) ? 1 : 0) ^ ((sign == 0L) ? 1 : 0)) != 0)@b@ return Double.longBitsToDouble(sign | (bits & 0xFFFFFFFF) + 1L);@b@@b@ return Double.longBitsToDouble(sign | (bits & 0xFFFFFFFF) - 1L);@b@ }@b@@b@ public static float nextAfter(float f, double direction)@b@ {@b@ if ((Double.isNaN(f)) || (Double.isNaN(direction)))@b@ return (0.0F / 0.0F);@b@ if (f == direction)@b@ return (float)direction;@b@ if (Float.isInfinite(f))@b@ return ((f < 0.0F) ? -3.402824E+038F : 3.4028235E+38F);@b@ if (f == 0.0F) {@b@ return ((direction < 0.0D) ? -1.401299E-045F : 1.4E-45F);@b@ }@b@@b@ int bits = Float.floatToIntBits(f);@b@ int sign = bits & 0x80000000;@b@ if ((((direction < f) ? 1 : 0) ^ ((sign == 0) ? 1 : 0)) != 0)@b@ return Float.intBitsToFloat(sign | (bits & 0x7FFFFFFF) + 1);@b@@b@ return Float.intBitsToFloat(sign | (bits & 0x7FFFFFFF) - 1);@b@ }@b@@b@ public static double floor(double x)@b@ {@b@ if (x != x) {@b@ return x;@b@ }@b@@b@ if ((x >= 4503599627370496.0D) || (x <= -4503599627370496.0D)) {@b@ return x;@b@ }@b@@b@ long y = ()x;@b@ if ((x < 0.0D) && (y != x)) {@b@ y -= 1L;@b@ }@b@@b@ if (y == 0L) {@b@ return (x * y);@b@ }@b@@b@ return y;@b@ }@b@@b@ public static double ceil(double x)@b@ {@b@ if (x != x) {@b@ return x;@b@ }@b@@b@ double y = floor(x);@b@ if (y == x) {@b@ return y;@b@ }@b@@b@ y += 1.0D;@b@@b@ if (y == 0.0D) {@b@ return (x * y);@b@ }@b@@b@ return y;@b@ }@b@@b@ public static double rint(double x)@b@ {@b@ double y = floor(x);@b@ double d = x - y;@b@@b@ if (d > 0.5D) {@b@ if (y == -1.0D)@b@ return -0.0D;@b@@b@ return (y + 1.0D);@b@ }@b@ if (d < 0.5D) {@b@ return y;@b@ }@b@@b@ long z = ()y;@b@ return (((z & 1L) == 0L) ? y : y + 1.0D);@b@ }@b@@b@ public static long round(double x)@b@ {@b@ return ()floor(x + 0.5D);@b@ }@b@@b@ public static int round(float x)@b@ {@b@ return (int)floor(x + 0.5F);@b@ }@b@@b@ public static int min(int a, int b)@b@ {@b@ return ((a <= b) ? a : b);@b@ }@b@@b@ public static long min(long a, long b)@b@ {@b@ return ((a <= b) ? a : b);@b@ }@b@@b@ public static float min(float a, float b)@b@ {@b@ if (a > b)@b@ return b;@b@@b@ if (a < b) {@b@ return a;@b@ }@b@@b@ if (a != b) {@b@ return (0.0F / 0.0F);@b@ }@b@@b@ int bits = Float.floatToRawIntBits(a);@b@ if (bits == -2147483648)@b@ return a;@b@@b@ return b;@b@ }@b@@b@ public static double min(double a, double b)@b@ {@b@ if (a > b)@b@ return b;@b@@b@ if (a < b) {@b@ return a;@b@ }@b@@b@ if (a != b) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ long bits = Double.doubleToRawLongBits(a);@b@ if (bits == -9223372036854775808L)@b@ return a;@b@@b@ return b;@b@ }@b@@b@ public static int max(int a, int b)@b@ {@b@ return ((a <= b) ? b : a);@b@ }@b@@b@ public static long max(long a, long b)@b@ {@b@ return ((a <= b) ? b : a);@b@ }@b@@b@ public static float max(float a, float b)@b@ {@b@ if (a > b)@b@ return a;@b@@b@ if (a < b) {@b@ return b;@b@ }@b@@b@ if (a != b) {@b@ return (0.0F / 0.0F);@b@ }@b@@b@ int bits = Float.floatToRawIntBits(a);@b@ if (bits == -2147483648)@b@ return b;@b@@b@ return a;@b@ }@b@@b@ public static double max(double a, double b)@b@ {@b@ if (a > b)@b@ return a;@b@@b@ if (a < b) {@b@ return b;@b@ }@b@@b@ if (a != b) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ long bits = Double.doubleToRawLongBits(a);@b@ if (bits == -9223372036854775808L)@b@ return b;@b@@b@ return a;@b@ }@b@@b@ public static double hypot(double x, double y)@b@ {@b@ if ((Double.isInfinite(x)) || (Double.isInfinite(y)))@b@ return (1.0D / 0.0D);@b@ if ((Double.isNaN(x)) || (Double.isNaN(y))) {@b@ return (0.0D / 0.0D);@b@ }@b@@b@ int expX = getExponent(x);@b@ int expY = getExponent(y);@b@ if (expX > expY + 27)@b@ {@b@ return abs(x); }@b@ if (expY > expX + 27)@b@ {@b@ return abs(y);@b@ }@b@@b@ int middleExp = (expX + expY) / 2;@b@@b@ double scaledX = scalb(x, -middleExp);@b@ double scaledY = scalb(y, -middleExp);@b@@b@ double scaledH = sqrt(scaledX * scaledX + scaledY * scaledY);@b@@b@ return scalb(scaledH, middleExp);@b@ }@b@@b@ public static double IEEEremainder(double dividend, double divisor)@b@ {@b@ return StrictMath.IEEEremainder(dividend, divisor);@b@ }@b@@b@ public static double copySign(double magnitude, double sign)@b@ {@b@ long m = Double.doubleToRawLongBits(magnitude);@b@ long s = Double.doubleToRawLongBits(sign);@b@ if ((m ^ s) >= 0L)@b@ return magnitude;@b@@b@ return (-magnitude);@b@ }@b@@b@ public static float copySign(float magnitude, float sign)@b@ {@b@ int m = Float.floatToRawIntBits(magnitude);@b@ int s = Float.floatToRawIntBits(sign);@b@ if ((m ^ s) >= 0)@b@ return magnitude;@b@@b@ return (-magnitude);@b@ }@b@@b@ public static int getExponent(double d)@b@ {@b@ return ((int)(Double.doubleToRawLongBits(d) >>> 52 & 0x7FF) - 1023);@b@ }@b@@b@ public static int getExponent(float f)@b@ {@b@ return ((Float.floatToRawIntBits(f) >>> 23 & 0xFF) - 127);@b@ }@b@@b@ public static void main(String[] a)@b@ {@b@ PrintStream out = System.out;@b@ FastMathCalc.printarray(out, "EXP_INT_TABLE_A", 1500, ExpIntTable.access$000());@b@ FastMathCalc.printarray(out, "EXP_INT_TABLE_B", 1500, ExpIntTable.access$100());@b@ FastMathCalc.printarray(out, "EXP_FRAC_TABLE_A", 1025, ExpFracTable.access$200());@b@ FastMathCalc.printarray(out, "EXP_FRAC_TABLE_B", 1025, ExpFracTable.access$300());@b@ FastMathCalc.printarray(out, "LN_MANT", 1024, lnMant.access$400());@b@ FastMathCalc.printarray(out, "SINE_TABLE_A", 14, SINE_TABLE_A);@b@ FastMathCalc.printarray(out, "SINE_TABLE_B", 14, SINE_TABLE_B);@b@ FastMathCalc.printarray(out, "COSINE_TABLE_A", 14, COSINE_TABLE_A);@b@ FastMathCalc.printarray(out, "COSINE_TABLE_B", 14, COSINE_TABLE_B);@b@ FastMathCalc.printarray(out, "TANGENT_TABLE_A", 14, TANGENT_TABLE_A);@b@ FastMathCalc.printarray(out, "TANGENT_TABLE_B", 14, TANGENT_TABLE_B);@b@ }@b@@b@ private static class CodyWaite@b@ {@b@ private final int finalK;@b@ private final double finalRemA;@b@ private final double finalRemB;@b@@b@ CodyWaite(double xa)@b@ {@b@ int k = (int)(xa * 0.6366197723675814D);@b@ while (true)@b@ {@b@ double a = -k * 1.570796251296997D;@b@ remA = xa + a;@b@ remB = -(remA - xa - a);@b@@b@ a = -k * 7.549789948768648E-008D;@b@ double b = remA;@b@ remA = a + b;@b@ remB += -(remA - b - a);@b@@b@ a = -k * 6.123233995736766E-017D;@b@ b = remA;@b@ remA = a + b;@b@ remB += -(remA - b - a);@b@@b@ if (remA > 0.0D) {@b@ break;@b@ }@b@@b@ --k;@b@ }@b@@b@ this.finalK = k;@b@ this.finalRemA = remA;@b@ this.finalRemB = remB;@b@ }@b@@b@ int getK()@b@ {@b@ return this.finalK;@b@ }@b@@b@ double getRemA()@b@ {@b@ return this.finalRemA;@b@ }@b@@b@ double getRemB()@b@ {@b@ return this.finalRemB;@b@ }@b@ }@b@@b@ private static class lnMant@b@ {@b@ private static final double[][] LN_MANT = FastMathLiteralArrays.loadLnMant();@b@ }@b@@b@ private static class ExpFracTable@b@ {@b@ private static final double[] EXP_FRAC_TABLE_A = FastMathLiteralArrays.loadExpFracA();@b@ private static final double[] EXP_FRAC_TABLE_B = FastMathLiteralArrays.loadExpFracB();@b@ }@b@@b@ private static class ExpIntTable@b@ {@b@ private static final double[] EXP_INT_TABLE_A = FastMathLiteralArrays.loadExpIntA();@b@ private static final double[] EXP_INT_TABLE_B = FastMathLiteralArrays.loadExpIntB();@b@ }@b@}