diff --git a/src/riemanniantrsolver.cc b/src/riemanniantrsolver.cc
index 124d095dc3a5a5776af331f0c68d0ed41a2cde89..df590b34365a416912cc20994758a99359735176 100644
--- a/src/riemanniantrsolver.cc
+++ b/src/riemanniantrsolver.cc
@@ -115,10 +115,7 @@ setup(const GridType& grid,
     //    Create Hessian matrix and its occupation structure
     // ////////////////////////////////////////////////////////////
     
-    if (hessianMatrix_)
-        delete hessianMatrix_;
-
-    hessianMatrix_ = new MatrixType;
+    hessianMatrix_ = std::auto_ptr<MatrixType>(new MatrixType);
     Dune::MatrixIndexSet indices(grid_->size(1), grid_->size(1));
     assembler_->getNeighborsPerVertex(indices);
     indices.exportIdx(*hessianMatrix_);
diff --git a/src/riemanniantrsolver.hh b/src/riemanniantrsolver.hh
index 6ee746168be73974be08931d54f6b8daae9c8fc6..fe0b26d18d9e288f46fb792ea29c6694b82ecea3 100644
--- a/src/riemanniantrsolver.hh
+++ b/src/riemanniantrsolver.hh
@@ -38,7 +38,7 @@ public:
 
     RiemannianTrustRegionSolver()
         : IterativeSolver<std::vector<TargetSpace>, Dune::BitSetVector<blocksize> >(0,100,NumProc::FULL),
-          hessianMatrix_(NULL), h1SemiNorm_(NULL)
+          hessianMatrix_(std::auto_ptr<MatrixType>(NULL)), h1SemiNorm_(NULL)
     {}
 
     void setup(const GridType& grid, 
@@ -100,7 +100,7 @@ protected:
     double qpTolerance_;
 
     /** \brief Hessian matrix */
-    MatrixType* hessianMatrix_;
+    std::auto_ptr<MatrixType> hessianMatrix_;
 
     /** \brief The assembler for the material law */
     const GeodesicFEAssembler<typename GridType::LeafGridView, TargetSpace>* assembler_;