Fix epsilon and dummy_precision values in long double for double doubles. Prevented some algorithms from converging on PPC.

This commit is contained in:
Chip Kerchner 2023-02-16 23:35:42 +00:00 committed by Rasmus Munk Larsen
parent a16fb889dd
commit 54459214a1
4 changed files with 17 additions and 7 deletions

View File

@ -246,8 +246,18 @@ template<> struct NumTraits<double> : GenericNumTraits<double>
template<> struct NumTraits<long double>
: GenericNumTraits<long double>
{
EIGEN_CONSTEXPR
static inline long double dummy_precision() { return 1e-15l; }
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR
static inline long double dummy_precision() { return static_cast<long double>(1e-15l); }
#if defined(EIGEN_ARCH_PPC) && (__LDBL_MANT_DIG__ == 106)
// PowerPC double double causes issues with some values
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR
static inline long double epsilon()
{
// 2^(-(__LDBL_MANT_DIG__)+1)
return static_cast<long double>(2.4651903288156618919116517665087e-32l);
}
#endif
};
template<typename Real_> struct NumTraits<std::complex<Real_> >

View File

@ -1109,7 +1109,7 @@ void set_from_triplets(const InputIterator& begin, const InputIterator& end, Spa
IndexMap outerIndexMap(mat.outerIndexPtr(), mat.outerSize() + 1);
for (InputIterator it(begin); it != end; ++it) {
eigen_assert(it->row() >= 0 && it->row() < mat.rows() && it->col() >= 0 && it->col() < mat.cols());
StorageIndex j = IsRowMajor ? it->row() : it->col();
StorageIndex j = static_cast<StorageIndex>(IsRowMajor ? it->row() : it->col());
outerIndexMap.coeffRef(j + 1)++;
}
@ -1124,8 +1124,8 @@ void set_from_triplets(const InputIterator& begin, const InputIterator& end, Spa
// push triplets to back of each inner vector
for (InputIterator it(begin); it != end; ++it) {
StorageIndex j = IsRowMajor ? it->row() : it->col();
StorageIndex i = IsRowMajor ? it->col() : it->row();
StorageIndex j = static_cast<StorageIndex>(IsRowMajor ? it->row() : it->col());
StorageIndex i = static_cast<StorageIndex>(IsRowMajor ? it->col() : it->row());
mat.data().index(back.coeff(j)) = i;
mat.data().value(back.coeff(j)) = it->value();
back.coeffRef(j)++;

View File

@ -283,7 +283,7 @@ inline int MatrixPowerAtomic<MatrixType>::getPadeDegree(long double normIminusT)
#endif
int degree = 3;
for (; degree <= maxPadeDegree; ++degree)
if (normIminusT <= maxNormForPade[degree - 3])
if (normIminusT <= static_cast<long double>(maxNormForPade[degree - 3]))
break;
return degree;
}

View File

@ -485,7 +485,7 @@ void test_sum_accuracy() {
// Test against probabilistic forward error bound. In reality, the error is much smaller
// when we use tree summation.
double err = Eigen::numext::abs(static_cast<double>(sum()) - expected_sum);
double tol = numext::sqrt(num_elements) * NumTraits<ScalarType>::epsilon() * static_cast<ScalarType>(abs_sum);
double tol = numext::sqrt(static_cast<double>(num_elements)) * NumTraits<ScalarType>::epsilon() * static_cast<ScalarType>(abs_sum);
VERIFY_LE(err, tol);
}
}