1 // SPDX-License-Identifier: GPL-2.0 2 /*---------------------------------------------------------------------------+ 3 | fpu_trig.c | 4 | | 5 | Implementation of the FPU "transcendental" functions. | 6 | | 7 | Copyright (C) 1992,1993,1994,1997,1999 | 8 | W. Metzenthen, 22 Parker St, Ormond, Vic 3163, | 9 | Australia. E-mail billm@melbpc.org.au | 10 | | 11 | | 12 +---------------------------------------------------------------------------*/ 13 14 #include "fpu_system.h" 15 #include "exception.h" 16 #include "fpu_emu.h" 17 #include "status_w.h" 18 #include "control_w.h" 19 #include "reg_constant.h" 20 21 static void rem_kernel(unsigned long long st0, unsigned long long *y, 22 unsigned long long st1, unsigned long long q, int n); 23 24 #define BETTER_THAN_486 25 26 #define FCOS 4 27 28 /* Used only by fptan, fsin, fcos, and fsincos. */ 29 /* This routine produces very accurate results, similar to 30 using a value of pi with more than 128 bits precision. */ 31 /* Limited measurements show no results worse than 64 bit precision 32 except for the results for arguments close to 2^63, where the 33 precision of the result sometimes degrades to about 63.9 bits */ 34 static int trig_arg(FPU_REG *st0_ptr, int even) 35 { 36 FPU_REG tmp; 37 u_char tmptag; 38 unsigned long long q; 39 int old_cw = control_word, saved_status = partial_status; 40 int tag, st0_tag = TAG_Valid; 41 42 if (exponent(st0_ptr) >= 63) { 43 partial_status |= SW_C2; /* Reduction incomplete. */ 44 return -1; 45 } 46 47 control_word &= ~CW_RC; 48 control_word |= RC_CHOP; 49 50 setpositive(st0_ptr); 51 tag = FPU_u_div(st0_ptr, &CONST_PI2, &tmp, PR_64_BITS | RC_CHOP | 0x3f, 52 SIGN_POS); 53 54 FPU_round_to_int(&tmp, tag); /* Fortunately, this can't overflow 55 to 2^64 */ 56 q = significand(&tmp); 57 if (q) { 58 rem_kernel(significand(st0_ptr), 59 &significand(&tmp), 60 significand(&CONST_PI2), 61 q, exponent(st0_ptr) - exponent(&CONST_PI2)); 62 setexponent16(&tmp, exponent(&CONST_PI2)); 63 st0_tag = FPU_normalize(&tmp); 64 FPU_copy_to_reg0(&tmp, st0_tag); 65 } 66 67 if ((even && !(q & 1)) || (!even && (q & 1))) { 68 st0_tag = 69 FPU_sub(REV | LOADED | TAG_Valid, (int)&CONST_PI2, 70 FULL_PRECISION); 71 72 #ifdef BETTER_THAN_486 73 /* So far, the results are exact but based upon a 64 bit 74 precision approximation to pi/2. The technique used 75 now is equivalent to using an approximation to pi/2 which 76 is accurate to about 128 bits. */ 77 if ((exponent(st0_ptr) <= exponent(&CONST_PI2extra) + 64) 78 || (q > 1)) { 79 /* This code gives the effect of having pi/2 to better than 80 128 bits precision. */ 81 82 significand(&tmp) = q + 1; 83 setexponent16(&tmp, 63); 84 FPU_normalize(&tmp); 85 tmptag = 86 FPU_u_mul(&CONST_PI2extra, &tmp, &tmp, 87 FULL_PRECISION, SIGN_POS, 88 exponent(&CONST_PI2extra) + 89 exponent(&tmp)); 90 setsign(&tmp, getsign(&CONST_PI2extra)); 91 st0_tag = FPU_add(&tmp, tmptag, 0, FULL_PRECISION); 92 if (signnegative(st0_ptr)) { 93 /* CONST_PI2extra is negative, so the result of the addition 94 can be negative. This means that the argument is actually 95 in a different quadrant. The correction is always < pi/2, 96 so it can't overflow into yet another quadrant. */ 97 setpositive(st0_ptr); 98 q++; 99 } 100 } 101 #endif /* BETTER_THAN_486 */ 102 } 103 #ifdef BETTER_THAN_486 104 else { 105 /* So far, the results are exact but based upon a 64 bit 106 precision approximation to pi/2. The technique used 107 now is equivalent to using an approximation to pi/2 which 108 is accurate to about 128 bits. */ 109 if (((q > 0) 110 && (exponent(st0_ptr) <= exponent(&CONST_PI2extra) + 64)) 111 || (q > 1)) { 112 /* This code gives the effect of having p/2 to better than 113 128 bits precision. */ 114 115 significand(&tmp) = q; 116 setexponent16(&tmp, 63); 117 FPU_normalize(&tmp); /* This must return TAG_Valid */ 118 tmptag = 119 FPU_u_mul(&CONST_PI2extra, &tmp, &tmp, 120 FULL_PRECISION, SIGN_POS, 121 exponent(&CONST_PI2extra) + 122 exponent(&tmp)); 123 setsign(&tmp, getsign(&CONST_PI2extra)); 124 st0_tag = FPU_sub(LOADED | (tmptag & 0x0f), (int)&tmp, 125 FULL_PRECISION); 126 if ((exponent(st0_ptr) == exponent(&CONST_PI2)) && 127 ((st0_ptr->sigh > CONST_PI2.sigh) 128 || ((st0_ptr->sigh == CONST_PI2.sigh) 129 && (st0_ptr->sigl > CONST_PI2.sigl)))) { 130 /* CONST_PI2extra is negative, so the result of the 131 subtraction can be larger than pi/2. This means 132 that the argument is actually in a different quadrant. 133 The correction is always < pi/2, so it can't overflow 134 into yet another quadrant. */ 135 st0_tag = 136 FPU_sub(REV | LOADED | TAG_Valid, 137 (int)&CONST_PI2, FULL_PRECISION); 138 q++; 139 } 140 } 141 } 142 #endif /* BETTER_THAN_486 */ 143 144 FPU_settag0(st0_tag); 145 control_word = old_cw; 146 partial_status = saved_status & ~SW_C2; /* Reduction complete. */ 147 148 return (q & 3) | even; 149 } 150 151 /* Convert a long to register */ 152 static void convert_l2reg(long const *arg, int deststnr) 153 { 154 int tag; 155 long num = *arg; 156 u_char sign; 157 FPU_REG *dest = &st(deststnr); 158 159 if (num == 0) { 160 FPU_copy_to_regi(&CONST_Z, TAG_Zero, deststnr); 161 return; 162 } 163 164 if (num > 0) { 165 sign = SIGN_POS; 166 } else { 167 num = -num; 168 sign = SIGN_NEG; 169 } 170 171 dest->sigh = num; 172 dest->sigl = 0; 173 setexponent16(dest, 31); 174 tag = FPU_normalize(dest); 175 FPU_settagi(deststnr, tag); 176 setsign(dest, sign); 177 return; 178 } 179 180 static void single_arg_error(FPU_REG *st0_ptr, u_char st0_tag) 181 { 182 if (st0_tag == TAG_Empty) 183 FPU_stack_underflow(); /* Puts a QNaN in st(0) */ 184 else if (st0_tag == TW_NaN) 185 real_1op_NaN(st0_ptr); /* return with a NaN in st(0) */ 186 #ifdef PARANOID 187 else 188 EXCEPTION(EX_INTERNAL | 0x0112); 189 #endif /* PARANOID */ 190 } 191 192 static void single_arg_2_error(FPU_REG *st0_ptr, u_char st0_tag) 193 { 194 int isNaN; 195 196 switch (st0_tag) { 197 case TW_NaN: 198 isNaN = (exponent(st0_ptr) == EXP_OVER) 199 && (st0_ptr->sigh & 0x80000000); 200 if (isNaN && !(st0_ptr->sigh & 0x40000000)) { /* Signaling ? */ 201 EXCEPTION(EX_Invalid); 202 if (control_word & CW_Invalid) { 203 /* The masked response */ 204 /* Convert to a QNaN */ 205 st0_ptr->sigh |= 0x40000000; 206 push(); 207 FPU_copy_to_reg0(st0_ptr, TAG_Special); 208 } 209 } else if (isNaN) { 210 /* A QNaN */ 211 push(); 212 FPU_copy_to_reg0(st0_ptr, TAG_Special); 213 } else { 214 /* pseudoNaN or other unsupported */ 215 EXCEPTION(EX_Invalid); 216 if (control_word & CW_Invalid) { 217 /* The masked response */ 218 FPU_copy_to_reg0(&CONST_QNaN, TAG_Special); 219 push(); 220 FPU_copy_to_reg0(&CONST_QNaN, TAG_Special); 221 } 222 } 223 break; /* return with a NaN in st(0) */ 224 #ifdef PARANOID 225 default: 226 EXCEPTION(EX_INTERNAL | 0x0112); 227 #endif /* PARANOID */ 228 } 229 } 230 231 /*---------------------------------------------------------------------------*/ 232 233 static void f2xm1(FPU_REG *st0_ptr, u_char tag) 234 { 235 FPU_REG a; 236 237 clear_C1(); 238 239 if (tag == TAG_Valid) { 240 /* For an 80486 FPU, the result is undefined if the arg is >= 1.0 */ 241 if (exponent(st0_ptr) < 0) { 242 denormal_arg: 243 244 FPU_to_exp16(st0_ptr, &a); 245 246 /* poly_2xm1(x) requires 0 < st(0) < 1. */ 247 poly_2xm1(getsign(st0_ptr), &a, st0_ptr); 248 } 249 set_precision_flag_up(); /* 80486 appears to always do this */ 250 return; 251 } 252 253 if (tag == TAG_Zero) 254 return; 255 256 if (tag == TAG_Special) 257 tag = FPU_Special(st0_ptr); 258 259 switch (tag) { 260 case TW_Denormal: 261 if (denormal_operand() < 0) 262 return; 263 goto denormal_arg; 264 case TW_Infinity: 265 if (signnegative(st0_ptr)) { 266 /* -infinity gives -1 (p16-10) */ 267 FPU_copy_to_reg0(&CONST_1, TAG_Valid); 268 setnegative(st0_ptr); 269 } 270 return; 271 default: 272 single_arg_error(st0_ptr, tag); 273 } 274 } 275 276 static void fptan(FPU_REG *st0_ptr, u_char st0_tag) 277 { 278 FPU_REG *st_new_ptr; 279 int q; 280 u_char arg_sign = getsign(st0_ptr); 281 282 /* Stack underflow has higher priority */ 283 if (st0_tag == TAG_Empty) { 284 FPU_stack_underflow(); /* Puts a QNaN in st(0) */ 285 if (control_word & CW_Invalid) { 286 st_new_ptr = &st(-1); 287 push(); 288 FPU_stack_underflow(); /* Puts a QNaN in the new st(0) */ 289 } 290 return; 291 } 292 293 if (STACK_OVERFLOW) { 294 FPU_stack_overflow(); 295 return; 296 } 297 298 if (st0_tag == TAG_Valid) { 299 if (exponent(st0_ptr) > -40) { 300 if ((q = trig_arg(st0_ptr, 0)) == -1) { 301 /* Operand is out of range */ 302 return; 303 } 304 305 poly_tan(st0_ptr); 306 setsign(st0_ptr, (q & 1) ^ (arg_sign != 0)); 307 set_precision_flag_up(); /* We do not really know if up or down */ 308 } else { 309 /* For a small arg, the result == the argument */ 310 /* Underflow may happen */ 311 312 denormal_arg: 313 314 FPU_to_exp16(st0_ptr, st0_ptr); 315 316 st0_tag = 317 FPU_round(st0_ptr, 1, 0, FULL_PRECISION, arg_sign); 318 FPU_settag0(st0_tag); 319 } 320 push(); 321 FPU_copy_to_reg0(&CONST_1, TAG_Valid); 322 return; 323 } 324 325 if (st0_tag == TAG_Zero) { 326 push(); 327 FPU_copy_to_reg0(&CONST_1, TAG_Valid); 328 setcc(0); 329 return; 330 } 331 332 if (st0_tag == TAG_Special) 333 st0_tag = FPU_Special(st0_ptr); 334 335 if (st0_tag == TW_Denormal) { 336 if (denormal_operand() < 0) 337 return; 338 339 goto denormal_arg; 340 } 341 342 if (st0_tag == TW_Infinity) { 343 /* The 80486 treats infinity as an invalid operand */ 344 if (arith_invalid(0) >= 0) { 345 st_new_ptr = &st(-1); 346 push(); 347 arith_invalid(0); 348 } 349 return; 350 } 351 352 single_arg_2_error(st0_ptr, st0_tag); 353 } 354 355 static void fxtract(FPU_REG *st0_ptr, u_char st0_tag) 356 { 357 FPU_REG *st_new_ptr; 358 u_char sign; 359 register FPU_REG *st1_ptr = st0_ptr; /* anticipate */ 360 361 if (STACK_OVERFLOW) { 362 FPU_stack_overflow(); 363 return; 364 } 365 366 clear_C1(); 367 368 if (st0_tag == TAG_Valid) { 369 long e; 370 371 push(); 372 sign = getsign(st1_ptr); 373 reg_copy(st1_ptr, st_new_ptr); 374 setexponent16(st_new_ptr, exponent(st_new_ptr)); 375 376 denormal_arg: 377 378 e = exponent16(st_new_ptr); 379 convert_l2reg(&e, 1); 380 setexponentpos(st_new_ptr, 0); 381 setsign(st_new_ptr, sign); 382 FPU_settag0(TAG_Valid); /* Needed if arg was a denormal */ 383 return; 384 } else if (st0_tag == TAG_Zero) { 385 sign = getsign(st0_ptr); 386 387 if (FPU_divide_by_zero(0, SIGN_NEG) < 0) 388 return; 389 390 push(); 391 FPU_copy_to_reg0(&CONST_Z, TAG_Zero); 392 setsign(st_new_ptr, sign); 393 return; 394 } 395 396 if (st0_tag == TAG_Special) 397 st0_tag = FPU_Special(st0_ptr); 398 399 if (st0_tag == TW_Denormal) { 400 if (denormal_operand() < 0) 401 return; 402 403 push(); 404 sign = getsign(st1_ptr); 405 FPU_to_exp16(st1_ptr, st_new_ptr); 406 goto denormal_arg; 407 } else if (st0_tag == TW_Infinity) { 408 sign = getsign(st0_ptr); 409 setpositive(st0_ptr); 410 push(); 411 FPU_copy_to_reg0(&CONST_INF, TAG_Special); 412 setsign(st_new_ptr, sign); 413 return; 414 } else if (st0_tag == TW_NaN) { 415 if (real_1op_NaN(st0_ptr) < 0) 416 return; 417 418 push(); 419 FPU_copy_to_reg0(st0_ptr, TAG_Special); 420 return; 421 } else if (st0_tag == TAG_Empty) { 422 /* Is this the correct behaviour? */ 423 if (control_word & EX_Invalid) { 424 FPU_stack_underflow(); 425 push(); 426 FPU_stack_underflow(); 427 } else 428 EXCEPTION(EX_StackUnder); 429 } 430 #ifdef PARANOID 431 else 432 EXCEPTION(EX_INTERNAL | 0x119); 433 #endif /* PARANOID */ 434 } 435 436 static void fdecstp(void) 437 { 438 clear_C1(); 439 top--; 440 } 441 442 static void fincstp(void) 443 { 444 clear_C1(); 445 top++; 446 } 447 448 static void fsqrt_(FPU_REG *st0_ptr, u_char st0_tag) 449 { 450 int expon; 451 452 clear_C1(); 453 454 if (st0_tag == TAG_Valid) { 455 u_char tag; 456 457 if (signnegative(st0_ptr)) { 458 arith_invalid(0); /* sqrt(negative) is invalid */ 459 return; 460 } 461 462 /* make st(0) in [1.0 .. 4.0) */ 463 expon = exponent(st0_ptr); 464 465 denormal_arg: 466 467 setexponent16(st0_ptr, (expon & 1)); 468 469 /* Do the computation, the sign of the result will be positive. */ 470 tag = wm_sqrt(st0_ptr, 0, 0, control_word, SIGN_POS); 471 addexponent(st0_ptr, expon >> 1); 472 FPU_settag0(tag); 473 return; 474 } 475 476 if (st0_tag == TAG_Zero) 477 return; 478 479 if (st0_tag == TAG_Special) 480 st0_tag = FPU_Special(st0_ptr); 481 482 if (st0_tag == TW_Infinity) { 483 if (signnegative(st0_ptr)) 484 arith_invalid(0); /* sqrt(-Infinity) is invalid */ 485 return; 486 } else if (st0_tag == TW_Denormal) { 487 if (signnegative(st0_ptr)) { 488 arith_invalid(0); /* sqrt(negative) is invalid */ 489 return; 490 } 491 492 if (denormal_operand() < 0) 493 return; 494 495 FPU_to_exp16(st0_ptr, st0_ptr); 496 497 expon = exponent16(st0_ptr); 498 499 goto denormal_arg; 500 } 501 502 single_arg_error(st0_ptr, st0_tag); 503 504 } 505 506 static void frndint_(FPU_REG *st0_ptr, u_char st0_tag) 507 { 508 int flags, tag; 509 510 if (st0_tag == TAG_Valid) { 511 u_char sign; 512 513 denormal_arg: 514 515 sign = getsign(st0_ptr); 516 517 if (exponent(st0_ptr) > 63) 518 return; 519 520 if (st0_tag == TW_Denormal) { 521 if (denormal_operand() < 0) 522 return; 523 } 524 525 /* Fortunately, this can't overflow to 2^64 */ 526 if ((flags = FPU_round_to_int(st0_ptr, st0_tag))) 527 set_precision_flag(flags); 528 529 setexponent16(st0_ptr, 63); 530 tag = FPU_normalize(st0_ptr); 531 setsign(st0_ptr, sign); 532 FPU_settag0(tag); 533 return; 534 } 535 536 if (st0_tag == TAG_Zero) 537 return; 538 539 if (st0_tag == TAG_Special) 540 st0_tag = FPU_Special(st0_ptr); 541 542 if (st0_tag == TW_Denormal) 543 goto denormal_arg; 544 else if (st0_tag == TW_Infinity) 545 return; 546 else 547 single_arg_error(st0_ptr, st0_tag); 548 } 549 550 static int f_sin(FPU_REG *st0_ptr, u_char tag) 551 { 552 u_char arg_sign = getsign(st0_ptr); 553 554 if (tag == TAG_Valid) { 555 int q; 556 557 if (exponent(st0_ptr) > -40) { 558 if ((q = trig_arg(st0_ptr, 0)) == -1) { 559 /* Operand is out of range */ 560 return 1; 561 } 562 563 poly_sine(st0_ptr); 564 565 if (q & 2) 566 changesign(st0_ptr); 567 568 setsign(st0_ptr, getsign(st0_ptr) ^ arg_sign); 569 570 /* We do not really know if up or down */ 571 set_precision_flag_up(); 572 return 0; 573 } else { 574 /* For a small arg, the result == the argument */ 575 set_precision_flag_up(); /* Must be up. */ 576 return 0; 577 } 578 } 579 580 if (tag == TAG_Zero) { 581 setcc(0); 582 return 0; 583 } 584 585 if (tag == TAG_Special) 586 tag = FPU_Special(st0_ptr); 587 588 if (tag == TW_Denormal) { 589 if (denormal_operand() < 0) 590 return 1; 591 592 /* For a small arg, the result == the argument */ 593 /* Underflow may happen */ 594 FPU_to_exp16(st0_ptr, st0_ptr); 595 596 tag = FPU_round(st0_ptr, 1, 0, FULL_PRECISION, arg_sign); 597 598 FPU_settag0(tag); 599 600 return 0; 601 } else if (tag == TW_Infinity) { 602 /* The 80486 treats infinity as an invalid operand */ 603 arith_invalid(0); 604 return 1; 605 } else { 606 single_arg_error(st0_ptr, tag); 607 return 1; 608 } 609 } 610 611 static void fsin(FPU_REG *st0_ptr, u_char tag) 612 { 613 f_sin(st0_ptr, tag); 614 } 615 616 static int f_cos(FPU_REG *st0_ptr, u_char tag) 617 { 618 u_char st0_sign; 619 620 st0_sign = getsign(st0_ptr); 621 622 if (tag == TAG_Valid) { 623 int q; 624 625 if (exponent(st0_ptr) > -40) { 626 if ((exponent(st0_ptr) < 0) 627 || ((exponent(st0_ptr) == 0) 628 && (significand(st0_ptr) <= 629 0xc90fdaa22168c234LL))) { 630 poly_cos(st0_ptr); 631 632 /* We do not really know if up or down */ 633 set_precision_flag_down(); 634 635 return 0; 636 } else if ((q = trig_arg(st0_ptr, FCOS)) != -1) { 637 poly_sine(st0_ptr); 638 639 if ((q + 1) & 2) 640 changesign(st0_ptr); 641 642 /* We do not really know if up or down */ 643 set_precision_flag_down(); 644 645 return 0; 646 } else { 647 /* Operand is out of range */ 648 return 1; 649 } 650 } else { 651 denormal_arg: 652 653 setcc(0); 654 FPU_copy_to_reg0(&CONST_1, TAG_Valid); 655 #ifdef PECULIAR_486 656 set_precision_flag_down(); /* 80486 appears to do this. */ 657 #else 658 set_precision_flag_up(); /* Must be up. */ 659 #endif /* PECULIAR_486 */ 660 return 0; 661 } 662 } else if (tag == TAG_Zero) { 663 FPU_copy_to_reg0(&CONST_1, TAG_Valid); 664 setcc(0); 665 return 0; 666 } 667 668 if (tag == TAG_Special) 669 tag = FPU_Special(st0_ptr); 670 671 if (tag == TW_Denormal) { 672 if (denormal_operand() < 0) 673 return 1; 674 675 goto denormal_arg; 676 } else if (tag == TW_Infinity) { 677 /* The 80486 treats infinity as an invalid operand */ 678 arith_invalid(0); 679 return 1; 680 } else { 681 single_arg_error(st0_ptr, tag); /* requires st0_ptr == &st(0) */ 682 return 1; 683 } 684 } 685 686 static void fcos(FPU_REG *st0_ptr, u_char st0_tag) 687 { 688 f_cos(st0_ptr, st0_tag); 689 } 690 691 static void fsincos(FPU_REG *st0_ptr, u_char st0_tag) 692 { 693 FPU_REG *st_new_ptr; 694 FPU_REG arg; 695 u_char tag; 696 697 /* Stack underflow has higher priority */ 698 if (st0_tag == TAG_Empty) { 699 FPU_stack_underflow(); /* Puts a QNaN in st(0) */ 700 if (control_word & CW_Invalid) { 701 st_new_ptr = &st(-1); 702 push(); 703 FPU_stack_underflow(); /* Puts a QNaN in the new st(0) */ 704 } 705 return; 706 } 707 708 if (STACK_OVERFLOW) { 709 FPU_stack_overflow(); 710 return; 711 } 712 713 if (st0_tag == TAG_Special) 714 tag = FPU_Special(st0_ptr); 715 else 716 tag = st0_tag; 717 718 if (tag == TW_NaN) { 719 single_arg_2_error(st0_ptr, TW_NaN); 720 return; 721 } else if (tag == TW_Infinity) { 722 /* The 80486 treats infinity as an invalid operand */ 723 if (arith_invalid(0) >= 0) { 724 /* Masked response */ 725 push(); 726 arith_invalid(0); 727 } 728 return; 729 } 730 731 reg_copy(st0_ptr, &arg); 732 if (!f_sin(st0_ptr, st0_tag)) { 733 push(); 734 FPU_copy_to_reg0(&arg, st0_tag); 735 f_cos(&st(0), st0_tag); 736 } else { 737 /* An error, so restore st(0) */ 738 FPU_copy_to_reg0(&arg, st0_tag); 739 } 740 } 741 742 /*---------------------------------------------------------------------------*/ 743 /* The following all require two arguments: st(0) and st(1) */ 744 745 /* A lean, mean kernel for the fprem instructions. This relies upon 746 the division and rounding to an integer in do_fprem giving an 747 exact result. Because of this, rem_kernel() needs to deal only with 748 the least significant 64 bits, the more significant bits of the 749 result must be zero. 750 */ 751 static void rem_kernel(unsigned long long st0, unsigned long long *y, 752 unsigned long long st1, unsigned long long q, int n) 753 { 754 int dummy; 755 unsigned long long x; 756 757 x = st0 << n; 758 759 /* Do the required multiplication and subtraction in the one operation */ 760 761 /* lsw x -= lsw st1 * lsw q */ 762 asm volatile ("mull %4; subl %%eax,%0; sbbl %%edx,%1":"=m" 763 (((unsigned *)&x)[0]), "=m"(((unsigned *)&x)[1]), 764 "=a"(dummy) 765 :"2"(((unsigned *)&st1)[0]), "m"(((unsigned *)&q)[0]) 766 :"%dx"); 767 /* msw x -= msw st1 * lsw q */ 768 asm volatile ("mull %3; subl %%eax,%0":"=m" (((unsigned *)&x)[1]), 769 "=a"(dummy) 770 :"1"(((unsigned *)&st1)[1]), "m"(((unsigned *)&q)[0]) 771 :"%dx"); 772 /* msw x -= lsw st1 * msw q */ 773 asm volatile ("mull %3; subl %%eax,%0":"=m" (((unsigned *)&x)[1]), 774 "=a"(dummy) 775 :"1"(((unsigned *)&st1)[0]), "m"(((unsigned *)&q)[1]) 776 :"%dx"); 777 778 *y = x; 779 } 780 781 /* Remainder of st(0) / st(1) */ 782 /* This routine produces exact results, i.e. there is never any 783 rounding or truncation, etc of the result. */ 784 static void do_fprem(FPU_REG *st0_ptr, u_char st0_tag, int round) 785 { 786 FPU_REG *st1_ptr = &st(1); 787 u_char st1_tag = FPU_gettagi(1); 788 789 if (!((st0_tag ^ TAG_Valid) | (st1_tag ^ TAG_Valid))) { 790 FPU_REG tmp, st0, st1; 791 u_char st0_sign, st1_sign; 792 u_char tmptag; 793 int tag; 794 int old_cw; 795 int expdif; 796 long long q; 797 unsigned short saved_status; 798 int cc; 799 800 fprem_valid: 801 /* Convert registers for internal use. */ 802 st0_sign = FPU_to_exp16(st0_ptr, &st0); 803 st1_sign = FPU_to_exp16(st1_ptr, &st1); 804 expdif = exponent16(&st0) - exponent16(&st1); 805 806 old_cw = control_word; 807 cc = 0; 808 809 /* We want the status following the denorm tests, but don't want 810 the status changed by the arithmetic operations. */ 811 saved_status = partial_status; 812 control_word &= ~CW_RC; 813 control_word |= RC_CHOP; 814 815 if (expdif < 64) { 816 /* This should be the most common case */ 817 818 if (expdif > -2) { 819 u_char sign = st0_sign ^ st1_sign; 820 tag = FPU_u_div(&st0, &st1, &tmp, 821 PR_64_BITS | RC_CHOP | 0x3f, 822 sign); 823 setsign(&tmp, sign); 824 825 if (exponent(&tmp) >= 0) { 826 FPU_round_to_int(&tmp, tag); /* Fortunately, this can't 827 overflow to 2^64 */ 828 q = significand(&tmp); 829 830 rem_kernel(significand(&st0), 831 &significand(&tmp), 832 significand(&st1), 833 q, expdif); 834 835 setexponent16(&tmp, exponent16(&st1)); 836 } else { 837 reg_copy(&st0, &tmp); 838 q = 0; 839 } 840 841 if ((round == RC_RND) 842 && (tmp.sigh & 0xc0000000)) { 843 /* We may need to subtract st(1) once more, 844 to get a result <= 1/2 of st(1). */ 845 unsigned long long x; 846 expdif = 847 exponent16(&st1) - exponent16(&tmp); 848 if (expdif <= 1) { 849 if (expdif == 0) 850 x = significand(&st1) - 851 significand(&tmp); 852 else /* expdif is 1 */ 853 x = (significand(&st1) 854 << 1) - 855 significand(&tmp); 856 if ((x < significand(&tmp)) || 857 /* or equi-distant (from 0 & st(1)) and q is odd */ 858 ((x == significand(&tmp)) 859 && (q & 1))) { 860 st0_sign = !st0_sign; 861 significand(&tmp) = x; 862 q++; 863 } 864 } 865 } 866 867 if (q & 4) 868 cc |= SW_C0; 869 if (q & 2) 870 cc |= SW_C3; 871 if (q & 1) 872 cc |= SW_C1; 873 } else { 874 control_word = old_cw; 875 setcc(0); 876 return; 877 } 878 } else { 879 /* There is a large exponent difference ( >= 64 ) */ 880 /* To make much sense, the code in this section should 881 be done at high precision. */ 882 int exp_1, N; 883 u_char sign; 884 885 /* prevent overflow here */ 886 /* N is 'a number between 32 and 63' (p26-113) */ 887 reg_copy(&st0, &tmp); 888 tmptag = st0_tag; 889 N = (expdif & 0x0000001f) + 32; /* This choice gives results 890 identical to an AMD 486 */ 891 setexponent16(&tmp, N); 892 exp_1 = exponent16(&st1); 893 setexponent16(&st1, 0); 894 expdif -= N; 895 896 sign = getsign(&tmp) ^ st1_sign; 897 tag = 898 FPU_u_div(&tmp, &st1, &tmp, 899 PR_64_BITS | RC_CHOP | 0x3f, sign); 900 setsign(&tmp, sign); 901 902 FPU_round_to_int(&tmp, tag); /* Fortunately, this can't 903 overflow to 2^64 */ 904 905 rem_kernel(significand(&st0), 906 &significand(&tmp), 907 significand(&st1), 908 significand(&tmp), exponent(&tmp) 909 ); 910 setexponent16(&tmp, exp_1 + expdif); 911 912 /* It is possible for the operation to be complete here. 913 What does the IEEE standard say? The Intel 80486 manual 914 implies that the operation will never be completed at this 915 point, and the behaviour of a real 80486 confirms this. 916 */ 917 if (!(tmp.sigh | tmp.sigl)) { 918 /* The result is zero */ 919 control_word = old_cw; 920 partial_status = saved_status; 921 FPU_copy_to_reg0(&CONST_Z, TAG_Zero); 922 setsign(&st0, st0_sign); 923 #ifdef PECULIAR_486 924 setcc(SW_C2); 925 #else 926 setcc(0); 927 #endif /* PECULIAR_486 */ 928 return; 929 } 930 cc = SW_C2; 931 } 932 933 control_word = old_cw; 934 partial_status = saved_status; 935 tag = FPU_normalize_nuo(&tmp); 936 reg_copy(&tmp, st0_ptr); 937 938 /* The only condition to be looked for is underflow, 939 and it can occur here only if underflow is unmasked. */ 940 if ((exponent16(&tmp) <= EXP_UNDER) && (tag != TAG_Zero) 941 && !(control_word & CW_Underflow)) { 942 setcc(cc); 943 tag = arith_underflow(st0_ptr); 944 setsign(st0_ptr, st0_sign); 945 FPU_settag0(tag); 946 return; 947 } else if ((exponent16(&tmp) > EXP_UNDER) || (tag == TAG_Zero)) { 948 stdexp(st0_ptr); 949 setsign(st0_ptr, st0_sign); 950 } else { 951 tag = 952 FPU_round(st0_ptr, 0, 0, FULL_PRECISION, st0_sign); 953 } 954 FPU_settag0(tag); 955 setcc(cc); 956 957 return; 958 } 959 960 if (st0_tag == TAG_Special) 961 st0_tag = FPU_Special(st0_ptr); 962 if (st1_tag == TAG_Special) 963 st1_tag = FPU_Special(st1_ptr); 964 965 if (((st0_tag == TAG_Valid) && (st1_tag == TW_Denormal)) 966 || ((st0_tag == TW_Denormal) && (st1_tag == TAG_Valid)) 967 || ((st0_tag == TW_Denormal) && (st1_tag == TW_Denormal))) { 968 if (denormal_operand() < 0) 969 return; 970 goto fprem_valid; 971 } else if ((st0_tag == TAG_Empty) || (st1_tag == TAG_Empty)) { 972 FPU_stack_underflow(); 973 return; 974 } else if (st0_tag == TAG_Zero) { 975 if (st1_tag == TAG_Valid) { 976 setcc(0); 977 return; 978 } else if (st1_tag == TW_Denormal) { 979 if (denormal_operand() < 0) 980 return; 981 setcc(0); 982 return; 983 } else if (st1_tag == TAG_Zero) { 984 arith_invalid(0); 985 return; 986 } /* fprem(?,0) always invalid */ 987 else if (st1_tag == TW_Infinity) { 988 setcc(0); 989 return; 990 } 991 } else if ((st0_tag == TAG_Valid) || (st0_tag == TW_Denormal)) { 992 if (st1_tag == TAG_Zero) { 993 arith_invalid(0); /* fprem(Valid,Zero) is invalid */ 994 return; 995 } else if (st1_tag != TW_NaN) { 996 if (((st0_tag == TW_Denormal) 997 || (st1_tag == TW_Denormal)) 998 && (denormal_operand() < 0)) 999 return; 1000 1001 if (st1_tag == TW_Infinity) { 1002 /* fprem(Valid,Infinity) is o.k. */ 1003 setcc(0); 1004 return; 1005 } 1006 } 1007 } else if (st0_tag == TW_Infinity) { 1008 if (st1_tag != TW_NaN) { 1009 arith_invalid(0); /* fprem(Infinity,?) is invalid */ 1010 return; 1011 } 1012 } 1013 1014 /* One of the registers must contain a NaN if we got here. */ 1015 1016 #ifdef PARANOID 1017 if ((st0_tag != TW_NaN) && (st1_tag != TW_NaN)) 1018 EXCEPTION(EX_INTERNAL | 0x118); 1019 #endif /* PARANOID */ 1020 1021 real_2op_NaN(st1_ptr, st1_tag, 0, st1_ptr); 1022 1023 } 1024 1025 /* ST(1) <- ST(1) * log ST; pop ST */ 1026 static void fyl2x(FPU_REG *st0_ptr, u_char st0_tag) 1027 { 1028 FPU_REG *st1_ptr = &st(1), exponent; 1029 u_char st1_tag = FPU_gettagi(1); 1030 u_char sign; 1031 int e, tag; 1032 1033 clear_C1(); 1034 1035 if ((st0_tag == TAG_Valid) && (st1_tag == TAG_Valid)) { 1036 both_valid: 1037 /* Both regs are Valid or Denormal */ 1038 if (signpositive(st0_ptr)) { 1039 if (st0_tag == TW_Denormal) 1040 FPU_to_exp16(st0_ptr, st0_ptr); 1041 else 1042 /* Convert st(0) for internal use. */ 1043 setexponent16(st0_ptr, exponent(st0_ptr)); 1044 1045 if ((st0_ptr->sigh == 0x80000000) 1046 && (st0_ptr->sigl == 0)) { 1047 /* Special case. The result can be precise. */ 1048 u_char esign; 1049 e = exponent16(st0_ptr); 1050 if (e >= 0) { 1051 exponent.sigh = e; 1052 esign = SIGN_POS; 1053 } else { 1054 exponent.sigh = -e; 1055 esign = SIGN_NEG; 1056 } 1057 exponent.sigl = 0; 1058 setexponent16(&exponent, 31); 1059 tag = FPU_normalize_nuo(&exponent); 1060 stdexp(&exponent); 1061 setsign(&exponent, esign); 1062 tag = 1063 FPU_mul(&exponent, tag, 1, FULL_PRECISION); 1064 if (tag >= 0) 1065 FPU_settagi(1, tag); 1066 } else { 1067 /* The usual case */ 1068 sign = getsign(st1_ptr); 1069 if (st1_tag == TW_Denormal) 1070 FPU_to_exp16(st1_ptr, st1_ptr); 1071 else 1072 /* Convert st(1) for internal use. */ 1073 setexponent16(st1_ptr, 1074 exponent(st1_ptr)); 1075 poly_l2(st0_ptr, st1_ptr, sign); 1076 } 1077 } else { 1078 /* negative */ 1079 if (arith_invalid(1) < 0) 1080 return; 1081 } 1082 1083 FPU_pop(); 1084 1085 return; 1086 } 1087 1088 if (st0_tag == TAG_Special) 1089 st0_tag = FPU_Special(st0_ptr); 1090 if (st1_tag == TAG_Special) 1091 st1_tag = FPU_Special(st1_ptr); 1092 1093 if ((st0_tag == TAG_Empty) || (st1_tag == TAG_Empty)) { 1094 FPU_stack_underflow_pop(1); 1095 return; 1096 } else if ((st0_tag <= TW_Denormal) && (st1_tag <= TW_Denormal)) { 1097 if (st0_tag == TAG_Zero) { 1098 if (st1_tag == TAG_Zero) { 1099 /* Both args zero is invalid */ 1100 if (arith_invalid(1) < 0) 1101 return; 1102 } else { 1103 u_char sign; 1104 sign = getsign(st1_ptr) ^ SIGN_NEG; 1105 if (FPU_divide_by_zero(1, sign) < 0) 1106 return; 1107 1108 setsign(st1_ptr, sign); 1109 } 1110 } else if (st1_tag == TAG_Zero) { 1111 /* st(1) contains zero, st(0) valid <> 0 */ 1112 /* Zero is the valid answer */ 1113 sign = getsign(st1_ptr); 1114 1115 if (signnegative(st0_ptr)) { 1116 /* log(negative) */ 1117 if (arith_invalid(1) < 0) 1118 return; 1119 } else if ((st0_tag == TW_Denormal) 1120 && (denormal_operand() < 0)) 1121 return; 1122 else { 1123 if (exponent(st0_ptr) < 0) 1124 sign ^= SIGN_NEG; 1125 1126 FPU_copy_to_reg1(&CONST_Z, TAG_Zero); 1127 setsign(st1_ptr, sign); 1128 } 1129 } else { 1130 /* One or both operands are denormals. */ 1131 if (denormal_operand() < 0) 1132 return; 1133 goto both_valid; 1134 } 1135 } else if ((st0_tag == TW_NaN) || (st1_tag == TW_NaN)) { 1136 if (real_2op_NaN(st0_ptr, st0_tag, 1, st0_ptr) < 0) 1137 return; 1138 } 1139 /* One or both arg must be an infinity */ 1140 else if (st0_tag == TW_Infinity) { 1141 if ((signnegative(st0_ptr)) || (st1_tag == TAG_Zero)) { 1142 /* log(-infinity) or 0*log(infinity) */ 1143 if (arith_invalid(1) < 0) 1144 return; 1145 } else { 1146 u_char sign = getsign(st1_ptr); 1147 1148 if ((st1_tag == TW_Denormal) 1149 && (denormal_operand() < 0)) 1150 return; 1151 1152 FPU_copy_to_reg1(&CONST_INF, TAG_Special); 1153 setsign(st1_ptr, sign); 1154 } 1155 } 1156 /* st(1) must be infinity here */ 1157 else if (((st0_tag == TAG_Valid) || (st0_tag == TW_Denormal)) 1158 && (signpositive(st0_ptr))) { 1159 if (exponent(st0_ptr) >= 0) { 1160 if ((exponent(st0_ptr) == 0) && 1161 (st0_ptr->sigh == 0x80000000) && 1162 (st0_ptr->sigl == 0)) { 1163 /* st(0) holds 1.0 */ 1164 /* infinity*log(1) */ 1165 if (arith_invalid(1) < 0) 1166 return; 1167 } 1168 /* else st(0) is positive and > 1.0 */ 1169 } else { 1170 /* st(0) is positive and < 1.0 */ 1171 1172 if ((st0_tag == TW_Denormal) 1173 && (denormal_operand() < 0)) 1174 return; 1175 1176 changesign(st1_ptr); 1177 } 1178 } else { 1179 /* st(0) must be zero or negative */ 1180 if (st0_tag == TAG_Zero) { 1181 /* This should be invalid, but a real 80486 is happy with it. */ 1182 1183 #ifndef PECULIAR_486 1184 sign = getsign(st1_ptr); 1185 if (FPU_divide_by_zero(1, sign) < 0) 1186 return; 1187 #endif /* PECULIAR_486 */ 1188 1189 changesign(st1_ptr); 1190 } else if (arith_invalid(1) < 0) /* log(negative) */ 1191 return; 1192 } 1193 1194 FPU_pop(); 1195 } 1196 1197 static void fpatan(FPU_REG *st0_ptr, u_char st0_tag) 1198 { 1199 FPU_REG *st1_ptr = &st(1); 1200 u_char st1_tag = FPU_gettagi(1); 1201 int tag; 1202 1203 clear_C1(); 1204 if (!((st0_tag ^ TAG_Valid) | (st1_tag ^ TAG_Valid))) { 1205 valid_atan: 1206 1207 poly_atan(st0_ptr, st0_tag, st1_ptr, st1_tag); 1208 1209 FPU_pop(); 1210 1211 return; 1212 } 1213 1214 if (st0_tag == TAG_Special) 1215 st0_tag = FPU_Special(st0_ptr); 1216 if (st1_tag == TAG_Special) 1217 st1_tag = FPU_Special(st1_ptr); 1218 1219 if (((st0_tag == TAG_Valid) && (st1_tag == TW_Denormal)) 1220 || ((st0_tag == TW_Denormal) && (st1_tag == TAG_Valid)) 1221 || ((st0_tag == TW_Denormal) && (st1_tag == TW_Denormal))) { 1222 if (denormal_operand() < 0) 1223 return; 1224 1225 goto valid_atan; 1226 } else if ((st0_tag == TAG_Empty) || (st1_tag == TAG_Empty)) { 1227 FPU_stack_underflow_pop(1); 1228 return; 1229 } else if ((st0_tag == TW_NaN) || (st1_tag == TW_NaN)) { 1230 if (real_2op_NaN(st0_ptr, st0_tag, 1, st0_ptr) >= 0) 1231 FPU_pop(); 1232 return; 1233 } else if ((st0_tag == TW_Infinity) || (st1_tag == TW_Infinity)) { 1234 u_char sign = getsign(st1_ptr); 1235 if (st0_tag == TW_Infinity) { 1236 if (st1_tag == TW_Infinity) { 1237 if (signpositive(st0_ptr)) { 1238 FPU_copy_to_reg1(&CONST_PI4, TAG_Valid); 1239 } else { 1240 setpositive(st1_ptr); 1241 tag = 1242 FPU_u_add(&CONST_PI4, &CONST_PI2, 1243 st1_ptr, FULL_PRECISION, 1244 SIGN_POS, 1245 exponent(&CONST_PI4), 1246 exponent(&CONST_PI2)); 1247 if (tag >= 0) 1248 FPU_settagi(1, tag); 1249 } 1250 } else { 1251 if ((st1_tag == TW_Denormal) 1252 && (denormal_operand() < 0)) 1253 return; 1254 1255 if (signpositive(st0_ptr)) { 1256 FPU_copy_to_reg1(&CONST_Z, TAG_Zero); 1257 setsign(st1_ptr, sign); /* An 80486 preserves the sign */ 1258 FPU_pop(); 1259 return; 1260 } else { 1261 FPU_copy_to_reg1(&CONST_PI, TAG_Valid); 1262 } 1263 } 1264 } else { 1265 /* st(1) is infinity, st(0) not infinity */ 1266 if ((st0_tag == TW_Denormal) 1267 && (denormal_operand() < 0)) 1268 return; 1269 1270 FPU_copy_to_reg1(&CONST_PI2, TAG_Valid); 1271 } 1272 setsign(st1_ptr, sign); 1273 } else if (st1_tag == TAG_Zero) { 1274 /* st(0) must be valid or zero */ 1275 u_char sign = getsign(st1_ptr); 1276 1277 if ((st0_tag == TW_Denormal) && (denormal_operand() < 0)) 1278 return; 1279 1280 if (signpositive(st0_ptr)) { 1281 /* An 80486 preserves the sign */ 1282 FPU_pop(); 1283 return; 1284 } 1285 1286 FPU_copy_to_reg1(&CONST_PI, TAG_Valid); 1287 setsign(st1_ptr, sign); 1288 } else if (st0_tag == TAG_Zero) { 1289 /* st(1) must be TAG_Valid here */ 1290 u_char sign = getsign(st1_ptr); 1291 1292 if ((st1_tag == TW_Denormal) && (denormal_operand() < 0)) 1293 return; 1294 1295 FPU_copy_to_reg1(&CONST_PI2, TAG_Valid); 1296 setsign(st1_ptr, sign); 1297 } 1298 #ifdef PARANOID 1299 else 1300 EXCEPTION(EX_INTERNAL | 0x125); 1301 #endif /* PARANOID */ 1302 1303 FPU_pop(); 1304 set_precision_flag_up(); /* We do not really know if up or down */ 1305 } 1306 1307 static void fprem(FPU_REG *st0_ptr, u_char st0_tag) 1308 { 1309 do_fprem(st0_ptr, st0_tag, RC_CHOP); 1310 } 1311 1312 static void fprem1(FPU_REG *st0_ptr, u_char st0_tag) 1313 { 1314 do_fprem(st0_ptr, st0_tag, RC_RND); 1315 } 1316 1317 static void fyl2xp1(FPU_REG *st0_ptr, u_char st0_tag) 1318 { 1319 u_char sign, sign1; 1320 FPU_REG *st1_ptr = &st(1), a, b; 1321 u_char st1_tag = FPU_gettagi(1); 1322 1323 clear_C1(); 1324 if (!((st0_tag ^ TAG_Valid) | (st1_tag ^ TAG_Valid))) { 1325 valid_yl2xp1: 1326 1327 sign = getsign(st0_ptr); 1328 sign1 = getsign(st1_ptr); 1329 1330 FPU_to_exp16(st0_ptr, &a); 1331 FPU_to_exp16(st1_ptr, &b); 1332 1333 if (poly_l2p1(sign, sign1, &a, &b, st1_ptr)) 1334 return; 1335 1336 FPU_pop(); 1337 return; 1338 } 1339 1340 if (st0_tag == TAG_Special) 1341 st0_tag = FPU_Special(st0_ptr); 1342 if (st1_tag == TAG_Special) 1343 st1_tag = FPU_Special(st1_ptr); 1344 1345 if (((st0_tag == TAG_Valid) && (st1_tag == TW_Denormal)) 1346 || ((st0_tag == TW_Denormal) && (st1_tag == TAG_Valid)) 1347 || ((st0_tag == TW_Denormal) && (st1_tag == TW_Denormal))) { 1348 if (denormal_operand() < 0) 1349 return; 1350 1351 goto valid_yl2xp1; 1352 } else if ((st0_tag == TAG_Empty) | (st1_tag == TAG_Empty)) { 1353 FPU_stack_underflow_pop(1); 1354 return; 1355 } else if (st0_tag == TAG_Zero) { 1356 switch (st1_tag) { 1357 case TW_Denormal: 1358 if (denormal_operand() < 0) 1359 return; 1360 fallthrough; 1361 case TAG_Zero: 1362 case TAG_Valid: 1363 setsign(st0_ptr, getsign(st0_ptr) ^ getsign(st1_ptr)); 1364 FPU_copy_to_reg1(st0_ptr, st0_tag); 1365 break; 1366 1367 case TW_Infinity: 1368 /* Infinity*log(1) */ 1369 if (arith_invalid(1) < 0) 1370 return; 1371 break; 1372 1373 case TW_NaN: 1374 if (real_2op_NaN(st0_ptr, st0_tag, 1, st0_ptr) < 0) 1375 return; 1376 break; 1377 1378 default: 1379 #ifdef PARANOID 1380 EXCEPTION(EX_INTERNAL | 0x116); 1381 return; 1382 #endif /* PARANOID */ 1383 break; 1384 } 1385 } else if ((st0_tag == TAG_Valid) || (st0_tag == TW_Denormal)) { 1386 switch (st1_tag) { 1387 case TAG_Zero: 1388 if (signnegative(st0_ptr)) { 1389 if (exponent(st0_ptr) >= 0) { 1390 /* st(0) holds <= -1.0 */ 1391 #ifdef PECULIAR_486 /* Stupid 80486 doesn't worry about log(negative). */ 1392 changesign(st1_ptr); 1393 #else 1394 if (arith_invalid(1) < 0) 1395 return; 1396 #endif /* PECULIAR_486 */ 1397 } else if ((st0_tag == TW_Denormal) 1398 && (denormal_operand() < 0)) 1399 return; 1400 else 1401 changesign(st1_ptr); 1402 } else if ((st0_tag == TW_Denormal) 1403 && (denormal_operand() < 0)) 1404 return; 1405 break; 1406 1407 case TW_Infinity: 1408 if (signnegative(st0_ptr)) { 1409 if ((exponent(st0_ptr) >= 0) && 1410 !((st0_ptr->sigh == 0x80000000) && 1411 (st0_ptr->sigl == 0))) { 1412 /* st(0) holds < -1.0 */ 1413 #ifdef PECULIAR_486 /* Stupid 80486 doesn't worry about log(negative). */ 1414 changesign(st1_ptr); 1415 #else 1416 if (arith_invalid(1) < 0) 1417 return; 1418 #endif /* PECULIAR_486 */ 1419 } else if ((st0_tag == TW_Denormal) 1420 && (denormal_operand() < 0)) 1421 return; 1422 else 1423 changesign(st1_ptr); 1424 } else if ((st0_tag == TW_Denormal) 1425 && (denormal_operand() < 0)) 1426 return; 1427 break; 1428 1429 case TW_NaN: 1430 if (real_2op_NaN(st0_ptr, st0_tag, 1, st0_ptr) < 0) 1431 return; 1432 } 1433 1434 } else if (st0_tag == TW_NaN) { 1435 if (real_2op_NaN(st0_ptr, st0_tag, 1, st0_ptr) < 0) 1436 return; 1437 } else if (st0_tag == TW_Infinity) { 1438 if (st1_tag == TW_NaN) { 1439 if (real_2op_NaN(st0_ptr, st0_tag, 1, st0_ptr) < 0) 1440 return; 1441 } else if (signnegative(st0_ptr)) { 1442 #ifndef PECULIAR_486 1443 /* This should have higher priority than denormals, but... */ 1444 if (arith_invalid(1) < 0) /* log(-infinity) */ 1445 return; 1446 #endif /* PECULIAR_486 */ 1447 if ((st1_tag == TW_Denormal) 1448 && (denormal_operand() < 0)) 1449 return; 1450 #ifdef PECULIAR_486 1451 /* Denormal operands actually get higher priority */ 1452 if (arith_invalid(1) < 0) /* log(-infinity) */ 1453 return; 1454 #endif /* PECULIAR_486 */ 1455 } else if (st1_tag == TAG_Zero) { 1456 /* log(infinity) */ 1457 if (arith_invalid(1) < 0) 1458 return; 1459 } 1460 1461 /* st(1) must be valid here. */ 1462 1463 else if ((st1_tag == TW_Denormal) && (denormal_operand() < 0)) 1464 return; 1465 1466 /* The Manual says that log(Infinity) is invalid, but a real 1467 80486 sensibly says that it is o.k. */ 1468 else { 1469 u_char sign = getsign(st1_ptr); 1470 FPU_copy_to_reg1(&CONST_INF, TAG_Special); 1471 setsign(st1_ptr, sign); 1472 } 1473 } 1474 #ifdef PARANOID 1475 else { 1476 EXCEPTION(EX_INTERNAL | 0x117); 1477 return; 1478 } 1479 #endif /* PARANOID */ 1480 1481 FPU_pop(); 1482 return; 1483 1484 } 1485 1486 static void fscale(FPU_REG *st0_ptr, u_char st0_tag) 1487 { 1488 FPU_REG *st1_ptr = &st(1); 1489 u_char st1_tag = FPU_gettagi(1); 1490 int old_cw = control_word; 1491 u_char sign = getsign(st0_ptr); 1492 1493 clear_C1(); 1494 if (!((st0_tag ^ TAG_Valid) | (st1_tag ^ TAG_Valid))) { 1495 long scale; 1496 FPU_REG tmp; 1497 1498 /* Convert register for internal use. */ 1499 setexponent16(st0_ptr, exponent(st0_ptr)); 1500 1501 valid_scale: 1502 1503 if (exponent(st1_ptr) > 30) { 1504 /* 2^31 is far too large, would require 2^(2^30) or 2^(-2^30) */ 1505 1506 if (signpositive(st1_ptr)) { 1507 EXCEPTION(EX_Overflow); 1508 FPU_copy_to_reg0(&CONST_INF, TAG_Special); 1509 } else { 1510 EXCEPTION(EX_Underflow); 1511 FPU_copy_to_reg0(&CONST_Z, TAG_Zero); 1512 } 1513 setsign(st0_ptr, sign); 1514 return; 1515 } 1516 1517 control_word &= ~CW_RC; 1518 control_word |= RC_CHOP; 1519 reg_copy(st1_ptr, &tmp); 1520 FPU_round_to_int(&tmp, st1_tag); /* This can never overflow here */ 1521 control_word = old_cw; 1522 scale = signnegative(st1_ptr) ? -tmp.sigl : tmp.sigl; 1523 scale += exponent16(st0_ptr); 1524 1525 setexponent16(st0_ptr, scale); 1526 1527 /* Use FPU_round() to properly detect under/overflow etc */ 1528 FPU_round(st0_ptr, 0, 0, control_word, sign); 1529 1530 return; 1531 } 1532 1533 if (st0_tag == TAG_Special) 1534 st0_tag = FPU_Special(st0_ptr); 1535 if (st1_tag == TAG_Special) 1536 st1_tag = FPU_Special(st1_ptr); 1537 1538 if ((st0_tag == TAG_Valid) || (st0_tag == TW_Denormal)) { 1539 switch (st1_tag) { 1540 case TAG_Valid: 1541 /* st(0) must be a denormal */ 1542 if ((st0_tag == TW_Denormal) 1543 && (denormal_operand() < 0)) 1544 return; 1545 1546 FPU_to_exp16(st0_ptr, st0_ptr); /* Will not be left on stack */ 1547 goto valid_scale; 1548 1549 case TAG_Zero: 1550 if (st0_tag == TW_Denormal) 1551 denormal_operand(); 1552 return; 1553 1554 case TW_Denormal: 1555 denormal_operand(); 1556 return; 1557 1558 case TW_Infinity: 1559 if ((st0_tag == TW_Denormal) 1560 && (denormal_operand() < 0)) 1561 return; 1562 1563 if (signpositive(st1_ptr)) 1564 FPU_copy_to_reg0(&CONST_INF, TAG_Special); 1565 else 1566 FPU_copy_to_reg0(&CONST_Z, TAG_Zero); 1567 setsign(st0_ptr, sign); 1568 return; 1569 1570 case TW_NaN: 1571 real_2op_NaN(st1_ptr, st1_tag, 0, st0_ptr); 1572 return; 1573 } 1574 } else if (st0_tag == TAG_Zero) { 1575 switch (st1_tag) { 1576 case TAG_Valid: 1577 case TAG_Zero: 1578 return; 1579 1580 case TW_Denormal: 1581 denormal_operand(); 1582 return; 1583 1584 case TW_Infinity: 1585 if (signpositive(st1_ptr)) 1586 arith_invalid(0); /* Zero scaled by +Infinity */ 1587 return; 1588 1589 case TW_NaN: 1590 real_2op_NaN(st1_ptr, st1_tag, 0, st0_ptr); 1591 return; 1592 } 1593 } else if (st0_tag == TW_Infinity) { 1594 switch (st1_tag) { 1595 case TAG_Valid: 1596 case TAG_Zero: 1597 return; 1598 1599 case TW_Denormal: 1600 denormal_operand(); 1601 return; 1602 1603 case TW_Infinity: 1604 if (signnegative(st1_ptr)) 1605 arith_invalid(0); /* Infinity scaled by -Infinity */ 1606 return; 1607 1608 case TW_NaN: 1609 real_2op_NaN(st1_ptr, st1_tag, 0, st0_ptr); 1610 return; 1611 } 1612 } else if (st0_tag == TW_NaN) { 1613 if (st1_tag != TAG_Empty) { 1614 real_2op_NaN(st1_ptr, st1_tag, 0, st0_ptr); 1615 return; 1616 } 1617 } 1618 #ifdef PARANOID 1619 if (!((st0_tag == TAG_Empty) || (st1_tag == TAG_Empty))) { 1620 EXCEPTION(EX_INTERNAL | 0x115); 1621 return; 1622 } 1623 #endif 1624 1625 /* At least one of st(0), st(1) must be empty */ 1626 FPU_stack_underflow(); 1627 1628 } 1629 1630 /*---------------------------------------------------------------------------*/ 1631 1632 static FUNC_ST0 const trig_table_a[] = { 1633 f2xm1, fyl2x, fptan, fpatan, 1634 fxtract, fprem1, (FUNC_ST0) fdecstp, (FUNC_ST0) fincstp 1635 }; 1636 1637 void FPU_triga(void) 1638 { 1639 (trig_table_a[FPU_rm]) (&st(0), FPU_gettag0()); 1640 } 1641 1642 static FUNC_ST0 const trig_table_b[] = { 1643 fprem, fyl2xp1, fsqrt_, fsincos, frndint_, fscale, fsin, fcos 1644 }; 1645 1646 void FPU_trigb(void) 1647 { 1648 (trig_table_b[FPU_rm]) (&st(0), FPU_gettag0()); 1649 } 1650