From 151d6272fdab2801e891fdcd1ca8ba876ba613c2 Mon Sep 17 00:00:00 2001 From: Daniel Kroening Date: Sat, 21 Dec 2024 14:04:37 +0000 Subject: [PATCH] Introduce floatbv_round_to_integral_exprt This adds a new expression, floatbv_round_to_integral, which rounds an IEEE 754 floating-point number given as bit-vector to the nearest integer, considering the explicitly given rounding mode. --- .../main.c | 9 - .../test.desc | 8 - .../main.c | 9 - .../test.desc | 8 - .../main.c | 9 - .../test.desc | 8 - .../smt2_solver/fp/roundToIntegral1.smt2 | 76 +++++ src/ansi-c/c_typecheck_expr.cpp | 18 + src/ansi-c/cprover_builtin_headers.h | 3 + src/ansi-c/library/math.c | 320 ++++-------------- src/solvers/flattening/boolbv.cpp | 3 + src/solvers/flattening/boolbv.h | 3 + src/solvers/flattening/boolbv_floatbv_op.cpp | 18 + src/solvers/floatbv/float_utils.cpp | 49 ++- src/solvers/floatbv/float_utils.h | 8 +- src/solvers/smt2/smt2_conv.cpp | 21 ++ src/solvers/smt2/smt2_conv.h | 3 + src/solvers/smt2/smt2_parser.cpp | 13 + src/util/floatbv_expr.h | 64 ++++ src/util/ieee_float.cpp | 63 +++- src/util/ieee_float.h | 91 +++-- src/util/irep_ids.def | 1 + unit/solvers/floatbv/float_utils.cpp | 88 +++++ unit/util/ieee_float.cpp | 70 +++- 24 files changed, 622 insertions(+), 341 deletions(-) delete mode 100644 regression/cbmc-library/__sort_of_CPROVER_round_to_integral-01/main.c delete mode 100644 regression/cbmc-library/__sort_of_CPROVER_round_to_integral-01/test.desc delete mode 100644 regression/cbmc-library/__sort_of_CPROVER_round_to_integralf-01/main.c delete mode 100644 regression/cbmc-library/__sort_of_CPROVER_round_to_integralf-01/test.desc delete mode 100644 regression/cbmc-library/__sort_of_CPROVER_round_to_integrall-01/main.c delete mode 100644 regression/cbmc-library/__sort_of_CPROVER_round_to_integrall-01/test.desc create mode 100644 regression/smt2_solver/fp/roundToIntegral1.smt2 diff --git a/regression/cbmc-library/__sort_of_CPROVER_round_to_integral-01/main.c b/regression/cbmc-library/__sort_of_CPROVER_round_to_integral-01/main.c deleted file mode 100644 index fcf4ed3ca10..00000000000 --- a/regression/cbmc-library/__sort_of_CPROVER_round_to_integral-01/main.c +++ /dev/null @@ -1,9 +0,0 @@ -#include -#include - -int main() -{ - __sort_of_CPROVER_round_to_integral(); - assert(0); - return 0; -} diff --git a/regression/cbmc-library/__sort_of_CPROVER_round_to_integral-01/test.desc b/regression/cbmc-library/__sort_of_CPROVER_round_to_integral-01/test.desc deleted file mode 100644 index 9542d988e8d..00000000000 --- a/regression/cbmc-library/__sort_of_CPROVER_round_to_integral-01/test.desc +++ /dev/null @@ -1,8 +0,0 @@ -KNOWNBUG -main.c ---pointer-check --bounds-check -^EXIT=0$ -^SIGNAL=0$ -^VERIFICATION SUCCESSFUL$ --- -^warning: ignoring diff --git a/regression/cbmc-library/__sort_of_CPROVER_round_to_integralf-01/main.c b/regression/cbmc-library/__sort_of_CPROVER_round_to_integralf-01/main.c deleted file mode 100644 index c2886c33a89..00000000000 --- a/regression/cbmc-library/__sort_of_CPROVER_round_to_integralf-01/main.c +++ /dev/null @@ -1,9 +0,0 @@ -#include -#include - -int main() -{ - __sort_of_CPROVER_round_to_integralf(); - assert(0); - return 0; -} diff --git a/regression/cbmc-library/__sort_of_CPROVER_round_to_integralf-01/test.desc b/regression/cbmc-library/__sort_of_CPROVER_round_to_integralf-01/test.desc deleted file mode 100644 index 9542d988e8d..00000000000 --- a/regression/cbmc-library/__sort_of_CPROVER_round_to_integralf-01/test.desc +++ /dev/null @@ -1,8 +0,0 @@ -KNOWNBUG -main.c ---pointer-check --bounds-check -^EXIT=0$ -^SIGNAL=0$ -^VERIFICATION SUCCESSFUL$ --- -^warning: ignoring diff --git a/regression/cbmc-library/__sort_of_CPROVER_round_to_integrall-01/main.c b/regression/cbmc-library/__sort_of_CPROVER_round_to_integrall-01/main.c deleted file mode 100644 index 5a8db444a8e..00000000000 --- a/regression/cbmc-library/__sort_of_CPROVER_round_to_integrall-01/main.c +++ /dev/null @@ -1,9 +0,0 @@ -#include -#include - -int main() -{ - __sort_of_CPROVER_round_to_integrall(); - assert(0); - return 0; -} diff --git a/regression/cbmc-library/__sort_of_CPROVER_round_to_integrall-01/test.desc b/regression/cbmc-library/__sort_of_CPROVER_round_to_integrall-01/test.desc deleted file mode 100644 index 9542d988e8d..00000000000 --- a/regression/cbmc-library/__sort_of_CPROVER_round_to_integrall-01/test.desc +++ /dev/null @@ -1,8 +0,0 @@ -KNOWNBUG -main.c ---pointer-check --bounds-check -^EXIT=0$ -^SIGNAL=0$ -^VERIFICATION SUCCESSFUL$ --- -^warning: ignoring diff --git a/regression/smt2_solver/fp/roundToIntegral1.smt2 b/regression/smt2_solver/fp/roundToIntegral1.smt2 new file mode 100644 index 00000000000..c51496a5a38 --- /dev/null +++ b/regression/smt2_solver/fp/roundToIntegral1.smt2 @@ -0,0 +1,76 @@ +(set-logic FP) + +(define-fun zero () (_ FloatingPoint 11 53) ((_ to_fp 11 53) RNE 0)) +(define-fun minus-zero () (_ FloatingPoint 11 53) (- zero)) +(define-fun one () (_ FloatingPoint 11 53) ((_ to_fp 11 53) RNE 1)) +(define-fun minus-one () (_ FloatingPoint 11 53) (- one)) +(define-fun zero-point-one () (_ FloatingPoint 11 53) ((_ to_fp 11 53) RNE 0.1)) +(define-fun minus-zero-point-one () (_ FloatingPoint 11 53) (- zero-point-one)) +(define-fun ten-point-one () (_ FloatingPoint 11 53) ((_ to_fp 11 53) RNE 10.1)) +(define-fun minus-ten-point-one () (_ FloatingPoint 11 53) (- ten-point-one)) +(define-fun ten () (_ FloatingPoint 11 53) ((_ to_fp 11 53) RNE 10)) +(define-fun minus-ten () (_ FloatingPoint 11 53) (- ten)) +(define-fun eleven () (_ FloatingPoint 11 53) ((_ to_fp 11 53) RNE 11)) +(define-fun minus-eleven () (_ FloatingPoint 11 53) (- eleven)) + +(assert (not (and + + ; round up + (= (fp.roundToIntegral RTP (_ NaN 11 53)) (_ NaN 11 53)) + (= (fp.roundToIntegral RTP (_ +oo 11 53)) (_ +oo 11 53)) + (= (fp.roundToIntegral RTP (_ -oo 11 53)) (_ -oo 11 53)) + (= (fp.roundToIntegral RTP zero) zero) + (= (fp.roundToIntegral RTP minus-zero) minus-zero) + (= (fp.roundToIntegral RTP one) one) + (= (fp.roundToIntegral RTP zero-point-one) one) + (= (fp.roundToIntegral RTP minus-zero-point-one) minus-zero) + (= (fp.roundToIntegral RTP ten-point-one) eleven) + (= (fp.roundToIntegral RTP minus-ten-point-one) minus-ten) + ;(= (fp.roundToIntegral RTP ((_ to_fp 11 53) RTN 0x1.0p+52)) ((_ to_fp 11 53) RTN 0x1.0p+52)) + ;(= (fp.roundToIntegral RTP dmax) dmax) + + ; round down + (= (fp.roundToIntegral RTN (_ NaN 11 53)) (_ NaN 11 53)) + (= (fp.roundToIntegral RTN (_ +oo 11 53)) (_ +oo 11 53)) + (= (fp.roundToIntegral RTN (_ -oo 11 53)) (_ -oo 11 53)) + (= (fp.roundToIntegral RTN zero) zero) + (= (fp.roundToIntegral RTN minus-zero) minus-zero) + (= (fp.roundToIntegral RTN one) one) + (= (fp.roundToIntegral RTN zero-point-one) zero) + (= (fp.roundToIntegral RTN minus-zero-point-one) minus-one) + (= (fp.roundToIntegral RTN ten-point-one) ten) + (= (fp.roundToIntegral RTN minus-ten-point-one) minus-eleven) + ;(= (fp.roundToIntegral RTN ((_ to_fp 11 53) RTN 0x1.0p+52)) ((_ to_fp 11 53) RTN 0x1.0p+52)) + ;(= (fp.roundToIntegral RTN dmax)) dmax) + + ; round to nearest ties to even + (= (fp.roundToIntegral RNE (_ NaN 11 53)) (_ NaN 11 53)) + (= (fp.roundToIntegral RNE (_ +oo 11 53)) (_ +oo 11 53)) + (= (fp.roundToIntegral RNE (_ -oo 11 53)) (_ -oo 11 53)) + (= (fp.roundToIntegral RNE zero) zero) + (= (fp.roundToIntegral RNE minus-zero) minus-zero) + (= (fp.roundToIntegral RNE one) one) + (= (fp.roundToIntegral RNE zero-point-one) zero) + (= (fp.roundToIntegral RNE minus-zero-point-one) minus-zero) + (= (fp.roundToIntegral RNE ten-point-one) ten) + (= (fp.roundToIntegral RNE minus-ten-point-one) minus-ten) + ;(= (fp.roundToIntegral RNE ((_ to_fp 11 53) RTN 0x1.0p+52)) ((_ to_fp 11 53) RTN 0x1.0p+52)) + ;(= (fp.roundToIntegral RNE dmax)) dmax) + + ; round to zero + (= (fp.roundToIntegral RTZ (_ NaN 11 53)) (_ NaN 11 53)) + (= (fp.roundToIntegral RTZ (_ +oo 11 53)) (_ +oo 11 53)) + (= (fp.roundToIntegral RTZ (_ -oo 11 53)) (_ -oo 11 53)) + (= (fp.roundToIntegral RTZ zero) zero) + (= (fp.roundToIntegral RTZ minus-zero) minus-zero) + (= (fp.roundToIntegral RTZ one) one) + (= (fp.roundToIntegral RTZ zero-point-one) zero) + (= (fp.roundToIntegral RTZ minus-zero-point-one) minus-zero) + (= (fp.roundToIntegral RTZ ten-point-one) ten) + (= (fp.roundToIntegral RTZ minus-ten-point-one) minus-ten) + ;(= (fp.roundToIntegral RTZ ((_ to_fp 11 53) RTN 0x1.0p+52)) ((_ to_fp 11 53) RTN 0x1.0p+52)) + ;(= (fp.roundToIntegral RTZ dmax)) dmax) +))) + +; should be unsat +(check-sat) diff --git a/src/ansi-c/c_typecheck_expr.cpp b/src/ansi-c/c_typecheck_expr.cpp index 35bc45cb257..942c1ebbb70 100644 --- a/src/ansi-c/c_typecheck_expr.cpp +++ b/src/ansi-c/c_typecheck_expr.cpp @@ -3237,6 +3237,24 @@ exprt c_typecheck_baset::do_special_functions( return std::move(infl_expr); } + else if( + identifier == CPROVER_PREFIX "round_to_integralf" || + identifier == CPROVER_PREFIX "round_to_integrald" || + identifier == CPROVER_PREFIX "round_to_integralld") + { + if(expr.arguments().size() != 2) + { + error().source_location = f_op.source_location(); + error() << identifier << " expects two arguments" << eom; + throw 0; + } + + auto round_to_integral_expr = + floatbv_round_to_integral_exprt{expr.arguments()[0], expr.arguments()[1]}; + round_to_integral_expr.add_source_location() = source_location; + + return std::move(round_to_integral_expr); + } else if( identifier == CPROVER_PREFIX "abs" || identifier == CPROVER_PREFIX "labs" || identifier == CPROVER_PREFIX "llabs" || diff --git a/src/ansi-c/cprover_builtin_headers.h b/src/ansi-c/cprover_builtin_headers.h index f0ae5e74004..998e81f0445 100644 --- a/src/ansi-c/cprover_builtin_headers.h +++ b/src/ansi-c/cprover_builtin_headers.h @@ -93,6 +93,9 @@ int __CPROVER_islessgreaterf(float f, float g); int __CPROVER_islessgreaterd(double f, double g); int __CPROVER_isunorderedf(float f, float g); int __CPROVER_isunorderedd(double f, double g); +float __CPROVER_round_to_integralf(float, int); +double __CPROVER_round_to_integrald(double, int); +long double __CPROVER_round_to_integralld(long double, int); // absolute value int __CPROVER_abs(int x); diff --git a/src/ansi-c/library/math.c b/src/ansi-c/library/math.c index 6d4fc796b24..dcf0bae334b 100644 --- a/src/ansi-c/library/math.c +++ b/src/ansi-c/library/math.c @@ -1244,126 +1244,6 @@ float fdimf(float f, float g) { return ((f > g) ? f - g : +0.0f); } long double fdiml(long double f, long double g) { return ((f > g) ? f - g : +0.0); } - - -/* FUNCTION: __sort_of_CPROVER_round_to_integral */ -// TODO : Should be a real __CPROVER function to convert to SMT-LIB - -#ifndef __CPROVER_MATH_H_INCLUDED -#include -#define __CPROVER_MATH_H_INCLUDED -#endif - -#ifndef __CPROVER_FENV_H_INCLUDED -#include -#define __CPROVER_FENV_H_INCLUDED -#endif - -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d) -{ - double magicConst = 0x1.0p+52; - double return_value; - int saved_rounding_mode = fegetround(); - fesetround(rounding_mode); - - if (fabs(d) >= magicConst || d == 0.0) - { - return_value = d; - } - else - { - if (!signbit(d)) { - double tmp = d + magicConst; - return_value = tmp - magicConst; - } else { - double tmp = d - magicConst; - return_value = tmp + magicConst; - } - } - - fesetround(saved_rounding_mode); - return return_value; -} - -/* FUNCTION: __sort_of_CPROVER_round_to_integralf */ -// TODO : Should be a real __CPROVER function to convert to SMT-LIB - -#ifndef __CPROVER_MATH_H_INCLUDED -#include -#define __CPROVER_MATH_H_INCLUDED -#endif - -#ifndef __CPROVER_FENV_H_INCLUDED -#include -#define __CPROVER_FENV_H_INCLUDED -#endif - -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d) -{ - float magicConst = 0x1.0p+23f; // 23 is significant - float return_value; - int saved_rounding_mode = fegetround(); - fesetround(rounding_mode); - - if (fabsf(d) >= magicConst || d == 0.0) - { - return_value = d; - } - else - { - if (!signbit(d)) { - float tmp = d + magicConst; - return_value = tmp - magicConst; - } else { - float tmp = d - magicConst; - return_value = tmp + magicConst; - } - } - - fesetround(saved_rounding_mode); - return return_value; -} - - -/* FUNCTION: __sort_of_CPROVER_round_to_integrall */ -// TODO : Should be a real __CPROVER function to convert to SMT-LIB - -#ifndef __CPROVER_MATH_H_INCLUDED -#include -#define __CPROVER_MATH_H_INCLUDED -#endif - -#ifndef __CPROVER_FENV_H_INCLUDED -#include -#define __CPROVER_FENV_H_INCLUDED -#endif - -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d) -{ - long double magicConst = 0x1.0p+64; - long double return_value; - int saved_rounding_mode = fegetround(); - fesetround(rounding_mode); - - if (fabsl(d) >= magicConst || d == 0.0) - { - return_value = d; - } - else - { - if (!signbit(d)) { - long double tmp = d + magicConst; - return_value = tmp - magicConst; - } else { - long double tmp = d - magicConst; - return_value = tmp + magicConst; - } - } - - fesetround(saved_rounding_mode); - return return_value; -} - /* ISO 9899:2011 * * The ceil functions compute the smallest integer value not less than @@ -1382,11 +1262,9 @@ long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); - double ceil(double x) { - return __sort_of_CPROVER_round_to_integral(FE_UPWARD, x); + return __CPROVER_round_to_integrald(x, 2); // FE_UPWARD } /* FUNCTION: ceilf */ @@ -1401,11 +1279,9 @@ double ceil(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); - float ceilf(float x) { - return __sort_of_CPROVER_round_to_integralf(FE_UPWARD, x); + return __CPROVER_round_to_integralf(x, 2); // FE_UPWARD } @@ -1421,11 +1297,9 @@ float ceilf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); - long double ceill(long double x) { - return __sort_of_CPROVER_round_to_integrall(FE_UPWARD, x); + return __CPROVER_round_to_integralld(x, 2); // FE_UPWARD } @@ -1446,11 +1320,9 @@ long double ceill(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); - double floor(double x) { - return __sort_of_CPROVER_round_to_integral(FE_DOWNWARD, x); + return __CPROVER_round_to_integrald(x, 3); // FE_DOWNWARD } /* FUNCTION: floorf */ @@ -1465,11 +1337,9 @@ double floor(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); - float floorf(float x) { - return __sort_of_CPROVER_round_to_integralf(FE_DOWNWARD, x); + return __CPROVER_round_to_integralf(x, 3); // FE_DOWNWARD } @@ -1485,11 +1355,9 @@ float floorf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); - long double floorl(long double x) { - return __sort_of_CPROVER_round_to_integrall(FE_DOWNWARD, x); + return __CPROVER_round_to_integralld(x, 3); // FE_DOWNWARD } @@ -1511,11 +1379,9 @@ long double floorl(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); - double trunc(double x) { - return __sort_of_CPROVER_round_to_integral(FE_TOWARDZERO, x); + return __CPROVER_round_to_integrald(x, 0); // FE_TOWARDZERO } /* FUNCTION: truncf */ @@ -1530,11 +1396,9 @@ double trunc(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); - float truncf(float x) { - return __sort_of_CPROVER_round_to_integralf(FE_TOWARDZERO, x); + return __CPROVER_round_to_integralf(x, 0); // FE_TOWARDZERO } @@ -1550,11 +1414,9 @@ float truncf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); - long double truncl(long double x) { - return __sort_of_CPROVER_round_to_integrall(FE_TOWARDZERO, x); + return __CPROVER_round_to_integralld(x, 0); // FE_TOWARDZERO } @@ -1576,12 +1438,10 @@ long double truncl(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); - double round(double x) { // Tempting but RNE not RNA - // return __sort_of_CPROVER_round_to_integral(FE_TONEAREST, x); + // return __CPROVER_round_to_integrald(x, FE_TONEAREST); int saved_rounding_mode = fegetround(); fesetround(FE_TOWARDZERO); @@ -1596,8 +1456,8 @@ double round(double x) } fesetround(saved_rounding_mode); - - return __sort_of_CPROVER_round_to_integral(FE_TOWARDZERO, xp); + + return __CPROVER_round_to_integrald(xp, 0); // FE_TOWARDZERO } /* FUNCTION: roundf */ @@ -1612,12 +1472,10 @@ double round(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); - float roundf(float x) { // Tempting but RNE not RNA - // return __sort_of_CPROVER_round_to_integralf(FE_TONEAREST, x); + // return __CPROVER_round_to_integralf(x, FE_TONEAREST); int saved_rounding_mode = fegetround(); fesetround(FE_TOWARDZERO); @@ -1632,8 +1490,8 @@ float roundf(float x) } fesetround(saved_rounding_mode); - - return __sort_of_CPROVER_round_to_integralf(FE_TOWARDZERO, xp); + + return __CPROVER_round_to_integralf(xp, 0); // FE_TOWARDZERO } @@ -1649,12 +1507,10 @@ float roundf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); - long double roundl(long double x) { // Tempting but RNE not RNA - // return __sort_of_CPROVER_round_to_integrall(FE_TONEAREST, x); + // return __CPROVER_round_to_integralld(x, FE_TONEAREST); int saved_rounding_mode = fegetround(); fesetround(FE_TOWARDZERO); @@ -1669,8 +1525,8 @@ long double roundl(long double x) } fesetround(saved_rounding_mode); - - return __sort_of_CPROVER_round_to_integrall(FE_TOWARDZERO, xp); + + return __CPROVER_round_to_integralld(xp, 0); // FE_TOWARDZERO } @@ -1694,11 +1550,11 @@ long double roundl(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); +extern int __CPROVER_rounding_mode; double nearbyint(double x) { - return __sort_of_CPROVER_round_to_integral(fegetround(), x); + return __CPROVER_round_to_integrald(x, __CPROVER_rounding_mode); } /* FUNCTION: nearbyintf */ @@ -1713,11 +1569,11 @@ double nearbyint(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); +extern int __CPROVER_rounding_mode; float nearbyintf(float x) { - return __sort_of_CPROVER_round_to_integralf(fegetround(), x); + return __CPROVER_round_to_integralf(x, __CPROVER_rounding_mode); } @@ -1733,11 +1589,11 @@ float nearbyintf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); +extern int __CPROVER_rounding_mode; long double nearbyintl(long double x) { - return __sort_of_CPROVER_round_to_integrall(fegetround(), x); + return __CPROVER_round_to_integralld(x, __CPROVER_rounding_mode); } @@ -1761,11 +1617,11 @@ long double nearbyintl(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); +extern int __CPROVER_rounding_mode; double rint(double x) { - return __sort_of_CPROVER_round_to_integral(fegetround(), x); + return __CPROVER_round_to_integrald(x, __CPROVER_rounding_mode); } /* FUNCTION: rintf */ @@ -1780,11 +1636,11 @@ double rint(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); +extern int __CPROVER_rounding_mode; float rintf(float x) { - return __sort_of_CPROVER_round_to_integralf(fegetround(), x); + return __CPROVER_round_to_integralf(x, __CPROVER_rounding_mode); } /* FUNCTION: rintl */ @@ -1799,11 +1655,11 @@ float rintf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); +extern int __CPROVER_rounding_mode; long double rintl(long double x) { - return __sort_of_CPROVER_round_to_integrall(fegetround(), x); + return __CPROVER_round_to_integralld(x, __CPROVER_rounding_mode); } @@ -1829,13 +1685,11 @@ long double rintl(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); +extern int __CPROVER_rounding_mode; long int lrint(double x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT - double rti = __sort_of_CPROVER_round_to_integral(fegetround(), x); + double rti = __CPROVER_round_to_integrald(x, __CPROVER_rounding_mode); return (long int)rti; } @@ -1851,13 +1705,11 @@ long int lrint(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); +extern int __CPROVER_rounding_mode; long int lrintf(float x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT - float rti = __sort_of_CPROVER_round_to_integralf(fegetround(), x); + float rti = __CPROVER_round_to_integralf(x, __CPROVER_rounding_mode); return (long int)rti; } @@ -1874,13 +1726,11 @@ long int lrintf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); +extern int __CPROVER_rounding_mode; long int lrintl(long double x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT - long double rti = __sort_of_CPROVER_round_to_integrall(fegetround(), x); + long double rti = __CPROVER_round_to_integralld(x, __CPROVER_rounding_mode); return (long int)rti; } @@ -1897,13 +1747,11 @@ long int lrintl(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); +extern int __CPROVER_rounding_mode; long long int llrint(double x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT - double rti = __sort_of_CPROVER_round_to_integral(fegetround(), x); + double rti = __CPROVER_round_to_integrald(x, __CPROVER_rounding_mode); return (long long int)rti; } @@ -1919,13 +1767,11 @@ long long int llrint(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); +extern int __CPROVER_rounding_mode; long long int llrintf(float x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT - float rti = __sort_of_CPROVER_round_to_integralf(fegetround(), x); + float rti = __CPROVER_round_to_integralf(x, __CPROVER_rounding_mode); return (long long int)rti; } @@ -1942,13 +1788,11 @@ long long int llrintf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); +extern int __CPROVER_rounding_mode; long long int llrintl(long double x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT - long double rti = __sort_of_CPROVER_round_to_integrall(fegetround(), x); + long double rti = __CPROVER_round_to_integralld(x, __CPROVER_rounding_mode); return (long long int)rti; } @@ -1974,12 +1818,9 @@ long long int llrintl(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); - long int lround(double x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT, plus should use RNA + // TODO: should use RNA int saved_rounding_mode = fegetround(); fesetround(FE_TOWARDZERO); @@ -1994,8 +1835,8 @@ long int lround(double x) } fesetround(saved_rounding_mode); - - double rti = __sort_of_CPROVER_round_to_integral(FE_TOWARDZERO, xp); + + double rti = __CPROVER_round_to_integrald(xp, 0); // FE_TOWARDZERO return (long int)rti; } @@ -2011,12 +1852,9 @@ long int lround(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); - long int lroundf(float x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT, plus should use RNA + // should use RNA int saved_rounding_mode = fegetround(); fesetround(FE_TOWARDZERO); @@ -2030,8 +1868,8 @@ long int lroundf(float x) } fesetround(saved_rounding_mode); - - float rti = __sort_of_CPROVER_round_to_integralf(FE_TOWARDZERO, xp); + + float rti = __CPROVER_round_to_integralf(xp, 0); // FE_TOWARDZERO return (long int)rti; } @@ -2048,15 +1886,12 @@ long int lroundf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); - long int lroundl(long double x) { int saved_rounding_mode = fegetround(); fesetround(FE_TOWARDZERO); - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT, plus should use RNA + // TODO: should use RNA long double xp; if (x > 0.0) { xp = x + 0.5; @@ -2067,8 +1902,8 @@ long int lroundl(long double x) } fesetround(saved_rounding_mode); - - long double rti = __sort_of_CPROVER_round_to_integrall(FE_TOWARDZERO, xp); + + long double rti = __CPROVER_round_to_integralld(xp, 0); // FE_TOWARDZERO return (long int)rti; } @@ -2085,8 +1920,6 @@ long int lroundl(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); - long long int llround(double x) { // TODO : should be an all-in-one __CPROVER function to allow @@ -2104,8 +1937,8 @@ long long int llround(double x) } fesetround(saved_rounding_mode); - - double rti = __sort_of_CPROVER_round_to_integral(FE_TOWARDZERO, xp); + + double rti = __CPROVER_round_to_integrald(xp, 0); // FE_TOWARDZERO return (long long int)rti; } @@ -2121,12 +1954,9 @@ long long int llround(double x) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); - long long int llroundf(float x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT, plus should use RNA + // TODO: should use RNA int saved_rounding_mode = fegetround(); fesetround(FE_TOWARDZERO); @@ -2140,8 +1970,8 @@ long long int llroundf(float x) } fesetround(saved_rounding_mode); - - float rti = __sort_of_CPROVER_round_to_integralf(FE_TOWARDZERO, xp); + + float rti = __CPROVER_round_to_integralf(xp, 0); // FE_TOWARDZERO return (long long int)rti; } @@ -2158,12 +1988,9 @@ long long int llroundf(float x) #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); - long long int llroundl(long double x) { - // TODO : should be an all-in-one __CPROVER function to allow - // conversion to SMT, plus should use RNA + // TODO: should use RNA int saved_rounding_mode = fegetround(); fesetround(FE_TOWARDZERO); @@ -2177,8 +2004,8 @@ long long int llroundl(long double x) } fesetround(saved_rounding_mode); - - long double rti = __sort_of_CPROVER_round_to_integrall(FE_TOWARDZERO, xp); + + long double rti = __CPROVER_round_to_integralld(xp, 0); // FE_TOWARDZERO return (long long int)rti; } @@ -2203,11 +2030,9 @@ long long int llroundl(long double x) #define __CPROVER_FENV_H_INCLUDED #endif -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); - double modf(double x, double *iptr) { - *iptr = __sort_of_CPROVER_round_to_integral(FE_TOWARDZERO, x); + *iptr = __CPROVER_round_to_integrald(x, 0); // FE_TOWARDZERO return (x - *iptr); } @@ -2223,11 +2048,9 @@ double modf(double x, double *iptr) #define __CPROVER_FENV_H_INCLUDED #endif -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); - - float modff(float x, float *iptr) +float modff(float x, float *iptr) { - *iptr = __sort_of_CPROVER_round_to_integralf(FE_TOWARDZERO, x); + *iptr = __CPROVER_round_to_integralf(x, 0); // FE_TOWARDZERO return (x - *iptr); } @@ -2244,11 +2067,9 @@ float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); #define __CPROVER_FENV_H_INCLUDED #endif -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); - - long double modfl(long double x, long double *iptr) +long double modfl(long double x, long double *iptr) { - *iptr = __sort_of_CPROVER_round_to_integralf(FE_TOWARDZERO, x); + *iptr = __CPROVER_round_to_integralld(x, 0); // FE_TOWARDZERO return (x - *iptr); } @@ -2256,8 +2077,7 @@ long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double /* FUNCTION: __sort_of_CPROVER_remainder */ // TODO : Should be a real __CPROVER function to convert to SMT-LIB -double __sort_of_CPROVER_round_to_integral (int rounding_mode, double d); - + double __sort_of_CPROVER_remainder (int rounding_mode, double x, double y) { if (x == 0.0 || __CPROVER_isinfd(y)) @@ -2265,7 +2085,7 @@ double __sort_of_CPROVER_remainder (int rounding_mode, double x, double y) // Extended precision helps... a bit... long double div = x/y; - long double n = __sort_of_CPROVER_round_to_integral(rounding_mode,div); + long double n = __CPROVER_round_to_integrald(rounding_mode, div); long double res = (-y * n) + x; // TODO : FMA would be an improvement return res; } @@ -2273,8 +2093,6 @@ double __sort_of_CPROVER_remainder (int rounding_mode, double x, double y) /* FUNCTION: __sort_of_CPROVER_remainderf */ // TODO : Should be a real __CPROVER function to convert to SMT-LIB -float __sort_of_CPROVER_round_to_integralf (int rounding_mode, float d); - float __sort_of_CPROVER_remainderf (int rounding_mode, float x, float y) { if (x == 0.0f || __CPROVER_isinff(y)) @@ -2282,7 +2100,7 @@ float __sort_of_CPROVER_remainderf (int rounding_mode, float x, float y) // Extended precision helps... a bit... long double div = x/y; - long double n = __sort_of_CPROVER_round_to_integral(rounding_mode,div); + long double n = __CPROVER_round_to_integralf(rounding_mode, div); long double res = (-y * n) + x; // TODO : FMA would be an improvement return res; } @@ -2290,8 +2108,6 @@ float __sort_of_CPROVER_remainderf (int rounding_mode, float x, float y) /* FUNCTION: __sort_of_CPROVER_remainderl */ // TODO : Should be a real __CPROVER function to convert to SMT-LIB -long double __sort_of_CPROVER_round_to_integrall (int rounding_mode, long double d); - long double __sort_of_CPROVER_remainderl (int rounding_mode, long double x, long double y) { if (x == 0.0 || __CPROVER_isinfld(y)) @@ -2299,7 +2115,7 @@ long double __sort_of_CPROVER_remainderl (int rounding_mode, long double x, long // Extended precision helps... a bit... long double div = x/y; - long double n = __sort_of_CPROVER_round_to_integral(rounding_mode,div); + long double n = __CPROVER_round_to_integralld(rounding_mode, div); long double res = (-y * n) + x; // TODO : FMA would be an improvement return res; } diff --git a/src/solvers/flattening/boolbv.cpp b/src/solvers/flattening/boolbv.cpp index 38f856898be..b94cd6db0c3 100644 --- a/src/solvers/flattening/boolbv.cpp +++ b/src/solvers/flattening/boolbv.cpp @@ -159,6 +159,9 @@ bvt boolbvt::convert_bitvector(const exprt &expr) return convert_floatbv_mod_rem(to_binary_expr(expr)); else if(expr.id()==ID_floatbv_typecast) return convert_floatbv_typecast(to_floatbv_typecast_expr(expr)); + else if(expr.id() == ID_floatbv_round_to_integral) + return convert_floatbv_round_to_integral( + to_floatbv_round_to_integral_expr(expr)); else if(expr.id()==ID_concatenation) return convert_concatenation(to_concatenation_expr(expr)); else if(expr.id()==ID_replication) diff --git a/src/solvers/flattening/boolbv.h b/src/solvers/flattening/boolbv.h index 78987b734ec..62d2703ee10 100644 --- a/src/solvers/flattening/boolbv.h +++ b/src/solvers/flattening/boolbv.h @@ -33,6 +33,7 @@ class byte_update_exprt; class concatenation_exprt; class extractbit_exprt; class extractbits_exprt; +class floatbv_round_to_integral_exprt; class floatbv_typecast_exprt; class ieee_float_op_exprt; class overflow_result_exprt; @@ -176,6 +177,8 @@ class boolbvt:public arrayst virtual bvt convert_floatbv_op(const ieee_float_op_exprt &); virtual bvt convert_floatbv_mod_rem(const binary_exprt &); virtual bvt convert_floatbv_typecast(const floatbv_typecast_exprt &expr); + virtual bvt + convert_floatbv_round_to_integral(const floatbv_round_to_integral_exprt &); virtual bvt convert_member(const member_exprt &expr); virtual bvt convert_with(const with_exprt &expr); virtual bvt convert_update(const update_exprt &); diff --git a/src/solvers/flattening/boolbv_floatbv_op.cpp b/src/solvers/flattening/boolbv_floatbv_op.cpp index 87dd3a8cd31..eb6d7eec40d 100644 --- a/src/solvers/flattening/boolbv_floatbv_op.cpp +++ b/src/solvers/flattening/boolbv_floatbv_op.cpp @@ -80,6 +80,24 @@ bvt boolbvt::convert_floatbv_typecast(const floatbv_typecast_exprt &expr) return conversion_failed(expr); } +bvt boolbvt::convert_floatbv_round_to_integral( + const floatbv_round_to_integral_exprt &expr) +{ + if(expr.op().type().id() == ID_floatbv) + { + float_utilst float_utils(prop); + + float_utils.set_rounding_mode(convert_bv(expr.rounding_mode())); + float_utils.spec = ieee_float_spect{to_floatbv_type(expr.op().type())}; + + auto op_bv = convert_bv(expr.op()); + + return float_utils.round_to_integral(op_bv); + } + else + return conversion_failed(expr); +} + bvt boolbvt::convert_floatbv_op(const ieee_float_op_exprt &expr) { const exprt &lhs = expr.lhs(); diff --git a/src/solvers/floatbv/float_utils.cpp b/src/solvers/floatbv/float_utils.cpp index 13706f06bd5..af918f3ab9a 100644 --- a/src/solvers/floatbv/float_utils.cpp +++ b/src/solvers/floatbv/float_utils.cpp @@ -44,7 +44,7 @@ bvt float_utilst::from_signed_integer(const bvt &src) src.size()-1, address_bits(src.size() - 1) + 1); - return rounder(result); + return round_and_pack(result); } bvt float_utilst::from_unsigned_integer(const bvt &src) @@ -61,7 +61,7 @@ bvt float_utilst::from_unsigned_integer(const bvt &src) result.sign=const_literal(false); - return rounder(result); + return round_and_pack(result); } bvt float_utilst::to_signed_integer( @@ -136,6 +136,32 @@ bvt float_utilst::to_integer( return result; } +bvt float_utilst::round_to_integral(const bvt &src) +{ + PRECONDITION(src.size() == spec.width()); + + // Zero? NaN? Infinity? + auto unpacked = unpack(src); + auto is_special = prop.lor({unpacked.zero, unpacked.NaN, unpacked.infinity}); + + // add 2^f, where f is the number of fraction bits, + // by adding f to the exponent + auto magic_number = ieee_floatt{ + spec, ieee_floatt::rounding_modet::ROUND_TO_ZERO, power(2, spec.f)}; + + auto magic_number_bv = build_constant(magic_number); + magic_number_bv.back() = src.back(); // copy sign bit + + auto tmp1 = add_sub(src, magic_number_bv, false); + + auto tmp2 = add_sub(tmp1, magic_number_bv, true); + + // restore the original sign bit + tmp2.back() = src.back(); + + return bv_utils.select(is_special, src, tmp2); +} + bvt float_utilst::build_constant(const ieee_floatt &src) { unbiased_floatt result; @@ -216,7 +242,7 @@ bvt float_utilst::conversion( // we actually need to round unbiased_floatt result=unpack(src); spec=dest_spec; - return rounder(result); + return round_and_pack(result); } } @@ -382,7 +408,7 @@ bvt float_utilst::add_sub( return pack(bias(result)); #endif - return rounder(result); + return round_and_pack(result); } /// Limits the shift distance @@ -458,7 +484,7 @@ bvt float_utilst::mul(const bvt &src1, const bvt &src2) result.NaN=prop.lor(NaN_cond); } - return rounder(result); + return round_and_pack(result); } bvt float_utilst::div(const bvt &src1, const bvt &src2) @@ -537,7 +563,7 @@ bvt float_utilst::div(const bvt &src1, const bvt &src2) result.fraction=bv_utils.select(force_zero, bv_utils.zeros(result.fraction.size()), result.fraction); - return rounder(result); + return round_and_pack(result); } bvt float_utilst::rem(const bvt &src1, const bvt &src2) @@ -901,11 +927,11 @@ void float_utilst::denormalization_shift(bvt &fraction, bvt &exponent) exponent); } -bvt float_utilst::rounder(const unbiased_floatt &src) +float_utilst::unbiased_floatt float_utilst::rounder(const unbiased_floatt &src) { // incoming: some fraction (with explicit 1), // some exponent without bias - // outgoing: rounded, with right size, with hidden bit, bias + // outgoing: rounded, with right size, but still unpacked bvt aligned_fraction=src.fraction, aligned_exponent=src.exponent; @@ -936,7 +962,12 @@ bvt float_utilst::rounder(const unbiased_floatt &src) round_fraction(result); round_exponent(result); - return pack(bias(result)); + return result; +} + +bvt float_utilst::round_and_pack(const unbiased_floatt &src) +{ + return pack(bias(rounder(src))); } /// rounding decision for fraction using sticky bit diff --git a/src/solvers/floatbv/float_utils.h b/src/solvers/floatbv/float_utils.h index 989c8d53377..1bd20f58c71 100644 --- a/src/solvers/floatbv/float_utils.h +++ b/src/solvers/floatbv/float_utils.h @@ -133,6 +133,7 @@ class float_utilst bvt to_signed_integer(const bvt &src, std::size_t int_width); bvt to_unsigned_integer(const bvt &src, std::size_t int_width); bvt conversion(const bvt &src, const ieee_float_spect &dest_spec); + bvt round_to_integral(const bvt &); // relations enum class relt { LT, LE, EQ, GT, GE }; @@ -189,12 +190,15 @@ class float_utilst { }; + unbiased_floatt unpack(const bvt &); biased_floatt bias(const unbiased_floatt &); + // takes unpacked, returns unpacked + unbiased_floatt rounder(const unbiased_floatt &); + // this takes unpacked format, and returns packed - virtual bvt rounder(const unbiased_floatt &); + bvt round_and_pack(const unbiased_floatt &); bvt pack(const biased_floatt &); - unbiased_floatt unpack(const bvt &); void round_fraction(unbiased_floatt &result); void round_exponent(unbiased_floatt &result); diff --git a/src/solvers/smt2/smt2_conv.cpp b/src/solvers/smt2/smt2_conv.cpp index dba907efd9d..dd7483a5b24 100644 --- a/src/solvers/smt2/smt2_conv.cpp +++ b/src/solvers/smt2/smt2_conv.cpp @@ -1199,6 +1199,10 @@ void smt2_convt::convert_expr(const exprt &expr) { convert_floatbv_typecast(to_floatbv_typecast_expr(expr)); } + else if(expr.id() == ID_floatbv_round_to_integral) + { + convert_floatbv_round_to_integral(to_floatbv_round_to_integral_expr(expr)); + } else if(expr.id()==ID_struct) { convert_struct(to_struct_expr(expr)); @@ -3272,6 +3276,23 @@ void smt2_convt::convert_floatbv_typecast(const floatbv_typecast_exprt &expr) } } +void smt2_convt::convert_floatbv_round_to_integral( + const floatbv_round_to_integral_exprt &expr) +{ + PRECONDITION(expr.type().id() == ID_floatbv); + + if(use_FPA_theory) + { + out << "(fp.roundToIntegral "; + convert_rounding_mode_FPA(expr.rounding_mode()); + out << ' '; + convert_expr(expr.op()); + out << ")"; + } + else + UNEXPECTEDCASE("TODO floatbv_round_to_integral without FPA"); +} + void smt2_convt::convert_struct(const struct_exprt &expr) { const struct_typet &struct_type = diff --git a/src/solvers/smt2/smt2_conv.h b/src/solvers/smt2/smt2_conv.h index fa6d79aec05..a2a8168c8aa 100644 --- a/src/solvers/smt2/smt2_conv.h +++ b/src/solvers/smt2/smt2_conv.h @@ -32,6 +32,7 @@ Author: Daniel Kroening, kroening@kroening.com class floatbv_typecast_exprt; class ieee_float_op_exprt; +class floatbv_round_to_integral_exprt; class union_typet; class update_bit_exprt; class update_bits_exprt; @@ -145,6 +146,8 @@ class smt2_convt : public stack_decision_proceduret void convert_floatbv_div(const ieee_float_op_exprt &expr); void convert_floatbv_mult(const ieee_float_op_exprt &expr); void convert_floatbv_rem(const binary_exprt &expr); + void + convert_floatbv_round_to_integral(const floatbv_round_to_integral_exprt &); void convert_mod(const mod_exprt &expr); void convert_euclidean_mod(const euclidean_mod_exprt &expr); void convert_index(const index_exprt &expr); diff --git a/src/solvers/smt2/smt2_parser.cpp b/src/solvers/smt2/smt2_parser.cpp index c952506111f..f08f06821b7 100644 --- a/src/solvers/smt2/smt2_parser.cpp +++ b/src/solvers/smt2/smt2_parser.cpp @@ -1406,6 +1406,19 @@ void smt2_parsert::setup_expressions() return binary_exprt(op[0], ID_floatbv_rem, op[1]); }; + expressions["fp.roundToIntegral"] = [this] { + auto op = operands(); + + if(op.size() != 2) + throw error() << "fp.roundToIntegral takes two operands"; + + if(op[1].type().id() != ID_floatbv) + throw error() << "fp.roundToIntegral takes a FloatingPoint operand"; + + // Note swapped order. + return floatbv_round_to_integral_exprt(op[1], op[0]); + }; + expressions["fp.eq"] = [this] { return function_application_ieee_float_eq(operands()); }; diff --git a/src/util/floatbv_expr.h b/src/util/floatbv_expr.h index c34f1c15af5..f0de4051c0e 100644 --- a/src/util/floatbv_expr.h +++ b/src/util/floatbv_expr.h @@ -83,6 +83,70 @@ inline floatbv_typecast_exprt &to_floatbv_typecast_expr(exprt &expr) return ret; } +/// \brief Round a floating-point number to an integral value +/// considering the given rounding mode +class floatbv_round_to_integral_exprt : public binary_exprt +{ +public: + floatbv_round_to_integral_exprt(exprt op, exprt rounding) + : binary_exprt( + op, + ID_floatbv_round_to_integral, + std::move(rounding), + op.type()) + { + } + + exprt &op() + { + return op0(); + } + + const exprt &op() const + { + return op0(); + } + + exprt &rounding_mode() + { + return op1(); + } + + const exprt &rounding_mode() const + { + return op1(); + } +}; + +template <> +inline bool can_cast_expr(const exprt &base) +{ + return base.id() == ID_floatbv_round_to_integral; +} + +/// \brief Cast an exprt to a \ref floatbv_round_to_integral_exprt +/// +/// \a expr must be known to be \ref floatbv_round_to_integral_exprt. +/// +/// \param expr: Source expression +/// \return Object of type \ref floatbv_round_to_integral_exprt +inline const floatbv_round_to_integral_exprt & +to_floatbv_round_to_integral_expr(const exprt &expr) +{ + PRECONDITION(expr.id() == ID_floatbv_round_to_integral); + floatbv_round_to_integral_exprt::check(expr); + return static_cast(expr); +} + +/// \copydoc to_floatbv_round_to_integral_expr(const exprt &) +inline floatbv_round_to_integral_exprt & +to_floatbv_round_to_integral_expr(exprt &expr) +{ + PRECONDITION(expr.id() == ID_floatbv_round_to_integral); + floatbv_round_to_integral_exprt::check(expr); + return static_cast(expr); +} + /// \brief Evaluates to true if the operand is NaN class isnan_exprt : public unary_predicate_exprt { diff --git a/src/util/ieee_float.cpp b/src/util/ieee_float.cpp index 050e51d03f9..34b384ea5ce 100644 --- a/src/util/ieee_float.cpp +++ b/src/util/ieee_float.cpp @@ -57,6 +57,13 @@ void ieee_float_spect::from_type(const floatbv_typet &type) e=e-1; // no hidden bit } +ieee_floatt ieee_floatt::abs() const +{ + ieee_floatt result = *this; + result.sign_flag = false; + return result; +} + constant_exprt ieee_floatt::rounding_mode_expr(rounding_modet rm) { return floatbv_rounding_mode(static_cast(rm)); @@ -595,7 +602,7 @@ void ieee_floatt::align() if(biased_exponent>=spec.max_exponent()) { // we need to consider the rounding mode here - switch(rounding_mode) + switch(_rounding_mode) { case UNKNOWN: case NONDETERMINISTIC: @@ -671,7 +678,7 @@ void ieee_floatt::divide_and_round( if(remainder!=0) { - switch(rounding_mode) + switch(_rounding_mode) { case ROUND_TO_EVEN: { @@ -861,7 +868,7 @@ ieee_floatt &ieee_floatt::operator+=(const ieee_floatt &other) return *this; else { - if(rounding_mode==ROUND_TO_MINUS_INF) + if(_rounding_mode == ROUND_TO_MINUS_INF) { set_sign(true); return *this; @@ -902,7 +909,7 @@ ieee_floatt &ieee_floatt::operator+=(const ieee_floatt &other) // there is some set of rules to get the sign if(fraction==0) { - if(rounding_mode==ROUND_TO_MINUS_INF) + if(_rounding_mode == ROUND_TO_MINUS_INF) sign_flag=true; else sign_flag=false; @@ -1052,6 +1059,18 @@ bool ieee_floatt::operator==(int i) const return *this==other; } +bool ieee_floatt::operator==(double d) const +{ + PRECONDITION(spec == ieee_float_spect::double_precision()); + return *this == ieee_floatt::from_double(_rounding_mode, d); +} + +bool ieee_floatt::operator==(float f) const +{ + PRECONDITION(spec == ieee_float_spect::single_precision()); + return *this == ieee_floatt::from_float(_rounding_mode, f); +} + bool ieee_floatt::operator!=(const ieee_floatt &other) const { return !(*this==other); @@ -1115,6 +1134,42 @@ mp_integer ieee_floatt::to_integer() const return result; } +ieee_floatt ieee_floatt::round_to_integral() const +{ + if(NaN_flag || infinity_flag || is_zero()) + return *this; + + ieee_floatt magic_number{spec, _rounding_mode}; + magic_number.from_integer(power(2, spec.f)); + + if(this->abs() >= magic_number) + return *this; + else + { + ieee_floatt result = *this; + + // add and subtract 2^f + // Preserve the original sign bit + bool original_sign = result.get_sign(); + + if(result.is_negative()) + { + result -= magic_number; + result += magic_number; + } + else + { + result += magic_number; + result -= magic_number; + } + + // Restore the original sign bit + result.sign_flag = original_sign; + + return result; + } +} + void ieee_floatt::from_double(const double d) { spec.f=52; diff --git a/src/util/ieee_float.h b/src/util/ieee_float.h index 8f1696291b9..b86d5b8b9f5 100644 --- a/src/util/ieee_float.h +++ b/src/util/ieee_float.h @@ -129,11 +129,15 @@ class ieee_floatt // A helper to turn a rounding mode into a constant bitvector expression static constant_exprt rounding_mode_expr(rounding_modet); + rounding_modet rounding_mode() const + { + return _rounding_mode; + } + ieee_float_spect spec; - explicit ieee_floatt(const ieee_float_spect &_spec) - : spec(_spec), - rounding_mode(ROUND_TO_EVEN), + ieee_floatt() + : _rounding_mode(ROUND_TO_EVEN), sign_flag(false), exponent(0), fraction(0), @@ -142,9 +146,9 @@ class ieee_floatt { } - ieee_floatt(ieee_float_spect __spec, rounding_modet __rounding_mode) - : spec(std::move(__spec)), - rounding_mode(__rounding_mode), + explicit ieee_floatt(const ieee_float_spect &_spec) + : spec(_spec), + _rounding_mode(ROUND_TO_EVEN), sign_flag(false), exponent(0), fraction(0), @@ -155,7 +159,7 @@ class ieee_floatt explicit ieee_floatt(const floatbv_typet &type) : spec(ieee_float_spect(type)), - rounding_mode(ROUND_TO_EVEN), + _rounding_mode(ROUND_TO_EVEN), sign_flag(false), exponent(0), fraction(0), @@ -164,11 +168,15 @@ class ieee_floatt { } - explicit ieee_floatt( - const floatbv_typet &type, - rounding_modet __rounding_mode) - : spec(ieee_float_spect(type)), - rounding_mode(__rounding_mode), + explicit ieee_floatt(const constant_exprt &expr) + : _rounding_mode(ROUND_TO_EVEN) + { + from_expr(expr); + } + + ieee_floatt(ieee_float_spect __spec, rounding_modet __rounding_mode) + : spec(std::move(__spec)), + _rounding_mode(__rounding_mode), sign_flag(false), exponent(0), fraction(0), @@ -177,21 +185,34 @@ class ieee_floatt { } - ieee_floatt(): - rounding_mode(ROUND_TO_EVEN), - sign_flag(false), exponent(0), fraction(0), - NaN_flag(false), infinity_flag(false) + ieee_floatt( + ieee_float_spect __spec, + rounding_modet __rounding_mode, + const mp_integer &i) + : spec(std::move(__spec)), + _rounding_mode(__rounding_mode), + sign_flag(false), + exponent(0), + fraction(0), + NaN_flag(false), + infinity_flag(false) { + from_integer(i); } - explicit ieee_floatt(const constant_exprt &expr): - rounding_mode(ROUND_TO_EVEN) + ieee_floatt(const floatbv_typet &type, rounding_modet __rounding_mode) + : spec(ieee_float_spect(type)), + _rounding_mode(__rounding_mode), + sign_flag(false), + exponent(0), + fraction(0), + NaN_flag(false), + infinity_flag(false) { - from_expr(expr); } ieee_floatt(const constant_exprt &expr, rounding_modet __rounding_mode) - : rounding_mode(__rounding_mode) + : _rounding_mode(__rounding_mode) { from_expr(expr); } @@ -276,6 +297,10 @@ class ieee_floatt return !NaN_flag && !infinity_flag && fraction==0 && exponent==0; } bool get_sign() const { return sign_flag; } + bool is_negative() const + { + return sign_flag; + } bool is_NaN() const { return NaN_flag; } bool is_infinity() const { return !NaN_flag && infinity_flag; } bool is_normal() const; @@ -291,6 +316,22 @@ class ieee_floatt void from_double(const double d); void from_float(const float f); + static ieee_floatt from_double(rounding_modet __rounding_mode, double value) + { + ieee_floatt result; + result._rounding_mode = __rounding_mode; + result.from_double(value); + return result; + } + + static ieee_floatt from_float(rounding_modet __rounding_mode, float value) + { + ieee_floatt result; + result._rounding_mode = __rounding_mode; + result.from_float(value); + return result; + } + // performs conversions from IEEE float-point format // to something else double to_double() const; @@ -331,23 +372,29 @@ class ieee_floatt bool operator<=(const ieee_floatt &other) const; bool operator>(const ieee_floatt &other) const; bool operator>=(const ieee_floatt &other) const; + ieee_floatt abs() const; // warning: these do packed equality, not IEEE equality // e.g., NAN==NAN bool operator==(const ieee_floatt &other) const; bool operator!=(const ieee_floatt &other) const; - bool operator==(int i) const; + bool operator==(int) const; + bool operator==(double) const; + bool operator==(float) const; // these do IEEE equality, i.e., NAN!=NAN bool ieee_equal(const ieee_floatt &other) const; bool ieee_not_equal(const ieee_floatt &other) const; + // Rounds according to the configured rounding mode + ieee_floatt round_to_integral() const; + protected: void divide_and_round(mp_integer ÷nd, const mp_integer &divisor); void align(); void next_representable(bool greater); - rounding_modet rounding_mode; + rounding_modet _rounding_mode; // we store the number unpacked bool sign_flag; diff --git a/src/util/irep_ids.def b/src/util/irep_ids.def index 507043b0d7d..0fbb2a95263 100644 --- a/src/util/irep_ids.def +++ b/src/util/irep_ids.def @@ -560,6 +560,7 @@ IREP_ID_ONE(floatbv_div) IREP_ID_ONE(floatbv_mod) IREP_ID_ONE(floatbv_rem) IREP_ID_ONE(floatbv_typecast) +IREP_ID_ONE(floatbv_round_to_integral) IREP_ID_ONE(compound_literal) IREP_ID_ONE(custom_bv) IREP_ID_ONE(custom_unsignedbv) diff --git a/unit/solvers/floatbv/float_utils.cpp b/unit/solvers/floatbv/float_utils.cpp index b7cdcce9b6a..d70ea45c8a9 100644 --- a/unit/solvers/floatbv/float_utils.cpp +++ b/unit/solvers/floatbv/float_utils.cpp @@ -218,3 +218,91 @@ SCENARIO("float_approximation", "[core][solvers][floatbv][float_approximation]") } } } + +static ieee_floatt round_to_integral(ieee_floatt op) +{ + satcheckt satcheck(null_message_handler); + float_utilst float_utils(satcheck); + float_utils.rounding_mode_bits.set(op.rounding_mode()); + float_utils.spec = op.spec; + + // add constraint + auto op_bv = float_utils.build_constant(op); + bvt result_bv = float_utils.round_to_integral(op_bv); + + // solve + const auto result = satcheck.prop_solve(); + REQUIRE(result == satcheckt::resultt::P_SATISFIABLE); + + return float_utils.get(result_bv); +} + +SCENARIO( + "float_utils_round_to_integral", + "[core][solvers][floatbv][float_utils][round_to_integral]") +{ + const auto sp = ieee_float_spect::single_precision(); + + // special values + REQUIRE(round_to_integral(ieee_floatt::NaN(sp)) == ieee_floatt::NaN(sp)); + REQUIRE( + round_to_integral(ieee_floatt::plus_infinity(sp)) == + ieee_floatt::plus_infinity(sp)); + REQUIRE( + round_to_integral(ieee_floatt::minus_infinity(sp)) == + ieee_floatt::minus_infinity(sp)); + + const auto up = ieee_floatt::ROUND_TO_PLUS_INF; + + // const auto dmax = std::numeric_limits::max(); + + REQUIRE(round_to_integral(ieee_floatt::from_double(up, 0)) == 0); + REQUIRE(round_to_integral(ieee_floatt::from_double(up, -0.0)) == -0.0); + REQUIRE(round_to_integral(ieee_floatt::from_double(up, 1)) == 1); + REQUIRE(round_to_integral(ieee_floatt::from_double(up, 0.1)) == 1); + REQUIRE(round_to_integral(ieee_floatt::from_double(up, -0.1)) == -0.0); + REQUIRE(round_to_integral(ieee_floatt::from_double(up, 10.1)) == 11); + REQUIRE(round_to_integral(ieee_floatt::from_double(up, -10.1)) == -10); + REQUIRE( + round_to_integral(ieee_floatt::from_double(up, 0x1.0p+52)) == 0x1.0p+52); + // REQUIRE(round_to_integral(ieee_floatt::from_double(up, dmax)) == dmax); + + const auto down = ieee_floatt::ROUND_TO_MINUS_INF; + + REQUIRE(round_to_integral(ieee_floatt::from_double(down, 0)) == 0); + REQUIRE(round_to_integral(ieee_floatt::from_double(down, -0.0)) == -0.0); + REQUIRE(round_to_integral(ieee_floatt::from_double(down, 1)) == 1); + REQUIRE(round_to_integral(ieee_floatt::from_double(down, 0.1)) == 0); + REQUIRE(round_to_integral(ieee_floatt::from_double(down, -0.1)) == -1); + REQUIRE(round_to_integral(ieee_floatt::from_double(down, 10.1)) == 10); + REQUIRE(round_to_integral(ieee_floatt::from_double(down, -10.1)) == -11); + REQUIRE( + round_to_integral(ieee_floatt::from_double(down, 0x1.0p+52)) == 0x1.0p+52); + // REQUIRE(round_to_integral(ieee_floatt::from_double(down, dmax)) == dmax); + + const auto even = ieee_floatt::ROUND_TO_EVEN; + + REQUIRE(round_to_integral(ieee_floatt::from_double(even, 0)) == 0); + REQUIRE(round_to_integral(ieee_floatt::from_double(even, -0.0)) == -0.0); + REQUIRE(round_to_integral(ieee_floatt::from_double(even, 1)) == 1); + REQUIRE(round_to_integral(ieee_floatt::from_double(even, 0.1)) == 0); + REQUIRE(round_to_integral(ieee_floatt::from_double(even, -0.1)) == -0.0); + REQUIRE(round_to_integral(ieee_floatt::from_double(even, 10.1)) == 10); + REQUIRE(round_to_integral(ieee_floatt::from_double(even, -10.1)) == -10); + REQUIRE( + round_to_integral(ieee_floatt::from_double(even, 0x1.0p+52)) == 0x1.0p+52); + // REQUIRE(round_to_integral(ieee_floatt::from_double(even, dmax)) == dmax); + + const auto zero = ieee_floatt::ROUND_TO_ZERO; + + REQUIRE(round_to_integral(ieee_floatt::from_double(zero, 0)) == 0); + REQUIRE(round_to_integral(ieee_floatt::from_double(zero, -0.0)) == -0.0); + REQUIRE(round_to_integral(ieee_floatt::from_double(zero, 1)) == 1); + REQUIRE(round_to_integral(ieee_floatt::from_double(zero, 0.1)) == 0); + REQUIRE(round_to_integral(ieee_floatt::from_double(zero, -0.1)) == -0.0); + REQUIRE(round_to_integral(ieee_floatt::from_double(zero, 10.1)) == 10); + REQUIRE(round_to_integral(ieee_floatt::from_double(zero, -10.1)) == -10); + REQUIRE( + round_to_integral(ieee_floatt::from_double(zero, 0x1.0p+52)) == 0x1.0p+52); + // REQUIRE(round_to_integral(ieee_floatt::from_double(zero, dmax)) == dmax); +} diff --git a/unit/util/ieee_float.cpp b/unit/util/ieee_float.cpp index 2b4b6e7c2ca..dc008f1b166 100644 --- a/unit/util/ieee_float.cpp +++ b/unit/util/ieee_float.cpp @@ -1,6 +1,6 @@ /*******************************************************************\ -Module: Unit tests for ieee_floatt +Module: Unit test for util/ieee_float Author: Daniel Kroening, dkr@amazon.com @@ -10,8 +10,76 @@ Author: Daniel Kroening, dkr@amazon.com #include +#include + TEST_CASE("Make an IEEE 754 one", "[core][util][ieee_float]") { auto spec = ieee_float_spect::single_precision(); REQUIRE(ieee_floatt::one(spec) == 1); } + +TEST_CASE("round_to_integral", "[unit][util][ieee_float]") +{ + const auto sp = ieee_float_spect::single_precision(); + + // special values + REQUIRE(ieee_floatt::NaN(sp).round_to_integral() == ieee_floatt::NaN(sp)); + REQUIRE( + ieee_floatt::plus_infinity(sp).round_to_integral() == + ieee_floatt::plus_infinity(sp)); + REQUIRE( + ieee_floatt::minus_infinity(sp).round_to_integral() == + ieee_floatt::minus_infinity(sp)); + + const auto dmax = std::numeric_limits::max(); + + const auto up = ieee_floatt::ROUND_TO_PLUS_INF; + + REQUIRE(ieee_floatt::from_double(up, 0).round_to_integral() == 0); + REQUIRE(ieee_floatt::from_double(up, -0.0).round_to_integral() == -0.0); + REQUIRE(ieee_floatt::from_double(up, 1).round_to_integral() == 1); + REQUIRE(ieee_floatt::from_double(up, 0.1).round_to_integral() == 1); + REQUIRE(ieee_floatt::from_double(up, -0.1).round_to_integral() == -0.0); + REQUIRE(ieee_floatt::from_double(up, 10.1).round_to_integral() == 11); + REQUIRE(ieee_floatt::from_double(up, -10.1).round_to_integral() == -10); + REQUIRE(ieee_floatt::from_double(up, dmax).round_to_integral() == dmax); + + const auto down = ieee_floatt::ROUND_TO_MINUS_INF; + + REQUIRE(ieee_floatt::from_double(down, 0).round_to_integral() == 0); + REQUIRE(ieee_floatt::from_double(down, -0.0).round_to_integral() == -0.0); + REQUIRE(ieee_floatt::from_double(down, 1).round_to_integral() == 1); + REQUIRE(ieee_floatt::from_double(down, 0.1).round_to_integral() == 0); + REQUIRE(ieee_floatt::from_double(down, -0.1).round_to_integral() == -1); + REQUIRE(ieee_floatt::from_double(down, 10.1).round_to_integral() == 10); + REQUIRE(ieee_floatt::from_double(down, -10.1).round_to_integral() == -11); + REQUIRE( + ieee_floatt::from_double(down, 0x1.0p+52).round_to_integral() == 0x1.0p+52); + REQUIRE(ieee_floatt::from_double(down, dmax).round_to_integral() == dmax); + + const auto even = ieee_floatt::ROUND_TO_EVEN; + + REQUIRE(ieee_floatt::from_double(even, 0).round_to_integral() == 0); + REQUIRE(ieee_floatt::from_double(even, -0.0).round_to_integral() == -0.0); + REQUIRE(ieee_floatt::from_double(even, 1).round_to_integral() == 1); + REQUIRE(ieee_floatt::from_double(even, 0.1).round_to_integral() == 0); + REQUIRE(ieee_floatt::from_double(even, -0.1).round_to_integral() == -0.0); + REQUIRE(ieee_floatt::from_double(even, 10.1).round_to_integral() == 10); + REQUIRE(ieee_floatt::from_double(even, -10.1).round_to_integral() == -10); + REQUIRE( + ieee_floatt::from_double(even, 0x1.0p+52).round_to_integral() == 0x1.0p+52); + REQUIRE(ieee_floatt::from_double(even, dmax).round_to_integral() == dmax); + + const auto zero = ieee_floatt::ROUND_TO_ZERO; + + REQUIRE(ieee_floatt::from_double(zero, 0).round_to_integral() == 0); + REQUIRE(ieee_floatt::from_double(zero, -0.0).round_to_integral() == -0.0); + REQUIRE(ieee_floatt::from_double(zero, 1).round_to_integral() == 1); + REQUIRE(ieee_floatt::from_double(zero, 0.1).round_to_integral() == 0); + REQUIRE(ieee_floatt::from_double(zero, -0.1).round_to_integral() == -0.0); + REQUIRE(ieee_floatt::from_double(zero, 10.1).round_to_integral() == 10); + REQUIRE(ieee_floatt::from_double(zero, -10.1).round_to_integral() == -10); + REQUIRE( + ieee_floatt::from_double(zero, 0x1.0p+52).round_to_integral() == 0x1.0p+52); + REQUIRE(ieee_floatt::from_double(zero, dmax).round_to_integral() == dmax); +}