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