Theory of Floating-Points

examples/api/cpp/floating_point_arith.cpp

  1/******************************************************************************
  2 * Top contributors (to current version):
  3 *   Mudathir Mohamed, Mathias Preiner, Andres Noetzli
  4 *
  5 * This file is part of the cvc5 project.
  6 *
  7 * Copyright (c) 2009-2022 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 * An example of solving floating-point problems with cvc5's cpp API
 14 *
 15 * This example shows to create floating-point types, variables and expressions,
 16 * and how to create rounding mode constants by solving toy problems. The
 17 * example also shows making special values (such as NaN and +oo) and converting
 18 * an IEEE 754-2008 bit-vector to a floating-point number.
 19 */
 20
 21#include <cvc5/cvc5.h>
 22
 23#include <iostream>
 24#include <cassert>
 25
 26using namespace std;
 27using namespace cvc5;
 28
 29int main()
 30{
 31  Solver solver;
 32  solver.setOption("produce-models", "true");
 33
 34  // Make single precision floating-point variables
 35  Sort fpt32 = solver.mkFloatingPointSort(8, 24);
 36  Term a = solver.mkConst(fpt32, "a");
 37  Term b = solver.mkConst(fpt32, "b");
 38  Term c = solver.mkConst(fpt32, "c");
 39  Term d = solver.mkConst(fpt32, "d");
 40  Term e = solver.mkConst(fpt32, "e");
 41
 42  // Assert that floating-point addition is not associative:
 43  // (a + (b + c)) != ((a + b) + c)
 44  Term rm = solver.mkRoundingMode(RoundingMode::ROUND_NEAREST_TIES_TO_EVEN);
 45  Term lhs = solver.mkTerm(
 46      Kind::FLOATINGPOINT_ADD,
 47      {rm, a, solver.mkTerm(Kind::FLOATINGPOINT_ADD, {rm, b, c})});
 48  Term rhs = solver.mkTerm(
 49      Kind::FLOATINGPOINT_ADD,
 50      {rm, solver.mkTerm(Kind::FLOATINGPOINT_ADD, {rm, a, b}), c});
 51  solver.assertFormula(
 52      solver.mkTerm(Kind::NOT, {solver.mkTerm(Kind::EQUAL, {a, b})}));
 53
 54  Result r = solver.checkSat();  // result is sat
 55  assert (r.isSat());
 56
 57  cout << "a = " << solver.getValue(a) << endl;
 58  cout << "b = " << solver.getValue(b) << endl;
 59  cout << "c = " << solver.getValue(c) << endl;
 60
 61  // Now, let's restrict `a` to be either NaN or positive infinity
 62  Term nan = solver.mkFloatingPointNaN(8, 24);
 63  Term inf = solver.mkFloatingPointPosInf(8, 24);
 64  solver.assertFormula(solver.mkTerm(Kind::OR,
 65                                     {solver.mkTerm(Kind::EQUAL, {a, inf}),
 66                                      solver.mkTerm(Kind::EQUAL, {a, nan})}));
 67
 68  r = solver.checkSat();  // result is sat
 69  assert (r.isSat());
 70
 71  cout << "a = " << solver.getValue(a) << endl;
 72  cout << "b = " << solver.getValue(b) << endl;
 73  cout << "c = " << solver.getValue(c) << endl;
 74
 75  // And now for something completely different. Let's try to find a (normal)
 76  // floating-point number that rounds to different integer values for
 77  // different rounding modes.
 78  Term rtp = solver.mkRoundingMode(RoundingMode::ROUND_TOWARD_POSITIVE);
 79  Term rtn = solver.mkRoundingMode(RoundingMode::ROUND_TOWARD_NEGATIVE);
 80  Op op = solver.mkOp(Kind::FLOATINGPOINT_TO_SBV, {16});  // (_ fp.to_sbv 16)
 81  lhs = solver.mkTerm(op, {rtp, d});
 82  rhs = solver.mkTerm(op, {rtn, d});
 83  solver.assertFormula(solver.mkTerm(Kind::FLOATINGPOINT_IS_NORMAL, {d}));
 84  solver.assertFormula(
 85      solver.mkTerm(Kind::NOT, {solver.mkTerm(Kind::EQUAL, {lhs, rhs})}));
 86
 87  r = solver.checkSat();  // result is sat
 88  assert (r.isSat());
 89
 90  // Convert the result to a rational and print it
 91  Term val = solver.getValue(d);
 92  Term realVal =
 93      solver.getValue(solver.mkTerm(Kind::FLOATINGPOINT_TO_REAL, {val}));
 94  cout << "d = " << val << " = " << realVal << endl;
 95  cout << "((_ fp.to_sbv 16) RTP d) = " << solver.getValue(lhs) << endl;
 96  cout << "((_ fp.to_sbv 16) RTN d) = " << solver.getValue(rhs) << endl;
 97
 98  // For our final trick, let's try to find a floating-point number between
 99  // positive zero and the smallest positive floating-point number
100  Term zero = solver.mkFloatingPointPosZero(8, 24);
101  Term smallest = solver.mkFloatingPoint(8, 24, solver.mkBitVector(32, 0b001));
102  solver.assertFormula(
103      solver.mkTerm(Kind::AND,
104                    {solver.mkTerm(Kind::FLOATINGPOINT_LT, {zero, e}),
105                     solver.mkTerm(Kind::FLOATINGPOINT_LT, {e, smallest})}));
106
107  r = solver.checkSat();  // result is unsat
108  assert (!r.isSat());
109}