Implement and test method countNonZeros for BCRSMatrix<double>

dune/istl/matrixutils.hh

@@ -63,12 +63,12 @@ namespace Dune
     };
 
     template<>
-    struct NonZeroCounter<1>
+    struct NonZeroCounter<0>
     {
       template<class M>
-      static typename M::size_type count(const M& matrix)
+      static auto count(const M& matrix)
       {
-        return matrix.N()*matrix.M();
+        return 1;
       }
     };
 

dune/istl/test/laplacian.hh

@@ -43,15 +43,20 @@ void setupLaplacian(Dune::BCRSMatrix<B,Alloc>& A, int N)
 
   setupSparsityPattern(A,N);
 
-  B diagonal(static_cast<FieldType>(0)), bone(static_cast<FieldType>(0)),
-  beps(static_cast<FieldType>(0));
-  for(typename B::RowIterator b = diagonal.begin(); b !=  diagonal.end(); ++b)
-    b->operator[](b.index())=4;
-
-
-  for(typename B::RowIterator b=bone.begin(); b !=  bone.end(); ++b)
-    b->operator[](b.index())=-1.0;
-
+  B diagonal(static_cast<FieldType>(0)), bone(static_cast<FieldType>(0));
+
+  Dune::Hybrid::ifElse(Dune::IsNumber<B>(),
+    [&](auto id) {
+      diagonal = B(4.0);
+      bone = B(-1.0);
+    },
+    [&](auto id) {
+      for (auto b = id(diagonal).begin(); b != id(diagonal).end(); ++b)
+        b->operator[](b.index())=4;
+
+      for (auto b=id(bone).begin(); b != id(bone).end(); ++b)
+        b->operator[](b.index())=-1.0;
+    });
 
   for (typename Dune::BCRSMatrix<B,Alloc>::RowIterator i = A.begin(); i != A.end(); ++i) {
     int x = i.index()%N; // x coordinate in the 2d field

dune/istl/test/matrixutilstest.cc

@@ -20,6 +20,15 @@ int main(int argc, char** argv)
 
   const int N=4;
 
+  // Test sparse matrix with scalar entries
+  Dune::BCRSMatrix<double> slaplace;
+  setupLaplacian(slaplace, N);
+
+  if(N*N*5-4*2-(N-2)*4!=countNonZeros(slaplace)) {
+    ++ret;
+    Dune::derr<<"Counting nonzeros of BCRSMatrix<double> failed!"<<std::endl;
+  }
+
   typedef Dune::BCRSMatrix<Dune::FieldMatrix<double,1,1> > BMatrix;
 
   BMatrix laplace;