arch/parisc/math-emu/dfdiv.c

   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
  36 int
  37 dbl_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 }