From c54f10479d054d9173ad8a11a2e55e121b3ad452 Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Wed, 13 Jul 2016 12:30:08 +0200
Subject: [PATCH] Add hasIgnore() to CanIgnore; use it
In particular, a CGStep will now, like a few other steps/solvers,
treat a missing ignore the same way as an ignore whose every entry is
false.
---
dune/solvers/common/canignore.hh | 7 +++++++
dune/solvers/iterationsteps/blockgsstep.cc | 2 +-
dune/solvers/iterationsteps/blockgssteps.hh | 2 +-
dune/solvers/iterationsteps/cgstep.cc | 16 ++++++++--------
dune/solvers/iterationsteps/trustregiongsstep.cc | 2 +-
dune/solvers/solvers/umfpacksolver.hh | 2 +-
6 files changed, 19 insertions(+), 12 deletions(-)
diff --git a/dune/solvers/common/canignore.hh b/dune/solvers/common/canignore.hh
index 4464ac0..019e5f0 100644
--- a/dune/solvers/common/canignore.hh
+++ b/dune/solvers/common/canignore.hh
@@ -47,6 +47,13 @@ public:
ignoreNodes_ = &i;
}
+ /**
+ * \brief Only if this returns true can ignore() safely be called.
+ */
+ bool hasIgnore() const {
+ return ignoreNodes_ != nullptr;
+ }
+
/**
* \brief Const access to bit vector
*/
diff --git a/dune/solvers/iterationsteps/blockgsstep.cc b/dune/solvers/iterationsteps/blockgsstep.cc
index 284351e..d49c588 100644
--- a/dune/solvers/iterationsteps/blockgsstep.cc
+++ b/dune/solvers/iterationsteps/blockgsstep.cc
@@ -24,7 +24,7 @@ template<class MatrixType, class DiscFuncType, class BitVectorType>
inline
void BlockGSStep<MatrixType, DiscFuncType, BitVectorType>::iterate()
{
- assert(this->ignoreNodes_ != nullptr);
+ assert(this->hasIgnore());
if (gs_type_ != Direction::BACKWARD)
for (std::size_t i=0; i<this->x_->size(); i++)
diff --git a/dune/solvers/iterationsteps/blockgssteps.hh b/dune/solvers/iterationsteps/blockgssteps.hh
index aa995c0..c8db9a2 100644
--- a/dune/solvers/iterationsteps/blockgssteps.hh
+++ b/dune/solvers/iterationsteps/blockgssteps.hh
@@ -310,7 +310,7 @@ struct BlockGSStep : public LinearIterationStep<Matrix, Vector, BitVector> {
, direction_(direction) {}
void iterate() {
- assert(this->ignoreNodes_ != nullptr);
+ assert(this->hasIgnore());
assert(this->mat_ != nullptr);
assert(this->x_ != nullptr);
assert(this->rhs_ != nullptr);
diff --git a/dune/solvers/iterationsteps/cgstep.cc b/dune/solvers/iterationsteps/cgstep.cc
index 4e00a3b..3085b2b 100644
--- a/dune/solvers/iterationsteps/cgstep.cc
+++ b/dune/solvers/iterationsteps/cgstep.cc
@@ -4,8 +4,6 @@
template <class MatrixType, class VectorType>
void CGStep<MatrixType, VectorType>::check() const
{
- if (this->ignoreNodes_ == nullptr)
- DUNE_THROW(SolverError, "ignoreNodes_ member not set");
}
template <class MatrixType, class VectorType>
@@ -13,9 +11,10 @@ void CGStep<MatrixType, VectorType>::preprocess()
{
// Compute the residual (r starts out as the rhs)
matrix_->mmv(*x_,r_);
- for (size_t i=0; i < r_.size(); ++i)
- if (this->ignore()[i].any())
- r_[i] = 0;
+ if (this->hasIgnore())
+ for (size_t i=0; i < r_.size(); ++i)
+ if (this->ignore()[i].any())
+ r_[i] = 0;
if (preconditioner_) {
preconditioner_->setMatrix(*matrix_);
@@ -36,9 +35,10 @@ void CGStep<MatrixType, VectorType>::iterate()
x_->axpy(alpha, p_); // x_1 = x_0 + alpha_0 p_0
r_.axpy(-alpha, q); // r_1 = r_0 - alpha_0 Ap_0
- for (size_t i=0; i < r_.size(); ++i)
- if (this->ignore()[i].any())
- r_[i] = 0;
+ if (this->hasIgnore())
+ for (size_t i=0; i < r_.size(); ++i)
+ if (this->ignore()[i].any())
+ r_[i] = 0;
if (preconditioner_)
preconditioner_->apply(q, r_);
diff --git a/dune/solvers/iterationsteps/trustregiongsstep.cc b/dune/solvers/iterationsteps/trustregiongsstep.cc
index 8fd542e..f7f7208 100644
--- a/dune/solvers/iterationsteps/trustregiongsstep.cc
+++ b/dune/solvers/iterationsteps/trustregiongsstep.cc
@@ -13,7 +13,7 @@ void TrustRegionGSStep<MatrixType, VectorType>::iterate()
// Dirichlet nodes in this block
std::bitset<BlockSize> ignoreNodes(0);
- if (this->ignoreNodes_ != nullptr)
+ if (this->hasIgnore())
ignoreNodes = this->ignore()[i];
if (ignoreNodes.all())
diff --git a/dune/solvers/solvers/umfpacksolver.hh b/dune/solvers/solvers/umfpacksolver.hh
index f2a36ac..a89bfeb 100644
--- a/dune/solvers/solvers/umfpacksolver.hh
+++ b/dune/solvers/solvers/umfpacksolver.hh
@@ -70,7 +70,7 @@ public:
// We may use the original rhs, but ISTL modifies it, so we need a non-const type here
VectorType mutableRhs = *rhs_;
- if (this->ignoreNodes_ == nullptr)
+ if (not this->hasIgnore())
{
/////////////////////////////////////////////////////////////////
// Solve the system
--
GitLab