July 15, 2010: Added support for Eigen C++ template library for linear algebra.
June 15, 2010: Added support for the MPFR 3.0.0.
May 21, 2010: MPFR C++ has been included in MBLAS and MLAPACK.
full history…
Introduction
This is the home page of C++ wrapper for MPFR library. MPFR library allows user to conduct floating-point calculations with virtually any (restricted by available memory only) precision with correct rounding. Besides simple arithmetic operations like “+” and “/” the whole set of mathematical functions is supported: sin, sqrt, pow, log, etc.
MPFR defines custom C-language type to represent floating-point number – mpfr_t. Mathematical manipulations with mpfr_t – variables are done through assembler-look-like functions. For instance, to add two numbers x and y with result in z special function mpfr_add(z,x,y,GMP_RNDN) should be called.
To illustrate the situation, let’s consider two versions of the Schwefel function:
// double - version
double schwefel(double x)
{
return 418.9829-x*sin(sqrt(abs(x)));
}
//MPFR C - version
void mpfr_schwefel(mpfr_t y, mpfr_t x)
{
mpfr_t t;
mpfr_init(t);
mpfr_abs(t,x,GMP_RNDN);
mpfr_sqrt(t,t,GMP_RNDN);
mpfr_sin(t,t,GMP_RNDN);
mpfr_mul(t,t,x,GMP_RNDN);
mpfr_set_str(y,“418.9829“,10,GMP_RNDN);
mpfr_sub(y,y,t,GMP_RNDN);
mpfr_clear(t);
}
Obviously to write mathematical expressions in such cryptic form requires some time, practice and can be error prone for complex formulas. Moreover in order to use MPFR library existing software should be completely rewritten!
C++ wrapper for MPFR aims to alleviate these issues. It introduces new C++ type for floating point numbers – mpreal, which encapsulates native mpfr_t. All arithmetic and Boolean operators (+, -, *, /, >, !=, etc.) are implemented for mpreal numbers through operator overloading technique. Mathematical functions (sqrt, pow, sin, cos, etc.) are supported too. This makes possible to use MPFR calculations in the same simple way as calculations with numbers of built-in types double or float.
MPFR C++ version of the Schwefel function is:
// MPFR C++ - version
mpreal mpfr_schwefel(mpreal& x)
{
return "418.9829"-x*sin(sqrt(abs(x)));
}
Features
Main MPFR C++ features are:
- In contrast to other C++ interfaces MPFR C++ is up to date with the latest version of MPFR library – 3.0.0.
- MPFR C++ is the only C++ wrapper which natively supports both commonly used development environments: GNU G++ and Microsoft Visual C++.
- MPFR C++ allows usage of human-friendly notation for mathematical expressions. For example, to write
z = x + yinstead ofmpfr_add(z,x,y,…). All arithmetic and Boolean operators along with standard mathematical functions are supported. Precision and rounding mode can be easily controlled too. - MPFR C++ allows painless porting from built-in types to MPFR with minimal altering of already typed and tested mathematical expressions. In most cases only renaming of types is needed (e.g.
double->mpreal). - MPFR C++ is simple to use. It consists of two source files only: mpreal.h and mpreal.cpp. Just include header file mpreal.h and use
mprealnumbers the same way you usually usedoubleorfloatones. Example is included in distribution. - MPFR C++ provides high performance interface for plain MPFR. Most operations are designed to be inline which means native MPFR function is called directly without C++ overhead. In contrast to other wrappers MPFR C++ doesn’t use advanced C++ features like function objects, virtual functions, etc. because of significant decrease of speed they suffer.
- MPFR C++ keeps correct accuracy of intermediate calculations during complex expression evaluation in order to obtain precise final result. Other C++ wrappers use different strategies on handling intermediate calculations which could lead to significant accuracy decreasing of final result (see details in Internals section).
MPFR C++ is under LGPL v2.1 license, which permits usage in free and commercial projects.
To download MPFR C++ please visit download section.
Internals
One of the most important features that make MPFR C++ distinct from the others is treatment of intermediate results during calculation of whole mathematical expression.
C++ decomposes expressions into atomic operations like “+”, “*“ and stores its results in temporary variables for further usage. For example, expression x = a*b+c*d is calculated by C++ as following:
t1 = a*b
t2 = c*d
t3 = t1+t2
x = t3
where t1, t2, t3 are temporary variables, actually hidden from user. In order to get correct final value such intermediate results should be treated carefully with enough precision.
Main MPFR C++ rule for intermediate results is prec(rop) = max(prec(op1),prec(op2)): result of operation is stored with the max of the precisions of operands.
In other words intermediate operations are conducted with the reasonably maximum precision defined by precisions of arguments and doesn’t depend on precision of target variable (in our example – x). Only final result are rounded to the precision of target variable (x = t3).
This is the most natural rule dictated by practical applications, such as numerical differentiation, digital filter convolution, etc. where in order to obtain correct and accurate final result intermediate calculations have to be done with high accuracy.
Other libraries calculate intermediate results with the precision of final variable or round them to some not obvious to user precision which is independent from precision of arguments as well as from precision of final variable. This could lead to significant accuracy decreasing of final result.
Download
The MPFR C++ is available for download in the following archive formats:
- mpfrc++.tar.gz (20KB)
- mpfrc++.zip (30KB)
MPFR C++ consists of two files: mpreal.h and mpreal.cpp. In order to use MPFR C++ header file mpreal.h should be #included and mpreal.cpp should be compiled along with other source modules of the program. Besides this MPIR (or GMP) and MPFR libraries should be installed on PC and be available for linker.
MPFR C++ is under LGPL v2.1 license, which permits usage in free and commercial projects. If you use the library, I consider it appropriate to give clear reference on this project at some place. This can either be the About dialog of your application or the documentation of your software. If you like MPFR C++, I encourage you to support it by:
- posting a brief description of how you use it in your project (in comments section below).
- sending me bug reports, wishlists, suggestions.
Please remember that MPFR (>= 2.3.1) and MPIR (>=2.0.0) (or GMP >=4.2.1) are needed for MPFR C++. If target system is Windows then MPIR and MPFR can be compiled following to Brian Gladman’s guide.
Already compiled MPFR 3.0.0 and MPIR 2.0.0 for Windows can be found in this archive: win32_mpir_mpfr.zip (4.0MB)
Software using MPFR C++
Classes for Physics and Mathematics by Ulrich Mutze.
Eigen C++ template library for linear algebra: vectors, matrices, and related algorithms. Use additional header file mpreal_eigen.h to embed mpreal type into Eigen. Example of usage: mpreal_eigen_test.cpp. Contributed by Konstantin Holoborodko.-
The Template Numerical Toolkit (TNT) is a collection of interfaces and reference implementations of numerical objects useful for scientific computing in C++. - MPACK: Multiple precision arithmetic BLAS (MBLAS) and LAPACK (MLAPACK).
Let me know about your project so I can add it to the list.
History
July 15, 2010:
- Added functions: machine_epsilon(), mpreal_min(), mpreal_max().
- Added support for Eigen C++ template library by additional header file mpreal_eigen.h. Thanks to Konstantin Holoborodko!
June 15, 2010:
- Added support for the latest MPFR 3.0.0.
- Fixed several bugs. Added const_infinity() function. Thanks to Pere Constans!
May 21, 2010:
- MPFR C++ has been included in MBLAS and MLAPACK.
August 31, 2009:
- Added support for Linux Intel Compiler. Cause of minor warning removed. Thanks to Heinz van Saanen!
May 11, 2009:
- MPFR C++ is under LGPL v2.1 now.
March 9, 2009:
- Added support of the latest MPFR 2.4.1. Now MPFR C++ is the only C++ interface which is up to date with the newest MPFR version.
November 18, 2008:
- Added specialized std::swap() for mpreal numbers. Now standard generic algorithms will use efficient version of swap. Fixed some bugs. C++ standard conformance is improved. Thanks to Fokko Beekhof!
October 1, 2008:
- Added operators <<=, >>=. Improved output formatting. Fixed some bugs. Thanks to Helmut Jarausch!
September 30, 2008:
- Added functions frexp, ldexp, modf. Thanks to Helmut Jarausch!
September 29th, 2008:
- Fixed incompatibility with GNU C/C++. Updated output operator<< to take care of ostream.precision() setting. Many thanks to Helmut Jarausch!
- Minor changes to remove warnings in compilation by gcc.
- Test: gcc 4.2.3 – ok. Thanks to Dmitriy Gubanov!
September 4th, 2008:
- Add: check for inexact conversion from floating point to mpreal. Thanks to Brian Gladman!
August 15th, 2008:
- Add: operators for comparison of mpreal with built-in types.
August 14th, 2008:
- Add: division operators when first operand is double.
August 11th, 2008:
- Add: pre-/postfix operators of increment/decrement. Thanks to Brian Gladman!
- Update: processing of mpz_t and mpq_t types is improved.
August 8th, 2008:
- Update: more convenient copy – constructor for mpfr_t type.
August 5th, 2008:
- Update: majority of functions are changed to be inline.
- Test: gcc 4.2.3 – ok. Thanks to Dmitriy Gubanov!
August 4th, 2008:
- Test: Visual Studio 2005 – ok.
- Public release.
July 11th, 2008: Project start.
Acknowledgements
My sincere gratitude goes to the creators of MPFR for making their library available, Brian Gladman for making MPFR and GMP ports for Windows, authors of other C++ wrappers for MPFR for pioneering the subject.
Loading ...
31 Comments
Hi Pavel,
Thank you for acknowledging my Windows port of MPFR.
I run all the mpfr C++ wrappers and I tried yours out today.
The only problem I had is that you don’t overload the pre and post increment operations (++x, –x, x++, x–) for the mpreal type.
You might be interested in a special feature that I have added to mpfrcpp for ‘floating point to mpreal’ conversions. I have added a global value called ‘doublebits’ that sets the maximum number of bits allowed in a number when a floating point value is converted to an mpreal (doublebits has a default value and can be set like precision)
In any conversion from floating point to mpreal an exception is raised if the value being converted has more than ‘doublebits’ bits. An exception is hence generated when an inexact ‘double to mpfr’ conversion is attempted.
For example with doublebits = 20, the value 8.0 will convert but the value 3.1415926 will raise an exception since it has too many bits set and may hence create a low precision value . This makes it much easier to trap values that are not precise enough to be translated from floating point to mpreals.
I haqve not done a lot with your wrapper yet but it looks easy to use – thank you for your efforts on this.
best regards,
Brian Gladman
Hello Brian,
Thank you for testing mpreal and for valuable suggestions!
Post-/prefix increment/decrement operations are already added.
I’ll think about how to handle inexact “double to mpreal” attempts of conversion. Now user can use global MPFR’s functions (mpfr_inexflag_p, mpfr_overflow_p, etc.) to check for such abnormal situations. However, mpreal is focusing on processing expressions, not atomic operations (as MPFR doing). Complex expression could include many intermediate conversions from double to mpreal and how to inform user in which one we have a problem? Maybe using global functions (or exceptions) for whole expression is still the best and only one way.
Thank you!
Pavel.
Hello Mr. Holoborodko,
my g++ compiler (3.4.4 with cygwin) does not like your mreal.h. Your style to use a static class member as default value for a function argument as in
mpreal(const mpq_t u, mp_prec_t prec = default_prec, mp_rnd_t mode = default_rnd);
This causes the message
error: `default_rnd’ was not declared in this scope.
It is tempting and logical to use static data members for default initialization of function arguments but I found it always not to work with some compilers. So I eliminated it from the list of my personal idioms. Would be nice to have a gnu-compatible version of your mpreal – files.
Best wishes
Ulrich
Hello Mr. Mutze,
Thank you for your feedback!
MPFR C++ supports native compilers both for Unix (GNU G++) and Windows (Visual Studio).
Usually I do tests for GNU with g++ 4.1.2 and 4.2.3 using CentOS and Ubuntu. It works fine. Moreover many people are using MPFR C++ with GNU G++ on various platforms successfully.
I think that cygwin is using too old version of g++ compiler. It should be just updated to fix this problem.
Anyway I will check this issue. Or even better, other cygwin users would appreciate if you would provide fix for this problem. Personally I am not using cygwin, so contribution of its experienced user is very welcome.
Thank you!
Dear Mr. Holoborodko,
your wrapper seems to work like a charm together with eigen.
(See: http://eigen.tuxfamily.org/dox/CustomizingEigen.html#CustomScalarType)
Thanks,
Simon
Dear Simon,
I appreciate your feedback very much!
Thank you,
Pavel.
P.S. It would be fine if I would be just ‘Pavel’ for everyone.
Simon,
What do you mean “seems to work”? Sure, Eigen compiles fine with MPFR C++, but at least for me it crashes when the ~Matrix() destructor is called.
Pavel, unless you can point me to some working example of Eigen with MPFR C++, then I think it is false advertising to put Eigen in your “Software using MPFR C++”. If it really does work, I would appreciate a pointer to some more information.
Thanks.
As for July 15, 2010, Eigen developers have fixed these bug.
Hi Pavel.
I must say your wrapper is fantastic. I really don’t imagine myself changing all my code to those criptic MPFR functions.
However, I must ask you something too:
When I use some mpfr function, using your wrapper, the result always has the same precision as the argument supplied to the function, if the argument is a mpreal. I know this is a MPFR behavior, and in MPFR you can change this. In your wrapper, there is a way to select the precision of the result, without regard the precision of the argument? I haven’t find one yet. By the way, I tried to use mpreal::set_default_prec(), but the behavior remains.
I’d really appreciate your answer.
See you later!
Hi Abel.
Thank you for the good words!
You are right, in MPFR C++ precision of the result depends on the precision of function argument(s).
That’s because how C++ decomposes expressions into simple operations. For example, complex expression like y=(1+cos(x))/(1+sin(x)) C++ evaluates in several steps using temporal variables:
t1 = sin(x)
t2 = 1+t1
t3 = cos(x)
t4=1+t3
t5 = t4/t2
y = t5
On every step appropriate overloaded arithmetic operator is called. During these temporal evaluations we don’t know the precision of final variable – y. All we have is arguments.
Maybe it is the most natural choice to conduct operations with the information we only have – based on precision of arguments. Up to now MPFR C++ follows this strategy.
Generally speaking, by using any C++ interface to simplify usage of native MPFR user sacrifices flexibility of scrupulous precision control for every operation in the expression and have to rely on choices made in particular C++ wrapper.
MPFR C++ strategy works the best for the case when all variables (target, arguments) have the same precision, or target variable has lower precision than arguments. This is particularly useful for the applications where in order to obtain correct and accurate final result (even in low precision) intermediate calculations have to be done with high accuracy/precision. Examples are: numerical differentiation and digital filter convolution.
So the simplest workaround for your case is to use the same precision for all variables.
Another way is to redesign MPFR C++ interface so that all intermediate operations will be done with the precision adaptively selected based on precision of arguments and target variable. I’m planning on adding such feature in the future releases.
Hi Pavel!
Thank you for your wrapper, I find it extremely useful and intuitive.
Nevertheless, I have a problem. By running my code using the mpreal class and setting the precision to 64 bits via mpfr::mpreal::set_default_prec ( 64 ), I obtain different results from when I run it with standard C++ double precision (still 64 bits).
Is there an obvious reason for this behaviour, other than a possible mistake of mine?
Thank you for your attention,
Guido
P.S. The mpreal objects I use are created after calling mpfr::mpreal::set_default_prec ( 64 ).
Hi Guido!
Thank you for using MPFR C++!
There is no mistake in your code nor in MPFR C++.
Native MPFR uses different technique for floating number representation than double of IEEE-754 standard.
That is why MPFR cannot emulate standard double arithmetic by simply setting precision to 53 bits (mantissa of double occupies 53 bits, not 64 as you used).
Special not-obvious tricks are needed to do this: subnormalizing after each operation, etc. More detailed information you can find in MPFR manual (see pages 33-34): http://www.mpfr.org/mpfr-current/mpfr.pdf
Emulation of double arithmetic is not the goal of MPFR and these tricks require working with MPFR on low level, not through C++ wrapper.
Hi Pavel,
thank you very much for your prompt answer! I now understand why the results were different.
I would be grateful if you could answer another question: do you know if there is any arbitrary precision C++ wrapper/library that computes Sine and Cosine integrals – Si(x) & Ci(x) – efficiently, as defined in here http://en.wikipedia.org/wiki/Trigonometric_integral? I was thinking to just take GNU GSL source code for such functions and replace “double” with “mpreal”, but I guess it is not the best way to proceed. Another way could be to start from the Exponential Integral Ei(x) – which is included in MPFR – to obtain Si(x) and Ci(x) but the are complex numbers involved and I do not want to mess with them.
Thank you again,
Guido
Hi,
I don’t know any particular C++ multi-precision library which supports Si(x) and Ci(x).
Probably the best way is to write your own functions based on MPFR C++. As far as I see it is not very difficult. I would gladly add these functions to MPFR C++ if you will decide to share your source code then.
FYI GSL is already ported to native (not C++) MPFR: http://marcomaggi.github.com/docs/mpgsl.html
Maybe it supports what you are looking for.
Best regards,
Pavel.
Thank you for your answer.
Alas the MPFR port of GSL does not include the functions I am looking for. In the end I solved my problem in another way. I had to compute the difference between the cosine integral and sinc for very large arguments, and I needed high precision since Ci(x) ~ Sinc(x) when x is large.
However I found out that I could use series expansion. In fact it turned out that
CosIntegral[x] – Sinc[x] ~ (6/x^4 – 1/x^2) Cos[x] + (24/x^5 – 2/x^3 ) Sin[x].
Thank you anyway for the support!
Cheers,
Guido
Thanks for sharing the solution!
Dear Pavel,
thank you for writing this great wrapper and putting up the link to already compiled GMP and MPFR.
I would like to ask you one thing. I tried the exaple code, that comes with your wrapper and as soon as I move any variable out of the main(), I get an assertion error in ini2.c (about precision I guess).
So this code …
mpreal x = "420.968746359948568"; int main(int argc, char* argv[]){ mpreal::set_default_prec(128); ...generates the error, while this
int main(int argc, char* argv[]){ mpreal x = "420.968746359948568"; mpreal::set_default_prec(128); ...does not generate the error, even though in the second case the precision is set only after the init of X, which I though was be the problem in the first example. So, to sum it up, is there any way I can use global variables of mpreal type in my applications?
Thank you for a response in advance!
Daniel.
Hi Daniel!
Thank you for your comment!!
Current design of MPFR C++ doesn’t allow ‘mpreal’ number to be used as global variable. Generally, global variables are considered as bad solution in C++.
If you want them anyway you can to get around this restriction by using local static variables in global functions.
For details, please check Item 47 in Scott Meyers. Effective C++: 50 Specific Ways to Improve Your Programs and Design (2nd Edition)
.
Cheers,
Pavel.
Thank you for an answer, I know what statics are, but I will look into the book anyway:)
Hello Pavel,
me again. After trying to get the code rid of global variables, I stumbled upon a problem with linking static variables. The following code
#include "mpreal.h" using namespace mpfr; using namespace std; class X { public: static mpreal x; }; int main(int argc, char* argv[]){ X::x = "420.968746359948568"; return 0; };gives me a linker error telling me ‘error LNK2001: unresolved external symbol “public: static class mpfr::mpreal X::x” ‘. Can this be solved? And if so, could you please direct me to the topic or keyword, is used to solve this problem?
Thank you once again in advance,
Daniel B.
You should define static class member globally (not inside function):
using namespace mpfr; using namespace std; class X { public: static mpreal x; }; mpreal X::x = "420.968746359948568"; int main(int argc, char* argv[]){ return 0; };Please check this great site for explanation:
http://www.parashift.com/c++-faq-lite/ctors.html#faq-10.11
Dear Pavel,
thank you very much!
All the best,
Daniel.
Hello, thanks for your nice software!
I’d like to use it for MPACK (MBLAS/MLAPACK).
http://mplapack.sourceforge.net/
How about wrapping MPC as well?
I prepared general (not very well one) for GMP.
http://mplapack.cvs.sourceforge.net/viewvc/*checkout*/mplapack/mpack/include/mpc_class.h?revision=1.5
any comment is really appreciated.
Thanks!
Hi Maho!
Thank you for your comment!
I really like your ideas of MPACK and high precision optimization solvers.
I’ve checked your web site and I think you are doing very impressive stuff.
MPFR library offers very nice features over GMP, like correct rounding, high level functions, etc. I think it would suit MPACK needs much better than bare GMP.
Besides there are some problems with GMP usage on Windows platform: http://www.gladman.me.uk/. So, if you want MPACK to be used on Windows it is better to consider MPIR: http://www.mpir.org/. MPFR can use it instead of GMP.
C++ interface for GMP is already exist and distributed along with GMP source code (gmpxx.h). I didn’t check it personally, but maybe it worth considering before developing new one.
I have plans to make C++ interface for MPC but I cannot promise exact dates, sorry.
I would be honored if you would use my library MPFR C++ for MPACK. I am open for any questions and ready to provide you assistance on that matter.
Thank you!
Hi Pavel, thanks for your comment!
I have already been using GMP via gmpxx.h. Above mentioned complex C++ wrapper is
based on GMP. Problem here is not so many functions, as you suggest.
Thanks for MPIR. I didn’t know about it. MPACK infrastructure is designed to be easy to
port to a multiple precision library. What I need are four types; “REAL”, “COMPLEX”, “INTEGER” and “LOGICAL”. If someone make a wrapper to MPFR, I can port to MPFR immediately.
I’d like to use MPFR for better calculation evaluation.
Anyway please take your time. Hope I can do some contribution to you.
Best regards,
Nakata Maho
First I’d like to thank Pavel for the valuable wrapper he has written. It has saved many people lots of work! I’d also like to call everyone’s attention to another “templated” linear algebra library that I have been using with the mpreal type and with which I have been happy — with both MPFRC++ and the library:
http://math.nist.gov/tnt/index.html
It’s like a lighter weight version of Eigen, one more devoted to core linear algebra. Everything is defined in header files (hence no libraries to build) and it uses doxygen documentation. Note: I am using version 3.0.12 from the download page.
Jerry
Hi Pavel,
Really nice piece of work.
I have a couple of questions upon usage howerver. I am trying to assign a double value to an mp real and I get an exception.
mpreal::set_default_prec(128);
mpreal::set_double_bits(32); // I am assuming this is decimal places — please advice
mpreal d = 2.11 ;
terminate called after throwing an instance of ‘mpfr::conversion_overflow’
Is this correct? I would assume that I could assign a double to an mpreal.
Hi!
To convert double values to mpreal you can use two ways:
1. mpreal d = 2.11;
2. mpreal d = “2.11″; // preferable, since there is no conversion to double
Second way is preferable, since there is no transitive conversion to double. Also numbers way beyond double precision can be assigned, like
mpreal d = “2.1234567891234567891234567″;
Generally you do not need to invoke mpreal::set_double_bits(). It is needed only for specific purposes (read first comment for more information).
Thanks.
Hello
I’m trying to use mfprc++ to build a high precision code do a massive computation. I settled everything had the mpfrc++ test work successfully. When try to compile my code I bumbed into this error:
error: conversion from mpfr::mpreal to int is ambiguous
this happen because I’m trying to promote an integer to mpreal. Is there any function does this promotion? I’ll be very thankful for any help
Thank you!
Bassam
Hello
I’ve added type conversion operator for
mprealtointtransformation. Please download new version and try it – check if it solves your problem.Something like this should work fine
Pavel
Yes, it did work. Now the compiler (gcc) is accepting converting mpreal to int.
The package is very nice
Many thanks