From 44fc18c4a657f16f1a3d603ff9db0efa7194a795 Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Thu, 10 Dec 2015 19:56:22 +0100
Subject: [PATCH] MultiTypeBlockVector: Fix NaN behaviour of infinity_norm()
---
dune/istl/multitypeblockvector.hh | 37 ++++++++++++++++++++++++++++---
1 file changed, 34 insertions(+), 3 deletions(-)
diff --git a/dune/istl/multitypeblockvector.hh b/dune/istl/multitypeblockvector.hh
index dfc72d2fe..90fa3c874 100644
--- a/dune/istl/multitypeblockvector.hh
+++ b/dune/istl/multitypeblockvector.hh
@@ -233,12 +233,38 @@ namespace Dune {
Each element of the vector has to provide the method "infinity_norm()"
in order to calculate the whole vector's norm.
*/
- template<int count, typename T>
+ template<int count, typename T,
+ bool hasNaN = has_nan<typename T::field_type>::value>
class MultiTypeBlockVector_InfinityNorm {
public:
typedef typename T::field_type field_type;
typedef typename FieldTraits<field_type>::real_type real_type;
+ /**
+ * Take the maximum over all elements' norms
+ */
+ static real_type result (const T& a)
+ {
+ std::pair<real_type, real_type> const ret = resultHelper(a);
+ return ret.first * (ret.second / ret.second);
+ }
+
+ // Computes the maximum while keeping track of NaN values
+ static std::pair<real_type, real_type> resultHelper(const T& a) {
+ using std::max;
+ auto const rest =
+ MultiTypeBlockVector_InfinityNorm<count - 1, T>::resultHelper(a);
+ real_type const norm = std::get<count - 1>(a).infinity_norm();
+ return {max(norm, rest.first), norm + rest.second};
+ }
+ };
+
+ template <int count, typename T>
+ class MultiTypeBlockVector_InfinityNorm<count, T, false> {
+ public:
+ typedef typename T::field_type field_type;
+ typedef typename FieldTraits<field_type>::real_type real_type;
+
/**
* Take the maximum over all elements' norms
*/
@@ -249,8 +275,8 @@ namespace Dune {
}
};
- template<typename T> //recursion end
- class MultiTypeBlockVector_InfinityNorm<0,T> {
+ template<typename T, bool hasNaN> //recursion end
+ class MultiTypeBlockVector_InfinityNorm<0,T, hasNaN> {
public:
typedef typename T::field_type field_type;
typedef typename FieldTraits<field_type>::real_type real_type;
@@ -258,6 +284,11 @@ namespace Dune {
{
return 0.0;
}
+
+ static std::pair<real_type,real_type> resultHelper (const T&)
+ {
+ return {0.0, 1.0};
+ }
};
/**
--
GitLab