diff --git a/dune/common/densematrix.hh b/dune/common/densematrix.hh index d3ed73f489920949ed5624b04e5120a8f84a50d8..1725b8cc766e62ced1e1e2a499b1f4a9f347ff85 100644 --- a/dune/common/densematrix.hh +++ b/dune/common/densematrix.hh @@ -567,31 +567,69 @@ namespace Dune } //! infinity norm (row sum norm, how to generalize for blocks?) - typename FieldTraits<value_type>::real_type infinity_norm () const - { - if (size() == 0) - return 0.0; - - ConstIterator it = begin(); - typename remove_const< typename FieldTraits<value_type>::real_type >::type max = (*it).one_norm(); - for (it = it + 1; it != end(); ++it) - max = std::max(max, (*it).one_norm()); - - return max; + template <typename vt = value_type, + typename std::enable_if<!has_nan<vt>::value, int>::type = 0> + typename FieldTraits<vt>::real_type infinity_norm() const { + using real_type = typename FieldTraits<vt>::real_type; + using std::max; + + real_type norm = 0.0; + for (auto&& x : *this) { + real_type const a = x.one_norm(); + norm = max(a, norm); + } + return norm; } //! simplified infinity norm (uses Manhattan norm for complex values) - typename FieldTraits<value_type>::real_type infinity_norm_real () const - { - if (size() == 0) - return 0.0; + template <typename vt = value_type, + typename std::enable_if<!has_nan<vt>::value, int>::type = 0> + typename FieldTraits<vt>::real_type infinity_norm_real() const { + using real_type = typename FieldTraits<vt>::real_type; + using std::max; + + real_type norm = 0.0; + for (auto&& x : *this) { + real_type const a = x.one_norm_real(); + norm = max(a, norm); + } + return norm; + } - ConstIterator it = begin(); - typename remove_const< typename FieldTraits<value_type>::real_type >::type max = (*it).one_norm_real(); - for (it = it + 1; it != end(); ++it) - max = std::max(max, (*it).one_norm_real()); + //! infinity norm (row sum norm, how to generalize for blocks?) + template <typename vt = value_type, + typename std::enable_if<has_nan<vt>::value, int>::type = 0> + typename FieldTraits<vt>::real_type infinity_norm() const { + using real_type = typename FieldTraits<vt>::real_type; + using std::max; + + real_type norm = 0.0; + real_type isNaN = 1.0; + for (auto&& x : *this) { + real_type const a = x.one_norm(); + norm = max(a, norm); + isNaN += a; + } + isNaN /= isNaN; + return norm * isNaN; + } - return max; + //! simplified infinity norm (uses Manhattan norm for complex values) + template <typename vt = value_type, + typename std::enable_if<has_nan<vt>::value, int>::type = 0> + typename FieldTraits<vt>::real_type infinity_norm_real() const { + using real_type = typename FieldTraits<vt>::real_type; + using std::max; + + real_type norm = 0.0; + real_type isNaN = 1.0; + for (auto&& x : *this) { + real_type const a = x.one_norm_real(); + norm = max(a, norm); + isNaN += a; + } + isNaN /= isNaN; + return norm * isNaN; } //===== solve