root/kernel/FPU-emu/fpu_trig.c

/* [previous][next][first][last][top][bottom][index][help] */

DEFINITIONS

This source file includes following definitions.
  1. trig_arg
  2. convert_l2reg
  3. single_arg_error
  4. f2xm1
  5. fptan
  6. fxtract
  7. fdecstp
  8. fincstp
  9. fsqrt_
  10. frndint_
  11. fsin
  12. f_cos
  13. fcos
  14. fsincos
  15. fprem_kernel
  16. fyl2x
  17. fpatan
  18. fprem
  19. fprem1
  20. fyl2xp1
  21. fscale
  22. trig_a
  23. trig_b

   1 /*---------------------------------------------------------------------------+
   2  |  fpu_trig.c                                                               |
   3  |                                                                           |
   4  | Implementation of the FPU "transcendental" functions.                     |
   5  |                                                                           |
   6  | Copyright (C) 1992    W. Metzenthen, 22 Parker St, Ormond, Vic 3163,      |
   7  |                       Australia.  E-mail apm233m@vaxc.cc.monash.edu.au    |
   8  |                                                                           |
   9  |                                                                           |
  10  +---------------------------------------------------------------------------*/
  11 
  12 #include "fpu_system.h"
  13 #include "exception.h"
  14 #include "fpu_emu.h"
  15 #include "status_w.h"
  16 #include "control_w.h"
  17 #include "reg_constant.h"       
  18 
  19 
  20 
  21 static int trig_arg(FPU_REG *X)
     /* [previous][next][first][last][top][bottom][index][help] */
  22 {
  23   FPU_REG tmp, quot;
  24   int rv;
  25   long long q;
  26   int old_cw = control_word;
  27 
  28   control_word &= ~CW_RC;
  29   control_word |= RC_CHOP;
  30   
  31   reg_move(X, &quot);
  32   reg_div(&quot, &CONST_PI2, &quot);
  33 
  34   reg_move(&quot, &tmp);
  35   round_to_int(&tmp);
  36   if ( tmp.sigh & 0x80000000 )
  37     return -1;              /* |Arg| is >= 2^63 */
  38   tmp.exp = EXP_BIAS + 63;
  39   q = *(long long *)&(tmp.sigl);
  40   normalize(&tmp);
  41 
  42   reg_sub(&quot, &tmp, X);
  43   rv = q & 7;
  44   
  45   control_word = old_cw;
  46   return rv;;
  47 }
  48 
  49 
  50 /* Convert a long to register */
  51 void convert_l2reg(long *arg, FPU_REG *dest)
     /* [previous][next][first][last][top][bottom][index][help] */
  52 {
  53   long num = *arg;
  54   
  55   if (num == 0)
  56     { reg_move(&CONST_Z, dest); return; }
  57 
  58   if (num > 0)
  59     dest->sign = SIGN_POS;
  60   else
  61     { num = -num; dest->sign = SIGN_NEG; }
  62   
  63   dest->sigh = num;
  64   dest->sigl = 0;
  65   dest->exp = EXP_BIAS + 31;
  66   dest->tag = TW_Valid;
  67   normalize(dest);
  68 }
  69 
  70 
  71 static void single_arg_error(void)
     /* [previous][next][first][last][top][bottom][index][help] */
  72 {
  73   switch ( FPU_st0_tag )
  74     {
  75     case TW_NaN:
  76       if ( !(FPU_st0_ptr->sigh & 0x40000000) )   /* Signaling ? */
  77         {
  78           EXCEPTION(EX_Invalid);
  79           /* Convert to a QNaN */
  80           FPU_st0_ptr->sigh |= 0x40000000;
  81         }
  82       break;              /* return with a NaN in st(0) */
  83     case TW_Empty:
  84       stack_underflow();  /* Puts a QNaN in st(0) */
  85       break;
  86 #ifdef PARANOID
  87     default:
  88       EXCEPTION(EX_INTERNAL|0x0112);
  89 #endif PARANOID
  90     }
  91 }
  92 
  93 
  94 /*---------------------------------------------------------------------------*/
  95 
  96 static void f2xm1(void)
     /* [previous][next][first][last][top][bottom][index][help] */
  97 {
  98   switch ( FPU_st0_tag )
  99     {
 100     case TW_Valid:
 101       {
 102         FPU_REG rv, tmp;
 103         
 104         if ( FPU_st0_ptr->sign == SIGN_POS )
 105           {
 106             /* poly_2xm1(x) requires 0 < x < 1. */
 107             if ( poly_2xm1(FPU_st0_ptr, &rv) )
 108               return;
 109             reg_mul(&rv, FPU_st0_ptr, FPU_st0_ptr);
 110             return;
 111           }
 112         else
 113           {
 114 /* **** Should change poly_2xm1() to at least handle numbers near 0 */
 115             /* poly_2xm1(x) doesn't handle negative numbers. */
 116             /* So we compute (poly_2xm1(x+1)-1)/2, for -1 < x < 0 */
 117             reg_add(FPU_st0_ptr, &CONST_1, &tmp);
 118             poly_2xm1(&tmp, &rv);
 119             reg_mul(&rv, &tmp, &tmp);
 120             reg_sub(&tmp, &CONST_1, FPU_st0_ptr);
 121             FPU_st0_ptr->exp--;
 122           }
 123         if ( FPU_st0_ptr->exp <= EXP_UNDER )
 124           arith_underflow(FPU_st0_ptr);
 125         return;
 126       }
 127     case TW_Zero:
 128       return;
 129     case TW_Infinity:
 130       if ( FPU_st0_ptr->sign == SIGN_NEG )
 131         {
 132           /* -infinity gives -1 (p16-10) */
 133           reg_move(&CONST_1, FPU_st0_ptr);
 134           FPU_st0_ptr->sign = SIGN_NEG;
 135         }
 136       return;
 137     default:
 138       single_arg_error();
 139     }
 140 }
 141 
 142 static void fptan(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 143 {
 144   FPU_REG *st_new_ptr;
 145   int q;
 146   char arg_sign = FPU_st0_ptr->sign;
 147 
 148   if ( STACK_OVERFLOW )
 149     { stack_overflow(); return; }
 150 
 151   switch ( FPU_st0_tag )
 152     {
 153     case TW_Valid:
 154       FPU_st0_ptr->sign = SIGN_POS;
 155       if ( (q = trig_arg(FPU_st0_ptr)) != -1 )
 156         {
 157           if (q & 1)
 158             reg_sub(&CONST_1, FPU_st0_ptr, FPU_st0_ptr);
 159           
 160           poly_tan(FPU_st0_ptr, FPU_st0_ptr);
 161 
 162           FPU_st0_ptr->sign = (q & 1) ^ arg_sign;
 163           
 164           if ( FPU_st0_ptr->exp <= EXP_UNDER )
 165             arith_underflow(FPU_st0_ptr);
 166           
 167           push();
 168           reg_move(&CONST_1, FPU_st0_ptr);
 169           setcc(0);
 170         }
 171       else
 172         {
 173           /* Operand is out of range */
 174           setcc(SW_C2);
 175           FPU_st0_ptr->sign = arg_sign;         /* restore st(0) */
 176           return;
 177         }
 178       break;
 179     case TW_Infinity:
 180       arith_invalid(FPU_st0_ptr);
 181       setcc(0);
 182       return;
 183     case TW_Zero:
 184       push();
 185       reg_move(&CONST_1, FPU_st0_ptr);
 186       setcc(0);
 187       break;
 188     default:
 189       single_arg_error();
 190       break;
 191     }
 192 }
 193 
 194 
 195 static void fxtract(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 196 {
 197   FPU_REG *st_new_ptr;
 198   register FPU_REG *st1_ptr = FPU_st0_ptr;  /* anticipate */
 199 
 200   if ( STACK_OVERFLOW )
 201     {  stack_overflow(); return; }
 202 
 203   if ( !(FPU_st0_tag ^ TW_Valid) )
 204     {
 205       long e;
 206           
 207       push();
 208       reg_move(st1_ptr, FPU_st0_ptr);
 209       FPU_st0_ptr->exp = EXP_BIAS;
 210       e = st1_ptr->exp - EXP_BIAS;
 211       convert_l2reg(&e, st1_ptr);
 212       return;
 213     }
 214   else if ( FPU_st0_tag == TW_Zero )
 215     {
 216       char sign = FPU_st0_ptr->sign;
 217       divide_by_zero(SIGN_NEG, FPU_st0_ptr);
 218       push();
 219       reg_move(&CONST_Z, FPU_st0_ptr);
 220       FPU_st0_ptr->sign = sign;
 221       return;
 222     }
 223   else if ( FPU_st0_tag == TW_Infinity )
 224     {
 225       char sign = FPU_st0_ptr->sign;
 226       FPU_st0_ptr->sign = SIGN_POS;
 227       push();
 228       reg_move(&CONST_INF, FPU_st0_ptr);
 229       FPU_st0_ptr->sign = sign;
 230       return;
 231     }
 232   else if ( FPU_st0_tag == TW_NaN )
 233     {
 234       if ( !(FPU_st0_ptr->sigh & 0x40000000) )   /* Signaling ? */
 235         {
 236           EXCEPTION(EX_Invalid);
 237           /* Convert to a QNaN */
 238           FPU_st0_ptr->sigh |= 0x40000000;
 239         }
 240       push();
 241       reg_move(st1_ptr, FPU_st0_ptr);
 242       return;
 243     }
 244   else if ( FPU_st0_tag == TW_Empty )
 245     {
 246       /* Is this the correct behaviour? */
 247       if ( control_word & EX_Invalid )
 248         {
 249           stack_underflow();
 250           push();
 251           stack_underflow();
 252         }
 253       else
 254         EXCEPTION(EX_StackUnder);
 255     }
 256 #ifdef PARANOID
 257   else
 258     EXCEPTION(EX_INTERNAL | 0x119);
 259 #endif PARANOID
 260 }
 261 
 262 
 263 static void fdecstp(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 264 {
 265   top--;  /* FPU_st0_ptr will be fixed in math_emulate() before the next instr */
 266 }
 267 
 268 static void fincstp(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 269 {
 270   top++;  /* FPU_st0_ptr will be fixed in math_emulate() before the next instr */
 271 }
 272 
 273 
 274 static void fsqrt_(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 275 {
 276   if ( !(FPU_st0_tag ^ TW_Valid) )
 277     {
 278       int expon;
 279       
 280       if (FPU_st0_ptr->sign == SIGN_NEG)
 281         {
 282           arith_invalid(FPU_st0_ptr);
 283           return;
 284         }
 285 
 286       expon = FPU_st0_ptr->exp - EXP_BIAS;
 287       FPU_st0_ptr->exp = EXP_BIAS + (expon & 1);  /* make st(0) in  [1.0 .. 4.0) */
 288       
 289       wm_sqrt(FPU_st0_ptr);     /* Do the computation */
 290       
 291       FPU_st0_ptr->exp += expon >> 1;
 292       FPU_st0_ptr->tag = TW_Valid;
 293       FPU_st0_ptr->sign = SIGN_POS;
 294     }
 295   else if ( FPU_st0_tag == TW_Zero )
 296     return;
 297   else if ( FPU_st0_tag == TW_Infinity )
 298     {
 299       if ( FPU_st0_ptr->sign == SIGN_NEG )
 300         arith_invalid(FPU_st0_ptr);
 301       return;
 302     }
 303   else
 304     single_arg_error();
 305 }
 306 
 307 
 308 static void frndint_(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 309 {
 310   if ( !(FPU_st0_tag ^ TW_Valid) )
 311     {
 312       if (FPU_st0_ptr->exp > EXP_BIAS+63)
 313         return;
 314 
 315       round_to_int(FPU_st0_ptr);  /* Fortunately, this can't overflow to 2^64 */
 316       FPU_st0_ptr->exp = EXP_BIAS + 63;
 317       normalize(FPU_st0_ptr);
 318       return;
 319     }
 320   else if ( (FPU_st0_tag == TW_Zero) || (FPU_st0_tag == TW_Infinity) )
 321     return;
 322   else
 323     single_arg_error();
 324 }
 325 
 326 
 327 static void fsin(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 328 {
 329   if ( FPU_st0_tag == TW_Valid )
 330     {
 331       int q;
 332       char arg_sign = FPU_st0_ptr->sign;
 333       FPU_st0_ptr->sign = SIGN_POS;
 334       if ( (q = trig_arg(FPU_st0_ptr)) != -1 )
 335         {
 336           FPU_REG rv;
 337           
 338           if (q & 1)
 339             reg_sub(&CONST_1, FPU_st0_ptr, FPU_st0_ptr);
 340           
 341           poly_sine(FPU_st0_ptr, &rv);
 342           
 343           setcc(0);
 344           if (q & 2)
 345             rv.sign ^= SIGN_POS ^ SIGN_NEG;
 346           rv.sign ^= arg_sign;
 347           reg_move(&rv, FPU_st0_ptr);
 348 
 349           if ( FPU_st0_ptr->exp <= EXP_UNDER )
 350             arith_underflow(FPU_st0_ptr);
 351 
 352           return;
 353         }
 354       else
 355         {
 356           /* Operand is out of range */
 357           setcc(SW_C2);
 358           FPU_st0_ptr->sign = arg_sign;         /* restore st(0) */
 359           EXCEPTION(EX_Invalid);
 360           return;
 361         }
 362     }
 363   else if ( FPU_st0_tag == TW_Zero )
 364     {
 365       setcc(0);
 366       return;
 367     }
 368   else if ( FPU_st0_tag == TW_Infinity )
 369     {
 370       arith_invalid(FPU_st0_ptr);
 371       setcc(0);
 372       return;
 373     }
 374   else
 375     single_arg_error();
 376 }
 377 
 378 
 379 static int f_cos(FPU_REG *arg)
     /* [previous][next][first][last][top][bottom][index][help] */
 380 {
 381   if ( arg->tag == TW_Valid )
 382     {
 383       int q;
 384       char arg_sign = arg->sign;
 385       arg->sign = SIGN_POS;
 386       if ( (q = trig_arg(arg)) != -1 )
 387         {
 388           FPU_REG rv;
 389           
 390           if ( !(q & 1) )
 391             reg_sub(&CONST_1, arg, arg);
 392           
 393           poly_sine(arg, &rv);
 394           
 395           setcc(0);
 396           if ((q+1) & 2)
 397             rv.sign ^= SIGN_POS ^ SIGN_NEG;
 398           reg_move(&rv, arg);
 399           
 400           return 0;
 401         }
 402       else
 403         {
 404           /* Operand is out of range */
 405           setcc(SW_C2);
 406           arg->sign = arg_sign;         /* restore st(0) */
 407           EXCEPTION(EX_Invalid);
 408           return 1;
 409         }
 410     }
 411   else if ( arg->tag == TW_Zero )
 412     {
 413       reg_move(&CONST_1, arg);
 414       setcc(0);
 415       return 0;
 416     }
 417   else if ( FPU_st0_tag == TW_Infinity )
 418     {
 419       arith_invalid(FPU_st0_ptr);
 420       setcc(0);
 421       return 1;
 422     }
 423   else
 424     {
 425       single_arg_error();  /* requires arg == &st(0) */
 426       return 1;
 427     }
 428 }
 429 
 430 
 431 static void fcos(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 432 {
 433   f_cos(FPU_st0_ptr);
 434 }
 435 
 436 
 437 static void fsincos(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 438 {
 439   FPU_REG *st_new_ptr;
 440   FPU_REG arg;
 441           
 442   if ( STACK_OVERFLOW )
 443     { stack_overflow(); return; }
 444 
 445   reg_move(FPU_st0_ptr,&arg);
 446   if ( !f_cos(&arg) )
 447     {
 448       fsin();
 449       push();
 450       reg_move(&arg,FPU_st0_ptr);
 451     }
 452 
 453 }
 454 
 455 
 456 /*---------------------------------------------------------------------------*/
 457 /* The following all require two arguments: st(0) and st(1) */
 458 
 459 /* remainder of st(0) / st(1) */
 460 /* Assumes that st(0) and st(1) are both TW_Valid */
 461 static void fprem_kernel(int round)
     /* [previous][next][first][last][top][bottom][index][help] */
 462 {
 463   FPU_REG *st1_ptr = &st(1);
 464   char st1_tag = st1_ptr->tag;
 465 
 466   if ( !((FPU_st0_tag ^ TW_Valid) | (st1_tag ^ TW_Valid)) )
 467     {
 468       FPU_REG tmp;
 469       int old_cw = control_word;
 470       int expdif = FPU_st0_ptr->exp - (st1_ptr)->exp;
 471       
 472       control_word &= ~CW_RC;
 473       control_word |= round;
 474       
 475       if (expdif < 64)
 476         {
 477           /* This should be the most common case */
 478           long long q;
 479           int c = 0;
 480           reg_div(FPU_st0_ptr, st1_ptr, &tmp);
 481           
 482           round_to_int(&tmp);  /* Fortunately, this can't overflow to 2^64 */
 483           tmp.exp = EXP_BIAS + 63;
 484           q = *(long long *)&(tmp.sigl);
 485           normalize(&tmp);
 486           
 487           reg_mul(st1_ptr, &tmp, &tmp);
 488           reg_sub(FPU_st0_ptr, &tmp, FPU_st0_ptr);
 489           
 490           if (q&4) c |= SW_C3;
 491           if (q&2) c |= SW_C1;
 492           if (q&1) c |= SW_C0;
 493           
 494           setcc(c);
 495         }
 496       else
 497         {
 498           /* There is a large exponent difference ( >= 64 ) */
 499           int N_exp;
 500           
 501           reg_div(FPU_st0_ptr, st1_ptr, &tmp);
 502           /* N is 'a number between 32 and 63' (p26-113) */
 503           N_exp = (tmp.exp & 31) + 32;
 504           tmp.exp = EXP_BIAS + N_exp;
 505           
 506           round_to_int(&tmp);  /* Fortunately, this can't overflow to 2^64 */
 507           tmp.exp = EXP_BIAS + 63;
 508           normalize(&tmp);
 509           
 510           tmp.exp = EXP_BIAS + expdif - N_exp;
 511           
 512           reg_mul(st1_ptr, &tmp, &tmp);
 513           reg_sub(FPU_st0_ptr, &tmp, FPU_st0_ptr);
 514           
 515           setcc(SW_C2);
 516         }
 517       control_word = old_cw;
 518 
 519       if ( FPU_st0_ptr->exp <= EXP_UNDER )
 520         arith_underflow(FPU_st0_ptr);
 521       return;
 522     }
 523   else if ( (FPU_st0_tag == TW_Empty) | (st1_tag == TW_Empty) )
 524     { stack_underflow(); return; }
 525   else if ( FPU_st0_tag == TW_Zero )
 526     {
 527       if ( (st1_tag == TW_Valid) || (st1_tag == TW_Infinity) )
 528         { setcc(0); return; }
 529       if ( st1_tag == TW_Zero )
 530         { arith_invalid(FPU_st0_ptr); return; }
 531     }
 532 
 533   if ( (FPU_st0_tag == TW_NaN) | (st1_tag == TW_NaN) )
 534     { real_2op_NaN(FPU_st0_ptr, st1_ptr, FPU_st0_ptr); return; }
 535   else if ( FPU_st0_tag == TW_Infinity )
 536     { arith_invalid(FPU_st0_ptr); return; }
 537 #ifdef PARANOID
 538   else
 539     EXCEPTION(EX_INTERNAL | 0x118);
 540 #endif PARANOID
 541 
 542 }
 543 
 544 
 545 /* ST(1) <- ST(1) * log ST;  pop ST */
 546 static void fyl2x(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 547 {
 548   FPU_REG *st1_ptr = &st(1);
 549   char st1_tag = st1_ptr->tag;
 550 
 551   if ( !((FPU_st0_tag ^ TW_Valid) | (st1_tag ^ TW_Valid)) )
 552     {
 553       if ( FPU_st0_ptr->sign == SIGN_POS )
 554         {
 555           poly_l2(FPU_st0_ptr, FPU_st0_ptr);
 556           
 557           reg_mul(FPU_st0_ptr, st1_ptr, st1_ptr);
 558           pop(); FPU_st0_ptr = &st(0);
 559           if ( FPU_st0_ptr->exp <= EXP_UNDER )
 560             arith_underflow(FPU_st0_ptr);
 561           else if ( FPU_st0_ptr->exp >= EXP_OVER )
 562             arith_overflow(FPU_st0_ptr);
 563         }
 564       else
 565         {
 566           /* negative   */
 567           pop(); FPU_st0_ptr = &st(0);
 568           arith_invalid(FPU_st0_ptr);
 569         }
 570       return;
 571     }
 572 
 573   if ( (FPU_st0_tag == TW_Empty) || (st1_tag == TW_Empty) )
 574     { stack_underflow(); return; }
 575 
 576   if ( (FPU_st0_tag == TW_NaN) || (st1_tag == TW_NaN) )
 577     {
 578       real_2op_NaN(FPU_st0_ptr, st1_ptr, st1_ptr);
 579       pop();
 580       return;
 581     }
 582 
 583   if ( (FPU_st0_tag <= TW_Zero) && (st1_tag <= TW_Zero) )
 584     {
 585       /* one of the args is zero, the other valid, or both zero */
 586       if ( FPU_st0_tag == TW_Zero )
 587         {
 588           pop(); FPU_st0_ptr = &st(0);
 589           if ( FPU_st0_ptr->tag == TW_Zero )
 590             arith_invalid(FPU_st0_ptr);
 591           else
 592             divide_by_zero(st1_ptr->sign ^ SIGN_NEG, FPU_st0_ptr);
 593           return;
 594         }
 595       if ( st1_ptr->sign == SIGN_POS )
 596         {
 597           /* Zero is the valid answer */
 598           char sign = FPU_st0_ptr->sign;
 599           if ( FPU_st0_ptr->exp < EXP_BIAS ) sign ^= SIGN_NEG;
 600           pop(); FPU_st0_ptr = &st(0);
 601           reg_move(&CONST_Z, FPU_st0_ptr);
 602           FPU_st0_ptr->sign = sign;
 603           return;
 604         }
 605       pop(); FPU_st0_ptr = &st(0);
 606       arith_invalid(FPU_st0_ptr);
 607       return;
 608     }
 609 
 610   /* One or both arg must be an infinity */
 611   if ( FPU_st0_tag == TW_Infinity )
 612     {
 613       if ( (FPU_st0_ptr->sign == SIGN_NEG) || (st1_tag == TW_Zero) )
 614         { pop(); FPU_st0_ptr = &st(0); arith_invalid(FPU_st0_ptr); return; }
 615       else
 616         {
 617           char sign = st1_ptr->sign;
 618           pop(); FPU_st0_ptr = &st(0);
 619           reg_move(&CONST_INF, FPU_st0_ptr);
 620           FPU_st0_ptr->sign = sign;
 621           return;
 622         }
 623     }
 624 
 625   /* st(1) must be infinity here */
 626   if ( (FPU_st0_tag == TW_Valid) && (FPU_st0_ptr->sign == SIGN_POS) )
 627     {
 628       if ( FPU_st0_ptr->exp >= EXP_BIAS )
 629         {
 630           if ( (FPU_st0_ptr->exp == EXP_BIAS) &&
 631               (FPU_st0_ptr->sigh == 0x80000000) &&
 632               (FPU_st0_ptr->sigl == 0) )
 633             {
 634               pop(); FPU_st0_ptr = &st(0);
 635               arith_invalid(FPU_st0_ptr);
 636               return;
 637             }
 638           pop();
 639           return;
 640         }
 641       else
 642         {
 643           pop(); FPU_st0_ptr = &st(0);
 644           FPU_st0_ptr->sign ^= SIGN_NEG;
 645           return;
 646         }
 647     }
 648   /* st(0) must be zero or negative */
 649   pop(); FPU_st0_ptr = &st(0);
 650   arith_invalid(FPU_st0_ptr);
 651   return;
 652 }
 653 
 654 
 655 static void fpatan(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 656 {
 657   FPU_REG *st1_ptr = &st(1);
 658   char st1_tag = st1_ptr->tag;
 659 
 660   if ( !((FPU_st0_tag ^ TW_Valid) | (st1_tag ^ TW_Valid)) )
 661     {
 662       FPU_REG sum;
 663       int quadrant = st1_ptr->sign | ((FPU_st0_ptr->sign)<<1);
 664       st1_ptr->sign = FPU_st0_ptr->sign = SIGN_POS;
 665       if (compare(st1_ptr) == COMP_A_LT_B)
 666         {
 667           quadrant |= 4;
 668           reg_div(FPU_st0_ptr, st1_ptr, &sum);
 669         }
 670       else
 671         reg_div(st1_ptr, FPU_st0_ptr, &sum);
 672       
 673       poly_atan(&sum);
 674       
 675       if (quadrant & 4)
 676         {
 677           reg_sub(&CONST_PI2, &sum, &sum);
 678         }
 679       if (quadrant & 2)
 680         {
 681           reg_sub(&CONST_PI, &sum, &sum);
 682         }
 683       if (quadrant & 1)
 684         sum.sign ^= SIGN_POS^SIGN_NEG;
 685       
 686       reg_move(&sum, st1_ptr);
 687       pop(); FPU_st0_ptr = &st(0);
 688       if ( FPU_st0_ptr->exp <= EXP_UNDER )
 689         arith_underflow(FPU_st0_ptr);
 690       return;
 691     }
 692 
 693   if ( (FPU_st0_tag == TW_Empty) || (st1_tag == TW_Empty) )
 694     { stack_underflow(); return; }
 695 
 696   if ( (FPU_st0_tag == TW_NaN) || (st1_tag == TW_NaN) )
 697     {
 698       real_2op_NaN(FPU_st0_ptr, st1_ptr, st1_ptr);
 699       pop();
 700       return;
 701     }
 702 
 703   if ( (FPU_st0_tag == TW_Infinity) || (st1_tag == TW_Infinity) )
 704     {
 705       char sign = st1_ptr->sign;
 706       if ( FPU_st0_tag == TW_Infinity )
 707         {
 708           if ( st1_tag == TW_Infinity )
 709             {
 710               if ( FPU_st0_ptr->sign == SIGN_POS )
 711                 { reg_move(&CONST_PI4, st1_ptr); }
 712               else
 713                 reg_add(&CONST_PI4, &CONST_PI2, st1_ptr);
 714             }
 715           else
 716             {
 717               if ( FPU_st0_ptr->sign == SIGN_POS )
 718                 { reg_move(&CONST_Z, st1_ptr); }
 719               else
 720                 reg_move(&CONST_PI, st1_ptr);
 721             }
 722         }
 723       else
 724         {
 725           reg_move(&CONST_PI2, st1_ptr);
 726         }
 727       st1_ptr->sign = sign;
 728       pop();
 729       return;
 730     }
 731 
 732  if ( st1_tag == TW_Zero )
 733     {
 734       char sign = st1_ptr->sign;
 735       /* st(0) must be valid or zero */
 736       if ( FPU_st0_ptr->sign == SIGN_POS )
 737         { reg_move(&CONST_Z, st1_ptr); }
 738       else
 739         reg_move(&CONST_PI, st1_ptr);
 740       st1_tag = sign;
 741       pop();
 742       return;
 743     }
 744   else if ( FPU_st0_tag == TW_Zero )
 745     {
 746       char sign = st1_ptr->sign;
 747       /* st(1) must be TW_Valid here */
 748       reg_move(&CONST_PI2, st1_ptr);
 749       st1_tag = sign;
 750       pop();
 751       return;
 752     }
 753 #ifdef PARANOID
 754   EXCEPTION(EX_INTERNAL | 0x220);
 755 #endif PARANOID
 756 }
 757 
 758 
 759 static void fprem(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 760 {
 761   fprem_kernel(RC_CHOP);
 762 }
 763 
 764 
 765 static void fprem1(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 766 {
 767   fprem_kernel(RC_RND);
 768 }
 769 
 770 
 771 static void fyl2xp1(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 772 {
 773   FPU_REG *st1_ptr = &st(1);
 774   char st1_tag = st1_ptr->tag;
 775 
 776   if ( !((FPU_st0_tag ^ TW_Valid) | (st1_tag ^ TW_Valid)) )
 777     {
 778       if ( poly_l2p1(FPU_st0_ptr, FPU_st0_ptr) )
 779         {
 780           arith_invalid(st1_ptr); pop(); return;
 781         }
 782       
 783       reg_mul(FPU_st0_ptr, st1_ptr, st1_ptr);
 784       pop();
 785       return;
 786     }
 787   else if ( (FPU_st0_tag == TW_Empty) | (st1_tag == TW_Empty) )
 788     stack_underflow();
 789   else if ( FPU_st0_tag == TW_Zero )
 790     {
 791       if ( st1_tag <= TW_Zero )
 792         {
 793           st1_ptr->sign ^= FPU_st0_ptr->sign;
 794           reg_move(FPU_st0_ptr, st1_ptr);
 795         }
 796       else if ( st1_tag == TW_Infinity )
 797         {
 798           arith_invalid(st1_ptr);
 799         }
 800       else if ( st1_tag == TW_NaN )
 801         {
 802           if ( !(st1_ptr->sigh & 0x40000000) )
 803             EXCEPTION(EX_Invalid);            /* signaling NaN */
 804           st1_ptr->sigh |= 0x40000000;     /* QNaN */
 805         }
 806 #ifdef PARANOID
 807       else
 808         {
 809           EXCEPTION(EX_INTERNAL | 0x116);
 810         }
 811 #endif PARANOID
 812       pop();
 813       return;
 814     }
 815   else if ( FPU_st0_tag == TW_NaN )
 816     {
 817       real_2op_NaN(FPU_st0_ptr, st1_ptr, st1_ptr);
 818       pop();
 819       return;
 820     }
 821   else if ( FPU_st0_tag == TW_Infinity )
 822     {
 823       if ( st1_tag == TW_NaN )
 824         real_2op_NaN(FPU_st0_ptr, st1_ptr, st1_ptr);
 825       else
 826         arith_invalid(st1_ptr);
 827       pop();
 828       return;
 829     }
 830 #ifdef PARANOID
 831   else
 832     EXCEPTION(EX_INTERNAL | 0x117);
 833 #endif PARANOID
 834 }
 835 
 836 
 837 static void fscale(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 838 {
 839   FPU_REG *st1_ptr = &st(1);
 840   char st1_tag = st1_ptr->tag;
 841 
 842   if ( !((FPU_st0_tag ^ TW_Valid) | (st1_tag ^ TW_Valid)) )
 843     {
 844       long scale;
 845       FPU_REG tmp;
 846 
 847       /* 2^31 is far too large, 2^-31 is far too small */
 848       if ( st1_ptr->exp > EXP_BIAS + 30 )
 849         {
 850           char sign;
 851           EXCEPTION(EX_Overflow);
 852           sign = FPU_st0_ptr->sign;
 853           reg_move(&CONST_INF, FPU_st0_ptr);
 854           FPU_st0_ptr->sign = sign;
 855           return;
 856         }
 857       else if ( st1_ptr->exp < EXP_BIAS - 30 )
 858         {
 859           EXCEPTION(EX_Underflow);
 860           reg_move(&CONST_Z, FPU_st0_ptr);
 861           return;
 862         }
 863 
 864       reg_move(st1_ptr, &tmp);
 865       round_to_int(&tmp);               /* This can never overflow here */
 866       scale = st1_ptr->sign ? -tmp.sigl : tmp.sigl;
 867       scale += FPU_st0_ptr->exp;
 868       FPU_st0_ptr->exp = scale;
 869 
 870       if ( scale <= EXP_UNDER )
 871         arith_underflow(FPU_st0_ptr);
 872       else if ( scale >= EXP_OVER )
 873         arith_overflow(FPU_st0_ptr);
 874 
 875       return;
 876     }
 877   else if ( FPU_st0_tag == TW_Valid )
 878     {
 879       if ( st1_tag == TW_Zero )
 880         { return; }
 881       if ( st1_tag == TW_Infinity )
 882         {
 883           char sign = st1_ptr->sign;
 884           if ( sign == SIGN_POS )
 885             { reg_move(&CONST_INF, FPU_st0_ptr); }
 886           else
 887               reg_move(&CONST_Z, FPU_st0_ptr);
 888           FPU_st0_ptr->sign = sign;
 889           return;
 890         }
 891       if ( st1_tag == TW_NaN )
 892         { real_2op_NaN(FPU_st0_ptr, st1_ptr, FPU_st0_ptr); return; }
 893     }
 894   else if ( FPU_st0_tag == TW_Zero )
 895     {
 896       if ( st1_tag <= TW_Zero ) { return; }
 897       else if ( st1_tag == TW_Infinity )
 898         {
 899           if ( st1_ptr->sign == SIGN_NEG )
 900             return;
 901           else
 902             { arith_invalid(FPU_st0_ptr); return; }
 903         }
 904       else if ( st1_tag == TW_NaN )
 905         { real_2op_NaN(FPU_st0_ptr, st1_ptr, FPU_st0_ptr); return; }
 906     }
 907   else if ( FPU_st0_tag == TW_Infinity )
 908     {
 909       if ( ((st1_tag == TW_Infinity) && (st1_ptr->sign == SIGN_POS))
 910           || (st1_tag <= TW_Zero) )
 911         return;
 912       else if ( st1_tag == TW_Infinity )
 913         { arith_invalid(FPU_st0_ptr); return; }
 914       else if ( st1_tag == TW_NaN )
 915         { real_2op_NaN(FPU_st0_ptr, st1_ptr, FPU_st0_ptr); return; }
 916     }
 917   else if ( FPU_st0_tag == TW_NaN )
 918     {
 919       if ( st1_tag != TW_Empty )
 920         { real_2op_NaN(FPU_st0_ptr, st1_ptr, FPU_st0_ptr); return; }
 921     }
 922 
 923 #ifdef PARANOID
 924   if ( !((FPU_st0_tag == TW_Empty) || (st1_tag == TW_Empty)) )
 925     {
 926       EXCEPTION(EX_INTERNAL | 0x115);
 927       return;
 928     }
 929 #endif
 930 
 931   /* At least one of st(0), st(1) must be empty */
 932   stack_underflow();
 933 
 934 }
 935 
 936 
 937 /*---------------------------------------------------------------------------*/
 938 
 939 static FUNC trig_table_a[] = {
 940   f2xm1, fyl2x, fptan, fpatan, fxtract, fprem1, fdecstp, fincstp
 941 };
 942 
 943 void trig_a(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 944 {
 945   (trig_table_a[FPU_rm])();
 946 }
 947 
 948 
 949 static FUNC trig_table_b[] =
 950   {
 951     fprem, fyl2xp1, fsqrt_, fsincos, frndint_, fscale, fsin, fcos
 952   };
 953 
 954 void trig_b(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 955 {
 956   (trig_table_b[FPU_rm])();
 957 }

/* [previous][next][first][last][top][bottom][index][help] */