diff --git a/dune/solvers/common/staticmatrixtools.hh b/dune/solvers/common/staticmatrixtools.hh
index 733407ce05722122cfa3cd8a334873c41dc189d1..1d159a945c3acc74e474ff38ab2e533fcdf6a106 100644
--- a/dune/solvers/common/staticmatrixtools.hh
+++ b/dune/solvers/common/staticmatrixtools.hh
@@ -11,6 +11,7 @@
 
 #include <dune/matrix-vector/addtodiagonal.hh>
 #include <dune/matrix-vector/axpy.hh>
+#include <dune/matrix-vector/ldlt.hh>
 #include <dune/matrix-vector/scalartraits.hh>
 
 #include "arithmetic.hh"
@@ -458,19 +459,7 @@ namespace StaticMatrix
         template<class SymmetricMatrix, class LowerTriangularMatrix, class DiagonalMatrix>
         static void ldlt(const SymmetricMatrix& A, LowerTriangularMatrix& L, DiagonalMatrix& D)
         {
-          for(unsigned int i = 0; i < A.N(); ++i)
-          {
-            D[i][i] = A[i][i];
-            for(unsigned int j = 0; j < i; ++j)
-            {
-              L[i][j] = A[i][j];
-              for(unsigned int k = 0; k < j; ++k)
-                L[i][j] -= L[i][k] * L[j][k] * D[k][k];
-              L[i][j] /= D[j][j];
-            }
-            for(unsigned int k = 0; k < i; ++k)
-              D[i][i] -= L[i][k]*L[i][k] * D[k][k];
-          }
+          Dune::MatrixVector::ldlt(A,L,D);
         }
 
         /** \brief Solve linear system using a LDL^T decomposition.
@@ -492,24 +481,8 @@ namespace StaticMatrix
         template<class LowerTriangularMatrix, class DiagonalMatrix, class RhsVector, class SolVector>
         static void ldltSolve(const LowerTriangularMatrix& L, const DiagonalMatrix& D, const RhsVector& b, SolVector& x)
         {
-          for(unsigned int i = 0; i < x.size(); ++i)
-          {
-            x[i] = b[i];
-            for(unsigned int j = 0; j < i; ++j)
-              x[i] -= L[i][j] * x[j];
-          }
-          for(unsigned int i = 0; i < x.size(); ++i)
-            x[i] /= D[i][i];
-          for(int i = x.size()-1; i >=0; --i)
-          {
-            for(unsigned int j = i+1; j < x.size(); ++j)
-              x[i] -= L[j][i] * x[j];
-          }
+          Dune::MatrixVector::ldltSolve(L,D,b,x);
         }
-
-
-
-
 }
 
 #endif
diff --git a/dune/solvers/iterationsteps/blockgssteps.hh b/dune/solvers/iterationsteps/blockgssteps.hh
index 6835c835c3b1a18aa5bed1c437eb88175e8fd0f0..01d17c64456b2362cf1e14a9b6e470bc3adde5d0 100644
--- a/dune/solvers/iterationsteps/blockgssteps.hh
+++ b/dune/solvers/iterationsteps/blockgssteps.hh
@@ -6,7 +6,8 @@
 
 #include <dune/common/parametertree.hh>
 
-#include <dune/solvers/common/staticmatrixtools.hh>
+#include <dune/matrix-vector/ldlt.hh>
+
 #include <dune/solvers/common/defaultbitvector.hh>
 #include <dune/solvers/iterationsteps/lineariterationstep.hh>
 
@@ -92,8 +93,8 @@ VBlock direct(const MBlock& m, const VBlock& b) {
 template <class MBlock, class VBlock>
 VBlock ldlt(MBlock m, const VBlock& b) {
   VBlock x;
-  StaticMatrix::ldlt(m, m, m);
-  StaticMatrix::ldltSolve(m, m, b, x);
+  Dune::MatrixVector::ldlt(m, m, m);
+  Dune::MatrixVector::ldltSolve(m, m, b, x);
   return x;
 }