Remove BlockGaussSeidelStep

src/one-body-sample.cc

@@ -34,7 +34,6 @@
 #include <dune/fufem/functionspacebases/p0basis.hh>
 #include <dune/fufem/functionspacebases/p1nodalbasis.hh>
 #include <dune/solvers/common/numproc.hh> // Solver::FULL
-#include <dune/solvers/iterationsteps/blockgsstep.hh>
 #include <dune/solvers/iterationsteps/truncatedblockgsstep.hh>
 #include <dune/solvers/norms/energynorm.hh>
 #include <dune/solvers/solvers/loopsolver.hh>
@@ -188,12 +187,10 @@ int main() {
 
     VectorType u1(grid.size(grid.maxLevel(), dim));
     u1 = 0;
-    VectorType u2 = u1;
     VectorType u3 = u1;
     CellVectorType vonMisesStress;
 
     VectorType b1;
-    VectorType b2;
     VectorType b3;
 
     std::vector<Dune::FieldVector<double, 1>> nodalIntegrals;
@@ -207,11 +204,10 @@ int main() {
       // b = neumann
       assemble_neumann<GridType, GridType::LeafGridView, SmallVector, P1Basis>(
           leafView, p1Basis, neumannNodes, b1, run);
-      b3 = b2 = b1;
+      b3 = b1;
 
       // b += linear update
       stiffnessMatrix.umv(u1, b1);
-      stiffnessMatrix.umv(u2, b2);
       stiffnessMatrix.umv(u3, b3);
 
       // {{{ Assemble terms for the nonlinearity
@@ -265,18 +261,6 @@ int main() {
 
       // Use a linear solver for comparison; should return roughly the
       // same results if phi vanishes (e.g. because the normalstress is zero)
-      {
-        // TODO: Why does blockGSStep even provide a default constructor?
-        BlockGSStep<OperatorType, VectorType> blockGSStep(stiffnessMatrix, u2,
-                                                          b2);
-        blockGSStep.ignoreNodes_ = &ignoreNodes;
-
-        LoopSolver<VectorType> solver(&blockGSStep, solver_maxIterations,
-                                      solver_tolerance, &energyNorm,
-                                      Solver::QUIET); // Solver::QUIET);
-        solver.solve();
-      }
-
       {
         TruncatedBlockGSStep<OperatorType, VectorType> blockGSStep(
             stiffnessMatrix, u3, b3);
@@ -290,11 +274,6 @@ int main() {
     }
     std::cout << std::endl;
 
-    VectorType diff2 = u2;
-    diff2 -= u1;
-    std::cout << "sup |u1 - u2| = " << diff2.infinity_norm() << std::endl;
-    std::cout << "|u1 - u2| = " << diff2.two_norm() << std::endl;
-
     VectorType diff3 = u3;
     diff3 -= u1;
     std::cout << "sup |u1 - u3| = " << diff3.infinity_norm() << std::endl;
@@ -303,10 +282,9 @@ int main() {
     // Print displacement on frictional boundary
     for (size_t i = 0; i < frictionalNodes.size(); ++i)
       if (frictionalNodes[i][0])
-        std::cout << boost::format("u1[%02d] = %+3e, %|40t|u2[%02d] = %+3e, "
-                                   "%|80t|u3[%02d] = %+3e, %|120t|s[%02d] = "
-                                   "%+3e") %
-                         i % u1[i] % i % u2[i] % i % u3[i] % i %
+        std::cout << boost::format("u1[%02d] = %+3e, %|40t|u3[%02d] = %+3e "
+                                   "%|80t|s[%02d] = %+3e") %
+                         i % u1[i] % i % u3[i] % i %
                          vonMisesStress[i] << std::endl;
   }
   catch (Dune::Exception &e) {