1 /* Simple data type for real numbers for the GNU compiler.
2 Copyright (C) 2002-2014 Free Software Foundation, Inc.
4 This file is part of GCC.
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
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
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/>. */
20 /* This library supports real numbers;
21 inf and nan are NOT supported.
22 It is written to be simple and fast.
28 (for < 64-bit machines sig = sig_lo + sig_hi * 2 ^ SREAL_PART_BITS)
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.
36 Invariant: The numbers are normalized before and after each call of sreal_*.
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
43 If the number would be too large, it is set to upper bounds of these
46 If the number is zero or would be too small it meets following conditions:
47 sig == 0 && exp == -SREAL_MAX_EXP
52 #include "coretypes.h"
55 /* Print the content of struct sreal. */
58 sreal::dump (FILE *file
) const
60 fprintf (file
, "(%" PRIu64
" * 2^%d)", m_sig
, m_exp
);
75 fprintf (stderr
, "<nil>\n");
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).
83 sreal::shift_right (int s
)
85 gcc_checking_assert (s
> 0);
86 gcc_checking_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
90 gcc_checking_assert (m_exp
+ s
<= SREAL_MAX_EXP
);
94 m_sig
+= (uint64_t) 1 << (s
- 1);
98 /* Normalize *this. */
106 m_exp
= -SREAL_MAX_EXP
;
108 else if (m_sig
< SREAL_MIN_SIG
)
115 while (m_sig
< SREAL_MIN_SIG
);
117 /* Check underflow. */
118 if (m_exp
< -SREAL_MAX_EXP
)
120 m_exp
= -SREAL_MAX_EXP
;
124 else if (m_sig
> SREAL_MAX_SIG
)
129 last_bit
= m_sig
& 1;
133 while (m_sig
> SREAL_MAX_SIG
);
135 /* Round the number. */
137 if (m_sig
> SREAL_MAX_SIG
)
143 /* Check overflow. */
144 if (m_exp
> SREAL_MAX_EXP
)
146 m_exp
= SREAL_MAX_EXP
;
147 m_sig
= SREAL_MAX_SIG
;
152 /* Return integer value of *this. */
155 sreal::to_int () const
157 int64_t sign
= m_negative
? -1 : 1;
159 if (m_exp
<= -SREAL_BITS
)
161 if (m_exp
>= SREAL_PART_BITS
)
162 return sign
* INTTYPE_MAXIMUM (int64_t);
164 return sign
* (m_sig
<< m_exp
);
166 return sign
* (m_sig
>> -m_exp
);
170 /* Return *this + other. */
173 sreal::operator+ (const sreal
&other
) const
175 const sreal
*a_p
= this, *b_p
= &other
;
177 if (a_p
->m_negative
&& !b_p
->m_negative
)
178 std::swap (a_p
, b_p
);
180 /* a + -b => a - b. */
181 if (!a_p
->m_negative
&& b_p
->m_negative
)
185 return signedless_minus (tmp
, *a_p
, false);
187 return signedless_minus (*a_p
, tmp
, true);
190 gcc_checking_assert (a_p
->m_negative
== b_p
->m_negative
);
192 sreal r
= signedless_plus (*a_p
, *b_p
, a_p
->m_negative
);
198 sreal::signedless_plus (const sreal
&a
, const sreal
&b
, bool negative
)
203 const sreal
*a_p
= &a
;
204 const sreal
*b_p
= &b
;
207 std::swap (a_p
, b_p
);
209 dexp
= a_p
->m_exp
- b_p
->m_exp
;
210 r
.m_exp
= a_p
->m_exp
;
211 if (dexp
> SREAL_BITS
)
213 r
.m_sig
= a_p
->m_sig
;
222 tmp
.shift_right (dexp
);
226 r
.m_sig
= a_p
->m_sig
+ bb
->m_sig
;
229 r
.m_negative
= negative
;
233 /* Return *this - other. */
236 sreal::operator- (const sreal
&other
) const
238 /* -a - b => -a + (-b). */
239 if (m_negative
&& !other
.m_negative
)
240 return signedless_plus (*this, -other
, true);
242 /* a - (-b) => a + b. */
243 if (!m_negative
&& other
.m_negative
)
244 return signedless_plus (*this, -other
, false);
246 gcc_checking_assert (m_negative
== other
.m_negative
);
248 /* We want to substract a smaller number from bigger
249 for nonegative numbers. */
250 if (!m_negative
&& *this < other
)
251 return -signedless_minus (other
, *this, true);
253 /* Example: -2 - (-3) => 3 - 2 */
254 if (m_negative
&& *this > other
)
255 return signedless_minus (-other
, -(*this), true);
257 sreal r
= signedless_minus (*this, other
, m_negative
);
263 sreal::signedless_minus (const sreal
&a
, const sreal
&b
, bool negative
)
268 const sreal
*a_p
= &a
;
269 const sreal
*b_p
= &b
;
271 dexp
= a_p
->m_exp
- b_p
->m_exp
;
273 r
.m_exp
= a_p
->m_exp
;
274 if (dexp
> SREAL_BITS
)
276 r
.m_sig
= a_p
->m_sig
;
284 tmp
.shift_right (dexp
);
288 r
.m_sig
= a_p
->m_sig
- bb
->m_sig
;
291 r
.m_negative
= negative
;
295 /* Return *this * other. */
298 sreal::operator* (const sreal
&other
) const
301 if (m_sig
< SREAL_MIN_SIG
|| other
.m_sig
< SREAL_MIN_SIG
)
304 r
.m_exp
= -SREAL_MAX_EXP
;
308 r
.m_sig
= m_sig
* other
.m_sig
;
309 r
.m_exp
= m_exp
+ other
.m_exp
;
313 r
.m_negative
= m_negative
^ other
.m_negative
;
317 /* Return *this / other. */
320 sreal::operator/ (const sreal
&other
) const
322 gcc_checking_assert (other
.m_sig
!= 0);
324 r
.m_sig
= (m_sig
<< SREAL_PART_BITS
) / other
.m_sig
;
325 r
.m_exp
= m_exp
- other
.m_exp
- SREAL_PART_BITS
;
326 r
.m_negative
= m_negative
^ other
.m_negative
;