Branch data Line data Source code
1 : : /******************************************************************************
2 : : * Top contributors (to current version):
3 : : * Aina Niemetz, Martin Brain, Haniel Barbosa
4 : : *
5 : : * This file is part of the cvc5 project.
6 : : *
7 : : * Copyright (c) 2009-2025 by the authors listed in the file AUTHORS
8 : : * in the top-level source directory and their institutional affiliations.
9 : : * All rights reserved. See the file COPYING in the top-level source
10 : : * directory for licensing information.
11 : : * ****************************************************************************
12 : : *
13 : : * A floating-point value.
14 : : *
15 : : * This file contains the data structures used by the constant and parametric
16 : : * types of the floating point theory.
17 : : */
18 : :
19 : : #include "util/floatingpoint.h"
20 : :
21 : : #include <math.h>
22 : :
23 : : #include <limits>
24 : :
25 : : #include "base/check.h"
26 : : #include "util/floatingpoint_literal_symfpu.h"
27 : : #include "util/integer.h"
28 : :
29 : : /* -------------------------------------------------------------------------- */
30 : :
31 : : namespace cvc5::internal {
32 : :
33 : : /* -------------------------------------------------------------------------- */
34 : :
35 : 185 : uint32_t FloatingPoint::getUnpackedExponentWidth(FloatingPointSize& size)
36 : : {
37 : 185 : return FloatingPointLiteral::getUnpackedExponentWidth(size);
38 : : }
39 : :
40 : 185 : uint32_t FloatingPoint::getUnpackedSignificandWidth(FloatingPointSize& size)
41 : : {
42 : 185 : return FloatingPointLiteral::getUnpackedSignificandWidth(size);
43 : : }
44 : :
45 : 21742 : FloatingPoint::FloatingPoint(uint32_t d_exp_size,
46 : : uint32_t d_sig_size,
47 : 21742 : const BitVector& bv)
48 : 21742 : : d_fpl(new FloatingPointLiteral(d_exp_size, d_sig_size, bv))
49 : : {
50 : 21742 : }
51 : :
52 : 32 : FloatingPoint::FloatingPoint(const FloatingPointSize& size, const BitVector& bv)
53 : 32 : : d_fpl(new FloatingPointLiteral(size, bv))
54 : : {
55 : 32 : }
56 : :
57 : 192 : FloatingPoint::FloatingPoint(const FloatingPointSize& size,
58 : : const RoundingMode& rm,
59 : : const BitVector& bv,
60 : 192 : bool signedBV)
61 : 192 : : d_fpl(new FloatingPointLiteral(size, rm, bv, signedBV))
62 : : {
63 : 192 : }
64 : :
65 : 42862 : FloatingPoint::FloatingPoint(FloatingPointLiteral* fpl) { d_fpl.reset(fpl); }
66 : :
67 : 87161 : FloatingPoint::FloatingPoint(const FloatingPoint& fp)
68 : : {
69 : 87161 : d_fpl.reset(new FloatingPointLiteral(*fp.d_fpl));
70 : 87161 : }
71 : :
72 : 83 : FloatingPoint::FloatingPoint(const FloatingPointSize& size,
73 : : const RoundingMode& rm,
74 : 83 : const Rational& r)
75 : : {
76 : 166 : Rational two(2, 1);
77 : :
78 [ + + ]: 83 : if (r.isZero())
79 : : {
80 : : // In keeping with the SMT-LIB standard
81 : 16 : d_fpl.reset(new FloatingPointLiteral(
82 : 16 : size, FloatingPointLiteral::SpecialConstKind::FPZERO, false));
83 : : }
84 : : else
85 : : {
86 [ - + ]: 67 : uint32_t negative = (r.sgn() < 0) ? 1 : 0;
87 : 134 : Rational rabs(r.abs());
88 : :
89 : : // Compute the exponent
90 : 134 : Integer exp(0U);
91 : 134 : Integer inc(1U);
92 : 134 : Rational working(1, 1);
93 : :
94 [ + + ]: 67 : if (rabs != working)
95 : : {
96 [ + + ]: 54 : if (rabs < working)
97 : : {
98 [ + + ]: 94 : while (rabs < working)
99 : : {
100 : 70 : exp -= inc;
101 : 70 : working /= two;
102 : : }
103 : : }
104 : : else
105 : : {
106 [ + + ]: 130 : while (rabs >= working)
107 : : {
108 : 100 : exp += inc;
109 : 100 : working *= two;
110 : : }
111 : 30 : exp -= inc;
112 : 30 : working /= two;
113 : : }
114 : : }
115 : :
116 [ - + ][ - + ]: 67 : Assert(working <= rabs);
[ - - ]
117 [ - + ][ - + ]: 67 : Assert(rabs < working * two);
[ - - ]
118 : :
119 : : // Work out the number of bits required to represent the exponent for a
120 : : // normal number
121 : 67 : uint32_t expBits = 2; // No point starting with an invalid amount
122 : :
123 : 134 : Integer doubleInt(2);
124 [ + + ]: 67 : if (exp.strictlyPositive())
125 : : {
126 : : // 1 more than exactly representable with expBits
127 : 60 : Integer representable(4);
128 [ + + ]: 34 : while (representable <= exp)
129 : : { // hence <=
130 : 4 : representable *= doubleInt;
131 : 4 : ++expBits;
132 : : }
133 : : }
134 [ + + ]: 37 : else if (exp.strictlyNegative())
135 : : {
136 : 48 : Integer representable(-4); // Exactly representable with expBits + sign
137 : : // but -2^n and -(2^n - 1) are both subnormal
138 [ + + ]: 36 : while ((representable + doubleInt) > exp)
139 : : {
140 : 12 : representable *= doubleInt;
141 : 12 : ++expBits;
142 : : }
143 : : }
144 : 67 : ++expBits; // To allow for sign
145 : :
146 : 134 : BitVector exactExp(expBits, exp);
147 : :
148 : : // Compute the significand.
149 : 67 : uint32_t sigBits = size.significandWidth() + 2; // guard and sticky bits
150 : 134 : BitVector sig(sigBits, 0U);
151 : 134 : BitVector one(sigBits, 1U);
152 : 134 : Rational workingSig(0, 1);
153 [ + + ]: 2412 : for (uint32_t i = 0; i < sigBits - 1; ++i)
154 : : {
155 : 4690 : Rational mid(workingSig + working);
156 : :
157 [ + + ]: 2345 : if (mid <= rabs)
158 : : {
159 : 594 : sig = sig.setBit(0, true);
160 : 594 : workingSig = mid;
161 : : }
162 : :
163 : 2345 : sig = sig.leftShift(one);
164 : 2345 : working /= two;
165 : : }
166 : :
167 : : // Compute the sticky bit
168 : 134 : Rational remainder(rabs - workingSig);
169 [ - + ][ - + ]: 67 : Assert(Rational(0, 1) <= remainder);
[ - - ]
170 : :
171 [ + + ]: 67 : if (!remainder.isZero())
172 : : {
173 : 24 : sig = sig.setBit(0, true);
174 : : }
175 : :
176 : : // Build an exact float
177 : 67 : FloatingPointSize exactFormat(expBits, sigBits);
178 : :
179 : : // A small subtlety... if the format has expBits the unpacked format
180 : : // may have more to allow subnormals to be normalised.
181 : : // Thus...
182 : : uint32_t extension =
183 : 67 : FloatingPointLiteral::getUnpackedExponentWidth(exactFormat) - expBits;
184 : :
185 : : FloatingPointLiteral exactFloat(
186 : 134 : exactFormat, negative, exactExp.signExtend(extension), sig);
187 : :
188 : : // Then cast...
189 : 67 : d_fpl.reset(new FloatingPointLiteral(exactFloat.convert(size, rm)));
190 : : }
191 : 83 : }
192 : :
193 : 152034 : FloatingPoint::~FloatingPoint()
194 : : {
195 : 152034 : }
196 : :
197 : 215517 : const FloatingPointSize& FloatingPoint::getSize() const
198 : : {
199 : 215517 : return d_fpl->getSize();
200 : : }
201 : :
202 : 210 : FloatingPoint FloatingPoint::makeNaN(const FloatingPointSize& size)
203 : : {
204 : : return FloatingPoint(new FloatingPointLiteral(
205 : 210 : size, FloatingPointLiteral::SpecialConstKind::FPNAN));
206 : : }
207 : :
208 : 362 : FloatingPoint FloatingPoint::makeInf(const FloatingPointSize& size, bool sign)
209 : : {
210 : : return FloatingPoint(new FloatingPointLiteral(
211 : 362 : size, FloatingPointLiteral::SpecialConstKind::FPINF, sign));
212 : : }
213 : :
214 : 384 : FloatingPoint FloatingPoint::makeZero(const FloatingPointSize& size, bool sign)
215 : : {
216 : : return FloatingPoint(new FloatingPointLiteral(
217 : 384 : size, FloatingPointLiteral::SpecialConstKind::FPZERO, sign));
218 : : }
219 : :
220 : 8 : FloatingPoint FloatingPoint::makeMinSubnormal(const FloatingPointSize& size,
221 : : bool sign)
222 : : {
223 [ + + ]: 16 : BitVector bvsign = sign ? BitVector::mkOne(1) : BitVector::mkZero(1);
224 : 16 : BitVector bvexp = BitVector::mkZero(size.packedExponentWidth());
225 : 8 : BitVector bvsig = BitVector::mkOne(size.packedSignificandWidth());
226 : 24 : return FloatingPoint(size, bvsign.concat(bvexp).concat(bvsig));
227 : : }
228 : :
229 : 8 : FloatingPoint FloatingPoint::makeMaxSubnormal(const FloatingPointSize& size,
230 : : bool sign)
231 : : {
232 [ + + ]: 16 : BitVector bvsign = sign ? BitVector::mkOne(1) : BitVector::mkZero(1);
233 : 16 : BitVector bvexp = BitVector::mkZero(size.packedExponentWidth());
234 : 8 : BitVector bvsig = BitVector::mkOnes(size.packedSignificandWidth());
235 : 24 : return FloatingPoint(size, bvsign.concat(bvexp).concat(bvsig));
236 : : }
237 : :
238 : 8 : FloatingPoint FloatingPoint::makeMinNormal(const FloatingPointSize& size,
239 : : bool sign)
240 : : {
241 [ + + ]: 16 : BitVector bvsign = sign ? BitVector::mkOne(1) : BitVector::mkZero(1);
242 : 16 : BitVector bvexp = BitVector::mkOne(size.packedExponentWidth());
243 : 8 : BitVector bvsig = BitVector::mkZero(size.packedSignificandWidth());
244 : 24 : return FloatingPoint(size, bvsign.concat(bvexp).concat(bvsig));
245 : : }
246 : :
247 : 8 : FloatingPoint FloatingPoint::makeMaxNormal(const FloatingPointSize& size,
248 : : bool sign)
249 : : {
250 [ + + ]: 16 : BitVector bvsign = sign ? BitVector::mkOne(1) : BitVector::mkZero(1);
251 : 16 : BitVector bvexp = BitVector::mkOnes(size.packedExponentWidth());
252 : 8 : bvexp.setBit(0, false);
253 : 8 : BitVector bvsig = BitVector::mkOnes(size.packedSignificandWidth());
254 : 24 : return FloatingPoint(size, bvsign.concat(bvexp).concat(bvsig));
255 : : }
256 : :
257 : 890 : FloatingPoint FloatingPoint::absolute(void) const
258 : : {
259 : 890 : return FloatingPoint(new FloatingPointLiteral(d_fpl->absolute()));
260 : : }
261 : :
262 : 1636 : FloatingPoint FloatingPoint::negate(void) const
263 : : {
264 : 1636 : return FloatingPoint(new FloatingPointLiteral(d_fpl->negate()));
265 : : }
266 : :
267 : 15615 : FloatingPoint FloatingPoint::add(const RoundingMode& rm,
268 : : const FloatingPoint& arg) const
269 : : {
270 : 15615 : return FloatingPoint(new FloatingPointLiteral(d_fpl->add(rm, *arg.d_fpl)));
271 : : }
272 : :
273 : 0 : FloatingPoint FloatingPoint::sub(const RoundingMode& rm,
274 : : const FloatingPoint& arg) const
275 : : {
276 : 0 : return FloatingPoint(new FloatingPointLiteral(d_fpl->sub(rm, *arg.d_fpl)));
277 : : }
278 : :
279 : 8144 : FloatingPoint FloatingPoint::mult(const RoundingMode& rm,
280 : : const FloatingPoint& arg) const
281 : : {
282 : 8144 : return FloatingPoint(new FloatingPointLiteral(d_fpl->mult(rm, *arg.d_fpl)));
283 : : }
284 : :
285 : 0 : FloatingPoint FloatingPoint::fma(const RoundingMode& rm,
286 : : const FloatingPoint& arg1,
287 : : const FloatingPoint& arg2) const
288 : : {
289 : : return FloatingPoint(
290 : 0 : new FloatingPointLiteral(d_fpl->fma(rm, *arg1.d_fpl, *arg2.d_fpl)));
291 : : }
292 : :
293 : 9258 : FloatingPoint FloatingPoint::div(const RoundingMode& rm,
294 : : const FloatingPoint& arg) const
295 : : {
296 : 9258 : return FloatingPoint(new FloatingPointLiteral(d_fpl->div(rm, *arg.d_fpl)));
297 : : }
298 : :
299 : 3558 : FloatingPoint FloatingPoint::sqrt(const RoundingMode& rm) const
300 : : {
301 : 3558 : return FloatingPoint(new FloatingPointLiteral(d_fpl->sqrt(rm)));
302 : : }
303 : :
304 : 280 : FloatingPoint FloatingPoint::rti(const RoundingMode& rm) const
305 : : {
306 : 280 : return FloatingPoint(new FloatingPointLiteral(d_fpl->rti(rm)));
307 : : }
308 : :
309 : 2060 : FloatingPoint FloatingPoint::rem(const FloatingPoint& arg) const
310 : : {
311 : 2060 : return FloatingPoint(new FloatingPointLiteral(d_fpl->rem(*arg.d_fpl)));
312 : : }
313 : :
314 : 424 : FloatingPoint FloatingPoint::maxTotal(const FloatingPoint& arg,
315 : : bool zeroCaseLeft) const
316 : : {
317 : : return FloatingPoint(
318 : 424 : new FloatingPointLiteral(d_fpl->maxTotal(*arg.d_fpl, zeroCaseLeft)));
319 : : }
320 : :
321 : 18 : FloatingPoint FloatingPoint::minTotal(const FloatingPoint& arg,
322 : : bool zeroCaseLeft) const
323 : : {
324 : : return FloatingPoint(
325 : 18 : new FloatingPointLiteral(d_fpl->minTotal(*arg.d_fpl, zeroCaseLeft)));
326 : : }
327 : :
328 : 212 : FloatingPoint::PartialFloatingPoint FloatingPoint::max(
329 : : const FloatingPoint& arg) const
330 : : {
331 : 212 : FloatingPoint tmp(maxTotal(arg, true));
332 : 636 : return PartialFloatingPoint(tmp, tmp == maxTotal(arg, false));
333 : : }
334 : :
335 : 8 : FloatingPoint::PartialFloatingPoint FloatingPoint::min(
336 : : const FloatingPoint& arg) const
337 : : {
338 : 8 : FloatingPoint tmp(minTotal(arg, true));
339 : 24 : return PartialFloatingPoint(tmp, tmp == minTotal(arg, false));
340 : : }
341 : :
342 : 81908 : bool FloatingPoint::operator==(const FloatingPoint& fp) const
343 : : {
344 : 81908 : return *d_fpl == *fp.d_fpl;
345 : : }
346 : :
347 : 44 : bool FloatingPoint::operator<=(const FloatingPoint& fp) const
348 : : {
349 : 44 : return *d_fpl <= *fp.d_fpl;
350 : : }
351 : :
352 : 2 : bool FloatingPoint::operator<(const FloatingPoint& fp) const
353 : : {
354 : 2 : return *d_fpl < *fp.d_fpl;
355 : : }
356 : :
357 : 19 : BitVector FloatingPoint::getExponent() const { return d_fpl->getExponent(); }
358 : :
359 : 19 : BitVector FloatingPoint::getSignificand() const
360 : : {
361 : 19 : return d_fpl->getSignificand();
362 : : }
363 : :
364 : 19 : bool FloatingPoint::getSign() const { return d_fpl->getSign(); }
365 : :
366 : 66 : bool FloatingPoint::isNormal(void) const { return d_fpl->isNormal(); }
367 : :
368 : 73 : bool FloatingPoint::isSubnormal(void) const { return d_fpl->isSubnormal(); }
369 : :
370 : 530 : bool FloatingPoint::isZero(void) const { return d_fpl->isZero(); }
371 : :
372 : 273 : bool FloatingPoint::isInfinite(void) const { return d_fpl->isInfinite(); }
373 : :
374 : 285 : bool FloatingPoint::isNaN(void) const { return d_fpl->isNaN(); }
375 : :
376 : 150 : bool FloatingPoint::isNegative(void) const { return d_fpl->isNegative(); }
377 : :
378 : 141 : bool FloatingPoint::isPositive(void) const { return d_fpl->isPositive(); }
379 : :
380 : 23 : FloatingPoint FloatingPoint::convert(const FloatingPointSize& target,
381 : : const RoundingMode& rm) const
382 : : {
383 : 23 : return FloatingPoint(new FloatingPointLiteral(d_fpl->convert(target, rm)));
384 : : }
385 : :
386 : 12 : BitVector FloatingPoint::convertToBVTotal(BitVectorSize width,
387 : : const RoundingMode& rm,
388 : : bool signedBV,
389 : : BitVector undefinedCase) const
390 : : {
391 [ + + ]: 12 : if (signedBV)
392 : : {
393 : 4 : return d_fpl->convertToSBVTotal(width, rm, undefinedCase);
394 : : }
395 : 8 : return d_fpl->convertToUBVTotal(width, rm, undefinedCase);
396 : : }
397 : :
398 : 33 : Rational FloatingPoint::convertToRationalTotal(Rational undefinedCase) const
399 : : {
400 : 66 : PartialRational p(convertToRational());
401 : :
402 [ + + ]: 66 : return p.second ? p.first : undefinedCase;
403 : : }
404 : :
405 : 0 : FloatingPoint::PartialBitVector FloatingPoint::convertToBV(
406 : : BitVectorSize width, const RoundingMode& rm, bool signedBV) const
407 : : {
408 : 0 : BitVector tmp(convertToBVTotal(width, rm, signedBV, BitVector(width, 0U)));
409 : : BitVector confirm(
410 : 0 : convertToBVTotal(width, rm, signedBV, BitVector(width, 1U)));
411 : :
412 : 0 : return PartialBitVector(tmp, tmp == confirm);
413 : : }
414 : :
415 : 41 : FloatingPoint::PartialRational FloatingPoint::convertToRational(void) const
416 : : {
417 [ + + ][ - + ]: 41 : if (isNaN() || isInfinite())
[ + + ]
418 : : {
419 : 20 : return PartialRational(Rational(0U, 1U), false);
420 : : }
421 [ + + ]: 31 : if (isZero())
422 : : {
423 : 8 : return PartialRational(Rational(0U, 1U), true);
424 : : }
425 [ + + ]: 54 : Integer sign((d_fpl->getSign()) ? -1 : 1);
426 : : Integer exp(
427 : 54 : d_fpl->getExponent().toSignedInteger()
428 : 27 : - (Integer(d_fpl->getSize().significandWidth()
429 : 54 : - 1))); // -1 as forcibly normalised into the [1,2) range
430 : 54 : Integer significand(d_fpl->getSignificand().toInteger());
431 : 54 : Integer signedSignificand(sign * significand);
432 : :
433 : : // We only have multiplyByPow(uint32_t) so we can't convert all numbers.
434 : : // As we convert Integer -> unsigned int -> uint32_t we need that
435 : : // unsigned int is not smaller than uint32_t
436 : : static_assert(sizeof(unsigned int) >= sizeof(uint32_t),
437 : : "Conversion float -> real could loose data");
438 : : #ifdef CVC5_ASSERTIONS
439 : : // Note that multipling by 2^n requires n bits of space (worst case)
440 : : // so, in effect, these tests limit us to cases where the resultant
441 : : // number requires up to 2^32 bits = 512 megabyte to represent.
442 : 54 : Integer shiftLimit(std::numeric_limits<uint32_t>::max());
443 : : #endif
444 : :
445 [ + + ]: 27 : if (!(exp.strictlyNegative()))
446 : : {
447 [ - + ][ - + ]: 19 : Assert(exp <= shiftLimit);
[ - - ]
448 : 19 : Integer r(signedSignificand.multiplyByPow2(exp.toUnsignedInt()));
449 : 38 : return PartialRational(Rational(r), true);
450 : : }
451 : 16 : Integer one(1U);
452 [ - + ][ - + ]: 8 : Assert((-exp) <= shiftLimit);
[ - - ]
453 : 16 : Integer q(one.multiplyByPow2((-exp).toUnsignedInt()));
454 : 8 : Rational r(signedSignificand, q);
455 : 8 : return PartialRational(r, true);
456 : : }
457 : :
458 : 113332 : BitVector FloatingPoint::pack(void) const { return d_fpl->pack(); }
459 : :
460 : 519 : void FloatingPoint::getIEEEBitvectors(BitVector& sign,
461 : : BitVector& exp,
462 : : BitVector& sig) const
463 : : {
464 : : // retrieve BV value
465 : 519 : BitVector bv(pack());
466 : : uint32_t largestSignificandBit =
467 : 519 : getSize().significandWidth() - 2; // -1 for -inclusive, -1 for hidden
468 : : uint32_t largestExponentBit =
469 : 519 : (getSize().exponentWidth() - 1) + (largestSignificandBit + 1);
470 : 519 : sign = bv.extract(largestExponentBit + 1, largestExponentBit + 1);
471 : 519 : exp = bv.extract(largestExponentBit, largestSignificandBit + 1);
472 : 519 : sig = bv.extract(largestSignificandBit, 0);
473 : 519 : }
474 : :
475 : 500 : std::string FloatingPoint::toString(bool printAsIndexed) const
476 : : {
477 : 500 : std::string str;
478 : : // retrive BV value
479 [ + + ][ + + ]: 4000 : BitVector v[3];
[ - - ]
480 : 500 : getIEEEBitvectors(v[0], v[1], v[2]);
481 : 500 : str.append("(fp ");
482 [ + + ]: 2000 : for (uint32_t i = 0; i < 3; ++i)
483 : : {
484 [ - + ]: 1500 : if (printAsIndexed)
485 : : {
486 : 0 : str.append("(_ bv");
487 : 0 : str.append(v[i].getValue().toString());
488 : 0 : str.append(" ");
489 : 0 : str.append(std::to_string(v[i].getSize()));
490 : 0 : str.append(")");
491 : : }
492 : : else
493 : : {
494 : 1500 : str.append("#b");
495 : 1500 : str.append(v[i].toString());
496 : : }
497 [ + + ]: 1500 : if (i < 2)
498 : : {
499 : 1000 : str.append(" ");
500 : : }
501 : : }
502 : 500 : str.append(")");
503 : 1000 : return str;
504 : : }
505 : :
506 : 0 : std::ostream& operator<<(std::ostream& os, const FloatingPoint& fp)
507 : : {
508 : : // print in binary notation
509 : 0 : return os << fp.toString();
510 : : }
511 : :
512 : 0 : std::ostream& operator<<(std::ostream& os, const FloatingPointSize& fps)
513 : : {
514 : 0 : return os << "(_ FloatingPoint " << fps.exponentWidth() << " "
515 : 0 : << fps.significandWidth() << ")";
516 : : }
517 : :
518 : 0 : std::ostream& operator<<(std::ostream& os, const FloatingPointConvertSort& fpcs)
519 : : {
520 : 0 : return os << "(_ to_fp " << fpcs.getSize().exponentWidth() << " "
521 : 0 : << fpcs.getSize().significandWidth() << ")";
522 : : }
523 : : } // namespace cvc5::internal
|
