diff --git a/dune/istl/io.hh b/dune/istl/io.hh
index 966cea4adf55496cb46ac395e2ae856b1fdcd667..b01c049e97c184d00a037760ef429e9dec5a5247 100644
--- a/dune/istl/io.hh
+++ b/dune/istl/io.hh
@@ -246,6 +246,36 @@ namespace Dune {
     s.precision(oldprec);
   }
 
+  namespace Impl
+  {
+    template<class B>
+    void printInnerMatrixElement(std::ostream& s,
+                                 const B& innerMatrixElement,
+                                 int innerrow, int innercol)
+    {
+      s<<innerMatrixElement<<" ";
+    }
+
+    template<class B, int n>
+    void printInnerMatrixElement(std::ostream& s,
+                                 const ScaledIdentityMatrix<B,n> innerMatrixElement,
+                                 int innerrow, int innercol)
+    {
+      if (innerrow == innercol)
+        s<<innerMatrixElement.scalar()<<" ";
+      else
+        s<<"-";
+    }
+
+    template<class B, int n, int m>
+    void printInnerMatrixElement(std::ostream& s,
+                                 const FieldMatrix<B,n,m> innerMatrixElement,
+                                 int innerrow, int innercol)
+    {
+      s<<innerMatrixElement[innerrow][innercol]<<" ";
+    }
+  }
+
   /**
    * \brief Prints a BCRSMatrix with fixed sized blocks.
    *
@@ -290,8 +320,8 @@ namespace Dune {
 
     typedef typename Matrix::ConstRowIterator Row;
 
-    int n = InnerMatrixType::rows;
-    int m = InnerMatrixType::cols;
+    constexpr int n = std::decay_t<decltype(Impl::asMatrix(std::declval<InnerMatrixType>()))>::rows;
+    constexpr int m = std::decay_t<decltype(Impl::asMatrix(std::declval<InnerMatrixType>()))>::cols;
     for(Row row=mat.begin(); row != mat.end(); ++row) {
       int skipcols=0;
       bool reachedEnd=false;
@@ -325,7 +355,7 @@ namespace Dune {
             }
             for(int innercol=0; innercol < m; ++innercol) {
               s.width(9);
-              printInnerMatrixElement(s,*col,innerrow,innercol);
+              Impl::printInnerMatrixElement(s,*col,innerrow,innercol);
             }
 
             s<<"|";
@@ -346,25 +376,6 @@ namespace Dune {
     s.precision(oldprec);
   }
 
-  template<class B, int n>
-  void printInnerMatrixElement(std::ostream& s,
-                         const ScaledIdentityMatrix<B,n> innerMatrixElement,
-                         int innerrow, int innercol)
-  {
-    if (innerrow == innercol)
-      s<<innerMatrixElement.scalar()<<" ";
-    else
-      s<<"-";
-  }
-
-  template<class B, int n, int m, class A>
-  void printInnerMatrixElement(std::ostream& s,
-                         const FieldMatrix<B,n,m> innerMatrixElement,
-                         int innerrow, int innercol)
-  {
-    s<<innerMatrixElement[innerrow][innercol]<<" ";
-  }
-
   namespace
   {
     template<typename T>
diff --git a/dune/istl/test/iotest.cc b/dune/istl/test/iotest.cc
index 8704bb85a5f0bbcd57f5c33d0aa5e48a3f881019..a1a12fcf8ad4d19503253fd9819a4cf4b6186fa5 100644
--- a/dune/istl/test/iotest.cc
+++ b/dune/istl/test/iotest.cc
@@ -89,22 +89,30 @@ int main(int argc, char** argv)
   v3.push_back(2);
   Dune::writeVectorToMatlabHelper(v3, std::cout);
 
-  // Test the printmatrix method
+  // Test the printmatrix and printSparseMatrix methods
   // BCRSMatrix
   {
     Dune::BCRSMatrix<double> matrix;
     setupLaplacian(matrix, 3);
     Dune::printmatrix(std::cout, matrix, "BCRSMatrix<double>", "--");
+    Dune::printSparseMatrix(std::cout, matrix, "BCRSMatrix<double>", "--");
   }
   {
     Dune::BCRSMatrix<Dune::FieldMatrix<double,1,1> > matrix;
     setupLaplacian(matrix, 3);
     Dune::printmatrix(std::cout, matrix, "BCRSMatrix<FieldMatrix<double,1,1> >", "--");
+    Dune::printSparseMatrix(std::cout, matrix, "BCRSMatrix<FieldMatrix<double,1,1> >", "--");
+  }
+  {
+    Dune::BCRSMatrix<Dune::ScaledIdentityMatrix<double,1> > matrix;
+    setupLaplacian(matrix, 3);
+    Dune::printSparseMatrix(std::cout, matrix, "BCRSMatrix<ScaledIdentityMatrix<double,1> >", "--");
   }
   {
     Dune::BCRSMatrix<Dune::FieldMatrix<double,2,3> > matrix;
     setupLaplacian(matrix, 3);
     Dune::printmatrix(std::cout, matrix, "BCRSMatrix<FieldMatrix<double,2,3> >", "--");
+    Dune::printSparseMatrix(std::cout, matrix, "BCRSMatrix<FieldMatrix<double,2,3> >", "--");
   }
 
   // Matrix