Macaulay2 Engine
Loading...
Searching...
No Matches
mpreal.h
Go to the documentation of this file.
1/*
2 MPFR C++: Multi-precision floating point number class for C++.
3 Based on MPFR library: http://mpfr.org
4
5 Project homepage: http://www.holoborodko.com/pavel/mpfr
6 Contact e-mail: pavel@holoborodko.com
7
8 Copyright (c) 2008-2022 Pavel Holoborodko
9
10 Contributors:
11 Dmitriy Gubanov, Konstantin Holoborodko, Brian Gladman,
12 Helmut Jarausch, Fokko Beekhof, Ulrich Mutze, Heinz van Saanen,
13 Pere Constans, Peter van Hoof, Gael Guennebaud, Tsai Chia Cheng,
14 Alexei Zubanov, Jauhien Piatlicki, Victor Berger, John Westwood,
15 Petr Aleksandrov, Orion Poplawski, Charles Karney, Arash Partow,
16 Rodney James, Jorge Leitao, Jerome Benoit, Michal Maly.
17
18 Licensing:
19 (A) MPFR C++ is under GNU General Public License ("GPL").
20
21 (B) Non-free licenses may also be purchased from the author, for users who
22 do not want their programs protected by the GPL.
23
24 The non-free licenses are for users that wish to use MPFR C++ in
25 their products but are unwilling to release their software
26 under the GPL (which would require them to release source code
27 and allow free redistribution).
28
29 Such users can purchase an unlimited-use license from the author.
30 Contact us for more details.
31
32 GNU General Public License ("GPL") copyright permissions statement:
33 **************************************************************************
34 This program is free software: you can redistribute it and/or modify
35 it under the terms of the GNU General Public License as published by
36 the Free Software Foundation, either version 3 of the License, or
37 (at your option) any later version.
38
39 This program is distributed in the hope that it will be useful,
40 but WITHOUT ANY WARRANTY; without even the implied warranty of
41 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
42 GNU General Public License for more details.
43
44 You should have received a copy of the GNU General Public License
45 along with this program. If not, see <http://www.gnu.org/licenses/>.
46*/
47
48#ifndef __MPREAL_H__
49#define __MPREAL_H__
50
75
76#include <string>
77#include <iostream>
78#include <sstream>
79#include <stdexcept>
80#include <cfloat>
81#include <cmath>
82#include <cstring>
83#include <limits>
84#include <complex>
85#include <algorithm>
86#include <stdint.h>
87
88// Options
89#define MPREAL_HAVE_MSVC_DEBUGVIEW // Enable Debugger Visualizer for "Debug" builds in MSVC.
90#define MPREAL_HAVE_DYNAMIC_STD_NUMERIC_LIMITS // Enable extended std::numeric_limits<mpfr::mpreal> specialization.
91 // Meaning that "digits", "round_style" and similar members are defined as functions, not constants.
92 // See std::numeric_limits<mpfr::mpreal> at the end of the file for more information.
93
94// Library version
95#define MPREAL_VERSION_MAJOR 3
96#define MPREAL_VERSION_MINOR 6
97#define MPREAL_VERSION_PATCHLEVEL 9
98#define MPREAL_VERSION_STRING "3.6.9"
99
100// Detect compiler using signatures from http://predef.sourceforge.net/
101#if defined(__GNUC__) && defined(__INTEL_COMPILER)
102 #define IsInf(x) isinf(x) // Intel ICC compiler on Linux
103
104#elif defined(_MSC_VER) // Microsoft Visual C++
105 #define IsInf(x) (!_finite(x))
106
107#else
108 #define IsInf(x) std::isinf(x) // GNU C/C++ (and/or other compilers), just hope for C99 conformance
109#endif
110
111// A Clang feature extension to determine compiler features.
112#ifndef __has_feature
113 #define __has_feature(x) 0
114#endif
115
116// Detect support for r-value references (move semantic).
117// Move semantic should be enabled with great care in multi-threading environments,
118// especially if MPFR uses custom memory allocators.
119// Everything should be thread-safe and support passing ownership over thread boundary.
120#if (__has_feature(cxx_rvalue_references) || \
121 defined(__GXX_EXPERIMENTAL_CXX0X__) || __cplusplus >= 201103L || \
122 (defined(_MSC_VER) && _MSC_VER >= 1600) && !defined(MPREAL_DISABLE_MOVE_SEMANTIC))
123
124 #define MPREAL_HAVE_MOVE_SUPPORT
125
126 // Use fields in mpfr_t structure to check if it was initialized / set dummy initialization
127 #define mpfr_is_initialized(x) (0 != (x)->_mpfr_d)
128 #define mpfr_set_uninitialized(x) ((x)->_mpfr_d = 0 )
129#endif
130
131// Detect support for explicit converters.
132#if (__has_feature(cxx_explicit_conversions) || \
133 (defined(__GXX_EXPERIMENTAL_CXX0X__) && __GNUC_MINOR >= 5) || __cplusplus >= 201103L || \
134 (defined(_MSC_VER) && _MSC_VER >= 1800) || \
135 (defined(__INTEL_COMPILER) && __INTEL_COMPILER >= 1300))
136
137 #define MPREAL_HAVE_EXPLICIT_CONVERTERS
138#endif
139
140#define MPFR_USE_INTMAX_T // Enable 64-bit integer types - should be defined before mpfr.h
141
142#if defined(MPREAL_HAVE_MSVC_DEBUGVIEW) && defined(_MSC_VER) && defined(_DEBUG)
143 #define MPREAL_MSVC_DEBUGVIEW_CODE DebugView = toString();
144 #define MPREAL_MSVC_DEBUGVIEW_DATA std::string DebugView;
145#else
146 #define MPREAL_MSVC_DEBUGVIEW_CODE
147 #define MPREAL_MSVC_DEBUGVIEW_DATA
148#endif
149
150#define MPFR_USE_NO_MACRO
151#include <mpfr.h>
152
153#if (MPFR_VERSION < MPFR_VERSION_NUM(3,0,0))
154 #include <cstdlib> // Needed for random()
155#endif
156
157// Less important options
158#define MPREAL_DOUBLE_BITS_OVERFLOW -1 // Triggers overflow exception during conversion to double if mpreal
159 // cannot fit in MPREAL_DOUBLE_BITS_OVERFLOW bits
160 // = -1 disables overflow checks (default)
161
162// Fast replacement for mpfr_set_zero(x, +1):
163// (a) uses low-level data members, might not be forward compatible
164// (b) sign is not set, add (x)->_mpfr_sign = 1;
165#define mpfr_set_zero_fast(x) ((x)->_mpfr_exp = __MPFR_EXP_ZERO)
166
167#if defined(__GNUC__)
168 #define MPREAL_PERMISSIVE_EXPR __extension__
169#else
170 #define MPREAL_PERMISSIVE_EXPR
171#endif
172
173namespace mpfr {
174
175class mpreal {
176private:
177 mpfr_t mp;
178
179public:
180
181 // Get default rounding mode & precision
182 inline static mp_rnd_t get_default_rnd() { return (mp_rnd_t)(mpfr_get_default_rounding_mode()); }
183 inline static mp_prec_t get_default_prec() { return (mpfr_get_default_prec)(); }
184
185 // Constructors && type conversions
186 mpreal();
187 mpreal(const mpreal& u);
188 mpreal(const mpf_t u);
189 mpreal(const mpz_t u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
190 mpreal(const mpq_t u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
191 mpreal(const double u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
192 mpreal(const long double u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
193 mpreal(const unsigned long long int u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
194 mpreal(const long long int u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
195 mpreal(const unsigned long int u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
196 mpreal(const unsigned int u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
197 mpreal(const long int u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
198 mpreal(const int u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
199
200 // Construct mpreal from mpfr_t structure.
201 // shared = true allows to avoid deep copy, so that mpreal and 'u' share the same data & pointers.
202 mpreal(const mpfr_t u, bool shared = false);
203
204 mpreal(const char* s, mp_prec_t prec = mpreal::get_default_prec(), int base = 10, mp_rnd_t mode = mpreal::get_default_rnd());
205 mpreal(const std::string& s, mp_prec_t prec = mpreal::get_default_prec(), int base = 10, mp_rnd_t mode = mpreal::get_default_rnd());
206
207 ~mpreal();
208
209#ifdef MPREAL_HAVE_MOVE_SUPPORT
210 mpreal& operator=(mpreal&& v);
211 mpreal(mpreal&& u);
212#endif
213
214 // Operations
215 // =
216 // +, -, *, /, ++, --, <<, >>
217 // *=, +=, -=, /=,
218 // <, >, ==, <=, >=
219
220 // =
221 mpreal& operator=(const mpreal& v);
222 mpreal& operator=(const mpf_t v);
223 mpreal& operator=(const mpz_t v);
224 mpreal& operator=(const mpq_t v);
225 mpreal& operator=(const long double v);
226 mpreal& operator=(const double v);
227 mpreal& operator=(const unsigned long int v);
228 mpreal& operator=(const unsigned long long int v);
229 mpreal& operator=(const long long int v);
230 mpreal& operator=(const unsigned int v);
231 mpreal& operator=(const long int v);
232 mpreal& operator=(const int v);
233 mpreal& operator=(const char* s);
234 mpreal& operator=(const std::string& s);
235 template <typename real_t> mpreal& operator= (const std::complex<real_t>& z);
236
237 // +
238 mpreal& operator+=(const mpreal& v);
239 mpreal& operator+=(const mpf_t v);
240 mpreal& operator+=(const mpz_t v);
241 mpreal& operator+=(const mpq_t v);
242 mpreal& operator+=(const long double u);
243 mpreal& operator+=(const double u);
244 mpreal& operator+=(const unsigned long int u);
245 mpreal& operator+=(const unsigned int u);
246 mpreal& operator+=(const long int u);
247 mpreal& operator+=(const int u);
248
249 mpreal& operator+=(const long long int u);
250 mpreal& operator+=(const unsigned long long int u);
251 mpreal& operator-=(const long long int u);
252 mpreal& operator-=(const unsigned long long int u);
253 mpreal& operator*=(const long long int u);
254 mpreal& operator*=(const unsigned long long int u);
255 mpreal& operator/=(const long long int u);
256 mpreal& operator/=(const unsigned long long int u);
257
258 const mpreal operator+() const;
259 mpreal& operator++ ();
260 const mpreal operator++ (int);
261
262 // -
263 mpreal& operator-=(const mpreal& v);
264 mpreal& operator-=(const mpz_t v);
265 mpreal& operator-=(const mpq_t v);
266 mpreal& operator-=(const long double u);
267 mpreal& operator-=(const double u);
268 mpreal& operator-=(const unsigned long int u);
269 mpreal& operator-=(const unsigned int u);
270 mpreal& operator-=(const long int u);
271 mpreal& operator-=(const int u);
272 const mpreal operator-() const;
273 friend const mpreal operator-(const unsigned long int b, const mpreal& a);
274 friend const mpreal operator-(const unsigned int b, const mpreal& a);
275 friend const mpreal operator-(const long int b, const mpreal& a);
276 friend const mpreal operator-(const int b, const mpreal& a);
277 friend const mpreal operator-(const double b, const mpreal& a);
278 mpreal& operator-- ();
279 const mpreal operator-- (int);
280
281 // *
282 mpreal& operator*=(const mpreal& v);
283 mpreal& operator*=(const mpz_t v);
284 mpreal& operator*=(const mpq_t v);
285 mpreal& operator*=(const long double v);
286 mpreal& operator*=(const double v);
287 mpreal& operator*=(const unsigned long int v);
288 mpreal& operator*=(const unsigned int v);
289 mpreal& operator*=(const long int v);
290 mpreal& operator*=(const int v);
291
292 // /
293 mpreal& operator/=(const mpreal& v);
294 mpreal& operator/=(const mpz_t v);
295 mpreal& operator/=(const mpq_t v);
296 mpreal& operator/=(const long double v);
297 mpreal& operator/=(const double v);
298 mpreal& operator/=(const unsigned long int v);
299 mpreal& operator/=(const unsigned int v);
300 mpreal& operator/=(const long int v);
301 mpreal& operator/=(const int v);
302 friend const mpreal operator/(const unsigned long int b, const mpreal& a);
303 friend const mpreal operator/(const unsigned int b, const mpreal& a);
304 friend const mpreal operator/(const long int b, const mpreal& a);
305 friend const mpreal operator/(const int b, const mpreal& a);
306 friend const mpreal operator/(const double b, const mpreal& a);
307
308 //<<= Fast Multiplication by 2^u
309 mpreal& operator<<=(const unsigned long int u);
310 mpreal& operator<<=(const unsigned int u);
311 mpreal& operator<<=(const long int u);
312 mpreal& operator<<=(const int u);
313
314 //>>= Fast Division by 2^u
315 mpreal& operator>>=(const unsigned long int u);
316 mpreal& operator>>=(const unsigned int u);
317 mpreal& operator>>=(const long int u);
318 mpreal& operator>>=(const int u);
319
320 // Type Conversion operators
321 bool toBool ( ) const;
322 long toLong (mp_rnd_t mode = GMP_RNDZ) const;
323 unsigned long toULong (mp_rnd_t mode = GMP_RNDZ) const;
324 long long toLLong (mp_rnd_t mode = GMP_RNDZ) const;
325 unsigned long long toULLong (mp_rnd_t mode = GMP_RNDZ) const;
326 float toFloat (mp_rnd_t mode = GMP_RNDN) const;
327 double toDouble (mp_rnd_t mode = GMP_RNDN) const;
328 long double toLDouble (mp_rnd_t mode = GMP_RNDN) const;
329
330#if defined (MPREAL_HAVE_EXPLICIT_CONVERTERS)
331 explicit operator bool () const { return toBool(); }
332 explicit operator signed char () const { return (signed char)toLong(); }
333 explicit operator unsigned char () const { return (unsigned char)toULong(); }
334 explicit operator short () const { return (short)toLong(); }
335 explicit operator unsigned short () const { return (unsigned short)toULong();}
336 explicit operator int () const { return (int)toLong(); }
337 explicit operator unsigned int () const { return (unsigned int)toULong(); }
338 explicit operator long () const { return toLong(); }
339 explicit operator unsigned long () const { return toULong(); }
340 explicit operator long long () const { return toLLong(); }
341 explicit operator unsigned long long () const { return toULLong(); }
342 explicit operator float () const { return toFloat(); }
343 explicit operator double () const { return toDouble(); }
344 explicit operator long double () const { return toLDouble(); }
345#endif
346
347 // Get raw pointers so that mpreal can be directly used in raw mpfr_* functions
348 ::mpfr_ptr mpfr_ptr();
349 ::mpfr_srcptr mpfr_ptr() const;
350 ::mpfr_srcptr mpfr_srcptr() const;
351
352 // Convert mpreal to string with n significant digits in base b
353 // n = -1 -> convert with the maximum available digits
354 std::string toString(int n = -1, int b = 10, mp_rnd_t mode = mpreal::get_default_rnd()) const;
355
356#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
357 std::string toString(const std::string& format) const;
358#endif
359
360 std::ostream& output(std::ostream& os) const;
361
362 // Math Functions
363 friend const mpreal sqr (const mpreal& v, mp_rnd_t rnd_mode);
364 friend const mpreal sqrt(const mpreal& v, mp_rnd_t rnd_mode);
365 friend const mpreal sqrt(const unsigned long int v, mp_rnd_t rnd_mode);
366 friend const mpreal cbrt(const mpreal& v, mp_rnd_t rnd_mode);
367 friend const mpreal root(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode);
368 friend const mpreal pow (const mpreal& a, const mpreal& b, mp_rnd_t rnd_mode);
369 friend const mpreal pow (const mpreal& a, const mpz_t b, mp_rnd_t rnd_mode);
370 friend const mpreal pow (const mpreal& a, const unsigned long int b, mp_rnd_t rnd_mode);
371 friend const mpreal pow (const mpreal& a, const long int b, mp_rnd_t rnd_mode);
372 friend const mpreal pow (const unsigned long int a, const mpreal& b, mp_rnd_t rnd_mode);
373 friend const mpreal pow (const unsigned long int a, const unsigned long int b, mp_rnd_t rnd_mode);
374 friend const mpreal fabs(const mpreal& v, mp_rnd_t rnd_mode);
375
376 friend const mpreal abs(const mpreal& v, mp_rnd_t rnd_mode);
377 friend const mpreal dim(const mpreal& a, const mpreal& b, mp_rnd_t rnd_mode);
378 friend inline const mpreal mul_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode);
379 friend inline const mpreal mul_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode);
380 friend inline const mpreal div_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode);
381 friend inline const mpreal div_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode);
382 friend int cmpabs(const mpreal& a,const mpreal& b);
383
384 friend const mpreal log (const mpreal& v, mp_rnd_t rnd_mode);
385 friend const mpreal log2 (const mpreal& v, mp_rnd_t rnd_mode);
386 friend const mpreal logb (const mpreal& v, mp_rnd_t rnd_mode);
387 friend mp_exp_t ilogb(const mpreal& v);
388 friend const mpreal log10(const mpreal& v, mp_rnd_t rnd_mode);
389 friend const mpreal exp (const mpreal& v, mp_rnd_t rnd_mode);
390 friend const mpreal exp2 (const mpreal& v, mp_rnd_t rnd_mode);
391 friend const mpreal exp10(const mpreal& v, mp_rnd_t rnd_mode);
392 friend const mpreal log1p(const mpreal& v, mp_rnd_t rnd_mode);
393 friend const mpreal expm1(const mpreal& v, mp_rnd_t rnd_mode);
394
395 friend const mpreal nextpow2(const mpreal& v, mp_rnd_t rnd_mode);
396
397 friend const mpreal cos(const mpreal& v, mp_rnd_t rnd_mode);
398 friend const mpreal sin(const mpreal& v, mp_rnd_t rnd_mode);
399 friend const mpreal tan(const mpreal& v, mp_rnd_t rnd_mode);
400 friend const mpreal sec(const mpreal& v, mp_rnd_t rnd_mode);
401 friend const mpreal csc(const mpreal& v, mp_rnd_t rnd_mode);
402 friend const mpreal cot(const mpreal& v, mp_rnd_t rnd_mode);
403 friend int sin_cos(mpreal& s, mpreal& c, const mpreal& v, mp_rnd_t rnd_mode);
404
405 friend const mpreal acos (const mpreal& v, mp_rnd_t rnd_mode);
406 friend const mpreal asin (const mpreal& v, mp_rnd_t rnd_mode);
407 friend const mpreal atan (const mpreal& v, mp_rnd_t rnd_mode);
408 friend const mpreal atan2 (const mpreal& y, const mpreal& x, mp_rnd_t rnd_mode);
409 friend const mpreal acot (const mpreal& v, mp_rnd_t rnd_mode);
410 friend const mpreal asec (const mpreal& v, mp_rnd_t rnd_mode);
411 friend const mpreal acsc (const mpreal& v, mp_rnd_t rnd_mode);
412
413 friend const mpreal cosh (const mpreal& v, mp_rnd_t rnd_mode);
414 friend const mpreal sinh (const mpreal& v, mp_rnd_t rnd_mode);
415 friend const mpreal tanh (const mpreal& v, mp_rnd_t rnd_mode);
416 friend const mpreal sech (const mpreal& v, mp_rnd_t rnd_mode);
417 friend const mpreal csch (const mpreal& v, mp_rnd_t rnd_mode);
418 friend const mpreal coth (const mpreal& v, mp_rnd_t rnd_mode);
419 friend const mpreal acosh (const mpreal& v, mp_rnd_t rnd_mode);
420 friend const mpreal asinh (const mpreal& v, mp_rnd_t rnd_mode);
421 friend const mpreal atanh (const mpreal& v, mp_rnd_t rnd_mode);
422 friend const mpreal acoth (const mpreal& v, mp_rnd_t rnd_mode);
423 friend const mpreal asech (const mpreal& v, mp_rnd_t rnd_mode);
424 friend const mpreal acsch (const mpreal& v, mp_rnd_t rnd_mode);
425
426 friend const mpreal hypot (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode);
427
428 friend const mpreal fac_ui (unsigned long int v, mp_prec_t prec, mp_rnd_t rnd_mode);
429 friend const mpreal eint (const mpreal& v, mp_rnd_t rnd_mode);
430
431 friend const mpreal gamma (const mpreal& v, mp_rnd_t rnd_mode);
432 friend const mpreal tgamma (const mpreal& v, mp_rnd_t rnd_mode);
433 friend const mpreal lngamma (const mpreal& v, mp_rnd_t rnd_mode);
434 friend const mpreal lgamma (const mpreal& v, int *signp, mp_rnd_t rnd_mode);
435 friend const mpreal zeta (const mpreal& v, mp_rnd_t rnd_mode);
436 friend const mpreal erf (const mpreal& v, mp_rnd_t rnd_mode);
437 friend const mpreal erfc (const mpreal& v, mp_rnd_t rnd_mode);
438 friend const mpreal besselj0 (const mpreal& v, mp_rnd_t rnd_mode);
439 friend const mpreal besselj1 (const mpreal& v, mp_rnd_t rnd_mode);
440 friend const mpreal besseljn (long n, const mpreal& v, mp_rnd_t rnd_mode);
441 friend const mpreal bessely0 (const mpreal& v, mp_rnd_t rnd_mode);
442 friend const mpreal bessely1 (const mpreal& v, mp_rnd_t rnd_mode);
443 friend const mpreal besselyn (long n, const mpreal& v, mp_rnd_t rnd_mode);
444 friend const mpreal fma (const mpreal& v1, const mpreal& v2, const mpreal& v3, mp_rnd_t rnd_mode);
445 friend const mpreal fms (const mpreal& v1, const mpreal& v2, const mpreal& v3, mp_rnd_t rnd_mode);
446 friend const mpreal agm (const mpreal& v1, const mpreal& v2, mp_rnd_t rnd_mode);
447 friend const mpreal sum (const mpreal tab[], const unsigned long int n, int& status, mp_rnd_t rnd_mode);
448 friend int sgn (const mpreal& v);
449
450// MPFR 2.4.0 Specifics
451#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
452 friend int sinh_cosh (mpreal& s, mpreal& c, const mpreal& v, mp_rnd_t rnd_mode);
453 friend const mpreal li2 (const mpreal& v, mp_rnd_t rnd_mode);
454 friend const mpreal fmod (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode);
455 friend const mpreal rec_sqrt (const mpreal& v, mp_rnd_t rnd_mode);
456
457 // MATLAB's semantic equivalents
458 friend const mpreal rem (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode); // Remainder after division
459 friend const mpreal mod (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode); // Modulus after division
460#endif
461
462#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,0,0))
463 friend const mpreal digamma (const mpreal& v, mp_rnd_t rnd_mode);
464 friend const mpreal ai (const mpreal& v, mp_rnd_t rnd_mode);
465 friend const mpreal urandom (gmp_randstate_t& state, mp_rnd_t rnd_mode); // use gmp_randinit_default() to init state, gmp_randclear() to clear
466#endif
467
468#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,1,0))
469 friend const mpreal grandom (gmp_randstate_t& state, mp_rnd_t rnd_mode); // use gmp_randinit_default() to init state, gmp_randclear() to clear
470 friend const mpreal grandom (unsigned int seed);
471#endif
472
473 // Uniformly distributed random number generation in [0,1] using
474 // Mersenne-Twister algorithm by default.
475 // Use parameter to setup seed, e.g.: random((unsigned)time(NULL))
476 // Check urandom() for more precise control.
477 friend const mpreal random(unsigned int seed);
478
479 // Splits mpreal value into fractional and integer parts.
480 // Returns fractional part and stores integer part in n.
481 friend const mpreal modf(const mpreal& v, mpreal& n);
482
483 // Constants
484 // don't forget to call mpfr_free_cache() for every thread where you are using const-functions
485 friend const mpreal const_log2 (mp_prec_t prec, mp_rnd_t rnd_mode);
486 friend const mpreal const_pi (mp_prec_t prec, mp_rnd_t rnd_mode);
487 friend const mpreal const_euler (mp_prec_t prec, mp_rnd_t rnd_mode);
488 friend const mpreal const_catalan (mp_prec_t prec, mp_rnd_t rnd_mode);
489
490 // returns +inf iff sign>=0 otherwise -inf
491 friend const mpreal const_infinity(int sign, mp_prec_t prec);
492
493 // Output/ Input
494 friend std::ostream& operator<<(std::ostream& os, const mpreal& v);
495 friend std::istream& operator>>(std::istream& is, mpreal& v);
496
497 // Integer Related Functions
498 friend const mpreal rint (const mpreal& v, mp_rnd_t rnd_mode);
499 friend const mpreal ceil (const mpreal& v);
500 friend const mpreal floor(const mpreal& v);
501 friend const mpreal round(const mpreal& v);
502 friend long lround(const mpreal& v);
503 friend long long llround(const mpreal& v);
504 friend const mpreal trunc(const mpreal& v);
505 friend const mpreal rint_ceil (const mpreal& v, mp_rnd_t rnd_mode);
506 friend const mpreal rint_floor (const mpreal& v, mp_rnd_t rnd_mode);
507 friend const mpreal rint_round (const mpreal& v, mp_rnd_t rnd_mode);
508 friend const mpreal rint_trunc (const mpreal& v, mp_rnd_t rnd_mode);
509 friend const mpreal frac (const mpreal& v, mp_rnd_t rnd_mode);
510 friend const mpreal remainder (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode);
511 friend const mpreal remquo (const mpreal& x, const mpreal& y, int* q, mp_rnd_t rnd_mode);
512
513 // Miscellaneous Functions
514 friend const mpreal nexttoward (const mpreal& x, const mpreal& y);
515 friend const mpreal nextabove (const mpreal& x);
516 friend const mpreal nextbelow (const mpreal& x);
517
518 // use gmp_randinit_default() to init state, gmp_randclear() to clear
519 friend const mpreal urandomb (gmp_randstate_t& state);
520
521// MPFR < 2.4.2 Specifics
522#if (MPFR_VERSION <= MPFR_VERSION_NUM(2,4,2))
523 friend const mpreal random2 (mp_size_t size, mp_exp_t exp);
524#endif
525
526 // Instance Checkers
527 friend bool isnan (const mpreal& v);
528 friend bool isinf (const mpreal& v);
529 friend bool isfinite (const mpreal& v);
530
531 friend bool isnum (const mpreal& v);
532 friend bool iszero (const mpreal& v);
533 friend bool isint (const mpreal& v);
534
535#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,0,0))
536 friend bool isregular(const mpreal& v);
537#endif
538
539 // Set/Get instance properties
540 inline mp_prec_t get_prec() const;
541 inline void set_prec(mp_prec_t prec, mp_rnd_t rnd_mode = get_default_rnd()); // Change precision with rounding mode
542
543 // Aliases for get_prec(), set_prec() - needed for compatibility with std::complex<mpreal> interface
544 inline mpreal& setPrecision(int Precision, mp_rnd_t RoundingMode = get_default_rnd());
545 inline int getPrecision() const;
546
547 // Set mpreal to +/- inf, NaN, +/-0
548 mpreal& setInf (int Sign = +1);
549 mpreal& setNan ();
550 mpreal& setZero (int Sign = +1);
551 mpreal& setSign (int Sign, mp_rnd_t RoundingMode = get_default_rnd());
552
553 //Exponent
554 mp_exp_t get_exp() const;
555 int set_exp(mp_exp_t e);
556 int check_range (int t, mp_rnd_t rnd_mode = get_default_rnd());
557 int subnormalize (int t, mp_rnd_t rnd_mode = get_default_rnd());
558
559 // Inexact conversion from float
560 inline bool fits_in_bits(double x, int n);
561
562 // Set/Get global properties
563 static void set_default_prec(mp_prec_t prec);
564 static void set_default_rnd(mp_rnd_t rnd_mode);
565
566 static mp_exp_t get_emin (void);
567 static mp_exp_t get_emax (void);
568 static mp_exp_t get_emin_min (void);
569 static mp_exp_t get_emin_max (void);
570 static mp_exp_t get_emax_min (void);
571 static mp_exp_t get_emax_max (void);
572 static int set_emin (mp_exp_t exp);
573 static int set_emax (mp_exp_t exp);
574
575 // Efficient swapping of two mpreal values - needed for std algorithms
576 friend void swap(mpreal& x, mpreal& y);
577
578 friend const mpreal fmax(const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode);
579 friend const mpreal fmin(const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode);
580
581private:
582 // Human friendly Debug Preview in Visual Studio.
583 // Put one of these lines:
584 //
585 // mpfr::mpreal=<DebugView> ; Show value only
586 // mpfr::mpreal=<DebugView>, <mp[0]._mpfr_prec,u>bits ; Show value & precision
587 //
588 // at the beginning of
589 // [Visual Studio Installation Folder]\Common7\Packages\Debugger\autoexp.dat
591
592 // "Smart" resources deallocation. Checks if instance initialized before deletion.
593 void clear(::mpfr_ptr);
594};
595
597// Exceptions
598class conversion_overflow : public std::exception {
599public:
600 std::string why() { return "inexact conversion from floating point"; }
601};
602
604// Constructors & converters
605// Default constructor: creates mp number and initializes it to 0.
613
614inline mpreal::mpreal(const mpreal& u)
615{
616 mpfr_init2(mpfr_ptr(),mpfr_get_prec(u.mpfr_srcptr()));
618
620}
621
622#ifdef MPREAL_HAVE_MOVE_SUPPORT
623inline mpreal::mpreal(mpreal&& other)
624{
625 mpfr_set_uninitialized(mpfr_ptr()); // make sure "other" holds null-pointer (in uninitialized state)
626 mpfr_swap(mpfr_ptr(), other.mpfr_ptr());
627
629}
630
631inline mpreal& mpreal::operator=(mpreal&& other)
632{
633 if (this != &other)
634 {
635 mpfr_swap(mpfr_ptr(), other.mpfr_ptr()); // destructor for "other" will be called just afterwards
637 }
638 return *this;
639}
640#endif
641
642inline mpreal::mpreal(const mpfr_t u, bool shared)
643{
644 if(shared)
645 {
646 std::memcpy(mpfr_ptr(), u, sizeof(mpfr_t));
647 }
648 else
649 {
650 mpfr_init2(mpfr_ptr(), mpfr_get_prec(u));
651 mpfr_set (mpfr_ptr(), u, mpreal::get_default_rnd());
652 }
653
655}
656
657inline mpreal::mpreal(const mpf_t u)
658{
659 mpfr_init2(mpfr_ptr(),(mp_prec_t) mpf_get_prec(u)); // (gmp: mp_bitcnt_t) unsigned long -> long (mpfr: mp_prec_t)
660 mpfr_set_f(mpfr_ptr(),u,mpreal::get_default_rnd());
661
663}
664
665inline mpreal::mpreal(const mpz_t u, mp_prec_t prec, mp_rnd_t mode)
666{
667 mpfr_init2(mpfr_ptr(), prec);
668 mpfr_set_z(mpfr_ptr(), u, mode);
669
671}
672
673inline mpreal::mpreal(const mpq_t u, mp_prec_t prec, mp_rnd_t mode)
674{
675 mpfr_init2(mpfr_ptr(), prec);
676 mpfr_set_q(mpfr_ptr(), u, mode);
677
679}
680
681inline mpreal::mpreal(const double u, mp_prec_t prec, mp_rnd_t mode)
682{
683 mpfr_init2(mpfr_ptr(), prec);
684
685#if (MPREAL_DOUBLE_BITS_OVERFLOW > -1)
687 {
688 mpfr_set_d(mpfr_ptr(), u, mode);
689 }else
690 throw conversion_overflow();
691#else
692 mpfr_set_d(mpfr_ptr(), u, mode);
693#endif
694
696}
697
698inline mpreal::mpreal(const long double u, mp_prec_t prec, mp_rnd_t mode)
699{
700 mpfr_init2 (mpfr_ptr(), prec);
701 mpfr_set_ld(mpfr_ptr(), u, mode);
702
704}
705
706inline mpreal::mpreal(const unsigned long long int u, mp_prec_t prec, mp_rnd_t mode)
707{
708 mpfr_init2 (mpfr_ptr(), prec);
709 mpfr_set_uj(mpfr_ptr(), u, mode);
710
712}
713
714inline mpreal::mpreal(const long long int u, mp_prec_t prec, mp_rnd_t mode)
715{
716 mpfr_init2 (mpfr_ptr(), prec);
717 mpfr_set_sj(mpfr_ptr(), u, mode);
718
720}
721
722inline mpreal::mpreal(const unsigned long int u, mp_prec_t prec, mp_rnd_t mode)
723{
724 mpfr_init2 (mpfr_ptr(), prec);
725 mpfr_set_ui(mpfr_ptr(), u, mode);
726
728}
729
730inline mpreal::mpreal(const unsigned int u, mp_prec_t prec, mp_rnd_t mode)
731{
732 mpfr_init2 (mpfr_ptr(), prec);
733 mpfr_set_ui(mpfr_ptr(), u, mode);
734
736}
737
738inline mpreal::mpreal(const long int u, mp_prec_t prec, mp_rnd_t mode)
739{
740 mpfr_init2 (mpfr_ptr(), prec);
741 mpfr_set_si(mpfr_ptr(), u, mode);
742
744}
745
746inline mpreal::mpreal(const int u, mp_prec_t prec, mp_rnd_t mode)
747{
748 mpfr_init2 (mpfr_ptr(), prec);
749 mpfr_set_si(mpfr_ptr(), u, mode);
750
752}
753
754inline mpreal::mpreal(const char* s, mp_prec_t prec, int base, mp_rnd_t mode)
755{
756 mpfr_init2 (mpfr_ptr(), prec);
757 mpfr_set_str(mpfr_ptr(), s, base, mode);
758
760}
761
762inline mpreal::mpreal(const std::string& s, mp_prec_t prec, int base, mp_rnd_t mode)
763{
764 mpfr_init2 (mpfr_ptr(), prec);
765 mpfr_set_str(mpfr_ptr(), s.c_str(), base, mode);
766
768}
769
771{
772#ifdef MPREAL_HAVE_MOVE_SUPPORT
773 if(mpfr_is_initialized(x))
774#endif
775 mpfr_clear(x);
776}
777
779{
780 clear(mpfr_ptr());
781}
782
783// internal namespace needed for template magic
784namespace internal{
785
786 // Use SFINAE to restrict arithmetic operations instantiation only for numeric types
787 // This is needed for smooth integration with libraries based on expression templates, like Eigen.
788 // TODO: Do the same for boolean operators.
789 template <typename ArgumentType> struct result_type {};
790
791 template <> struct result_type<mpreal> {typedef mpreal type;};
792 template <> struct result_type<mpz_t> {typedef mpreal type;};
793 template <> struct result_type<mpq_t> {typedef mpreal type;};
794 template <> struct result_type<long double> {typedef mpreal type;};
795 template <> struct result_type<double> {typedef mpreal type;};
796 template <> struct result_type<unsigned long int> {typedef mpreal type;};
797 template <> struct result_type<unsigned int> {typedef mpreal type;};
798 template <> struct result_type<long int> {typedef mpreal type;};
799 template <> struct result_type<int> {typedef mpreal type;};
800 template <> struct result_type<long long> {typedef mpreal type;};
801 template <> struct result_type<unsigned long long> {typedef mpreal type;};
802}
803
804// + Addition
805template <typename Rhs>
806inline const typename internal::result_type<Rhs>::type
807 operator+(const mpreal& lhs, const Rhs& rhs){ return mpreal(lhs) += rhs; }
808
809template <typename Lhs>
810inline const typename internal::result_type<Lhs>::type
811 operator+(const Lhs& lhs, const mpreal& rhs){ return mpreal(rhs) += lhs; }
812
813// - Subtraction
814template <typename Rhs>
815inline const typename internal::result_type<Rhs>::type
816 operator-(const mpreal& lhs, const Rhs& rhs){ return mpreal(lhs) -= rhs; }
817
818template <typename Lhs>
819inline const typename internal::result_type<Lhs>::type
820 operator-(const Lhs& lhs, const mpreal& rhs){ return mpreal(lhs) -= rhs; }
821
822// * Multiplication
823template <typename Rhs>
824inline const typename internal::result_type<Rhs>::type
825 operator*(const mpreal& lhs, const Rhs& rhs){ return mpreal(lhs) *= rhs; }
826
827template <typename Lhs>
828inline const typename internal::result_type<Lhs>::type
829 operator*(const Lhs& lhs, const mpreal& rhs){ return mpreal(rhs) *= lhs; }
830
831// / Division
832template <typename Rhs>
833inline const typename internal::result_type<Rhs>::type
834 operator/(const mpreal& lhs, const Rhs& rhs){ return mpreal(lhs) /= rhs; }
835
836template <typename Lhs>
837inline const typename internal::result_type<Lhs>::type
838 operator/(const Lhs& lhs, const mpreal& rhs){ return mpreal(lhs) /= rhs; }
839
841// sqrt
842const mpreal sqrt(const unsigned int v, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
843const mpreal sqrt(const long int v, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
844const mpreal sqrt(const int v, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
845const mpreal sqrt(const long double v, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
846const mpreal sqrt(const double v, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
847
848// abs
849inline const mpreal abs(const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd());
850
852// pow
853const mpreal pow(const mpreal& a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
854const mpreal pow(const mpreal& a, const int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
855const mpreal pow(const mpreal& a, const long double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
856const mpreal pow(const mpreal& a, const double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
857
858const mpreal pow(const unsigned int a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
859const mpreal pow(const long int a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
860const mpreal pow(const int a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
861const mpreal pow(const long double a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
862const mpreal pow(const double a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
863
864const mpreal pow(const unsigned long int a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
865const mpreal pow(const unsigned long int a, const long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
866const mpreal pow(const unsigned long int a, const int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
867const mpreal pow(const unsigned long int a, const long double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
868const mpreal pow(const unsigned long int a, const double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
869
870const mpreal pow(const unsigned int a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
871const mpreal pow(const unsigned int a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
872const mpreal pow(const unsigned int a, const long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
873const mpreal pow(const unsigned int a, const int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
874const mpreal pow(const unsigned int a, const long double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
875const mpreal pow(const unsigned int a, const double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
876
877const mpreal pow(const long int a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
878const mpreal pow(const long int a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
879const mpreal pow(const long int a, const long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
880const mpreal pow(const long int a, const int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
881const mpreal pow(const long int a, const long double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
882const mpreal pow(const long int a, const double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
883
884const mpreal pow(const int a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
885const mpreal pow(const int a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
886const mpreal pow(const int a, const long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
887const mpreal pow(const int a, const int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
888const mpreal pow(const int a, const long double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
889const mpreal pow(const int a, const double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
890
891const mpreal pow(const long double a, const long double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
892const mpreal pow(const long double a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
893const mpreal pow(const long double a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
894const mpreal pow(const long double a, const long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
895const mpreal pow(const long double a, const int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
896
897const mpreal pow(const double a, const double b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
898const mpreal pow(const double a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
899const mpreal pow(const double a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
900const mpreal pow(const double a, const long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
901const mpreal pow(const double a, const int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
902
903inline const mpreal mul_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
904inline const mpreal mul_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
905inline const mpreal div_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
906inline const mpreal div_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode = mpreal::get_default_rnd());
907
909// Estimate machine epsilon for the given precision
910// Returns smallest eps such that 1.0 + eps != 1.0
911inline mpreal machine_epsilon(mp_prec_t prec = mpreal::get_default_prec());
912
913// Returns smallest eps such that x + eps != x (relative machine epsilon)
914inline mpreal machine_epsilon(const mpreal& x);
915
916// Gives max & min values for the required precision,
917// minval is 'safe' meaning 1 / minval does not overflow
918// maxval is 'safe' meaning 1 / maxval does not underflow
919inline mpreal minval(mp_prec_t prec = mpreal::get_default_prec());
920inline mpreal maxval(mp_prec_t prec = mpreal::get_default_prec());
921
922// 'Dirty' equality check 1: |a-b| < min{|a|,|b|} * eps
923inline bool isEqualFuzzy(const mpreal& a, const mpreal& b, const mpreal& eps);
924
925// 'Dirty' equality check 2: |a-b| < min{|a|,|b|} * eps( min{|a|,|b|} )
926inline bool isEqualFuzzy(const mpreal& a, const mpreal& b);
927
928// 'Bitwise' equality check
929// maxUlps - a and b can be apart by maxUlps binary numbers.
930inline bool isEqualUlps(const mpreal& a, const mpreal& b, int maxUlps);
931
933// Convert precision in 'bits' to decimal digits and vice versa.
934// bits = ceil(digits*log[2](10))
935// digits = floor(bits*log[10](2))
936
937inline mp_prec_t digits2bits(int d);
938inline int bits2digits(mp_prec_t b);
939
941// min, max
942const mpreal (max)(const mpreal& x, const mpreal& y);
943const mpreal (min)(const mpreal& x, const mpreal& y);
944
946// Implementation
948
950// Operators - Assignment
952{
953 if (this != &v)
954 {
955 mp_prec_t tp = mpfr_get_prec( mpfr_srcptr());
956 mp_prec_t vp = mpfr_get_prec(v.mpfr_srcptr());
957
958 if(tp != vp){
959 clear(mpfr_ptr());
960 mpfr_init2(mpfr_ptr(), vp);
961 }
962
963 mpfr_set(mpfr_ptr(), v.mpfr_srcptr(), mpreal::get_default_rnd());
964
966 }
967 return *this;
968}
969
970inline mpreal& mpreal::operator=(const mpf_t v)
971{
972 mpfr_set_f(mpfr_ptr(), v, mpreal::get_default_rnd());
973
975 return *this;
976}
977
978inline mpreal& mpreal::operator=(const mpz_t v)
979{
980 mpfr_set_z(mpfr_ptr(), v, mpreal::get_default_rnd());
981
983 return *this;
984}
985
986inline mpreal& mpreal::operator=(const mpq_t v)
987{
988 mpfr_set_q(mpfr_ptr(), v, mpreal::get_default_rnd());
989
991 return *this;
992}
993
994inline mpreal& mpreal::operator=(const long double v)
995{
996 mpfr_set_ld(mpfr_ptr(), v, mpreal::get_default_rnd());
997
999 return *this;
1000}
1001
1002inline mpreal& mpreal::operator=(const double v)
1003{
1004#if (MPREAL_DOUBLE_BITS_OVERFLOW > -1)
1006 {
1007 mpfr_set_d(mpfr_ptr(),v,mpreal::get_default_rnd());
1008 }else
1009 throw conversion_overflow();
1010#else
1011 mpfr_set_d(mpfr_ptr(),v,mpreal::get_default_rnd());
1012#endif
1013
1015 return *this;
1016}
1017
1018inline mpreal& mpreal::operator=(const unsigned long int v)
1019{
1020 mpfr_set_ui(mpfr_ptr(), v, mpreal::get_default_rnd());
1021
1023 return *this;
1024}
1025
1026inline mpreal& mpreal::operator=(const unsigned int v)
1027{
1028 mpfr_set_ui(mpfr_ptr(), v, mpreal::get_default_rnd());
1029
1031 return *this;
1032}
1033
1034inline mpreal& mpreal::operator=(const unsigned long long int v)
1035{
1036 mpfr_set_uj(mpfr_ptr(), v, mpreal::get_default_rnd());
1037
1039 return *this;
1040}
1041
1042inline mpreal& mpreal::operator=(const long long int v)
1043{
1044 mpfr_set_sj(mpfr_ptr(), v, mpreal::get_default_rnd());
1045
1047 return *this;
1048}
1049
1050inline mpreal& mpreal::operator=(const long int v)
1051{
1052 mpfr_set_si(mpfr_ptr(), v, mpreal::get_default_rnd());
1053
1055 return *this;
1056}
1057
1058inline mpreal& mpreal::operator=(const int v)
1059{
1060 mpfr_set_si(mpfr_ptr(), v, mpreal::get_default_rnd());
1061
1063 return *this;
1064}
1065
1066inline mpreal& mpreal::operator=(const char* s)
1067{
1068 // Use other converters for more precise control on base & precision & rounding:
1069 //
1070 // mpreal(const char* s, mp_prec_t prec, int base, mp_rnd_t mode)
1071 // mpreal(const std::string& s,mp_prec_t prec, int base, mp_rnd_t mode)
1072 //
1073 // Here we assume base = 10 and we use precision of target variable.
1074
1075 mpfr_t t;
1076
1077 mpfr_init2(t, mpfr_get_prec(mpfr_srcptr()));
1078
1079 if(0 == mpfr_set_str(t, s, 10, mpreal::get_default_rnd()))
1080 {
1081 mpfr_set(mpfr_ptr(), t, mpreal::get_default_rnd());
1083 }
1084
1085 clear(t);
1086 return *this;
1087}
1088
1089inline mpreal& mpreal::operator=(const std::string& s)
1090{
1091 // Use other converters for more precise control on base & precision & rounding:
1092 //
1093 // mpreal(const char* s, mp_prec_t prec, int base, mp_rnd_t mode)
1094 // mpreal(const std::string& s,mp_prec_t prec, int base, mp_rnd_t mode)
1095 //
1096 // Here we assume base = 10 and we use precision of target variable.
1097
1098 mpfr_t t;
1099
1100 mpfr_init2(t, mpfr_get_prec(mpfr_srcptr()));
1101
1102 if(0 == mpfr_set_str(t, s.c_str(), 10, mpreal::get_default_rnd()))
1103 {
1104 mpfr_set(mpfr_ptr(), t, mpreal::get_default_rnd());
1106 }
1107
1108 clear(t);
1109 return *this;
1110}
1111
1112template <typename real_t>
1113inline mpreal& mpreal::operator= (const std::complex<real_t>& z)
1114{
1115 return *this = z.real();
1116}
1117
1119// + Addition
1121{
1124 return *this;
1125}
1126
1127inline mpreal& mpreal::operator+=(const mpf_t u)
1128{
1129 *this += mpreal(u);
1131 return *this;
1132}
1133
1134inline mpreal& mpreal::operator+=(const mpz_t u)
1135{
1136 mpfr_add_z(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1138 return *this;
1139}
1140
1141inline mpreal& mpreal::operator+=(const mpq_t u)
1142{
1143 mpfr_add_q(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1145 return *this;
1146}
1147
1148inline mpreal& mpreal::operator+= (const long double u)
1149{
1150 *this += mpreal(u);
1152 return *this;
1153}
1154
1155inline mpreal& mpreal::operator+= (const double u)
1156{
1157#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1158 mpfr_add_d(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1159#else
1160 *this += mpreal(u);
1161#endif
1162
1164 return *this;
1165}
1166
1167inline mpreal& mpreal::operator+=(const unsigned long int u)
1168{
1169 mpfr_add_ui(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1171 return *this;
1172}
1173
1174inline mpreal& mpreal::operator+=(const unsigned int u)
1175{
1176 mpfr_add_ui(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1178 return *this;
1179}
1180
1181inline mpreal& mpreal::operator+=(const long int u)
1182{
1183 mpfr_add_si(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1185 return *this;
1186}
1187
1188inline mpreal& mpreal::operator+=(const int u)
1189{
1190 mpfr_add_si(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1192 return *this;
1193}
1194
1195inline mpreal& mpreal::operator+=(const long long int u) { *this += mpreal(u); MPREAL_MSVC_DEBUGVIEW_CODE; return *this; }
1196inline mpreal& mpreal::operator+=(const unsigned long long int u){ *this += mpreal(u); MPREAL_MSVC_DEBUGVIEW_CODE; return *this; }
1197inline mpreal& mpreal::operator-=(const long long int u) { *this -= mpreal(u); MPREAL_MSVC_DEBUGVIEW_CODE; return *this; }
1198inline mpreal& mpreal::operator-=(const unsigned long long int u){ *this -= mpreal(u); MPREAL_MSVC_DEBUGVIEW_CODE; return *this; }
1199inline mpreal& mpreal::operator*=(const long long int u) { *this *= mpreal(u); MPREAL_MSVC_DEBUGVIEW_CODE; return *this; }
1200inline mpreal& mpreal::operator*=(const unsigned long long int u){ *this *= mpreal(u); MPREAL_MSVC_DEBUGVIEW_CODE; return *this; }
1201inline mpreal& mpreal::operator/=(const long long int u) { *this /= mpreal(u); MPREAL_MSVC_DEBUGVIEW_CODE; return *this; }
1202inline mpreal& mpreal::operator/=(const unsigned long long int u){ *this /= mpreal(u); MPREAL_MSVC_DEBUGVIEW_CODE; return *this; }
1203
1204inline const mpreal mpreal::operator+()const { return mpreal(*this); }
1205
1206inline const mpreal operator+(const mpreal& a, const mpreal& b)
1207{
1208 mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_ptr()), mpfr_get_prec(b.mpfr_ptr())));
1209 mpfr_add(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd());
1210 return c;
1211}
1212
1214{
1215 return *this += 1;
1216}
1217
1218inline const mpreal mpreal::operator++ (int)
1219{
1220 mpreal x(*this);
1221 *this += 1;
1222 return x;
1223}
1224
1226{
1227 return *this -= 1;
1228}
1229
1230inline const mpreal mpreal::operator-- (int)
1231{
1232 mpreal x(*this);
1233 *this -= 1;
1234 return x;
1235}
1236
1238// - Subtraction
1240{
1243 return *this;
1244}
1245
1246inline mpreal& mpreal::operator-=(const mpz_t v)
1247{
1248 mpfr_sub_z(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1250 return *this;
1251}
1252
1253inline mpreal& mpreal::operator-=(const mpq_t v)
1254{
1255 mpfr_sub_q(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1257 return *this;
1258}
1259
1260inline mpreal& mpreal::operator-=(const long double v)
1261{
1262 *this -= mpreal(v);
1264 return *this;
1265}
1266
1267inline mpreal& mpreal::operator-=(const double v)
1268{
1269#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1270 mpfr_sub_d(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1271#else
1272 *this -= mpreal(v);
1273#endif
1274
1276 return *this;
1277}
1278
1279inline mpreal& mpreal::operator-=(const unsigned long int v)
1280{
1281 mpfr_sub_ui(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1283 return *this;
1284}
1285
1286inline mpreal& mpreal::operator-=(const unsigned int v)
1287{
1288 mpfr_sub_ui(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1290 return *this;
1291}
1292
1293inline mpreal& mpreal::operator-=(const long int v)
1294{
1295 mpfr_sub_si(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1297 return *this;
1298}
1299
1300inline mpreal& mpreal::operator-=(const int v)
1301{
1302 mpfr_sub_si(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1304 return *this;
1305}
1306
1307inline const mpreal mpreal::operator-()const
1308{
1309 mpreal u(*this);
1310 mpfr_neg(u.mpfr_ptr(),u.mpfr_srcptr(),mpreal::get_default_rnd());
1311 return u;
1312}
1313
1314inline const mpreal operator-(const mpreal& a, const mpreal& b)
1315{
1316 mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_ptr()), mpfr_get_prec(b.mpfr_ptr())));
1317 mpfr_sub(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd());
1318 return c;
1319}
1320
1321inline const mpreal operator-(const double b, const mpreal& a)
1322{
1323#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1324 mpreal x(0, mpfr_get_prec(a.mpfr_ptr()));
1325 mpfr_d_sub(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1326 return x;
1327#else
1328 mpreal x(b, mpfr_get_prec(a.mpfr_ptr()));
1329 x -= a;
1330 return x;
1331#endif
1332}
1333
1334inline const mpreal operator-(const unsigned long int b, const mpreal& a)
1335{
1336 mpreal x(0, mpfr_get_prec(a.mpfr_ptr()));
1337 mpfr_ui_sub(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1338 return x;
1339}
1340
1341inline const mpreal operator-(const unsigned int b, const mpreal& a)
1342{
1343 mpreal x(0, mpfr_get_prec(a.mpfr_ptr()));
1344 mpfr_ui_sub(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1345 return x;
1346}
1347
1348inline const mpreal operator-(const long int b, const mpreal& a)
1349{
1350 mpreal x(0, mpfr_get_prec(a.mpfr_ptr()));
1351 mpfr_si_sub(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1352 return x;
1353}
1354
1355inline const mpreal operator-(const int b, const mpreal& a)
1356{
1357 mpreal x(0, mpfr_get_prec(a.mpfr_ptr()));
1358 mpfr_si_sub(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1359 return x;
1360}
1361
1363// * Multiplication
1364inline mpreal& mpreal::operator*= (const mpreal& v)
1365{
1368 return *this;
1369}
1370
1371inline mpreal& mpreal::operator*=(const mpz_t v)
1372{
1373 mpfr_mul_z(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1375 return *this;
1376}
1377
1378inline mpreal& mpreal::operator*=(const mpq_t v)
1379{
1380 mpfr_mul_q(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1382 return *this;
1383}
1384
1385inline mpreal& mpreal::operator*=(const long double v)
1386{
1387 *this *= mpreal(v);
1389 return *this;
1390}
1391
1392inline mpreal& mpreal::operator*=(const double v)
1393{
1394#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1395 mpfr_mul_d(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1396#else
1397 *this *= mpreal(v);
1398#endif
1400 return *this;
1401}
1402
1403inline mpreal& mpreal::operator*=(const unsigned long int v)
1404{
1405 mpfr_mul_ui(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1407 return *this;
1408}
1409
1410inline mpreal& mpreal::operator*=(const unsigned int v)
1411{
1412 mpfr_mul_ui(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1414 return *this;
1415}
1416
1417inline mpreal& mpreal::operator*=(const long int v)
1418{
1419 mpfr_mul_si(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1421 return *this;
1422}
1423
1424inline mpreal& mpreal::operator*=(const int v)
1425{
1426 mpfr_mul_si(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1428 return *this;
1429}
1430
1431inline const mpreal operator*(const mpreal& a, const mpreal& b)
1432{
1433 mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_ptr()), mpfr_get_prec(b.mpfr_ptr())));
1434 mpfr_mul(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd());
1435 return c;
1436}
1437
1439// / Division
1441{
1444 return *this;
1445}
1446
1447inline mpreal& mpreal::operator/=(const mpz_t v)
1448{
1449 mpfr_div_z(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1451 return *this;
1452}
1453
1454inline mpreal& mpreal::operator/=(const mpq_t v)
1455{
1456 mpfr_div_q(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1458 return *this;
1459}
1460
1461inline mpreal& mpreal::operator/=(const long double v)
1462{
1463 *this /= mpreal(v);
1465 return *this;
1466}
1467
1468inline mpreal& mpreal::operator/=(const double v)
1469{
1470#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1471 mpfr_div_d(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1472#else
1473 *this /= mpreal(v);
1474#endif
1476 return *this;
1477}
1478
1479inline mpreal& mpreal::operator/=(const unsigned long int v)
1480{
1481 mpfr_div_ui(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1483 return *this;
1484}
1485
1486inline mpreal& mpreal::operator/=(const unsigned int v)
1487{
1488 mpfr_div_ui(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1490 return *this;
1491}
1492
1493inline mpreal& mpreal::operator/=(const long int v)
1494{
1495 mpfr_div_si(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1497 return *this;
1498}
1499
1500inline mpreal& mpreal::operator/=(const int v)
1501{
1502 mpfr_div_si(mpfr_ptr(),mpfr_srcptr(),v,mpreal::get_default_rnd());
1504 return *this;
1505}
1506
1507inline const mpreal operator/(const mpreal& a, const mpreal& b)
1508{
1509 mpreal c(0, (std::max)(mpfr_get_prec(a.mpfr_srcptr()), mpfr_get_prec(b.mpfr_srcptr())));
1510 mpfr_div(c.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), mpreal::get_default_rnd());
1511 return c;
1512}
1513
1514inline const mpreal operator/(const unsigned long int b, const mpreal& a)
1515{
1516 mpreal x(0, mpfr_get_prec(a.mpfr_srcptr()));
1517 mpfr_ui_div(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1518 return x;
1519}
1520
1521inline const mpreal operator/(const unsigned int b, const mpreal& a)
1522{
1523 mpreal x(0, mpfr_get_prec(a.mpfr_srcptr()));
1524 mpfr_ui_div(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1525 return x;
1526}
1527
1528inline const mpreal operator/(const long int b, const mpreal& a)
1529{
1530 mpreal x(0, mpfr_get_prec(a.mpfr_srcptr()));
1531 mpfr_si_div(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1532 return x;
1533}
1534
1535inline const mpreal operator/(const int b, const mpreal& a)
1536{
1537 mpreal x(0, mpfr_get_prec(a.mpfr_srcptr()));
1538 mpfr_si_div(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1539 return x;
1540}
1541
1542inline const mpreal operator/(const double b, const mpreal& a)
1543{
1544#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1545 mpreal x(0, mpfr_get_prec(a.mpfr_srcptr()));
1546 mpfr_d_div(x.mpfr_ptr(), b, a.mpfr_srcptr(), mpreal::get_default_rnd());
1547 return x;
1548#else
1549 mpreal x(b, mpfr_get_prec(a.mpfr_ptr()));
1550 x /= a;
1551 return x;
1552#endif
1553}
1554
1556// Shifts operators - Multiplication/Division by power of 2
1557inline mpreal& mpreal::operator<<=(const unsigned long int u)
1558{
1559 mpfr_mul_2ui(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1561 return *this;
1562}
1563
1564inline mpreal& mpreal::operator<<=(const unsigned int u)
1565{
1566 mpfr_mul_2ui(mpfr_ptr(),mpfr_srcptr(),static_cast<unsigned long int>(u),mpreal::get_default_rnd());
1568 return *this;
1569}
1570
1571inline mpreal& mpreal::operator<<=(const long int u)
1572{
1573 mpfr_mul_2si(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1575 return *this;
1576}
1577
1578inline mpreal& mpreal::operator<<=(const int u)
1579{
1580 mpfr_mul_2si(mpfr_ptr(),mpfr_srcptr(),static_cast<long int>(u),mpreal::get_default_rnd());
1582 return *this;
1583}
1584
1585inline mpreal& mpreal::operator>>=(const unsigned long int u)
1586{
1587 mpfr_div_2ui(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1589 return *this;
1590}
1591
1592inline mpreal& mpreal::operator>>=(const unsigned int u)
1593{
1594 mpfr_div_2ui(mpfr_ptr(),mpfr_srcptr(),static_cast<unsigned long int>(u),mpreal::get_default_rnd());
1596 return *this;
1597}
1598
1599inline mpreal& mpreal::operator>>=(const long int u)
1600{
1601 mpfr_div_2si(mpfr_ptr(),mpfr_srcptr(),u,mpreal::get_default_rnd());
1603 return *this;
1604}
1605
1606inline mpreal& mpreal::operator>>=(const int u)
1607{
1608 mpfr_div_2si(mpfr_ptr(),mpfr_srcptr(),static_cast<long int>(u),mpreal::get_default_rnd());
1610 return *this;
1611}
1612
1613inline const mpreal operator<<(const mpreal& v, const unsigned long int k)
1614{
1615 return mul_2ui(v,k);
1616}
1617
1618inline const mpreal operator<<(const mpreal& v, const unsigned int k)
1619{
1620 return mul_2ui(v,static_cast<unsigned long int>(k));
1621}
1622
1623inline const mpreal operator<<(const mpreal& v, const long int k)
1624{
1625 return mul_2si(v,k);
1626}
1627
1628inline const mpreal operator<<(const mpreal& v, const int k)
1629{
1630 return mul_2si(v,static_cast<long int>(k));
1631}
1632
1633inline const mpreal operator>>(const mpreal& v, const unsigned long int k)
1634{
1635 return div_2ui(v,k);
1636}
1637
1638inline const mpreal operator>>(const mpreal& v, const long int k)
1639{
1640 return div_2si(v,k);
1641}
1642
1643inline const mpreal operator>>(const mpreal& v, const unsigned int k)
1644{
1645 return div_2ui(v,static_cast<unsigned long int>(k));
1646}
1647
1648inline const mpreal operator>>(const mpreal& v, const int k)
1649{
1650 return div_2si(v,static_cast<long int>(k));
1651}
1652
1653// mul_2ui
1654inline const mpreal mul_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode)
1655{
1656 mpreal x(v);
1657 mpfr_mul_2ui(x.mpfr_ptr(),v.mpfr_srcptr(),k,rnd_mode);
1658 return x;
1659}
1660
1661// mul_2si
1662inline const mpreal mul_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode)
1663{
1664 mpreal x(v);
1665 mpfr_mul_2si(x.mpfr_ptr(),v.mpfr_srcptr(),k,rnd_mode);
1666 return x;
1667}
1668
1669inline const mpreal div_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode)
1670{
1671 mpreal x(v);
1672 mpfr_div_2ui(x.mpfr_ptr(),v.mpfr_srcptr(),k,rnd_mode);
1673 return x;
1674}
1675
1676inline const mpreal div_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode)
1677{
1678 mpreal x(v);
1679 mpfr_div_2si(x.mpfr_ptr(),v.mpfr_srcptr(),k,rnd_mode);
1680 return x;
1681}
1682
1684//Relational operators
1685
1686// WARNING:
1687//
1688// Please note that following checks for double-NaN are guaranteed to work only in IEEE math mode:
1689//
1690// isnan(b) = (b != b)
1691// isnan(b) = !(b == b) (we use in code below)
1692//
1693// Be cautions if you use compiler options which break strict IEEE compliance (e.g. -ffast-math in GCC).
1694// Use std::isnan instead (C++11).
1695
1696inline bool operator > (const mpreal& a, const mpreal& b ){ return (mpfr_greater_p(a.mpfr_srcptr(),b.mpfr_srcptr()) != 0 ); }
1697inline bool operator > (const mpreal& a, const unsigned long int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) > 0 ); }
1698inline bool operator > (const mpreal& a, const unsigned int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) > 0 ); }
1699inline bool operator > (const mpreal& a, const long int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) > 0 ); }
1700inline bool operator > (const mpreal& a, const int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) > 0 ); }
1701inline bool operator > (const mpreal& a, const long double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_ld(a.mpfr_srcptr(),b) > 0 ); }
1702inline bool operator > (const mpreal& a, const double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_d (a.mpfr_srcptr(),b) > 0 ); }
1703
1704inline bool operator >= (const mpreal& a, const mpreal& b ){ return (mpfr_greaterequal_p(a.mpfr_srcptr(),b.mpfr_srcptr()) != 0 ); }
1705inline bool operator >= (const mpreal& a, const unsigned long int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) >= 0 ); }
1706inline bool operator >= (const mpreal& a, const unsigned int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) >= 0 ); }
1707inline bool operator >= (const mpreal& a, const long int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) >= 0 ); }
1708inline bool operator >= (const mpreal& a, const int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) >= 0 ); }
1709inline bool operator >= (const mpreal& a, const long double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_ld(a.mpfr_srcptr(),b) >= 0 ); }
1710inline bool operator >= (const mpreal& a, const double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_d (a.mpfr_srcptr(),b) >= 0 ); }
1711
1712inline bool operator < (const mpreal& a, const mpreal& b ){ return (mpfr_less_p(a.mpfr_srcptr(),b.mpfr_srcptr()) != 0 ); }
1713inline bool operator < (const mpreal& a, const unsigned long int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) < 0 ); }
1714inline bool operator < (const mpreal& a, const unsigned int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) < 0 ); }
1715inline bool operator < (const mpreal& a, const long int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) < 0 ); }
1716inline bool operator < (const mpreal& a, const int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) < 0 ); }
1717inline bool operator < (const mpreal& a, const long double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_ld(a.mpfr_srcptr(),b) < 0 ); }
1718inline bool operator < (const mpreal& a, const double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_d (a.mpfr_srcptr(),b) < 0 ); }
1719
1720inline bool operator <= (const mpreal& a, const mpreal& b ){ return (mpfr_lessequal_p(a.mpfr_srcptr(),b.mpfr_srcptr()) != 0 ); }
1721inline bool operator <= (const mpreal& a, const unsigned long int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) <= 0 ); }
1722inline bool operator <= (const mpreal& a, const unsigned int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) <= 0 ); }
1723inline bool operator <= (const mpreal& a, const long int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) <= 0 ); }
1724inline bool operator <= (const mpreal& a, const int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) <= 0 ); }
1725inline bool operator <= (const mpreal& a, const long double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_ld(a.mpfr_srcptr(),b) <= 0 ); }
1726inline bool operator <= (const mpreal& a, const double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_d (a.mpfr_srcptr(),b) <= 0 ); }
1727
1728inline bool operator == (const mpreal& a, const mpreal& b ){ return (mpfr_equal_p(a.mpfr_srcptr(),b.mpfr_srcptr()) != 0 ); }
1729inline bool operator == (const mpreal& a, const unsigned long int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) == 0 ); }
1730inline bool operator == (const mpreal& a, const unsigned int b ){ return !isnan(a) && (mpfr_cmp_ui(a.mpfr_srcptr(),b) == 0 ); }
1731inline bool operator == (const mpreal& a, const long int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) == 0 ); }
1732inline bool operator == (const mpreal& a, const int b ){ return !isnan(a) && (mpfr_cmp_si(a.mpfr_srcptr(),b) == 0 ); }
1733inline bool operator == (const mpreal& a, const long double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_ld(a.mpfr_srcptr(),b) == 0 ); }
1734inline bool operator == (const mpreal& a, const double b ){ return !isnan(a) && (b == b) && (mpfr_cmp_d (a.mpfr_srcptr(),b) == 0 ); }
1735
1736inline bool operator != (const mpreal& a, const mpreal& b ){ return !(a == b); }
1737inline bool operator != (const mpreal& a, const unsigned long int b ){ return !(a == b); }
1738inline bool operator != (const mpreal& a, const unsigned int b ){ return !(a == b); }
1739inline bool operator != (const mpreal& a, const long int b ){ return !(a == b); }
1740inline bool operator != (const mpreal& a, const int b ){ return !(a == b); }
1741inline bool operator != (const mpreal& a, const long double b ){ return !(a == b); }
1742inline bool operator != (const mpreal& a, const double b ){ return !(a == b); }
1743
1744inline bool isnan (const mpreal& op){ return (mpfr_nan_p (op.mpfr_srcptr()) != 0 ); }
1745inline bool isinf (const mpreal& op){ return (mpfr_inf_p (op.mpfr_srcptr()) != 0 ); }
1746inline bool isfinite (const mpreal& op){ return (mpfr_number_p (op.mpfr_srcptr()) != 0 ); }
1747inline bool iszero (const mpreal& op){ return (mpfr_zero_p (op.mpfr_srcptr()) != 0 ); }
1748inline bool isint (const mpreal& op){ return (mpfr_integer_p(op.mpfr_srcptr()) != 0 ); }
1749
1750#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,0,0))
1751inline bool isregular(const mpreal& op){ return (mpfr_regular_p(op.mpfr_srcptr()));}
1752#endif
1753
1755// Type Converters
1756inline bool mpreal::toBool ( ) const { return mpfr_zero_p (mpfr_srcptr()) == 0; }
1757inline long mpreal::toLong (mp_rnd_t mode) const { return mpfr_get_si (mpfr_srcptr(), mode); }
1758inline unsigned long mpreal::toULong (mp_rnd_t mode) const { return mpfr_get_ui (mpfr_srcptr(), mode); }
1759inline float mpreal::toFloat (mp_rnd_t mode) const { return mpfr_get_flt(mpfr_srcptr(), mode); }
1760inline double mpreal::toDouble (mp_rnd_t mode) const { return mpfr_get_d (mpfr_srcptr(), mode); }
1761inline long double mpreal::toLDouble(mp_rnd_t mode) const { return mpfr_get_ld (mpfr_srcptr(), mode); }
1762inline long long mpreal::toLLong (mp_rnd_t mode) const { return mpfr_get_sj (mpfr_srcptr(), mode); }
1763inline unsigned long long mpreal::toULLong (mp_rnd_t mode) const { return mpfr_get_uj (mpfr_srcptr(), mode); }
1764
1765inline ::mpfr_ptr mpreal::mpfr_ptr() { return mp; }
1766inline ::mpfr_srcptr mpreal::mpfr_ptr() const { return mp; }
1767inline ::mpfr_srcptr mpreal::mpfr_srcptr() const { return mp; }
1768
1769template <class T>
1770inline std::string toString(T t, std::ios_base & (*f)(std::ios_base&))
1771{
1772 std::ostringstream oss;
1773 oss << f << t;
1774 return oss.str();
1775}
1776
1777#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1778
1779inline std::string mpreal::toString(const std::string& format) const
1780{
1781 char *s = NULL;
1782 std::string out;
1783
1784 if( !format.empty() )
1785 {
1786 if(!(mpfr_asprintf(&s, format.c_str(), mpfr_srcptr()) < 0))
1787 {
1788 out = std::string(s);
1789
1790 mpfr_free_str(s);
1791 }
1792 }
1793
1794 return out;
1795}
1796
1797#endif
1798
1799inline std::string mpreal::toString(int n, int b, mp_rnd_t mode) const
1800{
1801 // TODO: Add extended format specification (f, e, rounding mode) as it done in output operator
1802 (void)b;
1803 (void)mode;
1804
1805#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1806
1807 std::ostringstream format;
1808
1809 int digits = (n >= 0) ? n : 2 + bits2digits(mpfr_get_prec(mpfr_srcptr()));
1810
1811 format << "%." << digits << "RNg";
1812
1813 return toString(format.str());
1814
1815#else
1816
1817 char *s, *ns = NULL;
1818 size_t slen, nslen;
1819 mp_exp_t exp;
1820 std::string out;
1821
1822 if(mpfr_inf_p(mp))
1823 {
1824 if(mpfr_sgn(mp)>0) return "+Inf";
1825 else return "-Inf";
1826 }
1827
1828 if(mpfr_zero_p(mp)) return "0";
1829 if(mpfr_nan_p(mp)) return "NaN";
1830
1831 s = mpfr_get_str(NULL, &exp, b, 0, mp, mode);
1832 ns = mpfr_get_str(NULL, &exp, b, (std::max)(0,n), mp, mode);
1833
1834 if(s!=NULL && ns!=NULL)
1835 {
1836 slen = strlen(s);
1837 nslen = strlen(ns);
1838 if(nslen<=slen)
1839 {
1840 mpfr_free_str(s);
1841 s = ns;
1842 slen = nslen;
1843 }
1844 else {
1845 mpfr_free_str(ns);
1846 }
1847
1848 // Make human eye-friendly formatting if possible
1849 if (exp>0 && static_cast<size_t>(exp)<slen)
1850 {
1851 if(s[0]=='-')
1852 {
1853 // Remove zeros starting from right end
1854 char* ptr = s+slen-1;
1855 while (*ptr=='0' && ptr>s+exp) ptr--;
1856
1857 if(ptr==s+exp) out = std::string(s,exp+1);
1858 else out = std::string(s,exp+1)+'.'+std::string(s+exp+1,ptr-(s+exp+1)+1);
1859
1860 //out = string(s,exp+1)+'.'+string(s+exp+1);
1861 }
1862 else
1863 {
1864 // Remove zeros starting from right end
1865 char* ptr = s+slen-1;
1866 while (*ptr=='0' && ptr>s+exp-1) ptr--;
1867
1868 if(ptr==s+exp-1) out = std::string(s,exp);
1869 else out = std::string(s,exp)+'.'+std::string(s+exp,ptr-(s+exp)+1);
1870
1871 //out = string(s,exp)+'.'+string(s+exp);
1872 }
1873
1874 }else{ // exp<0 || exp>slen
1875 if(s[0]=='-')
1876 {
1877 // Remove zeros starting from right end
1878 char* ptr = s+slen-1;
1879 while (*ptr=='0' && ptr>s+1) ptr--;
1880
1881 if(ptr==s+1) out = std::string(s,2);
1882 else out = std::string(s,2)+'.'+std::string(s+2,ptr-(s+2)+1);
1883
1884 //out = string(s,2)+'.'+string(s+2);
1885 }
1886 else
1887 {
1888 // Remove zeros starting from right end
1889 char* ptr = s+slen-1;
1890 while (*ptr=='0' && ptr>s) ptr--;
1891
1892 if(ptr==s) out = std::string(s,1);
1893 else out = std::string(s,1)+'.'+std::string(s+1,ptr-(s+1)+1);
1894
1895 //out = string(s,1)+'.'+string(s+1);
1896 }
1897
1898 // Make final string
1899 if(--exp)
1900 {
1901 if(exp>0) out += "e+"+mpfr::toString<mp_exp_t>(exp,std::dec);
1902 else out += "e"+mpfr::toString<mp_exp_t>(exp,std::dec);
1903 }
1904 }
1905
1906 mpfr_free_str(s);
1907 return out;
1908 }else{
1909 return "conversion error!";
1910 }
1911#endif
1912}
1913
1914
1916// I/O
1917inline std::ostream& mpreal::output(std::ostream& os) const
1918{
1919 std::ostringstream format;
1920 const std::ios::fmtflags flags = os.flags();
1921
1922 format << ((flags & std::ios::showpos) ? "%+" : "%");
1923 if (os.precision() >= 0)
1924 format << '.' << os.precision() << "R*"
1925 << ((flags & std::ios::floatfield) == std::ios::fixed ? 'f' :
1926 (flags & std::ios::floatfield) == std::ios::scientific ? 'e' :
1927 'g');
1928 else
1929 format << "R*e";
1930
1931 char *s = NULL;
1932 if(!(mpfr_asprintf(&s, format.str().c_str(),
1934 mpfr_srcptr())
1935 < 0))
1936 {
1937 os << std::string(s);
1938 mpfr_free_str(s);
1939 }
1940 return os;
1941}
1942
1943inline std::ostream& operator<<(std::ostream& os, const mpreal& v)
1944{
1945 return v.output(os);
1946}
1947
1948inline std::istream& operator>>(std::istream &is, mpreal& v)
1949{
1950 // TODO: use cout::hexfloat and other flags to setup base
1951 std::string tmp;
1952 is >> tmp;
1953 mpfr_set_str(v.mpfr_ptr(), tmp.c_str(), 10, mpreal::get_default_rnd());
1954 return is;
1955}
1956
1958// Bits - decimal digits relation
1959// bits = ceil(digits*log[2](10))
1960// digits = floor(bits*log[10](2))
1961
1962inline mp_prec_t digits2bits(int d)
1963{
1964 const double LOG2_10 = 3.3219280948873624;
1965
1966 return mp_prec_t(std::ceil( d * LOG2_10 ));
1967}
1968
1969inline int bits2digits(mp_prec_t b)
1970{
1971 const double LOG10_2 = 0.30102999566398119;
1972
1973 return int(std::floor( b * LOG10_2 ));
1974}
1975
1977// Set/Get number properties
1978inline mpreal& mpreal::setSign(int sign, mp_rnd_t RoundingMode)
1979{
1980 mpfr_setsign(mpfr_ptr(), mpfr_srcptr(), sign < 0, RoundingMode);
1982 return *this;
1983}
1984
1985inline int mpreal::getPrecision() const
1986{
1987 return int(mpfr_get_prec(mpfr_srcptr()));
1988}
1989
1990inline mpreal& mpreal::setPrecision(int Precision, mp_rnd_t RoundingMode)
1991{
1992 mpfr_prec_round(mpfr_ptr(), Precision, RoundingMode);
1994 return *this;
1995}
1996
1997inline mpreal& mpreal::setInf(int sign)
1998{
1999 mpfr_set_inf(mpfr_ptr(), sign);
2001 return *this;
2002}
2003
2005{
2006 mpfr_set_nan(mpfr_ptr());
2008 return *this;
2009}
2010
2011inline mpreal& mpreal::setZero(int sign)
2012{
2013#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,0,0))
2014 mpfr_set_zero(mpfr_ptr(), sign);
2015#else
2016 mpfr_set_si(mpfr_ptr(), 0, (mpfr_get_default_rounding_mode)());
2017 setSign(sign);
2018#endif
2019
2021 return *this;
2022}
2023
2024inline mp_prec_t mpreal::get_prec() const
2025{
2026 return mpfr_get_prec(mpfr_srcptr());
2027}
2028
2029inline void mpreal::set_prec(mp_prec_t prec, mp_rnd_t rnd_mode)
2030{
2031 mpfr_prec_round(mpfr_ptr(),prec,rnd_mode);
2033}
2034
2035inline mp_exp_t mpreal::get_exp () const
2036{
2037 return mpfr_get_exp(mpfr_srcptr());
2038}
2039
2040inline int mpreal::set_exp (mp_exp_t e)
2041{
2042 int x = mpfr_set_exp(mpfr_ptr(), e);
2044 return x;
2045}
2046
2047inline mpreal& negate(mpreal& x) // -x in place
2048{
2049 mpfr_neg(x.mpfr_ptr(),x.mpfr_srcptr(),mpreal::get_default_rnd());
2050 return x;
2051}
2052
2053inline const mpreal frexp(const mpreal& x, mp_exp_t* exp, mp_rnd_t mode = mpreal::get_default_rnd())
2054{
2055 mpreal y(x);
2056#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,1,0))
2057 mpfr_frexp(exp,y.mpfr_ptr(),x.mpfr_srcptr(),mode);
2058#else
2059 *exp = mpfr_get_exp(y.mpfr_srcptr());
2060 mpfr_set_exp(y.mpfr_ptr(),0);
2061#endif
2062 return y;
2063}
2064
2065inline const mpreal frexp(const mpreal& x, int* exp, mp_rnd_t mode = mpreal::get_default_rnd())
2066{
2067 mp_exp_t expl;
2068 mpreal y = frexp(x, &expl, mode);
2069 *exp = int(expl);
2070 return y;
2071}
2072
2073inline const mpreal ldexp(const mpreal& v, mp_exp_t exp)
2074{
2075 mpreal x(v);
2076
2077 // rounding is not important since we are just increasing the exponent (= exact operation)
2078 mpfr_mul_2si(x.mpfr_ptr(), x.mpfr_srcptr(), exp, mpreal::get_default_rnd());
2079 return x;
2080}
2081
2082inline const mpreal scalbn(const mpreal& v, mp_exp_t exp)
2083{
2084 return ldexp(v, exp);
2085}
2086
2087inline mpreal machine_epsilon(mp_prec_t prec)
2088{
2089 /* the smallest eps such that 1 + eps != 1 */
2090 return machine_epsilon(mpreal(1, prec));
2091}
2092
2094{
2095 /* the smallest eps such that x + eps != x */
2096 if( x < 0)
2097 {
2098 return nextabove(-x) + x;
2099 }else{
2100 return nextabove( x) - x;
2101 }
2102}
2103
2104// minval is 'safe' meaning 1 / minval does not overflow
2105inline mpreal minval(mp_prec_t prec)
2106{
2107 /* min = 1/2 * 2^emin = 2^(emin - 1) */
2108 return mpreal(1, prec) << mpreal::get_emin()-1;
2109}
2110
2111// maxval is 'safe' meaning 1 / maxval does not underflow
2112inline mpreal maxval(mp_prec_t prec)
2113{
2114 /* max = (1 - eps) * 2^emax, eps is machine epsilon */
2115 return (mpreal(1, prec) - machine_epsilon(prec)) << mpreal::get_emax();
2116}
2117
2118inline bool isEqualUlps(const mpreal& a, const mpreal& b, int maxUlps)
2119{
2120 return abs(a - b) <= machine_epsilon((max)(abs(a), abs(b))) * maxUlps;
2121}
2122
2123inline bool isEqualFuzzy(const mpreal& a, const mpreal& b, const mpreal& eps)
2124{
2125 return abs(a - b) <= eps;
2126}
2127
2128inline bool isEqualFuzzy(const mpreal& a, const mpreal& b)
2129{
2130 return isEqualFuzzy(a, b, machine_epsilon((max)(1, (min)(abs(a), abs(b)))));
2131}
2132
2134// C++11 sign functions.
2135inline mpreal copysign(const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2136{
2137 mpreal rop(0, mpfr_get_prec(x.mpfr_ptr()));
2138 mpfr_setsign(rop.mpfr_ptr(), x.mpfr_srcptr(), mpfr_signbit(y.mpfr_srcptr()), rnd_mode);
2139 return rop;
2140}
2141
2142inline bool signbit(const mpreal& x)
2143{
2144 return mpfr_signbit(x.mpfr_srcptr());
2145}
2146
2147inline mpreal& setsignbit(mpreal& x, bool minus, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2148{
2149 mpfr_setsign(x.mpfr_ptr(), x.mpfr_srcptr(), minus, rnd_mode);
2150 return x;
2151}
2152
2153inline const mpreal modf(const mpreal& v, mpreal& n)
2154{
2155 mpreal f(v);
2156
2157 // rounding is not important since we are using the same number
2158 mpfr_frac (f.mpfr_ptr(),f.mpfr_srcptr(),mpreal::get_default_rnd());
2159 mpfr_trunc(n.mpfr_ptr(),v.mpfr_srcptr());
2160 return f;
2161}
2162
2163inline int mpreal::check_range (int t, mp_rnd_t rnd_mode)
2164{
2165 return mpfr_check_range(mpfr_ptr(),t,rnd_mode);
2166}
2167
2168inline int mpreal::subnormalize (int t,mp_rnd_t rnd_mode)
2169{
2170 int r = mpfr_subnormalize(mpfr_ptr(),t,rnd_mode);
2172 return r;
2173}
2174
2175inline mp_exp_t mpreal::get_emin (void)
2176{
2177 return mpfr_get_emin();
2178}
2179
2180inline int mpreal::set_emin (mp_exp_t exp)
2181{
2182 return mpfr_set_emin(exp);
2183}
2184
2185inline mp_exp_t mpreal::get_emax (void)
2186{
2187 return mpfr_get_emax();
2188}
2189
2190inline int mpreal::set_emax (mp_exp_t exp)
2191{
2192 return mpfr_set_emax(exp);
2193}
2194
2195inline mp_exp_t mpreal::get_emin_min (void)
2196{
2197 return mpfr_get_emin_min();
2198}
2199
2200inline mp_exp_t mpreal::get_emin_max (void)
2201{
2202 return mpfr_get_emin_max();
2203}
2204
2205inline mp_exp_t mpreal::get_emax_min (void)
2206{
2207 return mpfr_get_emax_min();
2208}
2209
2210inline mp_exp_t mpreal::get_emax_max (void)
2211{
2212 return mpfr_get_emax_max();
2213}
2214
2216// Mathematical Functions
2218
2219// Unary function template with single 'mpreal' argument
2220#define MPREAL_UNARY_MATH_FUNCTION_BODY(f) \
2221 mpreal y(0, mpfr_get_prec(x.mpfr_srcptr())); \
2222 mpfr_##f(y.mpfr_ptr(), x.mpfr_srcptr(), r); \
2223 return y;
2224
2225// Binary function template with 'mpreal' and 'unsigned long' arguments
2226#define MPREAL_BINARY_MATH_FUNCTION_UI_BODY(f, u) \
2227 mpreal y(0, mpfr_get_prec(x.mpfr_srcptr())); \
2228 mpfr_##f(y.mpfr_ptr(), x.mpfr_srcptr(), u, r); \
2229 return y;
2230
2231inline const mpreal sqr (const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd())
2233
2234inline const mpreal sqrt (const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd())
2236
2237inline const mpreal sqrt(const unsigned long int x, mp_rnd_t r)
2238{
2239 mpreal y;
2240 mpfr_sqrt_ui(y.mpfr_ptr(), x, r);
2241 return y;
2242}
2243
2244inline const mpreal sqrt(const unsigned int v, mp_rnd_t rnd_mode)
2245{
2246 return sqrt(static_cast<unsigned long int>(v),rnd_mode);
2247}
2248
2249inline const mpreal sqrt(const long int v, mp_rnd_t rnd_mode)
2250{
2251 if (v>=0) return sqrt(static_cast<unsigned long int>(v),rnd_mode);
2252 else return mpreal().setNan(); // NaN
2253}
2254
2255inline const mpreal sqrt(const int v, mp_rnd_t rnd_mode)
2256{
2257 if (v>=0) return sqrt(static_cast<unsigned long int>(v),rnd_mode);
2258 else return mpreal().setNan(); // NaN
2259}
2260
2261inline const mpreal root(const mpreal& x, unsigned long int k, mp_rnd_t r = mpreal::get_default_rnd())
2262{
2263 mpreal y(0, mpfr_get_prec(x.mpfr_srcptr()));
2264#if (MPFR_VERSION >= MPFR_VERSION_NUM(4,0,0))
2265 mpfr_rootn_ui(y.mpfr_ptr(), x.mpfr_srcptr(), k, r);
2266#else
2267 mpfr_root(y.mpfr_ptr(), x.mpfr_srcptr(), k, r);
2268#endif
2269 return y;
2270}
2271
2272inline const mpreal dim(const mpreal& a, const mpreal& b, mp_rnd_t r = mpreal::get_default_rnd())
2273{
2274 mpreal y(0, mpfr_get_prec(a.mpfr_srcptr()));
2275 mpfr_dim(y.mpfr_ptr(), a.mpfr_srcptr(), b.mpfr_srcptr(), r);
2276 return y;
2277}
2278
2279inline int cmpabs(const mpreal& a,const mpreal& b)
2280{
2281 return mpfr_cmpabs(a.mpfr_ptr(), b.mpfr_srcptr());
2282}
2283
2284inline int sin_cos(mpreal& s, mpreal& c, const mpreal& v, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2285{
2286 return mpfr_sin_cos(s.mpfr_ptr(), c.mpfr_ptr(), v.mpfr_srcptr(), rnd_mode);
2287}
2288
2289inline const mpreal sqrt (const long double v, mp_rnd_t rnd_mode) { return sqrt(mpreal(v),rnd_mode); }
2290inline const mpreal sqrt (const double v, mp_rnd_t rnd_mode) { return sqrt(mpreal(v),rnd_mode); }
2291
2294inline const mpreal abs (const mpreal& x, mp_rnd_t r) { MPREAL_UNARY_MATH_FUNCTION_BODY(abs ); }
2310
2311inline const mpreal logb (const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd()) { return log2 (abs(x),r); }
2312inline mp_exp_t ilogb (const mpreal& x) { return x.get_exp(); }
2313
2314inline const mpreal acot (const mpreal& v, mp_rnd_t r = mpreal::get_default_rnd()) { return atan (1/v, r); }
2315inline const mpreal asec (const mpreal& v, mp_rnd_t r = mpreal::get_default_rnd()) { return acos (1/v, r); }
2316inline const mpreal acsc (const mpreal& v, mp_rnd_t r = mpreal::get_default_rnd()) { return asin (1/v, r); }
2317inline const mpreal acoth (const mpreal& v, mp_rnd_t r = mpreal::get_default_rnd()) { return atanh(1/v, r); }
2318inline const mpreal asech (const mpreal& v, mp_rnd_t r = mpreal::get_default_rnd()) { return acosh(1/v, r); }
2319inline const mpreal acsch (const mpreal& v, mp_rnd_t r = mpreal::get_default_rnd()) { return asinh(1/v, r); }
2320
2330
2344
2345#if (MPFR_VERSION >= MPFR_VERSION_NUM(4,0,0))
2346inline const mpreal gammainc (const mpreal& a, const mpreal& x, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2347{
2348 /*
2349 The non-normalized (upper) incomplete gamma function of a and x:
2350 gammainc(a,x) := Gamma(a,x) = int(t^(a-1) * exp(-t), t=x..infinity)
2351 */
2352 mpreal y(0,(std::max)(a.getPrecision(), x.getPrecision()));
2353 mpfr_gamma_inc(y.mpfr_ptr(), a.mpfr_srcptr(), x.mpfr_srcptr(), rnd_mode);
2354 return y;
2355}
2356
2357inline const mpreal beta (const mpreal& z, const mpreal& w, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2358{
2359 /*
2360 Beta function, uses formula (6.2.2) from Abramowitz & Stegun:
2361 beta(z,w) = gamma(z)*gamma(w)/gamma(z+w)
2362 */
2363 mpreal y(0,(std::max)(z.getPrecision(), w.getPrecision()));
2364 mpfr_beta(y.mpfr_ptr(), z.mpfr_srcptr(), w.mpfr_srcptr(), rnd_mode);
2365 return y;
2366}
2367
2368inline const mpreal log_ui (unsigned long int n, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2369{
2370 /* Computes natural logarithm of an unsigned long */
2371 mpreal y(0, prec);
2372 mpfr_log_ui(y.mpfr_ptr(),n,rnd_mode);
2373 return y;
2374}
2375#endif
2376
2377#if (MPFR_VERSION >= MPFR_VERSION_NUM(4,2,0))
2378
2379/* f(x,u) = f(2*pi*x/u) */
2380inline const mpreal cosu (const mpreal& x, unsigned long u, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_BINARY_MATH_FUNCTION_UI_BODY(cosu, u); }
2381inline const mpreal sinu (const mpreal& x, unsigned long u, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_BINARY_MATH_FUNCTION_UI_BODY(sinu, u); }
2382inline const mpreal tanu (const mpreal& x, unsigned long u, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_BINARY_MATH_FUNCTION_UI_BODY(tanu, u); }
2383inline const mpreal acosu (const mpreal& x, unsigned long u, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_BINARY_MATH_FUNCTION_UI_BODY(acosu, u); }
2384inline const mpreal asinu (const mpreal& x, unsigned long u, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_BINARY_MATH_FUNCTION_UI_BODY(asinu, u); }
2385inline const mpreal atanu (const mpreal& x, unsigned long u, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_BINARY_MATH_FUNCTION_UI_BODY(atanu, u); }
2386
2387/* f(x) = f(pi*x) */
2394
2395inline const mpreal log2p1 (const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_UNARY_MATH_FUNCTION_BODY(log2p1 ); } /* log2 (1+x) */
2396inline const mpreal log10p1(const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_UNARY_MATH_FUNCTION_BODY(log10p1); } /* log10(1+x) */
2397inline const mpreal exp2m1 (const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_UNARY_MATH_FUNCTION_BODY(exp2m1 ); } /* 2^x-1 */
2398inline const mpreal exp10m1(const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd()) { MPREAL_UNARY_MATH_FUNCTION_BODY(exp10m1); } /* 10^x-1 */
2399
2400inline const mpreal atan2u(const mpreal& y, const mpreal& x, unsigned long u, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2401{
2402 /*
2403 atan2u(y,x,u) = atan(|y/x|)*u/(2*pi) for x > 0
2404 atan2u(y,x,u) = 1-atan(|y/x|)*u/(2*pi) for x < 0
2405 */
2406 mpreal a(0, (std::max)(y.getPrecision(), x.getPrecision()));
2407 mpfr_atan2u(a.mpfr_ptr(), y.mpfr_srcptr(), x.mpfr_srcptr(), u, rnd_mode);
2408 return a;
2409}
2410
2411inline const mpreal atan2pi(const mpreal& y, const mpreal& x, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2412{
2413 /* atan2pi(x) = atan2u(u=2) */
2414 mpreal a(0, (std::max)(y.getPrecision(), x.getPrecision()));
2415 mpfr_atan2pi(a.mpfr_ptr(), y.mpfr_srcptr(), x.mpfr_srcptr(), rnd_mode);
2416 return a;
2417}
2418
2419inline const mpreal powr(const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2420{
2421 /* powr(x,y) = exp(y*log(x)) */
2422 mpreal a(0, (std::max)(x.getPrecision(), y.getPrecision()));
2423 mpfr_powr(a.mpfr_ptr(), x.mpfr_srcptr(), y.mpfr_srcptr(), rnd_mode);
2424 return a;
2425}
2426
2427inline const mpreal compound(const mpreal& x, long n, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2428{
2429 /* compound(x,n) = (1+x)^n */
2430 mpreal y(0, x.getPrecision());
2431 mpfr_compound_si(y.mpfr_ptr(),x.mpfr_srcptr(),n,rnd_mode);
2432 return y;
2433}
2434
2435inline const mpreal fmod(const mpreal& x, unsigned long u, mp_rnd_t r = mpreal::get_default_rnd())
2436{
2437 /* x modulo a machine integer u */
2439}
2440#endif
2441
2442inline const mpreal nextpow2(const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd())
2443{
2444 mpreal y(0, x.getPrecision());
2445
2446 if(!iszero(x))
2447 y = ceil(log2(abs(x,r),r));
2448
2449 return y;
2450}
2451
2452inline const mpreal atan2 (const mpreal& y, const mpreal& x, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2453{
2454 mpreal a(0,(std::max)(y.getPrecision(), x.getPrecision()));
2455 mpfr_atan2(a.mpfr_ptr(), y.mpfr_srcptr(), x.mpfr_srcptr(), rnd_mode);
2456 return a;
2457}
2458
2459inline const mpreal hypot (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2460{
2461 mpreal a(0,(std::max)(y.getPrecision(), x.getPrecision()));
2462 mpfr_hypot(a.mpfr_ptr(), x.mpfr_srcptr(), y.mpfr_srcptr(), rnd_mode);
2463 return a;
2464}
2465
2466inline const mpreal hypot(const mpreal& a, const mpreal& b, const mpreal& c)
2467{
2468 if(isnan(a) || isnan(b) || isnan(c)) return mpreal().setNan();
2469 else
2470 {
2471 mpreal absa = abs(a), absb = abs(b), absc = abs(c);
2472 mpreal w = (std::max)(absa, (std::max)(absb, absc));
2473 mpreal r;
2474
2475 if (!iszero(w))
2476 {
2477 mpreal iw = 1/w;
2478 r = w * sqrt(sqr(absa*iw) + sqr(absb*iw) + sqr(absc*iw));
2479 }
2480
2481 return r;
2482 }
2483}
2484
2485inline const mpreal hypot(const mpreal& a, const mpreal& b, const mpreal& c, const mpreal& d)
2486{
2487 if(isnan(a) || isnan(b) || isnan(c) || isnan(d)) return mpreal().setNan();
2488 else
2489 {
2490 mpreal absa = abs(a), absb = abs(b), absc = abs(c), absd = abs(d);
2491 mpreal w = (std::max)(absa, (std::max)(absb, (std::max)(absc, absd)));
2492 mpreal r;
2493
2494 if (!iszero(w))
2495 {
2496 mpreal iw = 1/w;
2497 r = w * sqrt(sqr(absa*iw) + sqr(absb*iw) + sqr(absc*iw) + sqr(absd*iw));
2498 }
2499
2500 return r;
2501 }
2502}
2503
2504inline const mpreal remainder (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2505{
2506 mpreal a(0,(std::max)(y.getPrecision(), x.getPrecision()));
2507 mpfr_remainder(a.mpfr_ptr(), x.mpfr_srcptr(), y.mpfr_srcptr(), rnd_mode);
2508 return a;
2509}
2510
2511inline const mpreal remquo (const mpreal& x, const mpreal& y, int* q, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2512{
2513 long lq;
2514 mpreal a(0,(std::max)(y.getPrecision(), x.getPrecision()));
2515 mpfr_remquo(a.mpfr_ptr(), &lq, x.mpfr_srcptr(), y.mpfr_srcptr(), rnd_mode);
2516 if (q) *q = int(lq);
2517 return a;
2518}
2519
2520inline const mpreal fac_ui (unsigned long int v, mp_prec_t prec = mpreal::get_default_prec(),
2521 mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2522{
2523 mpreal x(0, prec);
2524 mpfr_fac_ui(x.mpfr_ptr(),v,rnd_mode);
2525 return x;
2526}
2527
2528
2529inline const mpreal lgamma (const mpreal& v, int *signp = 0, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2530{
2531 mpreal x(v);
2532 int tsignp;
2533
2534 if(signp) mpfr_lgamma(x.mpfr_ptr(), signp,v.mpfr_srcptr(),rnd_mode);
2535 else mpfr_lgamma(x.mpfr_ptr(),&tsignp,v.mpfr_srcptr(),rnd_mode);
2536
2537 return x;
2538}
2539
2540
2541inline const mpreal besseljn (long n, const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd())
2542{
2543 mpreal y(0, x.getPrecision());
2544 mpfr_jn(y.mpfr_ptr(), n, x.mpfr_srcptr(), r);
2545 return y;
2546}
2547
2548inline const mpreal besselyn (long n, const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd())
2549{
2550 mpreal y(0, x.getPrecision());
2551 mpfr_yn(y.mpfr_ptr(), n, x.mpfr_srcptr(), r);
2552 return y;
2553}
2554
2555inline const mpreal fma (const mpreal& v1, const mpreal& v2, const mpreal& v3, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2556{
2557 mpreal a;
2558 mp_prec_t p1, p2, p3;
2559
2560 p1 = v1.get_prec();
2561 p2 = v2.get_prec();
2562 p3 = v3.get_prec();
2563
2564 a.set_prec(p3>p2?(p3>p1?p3:p1):(p2>p1?p2:p1));
2565
2566 mpfr_fma(a.mp,v1.mp,v2.mp,v3.mp,rnd_mode);
2567 return a;
2568}
2569
2570inline const mpreal fms (const mpreal& v1, const mpreal& v2, const mpreal& v3, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2571{
2572 mpreal a;
2573 mp_prec_t p1, p2, p3;
2574
2575 p1 = v1.get_prec();
2576 p2 = v2.get_prec();
2577 p3 = v3.get_prec();
2578
2579 a.set_prec(p3>p2?(p3>p1?p3:p1):(p2>p1?p2:p1));
2580
2581 mpfr_fms(a.mp,v1.mp,v2.mp,v3.mp,rnd_mode);
2582 return a;
2583}
2584
2585inline const mpreal agm (const mpreal& v1, const mpreal& v2, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2586{
2587 mpreal a;
2588 mp_prec_t p1, p2;
2589
2590 p1 = v1.get_prec();
2591 p2 = v2.get_prec();
2592
2593 a.set_prec(p1>p2?p1:p2);
2594
2595 mpfr_agm(a.mp, v1.mp, v2.mp, rnd_mode);
2596
2597 return a;
2598}
2599
2600inline const mpreal sum (const mpreal tab[], const unsigned long int n, int& status, mp_rnd_t mode = mpreal::get_default_rnd())
2601{
2602 mpfr_srcptr *p = new mpfr_srcptr[n];
2603
2604 for (unsigned long int i = 0; i < n; i++)
2605 p[i] = tab[i].mpfr_srcptr();
2606
2607 mpreal x;
2608 status = mpfr_sum(x.mpfr_ptr(), (mpfr_ptr*)p, n, mode);
2609
2610 delete [] p;
2611 return x;
2612}
2613
2615// MPFR 2.4.0 Specifics
2616#if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
2617
2618inline int sinh_cosh(mpreal& s, mpreal& c, const mpreal& v, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2619{
2620 return mpfr_sinh_cosh(s.mp,c.mp,v.mp,rnd_mode);
2621}
2622
2623inline const mpreal li2 (const mpreal& x, mp_rnd_t r = mpreal::get_default_rnd())
2624{
2626}
2627
2628inline const mpreal rem (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2629{
2630 /* R = rem(X,Y) if Y != 0, returns X - n * Y where n = trunc(X/Y). */
2631 return fmod(x, y, rnd_mode);
2632}
2633
2634inline const mpreal mod (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2635{
2636 (void)rnd_mode;
2637
2638 /*
2639
2640 m = mod(x,y) if y != 0, returns x - n*y where n = floor(x/y)
2641
2642 The following are true by convention:
2643 - mod(x,0) is x
2644 - mod(x,x) is 0
2645 - mod(x,y) for x != y and y != 0 has the same sign as y.
2646
2647 */
2648
2649 if(iszero(y)) return x;
2650 if(x == y) return 0;
2651
2652 mpreal m = x - floor(x / y) * y;
2653
2654 return copysign(abs(m),y); // make sure result has the same sign as Y
2655}
2656
2657inline const mpreal fmod (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2658{
2659 mpreal a;
2660 mp_prec_t yp, xp;
2661
2662 yp = y.get_prec();
2663 xp = x.get_prec();
2664
2665 a.set_prec(yp>xp?yp:xp);
2666
2667 mpfr_fmod(a.mp, x.mp, y.mp, rnd_mode);
2668
2669 return a;
2670}
2671
2672inline const mpreal rec_sqrt(const mpreal& v, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2673{
2674 mpreal x(v);
2675 mpfr_rec_sqrt(x.mp,v.mp,rnd_mode);
2676 return x;
2677}
2678#endif // MPFR 2.4.0 Specifics
2679
2681// MPFR 3.0.0 Specifics
2682#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,0,0))
2685#endif // MPFR 3.0.0 Specifics
2686
2688// Constants
2689inline const mpreal const_log2 (mp_prec_t p = mpreal::get_default_prec(), mp_rnd_t r = mpreal::get_default_rnd())
2690{
2691 mpreal x(0, p);
2692 mpfr_const_log2(x.mpfr_ptr(), r);
2693 return x;
2694}
2695
2696inline const mpreal const_pi (mp_prec_t p = mpreal::get_default_prec(), mp_rnd_t r = mpreal::get_default_rnd())
2697{
2698 mpreal x(0, p);
2699 mpfr_const_pi(x.mpfr_ptr(), r);
2700 return x;
2701}
2702
2703inline const mpreal const_euler (mp_prec_t p = mpreal::get_default_prec(), mp_rnd_t r = mpreal::get_default_rnd())
2704{
2705 mpreal x(0, p);
2706 mpfr_const_euler(x.mpfr_ptr(), r);
2707 return x;
2708}
2709
2710inline const mpreal const_catalan (mp_prec_t p = mpreal::get_default_prec(), mp_rnd_t r = mpreal::get_default_rnd())
2711{
2712 mpreal x(0, p);
2713 mpfr_const_catalan(x.mpfr_ptr(), r);
2714 return x;
2715}
2716
2717inline const mpreal const_infinity (int sign = 1, mp_prec_t p = mpreal::get_default_prec())
2718{
2719 mpreal x(0, p);
2720 mpfr_set_inf(x.mpfr_ptr(), sign);
2721 return x;
2722}
2723
2725// Integer Related Functions
2726inline const mpreal ceil(const mpreal& v)
2727{
2728 mpreal x(v);
2729 mpfr_ceil(x.mp,v.mp);
2730 return x;
2731}
2732
2733inline const mpreal floor(const mpreal& v)
2734{
2735 mpreal x(v);
2736 mpfr_floor(x.mp,v.mp);
2737 return x;
2738}
2739
2740inline const mpreal round(const mpreal& v)
2741{
2742 mpreal x(v);
2743 mpfr_round(x.mp,v.mp);
2744 return x;
2745}
2746
2747inline long lround(const mpreal& v)
2748{
2749 long r = std::numeric_limits<long>::min();
2750 mpreal x = round(v);
2751 if (abs(x) < -mpreal(r)) // Assume mpreal(LONG_MIN) is exact
2752 r = x.toLong();
2753 return r;
2754}
2755
2756inline long long llround(const mpreal& v)
2757{
2758 long long r = std::numeric_limits<long long>::min();
2759 mpreal x = round(v);
2760 if (abs(x) < -mpreal(r)) // Assume mpreal(LLONG_MIN) is exact
2761 r = x.toLLong();
2762 return r;
2763}
2764
2765inline const mpreal trunc(const mpreal& v)
2766{
2767 mpreal x(v);
2768 mpfr_trunc(x.mp,v.mp);
2769 return x;
2770}
2771
2778
2780// Miscellaneous Functions
2781inline int sgn(const mpreal& op)
2782{
2783 // Please note, this is classic signum function which ignores sign of zero.
2784 // Use signbit if you need sign of zero.
2785 return mpfr_sgn(op.mpfr_srcptr());
2786}
2787
2789// Miscellaneous Functions
2790inline void swap (mpreal& a, mpreal& b) { mpfr_swap(a.mpfr_ptr(),b.mpfr_ptr()); }
2791inline const mpreal (max)(const mpreal& x, const mpreal& y){ return (x<y?y:x); }
2792inline const mpreal (min)(const mpreal& x, const mpreal& y){ return (y<x?y:x); }
2793
2794inline const mpreal fmax(const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2795{
2796 mpreal a;
2797 mpfr_max(a.mp,x.mp,y.mp,rnd_mode);
2798 return a;
2799}
2800
2801inline const mpreal fmin(const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2802{
2803 mpreal a;
2804 mpfr_min(a.mp,x.mp,y.mp,rnd_mode);
2805 return a;
2806}
2807
2808inline const mpreal nexttoward (const mpreal& x, const mpreal& y)
2809{
2810 mpreal a(x);
2811 mpfr_nexttoward(a.mp,y.mp);
2812 return a;
2813}
2814
2815inline const mpreal nextabove (const mpreal& x)
2816{
2817 mpreal a(x);
2818 mpfr_nextabove(a.mp);
2819 return a;
2820}
2821
2822inline const mpreal nextbelow (const mpreal& x)
2823{
2824 mpreal a(x);
2825 mpfr_nextbelow(a.mp);
2826 return a;
2827}
2828
2829inline const mpreal urandomb (gmp_randstate_t& state)
2830{
2831 mpreal x;
2832 mpfr_urandomb(x.mpfr_ptr(),state);
2833 return x;
2834}
2835
2836#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,0,0))
2837inline const mpreal urandom (gmp_randstate_t& state, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2838{
2839 mpreal x;
2840 mpfr_urandom(x.mpfr_ptr(), state, rnd_mode);
2841 return x;
2842}
2843#endif
2844
2845#if (MPFR_VERSION <= MPFR_VERSION_NUM(2,4,2))
2846inline const mpreal random2 (mp_size_t size, mp_exp_t exp)
2847{
2848 mpreal x;
2849 mpfr_random2(x.mpfr_ptr(),size,exp);
2850 return x;
2851}
2852#endif
2853
2854// Uniformly distributed random number generation
2855// a = random(seed); <- initialization & first random number generation
2856// a = random(); <- next random numbers generation
2857// seed != 0
2858inline const mpreal random(unsigned int seed = 0)
2859{
2860#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,0,0))
2861 static gmp_randstate_t state;
2862 static bool initialize = true;
2863
2864 if(initialize)
2865 {
2866 gmp_randinit_default(state);
2867 gmp_randseed_ui(state,0);
2868 initialize = false;
2869 }
2870
2871 if(seed != 0) gmp_randseed_ui(state,seed);
2872
2873 return mpfr::urandom(state);
2874#else
2875 if(seed != 0) std::srand(seed);
2876 return mpfr::mpreal(std::rand()/(double)RAND_MAX);
2877#endif
2878}
2879
2880#if (MPFR_VERSION >= MPFR_VERSION_NUM(3,1,0))
2881inline const mpreal grandom (gmp_randstate_t& state, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2882{
2883 mpreal x;
2884#if (MPFR_VERSION >= MPFR_VERSION_NUM(4,0,0))
2885 mpfr_nrandom(x.mpfr_ptr(), state, rnd_mode);
2886#else
2887 mpfr_grandom(x.mpfr_ptr(), NULL, state, rnd_mode);
2888#endif
2889 return x;
2890}
2891
2892inline const mpreal grandom(unsigned int seed = 0)
2893{
2894 static gmp_randstate_t state;
2895 static bool initialize = true;
2896
2897 if(initialize)
2898 {
2899 gmp_randinit_default(state);
2900 gmp_randseed_ui(state,0);
2901 initialize = false;
2902 }
2903
2904 if(seed != 0) gmp_randseed_ui(state,seed);
2905
2906 return mpfr::grandom(state);
2907}
2908#endif
2909
2911// Set/Get global properties
2912inline void mpreal::set_default_prec(mp_prec_t prec)
2913{
2914 mpfr_set_default_prec(prec);
2915}
2916
2917inline void mpreal::set_default_rnd(mp_rnd_t rnd_mode)
2918{
2919 mpfr_set_default_rounding_mode(rnd_mode);
2920}
2921
2922inline bool mpreal::fits_in_bits(double x, int n)
2923{
2924 int i;
2925 double t;
2926 return IsInf(x) || (std::modf ( std::ldexp ( std::frexp ( x, &i ), n ), &t ) == 0.0);
2927}
2928
2929inline const mpreal pow(const mpreal& a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2930{
2931 mpreal x(a);
2932 mpfr_pow(x.mp,x.mp,b.mp,rnd_mode);
2933 return x;
2934}
2935
2936inline const mpreal pow(const mpreal& a, const mpz_t b, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2937{
2938 mpreal x(a);
2939 mpfr_pow_z(x.mp,x.mp,b,rnd_mode);
2940 return x;
2941}
2942
2943inline const mpreal pow(const mpreal& a, const long long b, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2944{
2945 (void)rnd_mode;
2946 return pow(a,mpreal(b));
2947}
2948
2949inline const mpreal pow(const mpreal& a, const unsigned long long b, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2950{
2951 (void)rnd_mode;
2952 return pow(a,mpreal(b));
2953}
2954
2955inline const mpreal pow(const mpreal& a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2956{
2957 mpreal x(a);
2958 mpfr_pow_ui(x.mp,x.mp,b,rnd_mode);
2959 return x;
2960}
2961
2962inline const mpreal pow(const mpreal& a, const unsigned int b, mp_rnd_t rnd_mode)
2963{
2964 return pow(a,static_cast<unsigned long int>(b),rnd_mode);
2965}
2966
2967inline const mpreal pow(const mpreal& a, const long int b, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2968{
2969 mpreal x(a);
2970 mpfr_pow_si(x.mp,x.mp,b,rnd_mode);
2971 return x;
2972}
2973
2974inline const mpreal pow(const mpreal& a, const int b, mp_rnd_t rnd_mode)
2975{
2976 return pow(a,static_cast<long int>(b),rnd_mode);
2977}
2978
2979inline const mpreal pow(const mpreal& a, const long double b, mp_rnd_t rnd_mode)
2980{
2981 return pow(a,mpreal(b),rnd_mode);
2982}
2983
2984inline const mpreal pow(const mpreal& a, const double b, mp_rnd_t rnd_mode)
2985{
2986 return pow(a,mpreal(b),rnd_mode);
2987}
2988
2989inline const mpreal pow(const unsigned long int a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::get_default_rnd())
2990{
2991 mpreal x(a);
2992 mpfr_ui_pow(x.mp,a,b.mp,rnd_mode);
2993 return x;
2994}
2995
2996inline const mpreal pow(const unsigned int a, const mpreal& b, mp_rnd_t rnd_mode)
2997{
2998 return pow(static_cast<unsigned long int>(a),b,rnd_mode);
2999}
3000
3001inline const mpreal pow(const long int a, const mpreal& b, mp_rnd_t rnd_mode)
3002{
3003 if (a>=0) return pow(static_cast<unsigned long int>(a),b,rnd_mode);
3004 else return pow(mpreal(a),b,rnd_mode);
3005}
3006
3007inline const mpreal pow(const int a, const mpreal& b, mp_rnd_t rnd_mode)
3008{
3009 if (a>=0) return pow(static_cast<unsigned long int>(a),b,rnd_mode);
3010 else return pow(mpreal(a),b,rnd_mode);
3011}
3012
3013inline const mpreal pow(const long double a, const mpreal& b, mp_rnd_t rnd_mode)
3014{
3015 return pow(mpreal(a),b,rnd_mode);
3016}
3017
3018inline const mpreal pow(const double a, const mpreal& b, mp_rnd_t rnd_mode)
3019{
3020 return pow(mpreal(a),b,rnd_mode);
3021}
3022
3023// pow unsigned long int
3024inline const mpreal pow(const unsigned long int a, const unsigned long int b, mp_rnd_t rnd_mode)
3025{
3026 mpreal x(a);
3027 mpfr_ui_pow_ui(x.mp,a,b,rnd_mode);
3028 return x;
3029}
3030
3031inline const mpreal pow(const unsigned long int a, const unsigned int b, mp_rnd_t rnd_mode)
3032{
3033 return pow(a,static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3034}
3035
3036inline const mpreal pow(const unsigned long int a, const long int b, mp_rnd_t rnd_mode)
3037{
3038 if(b>0) return pow(a,static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3039 else return pow(a,mpreal(b),rnd_mode); //mpfr_ui_pow
3040}
3041
3042inline const mpreal pow(const unsigned long int a, const int b, mp_rnd_t rnd_mode)
3043{
3044 if(b>0) return pow(a,static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3045 else return pow(a,mpreal(b),rnd_mode); //mpfr_ui_pow
3046}
3047
3048inline const mpreal pow(const unsigned long int a, const long double b, mp_rnd_t rnd_mode)
3049{
3050 return pow(a,mpreal(b),rnd_mode); //mpfr_ui_pow
3051}
3052
3053inline const mpreal pow(const unsigned long int a, const double b, mp_rnd_t rnd_mode)
3054{
3055 return pow(a,mpreal(b),rnd_mode); //mpfr_ui_pow
3056}
3057
3058// pow unsigned int
3059inline const mpreal pow(const unsigned int a, const unsigned long int b, mp_rnd_t rnd_mode)
3060{
3061 return pow(static_cast<unsigned long int>(a),b,rnd_mode); //mpfr_ui_pow_ui
3062}
3063
3064inline const mpreal pow(const unsigned int a, const unsigned int b, mp_rnd_t rnd_mode)
3065{
3066 return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3067}
3068
3069inline const mpreal pow(const unsigned int a, const long int b, mp_rnd_t rnd_mode)
3070{
3071 if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3072 else return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3073}
3074
3075inline const mpreal pow(const unsigned int a, const int b, mp_rnd_t rnd_mode)
3076{
3077 if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3078 else return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3079}
3080
3081inline const mpreal pow(const unsigned int a, const long double b, mp_rnd_t rnd_mode)
3082{
3083 return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3084}
3085
3086inline const mpreal pow(const unsigned int a, const double b, mp_rnd_t rnd_mode)
3087{
3088 return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3089}
3090
3091// pow long int
3092inline const mpreal pow(const long int a, const unsigned long int b, mp_rnd_t rnd_mode)
3093{
3094 if (a>0) return pow(static_cast<unsigned long int>(a),b,rnd_mode); //mpfr_ui_pow_ui
3095 else return pow(mpreal(a),b,rnd_mode); //mpfr_pow_ui
3096}
3097
3098inline const mpreal pow(const long int a, const unsigned int b, mp_rnd_t rnd_mode)
3099{
3100 if (a>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3101 else return pow(mpreal(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_pow_ui
3102}
3103
3104inline const mpreal pow(const long int a, const long int b, mp_rnd_t rnd_mode)
3105{
3106 if (a>0)
3107 {
3108 if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3109 else return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3110 }else{
3111 return pow(mpreal(a),b,rnd_mode); // mpfr_pow_si
3112 }
3113}
3114
3115inline const mpreal pow(const long int a, const int b, mp_rnd_t rnd_mode)
3116{
3117 if (a>0)
3118 {
3119 if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3120 else return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3121 }else{
3122 return pow(mpreal(a),static_cast<long int>(b),rnd_mode); // mpfr_pow_si
3123 }
3124}
3125
3126inline const mpreal pow(const long int a, const long double b, mp_rnd_t rnd_mode)
3127{
3128 if (a>=0) return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3129 else return pow(mpreal(a),mpreal(b),rnd_mode); //mpfr_pow
3130}
3131
3132inline const mpreal pow(const long int a, const double b, mp_rnd_t rnd_mode)
3133{
3134 if (a>=0) return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3135 else return pow(mpreal(a),mpreal(b),rnd_mode); //mpfr_pow
3136}
3137
3138// pow int
3139inline const mpreal pow(const int a, const unsigned long int b, mp_rnd_t rnd_mode)
3140{
3141 if (a>0) return pow(static_cast<unsigned long int>(a),b,rnd_mode); //mpfr_ui_pow_ui
3142 else return pow(mpreal(a),b,rnd_mode); //mpfr_pow_ui
3143}
3144
3145inline const mpreal pow(const int a, const unsigned int b, mp_rnd_t rnd_mode)
3146{
3147 if (a>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3148 else return pow(mpreal(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_pow_ui
3149}
3150
3151inline const mpreal pow(const int a, const long int b, mp_rnd_t rnd_mode)
3152{
3153 if (a>0)
3154 {
3155 if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3156 else return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3157 }else{
3158 return pow(mpreal(a),b,rnd_mode); // mpfr_pow_si
3159 }
3160}
3161
3162inline const mpreal pow(const int a, const int b, mp_rnd_t rnd_mode)
3163{
3164 if (a>0)
3165 {
3166 if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
3167 else return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3168 }else{
3169 return pow(mpreal(a),static_cast<long int>(b),rnd_mode); // mpfr_pow_si
3170 }
3171}
3172
3173inline const mpreal pow(const int a, const long double b, mp_rnd_t rnd_mode)
3174{
3175 if (a>=0) return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3176 else return pow(mpreal(a),mpreal(b),rnd_mode); //mpfr_pow
3177}
3178
3179inline const mpreal pow(const int a, const double b, mp_rnd_t rnd_mode)
3180{
3181 if (a>=0) return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3182 else return pow(mpreal(a),mpreal(b),rnd_mode); //mpfr_pow
3183}
3184
3185// pow long double
3186inline const mpreal pow(const long double a, const long double b, mp_rnd_t rnd_mode)
3187{
3188 return pow(mpreal(a),mpreal(b),rnd_mode);
3189}
3190
3191inline const mpreal pow(const long double a, const unsigned long int b, mp_rnd_t rnd_mode)
3192{
3193 return pow(mpreal(a),b,rnd_mode); //mpfr_pow_ui
3194}
3195
3196inline const mpreal pow(const long double a, const unsigned int b, mp_rnd_t rnd_mode)
3197{
3198 return pow(mpreal(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_pow_ui
3199}
3200
3201inline const mpreal pow(const long double a, const long int b, mp_rnd_t rnd_mode)
3202{
3203 return pow(mpreal(a),b,rnd_mode); // mpfr_pow_si
3204}
3205
3206inline const mpreal pow(const long double a, const int b, mp_rnd_t rnd_mode)
3207{
3208 return pow(mpreal(a),static_cast<long int>(b),rnd_mode); // mpfr_pow_si
3209}
3210
3211inline const mpreal pow(const double a, const double b, mp_rnd_t rnd_mode)
3212{
3213 return pow(mpreal(a),mpreal(b),rnd_mode);
3214}
3215
3216inline const mpreal pow(const double a, const unsigned long int b, mp_rnd_t rnd_mode)
3217{
3218 return pow(mpreal(a),b,rnd_mode); // mpfr_pow_ui
3219}
3220
3221inline const mpreal pow(const double a, const unsigned int b, mp_rnd_t rnd_mode)
3222{
3223 return pow(mpreal(a),static_cast<unsigned long int>(b),rnd_mode); // mpfr_pow_ui
3224}
3225
3226inline const mpreal pow(const double a, const long int b, mp_rnd_t rnd_mode)
3227{
3228 return pow(mpreal(a),b,rnd_mode); // mpfr_pow_si
3229}
3230
3231inline const mpreal pow(const double a, const int b, mp_rnd_t rnd_mode)
3232{
3233 return pow(mpreal(a),static_cast<long int>(b),rnd_mode); // mpfr_pow_si
3234}
3235} // End of mpfr namespace
3236
3237// Explicit specialization of std::swap for mpreal numbers
3238// Thus standard algorithms will use efficient version of swap (due to Koenig lookup)
3239// Non-throwing swap C++ idiom: http://en.wikibooks.org/wiki/More_C%2B%2B_Idioms/Non-throwing_swap
3240namespace std
3241{
3242
3243 template <>
3244 inline void swap(mpfr::mpreal& x, mpfr::mpreal& y)
3245 {
3246 return mpfr::swap(x, y);
3247 }
3248
3249 template<>
3250 class numeric_limits<mpfr::mpreal>
3251 {
3252 public:
3253 static const bool is_specialized = true;
3254 static const bool is_signed = true;
3255 static const bool is_integer = false;
3256 static const bool is_exact = false;
3257 static const int radix = 2;
3258
3259 static const bool has_infinity = true;
3260 static const bool has_quiet_NaN = true;
3261 static const bool has_signaling_NaN = true;
3262
3263 static const bool is_iec559 = true; // = IEEE 754
3264 static const bool is_bounded = true;
3265 static const bool is_modulo = false;
3266 static const bool traps = true;
3267 static const bool tinyness_before = true;
3268
3269 static const float_denorm_style has_denorm = denorm_absent;
3270
3271 inline static mpfr::mpreal (min) (mp_prec_t precision = mpfr::mpreal::get_default_prec()) { return mpfr::minval(precision); }
3272 inline static mpfr::mpreal (max) (mp_prec_t precision = mpfr::mpreal::get_default_prec()) { return mpfr::maxval(precision); }
3273 inline static mpfr::mpreal lowest (mp_prec_t precision = mpfr::mpreal::get_default_prec()) { return -mpfr::maxval(precision); }
3274
3275 // Returns smallest eps such that 1 + eps != 1 (classic machine epsilon)
3276 inline static mpfr::mpreal epsilon(mp_prec_t precision = mpfr::mpreal::get_default_prec()) { return mpfr::machine_epsilon(precision); }
3277
3278 // Returns smallest eps such that x + eps != x (relative machine epsilon)
3279 inline static mpfr::mpreal epsilon(const mpfr::mpreal& x) { return mpfr::machine_epsilon(x); }
3280
3281 inline static mpfr::mpreal round_error(mp_prec_t precision = mpfr::mpreal::get_default_prec())
3282 {
3283 mp_rnd_t r = mpfr::mpreal::get_default_rnd();
3284
3285 if(r == GMP_RNDN) return mpfr::mpreal(0.5, precision);
3286 else return mpfr::mpreal(1.0, precision);
3287 }
3288
3289 inline static const mpfr::mpreal infinity() { return mpfr::const_infinity(); }
3290 inline static const mpfr::mpreal quiet_NaN() { return mpfr::mpreal().setNan(); }
3291 inline static const mpfr::mpreal signaling_NaN() { return mpfr::mpreal().setNan(); }
3292 inline static const mpfr::mpreal denorm_min() { return (min)(); }
3293
3294 // Please note, exponent range is not fixed in MPFR
3295 static const int min_exponent = MPFR_EMIN_DEFAULT;
3296 static const int max_exponent = MPFR_EMAX_DEFAULT;
3297 MPREAL_PERMISSIVE_EXPR static const int min_exponent10 = (int) (MPFR_EMIN_DEFAULT * 0.3010299956639811);
3298 MPREAL_PERMISSIVE_EXPR static const int max_exponent10 = (int) (MPFR_EMAX_DEFAULT * 0.3010299956639811);
3299
3300#ifdef MPREAL_HAVE_DYNAMIC_STD_NUMERIC_LIMITS
3301
3302 // Following members should be constant according to standard, but they can be variable in MPFR
3303 // So we define them as functions here.
3304 //
3305 // This is preferable way for std::numeric_limits<mpfr::mpreal> specialization.
3306 // But it is incompatible with standard std::numeric_limits and might not work with other libraries, e.g. boost.
3307 // See below for compatible implementation.
3308 inline static float_round_style round_style()
3309 {
3310 mp_rnd_t r = mpfr::mpreal::get_default_rnd();
3311
3312 switch (r)
3313 {
3314 case GMP_RNDN: return round_to_nearest;
3315 case GMP_RNDZ: return round_toward_zero;
3316 case GMP_RNDU: return round_toward_infinity;
3317 case GMP_RNDD: return round_toward_neg_infinity;
3318 default: return round_indeterminate;
3319 }
3320 }
3321
3322 inline static int digits() { return int(mpfr::mpreal::get_default_prec()); }
3323 inline static int digits(const mpfr::mpreal& x) { return x.getPrecision(); }
3324
3325 inline static int digits10(mp_prec_t precision = mpfr::mpreal::get_default_prec())
3326 {
3327 return mpfr::bits2digits(precision);
3328 }
3329
3330 inline static int digits10(const mpfr::mpreal& x)
3331 {
3332 return mpfr::bits2digits(x.getPrecision());
3333 }
3334
3335 inline static int max_digits10(mp_prec_t precision = mpfr::mpreal::get_default_prec())
3336 {
3337 return digits10(precision);
3338 }
3339#else
3340 // Digits and round_style are NOT constants when it comes to mpreal.
3341 // If possible, please use functions digits() and round_style() defined above.
3342 //
3343 // These (default) values are preserved for compatibility with existing libraries, e.g. boost.
3344 // Change them accordingly to your application.
3345 //
3346 // For example, if you use 256 bits of precision uniformly in your program, then:
3347 // digits = 256
3348 // digits10 = 77
3349 // max_digits10 = 78
3350 //
3351 // Approximate formula for decimal digits is: digits10 = floor(log10(2) * digits). See bits2digits() for more details.
3352
3353 static const std::float_round_style round_style = round_to_nearest;
3354 static const int digits = 53;
3355 static const int digits10 = 15;
3356 static const int max_digits10 = 16;
3357#endif
3358 };
3359
3360}
3361
3362#endif /* __MPREAL_H__ */
std::string why()
Definition mpreal.h:600
friend const mpreal rec_sqrt(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2672
friend const mpreal acot(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2314
friend const mpreal dim(const mpreal &a, const mpreal &b, mp_rnd_t rnd_mode)
Definition mpreal.h:2272
static mp_exp_t get_emax_min(void)
Definition mpreal.h:2205
friend const mpreal log1p(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2331
friend const mpreal acsc(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2316
friend const mpreal root(const mpreal &v, unsigned long int k, mp_rnd_t rnd_mode)
Definition mpreal.h:2261
friend const mpreal sum(const mpreal tab[], const unsigned long int n, int &status, mp_rnd_t rnd_mode)
Definition mpreal.h:2600
friend const mpreal remquo(const mpreal &x, const mpreal &y, int *q, mp_rnd_t rnd_mode)
Definition mpreal.h:2511
friend const mpreal log2(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2296
friend const mpreal fmax(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode)
Definition mpreal.h:2794
friend const mpreal sqr(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2231
friend const mpreal acos(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2307
friend const mpreal rint_trunc(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2776
friend const mpreal nextbelow(const mpreal &x)
Definition mpreal.h:2822
friend const mpreal acosh(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2327
friend const mpreal nextpow2(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2442
mpreal & operator=(const mpreal &v)
Definition mpreal.h:951
friend const mpreal sec(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2304
friend const mpreal li2(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2623
friend const mpreal cos(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2301
friend const mpreal mul_2ui(const mpreal &v, unsigned long int k, mp_rnd_t rnd_mode)
Definition mpreal.h:1654
friend const mpreal div_2si(const mpreal &v, long int k, mp_rnd_t rnd_mode)
Definition mpreal.h:1676
friend const mpreal agm(const mpreal &v1, const mpreal &v2, mp_rnd_t rnd_mode)
Definition mpreal.h:2585
friend const mpreal zeta(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2337
friend const mpreal rem(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode)
Definition mpreal.h:2628
friend const mpreal sin(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2302
mpreal & setNan()
Definition mpreal.h:2004
friend bool isfinite(const mpreal &v)
Definition mpreal.h:1746
const mpreal operator-() const
Definition mpreal.h:1307
mp_exp_t get_exp() const
Definition mpreal.h:2035
unsigned long toULong(mp_rnd_t mode=GMP_RNDZ) const
Definition mpreal.h:1758
friend const mpreal const_catalan(mp_prec_t prec, mp_rnd_t rnd_mode)
Definition mpreal.h:2710
friend const mpreal fmin(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode)
Definition mpreal.h:2801
static void set_default_rnd(mp_rnd_t rnd_mode)
Definition mpreal.h:2917
friend const mpreal coth(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2326
friend int sinh_cosh(mpreal &s, mpreal &c, const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2618
static mp_exp_t get_emin(void)
Definition mpreal.h:2175
friend const mpreal fac_ui(unsigned long int v, mp_prec_t prec, mp_rnd_t rnd_mode)
Definition mpreal.h:2520
friend const mpreal const_pi(mp_prec_t prec, mp_rnd_t rnd_mode)
Definition mpreal.h:2696
friend const mpreal abs(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2294
double toDouble(mp_rnd_t mode=GMP_RNDN) const
Definition mpreal.h:1760
std::ostream & output(std::ostream &os) const
Definition mpreal.h:1917
MPREAL_MSVC_DEBUGVIEW_DATA void clear(::mpfr_ptr)
Definition mpreal.h:770
friend const mpreal round(const mpreal &v)
Definition mpreal.h:2740
friend const mpreal exp2(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2299
mp_prec_t get_prec() const
Definition mpreal.h:2024
int subnormalize(int t, mp_rnd_t rnd_mode=get_default_rnd())
Definition mpreal.h:2168
const mpreal operator+() const
Definition mpreal.h:1204
friend const mpreal remainder(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode)
Definition mpreal.h:2504
mpreal & setSign(int Sign, mp_rnd_t RoundingMode=get_default_rnd())
Definition mpreal.h:1978
friend int sin_cos(mpreal &s, mpreal &c, const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2284
friend const mpreal atan(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2309
friend const mpreal besselj0(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2340
friend const mpreal operator/(const unsigned long int b, const mpreal &a)
Definition mpreal.h:1514
friend const mpreal acoth(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2317
friend const mpreal csc(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2305
friend const mpreal logb(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2311
mpreal & setInf(int Sign=+1)
Definition mpreal.h:1997
long double toLDouble(mp_rnd_t mode=GMP_RNDN) const
Definition mpreal.h:1761
friend const mpreal sqrt(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2234
long toLong(mp_rnd_t mode=GMP_RNDZ) const
Definition mpreal.h:1757
static int set_emin(mp_exp_t exp)
Definition mpreal.h:2180
friend const mpreal lgamma(const mpreal &v, int *signp, mp_rnd_t rnd_mode)
Definition mpreal.h:2529
friend const mpreal eint(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2333
::mpfr_srcptr mpfr_srcptr() const
Definition mpreal.h:1767
mpreal & operator+=(const mpreal &v)
Definition mpreal.h:1120
friend const mpreal cot(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2306
friend const mpreal asin(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2308
friend const mpreal pow(const mpreal &a, const mpreal &b, mp_rnd_t rnd_mode)
Definition mpreal.h:2929
friend const mpreal urandomb(gmp_randstate_t &state)
Definition mpreal.h:2829
friend const mpreal frac(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2777
mpreal & operator*=(const long long int u)
Definition mpreal.h:1199
mpreal & operator++()
Definition mpreal.h:1213
int set_exp(mp_exp_t e)
Definition mpreal.h:2040
friend const mpreal floor(const mpreal &v)
Definition mpreal.h:2733
friend const mpreal rint_floor(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2774
float toFloat(mp_rnd_t mode=GMP_RNDN) const
Definition mpreal.h:1759
static mp_exp_t get_emax(void)
Definition mpreal.h:2185
static mp_exp_t get_emax_max(void)
Definition mpreal.h:2210
friend const mpreal mod(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode)
Definition mpreal.h:2634
mpreal & operator--()
Definition mpreal.h:1225
friend long lround(const mpreal &v)
Definition mpreal.h:2747
mpreal & setZero(int Sign=+1)
Definition mpreal.h:2011
friend const mpreal exp10(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2300
friend const mpreal erf(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2338
friend const mpreal log(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2295
friend const mpreal rint_round(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2775
friend int sgn(const mpreal &v)
Definition mpreal.h:2781
friend const mpreal acsch(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2319
friend std::istream & operator>>(std::istream &is, mpreal &v)
Definition mpreal.h:1948
mpreal & operator>>=(const unsigned long int u)
Definition mpreal.h:1585
friend const mpreal fmod(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode)
Definition mpreal.h:2657
friend const mpreal ceil(const mpreal &v)
Definition mpreal.h:2726
friend const mpreal nextabove(const mpreal &x)
Definition mpreal.h:2815
friend const mpreal csch(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2325
friend const mpreal rint_ceil(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2773
friend const mpreal fms(const mpreal &v1, const mpreal &v2, const mpreal &v3, mp_rnd_t rnd_mode)
Definition mpreal.h:2570
friend const mpreal random(unsigned int seed)
Definition mpreal.h:2858
void set_prec(mp_prec_t prec, mp_rnd_t rnd_mode=get_default_rnd())
Definition mpreal.h:2029
friend const mpreal urandom(gmp_randstate_t &state, mp_rnd_t rnd_mode)
Definition mpreal.h:2837
friend const mpreal const_euler(mp_prec_t prec, mp_rnd_t rnd_mode)
Definition mpreal.h:2703
bool toBool() const
Definition mpreal.h:1756
friend const mpreal asec(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2315
friend bool isinf(const mpreal &v)
Definition mpreal.h:1745
mpreal & operator-=(const long long int u)
Definition mpreal.h:1197
friend const mpreal cosh(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2321
friend const mpreal gamma(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2334
friend const mpreal grandom(gmp_randstate_t &state, mp_rnd_t rnd_mode)
Definition mpreal.h:2881
friend const mpreal fabs(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2293
friend const mpreal rint(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2772
static int set_emax(mp_exp_t exp)
Definition mpreal.h:2190
friend const mpreal const_infinity(int sign, mp_prec_t prec)
Definition mpreal.h:2717
friend const mpreal besselyn(long n, const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2548
friend const mpreal erfc(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2339
static void set_default_prec(mp_prec_t prec)
Definition mpreal.h:2912
friend const mpreal const_log2(mp_prec_t prec, mp_rnd_t rnd_mode)
Definition mpreal.h:2689
mpreal & operator/=(const long long int u)
Definition mpreal.h:1201
friend const mpreal besselj1(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2341
int check_range(int t, mp_rnd_t rnd_mode=get_default_rnd())
Definition mpreal.h:2163
long long toLLong(mp_rnd_t mode=GMP_RNDZ) const
Definition mpreal.h:1762
friend bool isint(const mpreal &v)
Definition mpreal.h:1748
friend bool isregular(const mpreal &v)
Definition mpreal.h:1751
static mp_exp_t get_emin_min(void)
Definition mpreal.h:2195
friend const mpreal ai(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2684
static mp_prec_t get_default_prec()
Definition mpreal.h:183
int getPrecision() const
Definition mpreal.h:1985
friend const mpreal trunc(const mpreal &v)
Definition mpreal.h:2765
friend const mpreal modf(const mpreal &v, mpreal &n)
Definition mpreal.h:2153
friend bool isnan(const mpreal &v)
Definition mpreal.h:1744
friend int cmpabs(const mpreal &a, const mpreal &b)
Definition mpreal.h:2279
friend const mpreal bessely1(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2343
friend const mpreal fma(const mpreal &v1, const mpreal &v2, const mpreal &v3, mp_rnd_t rnd_mode)
Definition mpreal.h:2555
unsigned long long toULLong(mp_rnd_t mode=GMP_RNDZ) const
Definition mpreal.h:1763
friend bool iszero(const mpreal &v)
Definition mpreal.h:1747
mpreal & operator<<=(const unsigned long int u)
Definition mpreal.h:1557
friend const mpreal asinh(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2328
friend const mpreal cbrt(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2292
friend const mpreal lngamma(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2336
friend const mpreal digamma(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2683
static mp_exp_t get_emin_max(void)
Definition mpreal.h:2200
static mp_rnd_t get_default_rnd()
Definition mpreal.h:182
friend const mpreal atan2(const mpreal &y, const mpreal &x, mp_rnd_t rnd_mode)
Definition mpreal.h:2452
friend const mpreal hypot(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode)
Definition mpreal.h:2459
std::string toString(int n=-1, int b=10, mp_rnd_t mode=mpreal::get_default_rnd()) const
Definition mpreal.h:1799
mpfr_t mp
Definition mpreal.h:177
friend const mpreal besseljn(long n, const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2541
friend const mpreal sech(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2324
friend const mpreal div_2ui(const mpreal &v, unsigned long int k, mp_rnd_t rnd_mode)
Definition mpreal.h:1669
friend const mpreal bessely0(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2342
friend bool isnum(const mpreal &v)
friend const mpreal atanh(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2329
friend mp_exp_t ilogb(const mpreal &v)
Definition mpreal.h:2312
mpreal & setPrecision(int Precision, mp_rnd_t RoundingMode=get_default_rnd())
Definition mpreal.h:1990
friend const mpreal mul_2si(const mpreal &v, long int k, mp_rnd_t rnd_mode)
Definition mpreal.h:1662
friend const mpreal expm1(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2332
bool fits_in_bits(double x, int n)
Definition mpreal.h:2922
::mpfr_ptr mpfr_ptr()
Definition mpreal.h:1765
friend const mpreal tgamma(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2335
friend const mpreal random2(mp_size_t size, mp_exp_t exp)
Definition mpreal.h:2846
friend const mpreal exp(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2298
friend const mpreal asech(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2318
friend const mpreal tan(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2303
friend const mpreal nexttoward(const mpreal &x, const mpreal &y)
Definition mpreal.h:2808
friend long long llround(const mpreal &v)
Definition mpreal.h:2756
friend const mpreal sinh(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2322
friend const mpreal tanh(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2323
friend const mpreal log10(const mpreal &v, mp_rnd_t rnd_mode)
Definition mpreal.h:2297
static const mpfr::mpreal signaling_NaN()
Definition mpreal.h:3291
static const bool tinyness_before
Definition mpreal.h:3267
static const bool has_quiet_NaN
Definition mpreal.h:3260
static int digits10(mp_prec_t precision=mpfr::mpreal::get_default_prec())
Definition mpreal.h:3325
static const mpfr::mpreal quiet_NaN()
Definition mpreal.h:3290
static const bool has_signaling_NaN
Definition mpreal.h:3261
static const float_denorm_style has_denorm
Definition mpreal.h:3269
static mpfr::mpreal max(mp_prec_t precision=mpfr::mpreal::get_default_prec())
Definition mpreal.h:3272
static int max_digits10(mp_prec_t precision=mpfr::mpreal::get_default_prec())
Definition mpreal.h:3335
static mpfr::mpreal round_error(mp_prec_t precision=mpfr::mpreal::get_default_prec())
Definition mpreal.h:3281
static const mpfr::mpreal denorm_min()
Definition mpreal.h:3292
static const bool is_specialized
Definition mpreal.h:3253
static mpfr::mpreal epsilon(mp_prec_t precision=mpfr::mpreal::get_default_prec())
Definition mpreal.h:3276
static MPREAL_PERMISSIVE_EXPR const int min_exponent10
Definition mpreal.h:3297
static const mpfr::mpreal infinity()
Definition mpreal.h:3289
static int digits(const mpfr::mpreal &x)
Definition mpreal.h:3323
static MPREAL_PERMISSIVE_EXPR const int max_exponent10
Definition mpreal.h:3298
static int digits10(const mpfr::mpreal &x)
Definition mpreal.h:3330
static float_round_style round_style()
Definition mpreal.h:3308
static mpfr::mpreal min(mp_prec_t precision=mpfr::mpreal::get_default_prec())
Definition mpreal.h:3271
static mpfr::mpreal epsilon(const mpfr::mpreal &x)
Definition mpreal.h:3279
static mpfr::mpreal lowest(mp_prec_t precision=mpfr::mpreal::get_default_prec())
Definition mpreal.h:3273
static CanonicalForm base
Definition factory.cpp:289
int p
int p1
void size_t s
Definition m2-mem.cpp:271
std::ostream & operator<<(std::ostream &o, const std::pair< T1, T2 > &p)
std::ostringstream & toString(std::ostringstream &o, int len, int *p)
#define swap(a, b, t)
Definition monsort.hpp:127
#define MPREAL_UNARY_MATH_FUNCTION_BODY(f)
Definition mpreal.h:2220
#define MPREAL_PERMISSIVE_EXPR
Definition mpreal.h:170
#define mpfr_set_zero_fast(x)
Definition mpreal.h:165
#define MPREAL_BINARY_MATH_FUNCTION_UI_BODY(f, u)
Definition mpreal.h:2226
#define MPREAL_MSVC_DEBUGVIEW_DATA
Definition mpreal.h:147
#define IsInf(x)
Definition mpreal.h:108
#define MPREAL_MSVC_DEBUGVIEW_CODE
Definition mpreal.h:146
#define MPREAL_DOUBLE_BITS_OVERFLOW
Definition mpreal.h:158
bool operator!=(const StatsAllocator< T1 > &, const StatsAllocator< T2 > &) noexcept
Definition myalloc.hpp:134
bool operator==(const StatsAllocator< T1 > &, const StatsAllocator< T2 > &) noexcept
Definition myalloc.hpp:129
const mpreal ldexp(const mpreal &v, mp_exp_t exp)
Definition mpreal.h:2073
const mpreal sqrt(const unsigned int v, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2244
const mpreal div_2ui(const mpreal &v, unsigned long int k, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:1669
const mpreal grandom(gmp_randstate_t &state, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2881
const mpreal asinu(const mpreal &x, unsigned long u, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2384
const mpreal min(const mpreal &x, const mpreal &y)
Definition mpreal.h:2792
const mpreal cosu(const mpreal &x, unsigned long u, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2380
const mpreal sinpi(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2389
bool signbit(const mpreal &x)
Definition mpreal.h:2142
int bits2digits(mp_prec_t b)
Definition mpreal.h:1969
const mpreal urandom(gmp_randstate_t &state, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2837
const internal::result_type< Rhs >::type operator+(const mpreal &lhs, const Rhs &rhs)
Definition mpreal.h:807
const mpreal cospi(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2388
const mpreal mul_2ui(const mpreal &v, unsigned long int k, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:1654
const mpreal acosu(const mpreal &x, unsigned long u, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2383
const mpreal hypot(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2459
const mpreal scalbn(const mpreal &v, mp_exp_t exp)
Definition mpreal.h:2082
const internal::result_type< Rhs >::type operator-(const mpreal &lhs, const Rhs &rhs)
Definition mpreal.h:816
bool isnan(const mpreal &op)
Definition mpreal.h:1744
mpreal copysign(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2135
mpreal minval(mp_prec_t prec=mpreal::get_default_prec())
Definition mpreal.h:2105
const mpreal gammainc(const mpreal &a, const mpreal &x, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2346
const internal::result_type< Rhs >::type operator*(const mpreal &lhs, const Rhs &rhs)
Definition mpreal.h:825
const mpreal atan2pi(const mpreal &y, const mpreal &x, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2411
const mpreal atanpi(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2393
bool isEqualUlps(const mpreal &a, const mpreal &b, int maxUlps)
Definition mpreal.h:2118
const mpreal fmod(const mpreal &x, unsigned long u, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2435
const mpreal exp10m1(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2398
bool isEqualFuzzy(const mpreal &a, const mpreal &b, const mpreal &eps)
Definition mpreal.h:2123
mp_prec_t digits2bits(int d)
Definition mpreal.h:1962
const mpreal compound(const mpreal &x, long n, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2427
mpreal & negate(mpreal &x)
Definition mpreal.h:2047
const mpreal atan2u(const mpreal &y, const mpreal &x, unsigned long u, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2400
const mpreal log2p1(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2395
const mpreal beta(const mpreal &z, const mpreal &w, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2357
const mpreal tanpi(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2390
const internal::result_type< Rhs >::type operator/(const mpreal &lhs, const Rhs &rhs)
Definition mpreal.h:834
const mpreal mul_2si(const mpreal &v, long int k, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:1662
const mpreal acospi(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2391
const mpreal log10p1(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2396
const mpreal div_2si(const mpreal &v, long int k, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:1676
const mpreal exp2m1(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2397
bool iszero(const mpreal &op)
Definition mpreal.h:1747
mpreal machine_epsilon(mp_prec_t prec=mpreal::get_default_prec())
Definition mpreal.h:2087
const mpreal asinpi(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2392
const mpreal log_ui(unsigned long int n, mp_prec_t prec=mpreal::get_default_prec(), mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2368
const mpreal frexp(const mpreal &x, mp_exp_t *exp, mp_rnd_t mode=mpreal::get_default_rnd())
Definition mpreal.h:2053
const mpreal exp(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2298
mpreal & setsignbit(mpreal &x, bool minus, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2147
mpreal maxval(mp_prec_t prec=mpreal::get_default_prec())
Definition mpreal.h:2112
const mpreal tanu(const mpreal &x, unsigned long u, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2382
const mpreal sinu(const mpreal &x, unsigned long u, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2381
const mpreal sqr(const mpreal &x, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2231
const mpreal atanu(const mpreal &x, unsigned long u, mp_rnd_t r=mpreal::get_default_rnd())
Definition mpreal.h:2385
const mpreal max(const mpreal &x, const mpreal &y)
Definition mpreal.h:2791
const mpreal powr(const mpreal &x, const mpreal &y, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2419
const mpreal operator>>(const mpreal &v, const unsigned long int k)
Definition mpreal.h:1633
const mpreal nextabove(const mpreal &x)
Definition mpreal.h:2815
const mpreal pow(const mpreal &a, const unsigned int b, mp_rnd_t rnd_mode=mpreal::get_default_rnd())
Definition mpreal.h:2962
Definition mpreal.h:173
STL namespace.
volatile int x
#define abs(x)
Definition polyroots.cpp:51
#define max(a, b)
Definition polyroots.cpp:52
static gmp_randstate_t state
Definition random.cpp:18
#define T
Definition table.c:13