1 /* $OpenBSD: sfdiv.c,v 1.4 2001/03/29 03:58:19 mickey Exp $ */ 2 3 /* 4 * Copyright 1996 1995 by Open Software Foundation, Inc. 5 * All Rights Reserved 6 * 7 * Permission to use, copy, modify, and distribute this software and 8 * its documentation for any purpose and without fee is hereby granted, 9 * provided that the above copyright notice appears in all copies and 10 * that both the copyright notice and this permission notice appear in 11 * supporting documentation. 12 * 13 * OSF DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE 14 * INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS 15 * FOR A PARTICULAR PURPOSE. 16 * 17 * IN NO EVENT SHALL OSF BE LIABLE FOR ANY SPECIAL, INDIRECT, OR 18 * CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM 19 * LOSS OF USE, DATA OR PROFITS, WHETHER IN ACTION OF CONTRACT, 20 * NEGLIGENCE, OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION 21 * WITH THE USE OR PERFORMANCE OF THIS SOFTWARE. 22 * 23 */ 24 /* 25 * pmk1.1 26 */ 27 /* 28 * (c) Copyright 1986 HEWLETT-PACKARD COMPANY 29 * 30 * To anyone who acknowledges that this file is provided "AS IS" 31 * without any express or implied warranty: 32 * permission to use, copy, modify, and distribute this file 33 * for any purpose is hereby granted without fee, provided that 34 * the above copyright notice and this notice appears in all 35 * copies, and that the name of Hewlett-Packard Company not be 36 * used in advertising or publicity pertaining to distribution 37 * of the software without specific, written prior permission. 38 * Hewlett-Packard Company makes no representations about the 39 * suitability of this software for any purpose. 40 */ 41 42 #include "../spmath/float.h" 43 #include "../spmath/sgl_float.h" 44 45 /* 46 * Single Precision Floating-point Divide 47 */ 48 int 49 sgl_fdiv(srcptr1,srcptr2,dstptr,status) 50 51 sgl_floating_point *srcptr1, *srcptr2, *dstptr; 52 unsigned int *status; 53 { 54 register unsigned int opnd1, opnd2, opnd3, result; 55 register int dest_exponent, count; 56 register int inexact = FALSE, guardbit = FALSE, stickybit = FALSE; 57 int is_tiny; 58 59 opnd1 = *srcptr1; 60 opnd2 = *srcptr2; 61 /* 62 * set sign bit of result 63 */ 64 if (Sgl_sign(opnd1) ^ Sgl_sign(opnd2)) Sgl_setnegativezero(result); 65 else Sgl_setzero(result); 66 /* 67 * check first operand for NaN's or infinity 68 */ 69 if (Sgl_isinfinity_exponent(opnd1)) { 70 if (Sgl_iszero_mantissa(opnd1)) { 71 if (Sgl_isnotnan(opnd2)) { 72 if (Sgl_isinfinity(opnd2)) { 73 /* 74 * invalid since both operands 75 * are infinity 76 */ 77 if (Is_invalidtrap_enabled()) 78 return(INVALIDEXCEPTION); 79 Set_invalidflag(); 80 Sgl_makequietnan(result); 81 *dstptr = result; 82 return(NOEXCEPTION); 83 } 84 /* 85 * return infinity 86 */ 87 Sgl_setinfinity_exponentmantissa(result); 88 *dstptr = result; 89 return(NOEXCEPTION); 90 } 91 } 92 else { 93 /* 94 * is NaN; signaling or quiet? 95 */ 96 if (Sgl_isone_signaling(opnd1)) { 97 /* trap if INVALIDTRAP enabled */ 98 if (Is_invalidtrap_enabled()) 99 return(INVALIDEXCEPTION); 100 /* make NaN quiet */ 101 Set_invalidflag(); 102 Sgl_set_quiet(opnd1); 103 } 104 /* 105 * is second operand a signaling NaN? 106 */ 107 else if (Sgl_is_signalingnan(opnd2)) { 108 /* trap if INVALIDTRAP enabled */ 109 if (Is_invalidtrap_enabled()) 110 return(INVALIDEXCEPTION); 111 /* make NaN quiet */ 112 Set_invalidflag(); 113 Sgl_set_quiet(opnd2); 114 *dstptr = opnd2; 115 return(NOEXCEPTION); 116 } 117 /* 118 * return quiet NaN 119 */ 120 *dstptr = opnd1; 121 return(NOEXCEPTION); 122 } 123 } 124 /* 125 * check second operand for NaN's or infinity 126 */ 127 if (Sgl_isinfinity_exponent(opnd2)) { 128 if (Sgl_iszero_mantissa(opnd2)) { 129 /* 130 * return zero 131 */ 132 Sgl_setzero_exponentmantissa(result); 133 *dstptr = result; 134 return(NOEXCEPTION); 135 } 136 /* 137 * is NaN; signaling or quiet? 138 */ 139 if (Sgl_isone_signaling(opnd2)) { 140 /* trap if INVALIDTRAP enabled */ 141 if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION); 142 /* make NaN quiet */ 143 Set_invalidflag(); 144 Sgl_set_quiet(opnd2); 145 } 146 /* 147 * return quiet NaN 148 */ 149 *dstptr = opnd2; 150 return(NOEXCEPTION); 151 } 152 /* 153 * check for division by zero 154 */ 155 if (Sgl_iszero_exponentmantissa(opnd2)) { 156 if (Sgl_iszero_exponentmantissa(opnd1)) { 157 /* invalid since both operands are zero */ 158 if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION); 159 Set_invalidflag(); 160 Sgl_makequietnan(result); 161 *dstptr = result; 162 return(NOEXCEPTION); 163 } 164 if (Is_divisionbyzerotrap_enabled()) 165 return(DIVISIONBYZEROEXCEPTION); 166 Set_divisionbyzeroflag(); 167 Sgl_setinfinity_exponentmantissa(result); 168 *dstptr = result; 169 return(NOEXCEPTION); 170 } 171 /* 172 * Generate exponent 173 */ 174 dest_exponent = Sgl_exponent(opnd1) - Sgl_exponent(opnd2) + SGL_BIAS; 175 176 /* 177 * Generate mantissa 178 */ 179 if (Sgl_isnotzero_exponent(opnd1)) { 180 /* set hidden bit */ 181 Sgl_clear_signexponent_set_hidden(opnd1); 182 } 183 else { 184 /* check for zero */ 185 if (Sgl_iszero_mantissa(opnd1)) { 186 Sgl_setzero_exponentmantissa(result); 187 *dstptr = result; 188 return(NOEXCEPTION); 189 } 190 /* is denormalized; want to normalize */ 191 Sgl_clear_signexponent(opnd1); 192 Sgl_leftshiftby1(opnd1); 193 Sgl_normalize(opnd1,dest_exponent); 194 } 195 /* opnd2 needs to have hidden bit set with msb in hidden bit */ 196 if (Sgl_isnotzero_exponent(opnd2)) { 197 Sgl_clear_signexponent_set_hidden(opnd2); 198 } 199 else { 200 /* is denormalized; want to normalize */ 201 Sgl_clear_signexponent(opnd2); 202 Sgl_leftshiftby1(opnd2); 203 while(Sgl_iszero_hiddenhigh7mantissa(opnd2)) { 204 Sgl_leftshiftby8(opnd2); 205 dest_exponent += 8; 206 } 207 if(Sgl_iszero_hiddenhigh3mantissa(opnd2)) { 208 Sgl_leftshiftby4(opnd2); 209 dest_exponent += 4; 210 } 211 while(Sgl_iszero_hidden(opnd2)) { 212 Sgl_leftshiftby1(opnd2); 213 dest_exponent += 1; 214 } 215 } 216 217 /* Divide the source mantissas */ 218 219 /* 220 * A non_restoring divide algorithm is used. 221 */ 222 Sgl_subtract(opnd1,opnd2,opnd1); 223 Sgl_setzero(opnd3); 224 for (count=1;count<=SGL_P && Sgl_all(opnd1);count++) { 225 Sgl_leftshiftby1(opnd1); 226 Sgl_leftshiftby1(opnd3); 227 if (Sgl_iszero_sign(opnd1)) { 228 Sgl_setone_lowmantissa(opnd3); 229 Sgl_subtract(opnd1,opnd2,opnd1); 230 } 231 else Sgl_addition(opnd1,opnd2,opnd1); 232 } 233 if (count <= SGL_P) { 234 Sgl_leftshiftby1(opnd3); 235 Sgl_setone_lowmantissa(opnd3); 236 Sgl_leftshift(opnd3,SGL_P-count); 237 if (Sgl_iszero_hidden(opnd3)) { 238 Sgl_leftshiftby1(opnd3); 239 dest_exponent--; 240 } 241 } 242 else { 243 if (Sgl_iszero_hidden(opnd3)) { 244 /* need to get one more bit of result */ 245 Sgl_leftshiftby1(opnd1); 246 Sgl_leftshiftby1(opnd3); 247 if (Sgl_iszero_sign(opnd1)) { 248 Sgl_setone_lowmantissa(opnd3); 249 Sgl_subtract(opnd1,opnd2,opnd1); 250 } 251 else Sgl_addition(opnd1,opnd2,opnd1); 252 dest_exponent--; 253 } 254 if (Sgl_iszero_sign(opnd1)) guardbit = TRUE; 255 stickybit = Sgl_all(opnd1); 256 } 257 inexact = guardbit | stickybit; 258 259 /* 260 * round result 261 */ 262 if (inexact && (dest_exponent > 0 || Is_underflowtrap_enabled())) { 263 Sgl_clear_signexponent(opnd3); 264 switch (Rounding_mode()) { 265 case ROUNDPLUS: 266 if (Sgl_iszero_sign(result)) 267 Sgl_increment_mantissa(opnd3); 268 break; 269 case ROUNDMINUS: 270 if (Sgl_isone_sign(result)) 271 Sgl_increment_mantissa(opnd3); 272 break; 273 case ROUNDNEAREST: 274 if (guardbit && 275 (stickybit || Sgl_isone_lowmantissa(opnd3))) 276 Sgl_increment_mantissa(opnd3); 277 } 278 if (Sgl_isone_hidden(opnd3)) dest_exponent++; 279 } 280 Sgl_set_mantissa(result,opnd3); 281 282 /* 283 * Test for overflow 284 */ 285 if (dest_exponent >= SGL_INFINITY_EXPONENT) { 286 /* trap if OVERFLOWTRAP enabled */ 287 if (Is_overflowtrap_enabled()) { 288 /* 289 * Adjust bias of result 290 */ 291 Sgl_setwrapped_exponent(result,dest_exponent,ovfl); 292 *dstptr = result; 293 if (inexact) { 294 if (Is_inexacttrap_enabled()) 295 return(OVERFLOWEXCEPTION | INEXACTEXCEPTION); 296 else Set_inexactflag(); 297 } 298 return(OVERFLOWEXCEPTION); 299 } 300 Set_overflowflag(); 301 /* set result to infinity or largest number */ 302 Sgl_setoverflow(result); 303 inexact = TRUE; 304 } 305 /* 306 * Test for underflow 307 */ 308 else if (dest_exponent <= 0) { 309 /* trap if UNDERFLOWTRAP enabled */ 310 if (Is_underflowtrap_enabled()) { 311 /* 312 * Adjust bias of result 313 */ 314 Sgl_setwrapped_exponent(result,dest_exponent,unfl); 315 *dstptr = result; 316 if (inexact) { 317 if (Is_inexacttrap_enabled()) 318 return(UNDERFLOWEXCEPTION | INEXACTEXCEPTION); 319 else Set_inexactflag(); 320 } 321 return(UNDERFLOWEXCEPTION); 322 } 323 324 /* Determine if should set underflow flag */ 325 is_tiny = TRUE; 326 if (dest_exponent == 0 && inexact) { 327 switch (Rounding_mode()) { 328 case ROUNDPLUS: 329 if (Sgl_iszero_sign(result)) { 330 Sgl_increment(opnd3); 331 if (Sgl_isone_hiddenoverflow(opnd3)) 332 is_tiny = FALSE; 333 Sgl_decrement(opnd3); 334 } 335 break; 336 case ROUNDMINUS: 337 if (Sgl_isone_sign(result)) { 338 Sgl_increment(opnd3); 339 if (Sgl_isone_hiddenoverflow(opnd3)) 340 is_tiny = FALSE; 341 Sgl_decrement(opnd3); 342 } 343 break; 344 case ROUNDNEAREST: 345 if (guardbit && (stickybit || 346 Sgl_isone_lowmantissa(opnd3))) { 347 Sgl_increment(opnd3); 348 if (Sgl_isone_hiddenoverflow(opnd3)) 349 is_tiny = FALSE; 350 Sgl_decrement(opnd3); 351 } 352 break; 353 } 354 } 355 356 /* 357 * denormalize result or set to signed zero 358 */ 359 stickybit = inexact; 360 Sgl_denormalize(opnd3,dest_exponent,guardbit,stickybit,inexact); 361 362 /* return rounded number */ 363 if (inexact) { 364 switch (Rounding_mode()) { 365 case ROUNDPLUS: 366 if (Sgl_iszero_sign(result)) { 367 Sgl_increment(opnd3); 368 } 369 break; 370 case ROUNDMINUS: 371 if (Sgl_isone_sign(result)) { 372 Sgl_increment(opnd3); 373 } 374 break; 375 case ROUNDNEAREST: 376 if (guardbit && (stickybit || 377 Sgl_isone_lowmantissa(opnd3))) { 378 Sgl_increment(opnd3); 379 } 380 break; 381 } 382 if (is_tiny) Set_underflowflag(); 383 } 384 Sgl_set_exponentmantissa(result,opnd3); 385 } 386 else Sgl_set_exponent(result,dest_exponent); 387 *dstptr = result; 388 /* check for inexact */ 389 if (inexact) { 390 if (Is_inexacttrap_enabled()) return(INEXACTEXCEPTION); 391 else Set_inexactflag(); 392 } 393 return(NOEXCEPTION); 394 } 395