libgcc/config/libbid/bid128_mul.c

   1 /* Copyright (C) 2007  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 2, or (at your option) any later
   8 version.
   9
  10 In addition to the permissions in the GNU General Public License, the
  11 Free Software Foundation gives you unlimited permission to link the
  12 compiled version of this file into combinations with other programs,
  13 and to distribute those combinations without any restriction coming
  14 from the use of this file.  (The General Public License restrictions
  15 do apply in other respects; for example, they cover modification of
  16 the file, and distribution when not linked into a combine
  17 executable.)
  18
  19 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
  20 WARRANTY; without even the implied warranty of MERCHANTABILITY or
  21 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  22 for more details.
  23
  24 You should have received a copy of the GNU General Public License
  25 along with GCC; see the file COPYING.  If not, write to the Free
  26 Software Foundation, 51 Franklin Street, Fifth Floor, Boston, MA
  27 02110-1301, USA.  */
  28
  29 #include "bid_internal.h"
  30
  31 #if DECIMAL_CALL_BY_REFERENCE
  32 void
  33 bid64dq_mul (UINT64 * pres, UINT64 * px, UINT128 * py
  34              _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
  35              _EXC_INFO_PARAM) {
  36   UINT64 x = *px;
  37 #if !DECIMAL_GLOBAL_ROUNDING
  38   unsigned int rnd_mode = *prnd_mode;
  39 #endif
  40 #else
  41 UINT64
  42 bid64dq_mul (UINT64 x, UINT128 y
  43              _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
  44              _EXC_INFO_PARAM) {
  45 #endif
  46   UINT64 res = 0xbaddbaddbaddbaddull;
  47   UINT128 x1;
  48
  49 #if DECIMAL_CALL_BY_REFERENCE
  50   bid64_to_bid128 (&x1, &x _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
  51   bid64qq_mul (&res, &x1, py
  52                _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
  53                _EXC_INFO_ARG);
  54 #else
  55   x1 = bid64_to_bid128 (x _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
  56   res = bid64qq_mul (x1, y
  57                      _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
  58                      _EXC_INFO_ARG);
  59 #endif
  60   BID_RETURN (res);
  61 }
  62
  63
  64 #if DECIMAL_CALL_BY_REFERENCE
  65 void
  66 bid64qd_mul (UINT64 * pres, UINT128 * px, UINT64 * py
  67              _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
  68              _EXC_INFO_PARAM) {
  69   UINT64 y = *py;
  70 #if !DECIMAL_GLOBAL_ROUNDING
  71   unsigned int rnd_mode = *prnd_mode;
  72 #endif
  73 #else
  74 UINT64
  75 bid64qd_mul (UINT128 x, UINT64 y
  76              _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
  77              _EXC_INFO_PARAM) {
  78 #endif
  79   UINT64 res = 0xbaddbaddbaddbaddull;
  80   UINT128 y1;
  81
  82 #if DECIMAL_CALL_BY_REFERENCE
  83   bid64_to_bid128 (&y1, &y _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
  84   bid64qq_mul (&res, px, &y1
  85                _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
  86                _EXC_INFO_ARG);
  87 #else
  88   y1 = bid64_to_bid128 (y _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
  89   res = bid64qq_mul (x, y1
  90                      _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
  91                      _EXC_INFO_ARG);
  92 #endif
  93   BID_RETURN (res);
  94 }
  95
  96
  97 #if DECIMAL_CALL_BY_REFERENCE
  98 void
  99 bid64qq_mul (UINT64 * pres, UINT128 * px, UINT128 * py
 100              _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 101              _EXC_INFO_PARAM) {
 102   UINT128 x = *px, y = *py;
 103 #if !DECIMAL_GLOBAL_ROUNDING
 104   unsigned int rnd_mode = *prnd_mode;
 105 #endif
 106 #else
 107 UINT64
 108 bid64qq_mul (UINT128 x, UINT128 y
 109              _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 110              _EXC_INFO_PARAM) {
 111 #endif
 112
 113   UINT128 z = { {0x0000000000000000ull, 0x5ffe000000000000ull}
 114   };
 115   UINT64 res = 0xbaddbaddbaddbaddull;
 116   UINT64 x_sign, y_sign, p_sign;
 117   UINT64 x_exp, y_exp, p_exp;
 118   int true_p_exp;
 119   UINT128 C1, C2;
 120
 121   BID_SWAP128 (z);
 122   // skip cases where at least one operand is NaN or infinity
 123   if (!(((x.w[HIGH_128W] & MASK_NAN) == MASK_NAN) ||
 124         ((y.w[HIGH_128W] & MASK_NAN) == MASK_NAN) ||
 125         ((x.w[HIGH_128W] & MASK_ANY_INF) == MASK_INF) ||
 126         ((y.w[HIGH_128W] & MASK_ANY_INF) == MASK_INF))) {
 127     // x, y are 0 or f but not inf or NaN => unpack the arguments and check
 128     // for non-canonical values
 129
 130     x_sign = x.w[HIGH_128W] & MASK_SIGN;        // 0 for positive, MASK_SIGN for negative
 131     C1.w[1] = x.w[HIGH_128W] & MASK_COEFF;
 132     C1.w[0] = x.w[LOW_128W];
 133     // check for non-canonical values - treated as zero
 134     if ((x.w[HIGH_128W] & 0x6000000000000000ull) ==
 135         0x6000000000000000ull) {
 136       // G0_G1=11 => non-canonical
 137       x_exp = (x.w[HIGH_128W] << 2) & MASK_EXP; // biased and shifted left 49 bits
 138       C1.w[1] = 0;      // significand high
 139       C1.w[0] = 0;      // significand low
 140     } else {    // G0_G1 != 11
 141       x_exp = x.w[HIGH_128W] & MASK_EXP;        // biased and shifted left 49 bits
 142       if (C1.w[1] > 0x0001ed09bead87c0ull ||
 143           (C1.w[1] == 0x0001ed09bead87c0ull &&
 144            C1.w[0] > 0x378d8e63ffffffffull)) {
 145         // x is non-canonical if coefficient is larger than 10^34 -1
 146         C1.w[1] = 0;
 147         C1.w[0] = 0;
 148       } else {  // canonical
 149         ;
 150       }
 151     }
 152     y_sign = y.w[HIGH_128W] & MASK_SIGN;        // 0 for positive, MASK_SIGN for negative
 153     C2.w[1] = y.w[HIGH_128W] & MASK_COEFF;
 154     C2.w[0] = y.w[LOW_128W];
 155     // check for non-canonical values - treated as zero
 156     if ((y.w[HIGH_128W] & 0x6000000000000000ull) ==
 157         0x6000000000000000ull) {
 158       // G0_G1=11 => non-canonical
 159       y_exp = (y.w[HIGH_128W] << 2) & MASK_EXP; // biased and shifted left 49 bits
 160       C2.w[1] = 0;      // significand high
 161       C2.w[0] = 0;      // significand low
 162     } else {    // G0_G1 != 11
 163       y_exp = y.w[HIGH_128W] & MASK_EXP;        // biased and shifted left 49 bits
 164       if (C2.w[1] > 0x0001ed09bead87c0ull ||
 165           (C2.w[1] == 0x0001ed09bead87c0ull &&
 166            C2.w[0] > 0x378d8e63ffffffffull)) {
 167         // y is non-canonical if coefficient is larger than 10^34 -1
 168         C2.w[1] = 0;
 169         C2.w[0] = 0;
 170       } else {  // canonical
 171         ;
 172       }
 173     }
 174     p_sign = x_sign ^ y_sign;   // sign of the product
 175
 176     true_p_exp = (x_exp >> 49) - 6176 + (y_exp >> 49) - 6176;
 177     // true_p_exp, p_exp are used only for 0 * 0, 0 * f, or f * 0
 178     if (true_p_exp < -398)
 179       p_exp = 0;        // cannot be less than EXP_MIN
 180     else if (true_p_exp > 369)
 181       p_exp = (UINT64) (369 + 398) << 53;       // cannot be more than EXP_MAX
 182     else
 183       p_exp = (UINT64) (true_p_exp + 398) << 53;
 184
 185     if ((C1.w[1] == 0x0 && C1.w[0] == 0x0) ||
 186         (C2.w[1] == 0x0 && C2.w[0] == 0x0)) {
 187       // x = 0 or y = 0
 188       // the result is 0
 189       res = p_sign | p_exp;     // preferred exponent in [EXP_MIN, EXP_MAX]
 190       BID_RETURN (res)
 191     }   // else continue
 192   }
 193   // swap x and y - ensure that a NaN in x has 'higher precedence' than one in y
 194 #if DECIMAL_CALL_BY_REFERENCE
 195   bid64qqq_fma (&res, &y, &x, &z
 196                 _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 197                 _EXC_INFO_ARG);
 198 #else
 199   res = bid64qqq_fma (y, x, z
 200                       _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 201                       _EXC_INFO_ARG);
 202 #endif
 203   BID_RETURN (res);
 204 }
 205
 206
 207 #if DECIMAL_CALL_BY_REFERENCE
 208 void
 209 bid128dd_mul (UINT128 * pres, UINT64 * px, UINT64 * py
 210               _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 211               _EXC_INFO_PARAM) {
 212   UINT64 x = *px, y = *py;
 213 #if !DECIMAL_GLOBAL_ROUNDING
 214   unsigned int rnd_mode = *prnd_mode;
 215 #endif
 216 #else
 217 UINT128
 218 bid128dd_mul (UINT64 x, UINT64 y
 219               _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 220               _EXC_INFO_PARAM) {
 221 #endif
 222   UINT128 res = { {0xbaddbaddbaddbaddull, 0xbaddbaddbaddbaddull}
 223   };
 224   UINT128 x1, y1;
 225
 226 #if DECIMAL_CALL_BY_REFERENCE
 227   bid64_to_bid128 (&x1, &x _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
 228   bid64_to_bid128 (&y1, &y _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
 229   bid128_mul (&res, &x1, &y1
 230               _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 231               _EXC_INFO_ARG);
 232 #else
 233   x1 = bid64_to_bid128 (x _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
 234   y1 = bid64_to_bid128 (y _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
 235   res = bid128_mul (x1, y1
 236                     _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 237                     _EXC_INFO_ARG);
 238 #endif
 239   BID_RETURN (res);
 240 }
 241
 242
 243 #if DECIMAL_CALL_BY_REFERENCE
 244 void
 245 bid128dq_mul (UINT128 * pres, UINT64 * px, UINT128 * py
 246               _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 247               _EXC_INFO_PARAM) {
 248   UINT64 x = *px;
 249 #if !DECIMAL_GLOBAL_ROUNDING
 250   unsigned int rnd_mode = *prnd_mode;
 251 #endif
 252 #else
 253 UINT128
 254 bid128dq_mul (UINT64 x, UINT128 y
 255               _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 256               _EXC_INFO_PARAM) {
 257 #endif
 258   UINT128 res = { {0xbaddbaddbaddbaddull, 0xbaddbaddbaddbaddull}
 259   };
 260   UINT128 x1;
 261
 262 #if DECIMAL_CALL_BY_REFERENCE
 263   bid64_to_bid128 (&x1, &x _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
 264   bid128_mul (&res, &x1, py
 265               _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 266               _EXC_INFO_ARG);
 267 #else
 268   x1 = bid64_to_bid128 (x _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
 269   res = bid128_mul (x1, y
 270                     _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 271                     _EXC_INFO_ARG);
 272 #endif
 273   BID_RETURN (res);
 274 }
 275
 276
 277 #if DECIMAL_CALL_BY_REFERENCE
 278 void
 279 bid128qd_mul (UINT128 * pres, UINT128 * px, UINT64 * py
 280               _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 281               _EXC_INFO_PARAM) {
 282   UINT64 y = *py;
 283 #if !DECIMAL_GLOBAL_ROUNDING
 284   unsigned int rnd_mode = *prnd_mode;
 285 #endif
 286 #else
 287 UINT128
 288 bid128qd_mul (UINT128 x, UINT64 y
 289               _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 290               _EXC_INFO_PARAM) {
 291 #endif
 292   UINT128 res = { {0xbaddbaddbaddbaddull, 0xbaddbaddbaddbaddull}
 293   };
 294   UINT128 y1;
 295
 296 #if DECIMAL_CALL_BY_REFERENCE
 297   bid64_to_bid128 (&y1, &y _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
 298   bid128_mul (&res, px, &y1
 299               _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 300               _EXC_INFO_ARG);
 301 #else
 302   y1 = bid64_to_bid128 (y _EXC_FLAGS_ARG _EXC_MASKS_ARG _EXC_INFO_ARG);
 303   res = bid128_mul (x, y1
 304                     _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 305                     _EXC_INFO_ARG);
 306 #endif
 307   BID_RETURN (res);
 308 }
 309
 310
 311 // bid128_mul stands for bid128qq_mul
 312 #if DECIMAL_CALL_BY_REFERENCE
 313 void
 314 bid128_mul (UINT128 * pres, UINT128 * px,
 315             UINT128 *
 316             py _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 317             _EXC_INFO_PARAM) {
 318   UINT128 x = *px, y = *py;
 319
 320 #if !DECIMAL_GLOBAL_ROUNDING
 321   unsigned int rnd_mode = *prnd_mode;
 322
 323 #endif
 324 #else
 325 UINT128
 326 bid128_mul (UINT128 x,
 327             UINT128 y _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
 328             _EXC_INFO_PARAM) {
 329
 330 #endif
 331   UINT128 z = { {0x0000000000000000ull, 0x5ffe000000000000ull}
 332   };
 333   UINT128 res = { {0xbaddbaddbaddbaddull, 0xbaddbaddbaddbaddull}
 334   };
 335   UINT64 x_sign, y_sign, p_sign;
 336   UINT64 x_exp, y_exp, p_exp;
 337   int true_p_exp;
 338   UINT128 C1, C2;
 339
 340   BID_SWAP128 (x);
 341   BID_SWAP128 (y);
 342   // skip cases where at least one operand is NaN or infinity
 343   if (!(((x.w[1] & MASK_NAN) == MASK_NAN) ||
 344         ((y.w[1] & MASK_NAN) == MASK_NAN) ||
 345         ((x.w[1] & MASK_ANY_INF) == MASK_INF) ||
 346         ((y.w[1] & MASK_ANY_INF) == MASK_INF))) {
 347     // x, y are 0 or f but not inf or NaN => unpack the arguments and check
 348     // for non-canonical values
 349
 350     x_sign = x.w[1] & MASK_SIGN;        // 0 for positive, MASK_SIGN for negative
 351     C1.w[1] = x.w[1] & MASK_COEFF;
 352     C1.w[0] = x.w[0];
 353     // check for non-canonical values - treated as zero
 354     if ((x.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull) {
 355       // G0_G1=11 => non-canonical
 356       x_exp = (x.w[1] << 2) & MASK_EXP; // biased and shifted left 49 bits
 357       C1.w[1] = 0;      // significand high
 358       C1.w[0] = 0;      // significand low
 359     } else {    // G0_G1 != 11
 360       x_exp = x.w[1] & MASK_EXP;        // biased and shifted left 49 bits
 361       if (C1.w[1] > 0x0001ed09bead87c0ull ||
 362           (C1.w[1] == 0x0001ed09bead87c0ull &&
 363            C1.w[0] > 0x378d8e63ffffffffull)) {
 364         // x is non-canonical if coefficient is larger than 10^34 -1
 365         C1.w[1] = 0;
 366         C1.w[0] = 0;
 367       } else {  // canonical
 368         ;
 369       }
 370     }
 371     y_sign = y.w[1] & MASK_SIGN;        // 0 for positive, MASK_SIGN for negative
 372     C2.w[1] = y.w[1] & MASK_COEFF;
 373     C2.w[0] = y.w[0];
 374     // check for non-canonical values - treated as zero
 375     if ((y.w[1] & 0x6000000000000000ull) == 0x6000000000000000ull) {
 376       // G0_G1=11 => non-canonical
 377       y_exp = (y.w[1] << 2) & MASK_EXP; // biased and shifted left 49 bits
 378       C2.w[1] = 0;      // significand high
 379       C2.w[0] = 0;      // significand low
 380     } else {    // G0_G1 != 11
 381       y_exp = y.w[1] & MASK_EXP;        // biased and shifted left 49 bits
 382       if (C2.w[1] > 0x0001ed09bead87c0ull ||
 383           (C2.w[1] == 0x0001ed09bead87c0ull &&
 384            C2.w[0] > 0x378d8e63ffffffffull)) {
 385         // y is non-canonical if coefficient is larger than 10^34 -1
 386         C2.w[1] = 0;
 387         C2.w[0] = 0;
 388       } else {  // canonical
 389         ;
 390       }
 391     }
 392     p_sign = x_sign ^ y_sign;   // sign of the product
 393
 394     true_p_exp = (x_exp >> 49) - 6176 + (y_exp >> 49) - 6176;
 395     // true_p_exp, p_exp are used only for 0 * 0, 0 * f, or f * 0
 396     if (true_p_exp < -6176)
 397       p_exp = 0;        // cannot be less than EXP_MIN
 398     else if (true_p_exp > 6111)
 399       p_exp = (UINT64) (6111 + 6176) << 49;     // cannot be more than EXP_MAX
 400     else
 401       p_exp = (UINT64) (true_p_exp + 6176) << 49;
 402
 403     if ((C1.w[1] == 0x0 && C1.w[0] == 0x0) ||
 404         (C2.w[1] == 0x0 && C2.w[0] == 0x0)) {
 405       // x = 0 or y = 0
 406       // the result is 0
 407       res.w[1] = p_sign | p_exp;        // preferred exponent in [EXP_MIN, EXP_MAX]
 408       res.w[0] = 0x0;
 409       BID_SWAP128 (res);
 410       BID_RETURN (res)
 411     }   // else continue
 412   }
 413
 414   BID_SWAP128 (x);
 415   BID_SWAP128 (y);
 416   BID_SWAP128 (z);
 417   // swap x and y - ensure that a NaN in x has 'higher precedence' than one in y
 418 #if DECIMAL_CALL_BY_REFERENCE
 419   bid128_fma (&res, &y, &x, &z
 420               _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 421               _EXC_INFO_ARG);
 422 #else
 423   res = bid128_fma (y, x, z
 424                     _RND_MODE_ARG _EXC_FLAGS_ARG _EXC_MASKS_ARG
 425                     _EXC_INFO_ARG);
 426 #endif
 427   BID_RETURN (res);
 428 }