GitHub-Proxy/eigen
1
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2025-08-14 04:35:57 +08:00
Code Issues Packages Projects Releases Wiki Activity

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

Browse Source 
(cherry picked from commit 54459214a1b9c67df04bc529474fca1ec9f4c84f)
This commit is contained in:
Chip Kerchner 2023-02-16 23:35:42 +00:00 committed by Antonio Sanchez
parent 079de53fa5
commit 6a4a0b66bd
3 changed files with 43 additions and 2 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

11
Eigen/src/Core/NumTraits.h
View File

@ -166,7 +166,16 @@ template<> struct NumTraits<double> : GenericNumTraits<double>
template<> struct NumTraits<long double>
: GenericNumTraits<long double>
{
static inline long double dummy_precision() { return 1e-15l; }
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
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> >

2
unsupported/Eigen/src/MatrixFunctions/MatrixPower.h
View File

@ -282,7 +282,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;
}

32
unsupported/test/cxx11_tensor_reduction.cpp
View File

@ -484,7 +484,39 @@ static void test_reduce_middle_dims() {
}
}
<<<<<<< HEAD
void test_cxx11_tensor_reduction() {
=======
template <typename ScalarType, int num_elements, int max_mean>
void test_sum_accuracy() {
Tensor<double, 1> double_tensor(num_elements);
Tensor<ScalarType, 1> tensor(num_elements);
for (double prescribed_mean = 0; prescribed_mean <= max_mean; prescribed_mean = numext::maxi(1.0, prescribed_mean*3.99)) {
// FIXME: NormalRandomGenerator doesn't work in bfloat and half.
double_tensor.setRandom<Eigen::internal::NormalRandomGenerator<double>>();
double_tensor += double_tensor.constant(prescribed_mean);
tensor = double_tensor.cast<ScalarType>();
Tensor<ScalarType, 0> sum;
sum = tensor.sum();
// Compute the reference value in double precsion.
double expected_sum = 0.0;
double abs_sum = 0.0;
for (int i = 0; i < num_elements; ++i) {
expected_sum += static_cast<double>(tensor(i));
abs_sum += static_cast<double>(numext::abs(tensor(i)));
}
// 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(static_cast<double>(num_elements)) * NumTraits<ScalarType>::epsilon() * static_cast<ScalarType>(abs_sum);
VERIFY_LE(err, tol);
}
}
EIGEN_DECLARE_TEST(cxx11_tensor_reduction) {
>>>>>>> 54459214a (Fix epsilon and dummy_precision values in long double for double doubles. Prevented some algorithms from converging on PPC.)
CALL_SUBTEST(test_trivial_reductions<ColMajor>());
CALL_SUBTEST(test_trivial_reductions<RowMajor>());
CALL_SUBTEST(test_simple_reductions<ColMajor>());