How to define function with arbitrary precision (Eigen/MPRealSupport)

How to define Matrix<mpreal, Dynamic, Dynamic> with typedef?

We usually put the defined type in a header file, but it needs mpreal::set_default_prec(256); and that's the problem.

I am quite new to MPFR, so my apology if it seems easy. So what is your idea?

#include <Eigen/LU>
#include <iostream>
#include <Eigen/Dense>
#include "mylibrary.hpp"
#include <unsupported/Eigen/MPRealSupport>

using namespace mpfr;
using namespace Eigen;


MatrixXmp hilbert_matrix(const int size)
{
    MatrixXmp A = MatrixXmp::Zero(size, size);
    for (int i = 1; i < size + 1; ++i)
        for (int j = 1; j < size + 1; ++j)
            A(i - 1, j - 1) = 1. / (i + j - 1.);

    return A;
}


int main()
{
    // set precision to 256 bits (double has only 53 bits)

    mpreal::set_default_prec(256);
    typedef Matrix<mpreal, Dynamic, Dynamic> MatrixXmp;
    typedef Matrix<mpreal, Dynamic, 1> VectorXmp;
    

    return 0;
}

Error:

‘MatrixXmp’ does not name a type

Solution

mpreal (attention, it is located inside a namespace mpfr) is a class and set_default_prec is a static method. This means that setting the precision has no effect on the template parameter and therefore also not on the typedef. Calling the static function will change the default precision only internally without having effects on the template.

You can typedef outside of your main inside your header and only then set the precision inside the main of the program which calls it.

 using MatrixXmp = Matrix<mpfr::mpreal, Dynamic, Dynamic>;
 using VectorXmp = Matrix<mpfr::mpreal, Dynamic, 1>;

 int main() {
   mpfr::mpreal::set_default_prec(256);
   return EXIT_SUCCESS;
 }

Alternatively there is a constructor
```
mpreal(const mpz_t u, mp_prec_t prec = mpreal::get_default_prec(), mp_rnd_t mode = mpreal::get_default_rnd());
```
This means you could actually add a mp_prec_t prec as a (template) argument to your function and default it to 256. Then when constructing any sort of mpfr::mpreal pass the precision and potentially also the rounding policy as second and third argument to the constructor.