[!252] Generalize countNonZeros to scalar matrices

Merge branch 'generalize-nonzerocounter-to-scalar-matrices' into 'master'

See merge request [!252]

  [!252]: Nonecore/dune-istl/merge_requests/252

dune/istl/matrixutils.hh

@@ -37,43 +37,6 @@ namespace Dune
    * @brief Some handy generic functions for ISTL matrices.
    * @author Markus Blatt
    */
-  namespace
-  {
-
-    template<int i>
-    struct NonZeroCounter
-    {
-      template<class M>
-      static typename M::size_type count(const M& matrix)
-      {
-        typedef typename M::ConstRowIterator RowIterator;
-
-        RowIterator endRow = matrix.end();
-        typename M::size_type nonZeros = 0;
-
-        for(RowIterator row = matrix.begin(); row != endRow; ++row) {
-          typedef typename M::ConstColIterator Entry;
-          Entry endEntry = row->end();
-          for(Entry entry = row->begin(); entry != endEntry; ++entry) {
-            nonZeros += NonZeroCounter<i-1>::count(*entry);
-          }
-        }
-        return nonZeros;
-      }
-    };
-
-    template<>
-    struct NonZeroCounter<1>
-    {
-      template<class M>
-      static typename M::size_type count(const M& matrix)
-      {
-        return matrix.N()*matrix.M();
-      }
-    };
-
-  }
-
   /**
    * @brief Check whether the a matrix has diagonal values
    * on blocklevel recursion levels.
@@ -150,10 +113,23 @@ namespace Dune
    * number of nonzero blocks time n*m.
    */
   template<class M>
-  inline int countNonZeros(const M& matrix)
+  inline auto countNonZeros(const M& matrix,typename std::enable_if_t<Dune::IsNumber<M>::value>* sfinae = nullptr)
   {
-    return NonZeroCounter<M::blocklevel>::count(matrix);
+    return 1;
   }
+
+  template<class M>
+  inline auto countNonZeros(const M& matrix,typename std::enable_if_t<!Dune::IsNumber<M>::value>* sfinae = nullptr)
+  {
+    typename M::size_type nonZeros = 0;
+
+    for (auto row = matrix.begin(); row != matrix.end(); ++row)
+      for (auto entry = row->begin(); entry != row->end(); ++entry)
+        nonZeros += countNonZeros(*entry);
+
+    return nonZeros;
+  }
+
   /*
      template<class M>
      struct ProcessOnFieldsOfMatrix

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;