![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to softfloat_fpsp.c CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [Previous NeXT emulator] / previous / src / softfloat|
Return to softfloat_fpsp.c CVS log | Up to [Previous NeXT emulator] / previous / src / softfloat |
1.1 ! root 1: ! 2: /*============================================================================ ! 3: ! 4: This C source file is an extension to the SoftFloat IEC/IEEE Floating-point ! 5: Arithmetic Package, Release 2a. ! 6: ! 7: Written by Andreas Grabher for Previous, NeXT Computer Emulator. ! 8: ! 9: =============================================================================*/ ! 10: ! 11: #include "softfloat.h" ! 12: #include "softfloat_fpsp_tables.h" ! 13: ! 14: ! 15: /*---------------------------------------------------------------------------- ! 16: | Algorithms for transcendental functions supported by MC68881 and MC68882 ! 17: | mathematical coprocessors. The functions are derived from FPSP library. ! 18: *----------------------------------------------------------------------------*/ ! 19: ! 20: #define pi_sig LIT64(0xc90fdaa22168c235) ! 21: #define pi_sig0 LIT64(0xc90fdaa22168c234) ! 22: #define pi_sig1 LIT64(0xc4c6628b80dc1cd1) ! 23: ! 24: #define pi_exp 0x4000 ! 25: #define piby2_exp 0x3FFF ! 26: #define piby4_exp 0x3FFE ! 27: ! 28: #define one_exp 0x3FFF ! 29: #define one_sig LIT64(0x8000000000000000) ! 30: ! 31: ! 32: /*---------------------------------------------------------------------------- ! 33: | Function for compactifying extended double-precision floating point values. ! 34: *----------------------------------------------------------------------------*/ ! 35: ! 36: int32 floatx80_make_compact(int32 aExp, bits64 aSig) ! 37: { ! 38: return (aExp<<16)|(aSig>>48); ! 39: } ! 40: ! 41: ! 42: /*---------------------------------------------------------------------------- ! 43: | Arc cosine ! 44: *----------------------------------------------------------------------------*/ ! 45: ! 46: floatx80 floatx80_acos(floatx80 a) ! 47: { ! 48: flag aSign; ! 49: int32 aExp; ! 50: bits64 aSig; ! 51: ! 52: int8 user_rnd_mode, user_rnd_prec; ! 53: ! 54: int32 compact; ! 55: floatx80 fp0, fp1, one; ! 56: ! 57: aSig = extractFloatx80Frac(a); ! 58: aExp = extractFloatx80Exp(a); ! 59: aSign = extractFloatx80Sign(a); ! 60: ! 61: if (aExp == 0x7FFF && (bits64) (aSig<<1)) { ! 62: return propagateFloatx80NaNOneArg(a); ! 63: } ! 64: if (aExp == 0 && aSig == 0) { ! 65: float_raise(float_flag_inexact); ! 66: return roundAndPackFloatx80(floatx80_rounding_precision, 0, piby2_exp, pi_sig, 0); ! 67: } ! 68: ! 69: compact = floatx80_make_compact(aExp, aSig); ! 70: ! 71: if (compact >= 0x3FFF8000) { // |X| >= 1 ! 72: if (aExp == one_exp && aSig == one_sig) { // |X| == 1 ! 73: if (aSign) { // X == -1 ! 74: a = packFloatx80(0, pi_exp, pi_sig); ! 75: float_raise(float_flag_inexact); ! 76: return floatx80_move(a); ! 77: } else { // X == +1 ! 78: return packFloatx80(0, 0, 0); ! 79: } ! 80: } else { // |X| > 1 ! 81: float_raise(float_flag_invalid); ! 82: a.low = floatx80_default_nan_low; ! 83: a.high = floatx80_default_nan_high; ! 84: return a; ! 85: } ! 86: } // |X| < 1 ! 87: ! 88: user_rnd_mode = float_rounding_mode; ! 89: user_rnd_prec = floatx80_rounding_precision; ! 90: float_rounding_mode = float_round_nearest_even; ! 91: floatx80_rounding_precision = 80; ! 92: ! 93: one = packFloatx80(0, one_exp, one_sig); ! 94: fp0 = a; ! 95: ! 96: fp1 = floatx80_add(one, fp0); // 1 + X ! 97: fp0 = floatx80_sub(one, fp0); // 1 - X ! 98: fp0 = floatx80_div(fp0, fp1); // (1-X)/(1+X) ! 99: fp0 = floatx80_sqrt(fp0); // SQRT((1-X)/(1+X)) ! 100: fp0 = floatx80_atan(fp0); // ATAN(SQRT((1-X)/(1+X))) ! 101: ! 102: float_rounding_mode = user_rnd_mode; ! 103: floatx80_rounding_precision = user_rnd_prec; ! 104: ! 105: a = floatx80_add(fp0, fp0); // 2 * ATAN(SQRT((1-X)/(1+X))) ! 106: ! 107: float_raise(float_flag_inexact); ! 108: ! 109: return a; ! 110: } ! 111: ! 112: /*---------------------------------------------------------------------------- ! 113: | Arc sine ! 114: *----------------------------------------------------------------------------*/ ! 115: ! 116: floatx80 floatx80_asin(floatx80 a) ! 117: { ! 118: flag aSign; ! 119: int32 aExp; ! 120: bits64 aSig; ! 121: ! 122: int8 user_rnd_mode, user_rnd_prec; ! 123: ! 124: int32 compact; ! 125: floatx80 fp0, fp1, fp2, one; ! 126: ! 127: aSig = extractFloatx80Frac(a); ! 128: aExp = extractFloatx80Exp(a); ! 129: aSign = extractFloatx80Sign(a); ! 130: ! 131: if (aExp == 0x7FFF && (bits64) (aSig<<1)) { ! 132: return propagateFloatx80NaNOneArg(a); ! 133: } ! 134: ! 135: if (aExp == 0 && aSig == 0) { ! 136: return packFloatx80(aSign, 0, 0); ! 137: } ! 138: ! 139: compact = floatx80_make_compact(aExp, aSig); ! 140: ! 141: if (compact >= 0x3FFF8000) { // |X| >= 1 ! 142: if (aExp == one_exp && aSig == one_sig) { // |X| == 1 ! 143: float_raise(float_flag_inexact); ! 144: a = packFloatx80(aSign, piby2_exp, pi_sig); ! 145: return floatx80_move(a); ! 146: } else { // |X| > 1 ! 147: float_raise(float_flag_invalid); ! 148: a.low = floatx80_default_nan_low; ! 149: a.high = floatx80_default_nan_high; ! 150: return a; ! 151: } ! 152: ! 153: } // |X| < 1 ! 154: ! 155: user_rnd_mode = float_rounding_mode; ! 156: user_rnd_prec = floatx80_rounding_precision; ! 157: float_rounding_mode = float_round_nearest_even; ! 158: floatx80_rounding_precision = 80; ! 159: ! 160: one = packFloatx80(0, one_exp, one_sig); ! 161: fp0 = a; ! 162: ! 163: fp1 = floatx80_sub(one, fp0); // 1 - X ! 164: fp2 = floatx80_add(one, fp0); // 1 + X ! 165: fp1 = floatx80_mul(fp2, fp1); // (1+X)*(1-X) ! 166: fp1 = floatx80_sqrt(fp1); // SQRT((1+X)*(1-X)) ! 167: fp0 = floatx80_div(fp0, fp1); // X/SQRT((1+X)*(1-X)) ! 168: ! 169: float_rounding_mode = user_rnd_mode; ! 170: floatx80_rounding_precision = user_rnd_prec; ! 171: ! 172: a = floatx80_atan(fp0); // ATAN(X/SQRT((1+X)*(1-X))) ! 173: ! 174: float_raise(float_flag_inexact); ! 175: ! 176: return a; ! 177: } ! 178: ! 179: /*---------------------------------------------------------------------------- ! 180: | Arc tangent ! 181: *----------------------------------------------------------------------------*/ ! 182: ! 183: floatx80 floatx80_atan(floatx80 a) ! 184: { ! 185: flag aSign; ! 186: int32 aExp; ! 187: bits64 aSig; ! 188: ! 189: int8 user_rnd_mode, user_rnd_prec; ! 190: ! 191: int32 compact, tbl_index; ! 192: floatx80 fp0, fp1, fp2, fp3, xsave; ! 193: ! 194: aSig = extractFloatx80Frac(a); ! 195: aExp = extractFloatx80Exp(a); ! 196: aSign = extractFloatx80Sign(a); ! 197: ! 198: if (aExp == 0x7FFF) { ! 199: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 200: a = packFloatx80(aSign, piby2_exp, pi_sig); ! 201: float_raise(float_flag_inexact); ! 202: return floatx80_move(a); ! 203: } ! 204: ! 205: if (aExp == 0 && aSig == 0) { ! 206: return packFloatx80(aSign, 0, 0); ! 207: } ! 208: ! 209: compact = floatx80_make_compact(aExp, aSig); ! 210: ! 211: user_rnd_mode = float_rounding_mode; ! 212: user_rnd_prec = floatx80_rounding_precision; ! 213: float_rounding_mode = float_round_nearest_even; ! 214: floatx80_rounding_precision = 80; ! 215: ! 216: if (compact < 0x3FFB8000 || compact > 0x4002FFFF) { // |X| >= 16 or |X| < 1/16 ! 217: if (compact > 0x3FFF8000) { // |X| >= 16 ! 218: if (compact > 0x40638000) { // |X| > 2^(100) ! 219: fp0 = packFloatx80(aSign, piby2_exp, pi_sig); ! 220: fp1 = packFloatx80(aSign, 0x0001, one_sig); ! 221: ! 222: float_rounding_mode = user_rnd_mode; ! 223: floatx80_rounding_precision = user_rnd_prec; ! 224: ! 225: a = floatx80_sub(fp0, fp1); ! 226: ! 227: float_raise(float_flag_inexact); ! 228: ! 229: return a; ! 230: } else { ! 231: fp0 = a; ! 232: fp1 = packFloatx80(1, one_exp, one_sig); // -1 ! 233: fp1 = floatx80_div(fp1, fp0); // X' = -1/X ! 234: xsave = fp1; ! 235: fp0 = floatx80_mul(fp1, fp1); // Y = X'*X' ! 236: fp1 = floatx80_mul(fp0, fp0); // Z = Y*Y ! 237: fp3 = float64_to_floatx80(LIT64(0xBFB70BF398539E6A)); // C5 ! 238: fp2 = float64_to_floatx80(LIT64(0x3FBC7187962D1D7D)); // C4 ! 239: fp3 = floatx80_mul(fp3, fp1); // Z*C5 ! 240: fp2 = floatx80_mul(fp2, fp1); // Z*C4 ! 241: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBFC24924827107B8))); // C3+Z*C5 ! 242: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FC999999996263E))); // C2+Z*C4 ! 243: fp1 = floatx80_mul(fp1, fp3); // Z*(C3+Z*C5) ! 244: fp2 = floatx80_mul(fp2, fp0); // Y*(C2+Z*C4) ! 245: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0xBFD5555555555536))); // C1+Z*(C3+Z*C5) ! 246: fp0 = floatx80_mul(fp0, xsave); // X'*Y ! 247: fp1 = floatx80_add(fp1, fp2); // [Y*(C2+Z*C4)]+[C1+Z*(C3+Z*C5)] ! 248: fp0 = floatx80_mul(fp0, fp1); // X'*Y*([B1+Z*(B3+Z*B5)]+[Y*(B2+Z*(B4+Z*B6))]) ?? ! 249: fp0 = floatx80_add(fp0, xsave); ! 250: fp1 = packFloatx80(aSign, piby2_exp, pi_sig); ! 251: ! 252: float_rounding_mode = user_rnd_mode; ! 253: floatx80_rounding_precision = user_rnd_prec; ! 254: ! 255: a = floatx80_add(fp0, fp1); ! 256: ! 257: float_raise(float_flag_inexact); ! 258: ! 259: return a; ! 260: } ! 261: } else { // |X| < 1/16 ! 262: if (compact < 0x3FD78000) { // |X| < 2^(-40) ! 263: float_rounding_mode = user_rnd_mode; ! 264: floatx80_rounding_precision = user_rnd_prec; ! 265: ! 266: a = floatx80_move(a); ! 267: ! 268: float_raise(float_flag_inexact); ! 269: ! 270: return a; ! 271: } else { ! 272: fp0 = a; ! 273: xsave = a; ! 274: fp0 = floatx80_mul(fp0, fp0); // Y = X*X ! 275: fp1 = floatx80_mul(fp0, fp0); // Z = Y*Y ! 276: fp2 = float64_to_floatx80(LIT64(0x3FB344447F876989)); // B6 ! 277: fp3 = float64_to_floatx80(LIT64(0xBFB744EE7FAF45DB)); // B5 ! 278: fp2 = floatx80_mul(fp2, fp1); // Z*B6 ! 279: fp3 = floatx80_mul(fp3, fp1); // Z*B5 ! 280: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FBC71C646940220))); // B4+Z*B6 ! 281: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBFC24924921872F9))); // B3+Z*B5 ! 282: fp2 = floatx80_mul(fp2, fp1); // Z*(B4+Z*B6) ! 283: fp1 = floatx80_mul(fp1, fp3); // Z*(B3+Z*B5) ! 284: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FC9999999998FA9))); // B2+Z*(B4+Z*B6) ! 285: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0xBFD5555555555555))); // B1+Z*(B3+Z*B5) ! 286: fp2 = floatx80_mul(fp2, fp0); // Y*(B2+Z*(B4+Z*B6)) ! 287: fp0 = floatx80_mul(fp0, xsave); // X*Y ! 288: fp1 = floatx80_add(fp1, fp2); // [B1+Z*(B3+Z*B5)]+[Y*(B2+Z*(B4+Z*B6))] ! 289: fp0 = floatx80_mul(fp0, fp1); // X*Y*([B1+Z*(B3+Z*B5)]+[Y*(B2+Z*(B4+Z*B6))]) ! 290: ! 291: float_rounding_mode = user_rnd_mode; ! 292: floatx80_rounding_precision = user_rnd_prec; ! 293: ! 294: a = floatx80_add(fp0, xsave); ! 295: ! 296: float_raise(float_flag_inexact); ! 297: ! 298: return a; ! 299: } ! 300: } ! 301: } else { ! 302: aSig &= LIT64(0xF800000000000000); ! 303: aSig |= LIT64(0x0400000000000000); ! 304: xsave = packFloatx80(aSign, aExp, aSig); // F ! 305: fp0 = a; ! 306: fp1 = a; // X ! 307: fp2 = packFloatx80(0, one_exp, one_sig); // 1 ! 308: fp1 = floatx80_mul(fp1, xsave); // X*F ! 309: fp0 = floatx80_sub(fp0, xsave); // X-F ! 310: fp1 = floatx80_add(fp1, fp2); // 1 + X*F ! 311: fp0 = floatx80_div(fp0, fp1); // U = (X-F)/(1+X*F) ! 312: ! 313: tbl_index = compact; ! 314: ! 315: tbl_index &= 0x7FFF0000; ! 316: tbl_index -= 0x3FFB0000; ! 317: tbl_index >>= 1; ! 318: tbl_index += compact&0x00007800; ! 319: tbl_index >>= 11; ! 320: ! 321: fp3 = atan_tbl[tbl_index]; ! 322: ! 323: fp3.high |= aSign ? 0x8000 : 0; // ATAN(F) ! 324: ! 325: fp1 = floatx80_mul(fp0, fp0); // V = U*U ! 326: fp2 = float64_to_floatx80(LIT64(0xBFF6687E314987D8)); // A3 ! 327: fp2 = floatx80_add(fp2, fp1); // A3+V ! 328: fp2 = floatx80_mul(fp2, fp1); // V*(A3+V) ! 329: fp1 = floatx80_mul(fp1, fp0); // U*V ! 330: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x4002AC6934A26DB3))); // A2+V*(A3+V) ! 331: fp1 = floatx80_mul(fp1, float64_to_floatx80(LIT64(0xBFC2476F4E1DA28E))); // A1+U*V ! 332: fp1 = floatx80_mul(fp1, fp2); // A1*U*V*(A2+V*(A3+V)) ! 333: fp0 = floatx80_add(fp0, fp1); // ATAN(U) ! 334: ! 335: float_rounding_mode = user_rnd_mode; ! 336: floatx80_rounding_precision = user_rnd_prec; ! 337: ! 338: a = floatx80_add(fp0, fp3); // ATAN(X) ! 339: ! 340: float_raise(float_flag_inexact); ! 341: ! 342: return a; ! 343: } ! 344: } ! 345: ! 346: /*---------------------------------------------------------------------------- ! 347: | Hyperbolic arc tangent ! 348: *----------------------------------------------------------------------------*/ ! 349: ! 350: floatx80 floatx80_atanh(floatx80 a) ! 351: { ! 352: flag aSign; ! 353: int32 aExp; ! 354: bits64 aSig; ! 355: ! 356: int8 user_rnd_mode, user_rnd_prec; ! 357: ! 358: int32 compact; ! 359: floatx80 fp0, fp1, fp2, one; ! 360: ! 361: aSig = extractFloatx80Frac(a); ! 362: aExp = extractFloatx80Exp(a); ! 363: aSign = extractFloatx80Sign(a); ! 364: ! 365: if (aExp == 0x7FFF && (bits64) (aSig<<1)) { ! 366: return propagateFloatx80NaNOneArg(a); ! 367: } ! 368: ! 369: if (aExp == 0 && aSig == 0) { ! 370: return packFloatx80(aSign, 0, 0); ! 371: } ! 372: ! 373: compact = floatx80_make_compact(aExp, aSig); ! 374: ! 375: if (compact >= 0x3FFF8000) { // |X| >= 1 ! 376: if (aExp == one_exp && aSig == one_sig) { // |X| == 1 ! 377: float_raise(float_flag_divbyzero); ! 378: return packFloatx80(aSign, 0x7FFF, floatx80_default_infinity_low); ! 379: } else { // |X| > 1 ! 380: float_raise(float_flag_invalid); ! 381: a.low = floatx80_default_nan_low; ! 382: a.high = floatx80_default_nan_high; ! 383: return a; ! 384: } ! 385: } // |X| < 1 ! 386: ! 387: user_rnd_mode = float_rounding_mode; ! 388: user_rnd_prec = floatx80_rounding_precision; ! 389: float_rounding_mode = float_round_nearest_even; ! 390: floatx80_rounding_precision = 80; ! 391: ! 392: one = packFloatx80(0, one_exp, one_sig); ! 393: fp2 = packFloatx80(aSign, 0x3FFE, one_sig); // SIGN(X) * (1/2) ! 394: fp0 = packFloatx80(0, aExp, aSig); // Y = |X| ! 395: fp1 = packFloatx80(1, aExp, aSig); // -Y ! 396: fp0 = floatx80_add(fp0, fp0); // 2Y ! 397: fp1 = floatx80_add(fp1, one); // 1-Y ! 398: fp0 = floatx80_div(fp0, fp1); // Z = 2Y/(1-Y) ! 399: fp0 = floatx80_lognp1(fp0); // LOG1P(Z) ! 400: ! 401: float_rounding_mode = user_rnd_mode; ! 402: floatx80_rounding_precision = user_rnd_prec; ! 403: ! 404: a = floatx80_mul(fp0, fp2); // ATANH(X) = SIGN(X) * (1/2) * LOG1P(Z) ! 405: ! 406: float_raise(float_flag_inexact); ! 407: ! 408: return a; ! 409: } ! 410: ! 411: /*---------------------------------------------------------------------------- ! 412: | Cosine ! 413: *----------------------------------------------------------------------------*/ ! 414: ! 415: floatx80 floatx80_cos(floatx80 a) ! 416: { ! 417: flag aSign, xSign; ! 418: int32 aExp, xExp; ! 419: bits64 aSig, xSig; ! 420: ! 421: int8 user_rnd_mode, user_rnd_prec; ! 422: ! 423: int32 compact, l, n, j; ! 424: floatx80 fp0, fp1, fp2, fp3, fp4, fp5, x, invtwopi, twopi1, twopi2; ! 425: float32 posneg1, twoto63; ! 426: flag adjn, endflag; ! 427: ! 428: aSig = extractFloatx80Frac(a); ! 429: aExp = extractFloatx80Exp(a); ! 430: aSign = extractFloatx80Sign(a); ! 431: ! 432: if (aExp == 0x7FFF) { ! 433: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 434: float_raise(float_flag_invalid); ! 435: a.low = floatx80_default_nan_low; ! 436: a.high = floatx80_default_nan_high; ! 437: return a; ! 438: } ! 439: ! 440: if (aExp == 0 && aSig == 0) { ! 441: return packFloatx80(0, one_exp, one_sig); ! 442: } ! 443: ! 444: adjn = 1; ! 445: ! 446: user_rnd_mode = float_rounding_mode; ! 447: user_rnd_prec = floatx80_rounding_precision; ! 448: float_rounding_mode = float_round_nearest_even; ! 449: floatx80_rounding_precision = 80; ! 450: ! 451: compact = floatx80_make_compact(aExp, aSig); ! 452: ! 453: fp0 = a; ! 454: ! 455: if (compact < 0x3FD78000 || compact > 0x4004BC7E) { // 2^(-40) > |X| > 15 PI ! 456: if (compact > 0x3FFF8000) { // |X| >= 15 PI ! 457: // REDUCEX ! 458: fp1 = packFloatx80(0, 0, 0); ! 459: if (compact == 0x7FFEFFFF) { ! 460: twopi1 = packFloatx80(aSign ^ 1, 0x7FFE, LIT64(0xC90FDAA200000000)); ! 461: twopi2 = packFloatx80(aSign ^ 1, 0x7FDC, LIT64(0x85A308D300000000)); ! 462: fp0 = floatx80_add(fp0, twopi1); ! 463: fp1 = fp0; ! 464: fp0 = floatx80_add(fp0, twopi2); ! 465: fp1 = floatx80_sub(fp1, fp0); ! 466: fp1 = floatx80_add(fp1, twopi2); ! 467: } ! 468: loop: ! 469: xSign = extractFloatx80Sign(fp0); ! 470: xExp = extractFloatx80Exp(fp0); ! 471: xExp -= 0x3FFF; ! 472: if (xExp <= 28) { ! 473: l = 0; ! 474: endflag = 1; ! 475: } else { ! 476: l = xExp - 27; ! 477: endflag = 0; ! 478: } ! 479: invtwopi = packFloatx80(0, 0x3FFE - l, LIT64(0xA2F9836E4E44152A)); // INVTWOPI ! 480: twopi1 = packFloatx80(0, 0x3FFF + l, LIT64(0xC90FDAA200000000)); ! 481: twopi2 = packFloatx80(0, 0x3FDD + l, LIT64(0x85A308D300000000)); ! 482: ! 483: twoto63 = 0x5F000000; ! 484: twoto63 |= xSign ? 0x80000000 : 0x00000000; // SIGN(INARG)*2^63 IN SGL ! 485: ! 486: fp2 = floatx80_mul(fp0, invtwopi); ! 487: fp2 = floatx80_add(fp2, float32_to_floatx80(twoto63)); // THE FRACTIONAL PART OF FP2 IS ROUNDED ! 488: fp2 = floatx80_sub(fp2, float32_to_floatx80(twoto63)); // FP2 is N ! 489: fp4 = floatx80_mul(twopi1, fp2); // W = N*P1 ! 490: fp5 = floatx80_mul(twopi2, fp2); // w = N*P2 ! 491: fp3 = floatx80_add(fp4, fp5); // FP3 is P ! 492: fp4 = floatx80_sub(fp4, fp3); // W-P ! 493: fp0 = floatx80_sub(fp0, fp3); // FP0 is A := R - P ! 494: fp4 = floatx80_add(fp4, fp5); // FP4 is p = (W-P)+w ! 495: fp3 = fp0; // FP3 is A ! 496: fp1 = floatx80_sub(fp1, fp4); // FP1 is a := r - p ! 497: fp0 = floatx80_add(fp0, fp1); // FP0 is R := A+a ! 498: ! 499: if (endflag > 0) { ! 500: n = floatx80_to_int32(fp2); ! 501: goto sincont; ! 502: } ! 503: fp3 = floatx80_sub(fp3, fp0); // A-R ! 504: fp1 = floatx80_add(fp1, fp3); // FP1 is r := (A-R)+a ! 505: goto loop; ! 506: } else { ! 507: // SINSM ! 508: fp0 = float32_to_floatx80(0x3F800000); // 1 ! 509: ! 510: float_rounding_mode = user_rnd_mode; ! 511: floatx80_rounding_precision = user_rnd_prec; ! 512: ! 513: if (adjn) { ! 514: // COSTINY ! 515: a = floatx80_sub(fp0, float32_to_floatx80(0x00800000)); ! 516: } else { ! 517: // SINTINY ! 518: a = floatx80_move(a); ! 519: } ! 520: float_raise(float_flag_inexact); ! 521: ! 522: return a; ! 523: } ! 524: } else { ! 525: fp1 = floatx80_mul(fp0, float64_to_floatx80(LIT64(0x3FE45F306DC9C883))); // X*2/PI ! 526: ! 527: n = floatx80_to_int32(fp1); ! 528: j = 32 + n; ! 529: ! 530: fp0 = floatx80_sub(fp0, pi_tbl[j]); // X-Y1 ! 531: fp0 = floatx80_sub(fp0, float32_to_floatx80(pi_tbl2[j])); // FP0 IS R = (X-Y1)-Y2 ! 532: ! 533: sincont: ! 534: if ((n + adjn) & 1) { ! 535: // COSPOLY ! 536: fp0 = floatx80_mul(fp0, fp0); // FP0 IS S ! 537: fp1 = floatx80_mul(fp0, fp0); // FP1 IS T ! 538: fp2 = float64_to_floatx80(LIT64(0x3D2AC4D0D6011EE3)); // B8 ! 539: fp3 = float64_to_floatx80(LIT64(0xBDA9396F9F45AC19)); // B7 ! 540: ! 541: xSign = extractFloatx80Sign(fp0); // X IS S ! 542: xExp = extractFloatx80Exp(fp0); ! 543: xSig = extractFloatx80Frac(fp0); ! 544: ! 545: if (((n + adjn) >> 1) & 1) { ! 546: xSign ^= 1; ! 547: posneg1 = 0xBF800000; // -1 ! 548: } else { ! 549: xSign ^= 0; ! 550: posneg1 = 0x3F800000; // 1 ! 551: } // X IS NOW R'= SGN*R ! 552: ! 553: fp2 = floatx80_mul(fp2, fp1); // TB8 ! 554: fp3 = floatx80_mul(fp3, fp1); // TB7 ! 555: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3E21EED90612C972))); // B6+TB8 ! 556: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBE927E4FB79D9FCF))); // B5+TB7 ! 557: fp2 = floatx80_mul(fp2, fp1); // T(B6+TB8) ! 558: fp3 = floatx80_mul(fp3, fp1); // T(B5+TB7) ! 559: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3EFA01A01A01D423))); // B4+T(B6+TB8) ! 560: fp4 = packFloatx80(1, 0x3FF5, LIT64(0xB60B60B60B61D438)); ! 561: fp3 = floatx80_add(fp3, fp4); // B3+T(B5+TB7) ! 562: fp2 = floatx80_mul(fp2, fp1); // T(B4+T(B6+TB8)) ! 563: fp1 = floatx80_mul(fp1, fp3); // T(B3+T(B5+TB7)) ! 564: fp4 = packFloatx80(0, 0x3FFA, LIT64(0xAAAAAAAAAAAAAB5E)); ! 565: fp2 = floatx80_add(fp2, fp4); // B2+T(B4+T(B6+TB8)) ! 566: fp1 = floatx80_add(fp1, float32_to_floatx80(0xBF000000)); // B1+T(B3+T(B5+TB7)) ! 567: fp0 = floatx80_mul(fp0, fp2); // S(B2+T(B4+T(B6+TB8))) ! 568: fp0 = floatx80_add(fp0, fp1); // [B1+T(B3+T(B5+TB7))]+[S(B2+T(B4+T(B6+TB8)))] ! 569: ! 570: x = packFloatx80(xSign, xExp, xSig); ! 571: fp0 = floatx80_mul(fp0, x); ! 572: ! 573: float_rounding_mode = user_rnd_mode; ! 574: floatx80_rounding_precision = user_rnd_prec; ! 575: ! 576: a = floatx80_add(fp0, float32_to_floatx80(posneg1)); ! 577: ! 578: float_raise(float_flag_inexact); ! 579: ! 580: return a; ! 581: } else { ! 582: // SINPOLY ! 583: xSign = extractFloatx80Sign(fp0); // X IS R ! 584: xExp = extractFloatx80Exp(fp0); ! 585: xSig = extractFloatx80Frac(fp0); ! 586: ! 587: xSign ^= ((n + adjn) >> 1) & 1; // X IS NOW R'= SGN*R ! 588: ! 589: fp0 = floatx80_mul(fp0, fp0); // FP0 IS S ! 590: fp1 = floatx80_mul(fp0, fp0); // FP1 IS T ! 591: fp3 = float64_to_floatx80(LIT64(0xBD6AAA77CCC994F5)); // A7 ! 592: fp2 = float64_to_floatx80(LIT64(0x3DE612097AAE8DA1)); // A6 ! 593: fp3 = floatx80_mul(fp3, fp1); // T*A7 ! 594: fp2 = floatx80_mul(fp2, fp1); // T*A6 ! 595: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBE5AE6452A118AE4))); // A5+T*A7 ! 596: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3EC71DE3A5341531))); // A4+T*A6 ! 597: fp3 = floatx80_mul(fp3, fp1); // T(A5+TA7) ! 598: fp2 = floatx80_mul(fp2, fp1); // T(A4+TA6) ! 599: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBF2A01A01A018B59))); // A3+T(A5+TA7) ! 600: fp4 = packFloatx80(0, 0x3FF8, LIT64(0x88888888888859AF)); ! 601: fp2 = floatx80_add(fp2, fp4); // A2+T(A4+TA6) ! 602: fp1 = floatx80_mul(fp1, fp3); // T(A3+T(A5+TA7)) ! 603: fp2 = floatx80_mul(fp2, fp0); // S(A2+T(A4+TA6)) ! 604: fp4 = packFloatx80(1, 0x3FFC, LIT64(0xAAAAAAAAAAAAAA99)); ! 605: fp1 = floatx80_add(fp1, fp4); // A1+T(A3+T(A5+TA7)) ! 606: fp1 = floatx80_add(fp1, fp2); // [A1+T(A3+T(A5+TA7))]+[S(A2+T(A4+TA6))] ! 607: ! 608: x = packFloatx80(xSign, xExp, xSig); ! 609: fp0 = floatx80_mul(fp0, x); // R'*S ! 610: fp0 = floatx80_mul(fp0, fp1); // SIN(R')-R' ! 611: ! 612: float_rounding_mode = user_rnd_mode; ! 613: floatx80_rounding_precision = user_rnd_prec; ! 614: ! 615: a = floatx80_add(fp0, x); ! 616: ! 617: float_raise(float_flag_inexact); ! 618: ! 619: return a; ! 620: } ! 621: } ! 622: } ! 623: ! 624: /*---------------------------------------------------------------------------- ! 625: | Hyperbolic cosine ! 626: *----------------------------------------------------------------------------*/ ! 627: ! 628: floatx80 floatx80_cosh(floatx80 a) ! 629: { ! 630: flag aSign; ! 631: int32 aExp; ! 632: bits64 aSig; ! 633: ! 634: int8 user_rnd_mode, user_rnd_prec; ! 635: ! 636: int32 compact; ! 637: floatx80 fp0, fp1; ! 638: ! 639: aSig = extractFloatx80Frac(a); ! 640: aExp = extractFloatx80Exp(a); ! 641: aSign = extractFloatx80Sign(a); ! 642: ! 643: if (aExp == 0x7FFF) { ! 644: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 645: return packFloatx80(0, 0x7FFF, floatx80_default_infinity_low); ! 646: } ! 647: ! 648: if (aExp == 0 && aSig == 0) { ! 649: return packFloatx80(0, one_exp, one_sig); ! 650: } ! 651: ! 652: user_rnd_mode = float_rounding_mode; ! 653: user_rnd_prec = floatx80_rounding_precision; ! 654: float_rounding_mode = float_round_nearest_even; ! 655: floatx80_rounding_precision = 80; ! 656: ! 657: compact = floatx80_make_compact(aExp, aSig); ! 658: ! 659: if (compact > 0x400CB167) { ! 660: if (compact > 0x400CB2B3) { ! 661: float_rounding_mode = user_rnd_mode; ! 662: floatx80_rounding_precision = user_rnd_prec; ! 663: return roundAndPackFloatx80(floatx80_rounding_precision, 0, 0x8000, one_sig, 0); ! 664: } else { ! 665: fp0 = packFloatx80(0, aExp, aSig); ! 666: fp0 = floatx80_sub(fp0, float64_to_floatx80(LIT64(0x40C62D38D3D64634))); ! 667: fp0 = floatx80_sub(fp0, float64_to_floatx80(LIT64(0x3D6F90AEB1E75CC7))); ! 668: fp0 = floatx80_etox(fp0); ! 669: fp1 = packFloatx80(0, 0x7FFB, one_sig); ! 670: ! 671: float_rounding_mode = user_rnd_mode; ! 672: floatx80_rounding_precision = user_rnd_prec; ! 673: ! 674: a = floatx80_mul(fp0, fp1); ! 675: ! 676: float_raise(float_flag_inexact); ! 677: ! 678: return a; ! 679: } ! 680: } ! 681: ! 682: fp0 = packFloatx80(0, aExp, aSig); // |X| ! 683: fp0 = floatx80_etox(fp0); // EXP(|X|) ! 684: fp0 = floatx80_mul(fp0, float32_to_floatx80(0x3F000000)); // (1/2)*EXP(|X|) ! 685: fp1 = float32_to_floatx80(0x3E800000); // 1/4 ! 686: fp1 = floatx80_div(fp1, fp0); // 1/(2*EXP(|X|)) ! 687: ! 688: float_rounding_mode = user_rnd_mode; ! 689: floatx80_rounding_precision = user_rnd_prec; ! 690: ! 691: a = floatx80_add(fp0, fp1); ! 692: ! 693: float_raise(float_flag_inexact); ! 694: ! 695: return a; ! 696: } ! 697: ! 698: /*---------------------------------------------------------------------------- ! 699: | e to x ! 700: *----------------------------------------------------------------------------*/ ! 701: ! 702: floatx80 floatx80_etox(floatx80 a) ! 703: { ! 704: flag aSign; ! 705: int32 aExp; ! 706: bits64 aSig; ! 707: ! 708: int8 user_rnd_mode, user_rnd_prec; ! 709: ! 710: int32 compact, n, j, k, m, m1; ! 711: floatx80 fp0, fp1, fp2, fp3, l2, scale, adjscale; ! 712: flag adjflag; ! 713: ! 714: aSig = extractFloatx80Frac(a); ! 715: aExp = extractFloatx80Exp(a); ! 716: aSign = extractFloatx80Sign(a); ! 717: ! 718: if (aExp == 0x7FFF) { ! 719: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 720: if (aSign) return packFloatx80(0, 0, 0); ! 721: return packFloatx80(0, 0x7FFF, floatx80_default_infinity_low); ! 722: } ! 723: ! 724: if (aExp == 0 && aSig == 0) { ! 725: return packFloatx80(0, one_exp, one_sig); ! 726: } ! 727: ! 728: user_rnd_mode = float_rounding_mode; ! 729: user_rnd_prec = floatx80_rounding_precision; ! 730: float_rounding_mode = float_round_nearest_even; ! 731: floatx80_rounding_precision = 80; ! 732: ! 733: adjflag = 0; ! 734: ! 735: if (aExp >= 0x3FBE) { // |X| >= 2^(-65) ! 736: compact = floatx80_make_compact(aExp, aSig); ! 737: ! 738: if (compact < 0x400CB167) { // |X| < 16380 log2 ! 739: fp0 = a; ! 740: fp1 = a; ! 741: fp0 = floatx80_mul(fp0, float32_to_floatx80(0x42B8AA3B)); // 64/log2 * X ! 742: adjflag = 0; ! 743: n = floatx80_to_int32(fp0); // int(64/log2*X) ! 744: fp0 = int32_to_floatx80(n); ! 745: ! 746: j = n & 0x3F; // J = N mod 64 ! 747: m = n / 64; // NOTE: this is really arithmetic right shift by 6 ! 748: if (n < 0 && j) { // arithmetic right shift is division and round towards minus infinity ! 749: m--; ! 750: } ! 751: m += 0x3FFF; // biased exponent of 2^(M) ! 752: ! 753: expcont1: ! 754: fp2 = fp0; // N ! 755: fp0 = floatx80_mul(fp0, float32_to_floatx80(0xBC317218)); // N * L1, L1 = lead(-log2/64) ! 756: l2 = packFloatx80(0, 0x3FDC, LIT64(0x82E308654361C4C6)); ! 757: fp2 = floatx80_mul(fp2, l2); // N * L2, L1+L2 = -log2/64 ! 758: fp0 = floatx80_add(fp0, fp1); // X + N*L1 ! 759: fp0 = floatx80_add(fp0, fp2); // R ! 760: ! 761: fp1 = floatx80_mul(fp0, fp0); // S = R*R ! 762: fp2 = float32_to_floatx80(0x3AB60B70); // A5 ! 763: fp2 = floatx80_mul(fp2, fp1); // fp2 is S*A5 ! 764: fp3 = floatx80_mul(float32_to_floatx80(0x3C088895), fp1); // fp3 is S*A4 ! 765: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FA5555555554431))); // fp2 is A3+S*A5 ! 766: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0x3FC5555555554018))); // fp3 is A2+S*A4 ! 767: fp2 = floatx80_mul(fp2, fp1); // fp2 is S*(A3+S*A5) ! 768: fp3 = floatx80_mul(fp3, fp1); // fp3 is S*(A2+S*A4) ! 769: fp2 = floatx80_add(fp2, float32_to_floatx80(0x3F000000)); // fp2 is A1+S*(A3+S*A5) ! 770: fp3 = floatx80_mul(fp3, fp0); // fp3 IS R*S*(A2+S*A4) ! 771: fp2 = floatx80_mul(fp2, fp1); // fp2 IS S*(A1+S*(A3+S*A5)) ! 772: fp0 = floatx80_add(fp0, fp3); // fp0 IS R+R*S*(A2+S*A4) ! 773: fp0 = floatx80_add(fp0, fp2); // fp0 IS EXP(R) - 1 ! 774: ! 775: fp1 = exp_tbl[j]; ! 776: fp0 = floatx80_mul(fp0, fp1); // 2^(J/64)*(Exp(R)-1) ! 777: fp0 = floatx80_add(fp0, float32_to_floatx80(exp_tbl2[j])); // accurate 2^(J/64) ! 778: fp0 = floatx80_add(fp0, fp1); // 2^(J/64) + 2^(J/64)*(Exp(R)-1) ! 779: ! 780: scale = packFloatx80(0, m, one_sig); ! 781: if (adjflag) { ! 782: adjscale = packFloatx80(0, m1, one_sig); ! 783: fp0 = floatx80_mul(fp0, adjscale); ! 784: } ! 785: ! 786: float_rounding_mode = user_rnd_mode; ! 787: floatx80_rounding_precision = user_rnd_prec; ! 788: ! 789: a = floatx80_mul(fp0, scale); ! 790: ! 791: float_raise(float_flag_inexact); ! 792: ! 793: return a; ! 794: } else { // |X| >= 16380 log2 ! 795: if (compact > 0x400CB27C) { // |X| >= 16480 log2 ! 796: float_rounding_mode = user_rnd_mode; ! 797: floatx80_rounding_precision = user_rnd_prec; ! 798: if (aSign) { ! 799: a = roundAndPackFloatx80(floatx80_rounding_precision, 0, -0x1000, aSig, 0); ! 800: } else { ! 801: a = roundAndPackFloatx80(floatx80_rounding_precision, 0, 0x8000, aSig, 0); ! 802: } ! 803: float_raise(float_flag_inexact); ! 804: ! 805: return a; ! 806: } else { ! 807: fp0 = a; ! 808: fp1 = a; ! 809: fp0 = floatx80_mul(fp0, float32_to_floatx80(0x42B8AA3B)); // 64/log2 * X ! 810: adjflag = 1; ! 811: n = floatx80_to_int32(fp0); // int(64/log2*X) ! 812: fp0 = int32_to_floatx80(n); ! 813: ! 814: j = n & 0x3F; // J = N mod 64 ! 815: k = n / 64; // NOTE: this is really arithmetic right shift by 6 ! 816: if (n < 0 && j) { // arithmetic right shift is division and round towards minus infinity ! 817: k--; ! 818: } ! 819: m1 = k / 2; // NOTE: this is really arithmetic right shift by 1 ! 820: if (k < 0 && (k & 1)) { // arithmetic right shift is division and round towards minus infinity ! 821: m1--; ! 822: } ! 823: m = k - m1; ! 824: m1 += 0x3FFF; // biased exponent of 2^(M1) ! 825: m += 0x3FFF; // biased exponent of 2^(M) ! 826: ! 827: goto expcont1; ! 828: } ! 829: } ! 830: } else { // |X| < 2^(-65) ! 831: float_rounding_mode = user_rnd_mode; ! 832: floatx80_rounding_precision = user_rnd_prec; ! 833: ! 834: a = floatx80_add(a, float32_to_floatx80(0x3F800000)); // 1 + X ! 835: ! 836: float_raise(float_flag_inexact); ! 837: ! 838: return a; ! 839: } ! 840: } ! 841: ! 842: /*---------------------------------------------------------------------------- ! 843: | e to x minus 1 ! 844: *----------------------------------------------------------------------------*/ ! 845: ! 846: floatx80 floatx80_etoxm1(floatx80 a) ! 847: { ! 848: flag aSign; ! 849: int32 aExp; ! 850: bits64 aSig; ! 851: ! 852: int8 user_rnd_mode, user_rnd_prec; ! 853: ! 854: int32 compact, n, j, m, m1; ! 855: floatx80 fp0, fp1, fp2, fp3, l2, sc, onebysc; ! 856: ! 857: aSig = extractFloatx80Frac(a); ! 858: aExp = extractFloatx80Exp(a); ! 859: aSign = extractFloatx80Sign(a); ! 860: ! 861: if (aExp == 0x7FFF) { ! 862: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 863: if (aSign) return packFloatx80(aSign, one_exp, one_sig); ! 864: return packFloatx80(0, 0x7FFF, floatx80_default_infinity_low); ! 865: } ! 866: ! 867: if (aExp == 0 && aSig == 0) { ! 868: return packFloatx80(aSign, 0, 0); ! 869: } ! 870: ! 871: user_rnd_mode = float_rounding_mode; ! 872: user_rnd_prec = floatx80_rounding_precision; ! 873: float_rounding_mode = float_round_nearest_even; ! 874: floatx80_rounding_precision = 80; ! 875: ! 876: if (aExp >= 0x3FFD) { // |X| >= 1/4 ! 877: compact = floatx80_make_compact(aExp, aSig); ! 878: ! 879: if (compact <= 0x4004C215) { // |X| <= 70 log2 ! 880: fp0 = a; ! 881: fp1 = a; ! 882: fp0 = floatx80_mul(fp0, float32_to_floatx80(0x42B8AA3B)); // 64/log2 * X ! 883: n = floatx80_to_int32(fp0); // int(64/log2*X) ! 884: fp0 = int32_to_floatx80(n); ! 885: ! 886: j = n & 0x3F; // J = N mod 64 ! 887: m = n / 64; // NOTE: this is really arithmetic right shift by 6 ! 888: if (n < 0 && j) { // arithmetic right shift is division and round towards minus infinity ! 889: m--; ! 890: } ! 891: m1 = -m; ! 892: //m += 0x3FFF; // biased exponent of 2^(M) ! 893: //m1 += 0x3FFF; // biased exponent of -2^(-M) ! 894: ! 895: fp2 = fp0; // N ! 896: fp0 = floatx80_mul(fp0, float32_to_floatx80(0xBC317218)); // N * L1, L1 = lead(-log2/64) ! 897: l2 = packFloatx80(0, 0x3FDC, LIT64(0x82E308654361C4C6)); ! 898: fp2 = floatx80_mul(fp2, l2); // N * L2, L1+L2 = -log2/64 ! 899: fp0 = floatx80_add(fp0, fp1); // X + N*L1 ! 900: fp0 = floatx80_add(fp0, fp2); // R ! 901: ! 902: fp1 = floatx80_mul(fp0, fp0); // S = R*R ! 903: fp2 = float32_to_floatx80(0x3950097B); // A6 ! 904: fp2 = floatx80_mul(fp2, fp1); // fp2 is S*A6 ! 905: fp3 = floatx80_mul(float32_to_floatx80(0x3AB60B6A), fp1); // fp3 is S*A5 ! 906: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3F81111111174385))); // fp2 IS A4+S*A6 ! 907: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0x3FA5555555554F5A))); // fp3 is A3+S*A5 ! 908: fp2 = floatx80_mul(fp2, fp1); // fp2 IS S*(A4+S*A6) ! 909: fp3 = floatx80_mul(fp3, fp1); // fp3 IS S*(A3+S*A5) ! 910: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FC5555555555555))); // fp2 IS A2+S*(A4+S*A6) ! 911: fp3 = floatx80_add(fp3, float32_to_floatx80(0x3F000000)); // fp3 IS A1+S*(A3+S*A5) ! 912: fp2 = floatx80_mul(fp2, fp1); // fp2 IS S*(A2+S*(A4+S*A6)) ! 913: fp1 = floatx80_mul(fp1, fp3); // fp1 IS S*(A1+S*(A3+S*A5)) ! 914: fp2 = floatx80_mul(fp2, fp0); // fp2 IS R*S*(A2+S*(A4+S*A6)) ! 915: fp0 = floatx80_add(fp0, fp1); // fp0 IS R+S*(A1+S*(A3+S*A5)) ! 916: fp0 = floatx80_add(fp0, fp2); // fp0 IS EXP(R) - 1 ! 917: ! 918: fp0 = floatx80_mul(fp0, exp_tbl[j]); // 2^(J/64)*(Exp(R)-1) ! 919: ! 920: if (m >= 64) { ! 921: fp1 = float32_to_floatx80(exp_tbl2[j]); ! 922: onebysc = packFloatx80(1, m1 + 0x3FFF, one_sig); // -2^(-M) ! 923: fp1 = floatx80_add(fp1, onebysc); ! 924: fp0 = floatx80_add(fp0, fp1); ! 925: fp0 = floatx80_add(fp0, exp_tbl[j]); ! 926: } else if (m < -3) { ! 927: fp0 = floatx80_add(fp0, float32_to_floatx80(exp_tbl2[j])); ! 928: fp0 = floatx80_add(fp0, exp_tbl[j]); ! 929: onebysc = packFloatx80(1, m1 + 0x3FFF, one_sig); // -2^(-M) ! 930: fp0 = floatx80_add(fp0, onebysc); ! 931: } else { // -3 <= m <= 63 ! 932: fp1 = exp_tbl[j]; ! 933: fp0 = floatx80_add(fp0, float32_to_floatx80(exp_tbl2[j])); ! 934: onebysc = packFloatx80(1, m1 + 0x3FFF, one_sig); // -2^(-M) ! 935: fp1 = floatx80_add(fp1, onebysc); ! 936: fp0 = floatx80_add(fp0, fp1); ! 937: } ! 938: ! 939: sc = packFloatx80(0, m + 0x3FFF, one_sig); ! 940: ! 941: float_rounding_mode = user_rnd_mode; ! 942: floatx80_rounding_precision = user_rnd_prec; ! 943: ! 944: a = floatx80_mul(fp0, sc); ! 945: ! 946: float_raise(float_flag_inexact); ! 947: ! 948: return a; ! 949: } else { // |X| > 70 log2 ! 950: if (aSign) { ! 951: fp0 = float32_to_floatx80(0xBF800000); // -1 ! 952: ! 953: float_rounding_mode = user_rnd_mode; ! 954: floatx80_rounding_precision = user_rnd_prec; ! 955: ! 956: a = floatx80_add(fp0, float32_to_floatx80(0x00800000)); // -1 + 2^(-126) ! 957: ! 958: float_raise(float_flag_inexact); ! 959: ! 960: return a; ! 961: } else { ! 962: float_rounding_mode = user_rnd_mode; ! 963: floatx80_rounding_precision = user_rnd_prec; ! 964: ! 965: return floatx80_etox(a); ! 966: } ! 967: } ! 968: } else { // |X| < 1/4 ! 969: if (aExp >= 0x3FBE) { ! 970: fp0 = a; ! 971: fp0 = floatx80_mul(fp0, fp0); // S = X*X ! 972: fp1 = float32_to_floatx80(0x2F30CAA8); // B12 ! 973: fp1 = floatx80_mul(fp1, fp0); // S * B12 ! 974: fp2 = float32_to_floatx80(0x310F8290); // B11 ! 975: fp1 = floatx80_add(fp1, float32_to_floatx80(0x32D73220)); // B10 ! 976: fp2 = floatx80_mul(fp2, fp0); ! 977: fp1 = floatx80_mul(fp1, fp0); ! 978: fp2 = floatx80_add(fp2, float32_to_floatx80(0x3493F281)); // B9 ! 979: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0x3EC71DE3A5774682))); // B8 ! 980: fp2 = floatx80_mul(fp2, fp0); ! 981: fp1 = floatx80_mul(fp1, fp0); ! 982: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3EFA01A019D7CB68))); // B7 ! 983: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0x3F2A01A01A019DF3))); // B6 ! 984: fp2 = floatx80_mul(fp2, fp0); ! 985: fp1 = floatx80_mul(fp1, fp0); ! 986: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3F56C16C16C170E2))); // B5 ! 987: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0x3F81111111111111))); // B4 ! 988: fp2 = floatx80_mul(fp2, fp0); ! 989: fp1 = floatx80_mul(fp1, fp0); ! 990: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FA5555555555555))); // B3 ! 991: fp3 = packFloatx80(0, 0x3FFC, LIT64(0xAAAAAAAAAAAAAAAB)); ! 992: fp1 = floatx80_add(fp1, fp3); // B2 ! 993: fp2 = floatx80_mul(fp2, fp0); ! 994: fp1 = floatx80_mul(fp1, fp0); ! 995: ! 996: fp2 = floatx80_mul(fp2, fp0); ! 997: fp1 = floatx80_mul(fp1, a); ! 998: ! 999: fp0 = floatx80_mul(fp0, float32_to_floatx80(0x3F000000)); // S*B1 ! 1000: fp1 = floatx80_add(fp1, fp2); // Q ! 1001: fp0 = floatx80_add(fp0, fp1); // S*B1+Q ! 1002: ! 1003: float_rounding_mode = user_rnd_mode; ! 1004: floatx80_rounding_precision = user_rnd_prec; ! 1005: ! 1006: a = floatx80_add(fp0, a); ! 1007: ! 1008: float_raise(float_flag_inexact); ! 1009: ! 1010: return a; ! 1011: } else { // |X| < 2^(-65) ! 1012: sc = packFloatx80(1, 1, one_sig); ! 1013: fp0 = a; ! 1014: ! 1015: if (aExp < 0x0033) { // |X| < 2^(-16382) ! 1016: fp0 = floatx80_mul(fp0, float64_to_floatx80(LIT64(0x48B0000000000000))); ! 1017: fp0 = floatx80_add(fp0, sc); ! 1018: ! 1019: float_rounding_mode = user_rnd_mode; ! 1020: floatx80_rounding_precision = user_rnd_prec; ! 1021: ! 1022: a = floatx80_mul(fp0, float64_to_floatx80(LIT64(0x3730000000000000))); ! 1023: } else { ! 1024: float_rounding_mode = user_rnd_mode; ! 1025: floatx80_rounding_precision = user_rnd_prec; ! 1026: ! 1027: a = floatx80_add(fp0, sc); ! 1028: } ! 1029: ! 1030: float_raise(float_flag_inexact); ! 1031: ! 1032: return a; ! 1033: } ! 1034: } ! 1035: } ! 1036: ! 1037: /*---------------------------------------------------------------------------- ! 1038: | Log base 10 ! 1039: *----------------------------------------------------------------------------*/ ! 1040: ! 1041: floatx80 floatx80_log10(floatx80 a) ! 1042: { ! 1043: flag aSign; ! 1044: int32 aExp; ! 1045: bits64 aSig; ! 1046: ! 1047: int8 user_rnd_mode, user_rnd_prec; ! 1048: ! 1049: floatx80 fp0, fp1; ! 1050: ! 1051: aSig = extractFloatx80Frac(a); ! 1052: aExp = extractFloatx80Exp(a); ! 1053: aSign = extractFloatx80Sign(a); ! 1054: ! 1055: if (aExp == 0x7FFF) { ! 1056: if ((bits64) (aSig<<1)) propagateFloatx80NaNOneArg(a); ! 1057: if (aSign == 0) ! 1058: return packFloatx80(0, 0x7FFF, floatx80_default_infinity_low); ! 1059: } ! 1060: ! 1061: if (aExp == 0 && aSig == 0) { ! 1062: float_raise(float_flag_divbyzero); ! 1063: return packFloatx80(1, 0x7FFF, floatx80_default_infinity_low); ! 1064: } ! 1065: ! 1066: if (aSign) { ! 1067: float_raise(float_flag_invalid); ! 1068: a.low = floatx80_default_nan_low; ! 1069: a.high = floatx80_default_nan_high; ! 1070: return a; ! 1071: } ! 1072: ! 1073: user_rnd_mode = float_rounding_mode; ! 1074: user_rnd_prec = floatx80_rounding_precision; ! 1075: float_rounding_mode = float_round_nearest_even; ! 1076: floatx80_rounding_precision = 80; ! 1077: ! 1078: fp0 = floatx80_logn(a); ! 1079: fp1 = packFloatx80(0, 0x3FFD, LIT64(0xDE5BD8A937287195)); // INV_L10 ! 1080: ! 1081: float_rounding_mode = user_rnd_mode; ! 1082: floatx80_rounding_precision = user_rnd_prec; ! 1083: ! 1084: a = floatx80_mul(fp0, fp1); // LOGN(X)*INV_L10 ! 1085: ! 1086: float_raise(float_flag_inexact); ! 1087: ! 1088: return a; ! 1089: } ! 1090: ! 1091: /*---------------------------------------------------------------------------- ! 1092: | Log base 2 ! 1093: *----------------------------------------------------------------------------*/ ! 1094: ! 1095: floatx80 floatx80_log2(floatx80 a) ! 1096: { ! 1097: flag aSign; ! 1098: int32 aExp; ! 1099: bits64 aSig; ! 1100: ! 1101: int8 user_rnd_mode, user_rnd_prec; ! 1102: ! 1103: floatx80 fp0, fp1; ! 1104: ! 1105: aSig = extractFloatx80Frac(a); ! 1106: aExp = extractFloatx80Exp(a); ! 1107: aSign = extractFloatx80Sign(a); ! 1108: ! 1109: if (aExp == 0x7FFF) { ! 1110: if ((bits64) (aSig<<1)) propagateFloatx80NaNOneArg(a); ! 1111: if (aSign == 0) ! 1112: return packFloatx80(0, 0x7FFF, floatx80_default_infinity_low); ! 1113: } ! 1114: ! 1115: if (aExp == 0) { ! 1116: if (aSig == 0) { ! 1117: float_raise(float_flag_divbyzero); ! 1118: return packFloatx80(1, 0x7FFF, floatx80_default_infinity_low); ! 1119: } ! 1120: normalizeFloatx80Subnormal(aSig, &aExp, &aSig); ! 1121: } ! 1122: ! 1123: if (aSign) { ! 1124: float_raise(float_flag_invalid); ! 1125: a.low = floatx80_default_nan_low; ! 1126: a.high = floatx80_default_nan_high; ! 1127: return a; ! 1128: } ! 1129: ! 1130: user_rnd_mode = float_rounding_mode; ! 1131: user_rnd_prec = floatx80_rounding_precision; ! 1132: float_rounding_mode = float_round_nearest_even; ! 1133: floatx80_rounding_precision = 80; ! 1134: ! 1135: if (aSig == one_sig) { // X is 2^k ! 1136: float_rounding_mode = user_rnd_mode; ! 1137: floatx80_rounding_precision = user_rnd_prec; ! 1138: ! 1139: a = int32_to_floatx80(aExp-0x3FFF); ! 1140: } else { ! 1141: fp0 = floatx80_logn(a); ! 1142: fp1 = packFloatx80(0, 0x3FFF, LIT64(0xB8AA3B295C17F0BC)); // INV_L2 ! 1143: ! 1144: float_rounding_mode = user_rnd_mode; ! 1145: floatx80_rounding_precision = user_rnd_prec; ! 1146: ! 1147: a = floatx80_mul(fp0, fp1); // LOGN(X)*INV_L2 ! 1148: } ! 1149: ! 1150: float_raise(float_flag_inexact); ! 1151: ! 1152: return a; ! 1153: } ! 1154: ! 1155: /*---------------------------------------------------------------------------- ! 1156: | Log base e ! 1157: *----------------------------------------------------------------------------*/ ! 1158: ! 1159: floatx80 floatx80_logn(floatx80 a) ! 1160: { ! 1161: flag aSign; ! 1162: int32 aExp; ! 1163: bits64 aSig, fSig; ! 1164: ! 1165: int8 user_rnd_mode, user_rnd_prec; ! 1166: ! 1167: int32 compact, j, k, adjk; ! 1168: floatx80 fp0, fp1, fp2, fp3, f, logof2, klog2, saveu; ! 1169: ! 1170: aSig = extractFloatx80Frac(a); ! 1171: aExp = extractFloatx80Exp(a); ! 1172: aSign = extractFloatx80Sign(a); ! 1173: ! 1174: if (aExp == 0x7FFF) { ! 1175: if ((bits64) (aSig<<1)) propagateFloatx80NaNOneArg(a); ! 1176: if (aSign == 0) ! 1177: return packFloatx80(0, 0x7FFF, floatx80_default_infinity_low); ! 1178: } ! 1179: ! 1180: adjk = 0; ! 1181: ! 1182: if (aExp == 0) { ! 1183: if (aSig == 0) { // zero ! 1184: float_raise(float_flag_divbyzero); ! 1185: return packFloatx80(1, 0x7FFF, floatx80_default_infinity_low); ! 1186: } ! 1187: #if 1 ! 1188: if ((aSig & one_sig) == 0) { // denormal ! 1189: normalizeFloatx80Subnormal(aSig, &aExp, &aSig); ! 1190: adjk = -100; ! 1191: aExp += 100; ! 1192: a = packFloatx80(aSign, aExp, aSig); ! 1193: } ! 1194: #else ! 1195: normalizeFloatx80Subnormal(aSig, &aExp, &aSig); ! 1196: #endif ! 1197: } ! 1198: ! 1199: if (aSign) { ! 1200: float_raise(float_flag_invalid); ! 1201: a.low = floatx80_default_nan_low; ! 1202: a.high = floatx80_default_nan_high; ! 1203: return a; ! 1204: } ! 1205: ! 1206: user_rnd_mode = float_rounding_mode; ! 1207: user_rnd_prec = floatx80_rounding_precision; ! 1208: float_rounding_mode = float_round_nearest_even; ! 1209: floatx80_rounding_precision = 80; ! 1210: ! 1211: compact = floatx80_make_compact(aExp, aSig); ! 1212: ! 1213: if (compact < 0x3FFEF07D || compact > 0x3FFF8841) { // |X| < 15/16 or |X| > 17/16 ! 1214: k = aExp - 0x3FFF; ! 1215: k += adjk; ! 1216: fp1 = int32_to_floatx80(k); ! 1217: ! 1218: fSig = (aSig & LIT64(0xFE00000000000000)) | LIT64(0x0100000000000000); ! 1219: j = (fSig >> 56) & 0x7E; // DISPLACEMENT FOR 1/F ! 1220: ! 1221: f = packFloatx80(0, 0x3FFF, fSig); // F ! 1222: fp0 = packFloatx80(0, 0x3FFF, aSig); // Y ! 1223: ! 1224: fp0 = floatx80_sub(fp0, f); // Y-F ! 1225: ! 1226: // LP1CONT1 ! 1227: fp0 = floatx80_mul(fp0, log_tbl[j]); // FP0 IS U = (Y-F)/F ! 1228: logof2 = packFloatx80(0, 0x3FFE, LIT64(0xB17217F7D1CF79AC)); ! 1229: klog2 = floatx80_mul(fp1, logof2); // FP1 IS K*LOG2 ! 1230: fp2 = floatx80_mul(fp0, fp0); // FP2 IS V=U*U ! 1231: ! 1232: fp3 = fp2; ! 1233: fp1 = fp2; ! 1234: ! 1235: fp1 = floatx80_mul(fp1, float64_to_floatx80(LIT64(0x3FC2499AB5E4040B))); // V*A6 ! 1236: fp2 = floatx80_mul(fp2, float64_to_floatx80(LIT64(0xBFC555B5848CB7DB))); // V*A5 ! 1237: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0x3FC99999987D8730))); // A4+V*A6 ! 1238: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0xBFCFFFFFFF6F7E97))); // A3+V*A5 ! 1239: fp1 = floatx80_mul(fp1, fp3); // V*(A4+V*A6) ! 1240: fp2 = floatx80_mul(fp2, fp3); // V*(A3+V*A5) ! 1241: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0x3FD55555555555A4))); // A2+V*(A4+V*A6) ! 1242: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0xBFE0000000000008))); // A1+V*(A3+V*A5) ! 1243: fp1 = floatx80_mul(fp1, fp3); // V*(A2+V*(A4+V*A6)) ! 1244: fp2 = floatx80_mul(fp2, fp3); // V*(A1+V*(A3+V*A5)) ! 1245: fp1 = floatx80_mul(fp1, fp0); // U*V*(A2+V*(A4+V*A6)) ! 1246: fp0 = floatx80_add(fp0, fp2); // U+V*(A1+V*(A3+V*A5)) ! 1247: ! 1248: fp1 = floatx80_add(fp1, log_tbl[j+1]); // LOG(F)+U*V*(A2+V*(A4+V*A6)) ! 1249: fp0 = floatx80_add(fp0, fp1); // FP0 IS LOG(F) + LOG(1+U) ! 1250: ! 1251: float_rounding_mode = user_rnd_mode; ! 1252: floatx80_rounding_precision = user_rnd_prec; ! 1253: ! 1254: a = floatx80_add(fp0, klog2); ! 1255: ! 1256: float_raise(float_flag_inexact); ! 1257: ! 1258: return a; ! 1259: } else { // |X-1| >= 1/16 ! 1260: fp0 = a; ! 1261: fp1 = a; ! 1262: fp1 = floatx80_sub(fp1, float32_to_floatx80(0x3F800000)); // FP1 IS X-1 ! 1263: fp0 = floatx80_add(fp0, float32_to_floatx80(0x3F800000)); // FP0 IS X+1 ! 1264: fp1 = floatx80_add(fp1, fp1); // FP1 IS 2(X-1) ! 1265: ! 1266: // LP1CONT2 ! 1267: fp1 = floatx80_div(fp1, fp0); // U ! 1268: saveu = fp1; ! 1269: fp0 = floatx80_mul(fp1, fp1); // FP0 IS V = U*U ! 1270: fp1 = floatx80_mul(fp0, fp0); // FP1 IS W = V*V ! 1271: ! 1272: fp3 = float64_to_floatx80(LIT64(0x3F175496ADD7DAD6)); // B5 ! 1273: fp2 = float64_to_floatx80(LIT64(0x3F3C71C2FE80C7E0)); // B4 ! 1274: fp3 = floatx80_mul(fp3, fp1); // W*B5 ! 1275: fp2 = floatx80_mul(fp2, fp1); // W*B4 ! 1276: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0x3F624924928BCCFF))); // B3+W*B5 ! 1277: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3F899999999995EC))); // B2+W*B4 ! 1278: fp1 = floatx80_mul(fp1, fp3); // W*(B3+W*B5) ! 1279: fp2 = floatx80_mul(fp2, fp0); // V*(B2+W*B4) ! 1280: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0x3FB5555555555555))); // B1+W*(B3+W*B5) ! 1281: ! 1282: fp0 = floatx80_mul(fp0, saveu); // FP0 IS U*V ! 1283: fp1 = floatx80_add(fp1, fp2); // B1+W*(B3+W*B5) + V*(B2+W*B4) ! 1284: fp0 = floatx80_mul(fp0, fp1); // U*V*( [B1+W*(B3+W*B5)] + [V*(B2+W*B4)] ) ! 1285: ! 1286: float_rounding_mode = user_rnd_mode; ! 1287: floatx80_rounding_precision = user_rnd_prec; ! 1288: ! 1289: a = floatx80_add(fp0, saveu); ! 1290: ! 1291: //if (!floatx80_is_zero(a)) { ! 1292: float_raise(float_flag_inexact); ! 1293: //} ! 1294: ! 1295: return a; ! 1296: } ! 1297: } ! 1298: ! 1299: /*---------------------------------------------------------------------------- ! 1300: | Log base e of x plus 1 ! 1301: *----------------------------------------------------------------------------*/ ! 1302: ! 1303: floatx80 floatx80_lognp1(floatx80 a) ! 1304: { ! 1305: flag aSign; ! 1306: int32 aExp; ! 1307: bits64 aSig, fSig; ! 1308: ! 1309: int8 user_rnd_mode, user_rnd_prec; ! 1310: ! 1311: int32 compact, j, k; ! 1312: floatx80 fp0, fp1, fp2, fp3, f, logof2, klog2, saveu; ! 1313: ! 1314: aSig = extractFloatx80Frac(a); ! 1315: aExp = extractFloatx80Exp(a); ! 1316: aSign = extractFloatx80Sign(a); ! 1317: ! 1318: if (aExp == 0x7FFF) { ! 1319: if ((bits64) (aSig<<1)) propagateFloatx80NaNOneArg(a); ! 1320: if (aSign) { ! 1321: float_raise(float_flag_invalid); ! 1322: a.low = floatx80_default_nan_low; ! 1323: a.high = floatx80_default_nan_high; ! 1324: return a; ! 1325: } ! 1326: return packFloatx80(0, 0x7FFF, floatx80_default_infinity_low); ! 1327: } ! 1328: ! 1329: if (aExp == 0 && aSig == 0) { ! 1330: return packFloatx80(aSign, 0, 0); ! 1331: } ! 1332: ! 1333: if (aSign && aExp >= one_exp) { ! 1334: if (aExp == one_exp && aSig == one_sig) { ! 1335: float_raise(float_flag_divbyzero); ! 1336: packFloatx80(aSign, 0x7FFF, floatx80_default_infinity_low); ! 1337: } ! 1338: float_raise(float_flag_invalid); ! 1339: a.low = floatx80_default_nan_low; ! 1340: a.high = floatx80_default_nan_high; ! 1341: return a; ! 1342: } ! 1343: ! 1344: if (aExp < 0x3f99 || (aExp == 0x3f99 && aSig == one_sig)) { // <= min threshold ! 1345: float_raise(float_flag_inexact); ! 1346: return floatx80_move(a); ! 1347: } ! 1348: ! 1349: user_rnd_mode = float_rounding_mode; ! 1350: user_rnd_prec = floatx80_rounding_precision; ! 1351: float_rounding_mode = float_round_nearest_even; ! 1352: floatx80_rounding_precision = 80; ! 1353: ! 1354: compact = floatx80_make_compact(aExp, aSig); ! 1355: ! 1356: fp0 = a; // Z ! 1357: fp1 = a; ! 1358: ! 1359: fp0 = floatx80_add(fp0, float32_to_floatx80(0x3F800000)); // X = (1+Z) ! 1360: ! 1361: aExp = extractFloatx80Exp(fp0); ! 1362: aSig = extractFloatx80Frac(fp0); ! 1363: ! 1364: compact = floatx80_make_compact(aExp, aSig); ! 1365: ! 1366: if (compact < 0x3FFE8000 || compact > 0x3FFFC000) { // |X| < 1/2 or |X| > 3/2 ! 1367: k = aExp - 0x3FFF; ! 1368: fp1 = int32_to_floatx80(k); ! 1369: ! 1370: fSig = (aSig & LIT64(0xFE00000000000000)) | LIT64(0x0100000000000000); ! 1371: j = (fSig >> 56) & 0x7E; // DISPLACEMENT FOR 1/F ! 1372: ! 1373: f = packFloatx80(0, 0x3FFF, fSig); // F ! 1374: fp0 = packFloatx80(0, 0x3FFF, aSig); // Y ! 1375: ! 1376: fp0 = floatx80_sub(fp0, f); // Y-F ! 1377: ! 1378: lp1cont1: ! 1379: // LP1CONT1 ! 1380: fp0 = floatx80_mul(fp0, log_tbl[j]); // FP0 IS U = (Y-F)/F ! 1381: logof2 = packFloatx80(0, 0x3FFE, LIT64(0xB17217F7D1CF79AC)); ! 1382: klog2 = floatx80_mul(fp1, logof2); // FP1 IS K*LOG2 ! 1383: fp2 = floatx80_mul(fp0, fp0); // FP2 IS V=U*U ! 1384: ! 1385: fp3 = fp2; ! 1386: fp1 = fp2; ! 1387: ! 1388: fp1 = floatx80_mul(fp1, float64_to_floatx80(LIT64(0x3FC2499AB5E4040B))); // V*A6 ! 1389: fp2 = floatx80_mul(fp2, float64_to_floatx80(LIT64(0xBFC555B5848CB7DB))); // V*A5 ! 1390: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0x3FC99999987D8730))); // A4+V*A6 ! 1391: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0xBFCFFFFFFF6F7E97))); // A3+V*A5 ! 1392: fp1 = floatx80_mul(fp1, fp3); // V*(A4+V*A6) ! 1393: fp2 = floatx80_mul(fp2, fp3); // V*(A3+V*A5) ! 1394: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0x3FD55555555555A4))); // A2+V*(A4+V*A6) ! 1395: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0xBFE0000000000008))); // A1+V*(A3+V*A5) ! 1396: fp1 = floatx80_mul(fp1, fp3); // V*(A2+V*(A4+V*A6)) ! 1397: fp2 = floatx80_mul(fp2, fp3); // V*(A1+V*(A3+V*A5)) ! 1398: fp1 = floatx80_mul(fp1, fp0); // U*V*(A2+V*(A4+V*A6)) ! 1399: fp0 = floatx80_add(fp0, fp2); // U+V*(A1+V*(A3+V*A5)) ! 1400: ! 1401: fp1 = floatx80_add(fp1, log_tbl[j+1]); // LOG(F)+U*V*(A2+V*(A4+V*A6)) ! 1402: fp0 = floatx80_add(fp0, fp1); // FP0 IS LOG(F) + LOG(1+U) ! 1403: ! 1404: float_rounding_mode = user_rnd_mode; ! 1405: floatx80_rounding_precision = user_rnd_prec; ! 1406: ! 1407: a = floatx80_add(fp0, klog2); ! 1408: ! 1409: float_raise(float_flag_inexact); ! 1410: ! 1411: return a; ! 1412: } else if (compact < 0x3FFEF07D || compact > 0x3FFF8841) { // |X| < 1/16 or |X| > -1/16 ! 1413: // LP1CARE ! 1414: fSig = (aSig & LIT64(0xFE00000000000000)) | LIT64(0x0100000000000000); ! 1415: f = packFloatx80(0, 0x3FFF, fSig); // F ! 1416: j = (fSig >> 56) & 0x7E; // DISPLACEMENT FOR 1/F ! 1417: ! 1418: if (compact >= 0x3FFF8000) { // 1+Z >= 1 ! 1419: // KISZERO ! 1420: fp0 = floatx80_sub(float32_to_floatx80(0x3F800000), f); // 1-F ! 1421: fp0 = floatx80_add(fp0, fp1); // FP0 IS Y-F = (1-F)+Z ! 1422: fp1 = packFloatx80(0, 0, 0); // K = 0 ! 1423: } else { ! 1424: // KISNEG ! 1425: fp0 = floatx80_sub(float32_to_floatx80(0x40000000), f); // 2-F ! 1426: fp1 = floatx80_add(fp1, fp1); // 2Z ! 1427: fp0 = floatx80_add(fp0, fp1); // FP0 IS Y-F = (2-F)+2Z ! 1428: fp1 = packFloatx80(1, one_exp, one_sig); // K = -1 ! 1429: } ! 1430: goto lp1cont1; ! 1431: } else { ! 1432: // LP1ONE16 ! 1433: fp1 = floatx80_add(fp1, fp1); // FP1 IS 2Z ! 1434: fp0 = floatx80_add(fp0, float32_to_floatx80(0x3F800000)); // FP0 IS 1+X ! 1435: ! 1436: // LP1CONT2 ! 1437: fp1 = floatx80_div(fp1, fp0); // U ! 1438: saveu = fp1; ! 1439: fp0 = floatx80_mul(fp1, fp1); // FP0 IS V = U*U ! 1440: fp1 = floatx80_mul(fp0, fp0); // FP1 IS W = V*V ! 1441: ! 1442: fp3 = float64_to_floatx80(LIT64(0x3F175496ADD7DAD6)); // B5 ! 1443: fp2 = float64_to_floatx80(LIT64(0x3F3C71C2FE80C7E0)); // B4 ! 1444: fp3 = floatx80_mul(fp3, fp1); // W*B5 ! 1445: fp2 = floatx80_mul(fp2, fp1); // W*B4 ! 1446: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0x3F624924928BCCFF))); // B3+W*B5 ! 1447: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3F899999999995EC))); // B2+W*B4 ! 1448: fp1 = floatx80_mul(fp1, fp3); // W*(B3+W*B5) ! 1449: fp2 = floatx80_mul(fp2, fp0); // V*(B2+W*B4) ! 1450: fp1 = floatx80_add(fp1, float64_to_floatx80(LIT64(0x3FB5555555555555))); // B1+W*(B3+W*B5) ! 1451: ! 1452: fp0 = floatx80_mul(fp0, saveu); // FP0 IS U*V ! 1453: fp1 = floatx80_add(fp1, fp2); // B1+W*(B3+W*B5) + V*(B2+W*B4) ! 1454: fp0 = floatx80_mul(fp0, fp1); // U*V*( [B1+W*(B3+W*B5)] + [V*(B2+W*B4)] ) ! 1455: ! 1456: float_rounding_mode = user_rnd_mode; ! 1457: floatx80_rounding_precision = user_rnd_prec; ! 1458: ! 1459: a = floatx80_add(fp0, saveu); ! 1460: ! 1461: //if (!floatx80_is_zero(a)) { ! 1462: float_raise(float_flag_inexact); ! 1463: //} ! 1464: ! 1465: return a; ! 1466: } ! 1467: } ! 1468: ! 1469: /*---------------------------------------------------------------------------- ! 1470: | Sine ! 1471: *----------------------------------------------------------------------------*/ ! 1472: ! 1473: floatx80 floatx80_sin(floatx80 a) ! 1474: { ! 1475: flag aSign, xSign; ! 1476: int32 aExp, xExp; ! 1477: bits64 aSig, xSig; ! 1478: ! 1479: int8 user_rnd_mode, user_rnd_prec; ! 1480: ! 1481: int32 compact, l, n, j; ! 1482: floatx80 fp0, fp1, fp2, fp3, fp4, fp5, x, invtwopi, twopi1, twopi2; ! 1483: float32 posneg1, twoto63; ! 1484: flag adjn, endflag; ! 1485: ! 1486: aSig = extractFloatx80Frac(a); ! 1487: aExp = extractFloatx80Exp(a); ! 1488: aSign = extractFloatx80Sign(a); ! 1489: ! 1490: if (aExp == 0x7FFF) { ! 1491: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 1492: float_raise(float_flag_invalid); ! 1493: a.low = floatx80_default_nan_low; ! 1494: a.high = floatx80_default_nan_high; ! 1495: return a; ! 1496: } ! 1497: ! 1498: if (aExp == 0 && aSig == 0) { ! 1499: return packFloatx80(aSign, 0, 0); ! 1500: } ! 1501: ! 1502: adjn = 0; ! 1503: ! 1504: user_rnd_mode = float_rounding_mode; ! 1505: user_rnd_prec = floatx80_rounding_precision; ! 1506: float_rounding_mode = float_round_nearest_even; ! 1507: floatx80_rounding_precision = 80; ! 1508: ! 1509: compact = floatx80_make_compact(aExp, aSig); ! 1510: ! 1511: fp0 = a; ! 1512: ! 1513: if (compact < 0x3FD78000 || compact > 0x4004BC7E) { // 2^(-40) > |X| > 15 PI ! 1514: if (compact > 0x3FFF8000) { // |X| >= 15 PI ! 1515: // REDUCEX ! 1516: fp1 = packFloatx80(0, 0, 0); ! 1517: if (compact == 0x7FFEFFFF) { ! 1518: twopi1 = packFloatx80(aSign ^ 1, 0x7FFE, LIT64(0xC90FDAA200000000)); ! 1519: twopi2 = packFloatx80(aSign ^ 1, 0x7FDC, LIT64(0x85A308D300000000)); ! 1520: fp0 = floatx80_add(fp0, twopi1); ! 1521: fp1 = fp0; ! 1522: fp0 = floatx80_add(fp0, twopi2); ! 1523: fp1 = floatx80_sub(fp1, fp0); ! 1524: fp1 = floatx80_add(fp1, twopi2); ! 1525: } ! 1526: loop: ! 1527: xSign = extractFloatx80Sign(fp0); ! 1528: xExp = extractFloatx80Exp(fp0); ! 1529: xExp -= 0x3FFF; ! 1530: if (xExp <= 28) { ! 1531: l = 0; ! 1532: endflag = 1; ! 1533: } else { ! 1534: l = xExp - 27; ! 1535: endflag = 0; ! 1536: } ! 1537: invtwopi = packFloatx80(0, 0x3FFE - l, LIT64(0xA2F9836E4E44152A)); // INVTWOPI ! 1538: twopi1 = packFloatx80(0, 0x3FFF + l, LIT64(0xC90FDAA200000000)); ! 1539: twopi2 = packFloatx80(0, 0x3FDD + l, LIT64(0x85A308D300000000)); ! 1540: ! 1541: twoto63 = 0x5F000000; ! 1542: twoto63 |= xSign ? 0x80000000 : 0x00000000; // SIGN(INARG)*2^63 IN SGL ! 1543: ! 1544: fp2 = floatx80_mul(fp0, invtwopi); ! 1545: fp2 = floatx80_add(fp2, float32_to_floatx80(twoto63)); // THE FRACTIONAL PART OF FP2 IS ROUNDED ! 1546: fp2 = floatx80_sub(fp2, float32_to_floatx80(twoto63)); // FP2 is N ! 1547: fp4 = floatx80_mul(twopi1, fp2); // W = N*P1 ! 1548: fp5 = floatx80_mul(twopi2, fp2); // w = N*P2 ! 1549: fp3 = floatx80_add(fp4, fp5); // FP3 is P ! 1550: fp4 = floatx80_sub(fp4, fp3); // W-P ! 1551: fp0 = floatx80_sub(fp0, fp3); // FP0 is A := R - P ! 1552: fp4 = floatx80_add(fp4, fp5); // FP4 is p = (W-P)+w ! 1553: fp3 = fp0; // FP3 is A ! 1554: fp1 = floatx80_sub(fp1, fp4); // FP1 is a := r - p ! 1555: fp0 = floatx80_add(fp0, fp1); // FP0 is R := A+a ! 1556: ! 1557: if (endflag > 0) { ! 1558: n = floatx80_to_int32(fp2); ! 1559: goto sincont; ! 1560: } ! 1561: fp3 = floatx80_sub(fp3, fp0); // A-R ! 1562: fp1 = floatx80_add(fp1, fp3); // FP1 is r := (A-R)+a ! 1563: goto loop; ! 1564: } else { ! 1565: // SINSM ! 1566: fp0 = float32_to_floatx80(0x3F800000); // 1 ! 1567: ! 1568: float_rounding_mode = user_rnd_mode; ! 1569: floatx80_rounding_precision = user_rnd_prec; ! 1570: ! 1571: if (adjn) { ! 1572: // COSTINY ! 1573: a = floatx80_sub(fp0, float32_to_floatx80(0x00800000)); ! 1574: } else { ! 1575: // SINTINY ! 1576: a = floatx80_move(a); ! 1577: } ! 1578: float_raise(float_flag_inexact); ! 1579: ! 1580: return a; ! 1581: } ! 1582: } else { ! 1583: fp1 = floatx80_mul(fp0, float64_to_floatx80(LIT64(0x3FE45F306DC9C883))); // X*2/PI ! 1584: ! 1585: n = floatx80_to_int32(fp1); ! 1586: j = 32 + n; ! 1587: ! 1588: fp0 = floatx80_sub(fp0, pi_tbl[j]); // X-Y1 ! 1589: fp0 = floatx80_sub(fp0, float32_to_floatx80(pi_tbl2[j])); // FP0 IS R = (X-Y1)-Y2 ! 1590: ! 1591: sincont: ! 1592: if ((n + adjn) & 1) { ! 1593: // COSPOLY ! 1594: fp0 = floatx80_mul(fp0, fp0); // FP0 IS S ! 1595: fp1 = floatx80_mul(fp0, fp0); // FP1 IS T ! 1596: fp2 = float64_to_floatx80(LIT64(0x3D2AC4D0D6011EE3)); // B8 ! 1597: fp3 = float64_to_floatx80(LIT64(0xBDA9396F9F45AC19)); // B7 ! 1598: ! 1599: xSign = extractFloatx80Sign(fp0); // X IS S ! 1600: xExp = extractFloatx80Exp(fp0); ! 1601: xSig = extractFloatx80Frac(fp0); ! 1602: ! 1603: if (((n + adjn) >> 1) & 1) { ! 1604: xSign ^= 1; ! 1605: posneg1 = 0xBF800000; // -1 ! 1606: } else { ! 1607: xSign ^= 0; ! 1608: posneg1 = 0x3F800000; // 1 ! 1609: } // X IS NOW R'= SGN*R ! 1610: ! 1611: fp2 = floatx80_mul(fp2, fp1); // TB8 ! 1612: fp3 = floatx80_mul(fp3, fp1); // TB7 ! 1613: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3E21EED90612C972))); // B6+TB8 ! 1614: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBE927E4FB79D9FCF))); // B5+TB7 ! 1615: fp2 = floatx80_mul(fp2, fp1); // T(B6+TB8) ! 1616: fp3 = floatx80_mul(fp3, fp1); // T(B5+TB7) ! 1617: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3EFA01A01A01D423))); // B4+T(B6+TB8) ! 1618: fp4 = packFloatx80(1, 0x3FF5, LIT64(0xB60B60B60B61D438)); ! 1619: fp3 = floatx80_add(fp3, fp4); // B3+T(B5+TB7) ! 1620: fp2 = floatx80_mul(fp2, fp1); // T(B4+T(B6+TB8)) ! 1621: fp1 = floatx80_mul(fp1, fp3); // T(B3+T(B5+TB7)) ! 1622: fp4 = packFloatx80(0, 0x3FFA, LIT64(0xAAAAAAAAAAAAAB5E)); ! 1623: fp2 = floatx80_add(fp2, fp4); // B2+T(B4+T(B6+TB8)) ! 1624: fp1 = floatx80_add(fp1, float32_to_floatx80(0xBF000000)); // B1+T(B3+T(B5+TB7)) ! 1625: fp0 = floatx80_mul(fp0, fp2); // S(B2+T(B4+T(B6+TB8))) ! 1626: fp0 = floatx80_add(fp0, fp1); // [B1+T(B3+T(B5+TB7))]+[S(B2+T(B4+T(B6+TB8)))] ! 1627: ! 1628: x = packFloatx80(xSign, xExp, xSig); ! 1629: fp0 = floatx80_mul(fp0, x); ! 1630: ! 1631: float_rounding_mode = user_rnd_mode; ! 1632: floatx80_rounding_precision = user_rnd_prec; ! 1633: ! 1634: a = floatx80_add(fp0, float32_to_floatx80(posneg1)); ! 1635: ! 1636: float_raise(float_flag_inexact); ! 1637: ! 1638: return a; ! 1639: } else { ! 1640: // SINPOLY ! 1641: xSign = extractFloatx80Sign(fp0); // X IS R ! 1642: xExp = extractFloatx80Exp(fp0); ! 1643: xSig = extractFloatx80Frac(fp0); ! 1644: ! 1645: xSign ^= ((n + adjn) >> 1) & 1; // X IS NOW R'= SGN*R ! 1646: ! 1647: fp0 = floatx80_mul(fp0, fp0); // FP0 IS S ! 1648: fp1 = floatx80_mul(fp0, fp0); // FP1 IS T ! 1649: fp3 = float64_to_floatx80(LIT64(0xBD6AAA77CCC994F5)); // A7 ! 1650: fp2 = float64_to_floatx80(LIT64(0x3DE612097AAE8DA1)); // A6 ! 1651: fp3 = floatx80_mul(fp3, fp1); // T*A7 ! 1652: fp2 = floatx80_mul(fp2, fp1); // T*A6 ! 1653: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBE5AE6452A118AE4))); // A5+T*A7 ! 1654: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3EC71DE3A5341531))); // A4+T*A6 ! 1655: fp3 = floatx80_mul(fp3, fp1); // T(A5+TA7) ! 1656: fp2 = floatx80_mul(fp2, fp1); // T(A4+TA6) ! 1657: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBF2A01A01A018B59))); // A3+T(A5+TA7) ! 1658: fp4 = packFloatx80(0, 0x3FF8, LIT64(0x88888888888859AF)); ! 1659: fp2 = floatx80_add(fp2, fp4); // A2+T(A4+TA6) ! 1660: fp1 = floatx80_mul(fp1, fp3); // T(A3+T(A5+TA7)) ! 1661: fp2 = floatx80_mul(fp2, fp0); // S(A2+T(A4+TA6)) ! 1662: fp4 = packFloatx80(1, 0x3FFC, LIT64(0xAAAAAAAAAAAAAA99)); ! 1663: fp1 = floatx80_add(fp1, fp4); // A1+T(A3+T(A5+TA7)) ! 1664: fp1 = floatx80_add(fp1, fp2); // [A1+T(A3+T(A5+TA7))]+[S(A2+T(A4+TA6))] ! 1665: ! 1666: x = packFloatx80(xSign, xExp, xSig); ! 1667: fp0 = floatx80_mul(fp0, x); // R'*S ! 1668: fp0 = floatx80_mul(fp0, fp1); // SIN(R')-R' ! 1669: ! 1670: float_rounding_mode = user_rnd_mode; ! 1671: floatx80_rounding_precision = user_rnd_prec; ! 1672: ! 1673: a = floatx80_add(fp0, x); ! 1674: ! 1675: float_raise(float_flag_inexact); ! 1676: ! 1677: return a; ! 1678: } ! 1679: } ! 1680: } ! 1681: ! 1682: /*---------------------------------------------------------------------------- ! 1683: | Hyperbolic sine ! 1684: *----------------------------------------------------------------------------*/ ! 1685: ! 1686: floatx80 floatx80_sinh(floatx80 a) ! 1687: { ! 1688: flag aSign; ! 1689: int32 aExp; ! 1690: bits64 aSig; ! 1691: ! 1692: int8 user_rnd_mode, user_rnd_prec; ! 1693: ! 1694: int32 compact; ! 1695: floatx80 fp0, fp1, fp2; ! 1696: float32 fact; ! 1697: ! 1698: aSig = extractFloatx80Frac(a); ! 1699: aExp = extractFloatx80Exp(a); ! 1700: aSign = extractFloatx80Sign(a); ! 1701: ! 1702: if (aExp == 0x7FFF) { ! 1703: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 1704: return packFloatx80(aSign, 0x7FFF, floatx80_default_infinity_low); ! 1705: } ! 1706: ! 1707: if (aExp == 0 && aSig == 0) { ! 1708: return packFloatx80(aSign, 0, 0); ! 1709: } ! 1710: ! 1711: user_rnd_mode = float_rounding_mode; ! 1712: user_rnd_prec = floatx80_rounding_precision; ! 1713: float_rounding_mode = float_round_nearest_even; ! 1714: floatx80_rounding_precision = 80; ! 1715: ! 1716: compact = floatx80_make_compact(aExp, aSig); ! 1717: ! 1718: if (compact > 0x400CB167) { ! 1719: // SINHBIG ! 1720: if (compact > 0x400CB2B3) { ! 1721: float_rounding_mode = user_rnd_mode; ! 1722: floatx80_rounding_precision = user_rnd_prec; ! 1723: ! 1724: return roundAndPackFloatx80(floatx80_rounding_precision, aSign, 0x8000, aSig, 0); ! 1725: } else { ! 1726: fp0 = floatx80_abs(a); // Y = |X| ! 1727: fp0 = floatx80_sub(fp0, float64_to_floatx80(LIT64(0x40C62D38D3D64634))); // (|X|-16381LOG2_LEAD) ! 1728: fp0 = floatx80_sub(fp0, float64_to_floatx80(LIT64(0x3D6F90AEB1E75CC7))); // |X| - 16381 LOG2, ACCURATE ! 1729: fp0 = floatx80_etox(fp0); ! 1730: fp2 = packFloatx80(aSign, 0x7FFB, one_sig); ! 1731: ! 1732: float_rounding_mode = user_rnd_mode; ! 1733: floatx80_rounding_precision = user_rnd_prec; ! 1734: ! 1735: a = floatx80_mul(fp0, fp2); ! 1736: ! 1737: float_raise(float_flag_inexact); ! 1738: ! 1739: return a; ! 1740: } ! 1741: } else { // |X| < 16380 LOG2 ! 1742: fp0 = floatx80_abs(a); // Y = |X| ! 1743: fp0 = floatx80_etoxm1(fp0); // FP0 IS Z = EXPM1(Y) ! 1744: fp1 = floatx80_add(fp0, float32_to_floatx80(0x3F800000)); // 1+Z ! 1745: fp2 = fp0; ! 1746: fp0 = floatx80_div(fp0, fp1); // Z/(1+Z) ! 1747: fp0 = floatx80_add(fp0, fp2); ! 1748: ! 1749: fact = 0x3F000000; ! 1750: fact |= aSign ? 0x80000000 : 0x00000000; ! 1751: ! 1752: float_rounding_mode = user_rnd_mode; ! 1753: floatx80_rounding_precision = user_rnd_prec; ! 1754: ! 1755: a = floatx80_mul(fp0, float32_to_floatx80(fact)); ! 1756: ! 1757: float_raise(float_flag_inexact); ! 1758: ! 1759: return a; ! 1760: } ! 1761: } ! 1762: ! 1763: /*---------------------------------------------------------------------------- ! 1764: | Tangent ! 1765: *----------------------------------------------------------------------------*/ ! 1766: ! 1767: floatx80 floatx80_tan(floatx80 a) ! 1768: { ! 1769: flag aSign, xSign; ! 1770: int32 aExp, xExp; ! 1771: bits64 aSig, xSig; ! 1772: ! 1773: int8 user_rnd_mode, user_rnd_prec; ! 1774: ! 1775: int32 compact, l, n, j; ! 1776: floatx80 fp0, fp1, fp2, fp3, fp4, fp5, invtwopi, twopi1, twopi2; ! 1777: float32 twoto63; ! 1778: flag endflag; ! 1779: ! 1780: aSig = extractFloatx80Frac(a); ! 1781: aExp = extractFloatx80Exp(a); ! 1782: aSign = extractFloatx80Sign(a); ! 1783: ! 1784: if (aExp == 0x7FFF) { ! 1785: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 1786: float_raise(float_flag_invalid); ! 1787: a.low = floatx80_default_nan_low; ! 1788: a.high = floatx80_default_nan_high; ! 1789: return a; ! 1790: } ! 1791: ! 1792: if (aExp == 0 && aSig == 0) { ! 1793: return packFloatx80(aSign, 0, 0); ! 1794: } ! 1795: ! 1796: user_rnd_mode = float_rounding_mode; ! 1797: user_rnd_prec = floatx80_rounding_precision; ! 1798: float_rounding_mode = float_round_nearest_even; ! 1799: floatx80_rounding_precision = 80; ! 1800: ! 1801: compact = floatx80_make_compact(aExp, aSig); ! 1802: ! 1803: fp0 = a; ! 1804: ! 1805: if (compact < 0x3FD78000 || compact > 0x4004BC7E) { // 2^(-40) > |X| > 15 PI ! 1806: if (compact > 0x3FFF8000) { // |X| >= 15 PI ! 1807: // REDUCEX ! 1808: fp1 = packFloatx80(0, 0, 0); ! 1809: if (compact == 0x7FFEFFFF) { ! 1810: twopi1 = packFloatx80(aSign ^ 1, 0x7FFE, LIT64(0xC90FDAA200000000)); ! 1811: twopi2 = packFloatx80(aSign ^ 1, 0x7FDC, LIT64(0x85A308D300000000)); ! 1812: fp0 = floatx80_add(fp0, twopi1); ! 1813: fp1 = fp0; ! 1814: fp0 = floatx80_add(fp0, twopi2); ! 1815: fp1 = floatx80_sub(fp1, fp0); ! 1816: fp1 = floatx80_add(fp1, twopi2); ! 1817: } ! 1818: loop: ! 1819: xSign = extractFloatx80Sign(fp0); ! 1820: xExp = extractFloatx80Exp(fp0); ! 1821: xExp -= 0x3FFF; ! 1822: if (xExp <= 28) { ! 1823: l = 0; ! 1824: endflag = 1; ! 1825: } else { ! 1826: l = xExp - 27; ! 1827: endflag = 0; ! 1828: } ! 1829: invtwopi = packFloatx80(0, 0x3FFE - l, LIT64(0xA2F9836E4E44152A)); // INVTWOPI ! 1830: twopi1 = packFloatx80(0, 0x3FFF + l, LIT64(0xC90FDAA200000000)); ! 1831: twopi2 = packFloatx80(0, 0x3FDD + l, LIT64(0x85A308D300000000)); ! 1832: ! 1833: twoto63 = 0x5F000000; ! 1834: twoto63 |= xSign ? 0x80000000 : 0x00000000; // SIGN(INARG)*2^63 IN SGL ! 1835: ! 1836: fp2 = floatx80_mul(fp0, invtwopi); ! 1837: fp2 = floatx80_add(fp2, float32_to_floatx80(twoto63)); // THE FRACTIONAL PART OF FP2 IS ROUNDED ! 1838: fp2 = floatx80_sub(fp2, float32_to_floatx80(twoto63)); // FP2 is N ! 1839: fp4 = floatx80_mul(twopi1, fp2); // W = N*P1 ! 1840: fp5 = floatx80_mul(twopi2, fp2); // w = N*P2 ! 1841: fp3 = floatx80_add(fp4, fp5); // FP3 is P ! 1842: fp4 = floatx80_sub(fp4, fp3); // W-P ! 1843: fp0 = floatx80_sub(fp0, fp3); // FP0 is A := R - P ! 1844: fp4 = floatx80_add(fp4, fp5); // FP4 is p = (W-P)+w ! 1845: fp3 = fp0; // FP3 is A ! 1846: fp1 = floatx80_sub(fp1, fp4); // FP1 is a := r - p ! 1847: fp0 = floatx80_add(fp0, fp1); // FP0 is R := A+a ! 1848: ! 1849: if (endflag > 0) { ! 1850: n = floatx80_to_int32(fp2); ! 1851: goto tancont; ! 1852: } ! 1853: fp3 = floatx80_sub(fp3, fp0); // A-R ! 1854: fp1 = floatx80_add(fp1, fp3); // FP1 is r := (A-R)+a ! 1855: goto loop; ! 1856: } else { ! 1857: float_rounding_mode = user_rnd_mode; ! 1858: floatx80_rounding_precision = user_rnd_prec; ! 1859: ! 1860: a = floatx80_move(a); ! 1861: ! 1862: float_raise(float_flag_inexact); ! 1863: ! 1864: return a; ! 1865: } ! 1866: } else { ! 1867: fp1 = floatx80_mul(fp0, float64_to_floatx80(LIT64(0x3FE45F306DC9C883))); // X*2/PI ! 1868: ! 1869: n = floatx80_to_int32(fp1); ! 1870: j = 32 + n; ! 1871: ! 1872: fp0 = floatx80_sub(fp0, pi_tbl[j]); // X-Y1 ! 1873: fp0 = floatx80_sub(fp0, float32_to_floatx80(pi_tbl2[j])); // FP0 IS R = (X-Y1)-Y2 ! 1874: ! 1875: tancont: ! 1876: if (n & 1) { ! 1877: // NODD ! 1878: fp1 = fp0; // R ! 1879: fp0 = floatx80_mul(fp0, fp0); // S = R*R ! 1880: fp3 = float64_to_floatx80(LIT64(0x3EA0B759F50F8688)); // Q4 ! 1881: fp2 = float64_to_floatx80(LIT64(0xBEF2BAA5A8924F04)); // P3 ! 1882: fp3 = floatx80_mul(fp3, fp0); // SQ4 ! 1883: fp2 = floatx80_mul(fp2, fp0); // SP3 ! 1884: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBF346F59B39BA65F))); // Q3+SQ4 ! 1885: fp4 = packFloatx80(0, 0x3FF6, LIT64(0xE073D3FC199C4A00)); ! 1886: fp2 = floatx80_add(fp2, fp4); // P2+SP3 ! 1887: fp3 = floatx80_mul(fp3, fp0); // S(Q3+SQ4) ! 1888: fp2 = floatx80_mul(fp2, fp0); // S(P2+SP3) ! 1889: fp4 = packFloatx80(0, 0x3FF9, LIT64(0xD23CD68415D95FA1)); ! 1890: fp3 = floatx80_add(fp3, fp4); // Q2+S(Q3+SQ4) ! 1891: fp4 = packFloatx80(1, 0x3FFC, LIT64(0x8895A6C5FB423BCA)); ! 1892: fp2 = floatx80_add(fp2, fp4); // P1+S(P2+SP3) ! 1893: fp3 = floatx80_mul(fp3, fp0); // S(Q2+S(Q3+SQ4)) ! 1894: fp2 = floatx80_mul(fp2, fp0); // S(P1+S(P2+SP3)) ! 1895: fp4 = packFloatx80(1, 0x3FFD, LIT64(0xEEF57E0DA84BC8CE)); ! 1896: fp3 = floatx80_add(fp3, fp4); // Q1+S(Q2+S(Q3+SQ4)) ! 1897: fp2 = floatx80_mul(fp2, fp1); // RS(P1+S(P2+SP3)) ! 1898: fp0 = floatx80_mul(fp0, fp3); // S(Q1+S(Q2+S(Q3+SQ4))) ! 1899: fp1 = floatx80_add(fp1, fp2); // R+RS(P1+S(P2+SP3)) ! 1900: fp0 = floatx80_add(fp0, float32_to_floatx80(0x3F800000)); // 1+S(Q1+S(Q2+S(Q3+SQ4))) ! 1901: ! 1902: xSign = extractFloatx80Sign(fp1); ! 1903: xExp = extractFloatx80Exp(fp1); ! 1904: xSig = extractFloatx80Frac(fp1); ! 1905: xSign ^= 1; ! 1906: fp1 = packFloatx80(xSign, xExp, xSig); ! 1907: ! 1908: float_rounding_mode = user_rnd_mode; ! 1909: floatx80_rounding_precision = user_rnd_prec; ! 1910: ! 1911: a = floatx80_div(fp0, fp1); ! 1912: ! 1913: float_raise(float_flag_inexact); ! 1914: ! 1915: return a; ! 1916: } else { ! 1917: fp1 = floatx80_mul(fp0, fp0); // S = R*R ! 1918: fp3 = float64_to_floatx80(LIT64(0x3EA0B759F50F8688)); // Q4 ! 1919: fp2 = float64_to_floatx80(LIT64(0xBEF2BAA5A8924F04)); // P3 ! 1920: fp3 = floatx80_mul(fp3, fp1); // SQ4 ! 1921: fp2 = floatx80_mul(fp2, fp1); // SP3 ! 1922: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0xBF346F59B39BA65F))); // Q3+SQ4 ! 1923: fp4 = packFloatx80(0, 0x3FF6, LIT64(0xE073D3FC199C4A00)); ! 1924: fp2 = floatx80_add(fp2, fp4); // P2+SP3 ! 1925: fp3 = floatx80_mul(fp3, fp1); // S(Q3+SQ4) ! 1926: fp2 = floatx80_mul(fp2, fp1); // S(P2+SP3) ! 1927: fp4 = packFloatx80(0, 0x3FF9, LIT64(0xD23CD68415D95FA1)); ! 1928: fp3 = floatx80_add(fp3, fp4); // Q2+S(Q3+SQ4) ! 1929: fp4 = packFloatx80(1, 0x3FFC, LIT64(0x8895A6C5FB423BCA)); ! 1930: fp2 = floatx80_add(fp2, fp4); // P1+S(P2+SP3) ! 1931: fp3 = floatx80_mul(fp3, fp1); // S(Q2+S(Q3+SQ4)) ! 1932: fp2 = floatx80_mul(fp2, fp1); // S(P1+S(P2+SP3)) ! 1933: fp4 = packFloatx80(1, 0x3FFD, LIT64(0xEEF57E0DA84BC8CE)); ! 1934: fp3 = floatx80_add(fp3, fp4); // Q1+S(Q2+S(Q3+SQ4)) ! 1935: fp2 = floatx80_mul(fp2, fp0); // RS(P1+S(P2+SP3)) ! 1936: fp1 = floatx80_mul(fp1, fp3); // S(Q1+S(Q2+S(Q3+SQ4))) ! 1937: fp0 = floatx80_add(fp0, fp2); // R+RS(P1+S(P2+SP3)) ! 1938: fp1 = floatx80_add(fp1, float32_to_floatx80(0x3F800000)); // 1+S(Q1+S(Q2+S(Q3+SQ4))) ! 1939: ! 1940: float_rounding_mode = user_rnd_mode; ! 1941: floatx80_rounding_precision = user_rnd_prec; ! 1942: ! 1943: a = floatx80_div(fp0, fp1); ! 1944: ! 1945: float_raise(float_flag_inexact); ! 1946: ! 1947: return a; ! 1948: } ! 1949: } ! 1950: } ! 1951: ! 1952: /*---------------------------------------------------------------------------- ! 1953: | Hyperbolic tangent ! 1954: *----------------------------------------------------------------------------*/ ! 1955: ! 1956: floatx80 floatx80_tanh(floatx80 a) ! 1957: { ! 1958: flag aSign, vSign; ! 1959: int32 aExp, vExp; ! 1960: bits64 aSig, vSig; ! 1961: ! 1962: int8 user_rnd_mode, user_rnd_prec; ! 1963: ! 1964: int32 compact; ! 1965: floatx80 fp0, fp1; ! 1966: float32 sign; ! 1967: ! 1968: aSig = extractFloatx80Frac(a); ! 1969: aExp = extractFloatx80Exp(a); ! 1970: aSign = extractFloatx80Sign(a); ! 1971: ! 1972: if (aExp == 0x7FFF) { ! 1973: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 1974: return packFloatx80(aSign, one_exp, one_sig); ! 1975: } ! 1976: ! 1977: if (aExp == 0 && aSig == 0) { ! 1978: return packFloatx80(aSign, 0, 0); ! 1979: } ! 1980: ! 1981: user_rnd_mode = float_rounding_mode; ! 1982: user_rnd_prec = floatx80_rounding_precision; ! 1983: float_rounding_mode = float_round_nearest_even; ! 1984: floatx80_rounding_precision = 80; ! 1985: ! 1986: compact = floatx80_make_compact(aExp, aSig); ! 1987: ! 1988: if (compact < 0x3FD78000 || compact > 0x3FFFDDCE) { ! 1989: // TANHBORS ! 1990: if (compact < 0x3FFF8000) { ! 1991: // TANHSM ! 1992: float_rounding_mode = user_rnd_mode; ! 1993: floatx80_rounding_precision = user_rnd_prec; ! 1994: ! 1995: a = floatx80_move(a); ! 1996: ! 1997: float_raise(float_flag_inexact); ! 1998: ! 1999: return a; ! 2000: } else { ! 2001: if (compact > 0x40048AA1) { ! 2002: // TANHHUGE ! 2003: sign = 0x3F800000; ! 2004: sign |= aSign ? 0x80000000 : 0x00000000; ! 2005: fp0 = float32_to_floatx80(sign); ! 2006: sign &= 0x80000000; ! 2007: sign ^= 0x80800000; // -SIGN(X)*EPS ! 2008: ! 2009: float_rounding_mode = user_rnd_mode; ! 2010: floatx80_rounding_precision = user_rnd_prec; ! 2011: ! 2012: a = floatx80_add(fp0, float32_to_floatx80(sign)); ! 2013: ! 2014: float_raise(float_flag_inexact); ! 2015: ! 2016: return a; ! 2017: } else { ! 2018: fp0 = packFloatx80(0, aExp+1, aSig); // Y = 2|X| ! 2019: fp0 = floatx80_etox(fp0); // FP0 IS EXP(Y) ! 2020: fp0 = floatx80_add(fp0, float32_to_floatx80(0x3F800000)); // EXP(Y)+1 ! 2021: sign = aSign ? 0x80000000 : 0x00000000; ! 2022: fp1 = floatx80_div(float32_to_floatx80(sign^0xC0000000), fp0); // -SIGN(X)*2 / [EXP(Y)+1] ! 2023: fp0 = float32_to_floatx80(sign | 0x3F800000); // SIGN ! 2024: ! 2025: float_rounding_mode = user_rnd_mode; ! 2026: floatx80_rounding_precision = user_rnd_prec; ! 2027: ! 2028: a = floatx80_add(fp1, fp0); ! 2029: ! 2030: float_raise(float_flag_inexact); ! 2031: ! 2032: return a; ! 2033: } ! 2034: } ! 2035: } else { // 2**(-40) < |X| < (5/2)LOG2 ! 2036: fp0 = packFloatx80(0, aExp+1, aSig); // Y = 2|X| ! 2037: fp0 = floatx80_etoxm1(fp0); // FP0 IS Z = EXPM1(Y) ! 2038: fp1 = floatx80_add(fp0, float32_to_floatx80(0x40000000)); // Z+2 ! 2039: ! 2040: vSign = extractFloatx80Sign(fp1); ! 2041: vExp = extractFloatx80Exp(fp1); ! 2042: vSig = extractFloatx80Frac(fp1); ! 2043: ! 2044: fp1 = packFloatx80(vSign ^ aSign, vExp, vSig); ! 2045: ! 2046: float_rounding_mode = user_rnd_mode; ! 2047: floatx80_rounding_precision = user_rnd_prec; ! 2048: ! 2049: a = floatx80_div(fp0, fp1); ! 2050: ! 2051: float_raise(float_flag_inexact); ! 2052: ! 2053: return a; ! 2054: } ! 2055: } ! 2056: ! 2057: /*---------------------------------------------------------------------------- ! 2058: | 10 to x ! 2059: *----------------------------------------------------------------------------*/ ! 2060: ! 2061: floatx80 floatx80_tentox(floatx80 a) ! 2062: { ! 2063: flag aSign; ! 2064: int32 aExp; ! 2065: bits64 aSig; ! 2066: ! 2067: int8 user_rnd_mode, user_rnd_prec; ! 2068: ! 2069: int32 compact, n, j, l, m, m1; ! 2070: floatx80 fp0, fp1, fp2, fp3, adjfact, fact1, fact2; ! 2071: ! 2072: aSig = extractFloatx80Frac(a); ! 2073: aExp = extractFloatx80Exp(a); ! 2074: aSign = extractFloatx80Sign(a); ! 2075: ! 2076: if (aExp == 0x7FFF) { ! 2077: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 2078: if (aSign) return packFloatx80(0, 0, 0); ! 2079: return packFloatx80(0, 0x7FFF, floatx80_default_infinity_low); ! 2080: } ! 2081: ! 2082: if (aExp == 0 && aSig == 0) { ! 2083: return packFloatx80(0, one_exp, one_sig); ! 2084: } ! 2085: ! 2086: user_rnd_mode = float_rounding_mode; ! 2087: user_rnd_prec = floatx80_rounding_precision; ! 2088: float_rounding_mode = float_round_nearest_even; ! 2089: floatx80_rounding_precision = 80; ! 2090: ! 2091: fp0 = a; ! 2092: ! 2093: compact = floatx80_make_compact(aExp, aSig); ! 2094: ! 2095: if (compact < 0x3FB98000 || compact > 0x400B9B07) { // |X| > 16480 LOG2/LOG10 or |X| < 2^(-70) ! 2096: if (compact > 0x3FFF8000) { // |X| > 16480 ! 2097: float_rounding_mode = user_rnd_mode; ! 2098: floatx80_rounding_precision = user_rnd_prec; ! 2099: ! 2100: if (aSign) { ! 2101: return roundAndPackFloatx80(floatx80_rounding_precision, 0, -0x1000, aSig, 0); ! 2102: } else { ! 2103: return roundAndPackFloatx80(floatx80_rounding_precision, 0, 0x8000, aSig, 0); ! 2104: } ! 2105: } else { // |X| < 2^(-70) ! 2106: float_rounding_mode = user_rnd_mode; ! 2107: floatx80_rounding_precision = user_rnd_prec; ! 2108: ! 2109: a = floatx80_add(fp0, float32_to_floatx80(0x3F800000)); // 1 + X ! 2110: ! 2111: float_raise(float_flag_inexact); ! 2112: ! 2113: return a; ! 2114: } ! 2115: } else { // 2^(-70) <= |X| <= 16480 LOG 2 / LOG 10 ! 2116: fp1 = fp0; // X ! 2117: fp1 = floatx80_mul(fp1, float64_to_floatx80(LIT64(0x406A934F0979A371))); // X*64*LOG10/LOG2 ! 2118: n = floatx80_to_int32(fp1); // N=INT(X*64*LOG10/LOG2) ! 2119: fp1 = int32_to_floatx80(n); ! 2120: ! 2121: j = n & 0x3F; ! 2122: l = n / 64; // NOTE: this is really arithmetic right shift by 6 ! 2123: if (n < 0 && j) { // arithmetic right shift is division and round towards minus infinity ! 2124: l--; ! 2125: } ! 2126: m = l / 2; // NOTE: this is really arithmetic right shift by 1 ! 2127: if (l < 0 && (l & 1)) { // arithmetic right shift is division and round towards minus infinity ! 2128: m--; ! 2129: } ! 2130: m1 = l - m; ! 2131: m1 += 0x3FFF; // ADJFACT IS 2^(M') ! 2132: ! 2133: adjfact = packFloatx80(0, m1, one_sig); ! 2134: fact1 = exp2_tbl[j]; ! 2135: fact1.high += m; ! 2136: fact2.high = exp2_tbl2[j]>>16; ! 2137: fact2.high += m; ! 2138: fact2.low = (bits64)(exp2_tbl2[j] & 0xFFFF); ! 2139: fact2.low <<= 48; ! 2140: ! 2141: fp2 = fp1; // N ! 2142: fp1 = floatx80_mul(fp1, float64_to_floatx80(LIT64(0x3F734413509F8000))); // N*(LOG2/64LOG10)_LEAD ! 2143: fp3 = packFloatx80(1, 0x3FCD, LIT64(0xC0219DC1DA994FD2)); ! 2144: fp2 = floatx80_mul(fp2, fp3); // N*(LOG2/64LOG10)_TRAIL ! 2145: fp0 = floatx80_sub(fp0, fp1); // X - N L_LEAD ! 2146: fp0 = floatx80_sub(fp0, fp2); // X - N L_TRAIL ! 2147: fp2 = packFloatx80(0, 0x4000, LIT64(0x935D8DDDAAA8AC17)); // LOG10 ! 2148: fp0 = floatx80_mul(fp0, fp2); // R ! 2149: ! 2150: // EXPR ! 2151: fp1 = floatx80_mul(fp0, fp0); // S = R*R ! 2152: fp2 = float64_to_floatx80(LIT64(0x3F56C16D6F7BD0B2)); // A5 ! 2153: fp3 = float64_to_floatx80(LIT64(0x3F811112302C712C)); // A4 ! 2154: fp2 = floatx80_mul(fp2, fp1); // S*A5 ! 2155: fp3 = floatx80_mul(fp3, fp1); // S*A4 ! 2156: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FA5555555554CC1))); // A3+S*A5 ! 2157: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0x3FC5555555554A54))); // A2+S*A4 ! 2158: fp2 = floatx80_mul(fp2, fp1); // S*(A3+S*A5) ! 2159: fp3 = floatx80_mul(fp3, fp1); // S*(A2+S*A4) ! 2160: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FE0000000000000))); // A1+S*(A3+S*A5) ! 2161: fp3 = floatx80_mul(fp3, fp0); // R*S*(A2+S*A4) ! 2162: ! 2163: fp2 = floatx80_mul(fp2, fp1); // S*(A1+S*(A3+S*A5)) ! 2164: fp0 = floatx80_add(fp0, fp3); // R+R*S*(A2+S*A4) ! 2165: fp0 = floatx80_add(fp0, fp2); // EXP(R) - 1 ! 2166: ! 2167: fp0 = floatx80_mul(fp0, fact1); ! 2168: fp0 = floatx80_add(fp0, fact2); ! 2169: fp0 = floatx80_add(fp0, fact1); ! 2170: ! 2171: float_rounding_mode = user_rnd_mode; ! 2172: floatx80_rounding_precision = user_rnd_prec; ! 2173: ! 2174: a = floatx80_mul(fp0, adjfact); ! 2175: ! 2176: float_raise(float_flag_inexact); ! 2177: ! 2178: return a; ! 2179: } ! 2180: } ! 2181: ! 2182: /*---------------------------------------------------------------------------- ! 2183: | 2 to x ! 2184: *----------------------------------------------------------------------------*/ ! 2185: ! 2186: floatx80 floatx80_twotox(floatx80 a) ! 2187: { ! 2188: flag aSign; ! 2189: int32 aExp; ! 2190: bits64 aSig; ! 2191: ! 2192: int8 user_rnd_mode, user_rnd_prec; ! 2193: ! 2194: int32 compact, n, j, l, m, m1; ! 2195: floatx80 fp0, fp1, fp2, fp3, adjfact, fact1, fact2; ! 2196: ! 2197: aSig = extractFloatx80Frac(a); ! 2198: aExp = extractFloatx80Exp(a); ! 2199: aSign = extractFloatx80Sign(a); ! 2200: ! 2201: if (aExp == 0x7FFF) { ! 2202: if ((bits64) (aSig<<1)) return propagateFloatx80NaNOneArg(a); ! 2203: if (aSign) return packFloatx80(0, 0, 0); ! 2204: return packFloatx80(0, 0x7FFF, floatx80_default_infinity_low); ! 2205: } ! 2206: ! 2207: if (aExp == 0 && aSig == 0) { ! 2208: return packFloatx80(0, one_exp, one_sig); ! 2209: } ! 2210: ! 2211: user_rnd_mode = float_rounding_mode; ! 2212: user_rnd_prec = floatx80_rounding_precision; ! 2213: float_rounding_mode = float_round_nearest_even; ! 2214: floatx80_rounding_precision = 80; ! 2215: ! 2216: fp0 = a; ! 2217: ! 2218: compact = floatx80_make_compact(aExp, aSig); ! 2219: ! 2220: if (compact < 0x3FB98000 || compact > 0x400D80C0) { // |X| > 16480 or |X| < 2^(-70) ! 2221: if (compact > 0x3FFF8000) { // |X| > 16480 ! 2222: float_rounding_mode = user_rnd_mode; ! 2223: floatx80_rounding_precision = user_rnd_prec; ! 2224: ! 2225: if (aSign) { ! 2226: return roundAndPackFloatx80(floatx80_rounding_precision, 0, -0x1000, aSig, 0); ! 2227: } else { ! 2228: return roundAndPackFloatx80(floatx80_rounding_precision, 0, 0x8000, aSig, 0); ! 2229: } ! 2230: } else { // |X| < 2^(-70) ! 2231: float_rounding_mode = user_rnd_mode; ! 2232: floatx80_rounding_precision = user_rnd_prec; ! 2233: ! 2234: a = floatx80_add(fp0, float32_to_floatx80(0x3F800000)); // 1 + X ! 2235: ! 2236: float_raise(float_flag_inexact); ! 2237: ! 2238: return a; ! 2239: } ! 2240: } else { // 2^(-70) <= |X| <= 16480 ! 2241: fp1 = fp0; // X ! 2242: fp1 = floatx80_mul(fp1, float32_to_floatx80(0x42800000)); // X * 64 ! 2243: n = floatx80_to_int32(fp1); ! 2244: fp1 = int32_to_floatx80(n); ! 2245: j = n & 0x3F; ! 2246: l = n / 64; // NOTE: this is really arithmetic right shift by 6 ! 2247: if (n < 0 && j) { // arithmetic right shift is division and round towards minus infinity ! 2248: l--; ! 2249: } ! 2250: m = l / 2; // NOTE: this is really arithmetic right shift by 1 ! 2251: if (l < 0 && (l & 1)) { // arithmetic right shift is division and round towards minus infinity ! 2252: m--; ! 2253: } ! 2254: m1 = l - m; ! 2255: m1 += 0x3FFF; // ADJFACT IS 2^(M') ! 2256: ! 2257: adjfact = packFloatx80(0, m1, one_sig); ! 2258: fact1 = exp2_tbl[j]; ! 2259: fact1.high += m; ! 2260: fact2.high = exp2_tbl2[j]>>16; ! 2261: fact2.high += m; ! 2262: fact2.low = (bits64)(exp2_tbl2[j] & 0xFFFF); ! 2263: fact2.low <<= 48; ! 2264: ! 2265: fp1 = floatx80_mul(fp1, float32_to_floatx80(0x3C800000)); // (1/64)*N ! 2266: fp0 = floatx80_sub(fp0, fp1); // X - (1/64)*INT(64 X) ! 2267: fp2 = packFloatx80(0, 0x3FFE, LIT64(0xB17217F7D1CF79AC)); // LOG2 ! 2268: fp0 = floatx80_mul(fp0, fp2); // R ! 2269: ! 2270: // EXPR ! 2271: fp1 = floatx80_mul(fp0, fp0); // S = R*R ! 2272: fp2 = float64_to_floatx80(LIT64(0x3F56C16D6F7BD0B2)); // A5 ! 2273: fp3 = float64_to_floatx80(LIT64(0x3F811112302C712C)); // A4 ! 2274: fp2 = floatx80_mul(fp2, fp1); // S*A5 ! 2275: fp3 = floatx80_mul(fp3, fp1); // S*A4 ! 2276: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FA5555555554CC1))); // A3+S*A5 ! 2277: fp3 = floatx80_add(fp3, float64_to_floatx80(LIT64(0x3FC5555555554A54))); // A2+S*A4 ! 2278: fp2 = floatx80_mul(fp2, fp1); // S*(A3+S*A5) ! 2279: fp3 = floatx80_mul(fp3, fp1); // S*(A2+S*A4) ! 2280: fp2 = floatx80_add(fp2, float64_to_floatx80(LIT64(0x3FE0000000000000))); // A1+S*(A3+S*A5) ! 2281: fp3 = floatx80_mul(fp3, fp0); // R*S*(A2+S*A4) ! 2282: ! 2283: fp2 = floatx80_mul(fp2, fp1); // S*(A1+S*(A3+S*A5)) ! 2284: fp0 = floatx80_add(fp0, fp3); // R+R*S*(A2+S*A4) ! 2285: fp0 = floatx80_add(fp0, fp2); // EXP(R) - 1 ! 2286: ! 2287: fp0 = floatx80_mul(fp0, fact1); ! 2288: fp0 = floatx80_add(fp0, fact2); ! 2289: fp0 = floatx80_add(fp0, fact1); ! 2290: ! 2291: float_rounding_mode = user_rnd_mode; ! 2292: floatx80_rounding_precision = user_rnd_prec; ! 2293: ! 2294: a = floatx80_mul(fp0, adjfact); ! 2295: ! 2296: float_raise(float_flag_inexact); ! 2297: ! 2298: return a; ! 2299: } ! 2300: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.