diff --git a/include/boost/numeric/ublas/detail/matrix_assign.hpp b/include/boost/numeric/ublas/detail/matrix_assign.hpp index 551c1dca..4e28b092 100644 --- a/include/boost/numeric/ublas/detail/matrix_assign.hpp +++ b/include/boost/numeric/ublas/detail/matrix_assign.hpp @@ -23,16 +23,26 @@ // Iterators based on ideas of Jeremy Siek namespace boost { namespace numeric { namespace ublas { +namespace detail { + + // Weak equality check - useful to compare equality two arbitary matrix expression results. + // Since the actual expressions are unknown, we check for and arbitary error bound + // on the relative error. + // For a linear expression the infinity norm makes sense as we do not know how the elements will be + // combined in the expression. False positive results are inevitable for arbirary expressions! + template + BOOST_UBLAS_INLINE + bool equals (const matrix_expression &e1, const matrix_expression &e2, S epsilon, S min_norm) { + return norm_inf (e1 - e2) < epsilon * + std::max (std::max (norm_inf (e1), norm_inf (e2)), min_norm); + } template BOOST_UBLAS_INLINE - bool equals (const matrix_expression &e1, const matrix_expression &e2) { + bool expression_type_check (const matrix_expression &e1, const matrix_expression &e2) { typedef typename type_traits::promote_type>::real_type real_type; - return norm_inf (e1 - e2) < BOOST_UBLAS_TYPE_CHECK_EPSILON * - std::max (std::max (norm_inf (e1), - norm_inf (e2)), - BOOST_UBLAS_TYPE_CHECK_MIN); + return equals (e1, e2, BOOST_UBLAS_TYPE_CHECK_EPSILON, BOOST_UBLAS_TYPE_CHECK_MIN); } @@ -244,6 +254,9 @@ namespace boost { namespace numeric { namespace ublas { m (index [k].first, index [k].second) = value_type/*zero*/(); } +}//namespace detail + + // Explicitly iterating row major template