1 /* IEC 559 floating point format related functions. 2 Copyright (C) 2001-2010 Roberto Bagnara <bagnara@cs.unipr.it> 3 Copyright (C) 2010-2016 BUGSENG srl (http://bugseng.com) 4 5 This file is part of the Parma Polyhedra Library (PPL). 6 7 The PPL is free software; you can redistribute it and/or modify it 8 under the terms of the GNU General Public License as published by the 9 Free Software Foundation; either version 3 of the License, or (at your 10 option) any later version. 11 12 The PPL is distributed in the hope that it will be useful, but WITHOUT 13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 15 for more details. 16 17 You should have received a copy of the GNU General Public License 18 along with this program; if not, write to the Free Software Foundation, 19 Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02111-1307, USA. 20 21 For the most up-to-date information see the Parma Polyhedra Library 22 site: http://bugseng.com/products/ppl/ . */ 23 24 #ifndef PPL_Float_defs_hh 25 #define PPL_Float_defs_hh 1 26 27 #include "globals_types.hh" 28 #include "meta_programming.hh" 29 #include "compiler.hh" 30 #include "Concrete_Expression_types.hh" 31 #include "Variable_types.hh" 32 #include "Linear_Form_types.hh" 33 #include <set> 34 #include <cmath> 35 #include <map> 36 #include <gmp.h> 37 38 #ifdef NAN 39 #define PPL_NAN NAN 40 #else 41 #define PPL_NAN (HUGE_VAL - HUGE_VAL) 42 #endif 43 44 namespace Parma_Polyhedra_Library { 45 46 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS 47 /*! \ingroup PPL_CXX_interface */ 48 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS) 49 struct float_ieee754_half { 50 uint16_t word; 51 static const uint16_t SGN_MASK = 0x8000U; 52 static const uint16_t EXP_MASK = 0xfc00U; 53 static const uint16_t WRD_MAX = 0x7bffU; 54 static const uint16_t POS_INF = 0x7c00U; 55 static const uint16_t NEG_INF = 0xfc00U; 56 static const uint16_t POS_ZERO = 0x0000U; 57 static const uint16_t NEG_ZERO = 0x8000U; 58 static const unsigned int BASE = 2; 59 static const unsigned int EXPONENT_BITS = 5; 60 static const unsigned int MANTISSA_BITS = 10; 61 static const int EXPONENT_MAX = (1 << (EXPONENT_BITS - 1)) - 1; 62 static const int EXPONENT_BIAS = EXPONENT_MAX; 63 static const int EXPONENT_MIN = -EXPONENT_MAX + 1; 64 static const int EXPONENT_MIN_DENORM = EXPONENT_MIN 65 - static_cast<int>(MANTISSA_BITS); 66 static const Floating_Point_Format floating_point_format = IEEE754_HALF; 67 int inf_sign() const; 68 bool is_nan() const; 69 int zero_sign() const; 70 bool sign_bit() const; 71 void negate(); 72 void dec(); 73 void inc(); 74 void set_max(bool negative); 75 void build(bool negative, mpz_t mantissa, int exponent); 76 77 }; 78 79 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS 80 /*! \ingroup PPL_CXX_interface */ 81 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS) 82 struct float_ieee754_single { 83 uint32_t word; 84 static const uint32_t SGN_MASK = 0x80000000U; 85 static const uint32_t EXP_MASK = 0x7f800000U; 86 static const uint32_t WRD_MAX = 0x7f7fffffU; 87 static const uint32_t POS_INF = 0x7f800000U; 88 static const uint32_t NEG_INF = 0xff800000U; 89 static const uint32_t POS_ZERO = 0x00000000U; 90 static const uint32_t NEG_ZERO = 0x80000000U; 91 static const unsigned int BASE = 2; 92 static const unsigned int EXPONENT_BITS = 8; 93 static const unsigned int MANTISSA_BITS = 23; 94 static const int EXPONENT_MAX = (1 << (EXPONENT_BITS - 1)) - 1; 95 static const int EXPONENT_BIAS = EXPONENT_MAX; 96 static const int EXPONENT_MIN = -EXPONENT_MAX + 1; 97 static const int EXPONENT_MIN_DENORM = EXPONENT_MIN 98 - static_cast<int>(MANTISSA_BITS); 99 static const Floating_Point_Format floating_point_format = IEEE754_SINGLE; 100 int inf_sign() const; 101 bool is_nan() const; 102 int zero_sign() const; 103 bool sign_bit() const; 104 void negate(); 105 void dec(); 106 void inc(); 107 void set_max(bool negative); 108 void build(bool negative, mpz_t mantissa, int exponent); 109 }; 110 111 #ifdef WORDS_BIGENDIAN 112 #ifndef PPL_WORDS_BIGENDIAN 113 #define PPL_WORDS_BIGENDIAN 114 #endif 115 #endif 116 117 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS 118 /*! \ingroup PPL_CXX_interface */ 119 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS) 120 struct float_ieee754_double { 121 #ifdef PPL_WORDS_BIGENDIAN 122 uint32_t msp; 123 uint32_t lsp; 124 #else 125 uint32_t lsp; 126 uint32_t msp; 127 #endif 128 static const uint32_t MSP_SGN_MASK = 0x80000000U; 129 static const uint32_t MSP_POS_INF = 0x7ff00000U; 130 static const uint32_t MSP_NEG_INF = 0xfff00000U; 131 static const uint32_t MSP_POS_ZERO = 0x00000000U; 132 static const uint32_t MSP_NEG_ZERO = 0x80000000U; 133 static const uint32_t LSP_INF = 0; 134 static const uint32_t LSP_ZERO = 0; 135 static const uint32_t MSP_MAX = 0x7fefffffU; 136 static const uint32_t LSP_MAX = 0xffffffffU; 137 static const unsigned int BASE = 2; 138 static const unsigned int EXPONENT_BITS = 11; 139 static const unsigned int MANTISSA_BITS = 52; 140 static const int EXPONENT_MAX = (1 << (EXPONENT_BITS - 1)) - 1; 141 static const int EXPONENT_BIAS = EXPONENT_MAX; 142 static const int EXPONENT_MIN = -EXPONENT_MAX + 1; 143 static const int EXPONENT_MIN_DENORM = EXPONENT_MIN 144 - static_cast<int>(MANTISSA_BITS); 145 static const Floating_Point_Format floating_point_format = IEEE754_DOUBLE; 146 int inf_sign() const; 147 bool is_nan() const; 148 int zero_sign() const; 149 bool sign_bit() const; 150 void negate(); 151 void dec(); 152 void inc(); 153 void set_max(bool negative); 154 void build(bool negative, mpz_t mantissa, int exponent); 155 }; 156 157 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS 158 /*! \ingroup PPL_CXX_interface */ 159 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS) 160 struct float_ibm_single { 161 uint32_t word; 162 static const uint32_t SGN_MASK = 0x80000000U; 163 static const uint32_t EXP_MASK = 0x7f000000U; 164 static const uint32_t WRD_MAX = 0x7fffffffU; 165 static const uint32_t POS_INF = 0x7f000000U; 166 static const uint32_t NEG_INF = 0xff000000U; 167 static const uint32_t POS_ZERO = 0x00000000U; 168 static const uint32_t NEG_ZERO = 0x80000000U; 169 static const unsigned int BASE = 16; 170 static const unsigned int EXPONENT_BITS = 7; 171 static const unsigned int MANTISSA_BITS = 24; 172 static const int EXPONENT_BIAS = 64; 173 static const int EXPONENT_MAX = (1 << (EXPONENT_BITS - 1)) - 1; 174 static const int EXPONENT_MIN = -EXPONENT_MAX + 1; 175 static const int EXPONENT_MIN_DENORM = EXPONENT_MIN 176 - static_cast<int>(MANTISSA_BITS); 177 static const Floating_Point_Format floating_point_format = IBM_SINGLE; 178 int inf_sign() const; 179 bool is_nan() const; 180 int zero_sign() const; 181 bool sign_bit() const; 182 void negate(); 183 void dec(); 184 void inc(); 185 void set_max(bool negative); 186 void build(bool negative, mpz_t mantissa, int exponent); 187 }; 188 189 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS 190 /*! \ingroup PPL_CXX_interface */ 191 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS) 192 struct float_ibm_double { 193 static const unsigned int BASE = 16; 194 static const unsigned int EXPONENT_BITS = 7; 195 static const unsigned int MANTISSA_BITS = 56; 196 static const int EXPONENT_BIAS = 64; 197 }; 198 199 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS 200 /*! \ingroup PPL_CXX_interface */ 201 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS) 202 struct float_intel_double_extended { 203 #ifdef PPL_WORDS_BIGENDIAN 204 uint32_t msp; 205 uint64_t lsp; 206 #else 207 uint64_t lsp; 208 uint32_t msp; 209 #endif 210 static const uint32_t MSP_SGN_MASK = 0x00008000U; 211 static const uint32_t MSP_POS_INF = 0x00007fffU; 212 static const uint32_t MSP_NEG_INF = 0x0000ffffU; 213 static const uint32_t MSP_POS_ZERO = 0x00000000U; 214 static const uint32_t MSP_NEG_ZERO = 0x00008000U; 215 static const uint64_t LSP_INF = static_cast<uint64_t>(0x8000000000000000ULL); 216 static const uint64_t LSP_ZERO = 0; 217 static const uint32_t MSP_MAX = 0x00007ffeU; 218 static const uint64_t LSP_DMAX = static_cast<uint64_t>(0x7fffffffffffffffULL); 219 static const uint64_t LSP_NMAX = static_cast<uint64_t>(0xffffffffffffffffULL); 220 static const unsigned int BASE = 2; 221 static const unsigned int EXPONENT_BITS = 15; 222 static const unsigned int MANTISSA_BITS = 63; 223 static const int EXPONENT_MAX = (1 << (EXPONENT_BITS - 1)) - 1; 224 static const int EXPONENT_BIAS = EXPONENT_MAX; 225 static const int EXPONENT_MIN = -EXPONENT_MAX + 1; 226 static const int EXPONENT_MIN_DENORM = EXPONENT_MIN 227 - static_cast<int>(MANTISSA_BITS); 228 static const Floating_Point_Format floating_point_format = 229 INTEL_DOUBLE_EXTENDED; 230 int inf_sign() const; 231 bool is_nan() const; 232 int zero_sign() const; 233 bool sign_bit() const; 234 void negate(); 235 void dec(); 236 void inc(); 237 void set_max(bool negative); 238 void build(bool negative, mpz_t mantissa, int exponent); 239 }; 240 241 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS 242 /*! \ingroup PPL_CXX_interface */ 243 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS) 244 struct float_ieee754_quad { 245 #ifdef PPL_WORDS_BIGENDIAN 246 uint64_t msp; 247 uint64_t lsp; 248 #else 249 uint64_t lsp; 250 uint64_t msp; 251 #endif 252 static const uint64_t MSP_SGN_MASK = static_cast<uint64_t>(0x8000000000000000ULL); 253 static const uint64_t MSP_POS_INF = static_cast<uint64_t>(0x7fff000000000000ULL); 254 static const uint64_t MSP_NEG_INF = static_cast<uint64_t>(0xffff000000000000ULL); 255 static const uint64_t MSP_POS_ZERO = static_cast<uint64_t>(0x0000000000000000ULL); 256 static const uint64_t MSP_NEG_ZERO = static_cast<uint64_t>(0x8000000000000000ULL); 257 static const uint64_t LSP_INF = 0; 258 static const uint64_t LSP_ZERO = 0; 259 static const uint64_t MSP_MAX = static_cast<uint64_t>(0x7ffeffffffffffffULL); 260 static const uint64_t LSP_MAX = static_cast<uint64_t>(0xffffffffffffffffULL); 261 static const unsigned int BASE = 2; 262 static const unsigned int EXPONENT_BITS = 15; 263 static const unsigned int MANTISSA_BITS = 112; 264 static const int EXPONENT_MAX = (1 << (EXPONENT_BITS - 1)) - 1; 265 static const int EXPONENT_BIAS = EXPONENT_MAX; 266 static const int EXPONENT_MIN = -EXPONENT_MAX + 1; 267 static const int EXPONENT_MIN_DENORM = EXPONENT_MIN 268 - static_cast<int>(MANTISSA_BITS); 269 static const Floating_Point_Format floating_point_format = IEEE754_QUAD; 270 int inf_sign() const; 271 bool is_nan() const; 272 int zero_sign() const; 273 bool sign_bit() const; 274 void negate(); 275 void dec(); 276 void inc(); 277 void set_max(bool negative); 278 void build(bool negative, mpz_t mantissa, int exponent); 279 }; 280 281 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS 282 /*! \ingroup PPL_CXX_interface */ 283 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS) 284 template <typename T> 285 class Float : public False { }; 286 287 #if PPL_SUPPORTED_FLOAT 288 template <> 289 class Float<float> : public True { 290 public: 291 #if PPL_CXX_FLOAT_BINARY_FORMAT == PPL_FLOAT_IEEE754_HALF 292 typedef float_ieee754_half Binary; 293 #elif PPL_CXX_FLOAT_BINARY_FORMAT == PPL_FLOAT_IEEE754_SINGLE 294 typedef float_ieee754_single Binary; 295 #elif PPL_CXX_FLOAT_BINARY_FORMAT == PPL_FLOAT_IEEE754_DOUBLE 296 typedef float_ieee754_double Binary; 297 #elif PPL_CXX_FLOAT_BINARY_FORMAT == PPL_FLOAT_IBM_SINGLE 298 typedef float_ibm_single Binary; 299 #elif PPL_CXX_FLOAT_BINARY_FORMAT == PPL_FLOAT_IEEE754_QUAD 300 typedef float_ieee754_quad Binary; 301 #elif PPL_CXX_FLOAT_BINARY_FORMAT == PPL_FLOAT_INTEL_DOUBLE_EXTENDED 302 typedef float_intel_double_extended Binary; 303 #else 304 #error "Invalid value for PPL_CXX_FLOAT_BINARY_FORMAT" 305 #endif 306 union { 307 float number; 308 Binary binary; 309 } u; 310 Float(); 311 Float(float v); 312 float value(); 313 }; 314 #endif 315 316 #if PPL_SUPPORTED_DOUBLE 317 template <> 318 class Float<double> : public True { 319 public: 320 #if PPL_CXX_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IEEE754_HALF 321 typedef float_ieee754_half Binary; 322 #elif PPL_CXX_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IEEE754_SINGLE 323 typedef float_ieee754_single Binary; 324 #elif PPL_CXX_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IEEE754_DOUBLE 325 typedef float_ieee754_double Binary; 326 #elif PPL_CXX_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IBM_SINGLE 327 typedef float_ibm_single Binary; 328 #elif PPL_CXX_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IEEE754_QUAD 329 typedef float_ieee754_quad Binary; 330 #elif PPL_CXX_DOUBLE_BINARY_FORMAT == PPL_FLOAT_INTEL_DOUBLE_EXTENDED 331 typedef float_intel_double_extended Binary; 332 #else 333 #error "Invalid value for PPL_CXX_DOUBLE_BINARY_FORMAT" 334 #endif 335 union { 336 double number; 337 Binary binary; 338 } u; 339 Float(); 340 Float(double v); 341 double value(); 342 }; 343 #endif 344 345 #if PPL_SUPPORTED_LONG_DOUBLE 346 template <> 347 class Float<long double> : public True { 348 public: 349 #if PPL_CXX_LONG_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IEEE754_HALF 350 typedef float_ieee754_half Binary; 351 #elif PPL_CXX_LONG_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IEEE754_SINGLE 352 typedef float_ieee754_single Binary; 353 #elif PPL_CXX_LONG_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IEEE754_DOUBLE 354 typedef float_ieee754_double Binary; 355 #elif PPL_CXX_LONG_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IBM_SINGLE 356 typedef float_ibm_single Binary; 357 #elif PPL_CXX_LONG_DOUBLE_BINARY_FORMAT == PPL_FLOAT_IEEE754_QUAD 358 typedef float_ieee754_quad Binary; 359 #elif PPL_CXX_LONG_DOUBLE_BINARY_FORMAT == PPL_FLOAT_INTEL_DOUBLE_EXTENDED 360 typedef float_intel_double_extended Binary; 361 #else 362 #error "Invalid value for PPL_CXX_LONG_DOUBLE_BINARY_FORMAT" 363 #endif 364 union { 365 long double number; 366 Binary binary; 367 } u; 368 Float(); 369 Float(long double v); 370 long double value(); 371 }; 372 #endif 373 374 // FIXME: is this the right place for this function? 375 #ifdef PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS 376 /*! \brief 377 If \p v is nonzero, returns the position of the most significant bit 378 in \p a. 379 380 The behavior is undefined if \p v is zero. 381 */ 382 #endif // defined(PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS) 383 unsigned int msb_position(unsigned long long v); 384 385 /*! \brief 386 An abstract class to be implemented by an external analyzer such 387 as ECLAIR in order to provide to the PPL the necessary information 388 for performing the analysis of floating point computations. 389 390 \par Template type parameters 391 392 - The class template parameter \p Target specifies the implementation 393 of Concrete_Expression to be used. 394 - The class template parameter \p FP_Interval_Type represents the type 395 of the intervals used in the abstract domain. The interval bounds 396 should have a floating point type. 397 */ 398 template <typename Target, typename FP_Interval_Type> 399 class FP_Oracle { 400 public: 401 /* 402 FIXME: the const qualifiers on expressions may raise CLANG 403 compatibility issues. It may be necessary to omit them. 404 */ 405 406 /*! \brief 407 Asks the external analyzer for an interval that correctly 408 approximates the floating point entity referenced by \p dim. 409 Result is stored into \p result. 410 411 \return <CODE>true</CODE> if the analyzer was able to find a correct 412 approximation, <CODE>false</CODE> otherwise. 413 */ 414 virtual bool get_interval(dimension_type dim, FP_Interval_Type& result) const 415 = 0; 416 417 /*! \brief 418 Asks the external analyzer for an interval that correctly 419 approximates the value of floating point constant \p expr. 420 Result is stored into \p result. 421 422 \return <CODE>true</CODE> if the analyzer was able to find a correct 423 approximation, <CODE>false</CODE> otherwise. 424 */ 425 virtual bool get_fp_constant_value( 426 const Floating_Point_Constant<Target>& expr, 427 FP_Interval_Type& result) const = 0; 428 429 /*! \brief 430 Asks the external analyzer for an interval that correctly approximates 431 the value of \p expr, which must be of integer type. 432 Result is stored into \p result. 433 434 \return <CODE>true</CODE> if the analyzer was able to find a correct 435 approximation, <CODE>false</CODE> otherwise. 436 */ 437 virtual bool get_integer_expr_value(const Concrete_Expression<Target>& expr, 438 FP_Interval_Type& result) const = 0; 439 440 /*! \brief 441 Asks the external analyzer for the possible space dimensions that 442 are associated to the approximable reference \p expr. 443 Result is stored into \p result. 444 445 \return <CODE>true</CODE> if the analyzer was able to return 446 the (possibly empty!) set, <CODE>false</CODE> otherwise. 447 448 The resulting set MUST NOT contain <CODE>not_a_dimension()</CODE>. 449 */ 450 virtual bool get_associated_dimensions( 451 const Approximable_Reference<Target>& expr, 452 std::set<dimension_type>& result) const = 0; 453 454 }; 455 456 /* FIXME: some of the following documentation should probably be 457 under PPL_DOXYGEN_INCLUDE_IMPLEMENTATION_DETAILS */ 458 459 /*! \brief \relates Float 460 Returns <CODE>true</CODE> if and only if there is some floating point 461 number that is representable by \p f2 but not by \p f1. 462 */ 463 bool is_less_precise_than(Floating_Point_Format f1, Floating_Point_Format f2); 464 465 /*! \brief \relates Float 466 Computes the absolute error of floating point computations. 467 468 \par Template type parameters 469 470 - The class template parameter \p FP_Interval_Type represents the type 471 of the intervals used in the abstract domain. The interval bounds 472 should have a floating point type. 473 474 \param analyzed_format The floating point format used by the analyzed 475 program. 476 477 \return The interval \f$[-\omega, \omega]\f$ where \f$\omega\f$ is the 478 smallest non-zero positive number in the less precise floating point 479 format between the analyzer format and the analyzed format. 480 */ 481 template <typename FP_Interval_Type> 482 const FP_Interval_Type& 483 compute_absolute_error(Floating_Point_Format analyzed_format); 484 485 /*! \brief \relates Linear_Form 486 Discards all linear forms containing variable \p var from the 487 linear form abstract store \p lf_store. 488 */ 489 template <typename FP_Interval_Type> 490 void 491 discard_occurrences(std::map<dimension_type, 492 Linear_Form<FP_Interval_Type> >& lf_store, 493 Variable var); 494 495 /*! \brief \relates Linear_Form 496 Assigns the linear form \p lf to \p var in the linear form abstract 497 store \p lf_store, then discards all occurrences of \p var from it. 498 */ 499 template <typename FP_Interval_Type> 500 void 501 affine_form_image(std::map<dimension_type, 502 Linear_Form<FP_Interval_Type> >& lf_store, 503 Variable var, 504 const Linear_Form<FP_Interval_Type>& lf); 505 506 /*! \brief \relates Linear_Form 507 Discards from \p ls1 all linear forms but those that are associated 508 to the same variable in \p ls2. 509 */ 510 template <typename FP_Interval_Type> 511 void 512 upper_bound_assign(std::map<dimension_type, 513 Linear_Form<FP_Interval_Type> >& ls1, 514 const std::map<dimension_type, 515 Linear_Form<FP_Interval_Type> >& ls2); 516 517 } // namespace Parma_Polyhedra_Library 518 519 #include "Float_inlines.hh" 520 #include "Float_templates.hh" 521 522 #endif // !defined(PPL_Float_defs_hh) 523