Back to home page

OSCL-LXR
 
 

     

0001 // SPDX-License-Identifier: GPL-2.0
0002 /*---------------------------------------------------------------------------+
0003  |  fpu_trig.c                                                               |
0004  |                                                                           |
0005  | Implementation of the FPU "transcendental" functions.                     |
0006  |                                                                           |
0007  | Copyright (C) 1992,1993,1994,1997,1999                                    |
0008  |                       W. Metzenthen, 22 Parker St, Ormond, Vic 3163,      |
0009  |                       Australia.  E-mail   billm@melbpc.org.au            |
0010  |                                                                           |
0011  |                                                                           |
0012  +---------------------------------------------------------------------------*/
0013 
0014 #include "fpu_system.h"
0015 #include "exception.h"
0016 #include "fpu_emu.h"
0017 #include "status_w.h"
0018 #include "control_w.h"
0019 #include "reg_constant.h"
0020 
0021 static void rem_kernel(unsigned long long st0, unsigned long long *y,
0022                unsigned long long st1, unsigned long long q, int n);
0023 
0024 #define BETTER_THAN_486
0025 
0026 #define FCOS  4
0027 
0028 /* Used only by fptan, fsin, fcos, and fsincos. */
0029 /* This routine produces very accurate results, similar to
0030    using a value of pi with more than 128 bits precision. */
0031 /* Limited measurements show no results worse than 64 bit precision
0032    except for the results for arguments close to 2^63, where the
0033    precision of the result sometimes degrades to about 63.9 bits */
0034 static int trig_arg(FPU_REG *st0_ptr, int even)
0035 {
0036     FPU_REG tmp;
0037     u_char tmptag;
0038     unsigned long long q;
0039     int old_cw = control_word, saved_status = partial_status;
0040     int tag, st0_tag = TAG_Valid;
0041 
0042     if (exponent(st0_ptr) >= 63) {
0043         partial_status |= SW_C2;    /* Reduction incomplete. */
0044         return -1;
0045     }
0046 
0047     control_word &= ~CW_RC;
0048     control_word |= RC_CHOP;
0049 
0050     setpositive(st0_ptr);
0051     tag = FPU_u_div(st0_ptr, &CONST_PI2, &tmp, PR_64_BITS | RC_CHOP | 0x3f,
0052             SIGN_POS);
0053 
0054     FPU_round_to_int(&tmp, tag);    /* Fortunately, this can't overflow
0055                        to 2^64 */
0056     q = significand(&tmp);
0057     if (q) {
0058         rem_kernel(significand(st0_ptr),
0059                &significand(&tmp),
0060                significand(&CONST_PI2),
0061                q, exponent(st0_ptr) - exponent(&CONST_PI2));
0062         setexponent16(&tmp, exponent(&CONST_PI2));
0063         st0_tag = FPU_normalize(&tmp);
0064         FPU_copy_to_reg0(&tmp, st0_tag);
0065     }
0066 
0067     if ((even && !(q & 1)) || (!even && (q & 1))) {
0068         st0_tag =
0069             FPU_sub(REV | LOADED | TAG_Valid, (int)&CONST_PI2,
0070                 FULL_PRECISION);
0071 
0072 #ifdef BETTER_THAN_486
0073         /* So far, the results are exact but based upon a 64 bit
0074            precision approximation to pi/2. The technique used
0075            now is equivalent to using an approximation to pi/2 which
0076            is accurate to about 128 bits. */
0077         if ((exponent(st0_ptr) <= exponent(&CONST_PI2extra) + 64)
0078             || (q > 1)) {
0079             /* This code gives the effect of having pi/2 to better than
0080                128 bits precision. */
0081 
0082             significand(&tmp) = q + 1;
0083             setexponent16(&tmp, 63);
0084             FPU_normalize(&tmp);
0085             tmptag =
0086                 FPU_u_mul(&CONST_PI2extra, &tmp, &tmp,
0087                       FULL_PRECISION, SIGN_POS,
0088                       exponent(&CONST_PI2extra) +
0089                       exponent(&tmp));
0090             setsign(&tmp, getsign(&CONST_PI2extra));
0091             st0_tag = FPU_add(&tmp, tmptag, 0, FULL_PRECISION);
0092             if (signnegative(st0_ptr)) {
0093                 /* CONST_PI2extra is negative, so the result of the addition
0094                    can be negative. This means that the argument is actually
0095                    in a different quadrant. The correction is always < pi/2,
0096                    so it can't overflow into yet another quadrant. */
0097                 setpositive(st0_ptr);
0098                 q++;
0099             }
0100         }
0101 #endif /* BETTER_THAN_486 */
0102     }
0103 #ifdef BETTER_THAN_486
0104     else {
0105         /* So far, the results are exact but based upon a 64 bit
0106            precision approximation to pi/2. The technique used
0107            now is equivalent to using an approximation to pi/2 which
0108            is accurate to about 128 bits. */
0109         if (((q > 0)
0110              && (exponent(st0_ptr) <= exponent(&CONST_PI2extra) + 64))
0111             || (q > 1)) {
0112             /* This code gives the effect of having p/2 to better than
0113                128 bits precision. */
0114 
0115             significand(&tmp) = q;
0116             setexponent16(&tmp, 63);
0117             FPU_normalize(&tmp);    /* This must return TAG_Valid */
0118             tmptag =
0119                 FPU_u_mul(&CONST_PI2extra, &tmp, &tmp,
0120                       FULL_PRECISION, SIGN_POS,
0121                       exponent(&CONST_PI2extra) +
0122                       exponent(&tmp));
0123             setsign(&tmp, getsign(&CONST_PI2extra));
0124             st0_tag = FPU_sub(LOADED | (tmptag & 0x0f), (int)&tmp,
0125                       FULL_PRECISION);
0126             if ((exponent(st0_ptr) == exponent(&CONST_PI2)) &&
0127                 ((st0_ptr->sigh > CONST_PI2.sigh)
0128                  || ((st0_ptr->sigh == CONST_PI2.sigh)
0129                  && (st0_ptr->sigl > CONST_PI2.sigl)))) {
0130                 /* CONST_PI2extra is negative, so the result of the
0131                    subtraction can be larger than pi/2. This means
0132                    that the argument is actually in a different quadrant.
0133                    The correction is always < pi/2, so it can't overflow
0134                    into yet another quadrant. */
0135                 st0_tag =
0136                     FPU_sub(REV | LOADED | TAG_Valid,
0137                         (int)&CONST_PI2, FULL_PRECISION);
0138                 q++;
0139             }
0140         }
0141     }
0142 #endif /* BETTER_THAN_486 */
0143 
0144     FPU_settag0(st0_tag);
0145     control_word = old_cw;
0146     partial_status = saved_status & ~SW_C2; /* Reduction complete. */
0147 
0148     return (q & 3) | even;
0149 }
0150 
0151 /* Convert a long to register */
0152 static void convert_l2reg(long const *arg, int deststnr)
0153 {
0154     int tag;
0155     long num = *arg;
0156     u_char sign;
0157     FPU_REG *dest = &st(deststnr);
0158 
0159     if (num == 0) {
0160         FPU_copy_to_regi(&CONST_Z, TAG_Zero, deststnr);
0161         return;
0162     }
0163 
0164     if (num > 0) {
0165         sign = SIGN_POS;
0166     } else {
0167         num = -num;
0168         sign = SIGN_NEG;
0169     }
0170 
0171     dest->sigh = num;
0172     dest->sigl = 0;
0173     setexponent16(dest, 31);
0174     tag = FPU_normalize(dest);
0175     FPU_settagi(deststnr, tag);
0176     setsign(dest, sign);
0177     return;
0178 }
0179 
0180 static void single_arg_error(FPU_REG *st0_ptr, u_char st0_tag)
0181 {
0182     if (st0_tag == TAG_Empty)
0183         FPU_stack_underflow();  /* Puts a QNaN in st(0) */
0184     else if (st0_tag == TW_NaN)
0185         real_1op_NaN(st0_ptr);  /* return with a NaN in st(0) */
0186 #ifdef PARANOID
0187     else
0188         EXCEPTION(EX_INTERNAL | 0x0112);
0189 #endif /* PARANOID */
0190 }
0191 
0192 static void single_arg_2_error(FPU_REG *st0_ptr, u_char st0_tag)
0193 {
0194     int isNaN;
0195 
0196     switch (st0_tag) {
0197     case TW_NaN:
0198         isNaN = (exponent(st0_ptr) == EXP_OVER)
0199             && (st0_ptr->sigh & 0x80000000);
0200         if (isNaN && !(st0_ptr->sigh & 0x40000000)) {   /* Signaling ? */
0201             EXCEPTION(EX_Invalid);
0202             if (control_word & CW_Invalid) {
0203                 /* The masked response */
0204                 /* Convert to a QNaN */
0205                 st0_ptr->sigh |= 0x40000000;
0206                 push();
0207                 FPU_copy_to_reg0(st0_ptr, TAG_Special);
0208             }
0209         } else if (isNaN) {
0210             /* A QNaN */
0211             push();
0212             FPU_copy_to_reg0(st0_ptr, TAG_Special);
0213         } else {
0214             /* pseudoNaN or other unsupported */
0215             EXCEPTION(EX_Invalid);
0216             if (control_word & CW_Invalid) {
0217                 /* The masked response */
0218                 FPU_copy_to_reg0(&CONST_QNaN, TAG_Special);
0219                 push();
0220                 FPU_copy_to_reg0(&CONST_QNaN, TAG_Special);
0221             }
0222         }
0223         break;      /* return with a NaN in st(0) */
0224 #ifdef PARANOID
0225     default:
0226         EXCEPTION(EX_INTERNAL | 0x0112);
0227 #endif /* PARANOID */
0228     }
0229 }
0230 
0231 /*---------------------------------------------------------------------------*/
0232 
0233 static void f2xm1(FPU_REG *st0_ptr, u_char tag)
0234 {
0235     FPU_REG a;
0236 
0237     clear_C1();
0238 
0239     if (tag == TAG_Valid) {
0240         /* For an 80486 FPU, the result is undefined if the arg is >= 1.0 */
0241         if (exponent(st0_ptr) < 0) {
0242               denormal_arg:
0243 
0244             FPU_to_exp16(st0_ptr, &a);
0245 
0246             /* poly_2xm1(x) requires 0 < st(0) < 1. */
0247             poly_2xm1(getsign(st0_ptr), &a, st0_ptr);
0248         }
0249         set_precision_flag_up();    /* 80486 appears to always do this */
0250         return;
0251     }
0252 
0253     if (tag == TAG_Zero)
0254         return;
0255 
0256     if (tag == TAG_Special)
0257         tag = FPU_Special(st0_ptr);
0258 
0259     switch (tag) {
0260     case TW_Denormal:
0261         if (denormal_operand() < 0)
0262             return;
0263         goto denormal_arg;
0264     case TW_Infinity:
0265         if (signnegative(st0_ptr)) {
0266             /* -infinity gives -1 (p16-10) */
0267             FPU_copy_to_reg0(&CONST_1, TAG_Valid);
0268             setnegative(st0_ptr);
0269         }
0270         return;
0271     default:
0272         single_arg_error(st0_ptr, tag);
0273     }
0274 }
0275 
0276 static void fptan(FPU_REG *st0_ptr, u_char st0_tag)
0277 {
0278     FPU_REG *st_new_ptr;
0279     int q;
0280     u_char arg_sign = getsign(st0_ptr);
0281 
0282     /* Stack underflow has higher priority */
0283     if (st0_tag == TAG_Empty) {
0284         FPU_stack_underflow();  /* Puts a QNaN in st(0) */
0285         if (control_word & CW_Invalid) {
0286             st_new_ptr = &st(-1);
0287             push();
0288             FPU_stack_underflow();  /* Puts a QNaN in the new st(0) */
0289         }
0290         return;
0291     }
0292 
0293     if (STACK_OVERFLOW) {
0294         FPU_stack_overflow();
0295         return;
0296     }
0297 
0298     if (st0_tag == TAG_Valid) {
0299         if (exponent(st0_ptr) > -40) {
0300             if ((q = trig_arg(st0_ptr, 0)) == -1) {
0301                 /* Operand is out of range */
0302                 return;
0303             }
0304 
0305             poly_tan(st0_ptr);
0306             setsign(st0_ptr, (q & 1) ^ (arg_sign != 0));
0307             set_precision_flag_up();    /* We do not really know if up or down */
0308         } else {
0309             /* For a small arg, the result == the argument */
0310             /* Underflow may happen */
0311 
0312               denormal_arg:
0313 
0314             FPU_to_exp16(st0_ptr, st0_ptr);
0315 
0316             st0_tag =
0317                 FPU_round(st0_ptr, 1, 0, FULL_PRECISION, arg_sign);
0318             FPU_settag0(st0_tag);
0319         }
0320         push();
0321         FPU_copy_to_reg0(&CONST_1, TAG_Valid);
0322         return;
0323     }
0324 
0325     if (st0_tag == TAG_Zero) {
0326         push();
0327         FPU_copy_to_reg0(&CONST_1, TAG_Valid);
0328         setcc(0);
0329         return;
0330     }
0331 
0332     if (st0_tag == TAG_Special)
0333         st0_tag = FPU_Special(st0_ptr);
0334 
0335     if (st0_tag == TW_Denormal) {
0336         if (denormal_operand() < 0)
0337             return;
0338 
0339         goto denormal_arg;
0340     }
0341 
0342     if (st0_tag == TW_Infinity) {
0343         /* The 80486 treats infinity as an invalid operand */
0344         if (arith_invalid(0) >= 0) {
0345             st_new_ptr = &st(-1);
0346             push();
0347             arith_invalid(0);
0348         }
0349         return;
0350     }
0351 
0352     single_arg_2_error(st0_ptr, st0_tag);
0353 }
0354 
0355 static void fxtract(FPU_REG *st0_ptr, u_char st0_tag)
0356 {
0357     FPU_REG *st_new_ptr;
0358     u_char sign;
0359     register FPU_REG *st1_ptr = st0_ptr;    /* anticipate */
0360 
0361     if (STACK_OVERFLOW) {
0362         FPU_stack_overflow();
0363         return;
0364     }
0365 
0366     clear_C1();
0367 
0368     if (st0_tag == TAG_Valid) {
0369         long e;
0370 
0371         push();
0372         sign = getsign(st1_ptr);
0373         reg_copy(st1_ptr, st_new_ptr);
0374         setexponent16(st_new_ptr, exponent(st_new_ptr));
0375 
0376           denormal_arg:
0377 
0378         e = exponent16(st_new_ptr);
0379         convert_l2reg(&e, 1);
0380         setexponentpos(st_new_ptr, 0);
0381         setsign(st_new_ptr, sign);
0382         FPU_settag0(TAG_Valid); /* Needed if arg was a denormal */
0383         return;
0384     } else if (st0_tag == TAG_Zero) {
0385         sign = getsign(st0_ptr);
0386 
0387         if (FPU_divide_by_zero(0, SIGN_NEG) < 0)
0388             return;
0389 
0390         push();
0391         FPU_copy_to_reg0(&CONST_Z, TAG_Zero);
0392         setsign(st_new_ptr, sign);
0393         return;
0394     }
0395 
0396     if (st0_tag == TAG_Special)
0397         st0_tag = FPU_Special(st0_ptr);
0398 
0399     if (st0_tag == TW_Denormal) {
0400         if (denormal_operand() < 0)
0401             return;
0402 
0403         push();
0404         sign = getsign(st1_ptr);
0405         FPU_to_exp16(st1_ptr, st_new_ptr);
0406         goto denormal_arg;
0407     } else if (st0_tag == TW_Infinity) {
0408         sign = getsign(st0_ptr);
0409         setpositive(st0_ptr);
0410         push();
0411         FPU_copy_to_reg0(&CONST_INF, TAG_Special);
0412         setsign(st_new_ptr, sign);
0413         return;
0414     } else if (st0_tag == TW_NaN) {
0415         if (real_1op_NaN(st0_ptr) < 0)
0416             return;
0417 
0418         push();
0419         FPU_copy_to_reg0(st0_ptr, TAG_Special);
0420         return;
0421     } else if (st0_tag == TAG_Empty) {
0422         /* Is this the correct behaviour? */
0423         if (control_word & EX_Invalid) {
0424             FPU_stack_underflow();
0425             push();
0426             FPU_stack_underflow();
0427         } else
0428             EXCEPTION(EX_StackUnder);
0429     }
0430 #ifdef PARANOID
0431     else
0432         EXCEPTION(EX_INTERNAL | 0x119);
0433 #endif /* PARANOID */
0434 }
0435 
0436 static void fdecstp(void)
0437 {
0438     clear_C1();
0439     top--;
0440 }
0441 
0442 static void fincstp(void)
0443 {
0444     clear_C1();
0445     top++;
0446 }
0447 
0448 static void fsqrt_(FPU_REG *st0_ptr, u_char st0_tag)
0449 {
0450     int expon;
0451 
0452     clear_C1();
0453 
0454     if (st0_tag == TAG_Valid) {
0455         u_char tag;
0456 
0457         if (signnegative(st0_ptr)) {
0458             arith_invalid(0);   /* sqrt(negative) is invalid */
0459             return;
0460         }
0461 
0462         /* make st(0) in  [1.0 .. 4.0) */
0463         expon = exponent(st0_ptr);
0464 
0465           denormal_arg:
0466 
0467         setexponent16(st0_ptr, (expon & 1));
0468 
0469         /* Do the computation, the sign of the result will be positive. */
0470         tag = wm_sqrt(st0_ptr, 0, 0, control_word, SIGN_POS);
0471         addexponent(st0_ptr, expon >> 1);
0472         FPU_settag0(tag);
0473         return;
0474     }
0475 
0476     if (st0_tag == TAG_Zero)
0477         return;
0478 
0479     if (st0_tag == TAG_Special)
0480         st0_tag = FPU_Special(st0_ptr);
0481 
0482     if (st0_tag == TW_Infinity) {
0483         if (signnegative(st0_ptr))
0484             arith_invalid(0);   /* sqrt(-Infinity) is invalid */
0485         return;
0486     } else if (st0_tag == TW_Denormal) {
0487         if (signnegative(st0_ptr)) {
0488             arith_invalid(0);   /* sqrt(negative) is invalid */
0489             return;
0490         }
0491 
0492         if (denormal_operand() < 0)
0493             return;
0494 
0495         FPU_to_exp16(st0_ptr, st0_ptr);
0496 
0497         expon = exponent16(st0_ptr);
0498 
0499         goto denormal_arg;
0500     }
0501 
0502     single_arg_error(st0_ptr, st0_tag);
0503 
0504 }
0505 
0506 static void frndint_(FPU_REG *st0_ptr, u_char st0_tag)
0507 {
0508     int flags, tag;
0509 
0510     if (st0_tag == TAG_Valid) {
0511         u_char sign;
0512 
0513           denormal_arg:
0514 
0515         sign = getsign(st0_ptr);
0516 
0517         if (exponent(st0_ptr) > 63)
0518             return;
0519 
0520         if (st0_tag == TW_Denormal) {
0521             if (denormal_operand() < 0)
0522                 return;
0523         }
0524 
0525         /* Fortunately, this can't overflow to 2^64 */
0526         if ((flags = FPU_round_to_int(st0_ptr, st0_tag)))
0527             set_precision_flag(flags);
0528 
0529         setexponent16(st0_ptr, 63);
0530         tag = FPU_normalize(st0_ptr);
0531         setsign(st0_ptr, sign);
0532         FPU_settag0(tag);
0533         return;
0534     }
0535 
0536     if (st0_tag == TAG_Zero)
0537         return;
0538 
0539     if (st0_tag == TAG_Special)
0540         st0_tag = FPU_Special(st0_ptr);
0541 
0542     if (st0_tag == TW_Denormal)
0543         goto denormal_arg;
0544     else if (st0_tag == TW_Infinity)
0545         return;
0546     else
0547         single_arg_error(st0_ptr, st0_tag);
0548 }
0549 
0550 static int f_sin(FPU_REG *st0_ptr, u_char tag)
0551 {
0552     u_char arg_sign = getsign(st0_ptr);
0553 
0554     if (tag == TAG_Valid) {
0555         int q;
0556 
0557         if (exponent(st0_ptr) > -40) {
0558             if ((q = trig_arg(st0_ptr, 0)) == -1) {
0559                 /* Operand is out of range */
0560                 return 1;
0561             }
0562 
0563             poly_sine(st0_ptr);
0564 
0565             if (q & 2)
0566                 changesign(st0_ptr);
0567 
0568             setsign(st0_ptr, getsign(st0_ptr) ^ arg_sign);
0569 
0570             /* We do not really know if up or down */
0571             set_precision_flag_up();
0572             return 0;
0573         } else {
0574             /* For a small arg, the result == the argument */
0575             set_precision_flag_up();    /* Must be up. */
0576             return 0;
0577         }
0578     }
0579 
0580     if (tag == TAG_Zero) {
0581         setcc(0);
0582         return 0;
0583     }
0584 
0585     if (tag == TAG_Special)
0586         tag = FPU_Special(st0_ptr);
0587 
0588     if (tag == TW_Denormal) {
0589         if (denormal_operand() < 0)
0590             return 1;
0591 
0592         /* For a small arg, the result == the argument */
0593         /* Underflow may happen */
0594         FPU_to_exp16(st0_ptr, st0_ptr);
0595 
0596         tag = FPU_round(st0_ptr, 1, 0, FULL_PRECISION, arg_sign);
0597 
0598         FPU_settag0(tag);
0599 
0600         return 0;
0601     } else if (tag == TW_Infinity) {
0602         /* The 80486 treats infinity as an invalid operand */
0603         arith_invalid(0);
0604         return 1;
0605     } else {
0606         single_arg_error(st0_ptr, tag);
0607         return 1;
0608     }
0609 }
0610 
0611 static void fsin(FPU_REG *st0_ptr, u_char tag)
0612 {
0613     f_sin(st0_ptr, tag);
0614 }
0615 
0616 static int f_cos(FPU_REG *st0_ptr, u_char tag)
0617 {
0618     u_char st0_sign;
0619 
0620     st0_sign = getsign(st0_ptr);
0621 
0622     if (tag == TAG_Valid) {
0623         int q;
0624 
0625         if (exponent(st0_ptr) > -40) {
0626             if ((exponent(st0_ptr) < 0)
0627                 || ((exponent(st0_ptr) == 0)
0628                 && (significand(st0_ptr) <=
0629                     0xc90fdaa22168c234LL))) {
0630                 poly_cos(st0_ptr);
0631 
0632                 /* We do not really know if up or down */
0633                 set_precision_flag_down();
0634 
0635                 return 0;
0636             } else if ((q = trig_arg(st0_ptr, FCOS)) != -1) {
0637                 poly_sine(st0_ptr);
0638 
0639                 if ((q + 1) & 2)
0640                     changesign(st0_ptr);
0641 
0642                 /* We do not really know if up or down */
0643                 set_precision_flag_down();
0644 
0645                 return 0;
0646             } else {
0647                 /* Operand is out of range */
0648                 return 1;
0649             }
0650         } else {
0651               denormal_arg:
0652 
0653             setcc(0);
0654             FPU_copy_to_reg0(&CONST_1, TAG_Valid);
0655 #ifdef PECULIAR_486
0656             set_precision_flag_down();  /* 80486 appears to do this. */
0657 #else
0658             set_precision_flag_up();    /* Must be up. */
0659 #endif /* PECULIAR_486 */
0660             return 0;
0661         }
0662     } else if (tag == TAG_Zero) {
0663         FPU_copy_to_reg0(&CONST_1, TAG_Valid);
0664         setcc(0);
0665         return 0;
0666     }
0667 
0668     if (tag == TAG_Special)
0669         tag = FPU_Special(st0_ptr);
0670 
0671     if (tag == TW_Denormal) {
0672         if (denormal_operand() < 0)
0673             return 1;
0674 
0675         goto denormal_arg;
0676     } else if (tag == TW_Infinity) {
0677         /* The 80486 treats infinity as an invalid operand */
0678         arith_invalid(0);
0679         return 1;
0680     } else {
0681         single_arg_error(st0_ptr, tag); /* requires st0_ptr == &st(0) */
0682         return 1;
0683     }
0684 }
0685 
0686 static void fcos(FPU_REG *st0_ptr, u_char st0_tag)
0687 {
0688     f_cos(st0_ptr, st0_tag);
0689 }
0690 
0691 static void fsincos(FPU_REG *st0_ptr, u_char st0_tag)
0692 {
0693     FPU_REG *st_new_ptr;
0694     FPU_REG arg;
0695     u_char tag;
0696 
0697     /* Stack underflow has higher priority */
0698     if (st0_tag == TAG_Empty) {
0699         FPU_stack_underflow();  /* Puts a QNaN in st(0) */
0700         if (control_word & CW_Invalid) {
0701             st_new_ptr = &st(-1);
0702             push();
0703             FPU_stack_underflow();  /* Puts a QNaN in the new st(0) */
0704         }
0705         return;
0706     }
0707 
0708     if (STACK_OVERFLOW) {
0709         FPU_stack_overflow();
0710         return;
0711     }
0712 
0713     if (st0_tag == TAG_Special)
0714         tag = FPU_Special(st0_ptr);
0715     else
0716         tag = st0_tag;
0717 
0718     if (tag == TW_NaN) {
0719         single_arg_2_error(st0_ptr, TW_NaN);
0720         return;
0721     } else if (tag == TW_Infinity) {
0722         /* The 80486 treats infinity as an invalid operand */
0723         if (arith_invalid(0) >= 0) {
0724             /* Masked response */
0725             push();
0726             arith_invalid(0);
0727         }
0728         return;
0729     }
0730 
0731     reg_copy(st0_ptr, &arg);
0732     if (!f_sin(st0_ptr, st0_tag)) {
0733         push();
0734         FPU_copy_to_reg0(&arg, st0_tag);
0735         f_cos(&st(0), st0_tag);
0736     } else {
0737         /* An error, so restore st(0) */
0738         FPU_copy_to_reg0(&arg, st0_tag);
0739     }
0740 }
0741 
0742 /*---------------------------------------------------------------------------*/
0743 /* The following all require two arguments: st(0) and st(1) */
0744 
0745 /* A lean, mean kernel for the fprem instructions. This relies upon
0746    the division and rounding to an integer in do_fprem giving an
0747    exact result. Because of this, rem_kernel() needs to deal only with
0748    the least significant 64 bits, the more significant bits of the
0749    result must be zero.
0750  */
0751 static void rem_kernel(unsigned long long st0, unsigned long long *y,
0752                unsigned long long st1, unsigned long long q, int n)
0753 {
0754     int dummy;
0755     unsigned long long x;
0756 
0757     x = st0 << n;
0758 
0759     /* Do the required multiplication and subtraction in the one operation */
0760 
0761     /* lsw x -= lsw st1 * lsw q */
0762     asm volatile ("mull %4; subl %%eax,%0; sbbl %%edx,%1":"=m"
0763               (((unsigned *)&x)[0]), "=m"(((unsigned *)&x)[1]),
0764               "=a"(dummy)
0765               :"2"(((unsigned *)&st1)[0]), "m"(((unsigned *)&q)[0])
0766               :"%dx");
0767     /* msw x -= msw st1 * lsw q */
0768     asm volatile ("mull %3; subl %%eax,%0":"=m" (((unsigned *)&x)[1]),
0769               "=a"(dummy)
0770               :"1"(((unsigned *)&st1)[1]), "m"(((unsigned *)&q)[0])
0771               :"%dx");
0772     /* msw x -= lsw st1 * msw q */
0773     asm volatile ("mull %3; subl %%eax,%0":"=m" (((unsigned *)&x)[1]),
0774               "=a"(dummy)
0775               :"1"(((unsigned *)&st1)[0]), "m"(((unsigned *)&q)[1])
0776               :"%dx");
0777 
0778     *y = x;
0779 }
0780 
0781 /* Remainder of st(0) / st(1) */
0782 /* This routine produces exact results, i.e. there is never any
0783    rounding or truncation, etc of the result. */
0784 static void do_fprem(FPU_REG *st0_ptr, u_char st0_tag, int round)
0785 {
0786     FPU_REG *st1_ptr = &st(1);
0787     u_char st1_tag = FPU_gettagi(1);
0788 
0789     if (!((st0_tag ^ TAG_Valid) | (st1_tag ^ TAG_Valid))) {
0790         FPU_REG tmp, st0, st1;
0791         u_char st0_sign, st1_sign;
0792         u_char tmptag;
0793         int tag;
0794         int old_cw;
0795         int expdif;
0796         long long q;
0797         unsigned short saved_status;
0798         int cc;
0799 
0800           fprem_valid:
0801         /* Convert registers for internal use. */
0802         st0_sign = FPU_to_exp16(st0_ptr, &st0);
0803         st1_sign = FPU_to_exp16(st1_ptr, &st1);
0804         expdif = exponent16(&st0) - exponent16(&st1);
0805 
0806         old_cw = control_word;
0807         cc = 0;
0808 
0809         /* We want the status following the denorm tests, but don't want
0810            the status changed by the arithmetic operations. */
0811         saved_status = partial_status;
0812         control_word &= ~CW_RC;
0813         control_word |= RC_CHOP;
0814 
0815         if (expdif < 64) {
0816             /* This should be the most common case */
0817 
0818             if (expdif > -2) {
0819                 u_char sign = st0_sign ^ st1_sign;
0820                 tag = FPU_u_div(&st0, &st1, &tmp,
0821                         PR_64_BITS | RC_CHOP | 0x3f,
0822                         sign);
0823                 setsign(&tmp, sign);
0824 
0825                 if (exponent(&tmp) >= 0) {
0826                     FPU_round_to_int(&tmp, tag);    /* Fortunately, this can't
0827                                        overflow to 2^64 */
0828                     q = significand(&tmp);
0829 
0830                     rem_kernel(significand(&st0),
0831                            &significand(&tmp),
0832                            significand(&st1),
0833                            q, expdif);
0834 
0835                     setexponent16(&tmp, exponent16(&st1));
0836                 } else {
0837                     reg_copy(&st0, &tmp);
0838                     q = 0;
0839                 }
0840 
0841                 if ((round == RC_RND)
0842                     && (tmp.sigh & 0xc0000000)) {
0843                     /* We may need to subtract st(1) once more,
0844                        to get a result <= 1/2 of st(1). */
0845                     unsigned long long x;
0846                     expdif =
0847                         exponent16(&st1) - exponent16(&tmp);
0848                     if (expdif <= 1) {
0849                         if (expdif == 0)
0850                             x = significand(&st1) -
0851                                 significand(&tmp);
0852                         else    /* expdif is 1 */
0853                             x = (significand(&st1)
0854                                  << 1) -
0855                                 significand(&tmp);
0856                         if ((x < significand(&tmp)) ||
0857                             /* or equi-distant (from 0 & st(1)) and q is odd */
0858                             ((x == significand(&tmp))
0859                              && (q & 1))) {
0860                             st0_sign = !st0_sign;
0861                             significand(&tmp) = x;
0862                             q++;
0863                         }
0864                     }
0865                 }
0866 
0867                 if (q & 4)
0868                     cc |= SW_C0;
0869                 if (q & 2)
0870                     cc |= SW_C3;
0871                 if (q & 1)
0872                     cc |= SW_C1;
0873             } else {
0874                 control_word = old_cw;
0875                 setcc(0);
0876                 return;
0877             }
0878         } else {
0879             /* There is a large exponent difference ( >= 64 ) */
0880             /* To make much sense, the code in this section should
0881                be done at high precision. */
0882             int exp_1, N;
0883             u_char sign;
0884 
0885             /* prevent overflow here */
0886             /* N is 'a number between 32 and 63' (p26-113) */
0887             reg_copy(&st0, &tmp);
0888             tmptag = st0_tag;
0889             N = (expdif & 0x0000001f) + 32; /* This choice gives results
0890                                identical to an AMD 486 */
0891             setexponent16(&tmp, N);
0892             exp_1 = exponent16(&st1);
0893             setexponent16(&st1, 0);
0894             expdif -= N;
0895 
0896             sign = getsign(&tmp) ^ st1_sign;
0897             tag =
0898                 FPU_u_div(&tmp, &st1, &tmp,
0899                       PR_64_BITS | RC_CHOP | 0x3f, sign);
0900             setsign(&tmp, sign);
0901 
0902             FPU_round_to_int(&tmp, tag);    /* Fortunately, this can't
0903                                overflow to 2^64 */
0904 
0905             rem_kernel(significand(&st0),
0906                    &significand(&tmp),
0907                    significand(&st1),
0908                    significand(&tmp), exponent(&tmp)
0909                 );
0910             setexponent16(&tmp, exp_1 + expdif);
0911 
0912             /* It is possible for the operation to be complete here.
0913                What does the IEEE standard say? The Intel 80486 manual
0914                implies that the operation will never be completed at this
0915                point, and the behaviour of a real 80486 confirms this.
0916              */
0917             if (!(tmp.sigh | tmp.sigl)) {
0918                 /* The result is zero */
0919                 control_word = old_cw;
0920                 partial_status = saved_status;
0921                 FPU_copy_to_reg0(&CONST_Z, TAG_Zero);
0922                 setsign(&st0, st0_sign);
0923 #ifdef PECULIAR_486
0924                 setcc(SW_C2);
0925 #else
0926                 setcc(0);
0927 #endif /* PECULIAR_486 */
0928                 return;
0929             }
0930             cc = SW_C2;
0931         }
0932 
0933         control_word = old_cw;
0934         partial_status = saved_status;
0935         tag = FPU_normalize_nuo(&tmp);
0936         reg_copy(&tmp, st0_ptr);
0937 
0938         /* The only condition to be looked for is underflow,
0939            and it can occur here only if underflow is unmasked. */
0940         if ((exponent16(&tmp) <= EXP_UNDER) && (tag != TAG_Zero)
0941             && !(control_word & CW_Underflow)) {
0942             setcc(cc);
0943             tag = arith_underflow(st0_ptr);
0944             setsign(st0_ptr, st0_sign);
0945             FPU_settag0(tag);
0946             return;
0947         } else if ((exponent16(&tmp) > EXP_UNDER) || (tag == TAG_Zero)) {
0948             stdexp(st0_ptr);
0949             setsign(st0_ptr, st0_sign);
0950         } else {
0951             tag =
0952                 FPU_round(st0_ptr, 0, 0, FULL_PRECISION, st0_sign);
0953         }
0954         FPU_settag0(tag);
0955         setcc(cc);
0956 
0957         return;
0958     }
0959 
0960     if (st0_tag == TAG_Special)
0961         st0_tag = FPU_Special(st0_ptr);
0962     if (st1_tag == TAG_Special)
0963         st1_tag = FPU_Special(st1_ptr);
0964 
0965     if (((st0_tag == TAG_Valid) && (st1_tag == TW_Denormal))
0966         || ((st0_tag == TW_Denormal) && (st1_tag == TAG_Valid))
0967         || ((st0_tag == TW_Denormal) && (st1_tag == TW_Denormal))) {
0968         if (denormal_operand() < 0)
0969             return;
0970         goto fprem_valid;
0971     } else if ((st0_tag == TAG_Empty) || (st1_tag == TAG_Empty)) {
0972         FPU_stack_underflow();
0973         return;
0974     } else if (st0_tag == TAG_Zero) {
0975         if (st1_tag == TAG_Valid) {
0976             setcc(0);
0977             return;
0978         } else if (st1_tag == TW_Denormal) {
0979             if (denormal_operand() < 0)
0980                 return;
0981             setcc(0);
0982             return;
0983         } else if (st1_tag == TAG_Zero) {
0984             arith_invalid(0);
0985             return;
0986         } /* fprem(?,0) always invalid */
0987         else if (st1_tag == TW_Infinity) {
0988             setcc(0);
0989             return;
0990         }
0991     } else if ((st0_tag == TAG_Valid) || (st0_tag == TW_Denormal)) {
0992         if (st1_tag == TAG_Zero) {
0993             arith_invalid(0);   /* fprem(Valid,Zero) is invalid */
0994             return;
0995         } else if (st1_tag != TW_NaN) {
0996             if (((st0_tag == TW_Denormal)
0997                  || (st1_tag == TW_Denormal))
0998                 && (denormal_operand() < 0))
0999                 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 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.1.0 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!