gcc/sreal.c

   1 /* Simple data type for positive real numbers for the GNU compiler.
   2    Copyright (C) 2002-2014 Free Software Foundation, Inc.
   3
   4 This file is part of GCC.
   5
   6 GCC is free software; you can redistribute it and/or modify it under
   7 the terms of the GNU General Public License as published by the Free
   8 Software Foundation; either version 3, or (at your option) any later
   9 version.
  10
  11 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
  12 WARRANTY; without even the implied warranty of MERCHANTABILITY or
  13 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  14 for more details.
  15
  16 You should have received a copy of the GNU General Public License
  17 along with GCC; see the file COPYING3.  If not see
  18 <http://www.gnu.org/licenses/>.  */
  19
  20 /* This library supports positive real numbers and 0;
  21    inf and nan are NOT supported.
  22    It is written to be simple and fast.
  23
  24    Value of sreal is
  25         x = sig * 2 ^ exp
  26    where
  27         sig = significant
  28           (for < 64-bit machines sig = sig_lo + sig_hi * 2 ^ SREAL_PART_BITS)
  29         exp = exponent
  30
  31    One uint64_t is used for the significant.
  32    Only a half of significant bits is used (in normalized sreals) so that we do
  33    not have problems with overflow, for example when c->sig = a->sig * b->sig.
  34    So the precision is 32-bit.
  35
  36    Invariant: The numbers are normalized before and after each call of sreal_*.
  37
  38    Normalized sreals:
  39    All numbers (except zero) meet following conditions:
  40          SREAL_MIN_SIG <= sig && sig <= SREAL_MAX_SIG
  41         -SREAL_MAX_EXP <= exp && exp <= SREAL_MAX_EXP
  42
  43    If the number would be too large, it is set to upper bounds of these
  44    conditions.
  45
  46    If the number is zero or would be too small it meets following conditions:
  47         sig == 0 && exp == -SREAL_MAX_EXP
  48 */
  49
  50 #include "config.h"
  51 #include "system.h"
  52 #include "coretypes.h"
  53 #include "sreal.h"
  54
  55 /* Print the content of struct sreal.  */
  56
  57 void
  58 sreal::dump (FILE *file) const
  59 {
  60   fprintf (file, "(%" PRIu64 " * 2^%d)", m_sig, m_exp);
  61 }
  62
  63 DEBUG_FUNCTION void
  64 debug (sreal &ref)
  65 {
  66   ref.dump (stderr);
  67 }
  68
  69 DEBUG_FUNCTION void
  70 debug (sreal *ptr)
  71 {
  72   if (ptr)
  73     debug (*ptr);
  74   else
  75     fprintf (stderr, "<nil>\n");
  76 }
  77
  78 /* Shift this right by S bits.  Needed: 0 < S <= SREAL_BITS.
  79    When the most significant bit shifted out is 1, add 1 to this (rounding).
  80    */
  81
  82 void
  83 sreal::shift_right (int s)
  84 {
  85   gcc_assert (s > 0);
  86   gcc_assert (s <= SREAL_BITS);
  87   /* Exponent should never be so large because shift_right is used only by
  88      sreal_add and sreal_sub ant thus the number cannot be shifted out from
  89      exponent range.  */
  90   gcc_assert (m_exp + s <= SREAL_MAX_EXP);
  91
  92   m_exp += s;
  93
  94   m_sig += (uint64_t) 1 << (s - 1);
  95   m_sig >>= s;
  96 }
  97
  98 /* Normalize *this.  */
  99
 100 void
 101 sreal::normalize ()
 102 {
 103   if (m_sig == 0)
 104     {
 105       m_exp = -SREAL_MAX_EXP;
 106     }
 107   else if (m_sig < SREAL_MIN_SIG)
 108     {
 109       do
 110         {
 111           m_sig <<= 1;
 112           m_exp--;
 113         }
 114       while (m_sig < SREAL_MIN_SIG);
 115
 116       /* Check underflow.  */
 117       if (m_exp < -SREAL_MAX_EXP)
 118         {
 119           m_exp = -SREAL_MAX_EXP;
 120           m_sig = 0;
 121         }
 122     }
 123   else if (m_sig > SREAL_MAX_SIG)
 124     {
 125       int last_bit;
 126       do
 127         {
 128           last_bit = m_sig & 1;
 129           m_sig >>= 1;
 130           m_exp++;
 131         }
 132       while (m_sig > SREAL_MAX_SIG);
 133
 134       /* Round the number.  */
 135       m_sig += last_bit;
 136       if (m_sig > SREAL_MAX_SIG)
 137         {
 138           m_sig >>= 1;
 139           m_exp++;
 140         }
 141
 142       /* Check overflow.  */
 143       if (m_exp > SREAL_MAX_EXP)
 144         {
 145           m_exp = SREAL_MAX_EXP;
 146           m_sig = SREAL_MAX_SIG;
 147         }
 148     }
 149 }
 150
 151 /* Return integer value of *this.  */
 152
 153 int64_t
 154 sreal::to_int () const
 155 {
 156   if (m_exp <= -SREAL_BITS)
 157     return 0;
 158   if (m_exp >= SREAL_PART_BITS)
 159     return INTTYPE_MAXIMUM (int64_t);
 160   if (m_exp > 0)
 161     return m_sig << m_exp;
 162   if (m_exp < 0)
 163     return m_sig >> -m_exp;
 164   return m_sig;
 165 }
 166
 167 /* Return *this + other.  */
 168
 169 sreal
 170 sreal::operator+ (const sreal &other) const
 171 {
 172   int dexp;
 173   sreal tmp, r;
 174 const sreal *a_p = this, *b_p = &other, *bb;
 175
 176   if (*a_p < *b_p)
 177     {
 178       const sreal *swap;
 179       swap = a_p;
 180       a_p = b_p;
 181       b_p = swap;
 182     }
 183
 184   dexp = a_p->m_exp - b_p->m_exp;
 185   r.m_exp = a_p->m_exp;
 186   if (dexp > SREAL_BITS)
 187     {
 188       r.m_sig = a_p->m_sig;
 189       return r;
 190     }
 191
 192   if (dexp == 0)
 193     bb = b_p;
 194   else
 195     {
 196       tmp = *b_p;
 197       tmp.shift_right (dexp);
 198       bb = &tmp;
 199     }
 200
 201   r.m_sig = a_p->m_sig + bb->m_sig;
 202   r.normalize ();
 203   return r;
 204 }
 205
 206 /* Return *this - other.  */
 207
 208 sreal
 209 sreal::operator- (const sreal &other) const
 210 {
 211   int dexp;
 212   sreal tmp, r;
 213   const sreal *bb;
 214
 215   gcc_assert (*this >= other);
 216
 217   dexp = m_exp - other.m_exp;
 218   r.m_exp = m_exp;
 219   if (dexp > SREAL_BITS)
 220     {
 221       r.m_sig = m_sig;
 222       return r;
 223     }
 224   if (dexp == 0)
 225     bb = &other;
 226   else
 227     {
 228       tmp = other;
 229       tmp.shift_right (dexp);
 230       bb = &tmp;
 231     }
 232
 233   r.m_sig = m_sig - bb->m_sig;
 234   r.normalize ();
 235   return r;
 236 }
 237
 238 /* Return *this * other.  */
 239
 240 sreal
 241 sreal::operator* (const sreal &other) const
 242 {
 243 sreal r;
 244   if (m_sig < SREAL_MIN_SIG || other.m_sig < SREAL_MIN_SIG)
 245     {
 246       r.m_sig = 0;
 247       r.m_exp = -SREAL_MAX_EXP;
 248     }
 249   else
 250     {
 251       r.m_sig = m_sig * other.m_sig;
 252       r.m_exp = m_exp + other.m_exp;
 253       r.normalize ();
 254     }
 255   return r;
 256 }
 257
 258 /* Return *this / other.  */
 259
 260 sreal
 261 sreal::operator/ (const sreal &other) const
 262 {
 263   gcc_assert (other.m_sig != 0);
 264 sreal r;
 265   r.m_sig = (m_sig << SREAL_PART_BITS) / other.m_sig;
 266   r.m_exp = m_exp - other.m_exp - SREAL_PART_BITS;
 267   r.normalize ();
 268   return r;
 269 }