# How to perform an eigendecomposition of a general complex matrix with arbitrary precision in C/C++

I need to obtain the Eigenvectors of a general complex matrix, but with quadruple precision. Is anyone aware of a means to do this?

I currently use Tux Eigen, and I see that in their unsupported modules, they offer arbitrary precision using MPFR C++, however this appears to be for real numbers only. Perhaps this could be augmented to complex numbers in some way?

• You might try boost.multiprecision with Eigen. – user14717 Nov 1 '19 at 12:17
• @user14717 thanks for the suggestion, but boost::multiprecision appears not to play nicely with Eigen, unless there is a way that I am not aware of – AlexD Nov 1 '19 at 19:31
• I've never had a problem using boost::multiprecision with Eigen. – user14717 Nov 1 '19 at 19:48
• the simple lines; typedef Matrix<boost::multiprecision::complex128, Dynamic, Dynamic> MatrixXcmp; MatrixXcmp C = MatrixXcmp::Random(8,8); cout<<C<<endl; invoke several errors starting with; Eigen/src/Core/NumTraits.h:34:20: error: no matching function for call to ‘ceil(boost::multiprecision::number<boost::multiprecision::backends::complex_adaptor<boost::multiprecision::backends::float128_backend>, (boost::multiprecision::expression_template_option)0>)’ return int(ceil(-log10(NumTraits<Real>::epsilon()))); ~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ – AlexD Nov 1 '19 at 19:58
• using std::ceil? – user14717 Nov 1 '19 at 20:08

It transpires that it is reasonably straightforward to utilise MPFR C++ complex numbers within Tux Eigen;

#include <iostream>
#include <cmath>
#include <Eigen/Dense>
#include <Eigen/Eigenvalues>
#include <Eigen/LU>
#include <unsupported/Eigen/MPRealSupport>
#include <iomanip>

using namespace Eigen;
using namespace std;
using namespace mpfr;

int main(){
// set precision to 128 bits (double has only 53 bits)
mpreal::set_default_prec(128);
typedef Matrix<complex<mpreal>, Dynamic, Dynamic> MatrixXcmp;
MatrixXcmp C = MatrixXcmp::Random(8,8);
ComplexEigenSolver<MatrixXcmp> ces;
ces.compute(C);
return 0;
}