![[BACK]](/icons/back.gif){kind=icon}
Return to dfdiv.c CVS log ![[TXT]](/icons/text.gif){kind=icon}
![[DIR]](/icons/dir.gif){kind=icon}
Up to [local] / sys / arch / hppa / spmath|
Return to dfdiv.c CVS log | Up to [local] / sys / arch / hppa / spmath |
1.1 ! nbrk 1: /* $OpenBSD: dfdiv.c,v 1.5 2002/05/07 22:19:30 mickey Exp $ */
! 2: /*
! 3: (c) Copyright 1986 HEWLETT-PACKARD COMPANY
! 4: To anyone who acknowledges that this file is provided "AS IS"
! 5: without any express or implied warranty:
! 6: permission to use, copy, modify, and distribute this file
! 7: for any purpose is hereby granted without fee, provided that
! 8: the above copyright notice and this notice appears in all
! 9: copies, and that the name of Hewlett-Packard Company not be
! 10: used in advertising or publicity pertaining to distribution
! 11: of the software without specific, written prior permission.
! 12: Hewlett-Packard Company makes no representations about the
! 13: suitability of this software for any purpose.
! 14: */
! 15: /* @(#)dfdiv.c: Revision: 2.11.88.1 Date: 93/12/07 15:05:39 */
! 16:
! 17: #include "float.h"
! 18: #include "dbl_float.h"
! 19:
! 20: /*
! 21: * Double Precision Floating-point Divide
! 22: */
! 23:
! 24: int
! 25: dbl_fdiv(srcptr1,srcptr2,dstptr,status)
! 26: dbl_floating_point *srcptr1, *srcptr2, *dstptr;
! 27: unsigned int *status;
! 28: {
! 29: register unsigned int opnd1p1, opnd1p2, opnd2p1, opnd2p2;
! 30: register unsigned int opnd3p1, opnd3p2, resultp1, resultp2;
! 31: register int dest_exponent, count;
! 32: register int inexact = FALSE, guardbit = FALSE, stickybit = FALSE;
! 33: int is_tiny;
! 34:
! 35: Dbl_copyfromptr(srcptr1,opnd1p1,opnd1p2);
! 36: Dbl_copyfromptr(srcptr2,opnd2p1,opnd2p2);
! 37: /*
! 38: * set sign bit of result
! 39: */
! 40: if (Dbl_sign(opnd1p1) ^ Dbl_sign(opnd2p1))
! 41: Dbl_setnegativezerop1(resultp1);
! 42: else Dbl_setzerop1(resultp1);
! 43: /*
! 44: * check first operand for NaN's or infinity
! 45: */
! 46: if (Dbl_isinfinity_exponent(opnd1p1)) {
! 47: if (Dbl_iszero_mantissa(opnd1p1,opnd1p2)) {
! 48: if (Dbl_isnotnan(opnd2p1,opnd2p2)) {
! 49: if (Dbl_isinfinity(opnd2p1,opnd2p2)) {
! 50: /*
! 51: * invalid since both operands
! 52: * are infinity
! 53: */
! 54: if (Is_invalidtrap_enabled())
! 55: return(INVALIDEXCEPTION);
! 56: Set_invalidflag();
! 57: Dbl_makequietnan(resultp1,resultp2);
! 58: Dbl_copytoptr(resultp1,resultp2,dstptr);
! 59: return(NOEXCEPTION);
! 60: }
! 61: /*
! 62: * return infinity
! 63: */
! 64: Dbl_setinfinity_exponentmantissa(resultp1,resultp2);
! 65: Dbl_copytoptr(resultp1,resultp2,dstptr);
! 66: return(NOEXCEPTION);
! 67: }
! 68: }
! 69: else {
! 70: /*
! 71: * is NaN; signaling or quiet?
! 72: */
! 73: if (Dbl_isone_signaling(opnd1p1)) {
! 74: /* trap if INVALIDTRAP enabled */
! 75: if (Is_invalidtrap_enabled())
! 76: return(INVALIDEXCEPTION);
! 77: /* make NaN quiet */
! 78: Set_invalidflag();
! 79: Dbl_set_quiet(opnd1p1);
! 80: }
! 81: /*
! 82: * is second operand a signaling NaN?
! 83: */
! 84: else if (Dbl_is_signalingnan(opnd2p1)) {
! 85: /* trap if INVALIDTRAP enabled */
! 86: if (Is_invalidtrap_enabled())
! 87: return(INVALIDEXCEPTION);
! 88: /* make NaN quiet */
! 89: Set_invalidflag();
! 90: Dbl_set_quiet(opnd2p1);
! 91: Dbl_copytoptr(opnd2p1,opnd2p2,dstptr);
! 92: return(NOEXCEPTION);
! 93: }
! 94: /*
! 95: * return quiet NaN
! 96: */
! 97: Dbl_copytoptr(opnd1p1,opnd1p2,dstptr);
! 98: return(NOEXCEPTION);
! 99: }
! 100: }
! 101: /*
! 102: * check second operand for NaN's or infinity
! 103: */
! 104: if (Dbl_isinfinity_exponent(opnd2p1)) {
! 105: if (Dbl_iszero_mantissa(opnd2p1,opnd2p2)) {
! 106: /*
! 107: * return zero
! 108: */
! 109: Dbl_setzero_exponentmantissa(resultp1,resultp2);
! 110: Dbl_copytoptr(resultp1,resultp2,dstptr);
! 111: return(NOEXCEPTION);
! 112: }
! 113: /*
! 114: * is NaN; signaling or quiet?
! 115: */
! 116: if (Dbl_isone_signaling(opnd2p1)) {
! 117: /* trap if INVALIDTRAP enabled */
! 118: if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
! 119: /* make NaN quiet */
! 120: Set_invalidflag();
! 121: Dbl_set_quiet(opnd2p1);
! 122: }
! 123: /*
! 124: * return quiet NaN
! 125: */
! 126: Dbl_copytoptr(opnd2p1,opnd2p2,dstptr);
! 127: return(NOEXCEPTION);
! 128: }
! 129: /*
! 130: * check for division by zero
! 131: */
! 132: if (Dbl_iszero_exponentmantissa(opnd2p1,opnd2p2)) {
! 133: if (Dbl_iszero_exponentmantissa(opnd1p1,opnd1p2)) {
! 134: /* invalid since both operands are zero */
! 135: if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
! 136: Set_invalidflag();
! 137: Dbl_makequietnan(resultp1,resultp2);
! 138: Dbl_copytoptr(resultp1,resultp2,dstptr);
! 139: return(NOEXCEPTION);
! 140: }
! 141: if (Is_divisionbyzerotrap_enabled())
! 142: return(DIVISIONBYZEROEXCEPTION);
! 143: Set_divisionbyzeroflag();
! 144: Dbl_setinfinity_exponentmantissa(resultp1,resultp2);
! 145: Dbl_copytoptr(resultp1,resultp2,dstptr);
! 146: return(NOEXCEPTION);
! 147: }
! 148: /*
! 149: * Generate exponent
! 150: */
! 151: dest_exponent = Dbl_exponent(opnd1p1) - Dbl_exponent(opnd2p1) + DBL_BIAS;
! 152:
! 153: /*
! 154: * Generate mantissa
! 155: */
! 156: if (Dbl_isnotzero_exponent(opnd1p1)) {
! 157: /* set hidden bit */
! 158: Dbl_clear_signexponent_set_hidden(opnd1p1);
! 159: }
! 160: else {
! 161: /* check for zero */
! 162: if (Dbl_iszero_mantissa(opnd1p1,opnd1p2)) {
! 163: Dbl_setzero_exponentmantissa(resultp1,resultp2);
! 164: Dbl_copytoptr(resultp1,resultp2,dstptr);
! 165: return(NOEXCEPTION);
! 166: }
! 167: /* is denormalized, want to normalize */
! 168: Dbl_clear_signexponent(opnd1p1);
! 169: Dbl_leftshiftby1(opnd1p1,opnd1p2);
! 170: Dbl_normalize(opnd1p1,opnd1p2,dest_exponent);
! 171: }
! 172: /* opnd2 needs to have hidden bit set with msb in hidden bit */
! 173: if (Dbl_isnotzero_exponent(opnd2p1)) {
! 174: Dbl_clear_signexponent_set_hidden(opnd2p1);
! 175: }
! 176: else {
! 177: /* is denormalized; want to normalize */
! 178: Dbl_clear_signexponent(opnd2p1);
! 179: Dbl_leftshiftby1(opnd2p1,opnd2p2);
! 180: while (Dbl_iszero_hiddenhigh7mantissa(opnd2p1)) {
! 181: dest_exponent+=8;
! 182: Dbl_leftshiftby8(opnd2p1,opnd2p2);
! 183: }
! 184: if (Dbl_iszero_hiddenhigh3mantissa(opnd2p1)) {
! 185: dest_exponent+=4;
! 186: Dbl_leftshiftby4(opnd2p1,opnd2p2);
! 187: }
! 188: while (Dbl_iszero_hidden(opnd2p1)) {
! 189: dest_exponent++;
! 190: Dbl_leftshiftby1(opnd2p1,opnd2p2);
! 191: }
! 192: }
! 193:
! 194: /* Divide the source mantissas */
! 195:
! 196: /*
! 197: * A non-restoring divide algorithm is used.
! 198: */
! 199: Twoword_subtract(opnd1p1,opnd1p2,opnd2p1,opnd2p2);
! 200: Dbl_setzero(opnd3p1,opnd3p2);
! 201: for (count=1; count <= DBL_P && (opnd1p1 || opnd1p2); count++) {
! 202: Dbl_leftshiftby1(opnd1p1,opnd1p2);
! 203: Dbl_leftshiftby1(opnd3p1,opnd3p2);
! 204: if (Dbl_iszero_sign(opnd1p1)) {
! 205: Dbl_setone_lowmantissap2(opnd3p2);
! 206: Twoword_subtract(opnd1p1,opnd1p2,opnd2p1,opnd2p2);
! 207: }
! 208: else {
! 209: Twoword_add(opnd1p1, opnd1p2, opnd2p1, opnd2p2);
! 210: }
! 211: }
! 212: if (count <= DBL_P) {
! 213: Dbl_leftshiftby1(opnd3p1,opnd3p2);
! 214: Dbl_setone_lowmantissap2(opnd3p2);
! 215: Dbl_leftshift(opnd3p1,opnd3p2,(DBL_P-count));
! 216: if (Dbl_iszero_hidden(opnd3p1)) {
! 217: Dbl_leftshiftby1(opnd3p1,opnd3p2);
! 218: dest_exponent--;
! 219: }
! 220: }
! 221: else {
! 222: if (Dbl_iszero_hidden(opnd3p1)) {
! 223: /* need to get one more bit of result */
! 224: Dbl_leftshiftby1(opnd1p1,opnd1p2);
! 225: Dbl_leftshiftby1(opnd3p1,opnd3p2);
! 226: if (Dbl_iszero_sign(opnd1p1)) {
! 227: Dbl_setone_lowmantissap2(opnd3p2);
! 228: Twoword_subtract(opnd1p1,opnd1p2,opnd2p1,opnd2p2);
! 229: }
! 230: else {
! 231: Twoword_add(opnd1p1,opnd1p2,opnd2p1,opnd2p2);
! 232: }
! 233: dest_exponent--;
! 234: }
! 235: if (Dbl_iszero_sign(opnd1p1)) guardbit = TRUE;
! 236: stickybit = Dbl_allp1(opnd1p1) || Dbl_allp2(opnd1p2);
! 237: }
! 238: inexact = guardbit | stickybit;
! 239:
! 240: /*
! 241: * round result
! 242: */
! 243: if (inexact && (dest_exponent > 0 || Is_underflowtrap_enabled())) {
! 244: Dbl_clear_signexponent(opnd3p1);
! 245: switch (Rounding_mode()) {
! 246: case ROUNDPLUS:
! 247: if (Dbl_iszero_sign(resultp1))
! 248: Dbl_increment(opnd3p1,opnd3p2);
! 249: break;
! 250: case ROUNDMINUS:
! 251: if (Dbl_isone_sign(resultp1))
! 252: Dbl_increment(opnd3p1,opnd3p2);
! 253: break;
! 254: case ROUNDNEAREST:
! 255: if (guardbit && (stickybit ||
! 256: Dbl_isone_lowmantissap2(opnd3p2))) {
! 257: Dbl_increment(opnd3p1,opnd3p2);
! 258: }
! 259: }
! 260: if (Dbl_isone_hidden(opnd3p1)) dest_exponent++;
! 261: }
! 262: Dbl_set_mantissa(resultp1,resultp2,opnd3p1,opnd3p2);
! 263:
! 264: /*
! 265: * Test for overflow
! 266: */
! 267: if (dest_exponent >= DBL_INFINITY_EXPONENT) {
! 268: /* trap if OVERFLOWTRAP enabled */
! 269: if (Is_overflowtrap_enabled()) {
! 270: /*
! 271: * Adjust bias of result
! 272: */
! 273: Dbl_setwrapped_exponent(resultp1,dest_exponent,ovfl);
! 274: Dbl_copytoptr(resultp1,resultp2,dstptr);
! 275: if (inexact) {
! 276: if (Is_inexacttrap_enabled())
! 277: return(OVERFLOWEXCEPTION | INEXACTEXCEPTION);
! 278: else
! 279: Set_inexactflag();
! 280: }
! 281: return(OVERFLOWEXCEPTION);
! 282: }
! 283: Set_overflowflag();
! 284: /* set result to infinity or largest number */
! 285: Dbl_setoverflow(resultp1,resultp2);
! 286: inexact = TRUE;
! 287: }
! 288: /*
! 289: * Test for underflow
! 290: */
! 291: else if (dest_exponent <= 0) {
! 292: /* trap if UNDERFLOWTRAP enabled */
! 293: if (Is_underflowtrap_enabled()) {
! 294: /*
! 295: * Adjust bias of result
! 296: */
! 297: Dbl_setwrapped_exponent(resultp1,dest_exponent,unfl);
! 298: Dbl_copytoptr(resultp1,resultp2,dstptr);
! 299: if (inexact) {
! 300: if (Is_inexacttrap_enabled())
! 301: return(UNDERFLOWEXCEPTION | INEXACTEXCEPTION);
! 302: else
! 303: Set_inexactflag();
! 304: }
! 305: return(UNDERFLOWEXCEPTION);
! 306: }
! 307:
! 308: /* Determine if should set underflow flag */
! 309: is_tiny = TRUE;
! 310: if (dest_exponent == 0 && inexact) {
! 311: switch (Rounding_mode()) {
! 312: case ROUNDPLUS:
! 313: if (Dbl_iszero_sign(resultp1)) {
! 314: Dbl_increment(opnd3p1,opnd3p2);
! 315: if (Dbl_isone_hiddenoverflow(opnd3p1))
! 316: is_tiny = FALSE;
! 317: Dbl_decrement(opnd3p1,opnd3p2);
! 318: }
! 319: break;
! 320: case ROUNDMINUS:
! 321: if (Dbl_isone_sign(resultp1)) {
! 322: Dbl_increment(opnd3p1,opnd3p2);
! 323: if (Dbl_isone_hiddenoverflow(opnd3p1))
! 324: is_tiny = FALSE;
! 325: Dbl_decrement(opnd3p1,opnd3p2);
! 326: }
! 327: break;
! 328: case ROUNDNEAREST:
! 329: if (guardbit && (stickybit ||
! 330: Dbl_isone_lowmantissap2(opnd3p2))) {
! 331: Dbl_increment(opnd3p1,opnd3p2);
! 332: if (Dbl_isone_hiddenoverflow(opnd3p1))
! 333: is_tiny = FALSE;
! 334: Dbl_decrement(opnd3p1,opnd3p2);
! 335: }
! 336: break;
! 337: }
! 338: }
! 339:
! 340: /*
! 341: * denormalize result or set to signed zero
! 342: */
! 343: stickybit = inexact;
! 344: Dbl_denormalize(opnd3p1,opnd3p2,dest_exponent,guardbit,
! 345: stickybit,inexact);
! 346:
! 347: /* return rounded number */
! 348: if (inexact) {
! 349: switch (Rounding_mode()) {
! 350: case ROUNDPLUS:
! 351: if (Dbl_iszero_sign(resultp1)) {
! 352: Dbl_increment(opnd3p1,opnd3p2);
! 353: }
! 354: break;
! 355: case ROUNDMINUS:
! 356: if (Dbl_isone_sign(resultp1)) {
! 357: Dbl_increment(opnd3p1,opnd3p2);
! 358: }
! 359: break;
! 360: case ROUNDNEAREST:
! 361: if (guardbit && (stickybit ||
! 362: Dbl_isone_lowmantissap2(opnd3p2))) {
! 363: Dbl_increment(opnd3p1,opnd3p2);
! 364: }
! 365: break;
! 366: }
! 367: if (is_tiny)
! 368: Set_underflowflag();
! 369: }
! 370: Dbl_set_exponentmantissa(resultp1,resultp2,opnd3p1,opnd3p2);
! 371: }
! 372: else Dbl_set_exponent(resultp1,dest_exponent);
! 373: Dbl_copytoptr(resultp1,resultp2,dstptr);
! 374:
! 375: /* check for inexact */
! 376: if (inexact) {
! 377: if (Is_inexacttrap_enabled()) return(INEXACTEXCEPTION);
! 378: else Set_inexactflag();
! 379: }
! 380: return(NOEXCEPTION);
! 381: }