(root)/
gcc-13.2.0/
libgcc/
config/
libbid/
bid64_minmax.c
       1  /* Copyright (C) 2007-2023 Free Software Foundation, Inc.
       2  
       3  This file is part of GCC.
       4  
       5  GCC is free software; you can redistribute it and/or modify it under
       6  the terms of the GNU General Public License as published by the Free
       7  Software Foundation; either version 3, or (at your option) any later
       8  version.
       9  
      10  GCC is distributed in the hope that it will be useful, but WITHOUT ANY
      11  WARRANTY; without even the implied warranty of MERCHANTABILITY or
      12  FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
      13  for more details.
      14  
      15  Under Section 7 of GPL version 3, you are granted additional
      16  permissions described in the GCC Runtime Library Exception, version
      17  3.1, as published by the Free Software Foundation.
      18  
      19  You should have received a copy of the GNU General Public License and
      20  a copy of the GCC Runtime Library Exception along with this program;
      21  see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
      22  <http://www.gnu.org/licenses/>.  */
      23  
      24  #include "bid_internal.h"
      25  
      26  /*****************************************************************************
      27   *  BID64 minimum function - returns greater of two numbers
      28   *****************************************************************************/
      29  
      30  static const UINT64 mult_factor[16] = {
      31    1ull, 10ull, 100ull, 1000ull,
      32    10000ull, 100000ull, 1000000ull, 10000000ull,
      33    100000000ull, 1000000000ull, 10000000000ull, 100000000000ull,
      34    1000000000000ull, 10000000000000ull,
      35    100000000000000ull, 1000000000000000ull
      36  };
      37  
      38  #if DECIMAL_CALL_BY_REFERENCE
      39  void
      40  bid64_minnum (UINT64 * pres, UINT64 * px, UINT64 * py _EXC_FLAGS_PARAM) {
      41    UINT64 x = *px;
      42    UINT64 y = *py;
      43  #else
      44  UINT64
      45  bid64_minnum (UINT64 x, UINT64 y _EXC_FLAGS_PARAM) {
      46  #endif
      47  
      48    UINT64 res;
      49    int exp_x, exp_y;
      50    UINT64 sig_x, sig_y;
      51    UINT128 sig_n_prime;
      52    char x_is_zero = 0, y_is_zero = 0;
      53  
      54    // check for non-canonical x
      55    if ((x & MASK_NAN) == MASK_NAN) {	// x is NaN
      56      x = x & 0xfe03ffffffffffffull;	// clear G6-G12
      57      if ((x & 0x0003ffffffffffffull) > 999999999999999ull) {
      58        x = x & 0xfe00000000000000ull;	// clear G6-G12 and the payload bits
      59      }
      60    } else if ((x & MASK_INF) == MASK_INF) {	// check for Infinity
      61      x = x & (MASK_SIGN | MASK_INF);
      62    } else {	// x is not special
      63      // check for non-canonical values - treated as zero
      64      if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
      65        // if the steering bits are 11, then the exponent is G[0:w+1]
      66        if (((x & MASK_BINARY_SIG2) | MASK_BINARY_OR2) >
      67  	  9999999999999999ull) {
      68  	// non-canonical
      69  	x = (x & MASK_SIGN) | ((x & MASK_BINARY_EXPONENT2) << 2);
      70        }	// else canonical
      71      }	// else canonical 
      72    }
      73  
      74    // check for non-canonical y
      75    if ((y & MASK_NAN) == MASK_NAN) {	// y is NaN
      76      y = y & 0xfe03ffffffffffffull;	// clear G6-G12
      77      if ((y & 0x0003ffffffffffffull) > 999999999999999ull) {
      78        y = y & 0xfe00000000000000ull;	// clear G6-G12 and the payload bits
      79      }
      80    } else if ((y & MASK_INF) == MASK_INF) {	// check for Infinity
      81      y = y & (MASK_SIGN | MASK_INF);
      82    } else {	// y is not special
      83      // check for non-canonical values - treated as zero
      84      if ((y & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
      85        // if the steering bits are 11, then the exponent is G[0:w+1]
      86        if (((y & MASK_BINARY_SIG2) | MASK_BINARY_OR2) >
      87  	  9999999999999999ull) {
      88  	// non-canonical
      89  	y = (y & MASK_SIGN) | ((y & MASK_BINARY_EXPONENT2) << 2);
      90        }	// else canonical
      91      }	// else canonical
      92    }
      93  
      94    // NaN (CASE1)
      95    if ((x & MASK_NAN) == MASK_NAN) {	// x is NAN
      96      if ((x & MASK_SNAN) == MASK_SNAN) {	// x is SNaN
      97        // if x is SNAN, then return quiet (x)
      98        *pfpsf |= INVALID_EXCEPTION;	// set exception if SNaN
      99        x = x & 0xfdffffffffffffffull;	// quietize x
     100        res = x;
     101      } else {	// x is QNaN
     102        if ((y & MASK_NAN) == MASK_NAN) {	// y is NAN
     103  	if ((y & MASK_SNAN) == MASK_SNAN) {	// y is SNAN
     104  	  *pfpsf |= INVALID_EXCEPTION;	// set invalid flag
     105  	}
     106  	res = x;
     107        } else {
     108  	res = y;
     109        }
     110      }
     111      BID_RETURN (res);
     112    } else if ((y & MASK_NAN) == MASK_NAN) {	// y is NaN, but x is not
     113      if ((y & MASK_SNAN) == MASK_SNAN) {
     114        *pfpsf |= INVALID_EXCEPTION;	// set exception if SNaN
     115        y = y & 0xfdffffffffffffffull;	// quietize y
     116        res = y;
     117      } else {
     118        // will return x (which is not NaN)
     119        res = x;
     120      }
     121      BID_RETURN (res);
     122    }
     123    // SIMPLE (CASE2)
     124    // if all the bits are the same, these numbers are equal, return either number
     125    if (x == y) {
     126      res = x;
     127      BID_RETURN (res);
     128    }
     129    // INFINITY (CASE3)
     130    if ((x & MASK_INF) == MASK_INF) {
     131      // if x is neg infinity, there is no way it is greater than y, return x
     132      if (((x & MASK_SIGN) == MASK_SIGN)) {
     133        res = x;
     134        BID_RETURN (res);
     135      }
     136      // x is pos infinity, return y
     137      else {
     138        res = y;
     139        BID_RETURN (res);
     140      }
     141    } else if ((y & MASK_INF) == MASK_INF) {
     142      // x is finite, so if y is positive infinity, then x is less, return y
     143      //                 if y is negative infinity, then x is greater, return x
     144      res = ((y & MASK_SIGN) == MASK_SIGN) ? y : x;
     145      BID_RETURN (res);
     146    }
     147    // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
     148    if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     149      exp_x = (x & MASK_BINARY_EXPONENT2) >> 51;
     150      sig_x = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
     151    } else {
     152      exp_x = (x & MASK_BINARY_EXPONENT1) >> 53;
     153      sig_x = (x & MASK_BINARY_SIG1);
     154    }
     155  
     156    // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
     157    if ((y & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     158      exp_y = (y & MASK_BINARY_EXPONENT2) >> 51;
     159      sig_y = (y & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
     160    } else {
     161      exp_y = (y & MASK_BINARY_EXPONENT1) >> 53;
     162      sig_y = (y & MASK_BINARY_SIG1);
     163    }
     164  
     165    // ZERO (CASE4)
     166    // some properties:
     167    //    (+ZERO == -ZERO) => therefore 
     168    //        ignore the sign, and neither number is greater
     169    //    (ZERO x 10^A == ZERO x 10^B) for any valid A, B => 
     170    //        ignore the exponent field
     171    //    (Any non-canonical # is considered 0)
     172    if (sig_x == 0) {
     173      x_is_zero = 1;
     174    }
     175    if (sig_y == 0) {
     176      y_is_zero = 1;
     177    }
     178  
     179    if (x_is_zero && y_is_zero) {
     180      // if both numbers are zero, neither is greater => return either
     181      res = y;
     182      BID_RETURN (res);
     183    } else if (x_is_zero) {
     184      // is x is zero, it is greater if Y is negative
     185      res = ((y & MASK_SIGN) == MASK_SIGN) ? y : x;
     186      BID_RETURN (res);
     187    } else if (y_is_zero) {
     188      // is y is zero, X is greater if it is positive
     189      res = ((x & MASK_SIGN) != MASK_SIGN) ? y : x;;
     190      BID_RETURN (res);
     191    }
     192    // OPPOSITE SIGN (CASE5)
     193    // now, if the sign bits differ, x is greater if y is negative
     194    if (((x ^ y) & MASK_SIGN) == MASK_SIGN) {
     195      res = ((y & MASK_SIGN) == MASK_SIGN) ? y : x;
     196      BID_RETURN (res);
     197    }
     198    // REDUNDANT REPRESENTATIONS (CASE6)
     199  
     200    // if both components are either bigger or smaller, 
     201    // it is clear what needs to be done
     202    if (sig_x > sig_y && exp_x >= exp_y) {
     203      res = ((x & MASK_SIGN) != MASK_SIGN) ? y : x;
     204      BID_RETURN (res);
     205    }
     206    if (sig_x < sig_y && exp_x <= exp_y) {
     207      res = ((x & MASK_SIGN) == MASK_SIGN) ? y : x;
     208      BID_RETURN (res);
     209    }
     210    // if exp_x is 15 greater than exp_y, no need for compensation
     211    if (exp_x - exp_y > 15) {
     212      res = ((x & MASK_SIGN) != MASK_SIGN) ? y : x;	// difference cannot be >10^15
     213      BID_RETURN (res);
     214    }
     215    // if exp_x is 15 less than exp_y, no need for compensation
     216    if (exp_y - exp_x > 15) {
     217      res = ((x & MASK_SIGN) == MASK_SIGN) ? y : x;
     218      BID_RETURN (res);
     219    }
     220    // if |exp_x - exp_y| < 15, it comes down to the compensated significand
     221    if (exp_x > exp_y) {	// to simplify the loop below,
     222  
     223      // otherwise adjust the x significand upwards
     224      __mul_64x64_to_128MACH (sig_n_prime, sig_x,
     225  			    mult_factor[exp_x - exp_y]);
     226      // if postitive, return whichever significand is larger 
     227      // (converse if negative)
     228      if (sig_n_prime.w[1] == 0 && (sig_n_prime.w[0] == sig_y)) {
     229        res = y;
     230        BID_RETURN (res);
     231      }
     232  
     233      res = (((sig_n_prime.w[1] > 0)
     234  	    || sig_n_prime.w[0] > sig_y) ^ ((x & MASK_SIGN) ==
     235  					    MASK_SIGN)) ? y : x;
     236      BID_RETURN (res);
     237    }
     238    // adjust the y significand upwards
     239    __mul_64x64_to_128MACH (sig_n_prime, sig_y,
     240  			  mult_factor[exp_y - exp_x]);
     241  
     242    // if postitive, return whichever significand is larger (converse if negative)
     243    if (sig_n_prime.w[1] == 0 && (sig_n_prime.w[0] == sig_x)) {
     244      res = y;
     245      BID_RETURN (res);
     246    }
     247    res = (((sig_n_prime.w[1] == 0)
     248  	  && (sig_x > sig_n_prime.w[0])) ^ ((x & MASK_SIGN) ==
     249  					    MASK_SIGN)) ? y : x;
     250    BID_RETURN (res);
     251  }
     252  
     253  /*****************************************************************************
     254   *  BID64 minimum magnitude function - returns greater of two numbers
     255   *****************************************************************************/
     256  
     257  #if DECIMAL_CALL_BY_REFERENCE
     258  void
     259  bid64_minnum_mag (UINT64 * pres, UINT64 * px,
     260  		  UINT64 * py _EXC_FLAGS_PARAM) {
     261    UINT64 x = *px;
     262    UINT64 y = *py;
     263  #else
     264  UINT64
     265  bid64_minnum_mag (UINT64 x, UINT64 y _EXC_FLAGS_PARAM) {
     266  #endif
     267  
     268    UINT64 res;
     269    int exp_x, exp_y;
     270    UINT64 sig_x, sig_y;
     271    UINT128 sig_n_prime;
     272  
     273    // check for non-canonical x
     274    if ((x & MASK_NAN) == MASK_NAN) {	// x is NaN
     275      x = x & 0xfe03ffffffffffffull;	// clear G6-G12
     276      if ((x & 0x0003ffffffffffffull) > 999999999999999ull) {
     277        x = x & 0xfe00000000000000ull;	// clear G6-G12 and the payload bits
     278      }
     279    } else if ((x & MASK_INF) == MASK_INF) {	// check for Infinity
     280      x = x & (MASK_SIGN | MASK_INF);
     281    } else {	// x is not special
     282      // check for non-canonical values - treated as zero
     283      if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     284        // if the steering bits are 11, then the exponent is G[0:w+1]
     285        if (((x & MASK_BINARY_SIG2) | MASK_BINARY_OR2) >
     286  	  9999999999999999ull) {
     287  	// non-canonical
     288  	x = (x & MASK_SIGN) | ((x & MASK_BINARY_EXPONENT2) << 2);
     289        }	// else canonical
     290      }	// else canonical 
     291    }
     292  
     293    // check for non-canonical y
     294    if ((y & MASK_NAN) == MASK_NAN) {	// y is NaN
     295      y = y & 0xfe03ffffffffffffull;	// clear G6-G12
     296      if ((y & 0x0003ffffffffffffull) > 999999999999999ull) {
     297        y = y & 0xfe00000000000000ull;	// clear G6-G12 and the payload bits
     298      }
     299    } else if ((y & MASK_INF) == MASK_INF) {	// check for Infinity
     300      y = y & (MASK_SIGN | MASK_INF);
     301    } else {	// y is not special
     302      // check for non-canonical values - treated as zero
     303      if ((y & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     304        // if the steering bits are 11, then the exponent is G[0:w+1]
     305        if (((y & MASK_BINARY_SIG2) | MASK_BINARY_OR2) >
     306  	  9999999999999999ull) {
     307  	// non-canonical
     308  	y = (y & MASK_SIGN) | ((y & MASK_BINARY_EXPONENT2) << 2);
     309        }	// else canonical
     310      }	// else canonical
     311    }
     312  
     313    // NaN (CASE1)
     314    if ((x & MASK_NAN) == MASK_NAN) {	// x is NAN
     315      if ((x & MASK_SNAN) == MASK_SNAN) {	// x is SNaN
     316        // if x is SNAN, then return quiet (x)
     317        *pfpsf |= INVALID_EXCEPTION;	// set exception if SNaN
     318        x = x & 0xfdffffffffffffffull;	// quietize x
     319        res = x;
     320      } else {	// x is QNaN
     321        if ((y & MASK_NAN) == MASK_NAN) {	// y is NAN
     322  	if ((y & MASK_SNAN) == MASK_SNAN) {	// y is SNAN
     323  	  *pfpsf |= INVALID_EXCEPTION;	// set invalid flag
     324  	}
     325  	res = x;
     326        } else {
     327  	res = y;
     328        }
     329      }
     330      BID_RETURN (res);
     331    } else if ((y & MASK_NAN) == MASK_NAN) {	// y is NaN, but x is not
     332      if ((y & MASK_SNAN) == MASK_SNAN) {
     333        *pfpsf |= INVALID_EXCEPTION;	// set exception if SNaN
     334        y = y & 0xfdffffffffffffffull;	// quietize y
     335        res = y;
     336      } else {
     337        // will return x (which is not NaN)
     338        res = x;
     339      }
     340      BID_RETURN (res);
     341    }
     342    // SIMPLE (CASE2)
     343    // if all the bits are the same, these numbers are equal, return either number
     344    if (x == y) {
     345      res = x;
     346      BID_RETURN (res);
     347    }
     348    // INFINITY (CASE3)
     349    if ((x & MASK_INF) == MASK_INF) {
     350      // x is infinity, its magnitude is greater than or equal to y
     351      // return x only if y is infinity and x is negative
     352      res = ((x & MASK_SIGN) == MASK_SIGN
     353  	   && (y & MASK_INF) == MASK_INF) ? x : y;
     354      BID_RETURN (res);
     355    } else if ((y & MASK_INF) == MASK_INF) {
     356      // y is infinity, then it must be greater in magnitude, return x
     357      res = x;
     358      BID_RETURN (res);
     359    }
     360    // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
     361    if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     362      exp_x = (x & MASK_BINARY_EXPONENT2) >> 51;
     363      sig_x = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
     364    } else {
     365      exp_x = (x & MASK_BINARY_EXPONENT1) >> 53;
     366      sig_x = (x & MASK_BINARY_SIG1);
     367    }
     368  
     369    // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
     370    if ((y & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     371      exp_y = (y & MASK_BINARY_EXPONENT2) >> 51;
     372      sig_y = (y & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
     373    } else {
     374      exp_y = (y & MASK_BINARY_EXPONENT1) >> 53;
     375      sig_y = (y & MASK_BINARY_SIG1);
     376    }
     377  
     378    // ZERO (CASE4)
     379    // some properties:
     380    //    (+ZERO == -ZERO) => therefore 
     381    //        ignore the sign, and neither number is greater
     382    //    (ZERO x 10^A == ZERO x 10^B) for any valid A, B => 
     383    //        ignore the exponent field
     384    //    (Any non-canonical # is considered 0)
     385    if (sig_x == 0) {
     386      res = x;	// x_is_zero, its magnitude must be smaller than y
     387      BID_RETURN (res);
     388    }
     389    if (sig_y == 0) {
     390      res = y;	// y_is_zero, its magnitude must be smaller than x
     391      BID_RETURN (res);
     392    }
     393    // REDUNDANT REPRESENTATIONS (CASE6)
     394    // if both components are either bigger or smaller, 
     395    // it is clear what needs to be done
     396    if (sig_x > sig_y && exp_x >= exp_y) {
     397      res = y;
     398      BID_RETURN (res);
     399    }
     400    if (sig_x < sig_y && exp_x <= exp_y) {
     401      res = x;
     402      BID_RETURN (res);
     403    }
     404    // if exp_x is 15 greater than exp_y, no need for compensation
     405    if (exp_x - exp_y > 15) {
     406      res = y;	// difference cannot be greater than 10^15
     407      BID_RETURN (res);
     408    }
     409    // if exp_x is 15 less than exp_y, no need for compensation
     410    if (exp_y - exp_x > 15) {
     411      res = x;
     412      BID_RETURN (res);
     413    }
     414    // if |exp_x - exp_y| < 15, it comes down to the compensated significand
     415    if (exp_x > exp_y) {	// to simplify the loop below,
     416      // otherwise adjust the x significand upwards
     417      __mul_64x64_to_128MACH (sig_n_prime, sig_x,
     418  			    mult_factor[exp_x - exp_y]);
     419      // now, sig_n_prime has: sig_x * 10^(exp_x-exp_y), this is 
     420      // the compensated signif.
     421      if (sig_n_prime.w[1] == 0 && (sig_n_prime.w[0] == sig_y)) {
     422        // two numbers are equal, return minNum(x,y)
     423        res = ((y & MASK_SIGN) == MASK_SIGN) ? y : x;
     424        BID_RETURN (res);
     425      }
     426      // now, if compensated_x (sig_n_prime) is greater than y, return y,  
     427      // otherwise return x
     428      res = ((sig_n_prime.w[1] != 0) || sig_n_prime.w[0] > sig_y) ? y : x;
     429      BID_RETURN (res);
     430    }
     431    // exp_y must be greater than exp_x, thus adjust the y significand upwards
     432    __mul_64x64_to_128MACH (sig_n_prime, sig_y,
     433  			  mult_factor[exp_y - exp_x]);
     434  
     435    if (sig_n_prime.w[1] == 0 && (sig_n_prime.w[0] == sig_x)) {
     436      res = ((y & MASK_SIGN) == MASK_SIGN) ? y : x;
     437      // two numbers are equal, return either
     438      BID_RETURN (res);
     439    }
     440  
     441    res = ((sig_n_prime.w[1] == 0) && (sig_x > sig_n_prime.w[0])) ? y : x;
     442    BID_RETURN (res);
     443  }
     444  
     445  /*****************************************************************************
     446   *  BID64 maximum function - returns greater of two numbers
     447   *****************************************************************************/
     448  
     449  #if DECIMAL_CALL_BY_REFERENCE
     450  void
     451  bid64_maxnum (UINT64 * pres, UINT64 * px, UINT64 * py _EXC_FLAGS_PARAM) {
     452    UINT64 x = *px;
     453    UINT64 y = *py;
     454  #else
     455  UINT64
     456  bid64_maxnum (UINT64 x, UINT64 y _EXC_FLAGS_PARAM) {
     457  #endif
     458  
     459    UINT64 res;
     460    int exp_x, exp_y;
     461    UINT64 sig_x, sig_y;
     462    UINT128 sig_n_prime;
     463    char x_is_zero = 0, y_is_zero = 0;
     464  
     465    // check for non-canonical x
     466    if ((x & MASK_NAN) == MASK_NAN) {	// x is NaN
     467      x = x & 0xfe03ffffffffffffull;	// clear G6-G12
     468      if ((x & 0x0003ffffffffffffull) > 999999999999999ull) {
     469        x = x & 0xfe00000000000000ull;	// clear G6-G12 and the payload bits
     470      }
     471    } else if ((x & MASK_INF) == MASK_INF) {	// check for Infinity
     472      x = x & (MASK_SIGN | MASK_INF);
     473    } else {	// x is not special
     474      // check for non-canonical values - treated as zero
     475      if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     476        // if the steering bits are 11, then the exponent is G[0:w+1]
     477        if (((x & MASK_BINARY_SIG2) | MASK_BINARY_OR2) >
     478  	  9999999999999999ull) {
     479  	// non-canonical
     480  	x = (x & MASK_SIGN) | ((x & MASK_BINARY_EXPONENT2) << 2);
     481        }	// else canonical
     482      }	// else canonical 
     483    }
     484  
     485    // check for non-canonical y
     486    if ((y & MASK_NAN) == MASK_NAN) {	// y is NaN
     487      y = y & 0xfe03ffffffffffffull;	// clear G6-G12
     488      if ((y & 0x0003ffffffffffffull) > 999999999999999ull) {
     489        y = y & 0xfe00000000000000ull;	// clear G6-G12 and the payload bits
     490      }
     491    } else if ((y & MASK_INF) == MASK_INF) {	// check for Infinity
     492      y = y & (MASK_SIGN | MASK_INF);
     493    } else {	// y is not special
     494      // check for non-canonical values - treated as zero
     495      if ((y & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     496        // if the steering bits are 11, then the exponent is G[0:w+1]
     497        if (((y & MASK_BINARY_SIG2) | MASK_BINARY_OR2) >
     498  	  9999999999999999ull) {
     499  	// non-canonical
     500  	y = (y & MASK_SIGN) | ((y & MASK_BINARY_EXPONENT2) << 2);
     501        }	// else canonical
     502      }	// else canonical
     503    }
     504  
     505    // NaN (CASE1)
     506    if ((x & MASK_NAN) == MASK_NAN) {	// x is NAN
     507      if ((x & MASK_SNAN) == MASK_SNAN) {	// x is SNaN
     508        // if x is SNAN, then return quiet (x)
     509        *pfpsf |= INVALID_EXCEPTION;	// set exception if SNaN
     510        x = x & 0xfdffffffffffffffull;	// quietize x
     511        res = x;
     512      } else {	// x is QNaN
     513        if ((y & MASK_NAN) == MASK_NAN) {	// y is NAN
     514  	if ((y & MASK_SNAN) == MASK_SNAN) {	// y is SNAN
     515  	  *pfpsf |= INVALID_EXCEPTION;	// set invalid flag
     516  	}
     517  	res = x;
     518        } else {
     519  	res = y;
     520        }
     521      }
     522      BID_RETURN (res);
     523    } else if ((y & MASK_NAN) == MASK_NAN) {	// y is NaN, but x is not
     524      if ((y & MASK_SNAN) == MASK_SNAN) {
     525        *pfpsf |= INVALID_EXCEPTION;	// set exception if SNaN
     526        y = y & 0xfdffffffffffffffull;	// quietize y
     527        res = y;
     528      } else {
     529        // will return x (which is not NaN)
     530        res = x;
     531      }
     532      BID_RETURN (res);
     533    }
     534    // SIMPLE (CASE2)
     535    // if all the bits are the same, these numbers are equal (not Greater).
     536    if (x == y) {
     537      res = x;
     538      BID_RETURN (res);
     539    }
     540    // INFINITY (CASE3)
     541    if ((x & MASK_INF) == MASK_INF) {
     542      // if x is neg infinity, there is no way it is greater than y, return y
     543      // x is pos infinity, it is greater, unless y is positive infinity => 
     544      // return y!=pos_infinity
     545      if (((x & MASK_SIGN) == MASK_SIGN)) {
     546        res = y;
     547        BID_RETURN (res);
     548      } else {
     549        res = (((y & MASK_INF) != MASK_INF)
     550  	     || ((y & MASK_SIGN) == MASK_SIGN)) ? x : y;
     551        BID_RETURN (res);
     552      }
     553    } else if ((y & MASK_INF) == MASK_INF) {
     554      // x is finite, so if y is positive infinity, then x is less, return y
     555      //                 if y is negative infinity, then x is greater, return x
     556      res = ((y & MASK_SIGN) == MASK_SIGN) ? x : y;
     557      BID_RETURN (res);
     558    }
     559    // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
     560    if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     561      exp_x = (x & MASK_BINARY_EXPONENT2) >> 51;
     562      sig_x = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
     563    } else {
     564      exp_x = (x & MASK_BINARY_EXPONENT1) >> 53;
     565      sig_x = (x & MASK_BINARY_SIG1);
     566    }
     567  
     568    // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
     569    if ((y & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     570      exp_y = (y & MASK_BINARY_EXPONENT2) >> 51;
     571      sig_y = (y & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
     572    } else {
     573      exp_y = (y & MASK_BINARY_EXPONENT1) >> 53;
     574      sig_y = (y & MASK_BINARY_SIG1);
     575    }
     576  
     577    // ZERO (CASE4)
     578    // some properties:
     579    //    (+ZERO == -ZERO) => therefore 
     580    //        ignore the sign, and neither number is greater
     581    //    (ZERO x 10^A == ZERO x 10^B) for any valid A, B => 
     582    //        ignore the exponent field
     583    //    (Any non-canonical # is considered 0)
     584    if (sig_x == 0) {
     585      x_is_zero = 1;
     586    }
     587    if (sig_y == 0) {
     588      y_is_zero = 1;
     589    }
     590  
     591    if (x_is_zero && y_is_zero) {
     592      // if both numbers are zero, neither is greater => return NOTGREATERTHAN
     593      res = y;
     594      BID_RETURN (res);
     595    } else if (x_is_zero) {
     596      // is x is zero, it is greater if Y is negative
     597      res = ((y & MASK_SIGN) == MASK_SIGN) ? x : y;
     598      BID_RETURN (res);
     599    } else if (y_is_zero) {
     600      // is y is zero, X is greater if it is positive
     601      res = ((x & MASK_SIGN) != MASK_SIGN) ? x : y;;
     602      BID_RETURN (res);
     603    }
     604    // OPPOSITE SIGN (CASE5)
     605    // now, if the sign bits differ, x is greater if y is negative
     606    if (((x ^ y) & MASK_SIGN) == MASK_SIGN) {
     607      res = ((y & MASK_SIGN) == MASK_SIGN) ? x : y;
     608      BID_RETURN (res);
     609    }
     610    // REDUNDANT REPRESENTATIONS (CASE6)
     611  
     612    // if both components are either bigger or smaller, 
     613    //     it is clear what needs to be done
     614    if (sig_x > sig_y && exp_x >= exp_y) {
     615      res = ((x & MASK_SIGN) != MASK_SIGN) ? x : y;
     616      BID_RETURN (res);
     617    }
     618    if (sig_x < sig_y && exp_x <= exp_y) {
     619      res = ((x & MASK_SIGN) == MASK_SIGN) ? x : y;
     620      BID_RETURN (res);
     621    }
     622    // if exp_x is 15 greater than exp_y, no need for compensation
     623    if (exp_x - exp_y > 15) {
     624      res = ((x & MASK_SIGN) != MASK_SIGN) ? x : y;
     625      // difference cannot be > 10^15
     626      BID_RETURN (res);
     627    }
     628    // if exp_x is 15 less than exp_y, no need for compensation
     629    if (exp_y - exp_x > 15) {
     630      res = ((x & MASK_SIGN) == MASK_SIGN) ? x : y;
     631      BID_RETURN (res);
     632    }
     633    // if |exp_x - exp_y| < 15, it comes down to the compensated significand
     634    if (exp_x > exp_y) {	// to simplify the loop below,
     635      // otherwise adjust the x significand upwards
     636      __mul_64x64_to_128MACH (sig_n_prime, sig_x,
     637  			    mult_factor[exp_x - exp_y]);
     638      // if postitive, return whichever significand is larger 
     639      // (converse if negative)
     640      if (sig_n_prime.w[1] == 0 && (sig_n_prime.w[0] == sig_y)) {
     641        res = y;
     642        BID_RETURN (res);
     643      }
     644      res = (((sig_n_prime.w[1] > 0)
     645  	    || sig_n_prime.w[0] > sig_y) ^ ((x & MASK_SIGN) ==
     646  					    MASK_SIGN)) ? x : y;
     647      BID_RETURN (res);
     648    }
     649    // adjust the y significand upwards
     650    __mul_64x64_to_128MACH (sig_n_prime, sig_y,
     651  			  mult_factor[exp_y - exp_x]);
     652  
     653    // if postitive, return whichever significand is larger (converse if negative)
     654    if (sig_n_prime.w[1] == 0 && (sig_n_prime.w[0] == sig_x)) {
     655      res = y;
     656      BID_RETURN (res);
     657    }
     658    res = (((sig_n_prime.w[1] == 0)
     659  	  && (sig_x > sig_n_prime.w[0])) ^ ((x & MASK_SIGN) ==
     660  					    MASK_SIGN)) ? x : y;
     661    BID_RETURN (res);
     662  }
     663  
     664  /*****************************************************************************
     665   *  BID64 maximum magnitude function - returns greater of two numbers
     666   *****************************************************************************/
     667  
     668  #if DECIMAL_CALL_BY_REFERENCE
     669  void
     670  bid64_maxnum_mag (UINT64 * pres, UINT64 * px,
     671  		  UINT64 * py _EXC_FLAGS_PARAM) {
     672    UINT64 x = *px;
     673    UINT64 y = *py;
     674  #else
     675  UINT64
     676  bid64_maxnum_mag (UINT64 x, UINT64 y _EXC_FLAGS_PARAM) {
     677  #endif
     678  
     679    UINT64 res;
     680    int exp_x, exp_y;
     681    UINT64 sig_x, sig_y;
     682    UINT128 sig_n_prime;
     683  
     684    // check for non-canonical x
     685    if ((x & MASK_NAN) == MASK_NAN) {	// x is NaN
     686      x = x & 0xfe03ffffffffffffull;	// clear G6-G12
     687      if ((x & 0x0003ffffffffffffull) > 999999999999999ull) {
     688        x = x & 0xfe00000000000000ull;	// clear G6-G12 and the payload bits
     689      }
     690    } else if ((x & MASK_INF) == MASK_INF) {	// check for Infinity
     691      x = x & (MASK_SIGN | MASK_INF);
     692    } else {	// x is not special
     693      // check for non-canonical values - treated as zero
     694      if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     695        // if the steering bits are 11, then the exponent is G[0:w+1]
     696        if (((x & MASK_BINARY_SIG2) | MASK_BINARY_OR2) >
     697  	  9999999999999999ull) {
     698  	// non-canonical
     699  	x = (x & MASK_SIGN) | ((x & MASK_BINARY_EXPONENT2) << 2);
     700        }	// else canonical
     701      }	// else canonical 
     702    }
     703  
     704    // check for non-canonical y
     705    if ((y & MASK_NAN) == MASK_NAN) {	// y is NaN
     706      y = y & 0xfe03ffffffffffffull;	// clear G6-G12
     707      if ((y & 0x0003ffffffffffffull) > 999999999999999ull) {
     708        y = y & 0xfe00000000000000ull;	// clear G6-G12 and the payload bits
     709      }
     710    } else if ((y & MASK_INF) == MASK_INF) {	// check for Infinity
     711      y = y & (MASK_SIGN | MASK_INF);
     712    } else {	// y is not special
     713      // check for non-canonical values - treated as zero
     714      if ((y & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     715        // if the steering bits are 11, then the exponent is G[0:w+1]
     716        if (((y & MASK_BINARY_SIG2) | MASK_BINARY_OR2) >
     717  	  9999999999999999ull) {
     718  	// non-canonical
     719  	y = (y & MASK_SIGN) | ((y & MASK_BINARY_EXPONENT2) << 2);
     720        }	// else canonical
     721      }	// else canonical
     722    }
     723  
     724    // NaN (CASE1)
     725    if ((x & MASK_NAN) == MASK_NAN) {	// x is NAN
     726      if ((x & MASK_SNAN) == MASK_SNAN) {	// x is SNaN
     727        // if x is SNAN, then return quiet (x)
     728        *pfpsf |= INVALID_EXCEPTION;	// set exception if SNaN
     729        x = x & 0xfdffffffffffffffull;	// quietize x
     730        res = x;
     731      } else {	// x is QNaN
     732        if ((y & MASK_NAN) == MASK_NAN) {	// y is NAN
     733  	if ((y & MASK_SNAN) == MASK_SNAN) {	// y is SNAN
     734  	  *pfpsf |= INVALID_EXCEPTION;	// set invalid flag
     735  	}
     736  	res = x;
     737        } else {
     738  	res = y;
     739        }
     740      }
     741      BID_RETURN (res);
     742    } else if ((y & MASK_NAN) == MASK_NAN) {	// y is NaN, but x is not
     743      if ((y & MASK_SNAN) == MASK_SNAN) {
     744        *pfpsf |= INVALID_EXCEPTION;	// set exception if SNaN
     745        y = y & 0xfdffffffffffffffull;	// quietize y
     746        res = y;
     747      } else {
     748        // will return x (which is not NaN)
     749        res = x;
     750      }
     751      BID_RETURN (res);
     752    }
     753    // SIMPLE (CASE2)
     754    // if all the bits are the same, these numbers are equal, return either number
     755    if (x == y) {
     756      res = x;
     757      BID_RETURN (res);
     758    }
     759    // INFINITY (CASE3)
     760    if ((x & MASK_INF) == MASK_INF) {
     761      // x is infinity, its magnitude is greater than or equal to y
     762      // return y as long as x isn't negative infinity
     763      res = ((x & MASK_SIGN) == MASK_SIGN
     764  	   && (y & MASK_INF) == MASK_INF) ? y : x;
     765      BID_RETURN (res);
     766    } else if ((y & MASK_INF) == MASK_INF) {
     767      // y is infinity, then it must be greater in magnitude
     768      res = y;
     769      BID_RETURN (res);
     770    }
     771    // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
     772    if ((x & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     773      exp_x = (x & MASK_BINARY_EXPONENT2) >> 51;
     774      sig_x = (x & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
     775    } else {
     776      exp_x = (x & MASK_BINARY_EXPONENT1) >> 53;
     777      sig_x = (x & MASK_BINARY_SIG1);
     778    }
     779  
     780    // if steering bits are 11 (condition will be 0), then exponent is G[0:w+1] =>
     781    if ((y & MASK_STEERING_BITS) == MASK_STEERING_BITS) {
     782      exp_y = (y & MASK_BINARY_EXPONENT2) >> 51;
     783      sig_y = (y & MASK_BINARY_SIG2) | MASK_BINARY_OR2;
     784    } else {
     785      exp_y = (y & MASK_BINARY_EXPONENT1) >> 53;
     786      sig_y = (y & MASK_BINARY_SIG1);
     787    }
     788  
     789    // ZERO (CASE4)
     790    // some properties:
     791    //    (+ZERO == -ZERO) => therefore 
     792    //        ignore the sign, and neither number is greater
     793    //    (ZERO x 10^A == ZERO x 10^B) for any valid A, B => 
     794    //        ignore the exponent field
     795    //    (Any non-canonical # is considered 0)
     796    if (sig_x == 0) {
     797      res = y;	// x_is_zero, its magnitude must be smaller than y
     798      BID_RETURN (res);
     799    }
     800    if (sig_y == 0) {
     801      res = x;	// y_is_zero, its magnitude must be smaller than x
     802      BID_RETURN (res);
     803    }
     804    // REDUNDANT REPRESENTATIONS (CASE6)
     805    // if both components are either bigger or smaller, 
     806    // it is clear what needs to be done
     807    if (sig_x > sig_y && exp_x >= exp_y) {
     808      res = x;
     809      BID_RETURN (res);
     810    }
     811    if (sig_x < sig_y && exp_x <= exp_y) {
     812      res = y;
     813      BID_RETURN (res);
     814    }
     815    // if exp_x is 15 greater than exp_y, no need for compensation
     816    if (exp_x - exp_y > 15) {
     817      res = x;	// difference cannot be greater than 10^15
     818      BID_RETURN (res);
     819    }
     820    // if exp_x is 15 less than exp_y, no need for compensation
     821    if (exp_y - exp_x > 15) {
     822      res = y;
     823      BID_RETURN (res);
     824    }
     825    // if |exp_x - exp_y| < 15, it comes down to the compensated significand
     826    if (exp_x > exp_y) {	// to simplify the loop below,
     827      // otherwise adjust the x significand upwards
     828      __mul_64x64_to_128MACH (sig_n_prime, sig_x,
     829  			    mult_factor[exp_x - exp_y]);
     830      // now, sig_n_prime has: sig_x * 10^(exp_x-exp_y), 
     831      // this is the compensated signif.
     832      if (sig_n_prime.w[1] == 0 && (sig_n_prime.w[0] == sig_y)) {
     833        // two numbers are equal, return maxNum(x,y)
     834        res = ((y & MASK_SIGN) == MASK_SIGN) ? x : y;
     835        BID_RETURN (res);
     836      }
     837      // now, if compensated_x (sig_n_prime) is greater than y return y,  
     838      // otherwise return x
     839      res = ((sig_n_prime.w[1] != 0) || sig_n_prime.w[0] > sig_y) ? x : y;
     840      BID_RETURN (res);
     841    }
     842    // exp_y must be greater than exp_x, thus adjust the y significand upwards
     843    __mul_64x64_to_128MACH (sig_n_prime, sig_y,
     844  			  mult_factor[exp_y - exp_x]);
     845  
     846    if (sig_n_prime.w[1] == 0 && (sig_n_prime.w[0] == sig_x)) {
     847      res = ((y & MASK_SIGN) == MASK_SIGN) ? x : y;
     848      // two numbers are equal, return either
     849      BID_RETURN (res);
     850    }
     851  
     852    res = ((sig_n_prime.w[1] == 0) && (sig_x > sig_n_prime.w[0])) ? x : y;
     853    BID_RETURN (res);
     854  }