From 78386c0dd9d4c600f9beaf8716f983f353db59dc Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Mon, 26 Sep 2022 10:24:27 +0200 Subject: [PATCH] Fix exponential complexity of FloatingPointComparator with nesting --- .../deal.II/base/floating_point_comparator.h | 180 +++++++++++------- 1 file changed, 107 insertions(+), 73 deletions(-) diff --git a/include/deal.II/base/floating_point_comparator.h b/include/deal.II/base/floating_point_comparator.h index 2a7f64ca9f..eef3f96b46 100644 --- a/include/deal.II/base/floating_point_comparator.h +++ b/include/deal.II/base/floating_point_comparator.h @@ -48,6 +48,18 @@ struct FloatingPointComparator static constexpr std::size_t width = dealii::internal::VectorizedArrayTrait::width; + /** + * An enum to decode whether a particular comparison returns less, larger, + * or equal. This is needed beyond the regular operator() because of + * efficiency reasons in nested checks as supported by this class. + */ + enum class ComparisonResult + { + less, + equal, + greater + }; + /** * Constructor. */ @@ -60,49 +72,58 @@ struct FloatingPointComparator FloatingPointComparator(FloatingPointComparator &&rhs) noexcept = default; + /** + * Compare two objects of the same type and return `true` in case the first + * object is less than the second one. This function calls the compare() + * functions below, so only types supported by compare() will work with this + * function. + */ + template + bool + operator()(const T &object1, const T &object2) const; + /** * Compare two vectors of numbers (not necessarily of the same length), * where vectors of different lengths are first sorted by their length and * then by the entries. */ template - bool - operator()(const std::vector &v1, const std::vector &v2) const; + ComparisonResult + compare(const std::vector &v1, const std::vector &v2) const; /** * Compare two arrays. */ template - bool - operator()(const std::array &t1, const std::array &t2) const; + ComparisonResult + compare(const std::array &t1, const std::array &t2) const; /** * Compare two tensors. */ template - bool - operator()(const Tensor &t1, - const Tensor &t2) const; + ComparisonResult + compare(const Tensor &t1, const Tensor &t2) const; /** * Compare two tables. */ template - bool - operator()(const Table<2, T> &t1, const Table<2, T> &t2) const; + ComparisonResult + compare(const Table<2, T> &t1, const Table<2, T> &t2) const; /** * Compare two scalar numbers. */ - bool - operator()(const ScalarNumber &s1, const ScalarNumber &s2) const; + ComparisonResult + compare(const ScalarNumber &s1, const ScalarNumber &s2) const; /** - * Compare two VectorizedArray instances. + * Compare two VectorizedArray instances. */ - bool - operator()(const VectorizedArray &v1, - const VectorizedArray &v2) const; + ComparisonResult + compare(const VectorizedArray &v1, + const VectorizedArray &v2) const; private: const ScalarNumber tolerance; @@ -122,123 +143,136 @@ FloatingPointComparator::FloatingPointComparator( : tolerance(tolerance) , use_absolute_tolerance(use_absolute_tolerance) , mask(mask) -{} +{ + Assert(mask.count() > 0, ExcMessage("No component selected")); +} template template bool -FloatingPointComparator::operator()(const std::vector &v1, - const std::vector &v2) const +FloatingPointComparator::operator()(const T &object1, + const T &object2) const +{ + return compare(object1, object2) == ComparisonResult::less; +} + + + +template +template +typename FloatingPointComparator::ComparisonResult +FloatingPointComparator::compare(const std::vector &v1, + const std::vector &v2) const { const unsigned int s1 = v1.size(), s2 = v2.size(); if (s1 < s2) - return true; + return ComparisonResult::less; else if (s1 > s2) - return false; + return ComparisonResult::greater; else for (unsigned int i = 0; i < s1; ++i) - if (this->operator()(v1[i], v2[i])) - return true; - else if (this->operator()(v2[i], v1[i])) - return false; - return false; + { + const ComparisonResult result = compare(v1[i], v2[i]); + if (result != ComparisonResult::equal) + return result; + } + return ComparisonResult::equal; } template template -bool -FloatingPointComparator::operator()(const std::array &t1, - const std::array &t2) const +typename FloatingPointComparator::ComparisonResult +FloatingPointComparator::compare(const std::array &t1, + const std::array &t2) const { for (unsigned int i = 0; i < t1.size(); ++i) - if (this->operator()(t1[i], t2[i])) - return true; - else if (this->operator()(t2[i], t1[i])) - return false; - return false; + { + const ComparisonResult result = compare(t1[i], t2[i]); + if (result != ComparisonResult::equal) + return result; + } + return ComparisonResult::equal; } template template -bool -FloatingPointComparator::operator()( - const Tensor &t1, - const Tensor &t2) const +typename FloatingPointComparator::ComparisonResult +FloatingPointComparator::compare(const Tensor &t1, + const Tensor &t2) const { - for (unsigned int k = 0; k < dim; ++k) - if (this->operator()(t1[k], t2[k])) - return true; - else if (this->operator()(t2[k], t1[k])) - return false; - return false; + for (unsigned int i = 0; i < dim; ++i) + { + const ComparisonResult result = compare(t1[i], t2[i]); + if (result != ComparisonResult::equal) + return result; + } + return ComparisonResult::equal; } template template -bool -FloatingPointComparator::operator()(const Table<2, T> &t1, - const Table<2, T> &t2) const +typename FloatingPointComparator::ComparisonResult +FloatingPointComparator::compare(const Table<2, T> &t1, + const Table<2, T> &t2) const { AssertDimension(t1.size(0), t2.size(0)); AssertDimension(t1.size(1), t2.size(1)); for (unsigned int i = 0; i < t1.size(0); ++i) for (unsigned int j = 0; j < t1.size(1); ++j) - if (this->operator()(t1[i][j], t2[i][j])) - return true; - else if (this->operator()(t2[i][j], t1[i][j])) - return false; - return false; + { + const ComparisonResult result = compare(t1[i][j], t2[i][j]); + if (result != ComparisonResult::equal) + return result; + } + return ComparisonResult::equal; } + + template -bool -FloatingPointComparator::operator()(const ScalarNumber &s1, - const ScalarNumber &s2) const +typename FloatingPointComparator::ComparisonResult +FloatingPointComparator::compare(const ScalarNumber &s1, + const ScalarNumber &s2) const { - if (mask[0]) + if (width == 1 || mask[0]) { const ScalarNumber tolerance = use_absolute_tolerance ? this->tolerance : (std::abs(s1 + s2) * this->tolerance); - if ((s1 < s2 - tolerance)) - return true; - else - return false; + if (s1 < s2 - tolerance) + return ComparisonResult::less; + else if (s1 > s2 + tolerance) + return ComparisonResult::greater; } - return false; + return ComparisonResult::equal; } template -bool -FloatingPointComparator::operator()( +typename FloatingPointComparator::ComparisonResult +FloatingPointComparator::compare( const VectorizedArray &v1, const VectorizedArray &v2) const { - for (unsigned int v = 0; v < width; ++v) - if (mask[v]) + for (unsigned int i = 0; i < width; ++i) + if (mask[i]) { - const ScalarNumber tolerance = - use_absolute_tolerance ? this->tolerance : - (std::abs(v1[v] + v2[v]) * this->tolerance); - - if (v1[v] < v2[v] - tolerance) - return true; - if (v1[v] > v2[v] + tolerance) - return false; + const ComparisonResult result = compare(v1[i], v2[i]); + if (result != ComparisonResult::equal) + return result; } - return false; + return ComparisonResult::equal; } -- 2.39.5