From f0d01da2293c8893b33156352e89810c5d368809 Mon Sep 17 00:00:00 2001
From: Oliver Sander <sander@igpm.rwth-aachen.de>
Date: Sun, 1 Feb 2015 09:21:38 +0000
Subject: [PATCH] Use template technique to test for P1NodalBasis, not the
preprocessor
The multigrid hierarchy is built in two different ways, depending on whether
we use a P1NodalBasis for the discretization, or anything else. Previously,
this switch was made using the preprocessor. The patches changes the code
to use std::is_same instead.
[[Imported from SVN: r10046]]
---
dune/gfe/riemanniantrsolver.cc | 33 ++++++++++++++++++++-------------
1 file changed, 20 insertions(+), 13 deletions(-)
diff --git a/dune/gfe/riemanniantrsolver.cc b/dune/gfe/riemanniantrsolver.cc
index 43152330..bf10e030 100644
--- a/dune/gfe/riemanniantrsolver.cc
+++ b/dune/gfe/riemanniantrsolver.cc
@@ -177,13 +177,24 @@ setup(const GridType& grid,
for (size_t k=0; k<mmgStep->mgTransfer_.size(); k++)
delete(mmgStep->mgTransfer_[k]);
-#if defined THIRD_ORDER || defined SECOND_ORDER
- mmgStep->mgTransfer_.resize(numLevels);
-#else
- mmgStep->mgTransfer_.resize(numLevels-1);
-#endif
-
-#if defined THIRD_ORDER || defined SECOND_ORDER
+ ////////////////////////////////////////////////////////////////////////
+ // The P1 space (actually P1/Q1, depending on the grid) is special:
+ // If we work in such a space, then the multigrid hierarchy of spaces
+ // is constructed in the usual way. For all other space, there is
+ // an additional restriction operator on the top of the hierarchy, which
+ // restricts the FE space to the P1/Q1 space on the same grid.
+ // On the lower grid levels a hierarchy of P1/Q1 spaces is used again.
+ ////////////////////////////////////////////////////////////////////////
+
+ bool isP1Basis = std::is_same<Basis,P1NodalBasis<typename Basis::GridView> >::value;
+
+ if (isP1Basis)
+ mmgStep->mgTransfer_.resize(numLevels-1);
+ else
+ mmgStep->mgTransfer_.resize(numLevels);
+
+ // Here we set up the restriction of the leaf grid space into the leaf grid P1/Q1 space
+ if (not isP1Basis)
{
typedef typename TruncatedCompressedMGTransfer<CorrectionType>::TransferOperatorType TransferOperatorType;
P1NodalBasis<typename GridType::LeafGridView,double> p1Basis(grid_->leafGridView());
@@ -213,13 +224,9 @@ setup(const GridType& grid,
Dune::shared_ptr<TransferOperatorType> topTransferOperator = Dune::make_shared<TransferOperatorType>(matrixComm.reduceCopy(pkToP1TransferMatrix));
dynamic_cast<TruncatedCompressedMGTransfer<CorrectionType>*>(mmgStep->mgTransfer_.back())->setMatrix(topTransferOperator);
}
-#endif
-#if defined THIRD_ORDER || defined SECOND_ORDER
- for (int i=0; i<mmgStep->mgTransfer_.size()-1; i++) {
-#else
- for (size_t i=0; i<mmgStep->mgTransfer_.size(); i++) {
-#endif
+ // Now the P1/Q1 restriction operators for the remaining levels
+ for (int i=0; i<numLevels-1; i++) {
// Construct the local multigrid transfer matrix
TruncatedCompressedMGTransfer<CorrectionType>* newTransferOp = new TruncatedCompressedMGTransfer<CorrectionType>;
--
GitLab