Many resources are needed to download a project. Please understand that we have to compensate our server costs. Thank you in advance. Project price only 1 $
You can buy this project and download/modify it how often you want.
/*++
Copyright (c) 2012 Microsoft Corporation
Module Name:
mpff.h
Abstract:
Multi precision fast floating point numbers.
The implementation is not compliant with the IEEE standard.
For an IEEE compliant implementation, see mpf.h
There are only two rounding modes: towards plus or minus inf.
Author:
Leonardo de Moura (leonardo) 2012-09-12
Revision History:
--*/
#pragma once
#include "util/id_gen.h"
#include "util/util.h"
#include "util/vector.h"
#include "util/z3_exception.h"
#include "util/scoped_numeral.h"
#include "util/scoped_numeral_vector.h"
#include "util/mpn.h"
class mpff_manager;
class mpff {
friend class mpff_manager;
unsigned m_sign:1;
unsigned m_sig_idx:31; // position where the significand is stored in the mpff_manager.
int m_exponent;
public:
mpff():
m_sign(0),
m_sig_idx(0),
m_exponent(0) {
}
void swap(mpff & other) noexcept {
unsigned sign = m_sign; m_sign = other.m_sign; other.m_sign = sign;
unsigned sig_idx = m_sig_idx; m_sig_idx = other.m_sig_idx; other.m_sig_idx = sig_idx;
std::swap(m_exponent, other.m_exponent);
}
};
inline void swap(mpff & m1, mpff & m2) noexcept { m1.swap(m2); }
class mpz;
class mpq;
template class mpz_manager;
template class mpq_manager;
#ifndef SINGLE_THREAD
typedef mpz_manager synch_mpz_manager;
typedef mpq_manager synch_mpq_manager;
#else
typedef mpz_manager synch_mpz_manager;
typedef mpq_manager synch_mpq_manager;
#endif
typedef mpz_manager unsynch_mpz_manager;
typedef mpq_manager unsynch_mpq_manager;
class mpff_manager {
// Some restrictions on mpff numbers
//
// - The exponent is always a machine integer. The main point is that 2^(2^31) is a huge number,
// we will not even be able to convert the mpff into mpq. Formulas that need this kind of huge number
// are usually out-of-reach for Z3.
//
// - The significand size is measured in words of 32-bit. The number of words is always even.
// This decision makes sure that the size (in bits) of mpff numbers is always a multiple of 64.
// Thus mpff objs can be easily packed in 64-bit machines.
//
// - The smallest mpff numeral has 128-bits total. mpff structure has always 64-bits.
// The minimal size for the significand is 64-bits.
//
// - All mpff numerals in a given manager use the same number of words for storing the significand.
// This is different from the mpf_manager where the same manager can be used to manipulate floating point numbers
// of different precision.
//
// - In the encoding used for mpff numbers, the most significand bit of the most significand word is always 1.
// The only exception is the number zero.
// For example, assuming we are using 64-bits for the significand, the number 1 is encoded as
// (sign = 0, significand = 0x800..0, exponent = -63)
// Note that, in this representation, the smallest positive integer is:
// (sign = 0, significand = 0x800..0, exponent = INT_MIN)
// instead of
// (sign = 0, significand = 0x000..1, exponent = INT_MIN)
//
// Remarks:
//
// - All values of type int, unsigned, int64_t and uint64_t can be precisely represented as mpff numerals.
//
// - Hardware float and double values (corresponding to rationals) can also be precisely represented as mpff numerals.
// That is, NaN, +oo and -oo are not supported by this module.
//
// - An exception (mpff_manager::exception) is thrown if overflow occurs. This can happen because the exponent is
// represented as a machine integer.
//
// - There are only two rounding modes: towards plus infinity and towards minus infinity.
// The rounding mode can be dynamically modified.
//
// - The mpff numerals are stored in a dynamic array.
// Type mpff is just an index (unsigned) into this array.
unsigned m_precision; //!< Number of words in the significand. Must be an even number.
unsigned m_precision_bits; //!< Number of bits in the significand. Must be 32*m_precision.
mutable unsigned_vector m_significands; //!< Array containing all significands.
unsigned m_capacity; //!< Number of significands that can be stored in m_significands.
bool m_to_plus_inf; //!< If True, then round to plus infinity, otherwise to minus infinity
id_gen m_id_gen;
static const unsigned MPFF_NUM_BUFFERS = 4;
svector m_buffers[MPFF_NUM_BUFFERS];
svector m_set_buffer;
mpff m_one;
mpn_manager m_mpn_manager;
unsigned * sig(mpff const & n) const { return m_significands.data() + (n.m_sig_idx * m_precision); }
void ensure_capacity(unsigned sig_idx) {
while (sig_idx >= m_capacity)
expand();
}
void expand();
void allocate_if_needed(mpff & n) {
if (n.m_sig_idx == 0)
allocate(n);
}
void allocate(mpff & n);
// copy n to buffer idx.
void to_buffer(unsigned idx, mpff const & n) const;
// copy n to buffer idx and add m_precision zeros.
void to_buffer_ext(unsigned idx, mpff const & n) const;
// copy (and shift by m_precision_bits) n to buffer idx
void to_buffer_shifting(unsigned idx, mpff const & n) const;
void inc_significand(unsigned * s, int64_t & exp);
void inc_significand(mpff & a);
void dec_significand(mpff & a);
bool min_significand(mpff const & a) const;
void set_min_significand(mpff & a);
void set_max_significand(mpff & a);
void set_big_exponent(mpff & a, int64_t e);
void set_exponent(mpff & a, int64_t e) {
if (e > INT_MAX || e < INT_MIN)
set_big_exponent(a, e);
else
a.m_exponent = static_cast(e);
}
template
void set_core(mpff & n, mpz_manager & m, mpz const & v);
template
void set_core(mpff & n, mpq_manager & m, mpq const & v);
template
void to_mpz_core(mpff const & n, mpz_manager & m, mpz & t);
template
void to_mpq_core(mpff const & n, mpq_manager & m, mpq & t);
template
void significand_core(mpff const & n, mpz_manager & m, mpz & r);
void add_sub(bool is_sub, mpff const & a, mpff const & b, mpff & c);
public:
typedef mpff numeral;
static bool precise() { return false; }
static bool field() { return true; }
class exception : public z3_exception {
char const * msg() const override { return "multi-precision floating point (mpff) exception"; }
};
class overflow_exception : public exception {
char const * msg() const override { return "multi-precision floating point (mpff) overflow"; }
};
class div0_exception : public exception {
char const * msg() const override { return "multi-precision floating point (mpff) division by zero"; }
};
mpff_manager(unsigned prec = 2, unsigned initial_capacity = 1024);
~mpff_manager();
void round_to_plus_inf() { m_to_plus_inf = true; }
void round_to_minus_inf() { m_to_plus_inf = false; }
void set_rounding(bool to_plus_inf) { m_to_plus_inf = to_plus_inf; }
bool rounding_to_plus_inf() const { return m_to_plus_inf; }
/**
\brief Return the exponent of n.
*/
static int exponent(mpff const & n) { return n.m_exponent; }
/**
\brief Update the exponent of n.
\remark It is a NOOP if n is zero.
*/
void set_exponent(mpff & n, int exp) { if (is_zero(n)) return; n.m_exponent = exp; SASSERT(check(n)); }
/**
\brief Return the significand as a mpz numeral.
*/
void significand(mpff const & n, unsynch_mpz_manager & m, mpz & r);
#ifndef SINGLE_THREAD
void significand(mpff const & n, synch_mpz_manager & m, mpz & r);
#endif
/**
\brief Return true if n is negative
*/
static bool sign(mpff const & n) { return is_neg(n); }
/**
\brief Set n to zero.
*/
void reset(mpff & n);
/**
\brief Return true if n is an integer.
*/
bool is_int(mpff const & n) const;
/**
\brief Return true if n is zero.
*/
static bool is_zero(mpff const & n) { return n.m_sig_idx == 0; }
/**
\brief Return true if n is positive.
*/
static bool is_pos(mpff const & n) { return n.m_sign == 0 && !is_zero(n); }
/**
\brief Return true if n is negative.
*/
static bool is_neg(mpff const & n) { return n.m_sign != 0; }
/**
\brief Return true if n is non positive.
*/
static bool is_nonpos(mpff const & n) { return !is_pos(n); }
/**
\brief Return true if n is non negative.
*/
static bool is_nonneg(mpff const & n) { return !is_neg(n); }
/**
\brief Return true if the absolute value of n is 1.
*/
bool is_abs_one(mpff const & n) const;
/**
\brief Return true if n is one.
*/
bool is_one(mpff const & n) const { return is_pos(n) && is_abs_one(n); }
/**
\brief Return true if n is minus one.
*/
bool is_minus_one(mpff const & n) const { return is_neg(n) && is_abs_one(n); }
/**
\brief Return true if n is two.
*/
bool is_two(mpff const & n) const;
/**
\brief Return true if \c a is the smallest representable negative number.
*/
bool is_minus_epsilon(mpff const & a) const;
/**
\brief Return true if \c a is the smallest representable positive number.
*/
bool is_plus_epsilon(mpff const & a) const;
/**
\brief Return true if \c a is an integer and fits in an int64_t machine integer.
*/
bool is_int64(mpff const & a) const;
/**
\brief Return true if \c a is a non-negative integer and fits in an int64_t machine integer.
*/
bool is_uint64(mpff const & a) const;
/**
\brief Delete the resources associated with n.
*/
void del(mpff & n);
/**
\brief a <- -a
*/
static void neg(mpff & a) { if (!is_zero(a)) a.m_sign = !a.m_sign; }
/**
\brief a <- |a|
*/
static void abs(mpff & a) { a.m_sign = 0; }
static void swap(mpff & a, mpff & b) noexcept { a.swap(b); }
/**
\brief c <- a + b
*/
void add(mpff const & a, mpff const & b, mpff & c);
/**
\brief c <- a - b
*/
void sub(mpff const & a, mpff const & b, mpff & c);
/**
\brief a <- a + 1
*/
void inc(mpff & a) { add(a, m_one, a); }
/**
\brief a <- a - 1
*/
void dec(mpff & a) { sub(a, m_one, a); }
/**
\brief c <- a * b
*/
void mul(mpff const & a, mpff const & b, mpff & c);
/**
\brief c <- a / b
\pre !is_zero(b)
*/
void div(mpff const & a, mpff const & b, mpff & c);
/**
\brief a <- 1/a
\pre !is_zero(a);
*/
void inv(mpff & a) { div(m_one, a, a); }
void inv(mpff const & a, mpff & b) { set(b, a); inv(b); }
/**
\brief b <- a^k
*/
void power(mpff const & a, unsigned k, mpff & b);
/**
\brief Return true if \c a is a power of 2. That is, a is equal to 2^k for some k >= 0.
*/
bool is_power_of_two(mpff const & a, unsigned & k) const;
bool is_power_of_two(mpff const & a) const;
bool eq(mpff const & a, mpff const & b) const;
bool neq(mpff const & a, mpff const & b) const { return !eq(a, b); }
bool lt(mpff const & a, mpff const & b) const;
bool gt(mpff const & a, mpff const & b) const { return lt(b, a); }
bool le(mpff const & a, mpff const & b) const { return !lt(b, a); }
bool ge(mpff const & a, mpff const & b) const { return !lt(a, b); }
void set(mpff & n, int v);
void set(mpff & n, unsigned v);
void set(mpff & n, int64_t v);
void set(mpff & n, uint64_t v);
void set(mpff & n, int num, unsigned den);
void set(mpff & n, int64_t num, uint64_t den);
void set(mpff & n, mpff const & v);
void set(mpff & n, unsynch_mpz_manager & m, mpz const & v);
void set(mpff & n, unsynch_mpq_manager & m, mpq const & v);
#ifndef SINGLE_THREAD
void set(mpff & n, synch_mpq_manager & m, mpq const & v);
void set(mpff & n, synch_mpz_manager & m, mpz const & v);
#endif
void set_plus_epsilon(mpff & n);
void set_minus_epsilon(mpff & n);
void set_max(mpff & n);
void set_min(mpff & n);
/**
\brief n <- floor(n)
*/
void floor(mpff & n);
void floor(mpff const & n, mpff & o) { set(o, n); floor(o); }
/**
\brief n <- ceil(n)
*/
void ceil(mpff & n);
void ceil(mpff const & n, mpff & o) { set(o, n); ceil(o); }
/**
\brief Update \c a to the next representable float.
Throws an exception if \c a is the maximal representable float.
*/
void next(mpff & a);
/**
\brief Update \c a to the previous representable float.
Throws an exception if \c a is the minimal representable float.
*/
void prev(mpff & a);
/**
\brief Convert n into a mpz numeral.
\pre is_int(n)
\remark if exponent(n) is too big, we may run out of memory.
*/
void to_mpz(mpff const & n, unsynch_mpz_manager & m, mpz & t);
#ifndef SINGLE_THREAD
/**
\brief Convert n into a mpz numeral.
\pre is_int(n)
\remark if exponent(n) is too big, we may run out of memory.
*/
void to_mpz(mpff const & n, synch_mpz_manager & m, mpz & t);
#endif
/**
\brief Convert n into a mpq numeral.
\remark if exponent(n) is too big, we may run out of memory.
*/
void to_mpq(mpff const & n, unsynch_mpq_manager & m, mpq & t);
#ifndef SINGLE_THREAD
/**
\brief Convert n into a mpq numeral.
\remark if exponent(n) is too big, we may run out of memory.
*/
void to_mpq(mpff const & n, synch_mpq_manager & m, mpq & t);
#endif
/**
\brief Return n as an int64.
\pre is_int64(n)
*/
int64_t get_int64(mpff const & n) const;
/**
\brief Return n as an uint64.
\pre is_uint64(n)
*/
uint64_t get_uint64(mpff const & n) const;
/**
\brief Return the biggest k s.t. 2^k <= a.
\remark Return 0 if a is not positive.
*/
unsigned prev_power_of_two(mpff const & a);
void display_raw(std::ostream & out, mpff const & n) const;
void display(std::ostream & out, mpff const & n) const;
void display_pp(std::ostream & out, mpff const & n) const { display(out, n); }
void display_decimal(std::ostream & out, mpff const & n, unsigned prec=32, unsigned max_power=128);
void display_smt2(std::ostream & out, mpff const & n, bool decimal=true) const;
std::string to_string(mpff const & a) const;
std::string to_rational_string(mpff const & a) const;
bool check(mpff const & n) const;
};
typedef _scoped_numeral scoped_mpff;
typedef _scoped_numeral_vector scoped_mpff_vector;
