From 223fb8f67cc9d7394e74785de660f1e42e6c3cba Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Carsten=20Gr=C3=A4ser?= <graeser@mi.fu-berlin.de>
Date: Tue, 10 Apr 2012 21:38:17 +0000
Subject: [PATCH] Implement (algebraic) nested iteration
This purely algebraic variant of nested iteration
can be used to compute initial iterates.
[[Imported from SVN: r5920]]
---
.../iterationsteps/obstacletnnmgstep.hh | 78 +++++++++++++++++++
1 file changed, 78 insertions(+)
diff --git a/dune/solvers/iterationsteps/obstacletnnmgstep.hh b/dune/solvers/iterationsteps/obstacletnnmgstep.hh
index 79990283..2b36ca51 100644
--- a/dune/solvers/iterationsteps/obstacletnnmgstep.hh
+++ b/dune/solvers/iterationsteps/obstacletnnmgstep.hh
@@ -17,6 +17,7 @@
#include <dune/solvers/solvers/loopsolver.hh>
#include <dune/solvers/norms/energynorm.hh>
#include <dune/solvers/common/boxconstraint.hh>
+#include <dune/solvers/transferoperators/mandelobsrestrictor.hh>
@@ -335,6 +336,83 @@ class ObstacleTNNMGStep
return s;
}
+
+ /**
+ * \brief Compute initial iterate by nested iteration
+ *
+ * This method performs a purely algebraic variant of nested iteration.
+ * The coarse matrix and rhs are obtained using the transfer operators.
+ * Coarse obstacles are constructed using a local max/min over the
+ * support of coarse basis functions as suggested by J. Mandel.
+ * The ignore-marks are restricted such that a coarse dof
+ * is ignored if it is associated to an ignored fine dof
+ * in the sense that the corresponding transfer operators entry is 1.
+ *
+ * On each level a fixed number of v-cycles is performed.
+ *
+ * This method can be called before the preprocess() method.
+ * It does not change the ObstacleTNNMGStep state despite
+ * updating the solution vector.
+ *
+ * \param coarseIterationSteps Number of v-cycle performed on the corse levels.
+ */
+ void nestedIteration(int coarseIterationSteps=2)
+ {
+ int maxLevel = transfer_.size();
+
+ std::vector<Matrix> coarseMatrix(maxLevel);
+ std::vector<Vector> coarseSolution(maxLevel);
+ std::vector<Vector> coarseRhs(maxLevel);
+ std::vector<BitVector> coarseIgnore(maxLevel);
+ std::vector<BoxConstraintVector> coarseObstacle(maxLevel);
+ MandelObstacleRestrictor<Vector> obstacleRestrictor;
+ BitVector critical;
+
+ critical.resize(x_->size());
+ critical.unsetAll();
+
+ transfer_[maxLevel-1]->galerkinRestrictSetOccupation(mat_, coarseMatrix[maxLevel-1]);
+ transfer_[maxLevel-1]->galerkinRestrict(mat_, coarseMatrix[maxLevel-1]);
+ transfer_[maxLevel-1]->restrict(rhs_, coarseRhs[maxLevel-1]);
+ transfer_[maxLevel-1]->restrictToFathers(*ignoreNodes_, coarseIgnore[maxLevel-1]);
+ obstacleRestrictor.restrict(obstacles_, coarseObstacle[maxLevel-1], hasObstacle_, hasObstacle_, *(transfer_[maxLevel-1]), critical);
+ for (int i = maxLevel-2; i>=0; --i)
+ {
+ transfer_[i]->galerkinRestrictSetOccupation(coarseMatrix[i+1], coarseMatrix[i]);
+ transfer_[i]->galerkinRestrict(coarseMatrix[i+1], coarseMatrix[i]);
+ transfer_[i]->restrict(coarseRhs[i+1], coarseRhs[i]);
+ transfer_[i]->restrictToFathers(coarseIgnore[i+1], coarseIgnore[i]);
+ obstacleRestrictor.restrict(coarseObstacle[i+1], coarseObstacle[i], hasObstacle_, hasObstacle_, *(transfer_[i]), critical);
+ coarseSolution[i].resize(coarseMatrix[i].N());
+ }
+
+ coarseSolution[0] = 0;
+ for (int i = 0; i<maxLevel; ++i)
+ {
+
+ TransferPointerVector coarseTransfer(i);
+ for (int j = 0; j < i; ++j)
+ coarseTransfer[j] = transfer_[j];
+
+ ObstacleTNNMGStep coarseTNNMGStep(
+ coarseMatrix[i],
+ coarseSolution[i],
+ coarseRhs[i],
+ coarseIgnore[i],
+ coarseObstacle[i],
+ coarseTransfer,
+ truncationTol_);
+ coarseTNNMGStep.preprocess();
+ for (int j = 0; j < coarseIterationSteps; ++j)
+ coarseTNNMGStep.iterate();
+ coarseSolution[i] = coarseTNNMGStep.getSol();
+
+ if (i<maxLevel-1)
+ transfer_[i]->prolong(coarseSolution[i], coarseSolution[i+1]);
+ }
+ transfer_[maxLevel-1]->prolong(coarseSolution[maxLevel-1], *x_);
+ }
+
protected:
// problem description given by the user
--
GitLab