1 /*
2   This file is part of MADNESS.
3 
4   Copyright (C) 2007,2010 Oak Ridge National Laboratory
5 
6   This program is free software; you can redistribute it and/or modify
7   it under the terms of the GNU General Public License as published by
8   the Free Software Foundation; either version 2 of the License, or
9   (at your option) any later version.
10 
11   This program is distributed in the hope that it will be useful,
12   but WITHOUT ANY WARRANTY; without even the implied warranty of
13   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14   GNU General Public License for more details.
15 
16   You should have received a copy of the GNU General Public License
17   along with this program; if not, write to the Free Software
18   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 
20   For more information please contact:
21 
22   Robert J. Harrison
23   Oak Ridge National Laboratory
24   One Bethel Valley Road
25   P.O. Box 2008, MS-6367
26 
27   email: harrisonrj@ornl.gov
28   tel:   865-241-3937
29   fax:   865-572-0680
30 
31   $Id$
32 */
33 /*
34 	Multi-precision real number class. C++ wrapper fo MPFR library.
35 	Project homepage: http://www.holoborodko.com/pavel/
36 	Contact e-mail:   pavel@holoborodko.com
37 
38 	Copyright (c) 2008 Pavel Holoborodko
39 
40 	This library is free software; you can redistribute it and/or
41 	modify it under the terms of the GNU Lesser General Public
42 	License as published by the Free Software Foundation; either
43 	version 2.1 of the License, or (at your option) any later version.
44 
45 	This library is distributed in the hope that it will be useful,
46 	but WITHOUT ANY WARRANTY; without even the implied warranty of
47 	MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
48 	Lesser General Public License for more details.
49 
50 	You should have received a copy of the GNU Lesser General Public
51 	License along with this library; if not, write to the Free Software
52 	Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
53 
54 	Contributors:
55 	Brain Gladman, Helmut Jarausch, Fokko Beekhof, Ulrich Mutze.
56 */
57 
58 #ifndef __MP_REAL_H__
59 #define __MP_REAL_H__
60 
61 #include <string>
62 #include <iostream>
63 #include <stdexcept>
64 #include <cfloat>
65 #include <cmath>
66 
67 #include "mpfr.h"
68 
69 // Detect compiler using signatures from http://predef.sourceforge.net/
70 // GNU C/C++
71 #if defined(__GNUC__)
72 	#define IsInf(x) std::isinf(x)				// C99 supported by C++ 2003 standard
73 
74 // Microsoft Visual C++
75 #elif defined(_MSC_VER)
76 	#define IsInf(x) (!_finite(x))				// Microsoft specific
77 
78 // Others, will be added as needed
79 #else
80 	#define IsInf(x) std::isinf(x)				// now just hope for the C99 conformance
81 #endif
82 
83 namespace mpfr {
84 
85 class mpreal {
86 private:
87 	mpfr_t mp;
88 
89 public:
90 	static mp_rnd_t   default_rnd;
91 	static mp_prec_t  default_prec;
92 	static int		  default_base;
93     static int        double_bits;
94 
95 public:
96 	// Constructors && type conversion
97 	mpreal();
98 	mpreal(const mpreal& u);
99 
100 	mpreal(const mpfr_t u);
101 	mpreal(const mpf_t u);
102 
103 	mpreal(const mpz_t u, mp_prec_t prec = default_prec, mp_rnd_t mode = default_rnd);
104 	mpreal(const mpq_t u, mp_prec_t prec = default_prec, mp_rnd_t mode = default_rnd);
105 	mpreal(const double u, mp_prec_t prec = default_prec, mp_rnd_t mode = default_rnd);
106 	mpreal(const long double u, mp_prec_t prec = default_prec, mp_rnd_t mode = default_rnd);
107 	mpreal(const unsigned long int u, mp_prec_t prec = default_prec, mp_rnd_t mode = default_rnd);
108 	mpreal(const unsigned int u, mp_prec_t prec = default_prec, mp_rnd_t mode = default_rnd);
109 	mpreal(const long int u, mp_prec_t prec = default_prec, mp_rnd_t mode = default_rnd);
110 	mpreal(const int u, mp_prec_t prec = default_prec, mp_rnd_t mode = default_rnd);
111 	mpreal(const char* s, mp_prec_t prec = default_prec, int base = default_base, mp_rnd_t mode = default_rnd);
112 
113 	~mpreal();
114 
115 	// Operations
116 	// =
117 	// +, -, *, /, ++, --, <<, >>
118 	// *=, +=, -=, /=,
119 	// <, >, ==, <=, >=
120 
121 	// =
122 	mpreal& operator=(const mpreal& v);
123 	mpreal& operator=(const mpf_t v);
124 	mpreal& operator=(const mpz_t v);
125 	mpreal& operator=(const mpq_t v);
126 	mpreal& operator=(const long double v);
127 	mpreal& operator=(const double v);
128 	mpreal& operator=(const unsigned long int v);
129 	mpreal& operator=(const unsigned int v);
130 	mpreal& operator=(const long int v);
131 	mpreal& operator=(const int v);
132 	mpreal& operator=(const char* s);
133 
134 	// +
135 	mpreal& operator+=(const mpreal& v);
136 	mpreal& operator+=(const mpf_t v);
137 	mpreal& operator+=(const mpz_t v);
138 	mpreal& operator+=(const mpq_t v);
139 	mpreal& operator+=(const long double u);
140 	mpreal& operator+=(const double u);
141 	mpreal& operator+=(const unsigned long int u);
142 	mpreal& operator+=(const unsigned int u);
143 	mpreal& operator+=(const long int u);
144 	mpreal& operator+=(const int u);
145 	const mpreal operator+() const;
146 	mpreal& operator++ ();
147 	const mpreal  operator++ (int);
148 
149 	// -
150 	mpreal& operator-=(const mpreal& v);
151 	mpreal& operator-=(const mpz_t v);
152 	mpreal& operator-=(const mpq_t v);
153 	mpreal& operator-=(const long double u);
154 	mpreal& operator-=(const double u);
155 	mpreal& operator-=(const unsigned long int u);
156 	mpreal& operator-=(const unsigned int u);
157 	mpreal& operator-=(const long int u);
158 	mpreal& operator-=(const int u);
159 	const mpreal operator-() const;
160 	friend const mpreal operator-(const unsigned long int b, const mpreal& a);
161 	friend const mpreal operator-(const unsigned int b, const mpreal& a);
162 	friend const mpreal operator-(const long int b, const mpreal& a);
163 	friend const mpreal operator-(const int b, const mpreal& a);
164 	friend const mpreal operator-(const double b, const mpreal& a);
165 	mpreal& operator-- ();
166 	const mpreal  operator-- (int);
167 
168 	// *
169 	mpreal& operator*=(const mpreal& v);
170 	mpreal& operator*=(const mpz_t v);
171 	mpreal& operator*=(const mpq_t v);
172 	mpreal& operator*=(const long double v);
173 	mpreal& operator*=(const double v);
174 	mpreal& operator*=(const unsigned long int v);
175 	mpreal& operator*=(const unsigned int v);
176 	mpreal& operator*=(const long int v);
177 	mpreal& operator*=(const int v);
178 
179 	// /
180 	mpreal& operator/=(const mpreal& v);
181 	mpreal& operator/=(const mpz_t v);
182 	mpreal& operator/=(const mpq_t v);
183 	mpreal& operator/=(const long double v);
184 	mpreal& operator/=(const double v);
185 	mpreal& operator/=(const unsigned long int v);
186 	mpreal& operator/=(const unsigned int v);
187 	mpreal& operator/=(const long int v);
188 	mpreal& operator/=(const int v);
189 	friend const mpreal operator/(const unsigned long int b, const mpreal& a);
190 	friend const mpreal operator/(const unsigned int b, const mpreal& a);
191 	friend const mpreal operator/(const long int b, const mpreal& a);
192 	friend const mpreal operator/(const int b, const mpreal& a);
193 	friend const mpreal operator/(const double b, const mpreal& a);
194 
195 	//<<= Fast Multiplication by 2^u
196 	mpreal& operator<<=(const unsigned long int u);
197 	mpreal& operator<<=(const unsigned int u);
198 	mpreal& operator<<=(const long int u);
199 	mpreal& operator<<=(const int u);
200 
201 	//>>= Fast Division by 2^u
202 	mpreal& operator>>=(const unsigned long int u);
203 	mpreal& operator>>=(const unsigned int u);
204 	mpreal& operator>>=(const long int u);
205 	mpreal& operator>>=(const int u);
206 
207 	// Boolean Operators
208 	friend bool operator >  (const mpreal& a, const mpreal& b);
209 	friend bool operator >= (const mpreal& a, const mpreal& b);
210 	friend bool operator <  (const mpreal& a, const mpreal& b);
211 	friend bool operator <= (const mpreal& a, const mpreal& b);
212 	friend bool operator == (const mpreal& a, const mpreal& b);
213 	friend bool operator != (const mpreal& a, const mpreal& b);
214 
215 	// Type Conversion operators
216 	operator long double() const;
217 	operator double() const;
218 	operator float() const;
219 	operator unsigned long() const;
220 	operator unsigned int() const;
221 	operator long() const;
222 	operator std::string() const;
223 	operator mpfr_ptr();
224 
225 	// Math Functions
226 	friend const mpreal sqr(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
227 	friend const mpreal sqrt(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
228 	friend const mpreal sqrt(const unsigned long int v, mp_rnd_t rnd_mode = mpreal::default_rnd);
229 	friend const mpreal cbrt(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
230 	friend const mpreal root(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode = mpreal::default_rnd);
231 	friend const mpreal pow(const mpreal& a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::default_rnd);
232 	friend const mpreal pow(const mpreal& a, const mpz_t b, mp_rnd_t rnd_mode = mpreal::default_rnd);
233 	friend const mpreal pow(const mpreal& a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
234 	friend const mpreal pow(const mpreal& a, const long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
235 	friend const mpreal pow(const unsigned long int a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::default_rnd);
236 	friend const mpreal pow(const unsigned long int a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
237 	friend const mpreal fabs(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
238 	friend const mpreal abs(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
239 	friend const mpreal dim(const mpreal& a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::default_rnd);
240 	friend const mpreal mul_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode = mpreal::default_rnd);
241 	friend const mpreal mul_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode = mpreal::default_rnd);
242 	friend const mpreal div_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode = mpreal::default_rnd);
243 	friend const mpreal div_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode = mpreal::default_rnd);
244 	friend int cmpabs(const mpreal& a,const mpreal& b);
245 
246 	friend const mpreal log  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
247 	friend const mpreal log2 (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
248 	friend const mpreal log10(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
249 	friend const mpreal exp  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
250 	friend const mpreal exp2 (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
251 	friend const mpreal exp10(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
252 
253 	friend const mpreal cos(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
254 	friend const mpreal sin(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
255 	friend const mpreal tan(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
256 	friend const mpreal sec(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
257 	friend const mpreal csc(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
258 	friend const mpreal cot(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
259 	friend int sin_cos(mpreal& s, mpreal& c, const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
260 
261 	friend const mpreal acos  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
262 	friend const mpreal asin  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
263 	friend const mpreal atan  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
264 	friend const mpreal atan2 (const mpreal& y, const mpreal& x, mp_rnd_t rnd_mode = mpreal::default_rnd);
265 	friend const mpreal cosh  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
266 	friend const mpreal sinh  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
267 	friend const mpreal tanh  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
268 	friend const mpreal sech  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
269 	friend const mpreal csch  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
270 	friend const mpreal coth  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
271 
272 	friend const mpreal acosh  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
273 	friend const mpreal asinh  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
274 	friend const mpreal atanh  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
275 	friend const mpreal fac_ui (unsigned long int v, mp_rnd_t rnd_mode = mpreal::default_rnd);
276 	friend const mpreal log1p  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
277 	friend const mpreal expm1  (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
278 	friend const mpreal eint   (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
279 
280 	friend const mpreal gamma (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
281 	friend const mpreal lngamma (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
282 	friend const mpreal lgamma (const mpreal& v, int *signp, mp_rnd_t rnd_mode = mpreal::default_rnd);
283 	friend const mpreal zeta (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
284 	friend const mpreal erf (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
285 	friend const mpreal erfc (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
286 	friend const mpreal _j0 (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
287 	friend const mpreal _j1 (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
288 	friend const mpreal _jn (long n, const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
289 	friend const mpreal _y0 (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
290 	friend const mpreal _y1 (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
291 	friend const mpreal _yn (long n, const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
292 	friend const mpreal fma (const mpreal& v1, const mpreal& v2, const mpreal& v3, mp_rnd_t rnd_mode = mpreal::default_rnd);
293 	friend const mpreal fms (const mpreal& v1, const mpreal& v2, const mpreal& v3, mp_rnd_t rnd_mode = mpreal::default_rnd);
294 	friend const mpreal agm (const mpreal& v1, const mpreal& v2, mp_rnd_t rnd_mode = mpreal::default_rnd);
295 	friend const mpreal hypot (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::default_rnd);
296 	friend const mpreal sum (const mpreal tab[], unsigned long int n, mp_rnd_t rnd_mode = mpreal::default_rnd);
297 	friend int sgn(const mpreal& v); // -1 or +1
298 
299 // MPFR 2.4.0 Specifics
300 #if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
301 	friend int sinh_cosh(mpreal& s, mpreal& c, const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
302 	friend const mpreal li2(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
303 	friend const mpreal fmod (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::default_rnd);
304 	friend const mpreal rec_sqrt(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
305 #endif
306 
307 	// Exponent and mantissa manipulation
308 	friend const mpreal frexp(const mpreal& v, mp_exp_t* exp);
309 	friend const mpreal ldexp(const mpreal& v, mp_exp_t exp);
310 
311 	// Splits mpreal value into fractional and integer parts.
312 	// Returns fractional part and stores integer part in n.
313 	friend const mpreal modf(const mpreal& v, mpreal& n);
314 
315 	// Constants
316 	// don't forget to call mpfr_free_cache() for every thread where you are using const-functions
317 	friend const mpreal const_log2 (mp_prec_t prec = mpreal::default_prec, mp_rnd_t rnd_mode = mpreal::default_rnd);
318 	friend const mpreal const_pi (mp_prec_t prec = mpreal::default_prec, mp_rnd_t rnd_mode = mpreal::default_rnd);
319 	friend const mpreal const_euler (mp_prec_t prec = mpreal::default_prec, mp_rnd_t rnd_mode = mpreal::default_rnd);
320 	friend const mpreal const_catalan (mp_prec_t prec = mpreal::default_prec, mp_rnd_t rnd_mode = mpreal::default_rnd);
321 
322 	// Output/ Input
323 	friend std::ostream& operator<<(std::ostream& os, const mpreal& v);
324     friend std::istream& operator>>(std::istream& is, mpreal& v);
325 
326 	// Integer Related Functions
327 	friend const mpreal rint (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
328 	friend const mpreal ceil (const mpreal& v);
329 	friend const mpreal floor(const mpreal& v);
330 	friend const mpreal round(const mpreal& v);
331 	friend const mpreal trunc(const mpreal& v);
332 	friend const mpreal rint_ceil (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
333 	friend const mpreal rint_floor(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
334 	friend const mpreal rint_round(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
335 	friend const mpreal rint_trunc(const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
336 	friend const mpreal frac (const mpreal& v, mp_rnd_t rnd_mode = mpreal::default_rnd);
337 	friend const mpreal remainder (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::default_rnd);
338 	friend const mpreal remquo (long* q, const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode = mpreal::default_rnd);
339 
340 	// Miscellaneous Functions
341 	friend const mpreal nexttoward (const mpreal& x, const mpreal& y);
342 	friend const mpreal nextabove  (const mpreal& x);
343 	friend const mpreal nextbelow  (const mpreal& x);
344 
345 	// use gmp_randinit_default() to init state, gmp_randclear() to clear
346 	friend const mpreal urandomb (gmp_randstate_t& state);
347 	friend const mpreal random2 (mp_size_t size, mp_exp_t exp);
348 
349 	// Instance Checkers
350 	friend bool _isnan(const mpreal& v);
351 	friend bool _isinf(const mpreal& v);
352 	friend bool _isnum(const mpreal& v);
353 	friend bool _iszero(const mpreal& v);
354 	friend bool _isint(const mpreal& v);
355 
356 	// Set/Get instance properties
357 	mp_prec_t	get_prec() const;
358 	void		set_prec(mp_prec_t prec, mp_rnd_t rnd_mode = default_rnd);	// Change precision with rounding mode
359 
360 	// Set mpreal to +-inf, NaN
361 	void      set_inf(int sign = +1);
362 	void	  set_nan();
363 
364 	// sign = -1 or +1
365 	void set_sign(int sign, mp_rnd_t rnd_mode = default_rnd);
366 
367 	//Exponent
368 	mp_exp_t get_exp();
369 	int set_exp(mp_exp_t e);
370 	int check_range (int t, mp_rnd_t rnd_mode = default_rnd);
371 	int subnormalize (int t,mp_rnd_t rnd_mode = default_rnd);
372 
373 	// Inexact conversion from float
374 	bool fits_in_bits(double x, int n);
375 
376 	// Set/Get global properties
377 	static void			set_default_prec(mp_prec_t prec);
378 	static mp_prec_t	get_default_prec();
379 	static void			set_default_base(int base);
380 	static int			get_default_base();
381 	static void			set_double_bits(int dbits);
382 	static int			get_double_bits();
383 	static void			set_default_rnd(mp_rnd_t rnd_mode);
384 	static mp_rnd_t		get_default_rnd();
385 	static mp_exp_t get_emin (void);
386 	static mp_exp_t get_emax (void);
387 	static mp_exp_t get_emin_min (void);
388 	static mp_exp_t get_emin_max (void);
389 	static mp_exp_t get_emax_min (void);
390 	static mp_exp_t get_emax_max (void);
391 	static int set_emin (mp_exp_t exp);
392 	static int set_emax (mp_exp_t exp);
393 
394 	// Get/Set conversions
395 	// Convert mpreal to string with n significant digits in base b
396 	// n = 0 -> convert with the maximum available digits
397 	std::string to_string(size_t n = 0, int b = default_base, mp_rnd_t mode = default_rnd) const;
398 
399 	// Efficient swapping of two mpreal values
400 	friend void swap(mpreal& x, mpreal& y);
401 
402 	//Min Max - macros is evil. Needed for systems which defines max and min globally as macros (e.g. Windows)
403 	//Hope that globally defined macros uses > < operations only
404 	#ifndef max
405 		friend const mpreal max(const mpreal& x, const mpreal& y);
406 	#endif
407 
408 	#ifndef min
409 		friend const mpreal min(const mpreal& x, const mpreal& y);
410 	#endif
411 };
412 
413 //////////////////////////////////////////////////////////////////////////
414 // Exceptions
415 class conversion_overflow : public std::exception {
416 public:
why()417 	std::string why() { return "inexact conversion from floating point"; }
418 };
419 
420 //////////////////////////////////////////////////////////////////////////
421 // + Addition
422 const mpreal operator+(const mpreal& a, const mpreal& b);
423 
424 // + Fast specialized addition - implemented through fast += operations
425 const mpreal operator+(const mpreal& a, const mpz_t b);
426 const mpreal operator+(const mpreal& a, const mpq_t b);
427 const mpreal operator+(const mpreal& a, const long double b);
428 const mpreal operator+(const mpreal& a, const double b);
429 const mpreal operator+(const mpreal& a, const unsigned long int b);
430 const mpreal operator+(const mpreal& a, const unsigned int b);
431 const mpreal operator+(const mpreal& a, const long int b);
432 const mpreal operator+(const mpreal& a, const int b);
433 const mpreal operator+(const mpreal& a, const char* b);
434 const mpreal operator+(const char* a, const mpreal& b);
435 const std::string operator+(const mpreal& a, const std::string b);
436 const std::string operator+(const std::string a, const mpreal& b);
437 
438 const mpreal operator+(const mpz_t b, const mpreal& a);
439 const mpreal operator+(const mpq_t b, const mpreal& a);
440 const mpreal operator+(const long double b, const mpreal& a);
441 const mpreal operator+(const double  b, const mpreal& a);
442 const mpreal operator+(const unsigned long int b, const mpreal& a);
443 const mpreal operator+(const unsigned int b, const mpreal& a);
444 const mpreal operator+(const long int b, const mpreal& a);
445 const mpreal operator+(const int b, const mpreal& a);
446 
447 //////////////////////////////////////////////////////////////////////////
448 // - Subtraction
449 const mpreal operator-(const mpreal& a, const mpreal& b);
450 
451 // - Fast specialized subtraction - implemented through fast -= operations
452 const mpreal operator-(const mpreal& a, const mpz_t b);
453 const mpreal operator-(const mpreal& a, const mpq_t b);
454 const mpreal operator-(const mpreal& a, const long double b);
455 const mpreal operator-(const mpreal& a, const double b);
456 const mpreal operator-(const mpreal& a, const unsigned long int b);
457 const mpreal operator-(const mpreal& a, const unsigned int b);
458 const mpreal operator-(const mpreal& a, const long int b);
459 const mpreal operator-(const mpreal& a, const int b);
460 const mpreal operator-(const mpreal& a, const char* b);
461 const mpreal operator-(const char* a, const mpreal& b);
462 
463 const mpreal operator-(const mpz_t b, const mpreal& a);
464 const mpreal operator-(const mpq_t b, const mpreal& a);
465 const mpreal operator-(const long double b, const mpreal& a);
466 //const mpreal operator-(const double  b, const mpreal& a);
467 
468 //////////////////////////////////////////////////////////////////////////
469 // * Multiplication
470 const mpreal operator*(const mpreal& a, const mpreal& b);
471 
472 // * Fast specialized multiplication - implemented through fast *= operations
473 const mpreal operator*(const mpreal& a, const mpz_t b);
474 const mpreal operator*(const mpreal& a, const mpq_t b);
475 const mpreal operator*(const mpreal& a, const long double b);
476 const mpreal operator*(const mpreal& a, const double b);
477 const mpreal operator*(const mpreal& a, const unsigned long int b);
478 const mpreal operator*(const mpreal& a, const unsigned int b);
479 const mpreal operator*(const mpreal& a, const long int b);
480 const mpreal operator*(const mpreal& a, const int b);
481 
482 const mpreal operator*(const mpz_t b, const mpreal& a);
483 const mpreal operator*(const mpq_t b, const mpreal& a);
484 const mpreal operator*(const long double b, const mpreal& a);
485 const mpreal operator*(const double  b, const mpreal& a);
486 const mpreal operator*(const unsigned long int b, const mpreal& a);
487 const mpreal operator*(const unsigned int b, const mpreal& a);
488 const mpreal operator*(const long int b, const mpreal& a);
489 const mpreal operator*(const int b, const mpreal& a);
490 
491 //////////////////////////////////////////////////////////////////////////
492 // / Division
493 const mpreal operator/(const mpreal& a, const mpreal& b);
494 
495 // / Fast specialized division - implemented through fast /= operations
496 const mpreal operator/(const mpreal& a, const mpz_t b);
497 const mpreal operator/(const mpreal& a, const mpq_t b);
498 const mpreal operator/(const mpreal& a, const long double b);
499 const mpreal operator/(const mpreal& a, const double b);
500 const mpreal operator/(const mpreal& a, const unsigned long int b);
501 const mpreal operator/(const mpreal& a, const unsigned int b);
502 const mpreal operator/(const mpreal& a, const long int b);
503 const mpreal operator/(const mpreal& a, const int b);
504 
505 const mpreal operator/(const long double b, const mpreal& a);
506 
507 //////////////////////////////////////////////////////////////////////////
508 // Shifts operators - Multiplication/Division by a power of 2
509 const mpreal operator<<(const mpreal& v, const unsigned long int k);
510 const mpreal operator<<(const mpreal& v, const unsigned int k);
511 const mpreal operator<<(const mpreal& v, const long int k);
512 const mpreal operator<<(const mpreal& v, const int k);
513 
514 const mpreal operator>>(const mpreal& v, const unsigned long int k);
515 const mpreal operator>>(const mpreal& v, const unsigned int k);
516 const mpreal operator>>(const mpreal& v, const long int k);
517 const mpreal operator>>(const mpreal& v, const int k);
518 
519 //////////////////////////////////////////////////////////////////////////
520 // Boolean operators
521 bool operator <  (const mpreal& a, const unsigned long int b);
522 bool operator <  (const mpreal& a, const unsigned int b);
523 bool operator <  (const mpreal& a, const long int b);
524 bool operator <  (const mpreal& a, const int b);
525 bool operator <  (const mpreal& a, const long double b);
526 bool operator <  (const mpreal& a, const double b);
527 
528 bool operator <  (const unsigned long int a,const mpreal& b);
529 bool operator <  (const unsigned int a,		const mpreal& b);
530 bool operator <  (const long int a,			const mpreal& b);
531 bool operator <  (const int a,				const mpreal& b);
532 bool operator <  (const long double a,		const mpreal& b);
533 bool operator <  (const double a,			const mpreal& b);
534 
535 bool operator >  (const mpreal& a, const unsigned long int b);
536 bool operator >  (const mpreal& a, const unsigned int b);
537 bool operator >  (const mpreal& a, const long int b);
538 bool operator >  (const mpreal& a, const int b);
539 bool operator >  (const mpreal& a, const long double b);
540 bool operator >  (const mpreal& a, const double b);
541 
542 bool operator >  (const unsigned long int a,const mpreal& b);
543 bool operator >  (const unsigned int a,		const mpreal& b);
544 bool operator >  (const long int a,			const mpreal& b);
545 bool operator >  (const int a,				const mpreal& b);
546 bool operator >  (const long double a,		const mpreal& b);
547 bool operator >  (const double a,			const mpreal& b);
548 
549 bool operator >=  (const mpreal& a, const unsigned long int b);
550 bool operator >=  (const mpreal& a, const unsigned int b);
551 bool operator >=  (const mpreal& a, const long int b);
552 bool operator >=  (const mpreal& a, const int b);
553 bool operator >=  (const mpreal& a, const long double b);
554 bool operator >=  (const mpreal& a, const double b);
555 
556 bool operator >=  (const unsigned long int a,const mpreal& b);
557 bool operator >=  (const unsigned int a,		const mpreal& b);
558 bool operator >=  (const long int a,			const mpreal& b);
559 bool operator >=  (const int a,				const mpreal& b);
560 bool operator >=  (const long double a,		const mpreal& b);
561 bool operator >=  (const double a,			const mpreal& b);
562 
563 bool operator <=  (const mpreal& a, const unsigned long int b);
564 bool operator <=  (const mpreal& a, const unsigned int b);
565 bool operator <=  (const mpreal& a, const long int b);
566 bool operator <=  (const mpreal& a, const int b);
567 bool operator <=  (const mpreal& a, const long double b);
568 bool operator <=  (const mpreal& a, const double b);
569 
570 bool operator <=  (const unsigned long int a,const mpreal& b);
571 bool operator <=  (const unsigned int a,		const mpreal& b);
572 bool operator <=  (const long int a,			const mpreal& b);
573 bool operator <=  (const int a,				const mpreal& b);
574 bool operator <=  (const long double a,		const mpreal& b);
575 bool operator <=  (const double a,			const mpreal& b);
576 
577 bool operator ==  (const mpreal& a, const unsigned long int b);
578 bool operator ==  (const mpreal& a, const unsigned int b);
579 bool operator ==  (const mpreal& a, const long int b);
580 bool operator ==  (const mpreal& a, const int b);
581 bool operator ==  (const mpreal& a, const long double b);
582 bool operator ==  (const mpreal& a, const double b);
583 
584 bool operator ==  (const unsigned long int a,const mpreal& b);
585 bool operator ==  (const unsigned int a,		const mpreal& b);
586 bool operator ==  (const long int a,			const mpreal& b);
587 bool operator ==  (const int a,				const mpreal& b);
588 bool operator ==  (const long double a,		const mpreal& b);
589 bool operator ==  (const double a,			const mpreal& b);
590 
591 bool operator !=  (const mpreal& a, const unsigned long int b);
592 bool operator !=  (const mpreal& a, const unsigned int b);
593 bool operator !=  (const mpreal& a, const long int b);
594 bool operator !=  (const mpreal& a, const int b);
595 bool operator !=  (const mpreal& a, const long double b);
596 bool operator !=  (const mpreal& a, const double b);
597 
598 bool operator !=  (const unsigned long int a,const mpreal& b);
599 bool operator !=  (const unsigned int a,		const mpreal& b);
600 bool operator !=  (const long int a,			const mpreal& b);
601 bool operator !=  (const int a,				const mpreal& b);
602 bool operator !=  (const long double a,		const mpreal& b);
603 bool operator !=  (const double a,			const mpreal& b);
604 
605 //////////////////////////////////////////////////////////////////////////
606 // sqrt
607 const mpreal sqrt(const unsigned int v, mp_rnd_t rnd_mode = mpreal::default_rnd);
608 const mpreal sqrt(const long int v, mp_rnd_t rnd_mode = mpreal::default_rnd);
609 const mpreal sqrt(const int v, mp_rnd_t rnd_mode = mpreal::default_rnd);
610 const mpreal sqrt(const long double v, mp_rnd_t rnd_mode = mpreal::default_rnd);
611 const mpreal sqrt(const double v, mp_rnd_t rnd_mode = mpreal::default_rnd);
612 
613 //////////////////////////////////////////////////////////////////////////
614 // pow
615 const mpreal pow(const mpreal& a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
616 const mpreal pow(const mpreal& a, const int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
617 const mpreal pow(const mpreal& a, const long double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
618 const mpreal pow(const mpreal& a, const double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
619 
620 const mpreal pow(const unsigned int a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::default_rnd);
621 const mpreal pow(const long int a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::default_rnd);
622 const mpreal pow(const int a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::default_rnd);
623 const mpreal pow(const long double a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::default_rnd);
624 const mpreal pow(const double a, const mpreal& b, mp_rnd_t rnd_mode = mpreal::default_rnd);
625 
626 const mpreal pow(const unsigned long int a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
627 const mpreal pow(const unsigned long int a, const long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
628 const mpreal pow(const unsigned long int a, const int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
629 const mpreal pow(const unsigned long int a, const long double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
630 const mpreal pow(const unsigned long int a, const double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
631 
632 const mpreal pow(const unsigned int a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
633 const mpreal pow(const unsigned int a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
634 const mpreal pow(const unsigned int a, const long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
635 const mpreal pow(const unsigned int a, const int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
636 const mpreal pow(const unsigned int a, const long double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
637 const mpreal pow(const unsigned int a, const double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
638 
639 const mpreal pow(const long int a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
640 const mpreal pow(const long int a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
641 const mpreal pow(const long int a, const long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
642 const mpreal pow(const long int a, const int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
643 const mpreal pow(const long int a, const long double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
644 const mpreal pow(const long int a, const double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
645 
646 const mpreal pow(const int a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
647 const mpreal pow(const int a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
648 const mpreal pow(const int a, const long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
649 const mpreal pow(const int a, const int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
650 const mpreal pow(const int a, const long double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
651 const mpreal pow(const int a, const double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
652 
653 const mpreal pow(const long double a, const long double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
654 const mpreal pow(const long double a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
655 const mpreal pow(const long double a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
656 const mpreal pow(const long double a, const long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
657 const mpreal pow(const long double a, const int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
658 
659 const mpreal pow(const double a, const double b, mp_rnd_t rnd_mode = mpreal::default_rnd);
660 const mpreal pow(const double a, const unsigned long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
661 const mpreal pow(const double a, const unsigned int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
662 const mpreal pow(const double a, const long int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
663 const mpreal pow(const double a, const int b, mp_rnd_t rnd_mode = mpreal::default_rnd);
664 
665 //////////////////////////////////////////////////////////////////////////
666 // Implementation of inline functions
667 //////////////////////////////////////////////////////////////////////////
668 
669 //////////////////////////////////////////////////////////////////////////
670 // Operators - Assignment
671 inline mpreal& mpreal::operator=(const mpreal& v)
672 {
673 	if (this!= &v)	mpfr_set(mp,v.mp,default_rnd);
674 	return *this;
675 }
676 
677 inline mpreal& mpreal::operator=(const mpf_t v)
678 {
679 	mpfr_set_f(mp,v,default_rnd);
680 	return *this;
681 }
682 
683 inline mpreal& mpreal::operator=(const mpz_t v)
684 {
685 	mpfr_set_z(mp,v,default_rnd);
686 	return *this;
687 }
688 
689 inline mpreal& mpreal::operator=(const mpq_t v)
690 {
691 	mpfr_set_q(mp,v,default_rnd);
692 	return *this;
693 }
694 
695 inline mpreal& mpreal::operator=(const long double v)
696 {
697     mpfr_set_ld(mp,v,default_rnd);
698 	return *this;
699 }
700 
701 inline mpreal& mpreal::operator=(const double v)
702 {
703     if(double_bits == -1 || fits_in_bits(v, double_bits))
704     {
705     	mpfr_set_d(mp,v,default_rnd);
706     }
707     else
708         throw conversion_overflow();
709 
710 	return *this;
711 }
712 
713 inline mpreal& mpreal::operator=(const unsigned long int v)
714 {
715 	mpfr_set_ui(mp,v,default_rnd);
716 	return *this;
717 }
718 
719 inline mpreal& mpreal::operator=(const unsigned int v)
720 {
721 	mpfr_set_ui(mp,v,default_rnd);
722 	return *this;
723 }
724 
725 inline mpreal& mpreal::operator=(const long int v)
726 {
727 	mpfr_set_si(mp,v,default_rnd);
728 	return *this;
729 }
730 
731 inline mpreal& mpreal::operator=(const int v)
732 {
733 	mpfr_set_si(mp,v,default_rnd);
734 	return *this;
735 }
736 
737 //////////////////////////////////////////////////////////////////////////
738 // + Addition
739 inline mpreal& mpreal::operator+=(const mpreal& v)
740 {
741 	mpfr_add(mp,mp,v.mp,default_rnd);
742 	return *this;
743 }
744 
745 inline mpreal& mpreal::operator+=(const mpf_t u)
746 {
747 	*this += mpreal(u);
748 	return *this;
749 }
750 
751 inline mpreal& mpreal::operator+=(const mpz_t u)
752 {
753 	mpfr_add_z(mp,mp,u,default_rnd);
754 	return *this;
755 }
756 
757 inline mpreal& mpreal::operator+=(const mpq_t u)
758 {
759 	mpfr_add_q(mp,mp,u,default_rnd);
760 	return *this;
761 }
762 
763 inline mpreal& mpreal::operator+= (const long double u)
764 {
765 	return *this += mpreal(u);
766 }
767 
768 inline mpreal& mpreal::operator+= (const double u)
769 {
770 #if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
771 	mpfr_add_d(mp,mp,u,default_rnd);
772 	return *this;
773 #else
774 	return *this += mpreal(u);
775 #endif
776 }
777 
778 inline mpreal& mpreal::operator+=(const unsigned long int u)
779 {
780 	mpfr_add_ui(mp,mp,u,default_rnd);
781 	return *this;
782 }
783 
784 inline mpreal& mpreal::operator+=(const unsigned int u)
785 {
786 	mpfr_add_ui(mp,mp,u,default_rnd);
787 	return *this;
788 }
789 
790 inline mpreal& mpreal::operator+=(const long int u)
791 {
792 	mpfr_add_si(mp,mp,u,default_rnd);
793 	return *this;
794 }
795 
796 inline mpreal& mpreal::operator+=(const int u)
797 {
798 	mpfr_add_si(mp,mp,u,default_rnd);
799 	return *this;
800 }
801 
802 inline const mpreal mpreal::operator+()const
803 {
804 	return mpreal(*this);
805 }
806 
807 inline const mpreal operator+(const mpreal& a, const mpreal& b)
808 {
809 	// prec(a+b) = max(prec(a),prec(b))
810 	if(a.get_prec()>b.get_prec()) return mpreal(a) += b;
811 	else						  return mpreal(b) += a;
812 }
813 
814 inline const std::string operator+(const mpreal& a, const std::string b)
815 {
816 	return (std::string)a+b;
817 }
818 
819 inline const std::string operator+(const std::string a, const mpreal& b)
820 {
821 	return a+(std::string)b;
822 }
823 
824 inline const mpreal operator+(const mpreal& a, const mpz_t b)
825 {
826 	return mpreal(a) += b;
827 }
828 
829 inline const mpreal operator+(const mpreal& a, const char* b)
830 {
831 	return a+mpreal(b);
832 }
833 
834 inline const mpreal operator+(const char* a, const mpreal& b)
835 {
836 	return mpreal(a)+b;
837 
838 }
839 
840 inline const mpreal operator+(const mpreal& a, const mpq_t b)
841 {
842 	return mpreal(a) += b;
843 }
844 
845 inline const mpreal operator+(const mpreal& a, const long double b)
846 {
847 	return mpreal(a) += b;
848 }
849 
850 inline const mpreal operator+(const mpreal& a, const double b)
851 {
852 	return mpreal(a) += b;
853 }
854 
855 inline const mpreal operator+(const mpreal& a, const unsigned long int b)
856 {
857 	return mpreal(a) += b;
858 }
859 
860 inline const mpreal operator+(const mpreal& a, const unsigned int b)
861 {
862 	return mpreal(a) += b;
863 }
864 
865 inline const mpreal operator+(const mpreal& a, const long int b)
866 {
867 	return mpreal(a) += b;
868 }
869 
870 inline const mpreal operator+(const mpreal& a, const int b)
871 {
872 	return mpreal(a) += b;
873 }
874 
875 inline const mpreal operator+(const mpz_t b, const mpreal& a)
876 {
877 	return mpreal(a) += b;
878 }
879 
880 inline const mpreal operator+(const mpq_t b, const mpreal& a)
881 {
882 	return mpreal(a) += b;
883 }
884 
885 inline const mpreal operator+(const long double b, const mpreal& a)
886 {
887 	return mpreal(a) += b;
888 }
889 
890 inline const mpreal operator+(const double  b, const mpreal& a)
891 {
892 	return mpreal(a) += b;
893 }
894 
895 inline const mpreal operator+(const unsigned long int b, const mpreal& a)
896 {
897 	return mpreal(a) += b;
898 }
899 
900 inline const mpreal operator+(const unsigned int b, const mpreal& a)
901 {
902 	return mpreal(a) += b;
903 }
904 
905 inline const mpreal operator+(const long int b, const mpreal& a)
906 {
907 	return mpreal(a) += b;
908 }
909 
910 inline const mpreal operator+(const int b, const mpreal& a)
911 {
912 	return mpreal(a) += b;
913 }
914 
915 inline mpreal& mpreal::operator++()
916 {
917 	*this += 1;
918 	return *this;
919 }
920 
921 inline const mpreal mpreal::operator++ (int)
922 {
923 	mpreal x(*this);
924 	*this += 1;
925 	return x;
926 }
927 
928 inline mpreal& mpreal::operator--()
929 {
930 	*this -= 1;
931 	return *this;
932 }
933 
934 inline const mpreal mpreal::operator-- (int)
935 {
936 	mpreal x(*this);
937 	*this -= 1;
938 	return x;
939 }
940 
941 //////////////////////////////////////////////////////////////////////////
942 // - Subtraction
943 inline mpreal& mpreal::operator-= (const mpreal& v)
944 {
945 	mpfr_sub(mp,mp,v.mp,default_rnd);
946 	return *this;
947 }
948 
949 inline mpreal& mpreal::operator-=(const mpz_t v)
950 {
951 	mpfr_sub_z(mp,mp,v,default_rnd);
952 	return *this;
953 }
954 
955 inline mpreal& mpreal::operator-=(const mpq_t v)
956 {
957 	mpfr_sub_q(mp,mp,v,default_rnd);
958 	return *this;
959 }
960 
961 inline mpreal& mpreal::operator-=(const long double v)
962 {
963 	return *this -= mpreal(v);
964 }
965 
966 inline mpreal& mpreal::operator-=(const double v)
967 {
968 #if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
969 	mpfr_sub_d(mp,mp,v,default_rnd);
970 	return *this;
971 #else
972 	return *this -= mpreal(v);
973 #endif
974 }
975 
976 inline mpreal& mpreal::operator-=(const unsigned long int v)
977 {
978 	mpfr_sub_ui(mp,mp,v,default_rnd);
979 	return *this;
980 }
981 
982 inline mpreal& mpreal::operator-=(const unsigned int v)
983 {
984 	mpfr_sub_ui(mp,mp,v,default_rnd);
985 	return *this;
986 }
987 
988 inline mpreal& mpreal::operator-=(const long int v)
989 {
990 	mpfr_sub_si(mp,mp,v,default_rnd);
991 	return *this;
992 }
993 
994 inline mpreal& mpreal::operator-=(const int v)
995 {
996 	mpfr_sub_si(mp,mp,v,default_rnd);
997 	return *this;
998 }
999 
1000 inline const mpreal mpreal::operator-()const
1001 {
1002 	mpreal u(*this);
1003 	mpfr_neg(u.mp,u.mp,default_rnd);
1004 	return u;
1005 }
1006 
1007 inline const mpreal operator-(const mpreal& a, const mpreal& b)
1008 {
1009 	// prec(a-b) = max(prec(a),prec(b))
1010 	if(a.get_prec()>b.get_prec())	return   mpreal(a) -= b;
1011 	else							return -(mpreal(b) -= a);
1012 }
1013 
1014 inline const mpreal operator-(const mpreal& a, const mpz_t b)
1015 {
1016 	return mpreal(a) -= b;
1017 }
1018 
1019 inline const mpreal operator-(const mpreal& a, const mpq_t b)
1020 {
1021 	return mpreal(a) -= b;
1022 }
1023 
1024 inline const mpreal operator-(const mpreal& a, const long double b)
1025 {
1026 	return mpreal(a) -= b;
1027 }
1028 
1029 inline const mpreal operator-(const mpreal& a, const double b)
1030 {
1031 	return mpreal(a) -= b;
1032 }
1033 
1034 inline const mpreal operator-(const mpreal& a, const unsigned long int b)
1035 {
1036 	return mpreal(a) -= b;
1037 }
1038 
1039 inline const mpreal operator-(const mpreal& a, const unsigned int b)
1040 {
1041 	return mpreal(a) -= b;
1042 }
1043 
1044 inline const mpreal operator-(const mpreal& a, const long int b)
1045 {
1046 	return mpreal(a) -= b;
1047 }
1048 
1049 inline const mpreal operator-(const mpreal& a, const int b)
1050 {
1051 	return mpreal(a) -= b;
1052 }
1053 
1054 inline const mpreal operator-(const mpz_t b, const mpreal& a)
1055 {
1056 	return -(mpreal(a) -= b);
1057 }
1058 
1059 inline const mpreal operator-(const mpq_t b, const mpreal& a)
1060 {
1061 	return -(mpreal(a) -= b);
1062 }
1063 
1064 inline const mpreal operator-(const long double b, const mpreal& a)
1065 {
1066 	return -(mpreal(a) -= b);
1067 }
1068 
1069 inline const mpreal operator-(const double  b, const mpreal& a)
1070 {
1071 #if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1072 	mpreal x(a);
1073 	mpfr_d_sub(x.mp,b,a.mp,mpreal::default_rnd);
1074 	return x;
1075 #else
1076 	return -(mpreal(a) -= b);
1077 #endif
1078 }
1079 
1080 inline const mpreal operator-(const unsigned long int b, const mpreal& a)
1081 {
1082 	mpreal x(a);
1083 	mpfr_ui_sub(x.mp,b,a.mp,mpreal::default_rnd);
1084 	return x;
1085 }
1086 
1087 inline const mpreal operator-(const unsigned int b, const mpreal& a)
1088 {
1089 	mpreal x(a);
1090 	mpfr_ui_sub(x.mp,b,a.mp,mpreal::default_rnd);
1091 	return x;
1092 }
1093 
1094 inline const mpreal operator-(const long int b, const mpreal& a)
1095 {
1096 	mpreal x(a);
1097 	mpfr_si_sub(x.mp,b,a.mp,mpreal::default_rnd);
1098 	return x;
1099 }
1100 
1101 inline const mpreal operator-(const int b, const mpreal& a)
1102 {
1103 	mpreal x(a);
1104 	mpfr_si_sub(x.mp,b,a.mp,mpreal::default_rnd);
1105 	return x;
1106 }
1107 
1108 inline const mpreal operator-(const mpreal& a, const char* b)
1109 {
1110 	return a-mpreal(b);
1111 }
1112 
1113 inline const mpreal operator-(const char* a, const mpreal& b)
1114 {
1115 	return mpreal(a)-b;
1116 }
1117 
1118 //////////////////////////////////////////////////////////////////////////
1119 // * Multiplication
1120 inline mpreal& mpreal::operator*= (const mpreal& v)
1121 {
1122 	mpfr_mul(mp,mp,v.mp,default_rnd);
1123 	return *this;
1124 }
1125 
1126 inline mpreal& mpreal::operator*=(const mpz_t v)
1127 {
1128 	mpfr_mul_z(mp,mp,v,default_rnd);
1129 	return *this;
1130 }
1131 
1132 inline mpreal& mpreal::operator*=(const mpq_t v)
1133 {
1134 	mpfr_mul_q(mp,mp,v,default_rnd);
1135 	return *this;
1136 }
1137 
1138 inline mpreal& mpreal::operator*=(const long double v)
1139 {
1140 	return *this *= mpreal(v);
1141 }
1142 
1143 inline mpreal& mpreal::operator*=(const double v)
1144 {
1145 #if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1146 	mpfr_mul_d(mp,mp,v,default_rnd);
1147 	return *this;
1148 #else
1149 	return *this *= mpreal(v);
1150 #endif
1151 }
1152 
1153 inline mpreal& mpreal::operator*=(const unsigned long int v)
1154 {
1155 	mpfr_mul_ui(mp,mp,v,default_rnd);
1156 	return *this;
1157 }
1158 
1159 inline mpreal& mpreal::operator*=(const unsigned int v)
1160 {
1161 	mpfr_mul_ui(mp,mp,v,default_rnd);
1162 	return *this;
1163 }
1164 
1165 inline mpreal& mpreal::operator*=(const long int v)
1166 {
1167 	mpfr_mul_si(mp,mp,v,default_rnd);
1168 	return *this;
1169 }
1170 
1171 inline mpreal& mpreal::operator*=(const int v)
1172 {
1173 	mpfr_mul_si(mp,mp,v,default_rnd);
1174 	return *this;
1175 }
1176 
1177 inline const mpreal operator*(const mpreal& a, const mpreal& b)
1178 {
1179 	// prec(a*b) = max(prec(a),prec(b))
1180 	if(a.get_prec()>b.get_prec())	return   mpreal(a) *= b;
1181 	else							return   mpreal(b) *= a;
1182 }
1183 
1184 inline const mpreal operator*(const mpreal& a, const mpz_t b)
1185 {
1186 	return mpreal(a) *= b;
1187 }
1188 
1189 inline const mpreal operator*(const mpreal& a, const mpq_t b)
1190 {
1191 	return mpreal(a) *= b;
1192 }
1193 
1194 inline const mpreal operator*(const mpreal& a, const long double b)
1195 {
1196 	return mpreal(a) *= b;
1197 }
1198 
1199 inline const mpreal operator*(const mpreal& a, const double b)
1200 {
1201 	return mpreal(a) *= b;
1202 }
1203 
1204 inline const mpreal operator*(const mpreal& a, const unsigned long int b)
1205 {
1206 	return mpreal(a) *= b;
1207 }
1208 
1209 inline const mpreal operator*(const mpreal& a, const unsigned int b)
1210 {
1211 	return mpreal(a) *= b;
1212 }
1213 
1214 inline const mpreal operator*(const mpreal& a, const long int b)
1215 {
1216 	return mpreal(a) *= b;
1217 }
1218 
1219 inline const mpreal operator*(const mpreal& a, const int b)
1220 {
1221 	return mpreal(a) *= b;
1222 }
1223 
1224 inline const mpreal operator*(const mpz_t b, const mpreal& a)
1225 {
1226 	return mpreal(a) *= b;
1227 }
1228 
1229 inline const mpreal operator*(const mpq_t b, const mpreal& a)
1230 {
1231 	return mpreal(a) *= b;
1232 }
1233 
1234 inline const mpreal operator*(const long double b, const mpreal& a)
1235 {
1236 	return mpreal(a) *= b;
1237 }
1238 
1239 inline const mpreal operator*(const double  b, const mpreal& a)
1240 {
1241 	return mpreal(a) *= b;
1242 }
1243 
1244 inline const mpreal operator*(const unsigned long int b, const mpreal& a)
1245 {
1246 	return mpreal(a) *= b;
1247 }
1248 
1249 inline const mpreal operator*(const unsigned int b, const mpreal& a)
1250 {
1251 	return mpreal(a) *= b;
1252 }
1253 
1254 inline const mpreal operator*(const long int b, const mpreal& a)
1255 {
1256 	return mpreal(a) *= b;
1257 }
1258 
1259 inline const mpreal operator*(const int b, const mpreal& a)
1260 {
1261 	return mpreal(a) *= b;
1262 }
1263 
1264 //////////////////////////////////////////////////////////////////////////
1265 // / Division
1266 inline mpreal& mpreal::operator/=(const mpreal& v)
1267 {
1268 	mpfr_div(mp,mp,v.mp,default_rnd);
1269 	return *this;
1270 }
1271 
1272 inline mpreal& mpreal::operator/=(const mpz_t v)
1273 {
1274 	mpfr_div_z(mp,mp,v,default_rnd);
1275 	return *this;
1276 }
1277 
1278 inline mpreal& mpreal::operator/=(const mpq_t v)
1279 {
1280 	mpfr_div_q(mp,mp,v,default_rnd);
1281 	return *this;
1282 }
1283 
1284 inline mpreal& mpreal::operator/=(const long double v)
1285 {
1286 	return *this /= mpreal(v);
1287 }
1288 
1289 inline mpreal& mpreal::operator/=(const double v)
1290 {
1291 #if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1292 	mpfr_div_d(mp,mp,v,default_rnd);
1293 	return *this;
1294 #else
1295 	return *this /= mpreal(v);
1296 #endif
1297 }
1298 
1299 inline mpreal& mpreal::operator/=(const unsigned long int v)
1300 {
1301 	mpfr_div_ui(mp,mp,v,default_rnd);
1302 	return *this;
1303 }
1304 
1305 inline mpreal& mpreal::operator/=(const unsigned int v)
1306 {
1307 	mpfr_div_ui(mp,mp,v,default_rnd);
1308 	return *this;
1309 }
1310 
1311 inline mpreal& mpreal::operator/=(const long int v)
1312 {
1313 	mpfr_div_si(mp,mp,v,default_rnd);
1314 	return *this;
1315 }
1316 
1317 inline mpreal& mpreal::operator/=(const int v)
1318 {
1319 	mpfr_div_si(mp,mp,v,default_rnd);
1320 	return *this;
1321 }
1322 
1323 inline const mpreal operator/(const mpreal& a, const mpreal& b)
1324 {
1325 	mpreal x(a);
1326 	mp_prec_t pb;
1327 	mp_prec_t pa;
1328 
1329 	// prec(a/b) = max(prec(a),prec(b))
1330 	pa = a.get_prec();
1331 	pb = b.get_prec();
1332 	if(pb>pa) x.set_prec(pb);
1333 
1334 	return   x /= b;
1335 }
1336 
1337 inline const mpreal operator/(const mpreal& a, const mpz_t b)
1338 {
1339 	return mpreal(a) /= b;
1340 }
1341 
1342 inline const mpreal operator/(const mpreal& a, const mpq_t b)
1343 {
1344 	return mpreal(a) /= b;
1345 }
1346 
1347 inline const mpreal operator/(const mpreal& a, const long double b)
1348 {
1349 	return mpreal(a) /= b;
1350 }
1351 
1352 inline const mpreal operator/(const mpreal& a, const double b)
1353 {
1354 	return mpreal(a) /= b;
1355 }
1356 
1357 inline const mpreal operator/(const mpreal& a, const unsigned long int b)
1358 {
1359 	return mpreal(a) /= b;
1360 }
1361 
1362 inline const mpreal operator/(const mpreal& a, const unsigned int b)
1363 {
1364 	return mpreal(a) /= b;
1365 }
1366 
1367 inline const mpreal operator/(const mpreal& a, const long int b)
1368 {
1369 	return mpreal(a) /= b;
1370 }
1371 
1372 inline const mpreal operator/(const mpreal& a, const int b)
1373 {
1374 	return mpreal(a) /= b;
1375 }
1376 
1377 inline const mpreal operator/(const unsigned long int b, const mpreal& a)
1378 {
1379 	mpreal x(a);
1380 	mpfr_ui_div(x.mp,b,a.mp,mpreal::default_rnd);
1381 	return x;
1382 }
1383 
1384 inline const mpreal operator/(const unsigned int b, const mpreal& a)
1385 {
1386 	mpreal x(a);
1387 	mpfr_ui_div(x.mp,b,a.mp,mpreal::default_rnd);
1388 	return x;
1389 }
1390 
1391 inline const mpreal operator/(const long int b, const mpreal& a)
1392 {
1393 	mpreal x(a);
1394 	mpfr_si_div(x.mp,b,a.mp,mpreal::default_rnd);
1395 	return x;
1396 }
1397 
1398 inline const mpreal operator/(const int b, const mpreal& a)
1399 {
1400 	mpreal x(a);
1401 	mpfr_si_div(x.mp,b,a.mp,mpreal::default_rnd);
1402 	return x;
1403 }
1404 
1405 inline const mpreal operator/(const long double b, const mpreal& a)
1406 {
1407 	mpreal x(b);
1408 	return x/a;
1409 }
1410 
1411 inline const mpreal operator/(const double  b, const mpreal& a)
1412 {
1413 #if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
1414 	mpreal x(a);
1415 	mpfr_d_div(x.mp,b,a.mp,mpreal::default_rnd);
1416 	return x;
1417 #else
1418 	mpreal x(b);
1419 	return x/a;
1420 #endif
1421 }
1422 
1423 //////////////////////////////////////////////////////////////////////////
1424 // Shifts operators - Multiplication/Division by power of 2
1425 inline mpreal& mpreal::operator<<=(const unsigned long int u)
1426 {
1427 	mpfr_mul_2ui(mp,mp,u,default_rnd);
1428 	return *this;
1429 }
1430 
1431 inline mpreal& mpreal::operator<<=(const unsigned int u)
1432 {
1433 	mpfr_mul_2ui(mp,mp,static_cast<unsigned long int>(u),default_rnd);
1434 	return *this;
1435 }
1436 
1437 inline mpreal& mpreal::operator<<=(const long int u)
1438 {
1439 	mpfr_mul_2si(mp,mp,u,default_rnd);
1440 	return *this;
1441 }
1442 
1443 inline mpreal& mpreal::operator<<=(const int u)
1444 {
1445 	mpfr_mul_2si(mp,mp,static_cast<long int>(u),default_rnd);
1446 	return *this;
1447 }
1448 
1449 inline mpreal& mpreal::operator>>=(const unsigned long int u)
1450 {
1451 	mpfr_div_2ui(mp,mp,u,default_rnd);
1452 	return *this;
1453 }
1454 
1455 inline mpreal& mpreal::operator>>=(const unsigned int u)
1456 {
1457 	mpfr_div_2ui(mp,mp,static_cast<unsigned long int>(u),default_rnd);
1458 	return *this;
1459 }
1460 
1461 inline mpreal& mpreal::operator>>=(const long int u)
1462 {
1463 	mpfr_div_2si(mp,mp,u,default_rnd);
1464 	return *this;
1465 }
1466 
1467 inline mpreal& mpreal::operator>>=(const int u)
1468 {
1469 	mpfr_div_2si(mp,mp,static_cast<long int>(u),default_rnd);
1470 	return *this;
1471 }
1472 
1473 inline const mpreal operator<<(const mpreal& v, const unsigned long int k)
1474 {
1475 	return mul_2ui(v,k);
1476 }
1477 
1478 inline const mpreal operator<<(const mpreal& v, const unsigned int k)
1479 {
1480 	return mul_2ui(v,static_cast<unsigned long int>(k));
1481 }
1482 
1483 inline const mpreal operator<<(const mpreal& v, const long int k)
1484 {
1485 	return mul_2si(v,k);
1486 }
1487 
1488 inline const mpreal operator<<(const mpreal& v, const int k)
1489 {
1490 	return mul_2si(v,static_cast<long int>(k));
1491 }
1492 
1493 inline const mpreal operator>>(const mpreal& v, const unsigned long int k)
1494 {
1495 	return div_2ui(v,k);
1496 }
1497 
1498 inline const mpreal operator>>(const mpreal& v, const long int k)
1499 {
1500 	return div_2si(v,k);
1501 }
1502 
1503 inline const mpreal operator>>(const mpreal& v, const unsigned int k)
1504 {
1505 	return div_2ui(v,static_cast<unsigned long int>(k));
1506 }
1507 
1508 inline const mpreal operator>>(const mpreal& v, const int k)
1509 {
1510 	return div_2si(v,static_cast<long int>(k));
1511 }
1512 
1513 // mul_2ui
mul_2ui(const mpreal & v,unsigned long int k,mp_rnd_t rnd_mode)1514 inline const mpreal mul_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode)
1515 {
1516 	mpreal x(v);
1517 	mpfr_mul_2ui(x.mp,v.mp,k,rnd_mode);
1518 	return x;
1519 }
1520 
1521 // mul_2si
mul_2si(const mpreal & v,long int k,mp_rnd_t rnd_mode)1522 inline const mpreal mul_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode)
1523 {
1524 	mpreal x(v);
1525 	mpfr_mul_2si(x.mp,v.mp,k,rnd_mode);
1526 	return x;
1527 }
1528 
div_2ui(const mpreal & v,unsigned long int k,mp_rnd_t rnd_mode)1529 inline const mpreal div_2ui(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode)
1530 {
1531 	mpreal x(v);
1532 	mpfr_div_2ui(x.mp,v.mp,k,rnd_mode);
1533 	return x;
1534 }
1535 
div_2si(const mpreal & v,long int k,mp_rnd_t rnd_mode)1536 inline const mpreal div_2si(const mpreal& v, long int k, mp_rnd_t rnd_mode)
1537 {
1538 	mpreal x(v);
1539 	mpfr_div_2si(x.mp,v.mp,k,rnd_mode);
1540 	return x;
1541 }
1542 
1543 //////////////////////////////////////////////////////////////////////////
1544 //Boolean operators
1545 inline bool operator > (const mpreal& a, const mpreal& b)
1546 {
1547 	return (mpfr_greater_p(a.mp,b.mp)!=0);
1548 }
1549 
1550 inline bool operator >  (const mpreal& a, const unsigned long int b)
1551 {
1552 	return a>mpreal(b);
1553 }
1554 
1555 inline bool operator >  (const mpreal& a, const unsigned int b)
1556 {
1557 	return a>mpreal(b);
1558 }
1559 
1560 inline bool operator >  (const mpreal& a, const long int b)
1561 {
1562 	return a>mpreal(b);
1563 }
1564 
1565 inline bool operator >  (const mpreal& a, const int b)
1566 {
1567 	return a>mpreal(b);
1568 }
1569 
1570 inline bool operator >  (const mpreal& a, const long double b)
1571 {
1572 	return a>mpreal(b);
1573 }
1574 
1575 inline bool operator >  (const mpreal& a, const double b)
1576 {
1577 	return a>mpreal(b);
1578 }
1579 
1580 inline bool operator >  (const unsigned long int a,	const mpreal& b)
1581 {
1582 	return mpreal(a)>b;
1583 }
1584 
1585 inline bool operator >  (const unsigned int a,		const mpreal& b)
1586 {
1587 	return mpreal(a)>b;
1588 }
1589 
1590 inline bool operator >  (const long int a,			const mpreal& b)
1591 {
1592 	return mpreal(a)>b;
1593 }
1594 
1595 inline bool operator >  (const int a,				const mpreal& b)
1596 {
1597 	return mpreal(a)>b;
1598 }
1599 
1600 inline bool operator >  (const long double a,		const mpreal& b)
1601 {
1602 	return mpreal(a)>b;
1603 }
1604 
1605 inline bool operator >  (const double a,			const mpreal& b)
1606 {
1607 	return mpreal(a)>b;
1608 }
1609 
1610 inline bool operator >= (const mpreal& a, const mpreal& b)
1611 {
1612 	return (mpfr_greaterequal_p(a.mp,b.mp)!=0);
1613 }
1614 
1615 inline bool operator >=  (const mpreal& a, const unsigned long int b)
1616 {
1617 	return a>=mpreal(b);
1618 }
1619 
1620 inline bool operator >=  (const mpreal& a, const unsigned int b)
1621 {
1622 	return a>=mpreal(b);
1623 }
1624 
1625 inline bool operator >=  (const mpreal& a, const long int b)
1626 {
1627 	return a>=mpreal(b);
1628 }
1629 
1630 inline bool operator >=  (const mpreal& a, const int b)
1631 {
1632 	return a>=mpreal(b);
1633 }
1634 
1635 inline bool operator >=  (const mpreal& a, const long double b)
1636 {
1637 	return a>=mpreal(b);
1638 }
1639 
1640 inline bool operator >=  (const mpreal& a, const double b)
1641 {
1642 	return a>=mpreal(b);
1643 }
1644 
1645 inline bool operator >=  (const unsigned long int a,const mpreal& b)
1646 {
1647 	return mpreal(a)>=b;
1648 }
1649 
1650 inline bool operator >=  (const unsigned int a,		const mpreal& b)
1651 {
1652 	return mpreal(a)>=b;
1653 }
1654 
1655 inline bool operator >=  (const long int a,			const mpreal& b)
1656 {
1657 	return mpreal(a)>=b;
1658 }
1659 
1660 inline bool operator >=  (const int a,				const mpreal& b)
1661 {
1662 	return mpreal(a)>=b;
1663 }
1664 
1665 inline bool operator >=  (const long double a,		const mpreal& b)
1666 {
1667 	return mpreal(a)>=b;
1668 }
1669 
1670 inline bool operator >=  (const double a,			const mpreal& b)
1671 {
1672 	return mpreal(a)>=b;
1673 }
1674 
1675 inline bool operator <  (const mpreal& a, const mpreal& b)
1676 {
1677 	return (mpfr_less_p(a.mp,b.mp)!=0);
1678 }
1679 
1680 inline bool operator <  (const mpreal& a, const unsigned long int b)
1681 {
1682 	return a<mpreal(b);
1683 }
1684 
1685 inline bool operator <  (const mpreal& a, const unsigned int b)
1686 {
1687 	return a<mpreal(b);
1688 }
1689 
1690 inline bool operator <  (const mpreal& a, const long int b)
1691 {
1692 	return a<mpreal(b);
1693 }
1694 
1695 inline bool operator <  (const mpreal& a, const int b)
1696 {
1697 	return a<mpreal(b);
1698 }
1699 
1700 inline bool operator <  (const mpreal& a, const long double b)
1701 {
1702 	return a<mpreal(b);
1703 }
1704 
1705 inline bool operator <  (const mpreal& a, const double b)
1706 {
1707 	return a<mpreal(b);
1708 }
1709 
1710 inline bool operator <  (const unsigned long int a,	const mpreal& b)
1711 {
1712 	return mpreal(a)<b;
1713 }
1714 
1715 inline bool operator <  (const unsigned int a,const mpreal& b)
1716 {
1717 	return mpreal(a)<b;
1718 }
1719 
1720 inline bool operator <  (const long int a,const mpreal& b)
1721 {
1722 	return mpreal(a)<b;
1723 }
1724 
1725 inline bool operator <  (const int a,const mpreal& b)
1726 {
1727 	return mpreal(a)<b;
1728 }
1729 
1730 inline bool operator <  (const long double a,const mpreal& b)
1731 {
1732 	return mpreal(a)<b;
1733 }
1734 
1735 inline bool operator <  (const double a,const mpreal& b)
1736 {
1737 	return mpreal(a)<b;
1738 }
1739 
1740 inline bool operator <= (const mpreal& a, const mpreal& b)
1741 {
1742 	return (mpfr_lessequal_p(a.mp,b.mp)!=0);
1743 }
1744 
1745 inline bool operator <=  (const mpreal& a, const unsigned long int b)
1746 {
1747 	return a<=mpreal(b);
1748 }
1749 
1750 inline bool operator <=  (const mpreal& a, const unsigned int b)
1751 {
1752 	return a<=mpreal(b);
1753 }
1754 
1755 inline bool operator <=  (const mpreal& a, const long int b)
1756 {
1757 	return a<=mpreal(b);
1758 }
1759 
1760 inline bool operator <=  (const mpreal& a, const int b)
1761 {
1762 	return a<=mpreal(b);
1763 }
1764 
1765 inline bool operator <=  (const mpreal& a, const long double b)
1766 {
1767 	return a<=mpreal(b);
1768 }
1769 
1770 inline bool operator <=  (const mpreal& a, const double b)
1771 {
1772 	return a<=mpreal(b);
1773 }
1774 
1775 inline bool operator <=  (const unsigned long int a,const mpreal& b)
1776 {
1777 	return mpreal(a)<=b;
1778 }
1779 
1780 inline bool operator <=  (const unsigned int a,		const mpreal& b)
1781 {
1782 	return mpreal(a)<=b;
1783 }
1784 
1785 inline bool operator <=  (const long int a,			const mpreal& b)
1786 {
1787 	return mpreal(a)<=b;
1788 }
1789 
1790 inline bool operator <=  (const int a,				const mpreal& b)
1791 {
1792 	return mpreal(a)<=b;
1793 }
1794 
1795 inline bool operator <=  (const long double a,		const mpreal& b)
1796 {
1797 	return mpreal(a)<=b;
1798 }
1799 
1800 inline bool operator <=  (const double a,			const mpreal& b)
1801 {
1802 	return mpreal(a)<=b;
1803 }
1804 
1805 inline bool operator == (const mpreal& a, const mpreal& b)
1806 {
1807 	return (mpfr_equal_p(a.mp,b.mp)!=0);
1808 }
1809 
1810 inline bool operator ==  (const mpreal& a, const unsigned long int b)
1811 {
1812 	return a==mpreal(b);
1813 }
1814 
1815 inline bool operator ==  (const mpreal& a, const unsigned int b)
1816 {
1817 	return a==mpreal(b);
1818 }
1819 
1820 inline bool operator ==  (const mpreal& a, const long int b)
1821 {
1822 	return a==mpreal(b);
1823 }
1824 
1825 inline bool operator ==  (const mpreal& a, const int b)
1826 {
1827 	return a==mpreal(b);
1828 }
1829 
1830 inline bool operator ==  (const mpreal& a, const long double b)
1831 {
1832 	return a==mpreal(b);
1833 }
1834 
1835 inline bool operator ==  (const mpreal& a, const double b)
1836 {
1837 	return a==mpreal(b);
1838 }
1839 
1840 inline bool operator ==  (const unsigned long int a,const mpreal& b)
1841 {
1842 	return mpreal(a)==b;
1843 }
1844 
1845 inline bool operator ==  (const unsigned int a,		const mpreal& b)
1846 {
1847 	return mpreal(a)==b;
1848 }
1849 
1850 inline bool operator ==  (const long int a,			const mpreal& b)
1851 {
1852 	return mpreal(a)==b;
1853 }
1854 
1855 inline bool operator ==  (const int a,				const mpreal& b)
1856 {
1857 	return mpreal(a)==b;
1858 }
1859 
1860 inline bool operator ==  (const long double a,		const mpreal& b)
1861 {
1862 	return mpreal(a)==b;
1863 }
1864 
1865 inline bool operator ==  (const double a,			const mpreal& b)
1866 {
1867 	return mpreal(a)==b;
1868 }
1869 
1870 inline bool operator != (const mpreal& a, const mpreal& b)
1871 {
1872 	return (mpfr_lessgreater_p(a.mp,b.mp)!=0);
1873 }
1874 
1875 inline bool operator !=  (const mpreal& a, const unsigned long int b)
1876 {
1877 	return a!=mpreal(b);
1878 }
1879 
1880 inline bool operator !=  (const mpreal& a, const unsigned int b)
1881 {
1882 	return a!=mpreal(b);
1883 }
1884 
1885 inline bool operator !=  (const mpreal& a, const long int b)
1886 {
1887 	return a!=mpreal(b);
1888 }
1889 
1890 inline bool operator !=  (const mpreal& a, const int b)
1891 {
1892 	return a!=mpreal(b);
1893 }
1894 
1895 inline bool operator !=  (const mpreal& a, const long double b)
1896 {
1897 	return a!=mpreal(b);
1898 }
1899 
1900 inline bool operator !=  (const mpreal& a, const double b)
1901 {
1902 	return a!=mpreal(b);
1903 }
1904 
1905 inline bool operator !=  (const unsigned long int a,const mpreal& b)
1906 {
1907 	return mpreal(a)!=b;
1908 }
1909 
1910 inline bool operator !=  (const unsigned int a,		const mpreal& b)
1911 {
1912 	return mpreal(a)!=b;
1913 }
1914 
1915 inline bool operator !=  (const long int a,			const mpreal& b)
1916 {
1917 	return mpreal(a)!=b;
1918 }
1919 
1920 inline bool operator !=  (const int a,				const mpreal& b)
1921 {
1922 	return mpreal(a)!=b;
1923 }
1924 
1925 inline bool operator !=  (const long double a,		const mpreal& b)
1926 {
1927 	return mpreal(a)!=b;
1928 }
1929 
1930 inline bool operator !=  (const double a,			const mpreal& b)
1931 {
1932 	return mpreal(a)!=b;
1933 }
1934 
_isnan(const mpreal & v)1935 inline bool _isnan(const mpreal& v)
1936 {
1937 	return (mpfr_nan_p(v.mp)!=0);
1938 }
1939 
_isinf(const mpreal & v)1940 inline bool _isinf(const mpreal& v)
1941 {
1942 	return (mpfr_inf_p(v.mp)!=0);
1943 }
1944 
_isnum(const mpreal & v)1945 inline bool _isnum(const mpreal& v)
1946 {
1947 	return (mpfr_number_p(v.mp)!=0);
1948 }
1949 
_iszero(const mpreal & v)1950 inline bool _iszero(const mpreal& v)
1951 {
1952 	return (mpfr_zero_p(v.mp)!=0);
1953 }
1954 
_isint(const mpreal & v)1955 inline bool _isint(const mpreal& v)
1956 {
1957 	return (mpfr_integer_p(v.mp)!=0);
1958 }
1959 
1960 //////////////////////////////////////////////////////////////////////////
1961 // Type Converters
1962 inline mpreal::operator double() const
1963 {
1964 	return mpfr_get_d(mp,default_rnd);
1965 }
1966 
1967 inline mpreal::operator float() const
1968 {
1969 	return (float)mpfr_get_d(mp,default_rnd);
1970 }
1971 
1972 inline mpreal::operator long double() const
1973 {
1974 	return mpfr_get_ld(mp,default_rnd);
1975 }
1976 
1977 inline mpreal::operator unsigned long() const
1978 {
1979 	return mpfr_get_ui(mp,default_rnd);
1980 }
1981 
1982 inline mpreal::operator unsigned int() const
1983 {
1984 	return mpfr_get_ui(mp,default_rnd);
1985 }
1986 
1987 inline mpreal::operator long() const
1988 {
1989 	return mpfr_get_si(mp,default_rnd);
1990 }
1991 
string()1992 inline mpreal::operator std::string() const
1993 {
1994 	return to_string();
1995 }
1996 
mpfr_ptr()1997 inline mpreal::operator mpfr_ptr()
1998 {
1999 	return mp;
2000 }
2001 
2002 //////////////////////////////////////////////////////////////////////////
2003 // Set/Get number properties
sgn(const mpreal & v)2004 inline int sgn(const mpreal& v)
2005 {
2006 	int r = mpfr_signbit(v.mp);
2007 	return (r>0?-1:1);
2008 }
2009 
set_sign(int sign,mp_rnd_t rnd_mode)2010 inline void mpreal::set_sign(int sign, mp_rnd_t rnd_mode)
2011 {
2012 	mpfr_setsign(mp,mp,(sign<0?1:0),rnd_mode);
2013 }
2014 
get_prec()2015 inline mp_prec_t mpreal::get_prec() const
2016 {
2017 	return mpfr_get_prec(mp);
2018 }
2019 
set_prec(mp_prec_t prec,mp_rnd_t rnd_mode)2020 inline void mpreal::set_prec(mp_prec_t prec, mp_rnd_t rnd_mode)
2021 {
2022 	mpfr_prec_round(mp,prec,rnd_mode);
2023 }
2024 
set_inf(int sign)2025 inline void mpreal::set_inf(int sign)
2026 {
2027 	mpfr_set_inf(mp,sign);
2028 }
2029 
set_nan()2030 inline void mpreal::set_nan()
2031 {
2032 	mpfr_set_nan(mp);
2033 }
2034 
get_exp()2035 inline mp_exp_t mpreal::get_exp ()
2036 {
2037 	return mpfr_get_exp(mp);
2038 }
2039 
set_exp(mp_exp_t e)2040 inline int mpreal::set_exp (mp_exp_t e)
2041 {
2042 	return mpfr_set_exp(mp,e);
2043 }
2044 
frexp(const mpreal & v,mp_exp_t * exp)2045 inline const mpreal frexp(const mpreal& v, mp_exp_t* exp)
2046 {
2047 	mpreal x(v);
2048 	*exp = x.get_exp();
2049 	x.set_exp(0);
2050 	return x;
2051 }
2052 
ldexp(const mpreal & v,mp_exp_t exp)2053 inline const mpreal ldexp(const mpreal& v, mp_exp_t exp)
2054 {
2055 	mpreal x(v);
2056 
2057 	// rounding is not important since we just increasing the exponent
2058 	mpfr_mul_2si(x.mp,x.mp,exp,mpreal::default_rnd);
2059 	return x;
2060 }
2061 
modf(const mpreal & v,mpreal & n)2062 inline const mpreal modf(const mpreal& v, mpreal& n)
2063 {
2064 	mpreal frac(v);
2065 
2066 	// rounding is not important since we are using the same number
2067 	mpfr_frac(frac.mp,frac.mp,mpreal::default_rnd);
2068 	mpfr_trunc(n.mp,v.mp);
2069 	return frac;
2070 }
2071 
check_range(int t,mp_rnd_t rnd_mode)2072 inline int mpreal::check_range (int t, mp_rnd_t rnd_mode)
2073 {
2074 	return mpfr_check_range(mp,t,rnd_mode);
2075 }
2076 
subnormalize(int t,mp_rnd_t rnd_mode)2077 inline int mpreal::subnormalize (int t,mp_rnd_t rnd_mode)
2078 {
2079 	return mpfr_subnormalize(mp,t,rnd_mode);
2080 }
2081 
get_emin(void)2082 inline mp_exp_t mpreal::get_emin (void)
2083 {
2084 	return mpfr_get_emin();
2085 }
2086 
set_emin(mp_exp_t exp)2087 inline int mpreal::set_emin (mp_exp_t exp)
2088 {
2089 	return mpfr_set_emin(exp);
2090 }
2091 
get_emax(void)2092 inline mp_exp_t mpreal::get_emax (void)
2093 {
2094 	return mpfr_get_emax();
2095 }
2096 
set_emax(mp_exp_t exp)2097 inline int mpreal::set_emax (mp_exp_t exp)
2098 {
2099 	return mpfr_set_emax(exp);
2100 }
2101 
get_emin_min(void)2102 inline mp_exp_t mpreal::get_emin_min (void)
2103 {
2104 	return mpfr_get_emin_min();
2105 }
2106 
get_emin_max(void)2107 inline mp_exp_t mpreal::get_emin_max (void)
2108 {
2109 	return mpfr_get_emin_max();
2110 }
2111 
get_emax_min(void)2112 inline mp_exp_t mpreal::get_emax_min (void)
2113 {
2114 	return mpfr_get_emax_min();
2115 }
2116 
get_emax_max(void)2117 inline mp_exp_t mpreal::get_emax_max (void)
2118 {
2119 	return mpfr_get_emax_max();
2120 }
2121 
2122 //////////////////////////////////////////////////////////////////////////
2123 // Mathematical Functions
2124 //////////////////////////////////////////////////////////////////////////
sqr(const mpreal & v,mp_rnd_t rnd_mode)2125 inline const mpreal sqr(const mpreal& v, mp_rnd_t rnd_mode)
2126 {
2127 	mpreal x(v);
2128 	mpfr_sqr(x.mp,x.mp,rnd_mode);
2129 	return x;
2130 }
2131 
sqrt(const mpreal & v,mp_rnd_t rnd_mode)2132 inline const mpreal sqrt(const mpreal& v, mp_rnd_t rnd_mode)
2133 {
2134 	mpreal x(v);
2135 	mpfr_sqrt(x.mp,x.mp,rnd_mode);
2136 	return x;
2137 }
2138 
sqrt(const unsigned long int v,mp_rnd_t rnd_mode)2139 inline const mpreal sqrt(const unsigned long int v, mp_rnd_t rnd_mode)
2140 {
2141 	mpreal x;
2142 	mpfr_sqrt_ui(x.mp,v,rnd_mode);
2143 	return x;
2144 }
2145 
sqrt(const unsigned int v,mp_rnd_t rnd_mode)2146 inline const mpreal sqrt(const unsigned int v, mp_rnd_t rnd_mode)
2147 {
2148 	return sqrt(static_cast<unsigned long int>(v),rnd_mode);
2149 }
2150 
sqrt(const long int v,mp_rnd_t rnd_mode)2151 inline const mpreal sqrt(const long int v, mp_rnd_t rnd_mode)
2152 {
2153 	if (v>=0)	return sqrt(static_cast<unsigned long int>(v),rnd_mode);
2154 	else		return mpreal(); // NaN
2155 }
2156 
sqrt(const int v,mp_rnd_t rnd_mode)2157 inline const mpreal sqrt(const int v, mp_rnd_t rnd_mode)
2158 {
2159 	if (v>=0)	return sqrt(static_cast<unsigned long int>(v),rnd_mode);
2160 	else		return mpreal(); // NaN
2161 }
2162 
sqrt(const long double v,mp_rnd_t rnd_mode)2163 inline const mpreal sqrt(const long double v, mp_rnd_t rnd_mode)
2164 {
2165 	return sqrt(mpreal(v),rnd_mode);
2166 }
2167 
sqrt(const double v,mp_rnd_t rnd_mode)2168 inline const mpreal sqrt(const double v, mp_rnd_t rnd_mode)
2169 {
2170 	return sqrt(mpreal(v),rnd_mode);
2171 }
2172 
cbrt(const mpreal & v,mp_rnd_t rnd_mode)2173 inline const mpreal cbrt(const mpreal& v, mp_rnd_t rnd_mode)
2174 {
2175 	mpreal x(v);
2176 	mpfr_cbrt(x.mp,x.mp,rnd_mode);
2177 	return x;
2178 }
2179 
root(const mpreal & v,unsigned long int k,mp_rnd_t rnd_mode)2180 inline const mpreal root(const mpreal& v, unsigned long int k, mp_rnd_t rnd_mode)
2181 {
2182 	mpreal x(v);
2183 	mpfr_root(x.mp,x.mp,k,rnd_mode);
2184 	return x;
2185 }
2186 
fabs(const mpreal & v,mp_rnd_t rnd_mode)2187 inline const mpreal fabs(const mpreal& v, mp_rnd_t rnd_mode)
2188 {
2189 	mpreal x(v);
2190 	mpfr_abs(x.mp,x.mp,rnd_mode);
2191 	return x;
2192 }
2193 
abs(const mpreal & v,mp_rnd_t rnd_mode)2194 inline const mpreal abs(const mpreal& v, mp_rnd_t rnd_mode)
2195 {
2196 	mpreal x(v);
2197 	mpfr_abs(x.mp,x.mp,rnd_mode);
2198 	return x;
2199 }
2200 
dim(const mpreal & a,const mpreal & b,mp_rnd_t rnd_mode)2201 inline const mpreal dim(const mpreal& a, const mpreal& b, mp_rnd_t rnd_mode)
2202 {
2203 	mpreal x(a);
2204 	mpfr_dim(x.mp,a.mp,b.mp,rnd_mode);
2205 	return x;
2206 }
2207 
cmpabs(const mpreal & a,const mpreal & b)2208 inline int cmpabs(const mpreal& a,const mpreal& b)
2209 {
2210 	return mpfr_cmpabs(a.mp,b.mp);
2211 }
2212 
log(const mpreal & v,mp_rnd_t rnd_mode)2213 inline const mpreal log  (const mpreal& v, mp_rnd_t rnd_mode)
2214 {
2215 	mpreal x(v);
2216 	mpfr_log(x.mp,v.mp,rnd_mode);
2217 	return x;
2218 }
2219 
log2(const mpreal & v,mp_rnd_t rnd_mode)2220 inline const mpreal log2(const mpreal& v, mp_rnd_t rnd_mode)
2221 {
2222 	mpreal x(v);
2223 	mpfr_log2(x.mp,v.mp,rnd_mode);
2224 	return x;
2225 }
2226 
log10(const mpreal & v,mp_rnd_t rnd_mode)2227 inline const mpreal log10(const mpreal& v, mp_rnd_t rnd_mode)
2228 {
2229 	mpreal x(v);
2230 	mpfr_log10(x.mp,v.mp,rnd_mode);
2231 	return x;
2232 }
2233 
exp(const mpreal & v,mp_rnd_t rnd_mode)2234 inline const mpreal exp(const mpreal& v, mp_rnd_t rnd_mode)
2235 {
2236 	mpreal x(v);
2237 	mpfr_exp(x.mp,v.mp,rnd_mode);
2238 	return x;
2239 }
2240 
exp2(const mpreal & v,mp_rnd_t rnd_mode)2241 inline const mpreal exp2(const mpreal& v, mp_rnd_t rnd_mode)
2242 {
2243 	mpreal x(v);
2244 	mpfr_exp2(x.mp,v.mp,rnd_mode);
2245 	return x;
2246 }
2247 
exp10(const mpreal & v,mp_rnd_t rnd_mode)2248 inline const mpreal exp10(const mpreal& v, mp_rnd_t rnd_mode)
2249 {
2250 	mpreal x(v);
2251 	mpfr_exp10(x.mp,v.mp,rnd_mode);
2252 	return x;
2253 }
2254 
cos(const mpreal & v,mp_rnd_t rnd_mode)2255 inline const mpreal cos(const mpreal& v, mp_rnd_t rnd_mode)
2256 {
2257 	mpreal x(v);
2258 	mpfr_cos(x.mp,v.mp,rnd_mode);
2259 	return x;
2260 }
2261 
sin(const mpreal & v,mp_rnd_t rnd_mode)2262 inline const mpreal sin(const mpreal& v, mp_rnd_t rnd_mode)
2263 {
2264 	mpreal x(v);
2265 	mpfr_sin(x.mp,v.mp,rnd_mode);
2266 	return x;
2267 }
2268 
tan(const mpreal & v,mp_rnd_t rnd_mode)2269 inline const mpreal tan(const mpreal& v, mp_rnd_t rnd_mode)
2270 {
2271 	mpreal x(v);
2272 	mpfr_tan(x.mp,v.mp,rnd_mode);
2273 	return x;
2274 }
2275 
sec(const mpreal & v,mp_rnd_t rnd_mode)2276 inline const mpreal sec(const mpreal& v, mp_rnd_t rnd_mode)
2277 {
2278 	mpreal x(v);
2279 	mpfr_sec(x.mp,v.mp,rnd_mode);
2280 	return x;
2281 }
2282 
csc(const mpreal & v,mp_rnd_t rnd_mode)2283 inline const mpreal csc(const mpreal& v, mp_rnd_t rnd_mode)
2284 {
2285 	mpreal x(v);
2286 	mpfr_csc(x.mp,v.mp,rnd_mode);
2287 	return x;
2288 }
2289 
cot(const mpreal & v,mp_rnd_t rnd_mode)2290 inline const mpreal cot(const mpreal& v, mp_rnd_t rnd_mode)
2291 {
2292 	mpreal x(v);
2293 	mpfr_cot(x.mp,v.mp,rnd_mode);
2294 	return x;
2295 }
2296 
sin_cos(mpreal & s,mpreal & c,const mpreal & v,mp_rnd_t rnd_mode)2297 inline int sin_cos(mpreal& s, mpreal& c, const mpreal& v, mp_rnd_t rnd_mode)
2298 {
2299 	return mpfr_sin_cos(s.mp,c.mp,v.mp,rnd_mode);
2300 }
2301 
acos(const mpreal & v,mp_rnd_t rnd_mode)2302 inline const mpreal acos (const mpreal& v, mp_rnd_t rnd_mode)
2303 {
2304 	mpreal x(v);
2305 	mpfr_acos(x.mp,v.mp,rnd_mode);
2306 	return x;
2307 }
2308 
asin(const mpreal & v,mp_rnd_t rnd_mode)2309 inline const mpreal asin (const mpreal& v, mp_rnd_t rnd_mode)
2310 {
2311 	mpreal x(v);
2312 	mpfr_asin(x.mp,v.mp,rnd_mode);
2313 	return x;
2314 }
2315 
atan(const mpreal & v,mp_rnd_t rnd_mode)2316 inline const mpreal atan (const mpreal& v, mp_rnd_t rnd_mode)
2317 {
2318 	mpreal x(v);
2319 	mpfr_atan(x.mp,v.mp,rnd_mode);
2320 	return x;
2321 }
2322 
atan2(const mpreal & y,const mpreal & x,mp_rnd_t rnd_mode)2323 inline const mpreal atan2 (const mpreal& y, const mpreal& x, mp_rnd_t rnd_mode)
2324 {
2325 	mpreal a;
2326 	mp_prec_t yp, xp;
2327 
2328 	yp = y.get_prec();
2329 	xp = x.get_prec();
2330 
2331 	a.set_prec(yp>xp?yp:xp);
2332 
2333 	mpfr_atan2(a.mp, y.mp, x.mp, rnd_mode);
2334 
2335 	return a;
2336 }
2337 
cosh(const mpreal & v,mp_rnd_t rnd_mode)2338 inline const mpreal cosh (const mpreal& v, mp_rnd_t rnd_mode)
2339 {
2340 	mpreal x(v);
2341 	mpfr_cosh(x.mp,v.mp,rnd_mode);
2342 	return x;
2343 }
2344 
sinh(const mpreal & v,mp_rnd_t rnd_mode)2345 inline const mpreal sinh (const mpreal& v, mp_rnd_t rnd_mode)
2346 {
2347 	mpreal x(v);
2348 	mpfr_sinh(x.mp,v.mp,rnd_mode);
2349 	return x;
2350 }
2351 
tanh(const mpreal & v,mp_rnd_t rnd_mode)2352 inline const mpreal tanh (const mpreal& v, mp_rnd_t rnd_mode)
2353 {
2354 	mpreal x(v);
2355 	mpfr_tanh(x.mp,v.mp,rnd_mode);
2356 	return x;
2357 }
2358 
sech(const mpreal & v,mp_rnd_t rnd_mode)2359 inline const mpreal sech (const mpreal& v, mp_rnd_t rnd_mode)
2360 {
2361 	mpreal x(v);
2362 	mpfr_sech(x.mp,v.mp,rnd_mode);
2363 	return x;
2364 }
2365 
csch(const mpreal & v,mp_rnd_t rnd_mode)2366 inline const mpreal csch (const mpreal& v, mp_rnd_t rnd_mode)
2367 {
2368 	mpreal x(v);
2369 	mpfr_csch(x.mp,v.mp,rnd_mode);
2370 	return x;
2371 }
2372 
coth(const mpreal & v,mp_rnd_t rnd_mode)2373 inline const mpreal coth (const mpreal& v, mp_rnd_t rnd_mode)
2374 {
2375 	mpreal x(v);
2376 	mpfr_coth(x.mp,v.mp,rnd_mode);
2377 	return x;
2378 }
2379 
acosh(const mpreal & v,mp_rnd_t rnd_mode)2380 inline const mpreal acosh  (const mpreal& v, mp_rnd_t rnd_mode)
2381 {
2382 	mpreal x(v);
2383 	mpfr_acosh(x.mp,v.mp,rnd_mode);
2384 	return x;
2385 }
2386 
asinh(const mpreal & v,mp_rnd_t rnd_mode)2387 inline const mpreal asinh  (const mpreal& v, mp_rnd_t rnd_mode)
2388 {
2389 	mpreal x(v);
2390 	mpfr_asinh(x.mp,v.mp,rnd_mode);
2391 	return x;
2392 }
2393 
atanh(const mpreal & v,mp_rnd_t rnd_mode)2394 inline const mpreal atanh  (const mpreal& v, mp_rnd_t rnd_mode)
2395 {
2396 	mpreal x(v);
2397 	mpfr_atanh(x.mp,v.mp,rnd_mode);
2398 	return x;
2399 }
2400 
fac_ui(unsigned long int v,mp_rnd_t rnd_mode)2401 inline const mpreal fac_ui (unsigned long int v, mp_rnd_t rnd_mode)
2402 {
2403 	mpreal x;
2404 	mpfr_fac_ui(x.mp,v,rnd_mode);
2405 	return x;
2406 }
2407 
log1p(const mpreal & v,mp_rnd_t rnd_mode)2408 inline const mpreal log1p  (const mpreal& v, mp_rnd_t rnd_mode)
2409 {
2410 	mpreal x(v);
2411 	mpfr_log1p(x.mp,v.mp,rnd_mode);
2412 	return x;
2413 }
2414 
expm1(const mpreal & v,mp_rnd_t rnd_mode)2415 inline const mpreal expm1  (const mpreal& v, mp_rnd_t rnd_mode)
2416 {
2417 	mpreal x(v);
2418 	mpfr_expm1(x.mp,v.mp,rnd_mode);
2419 	return x;
2420 }
2421 
eint(const mpreal & v,mp_rnd_t rnd_mode)2422 inline const mpreal eint   (const mpreal& v, mp_rnd_t rnd_mode)
2423 {
2424 	mpreal x(v);
2425 	mpfr_eint(x.mp,v.mp,rnd_mode);
2426 	return x;
2427 }
2428 
gamma(const mpreal & v,mp_rnd_t rnd_mode)2429 inline const mpreal gamma (const mpreal& v, mp_rnd_t rnd_mode)
2430 {
2431 	mpreal x(v);
2432 	mpfr_gamma(x.mp,v.mp,rnd_mode);
2433 	return x;
2434 }
2435 
lngamma(const mpreal & v,mp_rnd_t rnd_mode)2436 inline const mpreal lngamma (const mpreal& v, mp_rnd_t rnd_mode)
2437 {
2438 	mpreal x(v);
2439 	mpfr_lngamma(x.mp,v.mp,rnd_mode);
2440 	return x;
2441 }
2442 
lgamma(const mpreal & v,int * signp,mp_rnd_t rnd_mode)2443 inline const mpreal lgamma (const mpreal& v, int *signp, mp_rnd_t rnd_mode)
2444 {
2445 	mpreal x(v);
2446 	mpfr_lgamma(x.mp,signp,v.mp,rnd_mode);
2447 	return x;
2448 }
2449 
zeta(const mpreal & v,mp_rnd_t rnd_mode)2450 inline const mpreal zeta (const mpreal& v, mp_rnd_t rnd_mode)
2451 {
2452 	mpreal x(v);
2453 	mpfr_zeta(x.mp,v.mp,rnd_mode);
2454 	return x;
2455 }
2456 
erf(const mpreal & v,mp_rnd_t rnd_mode)2457 inline const mpreal erf (const mpreal& v, mp_rnd_t rnd_mode)
2458 {
2459 	mpreal x(v);
2460 	mpfr_erf(x.mp,v.mp,rnd_mode);
2461 	return x;
2462 }
2463 
erfc(const mpreal & v,mp_rnd_t rnd_mode)2464 inline const mpreal erfc (const mpreal& v, mp_rnd_t rnd_mode)
2465 {
2466 	mpreal x(v);
2467 	mpfr_erfc(x.mp,v.mp,rnd_mode);
2468 	return x;
2469 }
2470 
_j0(const mpreal & v,mp_rnd_t rnd_mode)2471 inline const mpreal _j0 (const mpreal& v, mp_rnd_t rnd_mode)
2472 {
2473 	mpreal x(v);
2474 	mpfr_j0(x.mp,v.mp,rnd_mode);
2475 	return x;
2476 }
2477 
_j1(const mpreal & v,mp_rnd_t rnd_mode)2478 inline const mpreal _j1 (const mpreal& v, mp_rnd_t rnd_mode)
2479 {
2480 	mpreal x(v);
2481 	mpfr_j1(x.mp,v.mp,rnd_mode);
2482 	return x;
2483 }
2484 
_jn(long n,const mpreal & v,mp_rnd_t rnd_mode)2485 inline const mpreal _jn (long n, const mpreal& v, mp_rnd_t rnd_mode)
2486 {
2487 	mpreal x(v);
2488 	mpfr_jn(x.mp,n,v.mp,rnd_mode);
2489 	return x;
2490 }
2491 
_y0(const mpreal & v,mp_rnd_t rnd_mode)2492 inline const mpreal _y0 (const mpreal& v, mp_rnd_t rnd_mode)
2493 {
2494 	mpreal x(v);
2495 	mpfr_y0(x.mp,v.mp,rnd_mode);
2496 	return x;
2497 }
2498 
_y1(const mpreal & v,mp_rnd_t rnd_mode)2499 inline const mpreal _y1 (const mpreal& v, mp_rnd_t rnd_mode)
2500 {
2501 	mpreal x(v);
2502 	mpfr_y1(x.mp,v.mp,rnd_mode);
2503 	return x;
2504 }
2505 
_yn(long n,const mpreal & v,mp_rnd_t rnd_mode)2506 inline const mpreal _yn (long n, const mpreal& v, mp_rnd_t rnd_mode)
2507 {
2508 	mpreal x(v);
2509 	mpfr_yn(x.mp,n,v.mp,rnd_mode);
2510 	return x;
2511 }
2512 
2513 
2514 // MPFR 2.4.0 Specifics
2515 #if (MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0))
2516 
sinh_cosh(mpreal & s,mpreal & c,const mpreal & v,mp_rnd_t rnd_mode)2517 inline int sinh_cosh(mpreal& s, mpreal& c, const mpreal& v, mp_rnd_t rnd_mode)
2518 {
2519 	return mpfr_sinh_cosh(s.mp,c.mp,v.mp,rnd_mode);
2520 }
2521 
li2(const mpreal & v,mp_rnd_t rnd_mode)2522 inline const mpreal li2(const mpreal& v, mp_rnd_t rnd_mode)
2523 {
2524 	mpreal x(v);
2525 	mpfr_li2(x.mp,v.mp,rnd_mode);
2526 	return x;
2527 }
2528 
fmod(const mpreal & x,const mpreal & y,mp_rnd_t rnd_mode)2529 inline const mpreal fmod (const mpreal& x, const mpreal& y, mp_rnd_t rnd_mode)
2530 {
2531 	mpreal a;
2532 	mp_prec_t yp, xp;
2533 
2534 	yp = y.get_prec();
2535 	xp = x.get_prec();
2536 
2537 	a.set_prec(yp>xp?yp:xp);
2538 
2539 	mpfr_fmod(a.mp, x.mp, y.mp, rnd_mode);
2540 
2541 	return a;
2542 }
2543 
rec_sqrt(const mpreal & v,mp_rnd_t rnd_mode)2544 inline const mpreal rec_sqrt(const mpreal& v, mp_rnd_t rnd_mode)
2545 {
2546 	mpreal x(v);
2547 	mpfr_rec_sqrt(x.mp,v.mp,rnd_mode);
2548 	return x;
2549 }
2550 
2551 #endif
2552 
2553 //////////////////////////////////////////////////////////////////////////
2554 // Constants
const_log2(mp_prec_t prec,mp_rnd_t rnd_mode)2555 inline const mpreal const_log2 (mp_prec_t prec, mp_rnd_t rnd_mode)
2556 {
2557 	mpreal x;
2558 	x.set_prec(prec);
2559 	mpfr_const_log2(x.mp,rnd_mode);
2560 	return x;
2561 }
2562 
const_pi(mp_prec_t prec,mp_rnd_t rnd_mode)2563 inline const mpreal const_pi (mp_prec_t prec, mp_rnd_t rnd_mode)
2564 {
2565 	mpreal x;
2566 	x.set_prec(prec);
2567 	mpfr_const_pi(x.mp,rnd_mode);
2568 	return x;
2569 }
2570 
const_euler(mp_prec_t prec,mp_rnd_t rnd_mode)2571 inline const mpreal const_euler (mp_prec_t prec, mp_rnd_t rnd_mode)
2572 {
2573 	mpreal x;
2574 	x.set_prec(prec);
2575 	mpfr_const_euler(x.mp,rnd_mode);
2576 	return x;
2577 }
2578 
const_catalan(mp_prec_t prec,mp_rnd_t rnd_mode)2579 inline const mpreal const_catalan (mp_prec_t prec, mp_rnd_t rnd_mode)
2580 {
2581 	mpreal x;
2582 	x.set_prec(prec);
2583 	mpfr_const_catalan(x.mp,rnd_mode);
2584 	return x;
2585 }
2586 
2587 //////////////////////////////////////////////////////////////////////////
2588 // Integer Related Functions
rint(const mpreal & v,mp_rnd_t rnd_mode)2589 inline const mpreal rint(const mpreal& v, mp_rnd_t rnd_mode)
2590 {
2591 	mpreal x(v);
2592 	mpfr_rint(x.mp,v.mp,rnd_mode);
2593 	return x;
2594 }
2595 
ceil(const mpreal & v)2596 inline const mpreal ceil(const mpreal& v)
2597 {
2598 	mpreal x(v);
2599 	mpfr_ceil(x.mp,v.mp);
2600 	return x;
2601 
2602 }
2603 
floor(const mpreal & v)2604 inline const mpreal floor(const mpreal& v)
2605 {
2606 	mpreal x(v);
2607 	mpfr_floor(x.mp,v.mp);
2608 	return x;
2609 }
2610 
round(const mpreal & v)2611 inline const mpreal round(const mpreal& v)
2612 {
2613 	mpreal x(v);
2614 	mpfr_round(x.mp,v.mp);
2615 	return x;
2616 }
2617 
trunc(const mpreal & v)2618 inline const mpreal trunc(const mpreal& v)
2619 {
2620 	mpreal x(v);
2621 	mpfr_trunc(x.mp,v.mp);
2622 	return x;
2623 }
2624 
rint_ceil(const mpreal & v,mp_rnd_t rnd_mode)2625 inline const mpreal rint_ceil (const mpreal& v, mp_rnd_t rnd_mode)
2626 {
2627 	mpreal x(v);
2628 	mpfr_rint_ceil(x.mp,v.mp,rnd_mode);
2629 	return x;
2630 }
2631 
rint_floor(const mpreal & v,mp_rnd_t rnd_mode)2632 inline const mpreal rint_floor(const mpreal& v, mp_rnd_t rnd_mode)
2633 {
2634 	mpreal x(v);
2635 	mpfr_rint_floor(x.mp,v.mp,rnd_mode);
2636 	return x;
2637 }
2638 
rint_round(const mpreal & v,mp_rnd_t rnd_mode)2639 inline const mpreal rint_round(const mpreal& v, mp_rnd_t rnd_mode)
2640 {
2641 	mpreal x(v);
2642 	mpfr_rint_round(x.mp,v.mp,rnd_mode);
2643 	return x;
2644 }
2645 
rint_trunc(const mpreal & v,mp_rnd_t rnd_mode)2646 inline const mpreal rint_trunc(const mpreal& v, mp_rnd_t rnd_mode)
2647 {
2648 	mpreal x(v);
2649 	mpfr_rint_trunc(x.mp,v.mp,rnd_mode);
2650 	return x;
2651 }
2652 
frac(const mpreal & v,mp_rnd_t rnd_mode)2653 inline const mpreal frac (const mpreal& v, mp_rnd_t rnd_mode)
2654 {
2655 	mpreal x(v);
2656 	mpfr_frac(x.mp,v.mp,rnd_mode);
2657 	return x;
2658 }
2659 
2660 //////////////////////////////////////////////////////////////////////////
2661 // Miscellaneous Functions
swap(mpreal & a,mpreal & b)2662 inline void swap(mpreal& a, mpreal& b)
2663 {
2664 	mpfr_swap(a.mp,b.mp);
2665 }
2666 
2667 #ifndef max
max(const mpreal & x,const mpreal & y)2668 inline const mpreal max(const mpreal& x, const mpreal& y)
2669 {
2670 	return (x>y?x:y);
2671 }
2672 #endif
2673 
2674 #ifndef min
min(const mpreal & x,const mpreal & y)2675 inline const mpreal min(const mpreal& x, const mpreal& y)
2676 {
2677 	return (x<y?x:y);
2678 }
2679 #endif
2680 
nexttoward(const mpreal & x,const mpreal & y)2681 inline const mpreal nexttoward (const mpreal& x, const mpreal& y)
2682 {
2683 	mpreal a(x);
2684 	mpfr_nexttoward(a.mp,y.mp);
2685 	return a;
2686 }
2687 
nextabove(const mpreal & x)2688 inline const mpreal nextabove  (const mpreal& x)
2689 {
2690 	mpreal a(x);
2691 	mpfr_nextabove(a.mp);
2692 	return a;
2693 }
2694 
nextbelow(const mpreal & x)2695 inline const mpreal nextbelow  (const mpreal& x)
2696 {
2697 	mpreal a(x);
2698 	mpfr_nextbelow(a.mp);
2699 	return a;
2700 }
2701 
urandomb(gmp_randstate_t & state)2702 inline const mpreal urandomb (gmp_randstate_t& state)
2703 {
2704 	mpreal x;
2705 	mpfr_urandomb(x.mp,state);
2706 	return x;
2707 }
2708 
2709 // inline const mpreal random2 (mp_size_t size, mp_exp_t exp)
2710 // {
2711 // 	mpreal x;
2712 // 	mpfr_random2(x.mp,size,exp);
2713 // 	return x;
2714 // }
2715 
2716 //////////////////////////////////////////////////////////////////////////
2717 // Set/Get global properties
set_default_prec(mp_prec_t prec)2718 inline void mpreal::set_default_prec(mp_prec_t prec)
2719 {
2720 	default_prec = prec;
2721 	mpfr_set_default_prec(prec);
2722 }
2723 
get_default_prec()2724 inline mp_prec_t mpreal::get_default_prec()
2725 {
2726 	return mpfr_get_default_prec();
2727 }
2728 
set_default_base(int base)2729 inline void mpreal::set_default_base(int base)
2730 {
2731 	default_base = base;
2732 }
2733 
get_default_base()2734 inline int mpreal::get_default_base()
2735 {
2736 	return default_base;
2737 }
2738 
set_default_rnd(mp_rnd_t rnd_mode)2739 inline void mpreal::set_default_rnd(mp_rnd_t rnd_mode)
2740 {
2741 	default_rnd =  rnd_mode;
2742 	mpfr_set_default_rounding_mode(rnd_mode);
2743 }
2744 
get_default_rnd()2745 inline mp_rnd_t mpreal::get_default_rnd()
2746 {
2747 	return mpfr_get_default_rounding_mode();
2748 }
2749 
set_double_bits(int dbits)2750 inline void mpreal::set_double_bits(int dbits)
2751 {
2752 	double_bits = dbits;
2753 }
2754 
get_double_bits()2755 inline int mpreal::get_double_bits()
2756 {
2757 	return double_bits;
2758 }
2759 
fits_in_bits(double x,int n)2760 inline bool mpreal::fits_in_bits(double x, int n)
2761 {
2762 	int i;
2763 	double t;
2764 	return IsInf(x) || (std::modf ( std::ldexp ( std::frexp ( x, &i ), n ), &t ) == 0.0);
2765 }
2766 
pow(const mpreal & a,const mpreal & b,mp_rnd_t rnd_mode)2767 inline const mpreal pow(const mpreal& a, const mpreal& b, mp_rnd_t rnd_mode)
2768 {
2769 	mpreal x(a);
2770 	mpfr_pow(x.mp,x.mp,b.mp,rnd_mode);
2771 	return x;
2772 }
2773 
pow(const mpreal & a,const mpz_t b,mp_rnd_t rnd_mode)2774 inline const mpreal pow(const mpreal& a, const mpz_t b, mp_rnd_t rnd_mode)
2775 {
2776 	mpreal x(a);
2777 	mpfr_pow_z(x.mp,x.mp,b,rnd_mode);
2778 	return x;
2779 }
2780 
pow(const mpreal & a,const unsigned long int b,mp_rnd_t rnd_mode)2781 inline const mpreal pow(const mpreal& a, const unsigned long int b, mp_rnd_t rnd_mode)
2782 {
2783 	mpreal x(a);
2784 	mpfr_pow_ui(x.mp,x.mp,b,rnd_mode);
2785 	return x;
2786 }
2787 
pow(const mpreal & a,const unsigned int b,mp_rnd_t rnd_mode)2788 inline const mpreal pow(const mpreal& a, const unsigned int b, mp_rnd_t rnd_mode)
2789 {
2790 	return pow(a,static_cast<unsigned long int>(b),rnd_mode);
2791 }
2792 
pow(const mpreal & a,const long int b,mp_rnd_t rnd_mode)2793 inline const mpreal pow(const mpreal& a, const long int b, mp_rnd_t rnd_mode)
2794 {
2795 	mpreal x(a);
2796 	mpfr_pow_si(x.mp,x.mp,b,rnd_mode);
2797 	return x;
2798 }
2799 
pow(const mpreal & a,const int b,mp_rnd_t rnd_mode)2800 inline const mpreal pow(const mpreal& a, const int b, mp_rnd_t rnd_mode)
2801 {
2802 	return pow(a,static_cast<long int>(b),rnd_mode);
2803 }
2804 
pow(const mpreal & a,const long double b,mp_rnd_t rnd_mode)2805 inline const mpreal pow(const mpreal& a, const long double b, mp_rnd_t rnd_mode)
2806 {
2807 	return pow(a,mpreal(b),rnd_mode);
2808 }
2809 
pow(const mpreal & a,const double b,mp_rnd_t rnd_mode)2810 inline const mpreal pow(const mpreal& a, const double b, mp_rnd_t rnd_mode)
2811 {
2812 	return pow(a,mpreal(b),rnd_mode);
2813 }
2814 
pow(const unsigned long int a,const mpreal & b,mp_rnd_t rnd_mode)2815 inline const mpreal pow(const unsigned long int a, const mpreal& b, mp_rnd_t rnd_mode)
2816 {
2817 	mpreal x(a);
2818 	mpfr_ui_pow(x.mp,a,b.mp,rnd_mode);
2819 	return x;
2820 }
2821 
pow(const unsigned int a,const mpreal & b,mp_rnd_t rnd_mode)2822 inline const mpreal pow(const unsigned int a, const mpreal& b, mp_rnd_t rnd_mode)
2823 {
2824 	return pow(static_cast<unsigned long int>(a),b,rnd_mode);
2825 }
2826 
pow(const long int a,const mpreal & b,mp_rnd_t rnd_mode)2827 inline const mpreal pow(const long int a, const mpreal& b, mp_rnd_t rnd_mode)
2828 {
2829 	if (a>=0) 	return pow(static_cast<unsigned long int>(a),b,rnd_mode);
2830 	else		return pow(mpreal(a),b,rnd_mode);
2831 }
2832 
pow(const int a,const mpreal & b,mp_rnd_t rnd_mode)2833 inline const mpreal pow(const int a, const mpreal& b, mp_rnd_t rnd_mode)
2834 {
2835 	if (a>=0) 	return pow(static_cast<unsigned long int>(a),b,rnd_mode);
2836 	else		return pow(mpreal(a),b,rnd_mode);
2837 }
2838 
pow(const long double a,const mpreal & b,mp_rnd_t rnd_mode)2839 inline const mpreal pow(const long double a, const mpreal& b, mp_rnd_t rnd_mode)
2840 {
2841 	return pow(mpreal(a),b,rnd_mode);
2842 }
2843 
pow(const double a,const mpreal & b,mp_rnd_t rnd_mode)2844 inline const mpreal pow(const double a, const mpreal& b, mp_rnd_t rnd_mode)
2845 {
2846 	return pow(mpreal(a),b,rnd_mode);
2847 }
2848 
2849 // pow unsigned long int
pow(const unsigned long int a,const unsigned long int b,mp_rnd_t rnd_mode)2850 inline const mpreal pow(const unsigned long int a, const unsigned long int b, mp_rnd_t rnd_mode)
2851 {
2852 	mpreal x(a);
2853 	mpfr_ui_pow_ui(x.mp,a,b,rnd_mode);
2854 	return x;
2855 }
2856 
pow(const unsigned long int a,const unsigned int b,mp_rnd_t rnd_mode)2857 inline const mpreal pow(const unsigned long int a, const unsigned int b, mp_rnd_t rnd_mode)
2858 {
2859 	return pow(a,static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2860 }
2861 
pow(const unsigned long int a,const long int b,mp_rnd_t rnd_mode)2862 inline const mpreal pow(const unsigned long int a, const long int b, mp_rnd_t rnd_mode)
2863 {
2864 	if(b>0)	return pow(a,static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2865 	else	return pow(a,mpreal(b),rnd_mode); //mpfr_ui_pow
2866 }
2867 
pow(const unsigned long int a,const int b,mp_rnd_t rnd_mode)2868 inline const mpreal pow(const unsigned long int a, const int b, mp_rnd_t rnd_mode)
2869 {
2870 	if(b>0)	return pow(a,static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2871 	else	return pow(a,mpreal(b),rnd_mode); //mpfr_ui_pow
2872 }
2873 
pow(const unsigned long int a,const long double b,mp_rnd_t rnd_mode)2874 inline const mpreal pow(const unsigned long int a, const long double b, mp_rnd_t rnd_mode)
2875 {
2876 	return pow(a,mpreal(b),rnd_mode); //mpfr_ui_pow
2877 }
2878 
pow(const unsigned long int a,const double b,mp_rnd_t rnd_mode)2879 inline const mpreal pow(const unsigned long int a, const double b, mp_rnd_t rnd_mode)
2880 {
2881 	return pow(a,mpreal(b),rnd_mode); //mpfr_ui_pow
2882 }
2883 
2884 // pow unsigned int
pow(const unsigned int a,const unsigned long int b,mp_rnd_t rnd_mode)2885 inline const mpreal pow(const unsigned int a, const unsigned long int b, mp_rnd_t rnd_mode)
2886 {
2887 	return pow(static_cast<unsigned long int>(a),b,rnd_mode); //mpfr_ui_pow_ui
2888 }
2889 
pow(const unsigned int a,const unsigned int b,mp_rnd_t rnd_mode)2890 inline const mpreal pow(const unsigned int a, const unsigned int b, mp_rnd_t rnd_mode)
2891 {
2892 	return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2893 }
2894 
pow(const unsigned int a,const long int b,mp_rnd_t rnd_mode)2895 inline const mpreal pow(const unsigned int a, const long int b, mp_rnd_t rnd_mode)
2896 {
2897 	if(b>0)	return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2898 	else	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2899 }
2900 
pow(const unsigned int a,const int b,mp_rnd_t rnd_mode)2901 inline const mpreal pow(const unsigned int a, const int b, mp_rnd_t rnd_mode)
2902 {
2903 	if(b>0)	return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2904 	else	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2905 }
2906 
pow(const unsigned int a,const long double b,mp_rnd_t rnd_mode)2907 inline const mpreal pow(const unsigned int a, const long double b, mp_rnd_t rnd_mode)
2908 {
2909 	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2910 }
2911 
pow(const unsigned int a,const double b,mp_rnd_t rnd_mode)2912 inline const mpreal pow(const unsigned int a, const double b, mp_rnd_t rnd_mode)
2913 {
2914 	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2915 }
2916 
2917 // pow long int
pow(const long int a,const unsigned long int b,mp_rnd_t rnd_mode)2918 inline const mpreal pow(const long int a, const unsigned long int b, mp_rnd_t rnd_mode)
2919 {
2920 	if (a>0) return pow(static_cast<unsigned long int>(a),b,rnd_mode); //mpfr_ui_pow_ui
2921 	else	 return pow(mpreal(a),b,rnd_mode); //mpfr_pow_ui
2922 }
2923 
pow(const long int a,const unsigned int b,mp_rnd_t rnd_mode)2924 inline const mpreal pow(const long int a, const unsigned int b, mp_rnd_t rnd_mode)
2925 {
2926 	if (a>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode);  //mpfr_ui_pow_ui
2927 	else	 return pow(mpreal(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_pow_ui
2928 }
2929 
pow(const long int a,const long int b,mp_rnd_t rnd_mode)2930 inline const mpreal pow(const long int a, const long int b, mp_rnd_t rnd_mode)
2931 {
2932 	if (a>0)
2933 	{
2934 		if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2935 		else	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2936 	}else{
2937 		return pow(mpreal(a),b,rnd_mode); // mpfr_pow_si
2938 	}
2939 }
2940 
pow(const long int a,const int b,mp_rnd_t rnd_mode)2941 inline const mpreal pow(const long int a, const int b, mp_rnd_t rnd_mode)
2942 {
2943 	if (a>0)
2944 	{
2945 		if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2946 		else	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2947 	}else{
2948 		return pow(mpreal(a),static_cast<long int>(b),rnd_mode); // mpfr_pow_si
2949 	}
2950 }
2951 
pow(const long int a,const long double b,mp_rnd_t rnd_mode)2952 inline const mpreal pow(const long int a, const long double b, mp_rnd_t rnd_mode)
2953 {
2954 	if (a>=0) 	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2955 	else		return pow(mpreal(a),mpreal(b),rnd_mode); //mpfr_pow
2956 }
2957 
pow(const long int a,const double b,mp_rnd_t rnd_mode)2958 inline const mpreal pow(const long int a, const double b, mp_rnd_t rnd_mode)
2959 {
2960 	if (a>=0) 	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2961 	else		return pow(mpreal(a),mpreal(b),rnd_mode); //mpfr_pow
2962 }
2963 
2964 // pow int
pow(const int a,const unsigned long int b,mp_rnd_t rnd_mode)2965 inline const mpreal pow(const int a, const unsigned long int b, mp_rnd_t rnd_mode)
2966 {
2967 	if (a>0) return pow(static_cast<unsigned long int>(a),b,rnd_mode); //mpfr_ui_pow_ui
2968 	else	 return pow(mpreal(a),b,rnd_mode); //mpfr_pow_ui
2969 }
2970 
pow(const int a,const unsigned int b,mp_rnd_t rnd_mode)2971 inline const mpreal pow(const int a, const unsigned int b, mp_rnd_t rnd_mode)
2972 {
2973 	if (a>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode);  //mpfr_ui_pow_ui
2974 	else	 return pow(mpreal(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_pow_ui
2975 }
2976 
pow(const int a,const long int b,mp_rnd_t rnd_mode)2977 inline const mpreal pow(const int a, const long int b, mp_rnd_t rnd_mode)
2978 {
2979 	if (a>0)
2980 	{
2981 		if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2982 		else	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2983 	}else{
2984 		return pow(mpreal(a),b,rnd_mode); // mpfr_pow_si
2985 	}
2986 }
2987 
pow(const int a,const int b,mp_rnd_t rnd_mode)2988 inline const mpreal pow(const int a, const int b, mp_rnd_t rnd_mode)
2989 {
2990 	if (a>0)
2991 	{
2992 		if(b>0) return pow(static_cast<unsigned long int>(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_ui_pow_ui
2993 		else	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
2994 	}else{
2995 		return pow(mpreal(a),static_cast<long int>(b),rnd_mode); // mpfr_pow_si
2996 	}
2997 }
2998 
pow(const int a,const long double b,mp_rnd_t rnd_mode)2999 inline const mpreal pow(const int a, const long double b, mp_rnd_t rnd_mode)
3000 {
3001 	if (a>=0) 	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3002 	else		return pow(mpreal(a),mpreal(b),rnd_mode); //mpfr_pow
3003 }
3004 
pow(const int a,const double b,mp_rnd_t rnd_mode)3005 inline const mpreal pow(const int a, const double b, mp_rnd_t rnd_mode)
3006 {
3007 	if (a>=0) 	return pow(static_cast<unsigned long int>(a),mpreal(b),rnd_mode); //mpfr_ui_pow
3008 	else		return pow(mpreal(a),mpreal(b),rnd_mode); //mpfr_pow
3009 }
3010 
3011 // pow long double
pow(const long double a,const long double b,mp_rnd_t rnd_mode)3012 inline const mpreal pow(const long double a, const long double b, mp_rnd_t rnd_mode)
3013 {
3014 	return pow(mpreal(a),mpreal(b),rnd_mode);
3015 }
3016 
pow(const long double a,const unsigned long int b,mp_rnd_t rnd_mode)3017 inline const mpreal pow(const long double a, const unsigned long int b, mp_rnd_t rnd_mode)
3018 {
3019 	return pow(mpreal(a),b,rnd_mode); //mpfr_pow_ui
3020 }
3021 
pow(const long double a,const unsigned int b,mp_rnd_t rnd_mode)3022 inline const mpreal pow(const long double a, const unsigned int b, mp_rnd_t rnd_mode)
3023 {
3024 	return pow(mpreal(a),static_cast<unsigned long int>(b),rnd_mode); //mpfr_pow_ui
3025 }
3026 
pow(const long double a,const long int b,mp_rnd_t rnd_mode)3027 inline const mpreal pow(const long double a, const long int b, mp_rnd_t rnd_mode)
3028 {
3029 	return pow(mpreal(a),b,rnd_mode); // mpfr_pow_si
3030 }
3031 
pow(const long double a,const int b,mp_rnd_t rnd_mode)3032 inline const mpreal pow(const long double a, const int b, mp_rnd_t rnd_mode)
3033 {
3034 	return pow(mpreal(a),static_cast<long int>(b),rnd_mode); // mpfr_pow_si
3035 }
3036 
pow(const double a,const double b,mp_rnd_t rnd_mode)3037 inline const mpreal pow(const double a, const double b, mp_rnd_t rnd_mode)
3038 {
3039 	return pow(mpreal(a),mpreal(b),rnd_mode);
3040 }
3041 
pow(const double a,const unsigned long int b,mp_rnd_t rnd_mode)3042 inline const mpreal pow(const double a, const unsigned long int b, mp_rnd_t rnd_mode)
3043 {
3044 	return pow(mpreal(a),b,rnd_mode); // mpfr_pow_ui
3045 }
3046 
pow(const double a,const unsigned int b,mp_rnd_t rnd_mode)3047 inline const mpreal pow(const double a, const unsigned int b, mp_rnd_t rnd_mode)
3048 {
3049 	return pow(mpreal(a),static_cast<unsigned long int>(b),rnd_mode); // mpfr_pow_ui
3050 }
3051 
pow(const double a,const long int b,mp_rnd_t rnd_mode)3052 inline const mpreal pow(const double a, const long int b, mp_rnd_t rnd_mode)
3053 {
3054 	return pow(mpreal(a),b,rnd_mode); // mpfr_pow_si
3055 }
3056 
pow(const double a,const int b,mp_rnd_t rnd_mode)3057 inline const mpreal pow(const double a, const int b, mp_rnd_t rnd_mode)
3058 {
3059 	return pow(mpreal(a),static_cast<long int>(b),rnd_mode); // mpfr_pow_si
3060 }
3061 }
3062 
3063 // Explicit specialization of std::swap for mpreal numbers
3064 // Thus standard algorithms will use efficient version of swap (due to Koenig lookup)
3065 // Non-throwing swap C++ idiom: http://en.wikibooks.org/wiki/More_C%2B%2B_Idioms/Non-throwing_swap
3066 namespace std
3067 {
3068 	template <>
swap(mpfr::mpreal & x,mpfr::mpreal & y)3069 	inline void swap(mpfr::mpreal& x, mpfr::mpreal& y)
3070 	{
3071 		return mpfr::swap(x, y);
3072 	}
3073 }
3074 
3075 #endif /* __MP_REAL_H__ */
3076