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
);
64 debug (const sreal
&ref
)
70 debug (const sreal
*ptr
)
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
+= (int64_t) 1 << (s
- 1);
98 /* Normalize *this. */
103 int64_t s
= m_sig
< 0 ? -1 : 1;
104 unsigned HOST_WIDE_INT sig
= absu_hwi (m_sig
);
108 m_exp
= -SREAL_MAX_EXP
;
110 else if (sig
< SREAL_MIN_SIG
)
117 while (sig
< SREAL_MIN_SIG
);
119 /* Check underflow. */
120 if (m_exp
< -SREAL_MAX_EXP
)
122 m_exp
= -SREAL_MAX_EXP
;
126 else if (sig
> SREAL_MAX_SIG
)
135 while (sig
> SREAL_MAX_SIG
);
137 /* Round the number. */
139 if (sig
> SREAL_MAX_SIG
)
145 /* Check overflow. */
146 if (m_exp
> SREAL_MAX_EXP
)
148 m_exp
= SREAL_MAX_EXP
;
156 /* Return integer value of *this. */
159 sreal::to_int () const
161 int64_t sign
= m_sig
< 0 ? -1 : 1;
163 if (m_exp
<= -SREAL_BITS
)
165 if (m_exp
>= SREAL_PART_BITS
)
166 return sign
* INTTYPE_MAXIMUM (int64_t);
168 return m_sig
<< m_exp
;
170 return m_sig
>> -m_exp
;
174 /* Return *this + other. */
177 sreal::operator+ (const sreal
&other
) const
182 const sreal
*a_p
= this, *b_p
= &other
, *bb
;
184 if (a_p
->m_exp
< b_p
->m_exp
)
185 std::swap (a_p
, b_p
);
187 dexp
= a_p
->m_exp
- b_p
->m_exp
;
188 r
.m_exp
= a_p
->m_exp
;
189 if (dexp
> SREAL_BITS
)
191 r
.m_sig
= a_p
->m_sig
;
200 tmp
.shift_right (dexp
);
204 r
.m_sig
= a_p
->m_sig
+ bb
->m_sig
;
210 /* Return *this - other. */
213 sreal::operator- (const sreal
&other
) const
218 const sreal
*a_p
= this, *b_p
= &other
;
221 if (a_p
->m_exp
< b_p
->m_exp
)
224 std::swap (a_p
, b_p
);
227 dexp
= a_p
->m_exp
- b_p
->m_exp
;
228 r
.m_exp
= a_p
->m_exp
;
229 if (dexp
> SREAL_BITS
)
231 r
.m_sig
= sign
* a_p
->m_sig
;
239 tmp
.shift_right (dexp
);
243 r
.m_sig
= sign
* (a_p
->m_sig
- bb
->m_sig
);
248 /* Return *this * other. */
251 sreal::operator* (const sreal
&other
) const
254 if (std::abs (m_sig
) < SREAL_MIN_SIG
|| std::abs (other
.m_sig
) < SREAL_MIN_SIG
)
257 r
.m_exp
= -SREAL_MAX_EXP
;
261 r
.m_sig
= m_sig
* other
.m_sig
;
262 r
.m_exp
= m_exp
+ other
.m_exp
;
269 /* Return *this / other. */
272 sreal::operator/ (const sreal
&other
) const
274 gcc_checking_assert (other
.m_sig
!= 0);
276 r
.m_sig
= (m_sig
<< SREAL_PART_BITS
) / other
.m_sig
;
277 r
.m_exp
= m_exp
- other
.m_exp
- SREAL_PART_BITS
;