12#ifndef EIGEN_MPREALSUPPORT_MODULE_H
13#define EIGEN_MPREALSUPPORT_MODULE_H
15#include "../../Eigen/Core"
63 : GenericNumTraits<mpfr::mpreal>
69 RequireInitialization = 1,
75 typedef mpfr::mpreal Real;
76 typedef mpfr::mpreal NonInteger;
78 static inline Real highest (
long Precision = mpfr::mpreal::get_default_prec()) {
return mpfr::maxval(Precision); }
79 static inline Real lowest (
long Precision = mpfr::mpreal::get_default_prec()) {
return -mpfr::maxval(Precision); }
82 static inline Real Pi (
long Precision = mpfr::mpreal::get_default_prec()) {
return mpfr::const_pi(Precision); }
83 static inline Real Euler (
long Precision = mpfr::mpreal::get_default_prec()) {
return mpfr::const_euler(Precision); }
84 static inline Real Log2 (
long Precision = mpfr::mpreal::get_default_prec()) {
return mpfr::const_log2(Precision); }
85 static inline Real Catalan (
long Precision = mpfr::mpreal::get_default_prec()) {
return mpfr::const_catalan(Precision); }
87 static inline Real epsilon (
long Precision = mpfr::mpreal::get_default_prec()) {
return mpfr::machine_epsilon(Precision); }
88 static inline Real epsilon (
const Real& x) {
return mpfr::machine_epsilon(x); }
90#ifdef MPREAL_HAVE_DYNAMIC_STD_NUMERIC_LIMITS
91 static inline int digits10 (
long Precision = mpfr::mpreal::get_default_prec()) {
return std::numeric_limits<Real>::digits10(Precision); }
92 static inline int digits10 (
const Real& x) {
return std::numeric_limits<Real>::digits10(x); }
94 static inline int digits () {
return std::numeric_limits<Real>::digits(); }
95 static inline int digits (
const Real& x) {
return std::numeric_limits<Real>::digits(x); }
98 static inline Real dummy_precision()
100 mpfr_prec_t weak_prec = ((mpfr::mpreal::get_default_prec()-1) * 90) / 100;
101 return mpfr::machine_epsilon(weak_prec);
107 template<>
inline mpfr::mpreal random<mpfr::mpreal>()
109 return mpfr::random();
112 template<>
inline mpfr::mpreal random<mpfr::mpreal>(
const mpfr::mpreal& a,
const mpfr::mpreal& b)
114 return a + (b-a) * random<mpfr::mpreal>();
117 inline bool isMuchSmallerThan(
const mpfr::mpreal& a,
const mpfr::mpreal& b,
const mpfr::mpreal& eps)
119 return mpfr::abs(a) <= mpfr::abs(b) * eps;
122 inline bool isApprox(
const mpfr::mpreal& a,
const mpfr::mpreal& b,
const mpfr::mpreal& eps)
124 return mpfr::isEqualFuzzy(a,b,eps);
127 inline bool isApproxOrLessThan(
const mpfr::mpreal& a,
const mpfr::mpreal& b,
const mpfr::mpreal& eps)
129 return a <= b || mpfr::isEqualFuzzy(a,b,eps);
132 template<>
inline long double cast<mpfr::mpreal,long double>(
const mpfr::mpreal& x)
133 {
return x.toLDouble(); }
135 template<>
inline double cast<mpfr::mpreal,double>(
const mpfr::mpreal& x)
136 {
return x.toDouble(); }
138 template<>
inline long cast<mpfr::mpreal,long>(
const mpfr::mpreal& x)
139 {
return x.toLong(); }
141 template<>
inline int cast<mpfr::mpreal,int>(
const mpfr::mpreal& x)
142 {
return int(x.toLong()); }
147 class gebp_traits<mpfr::mpreal, mpfr::mpreal, false, false>
150 typedef mpfr::mpreal ResScalar;
152 Vectorizable =
false,
156 NumberOfRegisters = 1,
162 typedef ResScalar LhsPacket;
163 typedef ResScalar RhsPacket;
164 typedef ResScalar ResPacket;
165 typedef LhsPacket LhsPacket4Packing;
171 template<
typename Index,
typename DataMapper,
bool ConjugateLhs,
bool ConjugateRhs>
172 struct gebp_kernel<mpfr::mpreal,mpfr::mpreal,
Index,DataMapper,1,1,ConjugateLhs,ConjugateRhs>
174 typedef mpfr::mpreal mpreal;
177 void operator()(
const DataMapper& res,
const mpreal* blockA,
const mpreal* blockB,
178 Index rows, Index depth, Index cols,
const mpreal& alpha,
179 Index strideA=-1, Index strideB=-1, Index offsetA=0, Index offsetB=0)
181 if(rows==0 || cols==0 || depth==0)
184 mpreal acc1(0,mpfr_get_prec(blockA[0].mpfr_srcptr())),
185 tmp (0,mpfr_get_prec(blockA[0].mpfr_srcptr()));
187 if(strideA==-1) strideA = depth;
188 if(strideB==-1) strideB = depth;
190 for(Index i=0; i<rows; ++i)
192 for(Index j=0; j<cols; ++j)
194 const mpreal *A = blockA + i*strideA + offsetA;
195 const mpreal *B = blockB + j*strideB + offsetB;
198 for(Index k=0; k<depth; k++)
200 mpfr_mul(tmp.mpfr_ptr(), A[k].mpfr_srcptr(), B[k].mpfr_srcptr(), mpreal::get_default_rnd());
201 mpfr_add(acc1.mpfr_ptr(), acc1.mpfr_ptr(), tmp.mpfr_ptr(), mpreal::get_default_rnd());
204 mpfr_mul(acc1.mpfr_ptr(), acc1.mpfr_srcptr(), alpha.mpfr_srcptr(), mpreal::get_default_rnd());
205 mpfr_add(res(i,j).mpfr_ptr(), res(i,j).mpfr_srcptr(), acc1.mpfr_srcptr(), mpreal::get_default_rnd());
Namespace containing all symbols from the Eigen library.
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index