Home | History | Annotate | Line # | Download | only in builtins
      1 //===----- lib/fp_add_impl.inc - floaing point addition -----------*- C -*-===//
      2 //
      3 //                     The LLVM Compiler Infrastructure
      4 //
      5 // This file is dual licensed under the MIT and the University of Illinois Open
      6 // Source Licenses. See LICENSE.TXT for details.
      7 //
      8 //===----------------------------------------------------------------------===//
      9 //
     10 // This file implements soft-float addition with the IEEE-754 default rounding
     11 // (to nearest, ties to even).
     12 //
     13 //===----------------------------------------------------------------------===//
     14 
     15 #include "fp_lib.h"
     16 
     17 static __inline fp_t __addXf3__(fp_t a, fp_t b) {
     18     rep_t aRep = toRep(a);
     19     rep_t bRep = toRep(b);
     20     const rep_t aAbs = aRep & absMask;
     21     const rep_t bAbs = bRep & absMask;
     22 
     23     // Detect if a or b is zero, infinity, or NaN.
     24     if (aAbs - REP_C(1) >= infRep - REP_C(1) ||
     25         bAbs - REP_C(1) >= infRep - REP_C(1)) {
     26         // NaN + anything = qNaN
     27         if (aAbs > infRep) return fromRep(toRep(a) | quietBit);
     28         // anything + NaN = qNaN
     29         if (bAbs > infRep) return fromRep(toRep(b) | quietBit);
     30 
     31         if (aAbs == infRep) {
     32             // +/-infinity + -/+infinity = qNaN
     33             if ((toRep(a) ^ toRep(b)) == signBit) return fromRep(qnanRep);
     34             // +/-infinity + anything remaining = +/- infinity
     35             else return a;
     36         }
     37 
     38         // anything remaining + +/-infinity = +/-infinity
     39         if (bAbs == infRep) return b;
     40 
     41         // zero + anything = anything
     42         if (!aAbs) {
     43             // but we need to get the sign right for zero + zero
     44             if (!bAbs) return fromRep(toRep(a) & toRep(b));
     45             else return b;
     46         }
     47 
     48         // anything + zero = anything
     49         if (!bAbs) return a;
     50     }
     51 
     52     // Swap a and b if necessary so that a has the larger absolute value.
     53     if (bAbs > aAbs) {
     54         const rep_t temp = aRep;
     55         aRep = bRep;
     56         bRep = temp;
     57     }
     58 
     59     // Extract the exponent and significand from the (possibly swapped) a and b.
     60     int aExponent = aRep >> significandBits & maxExponent;
     61     int bExponent = bRep >> significandBits & maxExponent;
     62     rep_t aSignificand = aRep & significandMask;
     63     rep_t bSignificand = bRep & significandMask;
     64 
     65     // Normalize any denormals, and adjust the exponent accordingly.
     66     if (aExponent == 0) aExponent = normalize(&aSignificand);
     67     if (bExponent == 0) bExponent = normalize(&bSignificand);
     68 
     69     // The sign of the result is the sign of the larger operand, a.  If they
     70     // have opposite signs, we are performing a subtraction; otherwise addition.
     71     const rep_t resultSign = aRep & signBit;
     72     const bool subtraction = (aRep ^ bRep) & signBit;
     73 
     74     // Shift the significands to give us round, guard and sticky, and or in the
     75     // implicit significand bit.  (If we fell through from the denormal path it
     76     // was already set by normalize( ), but setting it twice won't hurt
     77     // anything.)
     78     aSignificand = (aSignificand | implicitBit) << 3;
     79     bSignificand = (bSignificand | implicitBit) << 3;
     80 
     81     // Shift the significand of b by the difference in exponents, with a sticky
     82     // bottom bit to get rounding correct.
     83     const unsigned int align = aExponent - bExponent;
     84     if (align) {
     85         if (align < typeWidth) {
     86             const bool sticky = bSignificand << (typeWidth - align);
     87             bSignificand = bSignificand >> align | sticky;
     88         } else {
     89             bSignificand = 1; // sticky; b is known to be non-zero.
     90         }
     91     }
     92     if (subtraction) {
     93         aSignificand -= bSignificand;
     94         // If a == -b, return +zero.
     95         if (aSignificand == 0) return fromRep(0);
     96 
     97         // If partial cancellation occured, we need to left-shift the result
     98         // and adjust the exponent:
     99         if (aSignificand < implicitBit << 3) {
    100             const int shift = rep_clz(aSignificand) - rep_clz(implicitBit << 3);
    101             aSignificand <<= shift;
    102             aExponent -= shift;
    103         }
    104     }
    105     else /* addition */ {
    106         aSignificand += bSignificand;
    107 
    108         // If the addition carried up, we need to right-shift the result and
    109         // adjust the exponent:
    110         if (aSignificand & implicitBit << 4) {
    111             const bool sticky = aSignificand & 1;
    112             aSignificand = aSignificand >> 1 | sticky;
    113             aExponent += 1;
    114         }
    115     }
    116 
    117     // If we have overflowed the type, return +/- infinity:
    118     if (aExponent >= maxExponent) return fromRep(infRep | resultSign);
    119 
    120     if (aExponent <= 0) {
    121         // Result is denormal before rounding; the exponent is zero and we
    122         // need to shift the significand.
    123         const int shift = 1 - aExponent;
    124         const bool sticky = aSignificand << (typeWidth - shift);
    125         aSignificand = aSignificand >> shift | sticky;
    126         aExponent = 0;
    127     }
    128 
    129     // Low three bits are round, guard, and sticky.
    130     const int roundGuardSticky = aSignificand & 0x7;
    131 
    132     // Shift the significand into place, and mask off the implicit bit.
    133     rep_t result = aSignificand >> 3 & significandMask;
    134 
    135     // Insert the exponent and sign.
    136     result |= (rep_t)aExponent << significandBits;
    137     result |= resultSign;
    138 
    139     // Final rounding.  The result may overflow to infinity, but that is the
    140     // correct result in that case.
    141     if (roundGuardSticky > 0x4) result++;
    142     if (roundGuardSticky == 0x4) result += result & 1;
    143     return fromRep(result);
    144 }
    145