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,1993                                                   |
   7  |                       W. Metzenthen, 22 Parker St, Ormond, Vic 3163,      |
   8  |                       Australia.  E-mail apm233m@vaxc.cc.monash.edu.au    |
   9  |                                                                           |
  10  |                                                                           |
  11  +---------------------------------------------------------------------------*/
  12 
  13 #include "fpu_system.h"
  14 #include "exception.h"
  15 #include "fpu_emu.h"
  16 #include "status_w.h"
  17 #include "control_w.h"
  18 #include "reg_constant.h"       
  19 
  20 
  21 
  22 static int trig_arg(FPU_REG *X)
     /* [previous][next][first][last][top][bottom][index][help] */
  23 {
  24   FPU_REG tmp, quot;
  25   int rv;
  26   long long q;
  27   int old_cw = control_word;
  28 
  29   control_word &= ~CW_RC;
  30   control_word |= RC_CHOP;
  31   
  32   reg_move(X, &quot);
  33   reg_div(&quot, &CONST_PI2, &quot, FULL_PRECISION);
  34 
  35   reg_move(&quot, &tmp);
  36   round_to_int(&tmp);
  37   if ( tmp.sigh & 0x80000000 )
  38     return -1;              /* |Arg| is >= 2^63 */
  39   tmp.exp = EXP_BIAS + 63;
  40   q = *(long long *)&(tmp.sigl);
  41   normalize(&tmp);
  42 
  43   reg_sub(&quot, &tmp, X, FULL_PRECISION);
  44   rv = q & 7;
  45   
  46   control_word = old_cw;
  47   return rv;;
  48 }
  49 
  50 
  51 /* Convert a long to register */
  52 void convert_l2reg(long *arg, FPU_REG *dest)
     /* [previous][next][first][last][top][bottom][index][help] */
  53 {
  54   long num = *arg;
  55   
  56   if (num == 0)
  57     { reg_move(&CONST_Z, dest); return; }
  58 
  59   if (num > 0)
  60     dest->sign = SIGN_POS;
  61   else
  62     { num = -num; dest->sign = SIGN_NEG; }
  63   
  64   dest->sigh = num;
  65   dest->sigl = 0;
  66   dest->exp = EXP_BIAS + 31;
  67   dest->tag = TW_Valid;
  68   normalize(dest);
  69 }
  70 
  71 
  72 static void single_arg_error(void)
     /* [previous][next][first][last][top][bottom][index][help] */
  73 {
  74   switch ( FPU_st0_tag )
  75     {
  76     case TW_NaN:
  77       if ( !(FPU_st0_ptr->sigh & 0x40000000) )   /* Signaling ? */
  78         {
  79           EXCEPTION(EX_Invalid);
  80           /* Convert to a QNaN */
  81           FPU_st0_ptr->sigh |= 0x40000000;
  82         }
  83       break;              /* return with a NaN in st(0) */
  84     case TW_Empty:
  85       stack_underflow();  /* Puts a QNaN in st(0) */
  86       break;
  87 #ifdef PARANOID
  88     default:
  89       EXCEPTION(EX_INTERNAL|0x0112);
  90 #endif PARANOID
  91     }
  92 }
  93 
  94 
  95 /*---------------------------------------------------------------------------*/
  96 
  97 static void f2xm1(void)
     /* [previous][next][first][last][top][bottom][index][help] */
  98 {
  99   switch ( FPU_st0_tag )
 100     {
 101     case TW_Valid:
 102       {
 103         FPU_REG rv, tmp;
 104         
 105         if ( FPU_st0_ptr->sign == SIGN_POS )
 106           {
 107             /* poly_2xm1(x) requires 0 < x < 1. */
 108             if ( poly_2xm1(FPU_st0_ptr, &rv) )
 109               return;                  /* error */
 110             reg_mul(&rv, FPU_st0_ptr, FPU_st0_ptr, FULL_PRECISION);
 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, FULL_PRECISION);
 118             poly_2xm1(&tmp, &rv);
 119             reg_mul(&rv, &tmp, &tmp, FULL_PRECISION);
 120             reg_sub(&tmp, &CONST_1, FPU_st0_ptr, FULL_PRECISION);
 121             FPU_st0_ptr->exp--;
 122             if ( FPU_st0_ptr->exp <= EXP_UNDER )
 123               arith_underflow(FPU_st0_ptr);
 124           }
 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, FULL_PRECISION);
 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, control_word);       /* Do the computation */
 290       
 291       FPU_st0_ptr->exp += expon >> 1;
 292       FPU_st0_ptr->sign = SIGN_POS;
 293     }
 294   else if ( FPU_st0_tag == TW_Zero )
 295     return;
 296   else if ( FPU_st0_tag == TW_Infinity )
 297     {
 298       if ( FPU_st0_ptr->sign == SIGN_NEG )
 299         arith_invalid(FPU_st0_ptr);
 300       return;
 301     }
 302   else
 303     { single_arg_error(); return; }
 304 
 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, FULL_PRECISION);
 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, FULL_PRECISION);
 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, FULL_PRECISION);
 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, FULL_PRECISION);
 488           reg_sub(FPU_st0_ptr, &tmp, FPU_st0_ptr, FULL_PRECISION);
 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, FULL_PRECISION);
 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, FULL_PRECISION);
 513           reg_sub(FPU_st0_ptr, &tmp, FPU_st0_ptr, FULL_PRECISION);
 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, FULL_PRECISION);
 558           pop(); FPU_st0_ptr = &st(0);
 559           if ( FPU_st0_ptr->exp <= EXP_UNDER )
 560             { arith_underflow(FPU_st0_ptr); return; }
 561           else if ( FPU_st0_ptr->exp >= EXP_OVER )
 562             { arith_overflow(FPU_st0_ptr); return; }
 563         }
 564       else
 565         {
 566           /* negative   */
 567           pop(); FPU_st0_ptr = &st(0);
 568           arith_invalid(FPU_st0_ptr);
 569           return;
 570         }
 571     }
 572   else if ( (FPU_st0_tag == TW_Empty) || (st1_tag == TW_Empty) )
 573     { stack_underflow(); return; }
 574   else if ( (FPU_st0_tag == TW_NaN) || (st1_tag == TW_NaN) )
 575     {
 576       real_2op_NaN(FPU_st0_ptr, st1_ptr, st1_ptr);
 577       pop();
 578       return;
 579     }
 580   else if ( (FPU_st0_tag <= TW_Zero) && (st1_tag <= TW_Zero) )
 581     {
 582       /* one of the args is zero, the other valid, or both zero */
 583       if ( FPU_st0_tag == TW_Zero )
 584         {
 585           pop(); FPU_st0_ptr = &st(0);
 586           if ( FPU_st0_ptr->tag == TW_Zero )
 587             arith_invalid(FPU_st0_ptr);
 588           else
 589             divide_by_zero(st1_ptr->sign ^ SIGN_NEG, FPU_st0_ptr);
 590           return;
 591         }
 592       else if ( st1_ptr->sign == SIGN_POS )
 593         {
 594           /* Zero is the valid answer */
 595           char sign = FPU_st0_ptr->sign;
 596           if ( FPU_st0_ptr->exp < EXP_BIAS ) sign ^= SIGN_NEG;
 597           pop(); FPU_st0_ptr = &st(0);
 598           reg_move(&CONST_Z, FPU_st0_ptr);
 599           FPU_st0_ptr->sign = sign;
 600           return;
 601         }
 602       else
 603         {
 604           pop(); FPU_st0_ptr = &st(0);
 605           arith_invalid(FPU_st0_ptr);
 606           return;
 607         }
 608     }
 609   /* One or both arg must be an infinity */
 610   else if ( FPU_st0_tag == TW_Infinity )
 611     {
 612       if ( (FPU_st0_ptr->sign == SIGN_NEG) || (st1_tag == TW_Zero) )
 613         { pop(); FPU_st0_ptr = &st(0); arith_invalid(FPU_st0_ptr); return; }
 614       else
 615         {
 616           char sign = st1_ptr->sign;
 617           pop(); FPU_st0_ptr = &st(0);
 618           reg_move(&CONST_INF, FPU_st0_ptr);
 619           FPU_st0_ptr->sign = sign;
 620           return;
 621         }
 622     }
 623   /* st(1) must be infinity here */
 624   else if ( (FPU_st0_tag == TW_Valid) && (FPU_st0_ptr->sign == SIGN_POS) )
 625     {
 626       if ( FPU_st0_ptr->exp >= EXP_BIAS )
 627         {
 628           if ( (FPU_st0_ptr->exp == EXP_BIAS) &&
 629               (FPU_st0_ptr->sigh == 0x80000000) &&
 630               (FPU_st0_ptr->sigl == 0) )
 631             {
 632               pop(); FPU_st0_ptr = &st(0);
 633               arith_invalid(FPU_st0_ptr);
 634               return;
 635             }
 636           pop();
 637         }
 638       else
 639         {
 640           pop(); FPU_st0_ptr = &st(0);
 641           FPU_st0_ptr->sign ^= SIGN_NEG;
 642         }
 643       return;
 644     }
 645   else
 646     {
 647       /* st(0) must be zero or negative */
 648       pop(); FPU_st0_ptr = &st(0);
 649       arith_invalid(FPU_st0_ptr);
 650       return;
 651     }
 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, FULL_PRECISION);
 669         }
 670       else
 671         reg_div(st1_ptr, FPU_st0_ptr, &sum, FULL_PRECISION);
 672       
 673       poly_atan(&sum);
 674       
 675       if (quadrant & 4)
 676         {
 677           reg_sub(&CONST_PI2, &sum, &sum, FULL_PRECISION);
 678         }
 679       if (quadrant & 2)
 680         {
 681           reg_sub(&CONST_PI, &sum, &sum, FULL_PRECISION);
 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); return; }
 690     }
 691   else if ( (FPU_st0_tag == TW_Empty) || (st1_tag == TW_Empty) )
 692     { stack_underflow(); return; }
 693   else if ( (FPU_st0_tag == TW_NaN) || (st1_tag == TW_NaN) )
 694     {
 695       real_2op_NaN(FPU_st0_ptr, st1_ptr, st1_ptr);
 696       pop();
 697       return;
 698     }
 699   else if ( (FPU_st0_tag == TW_Infinity) || (st1_tag == TW_Infinity) )
 700     {
 701       char sign = st1_ptr->sign;
 702       if ( FPU_st0_tag == TW_Infinity )
 703         {
 704           if ( st1_tag == TW_Infinity )
 705             {
 706               if ( FPU_st0_ptr->sign == SIGN_POS )
 707                 { reg_move(&CONST_PI4, st1_ptr); }
 708               else
 709                 reg_add(&CONST_PI4, &CONST_PI2, st1_ptr, FULL_PRECISION);
 710             }
 711           else
 712             {
 713               if ( FPU_st0_ptr->sign == SIGN_POS )
 714                 { reg_move(&CONST_Z, st1_ptr); }
 715               else
 716                 reg_move(&CONST_PI, st1_ptr);
 717             }
 718         }
 719       else
 720         {
 721           reg_move(&CONST_PI2, st1_ptr);
 722         }
 723       st1_ptr->sign = sign;
 724       pop();
 725     }
 726   else if ( st1_tag == TW_Zero )
 727     {
 728       char sign = st1_ptr->sign;
 729       /* st(0) must be valid or zero */
 730       if ( FPU_st0_ptr->sign == SIGN_POS )
 731         { reg_move(&CONST_Z, st1_ptr); }
 732       else
 733         reg_move(&CONST_PI, st1_ptr);
 734       st1_ptr->sign = sign;
 735       pop();
 736     }
 737   else if ( FPU_st0_tag == TW_Zero )
 738     {
 739       char sign = st1_ptr->sign;
 740       /* st(1) must be TW_Valid here */
 741       reg_move(&CONST_PI2, st1_ptr);
 742       st1_ptr->sign = sign;
 743       pop();
 744     }
 745 #ifdef PARANOID
 746   else
 747     EXCEPTION(EX_INTERNAL | 0x220);
 748 #endif PARANOID
 749 }
 750 
 751 
 752 static void fprem(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 753 {
 754   fprem_kernel(RC_CHOP);
 755 }
 756 
 757 
 758 static void fprem1(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 759 {
 760   fprem_kernel(RC_RND);
 761 }
 762 
 763 
 764 static void fyl2xp1(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 765 {
 766   FPU_REG *st1_ptr = &st(1);
 767   char st1_tag = st1_ptr->tag;
 768 
 769   if ( !((FPU_st0_tag ^ TW_Valid) | (st1_tag ^ TW_Valid)) )
 770     {
 771       if ( poly_l2p1(FPU_st0_ptr, FPU_st0_ptr) )
 772         {
 773           arith_invalid(st1_ptr); pop(); return;
 774         }
 775       
 776       reg_mul(FPU_st0_ptr, st1_ptr, st1_ptr, FULL_PRECISION);
 777       pop();
 778     }
 779   else if ( (FPU_st0_tag == TW_Empty) | (st1_tag == TW_Empty) )
 780     { stack_underflow(); return; }
 781   else if ( FPU_st0_tag == TW_Zero )
 782     {
 783       if ( st1_tag <= TW_Zero )
 784         {
 785           st1_ptr->sign ^= FPU_st0_ptr->sign;
 786           reg_move(FPU_st0_ptr, st1_ptr);
 787         }
 788       else if ( st1_tag == TW_Infinity )
 789         {
 790           arith_invalid(st1_ptr);
 791           return;
 792         }
 793       else if ( st1_tag == TW_NaN )
 794         {
 795           if ( !(st1_ptr->sigh & 0x40000000) )
 796             EXCEPTION(EX_Invalid);            /* signaling NaN */
 797           st1_ptr->sigh |= 0x40000000;     /* QNaN */
 798           return;
 799         }
 800 #ifdef PARANOID
 801       else
 802         {
 803           EXCEPTION(EX_INTERNAL | 0x116);
 804           return;
 805         }
 806 #endif PARANOID
 807       pop();
 808     }
 809   else if ( FPU_st0_tag == TW_NaN )
 810     {
 811       real_2op_NaN(FPU_st0_ptr, st1_ptr, st1_ptr);
 812       pop();
 813       return;
 814     }
 815   else if ( FPU_st0_tag == TW_Infinity )
 816     {
 817       if ( st1_tag == TW_NaN )
 818         real_2op_NaN(FPU_st0_ptr, st1_ptr, st1_ptr);
 819       else
 820         arith_invalid(st1_ptr);
 821       pop();
 822       return;
 823     }
 824 #ifdef PARANOID
 825   else
 826     EXCEPTION(EX_INTERNAL | 0x117);
 827 #endif PARANOID
 828 }
 829 
 830 
 831 static void fscale(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 832 {
 833   FPU_REG *st1_ptr = &st(1);
 834   char st1_tag = st1_ptr->tag;
 835 
 836   if ( !((FPU_st0_tag ^ TW_Valid) | (st1_tag ^ TW_Valid)) )
 837     {
 838       long scale;
 839       FPU_REG tmp;
 840 
 841       /* 2^31 is far too large, 2^-31 is far too small */
 842       if ( st1_ptr->exp > EXP_BIAS + 30 )
 843         {
 844           char sign;
 845           EXCEPTION(EX_Overflow);
 846           sign = FPU_st0_ptr->sign;
 847           reg_move(&CONST_INF, FPU_st0_ptr);
 848           FPU_st0_ptr->sign = sign;
 849           return;
 850         }
 851       else if ( st1_ptr->exp < EXP_BIAS - 30 )
 852         {
 853           EXCEPTION(EX_Underflow);
 854           reg_move(&CONST_Z, FPU_st0_ptr);
 855           return;
 856         }
 857 
 858       reg_move(st1_ptr, &tmp);
 859       round_to_int(&tmp);               /* This can never overflow here */
 860       scale = st1_ptr->sign ? -tmp.sigl : tmp.sigl;
 861       scale += FPU_st0_ptr->exp;
 862       FPU_st0_ptr->exp = scale;
 863 
 864       if ( scale <= EXP_UNDER )
 865         { arith_underflow(FPU_st0_ptr); return; }
 866       else if ( scale >= EXP_OVER )
 867         { arith_overflow(FPU_st0_ptr); return; }
 868 
 869       return;
 870     }
 871   else if ( FPU_st0_tag == TW_Valid )
 872     {
 873       if ( st1_tag == TW_Zero )
 874         { return; }
 875       if ( st1_tag == TW_Infinity )
 876         {
 877           char sign = st1_ptr->sign;
 878           if ( sign == SIGN_POS )
 879             { reg_move(&CONST_INF, FPU_st0_ptr); }
 880           else
 881               reg_move(&CONST_Z, FPU_st0_ptr);
 882           FPU_st0_ptr->sign = sign;
 883           return;
 884         }
 885       if ( st1_tag == TW_NaN )
 886         { real_2op_NaN(FPU_st0_ptr, st1_ptr, FPU_st0_ptr); return; }
 887     }
 888   else if ( FPU_st0_tag == TW_Zero )
 889     {
 890       if ( st1_tag <= TW_Zero ) { return; }
 891       else if ( st1_tag == TW_Infinity )
 892         {
 893           if ( st1_ptr->sign == SIGN_NEG )
 894             return;
 895           else
 896             { arith_invalid(FPU_st0_ptr); return; }
 897         }
 898       else if ( st1_tag == TW_NaN )
 899         { real_2op_NaN(FPU_st0_ptr, st1_ptr, FPU_st0_ptr); return; }
 900     }
 901   else if ( FPU_st0_tag == TW_Infinity )
 902     {
 903       if ( ((st1_tag == TW_Infinity) && (st1_ptr->sign == SIGN_POS))
 904           || (st1_tag <= TW_Zero) )
 905         return;
 906       else if ( st1_tag == TW_Infinity )
 907         { arith_invalid(FPU_st0_ptr); return; }
 908       else if ( st1_tag == TW_NaN )
 909         { real_2op_NaN(FPU_st0_ptr, st1_ptr, FPU_st0_ptr); return; }
 910     }
 911   else if ( FPU_st0_tag == TW_NaN )
 912     {
 913       if ( st1_tag != TW_Empty )
 914         { real_2op_NaN(FPU_st0_ptr, st1_ptr, FPU_st0_ptr); return; }
 915     }
 916 
 917 #ifdef PARANOID
 918   if ( !((FPU_st0_tag == TW_Empty) || (st1_tag == TW_Empty)) )
 919     {
 920       EXCEPTION(EX_INTERNAL | 0x115);
 921       return;
 922     }
 923 #endif
 924 
 925   /* At least one of st(0), st(1) must be empty */
 926   stack_underflow();
 927 
 928 }
 929 
 930 
 931 /*---------------------------------------------------------------------------*/
 932 
 933 static FUNC trig_table_a[] = {
 934   f2xm1, fyl2x, fptan, fpatan, fxtract, fprem1, fdecstp, fincstp
 935 };
 936 
 937 void trig_a(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 938 {
 939   (trig_table_a[FPU_rm])();
 940 }
 941 
 942 
 943 static FUNC trig_table_b[] =
 944   {
 945     fprem, fyl2xp1, fsqrt_, fsincos, frndint_, fscale, fsin, fcos
 946   };
 947 
 948 void trig_b(void)
     /* [previous][next][first][last][top][bottom][index][help] */
 949 {
 950   (trig_table_b[FPU_rm])();
 951 }

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