Tests: Cover the case of partial Dirichlet BC

dune/solvers/test/common.hh

@@ -357,7 +357,7 @@ void markBoundaryDOFs(const GridView& gridView, BitVector& isBoundary)
     }
 }
 
-template <size_t blocksize, class GridView>
+template <size_t blocksize, class GridView, bool trivialDirichletOnly = true>
 class SymmetricSampleProblem {
 public:
   typedef typename Dune::FieldMatrix<double, blocksize, blocksize> MatrixBlock;
@@ -376,7 +376,15 @@ public:
 
     ignore.resize(A.N());
     ignore.unsetAll();
-    markBoundaryDOFs(gridView, ignore);
+    if (trivialDirichletOnly)
+        markBoundaryDOFs(gridView, ignore);
+    else {
+        // Mark the first component only
+        Dune::BitSetVector<blocksize> boundaryNodes(A.N(), false);
+        markBoundaryDOFs(gridView, boundaryNodes);
+        for (size_t i=0; i < boundaryNodes.size(); ++i)
+            ignore[i][0] = boundaryNodes[i][0];
+    }
 
     u.resize(A.N());
     u_ex.resize(A.N());
@@ -389,9 +397,9 @@ public:
     u = 0;
     for (size_t i = 0; i < u.size(); ++i)
       for (size_t j = 0; j < blocksize; ++j) {
-        // in case of no dirichlet values you can make the
-        // matrix pd hereby
-        // A[i][i][j][j] += 0.5*std::abs(A[0][0][0][0]);
+        // Make the matrix positive-definitive
+        if (not trivialDirichletOnly)
+          A[i][i][j][j] += 0.5*std::abs(A[0][0][0][0]);
         u_ex[i][j] = (1.0 * rand()) / RAND_MAX;
       }
 

dune/solvers/test/multigridtest.cc

@@ -32,7 +32,7 @@
   *  It furthermore tests the functionality to solve correctly for a different rhs with otherwise unaltered data without calling preprocess().
   *  Setting and using different smoothers in different levels is NOT a tested feature.
   */
-template <size_t blocksize>
+template <size_t blocksize, bool trivialDirichletOnly = true>
 struct MultigridTestSuite
 {
     template <class GridType>
@@ -44,7 +44,7 @@ struct MultigridTestSuite
         bool passed = true;
 
         using Problem =
-            SymmetricSampleProblem<blocksize, typename GridType::LevelGridView>;
+            SymmetricSampleProblem<blocksize, typename GridType::LevelGridView, trivialDirichletOnly>;
         Problem p(grid.levelGridView(grid.maxLevel()));
 
         typedef typename Problem::Vector Vector;
@@ -127,9 +127,11 @@ int main(int argc, char** argv)
 
     MultigridTestSuite<1> testsuite1;
     MultigridTestSuite<2> testsuite2;
+    MultigridTestSuite<3, false> testsuite3f;
     MultigridTestSuite<8> testsuite8;
     passed = checkWithStandardGrids(testsuite1);
     passed = checkWithStandardGrids(testsuite2);
+    passed = checkWithStandardGrids(testsuite3f);
     passed = checkWithStandardGrids(testsuite8);
 
     return passed ? 0 : 1;