GitHub-Proxy/eigen
1
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2025-09-14 18:33:16 +08:00
Code Issues Packages Projects Releases Wiki Activity

Update included mpreal header to 3.6.5 and fix deprecated warnings.

Browse Source
This commit is contained in:
Gael Guennebaud 2018-10-08 17:09:23 +02:00
parent 64b1a15318
commit e29bfe8479
1 changed files with 349 additions and 273 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

214
unsupported/test/mpreal/mpreal.h
View File

@ -5,7 +5,7 @@
Project homepage: http://www.holoborodko.com/pavel/mpfr Project homepage: http://www.holoborodko.com/pavel/mpfr
Contact e-mail: pavel@holoborodko.com Contact e-mail: pavel@holoborodko.com
Copyright (c) 2008-2015 Pavel Holoborodko Copyright (c) 2008-2016 Pavel Holoborodko
Contributors: Contributors:
Dmitriy Gubanov, Konstantin Holoborodko, Brian Gladman, Dmitriy Gubanov, Konstantin Holoborodko, Brian Gladman,
@ -13,7 +13,7 @@
Pere Constans, Peter van Hoof, Gael Guennebaud, Tsai Chia Cheng, Pere Constans, Peter van Hoof, Gael Guennebaud, Tsai Chia Cheng,
Alexei Zubanov, Jauhien Piatlicki, Victor Berger, John Westwood, Alexei Zubanov, Jauhien Piatlicki, Victor Berger, John Westwood,
Petr Aleksandrov, Orion Poplawski, Charles Karney, Arash Partow, Petr Aleksandrov, Orion Poplawski, Charles Karney, Arash Partow,
Rodney James, Jorge Leitao. Rodney James, Jorge Leitao, Jerome Benoit.
Licensing: Licensing:
(A) MPFR C++ is under GNU General Public License ("GPL"). (A) MPFR C++ is under GNU General Public License ("GPL").
@ -58,6 +58,7 @@
#include <limits> #include <limits>
#include <complex> #include <complex>
#include <algorithm> #include <algorithm>
#include <stdint.h>
// Options // Options
#define MPREAL_HAVE_MSVC_DEBUGVIEW // Enable Debugger Visualizer for "Debug" builds in MSVC. #define MPREAL_HAVE_MSVC_DEBUGVIEW // Enable Debugger Visualizer for "Debug" builds in MSVC.
@ -68,16 +69,18 @@
// Library version // Library version
#define MPREAL_VERSION_MAJOR 3 #define MPREAL_VERSION_MAJOR 3
#define MPREAL_VERSION_MINOR 6 #define MPREAL_VERSION_MINOR 6
#define MPREAL_VERSION_PATCHLEVEL 2 #define MPREAL_VERSION_PATCHLEVEL 5
#define MPREAL_VERSION_STRING "3.6.2" #define MPREAL_VERSION_STRING "3.6.5"
// Detect compiler using signatures from http://predef.sourceforge.net/ // Detect compiler using signatures from http://predef.sourceforge.net/
#if defined(__GNUC__) #if defined(__GNUC__) && defined(__INTEL_COMPILER)
#define IsInf(x) (isinf)(x) // GNU C++/Intel ICC compiler on Linux #define IsInf(x) isinf EIGEN_NOT_A_MACRO (x) // Intel ICC compiler on Linux
#elif defined(_MSC_VER) // Microsoft Visual C++ #elif defined(_MSC_VER) // Microsoft Visual C++
#define IsInf(x) (!_finite(x)) #define IsInf(x) (!_finite(x))
#else #else
#define IsInf(x) (std::isinf)(x) // GNU C/C++ (and/or other compilers), just hope for C99 conformance #define IsInf(x) std::isinf EIGEN_NOT_A_MACRO (x) // GNU C/C++ (and/or other compilers), just hope for C99 conformance
#endif #endif
// A Clang feature extension to determine compiler features. // A Clang feature extension to determine compiler features.
@ -85,10 +88,13 @@
#define __has_feature(x) 0 #define __has_feature(x) 0
#endif #endif
// Detect support for r-value references (move semantic). Borrowed from Eigen. // Detect support for r-value references (move semantic).
// Move semantic should be enabled with great care in multi-threading environments,
// especially if MPFR uses custom memory allocators.
// Everything should be thread-safe and support passing ownership over thread boundary.
#if (__has_feature(cxx_rvalue_references) || \ #if (__has_feature(cxx_rvalue_references) || \
defined(__GXX_EXPERIMENTAL_CXX0X__) || __cplusplus >= 201103L || \ defined(__GXX_EXPERIMENTAL_CXX0X__) || __cplusplus >= 201103L || \
(defined(_MSC_VER) && _MSC_VER >= 1600)) (defined(_MSC_VER) && _MSC_VER >= 1600) && !defined(MPREAL_DISABLE_MOVE_SEMANTIC))
#define MPREAL_HAVE_MOVE_SUPPORT #define MPREAL_HAVE_MOVE_SUPPORT
@ -99,8 +105,9 @@
// Detect support for explicit converters. // Detect support for explicit converters.
#if (__has_feature(cxx_explicit_conversions) || \ #if (__has_feature(cxx_explicit_conversions) || \
(defined(__GXX_EXPERIMENTAL_CXX0X__) && __GNUC_MINOR__ >= 5) || __cplusplus >= 201103L || \ (defined(__GXX_EXPERIMENTAL_CXX0X__) && __GNUC_MINOR >= 5) || __cplusplus >= 201103L || \
(defined(_MSC_VER) && _MSC_VER >= 1800)) (defined(_MSC_VER) && _MSC_VER >= 1800) || \
(defined(__INTEL_COMPILER) && __INTEL_COMPILER >= 1300))
#define MPREAL_HAVE_EXPLICIT_CONVERTERS #define MPREAL_HAVE_EXPLICIT_CONVERTERS
#endif #endif
@ -127,7 +134,7 @@
// = -1 disables overflow checks (default) // = -1 disables overflow checks (default)
// Fast replacement for mpfr_set_zero(x, +1): // Fast replacement for mpfr_set_zero(x, +1):
// (a) uses low-level data members, might not be compatible with new versions of MPFR // (a) uses low-level data members, might not be forward compatible
// (b) sign is not set, add (x)->_mpfr_sign = 1; // (b) sign is not set, add (x)->_mpfr_sign = 1;
#define mpfr_set_zero_fast(x) ((x)->_mpfr_exp = __MPFR_EXP_ZERO) #define mpfr_set_zero_fast(x) ((x)->_mpfr_exp = __MPFR_EXP_ZERO)
@ -147,7 +154,7 @@ public:
// Get default rounding mode & precision // Get default rounding mode & precision
inline static mp_rnd_t get_default_rnd() { return (mp_rnd_t)(mpfr_get_default_rounding_mode()); } inline static mp_rnd_t get_default_rnd() { return (mp_rnd_t)(mpfr_get_default_rounding_mode()); }
inline static mp_prec_t get_default_prec() { return mpfr_get_default_prec(); } inline static mp_prec_t get_default_prec() { return (mpfr_get_default_prec)(); }
// Constructors && type conversions // Constructors && type conversions
mpreal(); mpreal();
@ -354,6 +361,8 @@ public:
friend const mpreal log1p(const mpreal& v, mp_rnd_t rnd_mode); friend const mpreal log1p(const mpreal& v, mp_rnd_t rnd_mode);
friend const mpreal expm1(const mpreal& v, mp_rnd_t rnd_mode); friend const mpreal expm1(const mpreal& v, mp_rnd_t rnd_mode);
friend const mpreal nextpow2(const mpreal& v, mp_rnd_t rnd_mode);
friend const mpreal cos(const mpreal& v, mp_rnd_t rnd_mode); friend const mpreal cos(const mpreal& v, mp_rnd_t rnd_mode);
friend const mpreal sin(const mpreal& v, mp_rnd_t rnd_mode); friend const mpreal sin(const mpreal& v, mp_rnd_t rnd_mode);
friend const mpreal tan(const mpreal& v, mp_rnd_t rnd_mode); friend const mpreal tan(const mpreal& v, mp_rnd_t rnd_mode);
@ -405,7 +414,7 @@ public:
friend const mpreal fms (const mpreal& v1, const mpreal& v2, const mpreal& v3, mp_rnd_t rnd_mode); friend const mpreal fms (const mpreal& v1, const mpreal& v2, const mpreal& v3, mp_rnd_t rnd_mode);
friend const mpreal agm (const mpreal& v1, const mpreal& v2, mp_rnd_t rnd_mode); friend const mpreal agm (const mpreal& v1, const mpreal& v2, mp_rnd_t rnd_mode);
friend const mpreal sum (const mpreal tab[], const unsigned long int n, int& status, mp_rnd_t rnd_mode); friend const mpreal sum (const mpreal tab[], const unsigned long int n, int& status, mp_rnd_t rnd_mode);
friend int sgn(const mpreal& v); // returns -1 or +1 friend int sgn (const mpreal& v);
// MPFR 2.4.0 Specifics // MPFR 2.4.0 Specifics
#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0)) #if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
@ -482,7 +491,7 @@ public:
#endif #endif
// Instance Checkers // Instance Checkers
friend bool (isnan) (const mpreal& v); friend bool isnan EIGEN_NOT_A_MACRO (const mpreal& v);
friend bool (isinf) (const mpreal& v); friend bool (isinf) (const mpreal& v);
friend bool (isfinite) (const mpreal& v); friend bool (isfinite) (const mpreal& v);
@ -509,7 +518,7 @@ public:
mpreal& setSign (int Sign, mp_rnd_t RoundingMode = get_default_rnd()); mpreal& setSign (int Sign, mp_rnd_t RoundingMode = get_default_rnd());
//Exponent //Exponent
mp_exp_t get_exp(); mp_exp_t get_exp() const;
int set_exp(mp_exp_t e); int set_exp(mp_exp_t e);
int check_range (int t, mp_rnd_t rnd_mode = get_default_rnd()); int check_range (int t, mp_rnd_t rnd_mode = get_default_rnd());
int subnormalize (int t, mp_rnd_t rnd_mode = get_default_rnd()); int subnormalize (int t, mp_rnd_t rnd_mode = get_default_rnd());
@ -580,7 +589,7 @@ inline mpreal::mpreal(const mpreal& u)
#ifdef MPREAL_HAVE_MOVE_SUPPORT #ifdef MPREAL_HAVE_MOVE_SUPPORT
inline mpreal::mpreal(mpreal&& other) inline mpreal::mpreal(mpreal&& other)
{ {
mpfr_set_uninitialized(mpfr_ptr()); // make sure "other" holds no pointer to actual data mpfr_set_uninitialized(mpfr_ptr()); // make sure "other" holds null-pointer (in uninitialized state)
mpfr_swap(mpfr_ptr(), other.mpfr_ptr()); mpfr_swap(mpfr_ptr(), other.mpfr_ptr());
MPREAL_MSVC_DEBUGVIEW_CODE; MPREAL_MSVC_DEBUGVIEW_CODE;
@ -588,9 +597,11 @@ inline mpreal::mpreal(mpreal&& other)
inline mpreal& mpreal::operator=(mpreal&& other) inline mpreal& mpreal::operator=(mpreal&& other)
{ {
mpfr_swap(mpfr_ptr(), other.mpfr_ptr()); if (this != &other)
{
MPREAL_MSVC_DEBUGVIEW_CODE; mpfr_swap(mpfr_ptr(), other.mpfr_ptr()); // destructor for "other" will be called just afterwards
MPREAL_MSVC_DEBUGVIEW_CODE;
}
return *this; return *this;
} }
#endif #endif
@ -639,13 +650,13 @@ inline mpreal::mpreal(const double u, mp_prec_t prec, mp_rnd_t mode)
mpfr_init2(mpfr_ptr(), prec); mpfr_init2(mpfr_ptr(), prec);
#if (MPREAL_DOUBLE_BITS_OVERFLOW > -1) #if (MPREAL_DOUBLE_BITS_OVERFLOW > -1)
if(fits_in_bits(u, MPREAL_DOUBLE_BITS_OVERFLOW)) if(fits_in_bits(u, MPREAL_DOUBLE_BITS_OVERFLOW))
{ {
mpfr_set_d(mpfr_ptr(), u, mode); mpfr_set_d(mpfr_ptr(), u, mode);
}else }else
throw conversion_overflow(); throw conversion_overflow();
#else #else
mpfr_set_d(mpfr_ptr(), u, mode); mpfr_set_d(mpfr_ptr(), u, mode);
#endif #endif
MPREAL_MSVC_DEBUGVIEW_CODE; MPREAL_MSVC_DEBUGVIEW_CODE;
@ -908,13 +919,13 @@ inline mpreal& mpreal::operator=(const mpreal& v)
{ {
if (this != &v) if (this != &v)
{ {
mp_prec_t tp = mpfr_get_prec( mpfr_srcptr()); mp_prec_t tp = mpfr_get_prec( mpfr_srcptr());
mp_prec_t vp = mpfr_get_prec(v.mpfr_srcptr()); mp_prec_t vp = mpfr_get_prec(v.mpfr_srcptr());
if(tp != vp){ if(tp != vp){
clear(mpfr_ptr()); clear(mpfr_ptr());
mpfr_init2(mpfr_ptr(), vp); mpfr_init2(mpfr_ptr(), vp);
} }
mpfr_set(mpfr_ptr(), v.mpfr_srcptr(), mpreal::get_default_rnd()); mpfr_set(mpfr_ptr(), v.mpfr_srcptr(), mpreal::get_default_rnd());
@ -958,16 +969,16 @@ inline mpreal& mpreal::operator=(const long double v)
inline mpreal& mpreal::operator=(const double v) inline mpreal& mpreal::operator=(const double v)
{ {
#if (MPREAL_DOUBLE_BITS_OVERFLOW > -1) #if (MPREAL_DOUBLE_BITS_OVERFLOW > -1)
if(fits_in_bits(v, MPREAL_DOUBLE_BITS_OVERFLOW)) if(fits_in_bits(v, MPREAL_DOUBLE_BITS_OVERFLOW))
{ {
mpfr_set_d(mpfr_ptr(),v,mpreal::get_default_rnd()); mpfr_set_d(mpfr_ptr(),v,mpreal::get_default_rnd());
}else }else
throw conversion_overflow(); throw conversion_overflow();
#else #else
mpfr_set_d(mpfr_ptr(),v,mpreal::get_default_rnd()); mpfr_set_d(mpfr_ptr(),v,mpreal::get_default_rnd());
#endif #endif
MPREAL_MSVC_DEBUGVIEW_CODE; MPREAL_MSVC_DEBUGVIEW_CODE;
return *this; return *this;
} }
@ -1161,9 +1172,9 @@ inline const mpreal mpreal::operator+()const { return mpreal(*this); }
inline const mpreal operator+(const mpreal& a, const mpreal& b) inline const mpreal operator+(const mpreal& a, const mpreal& b)
{ {
mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_ptr()), mpfr_get_prec(b.mpfr_ptr()))); mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_ptr()), mpfr_get_prec(b.mpfr_ptr())));
mpfr_add(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd()); mpfr_add(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd());
return c; return c;
} }
inline mpreal& mpreal::operator++() inline mpreal& mpreal::operator++()
@ -1269,9 +1280,9 @@ inline const mpreal mpreal::operator-()const
inline const mpreal operator-(const mpreal& a, const mpreal& b) inline const mpreal operator-(const mpreal& a, const mpreal& b)
{ {
mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_ptr()), mpfr_get_prec(b.mpfr_ptr()))); mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_ptr()), mpfr_get_prec(b.mpfr_ptr())));
mpfr_sub(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd()); mpfr_sub(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd());
return c; return c;
} }
inline const mpreal operator-(const double b, const mpreal& a) inline const mpreal operator-(const double b, const mpreal& a)
@ -1386,9 +1397,9 @@ inline mpreal& mpreal::operator*=(const int v)
inline const mpreal operator*(const mpreal& a, const mpreal& b) inline const mpreal operator*(const mpreal& a, const mpreal& b)
{ {
mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_ptr()), mpfr_get_prec(b.mpfr_ptr()))); mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_ptr()), mpfr_get_prec(b.mpfr_ptr())));
mpfr_mul(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd()); mpfr_mul(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd());
return c; return c;
} }
////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////
@ -1462,9 +1473,9 @@ inline mpreal& mpreal::operator/=(const int v)
inline const mpreal operator/(const mpreal& a, const mpreal& b) inline const mpreal operator/(const mpreal& a, const mpreal& b)
{ {
mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_srcptr()), mpfr_get_prec(b.mpfr_srcptr()))); mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_srcptr()), mpfr_get_prec(b.mpfr_srcptr())));
mpfr_div(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd()); mpfr_div(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd());
return c; return c;
} }
inline const mpreal operator/(const unsigned long int b, const mpreal& a) inline const mpreal operator/(const unsigned long int b, const mpreal& a)
@ -1659,7 +1670,7 @@ inline bool operator > (const mpreal& a, const double b ){ return !
inline bool operator >= (const mpreal& a, const mpreal& b ){ return (mpfr_greaterequal_p(a.mpfr_srcptr(),b.mpfr_srcptr()) != 0 ); } inline bool operator >= (const mpreal& a, const mpreal& b ){ return (mpfr_greaterequal_p(a.mpfr_srcptr(),b.mpfr_srcptr()) != 0 ); }
inline bool operator >= (const mpreal& a, const unsigned long int b ){ return !isnan EIGEN_NOT_A_MACRO (a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) >= 0 ); } inline bool operator >= (const mpreal& a, const unsigned long int b ){ return !isnan EIGEN_NOT_A_MACRO (a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) >= 0 ); }
// inline bool operator >= (const mpreal& a, const unsigned int b ){ return !isnan EIGEN_NOT_A_MACRO (isnan()a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) >= 0 ); } inline bool operator >= (const mpreal& a, const unsigned int b ){ return !isnan EIGEN_NOT_A_MACRO (a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) >= 0 ); }
inline bool operator >= (const mpreal& a, const long int b ){ return !isnan EIGEN_NOT_A_MACRO (a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) >= 0 ); } inline bool operator >= (const mpreal& a, const long int b ){ return !isnan EIGEN_NOT_A_MACRO (a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) >= 0 ); }
inline bool operator >= (const mpreal& a, const int b ){ return !isnan EIGEN_NOT_A_MACRO (a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) >= 0 ); } inline bool operator >= (const mpreal& a, const int b ){ return !isnan EIGEN_NOT_A_MACRO (a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) >= 0 ); }
inline bool operator >= (const mpreal& a, const long double b ){ return !isnan EIGEN_NOT_A_MACRO (a) && (b == b) && (mpfr_cmp_ld(a.mpfr_srcptr(),b) >= 0 ); } inline bool operator >= (const mpreal& a, const long double b ){ return !isnan EIGEN_NOT_A_MACRO (a) && (b == b) && (mpfr_cmp_ld(a.mpfr_srcptr(),b) >= 0 ); }
@ -1697,7 +1708,7 @@ inline bool operator != (const mpreal& a, const int b ){ return !
inline bool operator != (const mpreal& a, const long double b ){ return !(a == b); } inline bool operator != (const mpreal& a, const long double b ){ return !(a == b); }
inline bool operator != (const mpreal& a, const double b ){ return !(a == b); } inline bool operator != (const mpreal& a, const double b ){ return !(a == b); }
inline bool (isnan) (const mpreal& op){ return (mpfr_nan_p (op.mpfr_srcptr()) != 0 ); } inline bool isnan EIGEN_NOT_A_MACRO (const mpreal& op){ return (mpfr_nan_p (op.mpfr_srcptr()) != 0 ); }
inline bool (isinf) (const mpreal& op){ return (mpfr_inf_p (op.mpfr_srcptr()) != 0 ); } inline bool (isinf) (const mpreal& op){ return (mpfr_inf_p (op.mpfr_srcptr()) != 0 ); }
inline bool (isfinite) (const mpreal& op){ return (mpfr_number_p (op.mpfr_srcptr()) != 0 ); } inline bool (isfinite) (const mpreal& op){ return (mpfr_number_p (op.mpfr_srcptr()) != 0 ); }
inline bool iszero (const mpreal& op){ return (mpfr_zero_p (op.mpfr_srcptr()) != 0 ); } inline bool iszero (const mpreal& op){ return (mpfr_zero_p (op.mpfr_srcptr()) != 0 ); }
@ -1762,7 +1773,7 @@ inline std::string mpreal::toString(int n, int b, mp_rnd_t mode) const
std::ostringstream format; std::ostringstream format;
int digits = (n >= 0) ? n : 1 + bits2digits(mpfr_get_prec(mpfr_srcptr())); int digits = (n >= 0) ? n : 2 + bits2digits(mpfr_get_prec(mpfr_srcptr()));
format << "%." << digits << "RNg"; format << "%." << digits << "RNg";
@ -1931,14 +1942,9 @@ inline int bits2digits(mp_prec_t b)
////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////
// Set/Get number properties // Set/Get number properties
inline int sgn(const mpreal& op)
{
return mpfr_sgn(op.mpfr_srcptr());
}
inline mpreal& mpreal::setSign(int sign, mp_rnd_t RoundingMode) inline mpreal& mpreal::setSign(int sign, mp_rnd_t RoundingMode)
{ {
mpfr_setsign(mpfr_ptr(), mpfr_srcptr(), (sign < 0 ? 1 : 0), RoundingMode); mpfr_setsign(mpfr_ptr(), mpfr_srcptr(), sign < 0, RoundingMode);
MPREAL_MSVC_DEBUGVIEW_CODE; MPREAL_MSVC_DEBUGVIEW_CODE;
return *this; return *this;
} }
@ -1993,7 +1999,7 @@ inline void mpreal::set_prec(mp_prec_t prec, mp_rnd_t rnd_mode)
MPREAL_MSVC_DEBUGVIEW_CODE; MPREAL_MSVC_DEBUGVIEW_CODE;
} }
inline mp_exp_t mpreal::get_exp () inline mp_exp_t mpreal::get_exp () const
{ {
return mpfr_get_exp(mpfr_srcptr()); return mpfr_get_exp(mpfr_srcptr());
} }
@ -2091,6 +2097,12 @@ inline bool signbit(const mpreal& x)
return mpfr_signbit(x.mpfr_srcptr()); return mpfr_signbit(x.mpfr_srcptr());
} }
inline mpreal& setsignbit(mpreal& x, bool minus, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
{
mpfr_setsign(x.mpfr_ptr(), x.mpfr_srcptr(), minus, rnd_mode);
return x;
}
inline const mpreal modf(const mpreal& v, mpreal& n) inline const mpreal modf(const mpreal& v, mpreal& n)
{ {
mpreal f(v); mpreal f(v);
@ -2194,7 +2206,7 @@ inline const mpreal sqrt(const int v, mp_rnd_t rnd_mode)
inline const mpreal root(const mpreal& x, unsigned long int k, mp_rnd_t r = mpreal::get_default_rnd()) inline const mpreal root(const mpreal& x, unsigned long int k, mp_rnd_t r = mpreal::get_default_rnd())
{ {
mpreal y(0, mpfr_get_prec(x.mpfr_srcptr())); mpreal y(0, mpfr_get_prec(x.mpfr_srcptr()));
mpfr_root(y.mpfr_ptr(), x.mpfr_srcptr(), k, r); mpfr_rootn_ui(y.mpfr_ptr(), x.mpfr_srcptr(), k, r);
return y; return y;
} }
@ -2270,6 +2282,16 @@ inline const mpreal besselj1(const mpreal& x, mp_rnd_t r = mpreal::get_default_r
inline const mpreal bessely0(const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_UNARY_MATH_FUNCTION_BODY(y0 ); } inline const mpreal bessely0(const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_UNARY_MATH_FUNCTION_BODY(y0 ); }
inline const mpreal bessely1(const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_UNARY_MATH_FUNCTION_BODY(y1 ); } inline const mpreal bessely1(const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_UNARY_MATH_FUNCTION_BODY(y1 ); }
inline const mpreal nextpow2(const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd())
{
mpreal y(0, x.getPrecision());
if(!iszero(x))
y = ceil(log2(abs(x,r),r));
return y;
}
inline const mpreal atan2 (const mpreal& y, const mpreal& x, mp_rnd_t rnd_mode = mpreal::get_default_rnd()) inline const mpreal atan2 (const mpreal& y, const mpreal& x, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
{ {
mpreal a(0,(std::max)(y.getPrecision(), x.getPrecision())); mpreal a(0,(std::max)(y.getPrecision(), x.getPrecision()));
@ -2284,6 +2306,44 @@ inline const mpreal hypot (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode =
return a; return a;
} }
inline const mpreal hypot(const mpreal& a, const mpreal& b, const mpreal& c)
{
if(isnan EIGEN_NOT_A_MACRO (a) || isnan EIGEN_NOT_A_MACRO (b) || isnan EIGEN_NOT_A_MACRO(c)) return mpreal().setNan();
else
{
mpreal absa = abs(a), absb = abs(b), absc = abs(c);
mpreal w = (std::max)(absa, (std::max)(absb, absc));
mpreal r;
if (!iszero(w))
{
mpreal iw = 1/w;
r = w * sqrt(sqr(absa*iw) + sqr(absb*iw) + sqr(absc*iw));
}
return r;
}
}
inline const mpreal hypot(const mpreal& a, const mpreal& b, const mpreal& c, const mpreal& d)
{
if(isnan EIGEN_NOT_A_MACRO (a) || isnan EIGEN_NOT_A_MACRO (b) || isnan EIGEN_NOT_A_MACRO (c) || isnan EIGEN_NOT_A_MACRO (d)) return mpreal().setNan();
else
{
mpreal absa = abs(a), absb = abs(b), absc = abs(c), absd = abs(d);
mpreal w = (std::max)(absa, (std::max)(absb, (std::max)(absc, absd)));
mpreal r;
if (!iszero(w))
{
mpreal iw = 1/w;
r = w * sqrt(sqr(absa*iw) + sqr(absb*iw) + sqr(absc*iw) + sqr(absd*iw));
}
return r;
}
}
inline const mpreal remainder (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd()) inline const mpreal remainder (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
{ {
mpreal a(0,(std::max)(y.getPrecision(), x.getPrecision())); mpreal a(0,(std::max)(y.getPrecision(), x.getPrecision()));
@ -2299,7 +2359,7 @@ inline const mpreal remquo (long* q, const mpreal& x, const mpreal& y, mp_rnd_t
} }
inline const mpreal fac_ui (unsigned long int v, mp_prec_t prec = mpreal::get_default_prec(), inline const mpreal fac_ui (unsigned long int v, mp_prec_t prec = mpreal::get_default_prec(),
mp_rnd_t rnd_mode = mpreal::get_default_rnd()) mp_rnd_t rnd_mode = mpreal::get_default_rnd())
{ {
mpreal x(0, prec); mpreal x(0, prec);
mpfr_fac_ui(x.mpfr_ptr(),v,rnd_mode); mpfr_fac_ui(x.mpfr_ptr(),v,rnd_mode);
@ -2432,9 +2492,7 @@ inline const mpreal mod (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = m
mpreal m = x - floor(x / y) * y; mpreal m = x - floor(x / y) * y;
m.setSign(sgn(y)); // make sure result has the same sign as Y return copysign(abs(m),y); // make sure result has the same sign as Y
return m;
} }
inline const mpreal fmod (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd()) inline const mpreal fmod (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
@ -2543,7 +2601,16 @@ inline const mpreal frac (const mpreal& x, mp_rnd_t r = mpreal::get_defaul
////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////
// Miscellaneous Functions // Miscellaneous Functions
inline void swap (mpreal& a, mpreal& b) { mpfr_swap(a.mp,b.mp); } inline int sgn(const mpreal& op)
{
// Please note, this is classic signum function which ignores sign of zero.
// Use signbit if you need sign of zero.
return mpfr_sgn(op.mpfr_srcptr());
}
//////////////////////////////////////////////////////////////////////////
// Miscellaneous Functions
inline void swap (mpreal& a, mpreal& b) { mpfr_swap(a.mpfr_ptr(),b.mpfr_ptr()); }
inline const mpreal (max)(const mpreal& x, const mpreal& y){ return (x>y?x:y); } inline const mpreal (max)(const mpreal& x, const mpreal& y){ return (x>y?x:y); }
inline const mpreal (min)(const mpreal& x, const mpreal& y){ return (x<y?x:y); } inline const mpreal (min)(const mpreal& x, const mpreal& y){ return (x<y?x:y); }
@ -2634,14 +2701,23 @@ inline const mpreal random(unsigned int seed = 0)
} }
#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,1,0)) #if (MPFR_VERSION >= MPFR_VERSION_NUM(3,1,0) && MPFR_VERSION < MPFR_VERSION_NUM(4,0,0))
// TODO:
// Use mpfr_nrandom since mpfr_grandom is deprecated
#if defined(_MSC_VER)
#pragma warning( push )
#pragma warning( disable : 1478)
#endif
inline const mpreal grandom (gmp_randstate_t& state, mp_rnd_t rnd_mode = mpreal::get_default_rnd()) inline const mpreal grandom (gmp_randstate_t& state, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
{ {
mpreal x; mpreal x;
mpfr_grandom(x.mpfr_ptr(), NULL, state, rnd_mode); mpfr_grandom(x.mpfr_ptr(), NULL, state, rnd_mode);
return x; return x;
} }
#if defined(_MSC_VER)
#pragma warning( pop )
#endif
inline const mpreal grandom(unsigned int seed = 0) inline const mpreal grandom(unsigned int seed = 0)
{ {
@ -2981,7 +3057,7 @@ inline const mpreal pow(const double a, const int b, mp_rnd_t rnd_mode)
// Non-throwing swap C++ idiom: http://en.wikibooks.org/wiki/More_C%2B%2B_Idioms/Non-throwing_swap // Non-throwing swap C++ idiom: http://en.wikibooks.org/wiki/More_C%2B%2B_Idioms/Non-throwing_swap
namespace std namespace std
{ {
// we are allowed to extend namespace std with specializations only // we are allowed to extend namespace std with specializations only
template <> template <>
inline void swap(mpfr::mpreal& x, mpfr::mpreal& y) inline void swap(mpfr::mpreal& x, mpfr::mpreal& y)
{ {