diff --git a/src/riemanniantrsolver.cc b/src/riemanniantrsolver.cc
index 67ac4d3b207de2a34ccaec3b5f520ec209c26ab5..fd442219342887b9e431d9ed0801fcde7545664c 100644
--- a/src/riemanniantrsolver.cc
+++ b/src/riemanniantrsolver.cc
@@ -103,16 +103,16 @@ setup(const GridType& grid,
 
     h1SemiNorm_ = new H1SemiNorm<CorrectionType>(*A);
 
-    innerSolver_ = new ::LoopSolver<CorrectionType>(mmgStep,
-                                                    innerIterations_,
-                                                    innerTolerance_,
-                                                    h1SemiNorm_,
-                                                    Solver::QUIET);
+    innerSolver_ = std::shared_ptr<LoopSolver<CorrectionType> >(new ::LoopSolver<CorrectionType>(mmgStep,
+                                                                                                   innerIterations_,
+                                                                                                   innerTolerance_,
+                                                                                                   h1SemiNorm_,
+                                                                                                 Solver::QUIET));
 
     // Write all intermediate solutions, if requested
     if (instrumented_
-        && dynamic_cast<IterativeSolver<CorrectionType>*>(innerSolver_))
-        dynamic_cast<IterativeSolver<CorrectionType>*>(innerSolver_)->historyBuffer_ = "tmp/mgHistory";
+        && dynamic_cast<IterativeSolver<CorrectionType>*>(innerSolver_.get()))
+        dynamic_cast<IterativeSolver<CorrectionType>*>(innerSolver_.get())->historyBuffer_ = "tmp/mgHistory";
 
     // ////////////////////////////////////////////////////////////
     //    Create Hessian matrix and its occupation structure
@@ -234,8 +234,8 @@ void RiemannianTrustRegionSolver<GridType,TargetSpace>::solve()
     MonotoneMGStep<MatrixType,CorrectionType>* mgStep = NULL;
 
     // if the inner solver is a monotone multigrid set up a max-norm trust-region
-    if (dynamic_cast<LoopSolver<CorrectionType>*>(innerSolver_)) {
-        mgStep = dynamic_cast<MonotoneMGStep<MatrixType,CorrectionType>*>(dynamic_cast<LoopSolver<CorrectionType>*>(innerSolver_)->iterationStep_);
+    if (dynamic_cast<LoopSolver<CorrectionType>*>(innerSolver_.get())) {
+        mgStep = dynamic_cast<MonotoneMGStep<MatrixType,CorrectionType>*>(dynamic_cast<LoopSolver<CorrectionType>*>(innerSolver_.get())->iterationStep_);
     
     }    
 
@@ -330,8 +330,8 @@ void RiemannianTrustRegionSolver<GridType,TargetSpace>::solve()
         
         } else {       // inner solver is a truncated cg
 
-            assert((dynamic_cast<TruncatedCGSolver<MatrixType,CorrectionType>*>(innerSolver_)));
-            dynamic_cast<TruncatedCGSolver<MatrixType,CorrectionType>*>(innerSolver_)->setProblem(*hessianMatrix_, 
+            assert((dynamic_cast<TruncatedCGSolver<MatrixType,CorrectionType>*>(innerSolver_.get())));
+            dynamic_cast<TruncatedCGSolver<MatrixType,CorrectionType>*>(innerSolver_.get())->setProblem(*hessianMatrix_, 
                                                                                                   &corr, 
                                                                                                   &backupRhs, 
                                                                                                   trustRegion.radius());
diff --git a/src/riemanniantrsolver.hh b/src/riemanniantrsolver.hh
index 7919ef3a2ca0206264598675619e25b0002e8ce4..75e40dec583ddea6d6812437cd4af8625b961b30 100644
--- a/src/riemanniantrsolver.hh
+++ b/src/riemanniantrsolver.hh
@@ -103,7 +103,7 @@ protected:
     const GeodesicFEAssembler<typename GridType::LeafGridView, TargetSpace>* assembler_;
 
     /** \brief The solver for the quadratic inner problems */
-    Solver* innerSolver_;
+    std::shared_ptr<Solver> innerSolver_;
 
     /** \brief Dummy fields containing 'true' everywhere.  The multigrid step
         expects them :-( */