Power function for mpf_class numbers in gmpxx

gmp

gmp library provides the function void mpf_pow_ui (mpf_t rop, const mpf_t op1, unsigned long int op2) to raise op1 to the power op2 (according to https://gmplib.org/manual/Float-Arithmetic#index-Powering-functions-1).

But the documentation seems to say nothing about it in the c++ interface. I've tried with names such as pow, pow_iu, power but none of them are defined.

Is there the way to raise a float to an exponent (either float or integer) using gmpxx?

gmpxx.h contains interfaces to some mathematical operations like sqrt (see line 3341)

__GMP_DEFINE_UNARY_FUNCTION_1(mpf_t, sqrt, __gmp_sqrt_function)

where the caps-named macro depends on a general evaluation class template named __gmp_expr<T, U>, which in turn allows for the eval of expressions with mixed native and arbitray-precision types, and that's what makes the C++ GMP interface much easier to use.

But there's no such definition for pow. As Marc Glisse pointed out, you must convert your C++ objects to the C types using the mpf_class::get_mpf_t function.

The following is an example code. Name it, say, test.cpp.

#include <iostream>
#include <gmpxx.h>
using namespace std;

int main (void) {
    mpf_class a = 12.3;
    unsigned long int b = 123UL;
    mpf_t c;
    mpf_set_default_prec(100000);
    mpf_init(c);
    mpf_pow_ui(c, a.get_mpf_t(), b);
    gmp_printf("c = %.50Ff\n", c);
    return 0;
}

And compile with

g++ test.cpp -o test -lgmpxx -lgmp

Which produces the output:

114374367934618002778643226182707594198913258409535335775583252201365538178632825702225459029661601216944929436371688246107986574246790.32099077871758646985223686110515186972735931183764

Unfortunately, even when increasing mpf precision, I can't make this agree further than 14 digits with the answer from WolframAlpha:

114374367934617190099880295228066276746218078451850229775887975052369504785666896446606568365201542169649974727730628842345343196581134.89591994282087444983721209947664895835902379607854

Note that if you want to raise an arbitrary-precision float to an arbitrary-precision float you'll need to use the function mpfr_pow from the MPFR library, which by the way is recommended library to handle arbitrary-precision float:

New projects should consider using the GMP extension library MPFR (http://mpfr.org) instead. MPFR provides well-defined precision and accurate rounding, and thereby naturally extends IEEE P754.

EDIT: Due to the discrepancy between GMP and WolframAlpha (which, by the way, uses GMP internally) I've posted this question.

EDIT2: As mentioned in this and in this comments, the discrepancy is expected, since when using the mpf_pow_ui function the compiler converts 12.3 into a double, which is not exactly representable in binary, while Mathematica uses arbitrary precision for that value, so its more accurate in this particular case.

EDIT3: GMP can actually match the result from WolframAlpha, as explained by John Bollinger in his answer.

The issue with my original code is in using the mpf_set_d function to set the value 12.3, because it converts to double and thus looses precision. John modified the code to use the mpf_set_str function instead, which, converted to C++, becomes:

#include <iostream>
#include <gmpxx.h>
using namespace std;

int main (void) {
    mpf_class a("12.3",2000);
    unsigned long int b = 123UL;
    mpf_t c;
    mpf_set_default_prec(2000);
    mpf_init(c);
    mpf_pow_ui(c, a.get_mpf_t(), b);
    gmp_printf("c = %.50Ff\n", c);
    return 0;
}

And this outputs the fully correct answer (with last digit rounded to '5' since that digit was a '4' followed by '9'):

114374367934617190099880295228066276746218078451850229775887975052369504785666896446606568365201542169649974727730628842345343196581134.89591994282087444983721209947664895835902379607855

Power function for mpf_class numbers in gmpxx

Related

Recent Posts