How to compare vectors approximately in Eigen?

Question

Is there a function in Eigen to compare vectors (matrices) using both relative and absolute tolerance aka numpy.allclose? Standard isApprox fails if one of the vectors is very close to zero.

Syconium

See eigen.tuxfamily.org/dox-2.0/TutorialCore.html – Highlight

Answer 1

There is no built-in function implementing numpy.allclose, but you easily write one yourself if that's really what you need. However, I'd rather suggest the use of isMuchSmallerThan with reference value:

(a-b).isMuchSmallerThan(ref)

where ref is a representative non zero for your problem.

EDIT: for reference here is a possible implementation of allclose:

template<typename DerivedA, typename DerivedB>
bool allclose(const Eigen::DenseBase<DerivedA>& a,
              const Eigen::DenseBase<DerivedB>& b,
              const typename DerivedA::RealScalar& rtol
                  = Eigen::NumTraits<typename DerivedA::RealScalar>::dummy_precision(),
              const typename DerivedA::RealScalar& atol
                  = Eigen::NumTraits<typename DerivedA::RealScalar>::epsilon())
{
  return ((a.derived() - b.derived()).array().abs()
          <= (atol + rtol * b.derived().array().abs())).all();
}

Orthman

isMuchSmallerThan can be used for absolute comparison, i.e. (a-b).isMuchSmallerThan(1.0, atol) is equivalent to np.allclose(a, b, 0.0, atol), so to mimic np.allclose we have to do something like this: (a-b).isMuchSmallerThan(1.0, atol) || a.isApprox(b, rtol). Am I correct? – Syconium

Not exactly because isMuchSmallerThan and isApprox are based on the L2 matrix norm and not element-wise comparisons (infinite norm). – Orthman

It's probably worth mentioning that this, like numpy's allclose, is not symmetric, and b is expected to be a reference value. If you need a symmetric comparison, multiply rtol by a.derived().array().abs().max(b.derived().array().abs()) instead. – Utah

For Quaternions which have a coeffs() method to get at the matrix and where negatives are equivelent: auto e{ atol + rtol * b.coeffs().derived().array().abs() }; return ((a.coeffs().derived() - b.coeffs().derived()).array().abs() <= e).all() || ((a.coeffs().derived() + b.coeffs().derived()).array().abs() <= e).all(); – Pressley

Answer 2

There is also isApprox function which was not working for me. I am just using ( expect - res).norm() < some small number.

Multipara

Answer 3

Here I figured out a version for eigen tensor:

template <typename DerivedA, typename DerivedB>
bool allclose(const Eigen::TensorBase<DerivedA, Eigen::WriteAccessors> &tensorA,
              const Eigen::TensorBase<DerivedB, Eigen::WriteAccessors> &tensorB,
              const typename DerivedA::Scalar &relativeTolerance = Eigen::NumTraits<typename DerivedA::Scalar>::dummy_precision(),
              const typename DerivedA::Scalar &absoluteTolerance = Eigen::NumTraits<typename DerivedA::Scalar>::epsilon())
{
    auto absoluteDifference = (tensorA.eval() - tensorB.eval()).abs();
    auto toleranceBound = (relativeTolerance * tensorB.abs() + absoluteTolerance);
    constexpr int layout = DerivedA::Layout == Eigen::RowMajor ? Eigen::RowMajor : Eigen::ColMajor;
    Eigen::Tensor<bool, 0, layout> isClose = (toleranceBound - absoluteDifference).unaryExpr([](auto x) { return x >= 0; }).all();
    return isClose(0);
}

Nacreous