diff --git a/src/riemanniantrsolver.cc b/src/riemanniantrsolver.cc
index 9c362fe1da6620ce6bc7a0c14c492713fa91e328..124d095dc3a5a5776af331f0c68d0ed41a2cde89 100644
--- a/src/riemanniantrsolver.cc
+++ b/src/riemanniantrsolver.cc
@@ -100,15 +100,16 @@ setup(const GridType& grid,
 
     h1SemiNorm_ = new H1SemiNorm<CorrectionType>(**A);
 
-    mmgSolver_ = new ::LoopSolver<CorrectionType>(mmgStep,
+    innerSolver_ = new ::LoopSolver<CorrectionType>(mmgStep,
                                                      multigridIterations_,
                                                      qpTolerance_,
                                                      h1SemiNorm_,
                                                      Solver::QUIET);
 
     // Write all intermediate solutions, if requested
-    if (instrumented_)
-        mmgSolver_->historyBuffer_ = "tmp/mgHistory";
+    if (instrumented_
+        && dynamic_cast<IterativeSolver<CorrectionType>*>(innerSolver_))
+        dynamic_cast<IterativeSolver<CorrectionType>*>(innerSolver_)->historyBuffer_ = "tmp/mgHistory";
 
     // ////////////////////////////////////////////////////////////
     //    Create Hessian matrix and its occupation structure
@@ -149,11 +150,19 @@ setup(const GridType& grid,
 template <class GridType, class TargetSpace>
 void RiemannianTrustRegionSolver<GridType,TargetSpace>::solve()
 {
-    MaxNormTrustRegion<blocksize> trustRegion(x_.size(), initialTrustRegionRadius_);
+    MonotoneMGStep<MatrixType,CorrectionType>* mgStep = NULL;
 
-    std::vector<std::vector<BoxConstraint<field_type,blocksize> > > trustRegionObstacles(dynamic_cast<MultigridStep<MatrixType,CorrectionType>*>(mmgSolver_->iterationStep_)->numLevels_);
+    // 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_);
     
-    // /////////////////////////////////////////////////////
+    }    
+
+    MaxNormTrustRegion<blocksize> trustRegion(x_.size(), initialTrustRegionRadius_);
+
+    std::vector<std::vector<BoxConstraint<field_type,blocksize> > > trustRegionObstacles(mgStep->numLevels_);
+
+   // /////////////////////////////////////////////////////
     //   Set up the log file, if requested
     // /////////////////////////////////////////////////////
     FILE* fp;
@@ -193,23 +202,23 @@ void RiemannianTrustRegionSolver<GridType,TargetSpace>::solve()
 
         rhs *= -1;
 
-        dynamic_cast<MultigridStep<MatrixType,CorrectionType>*>(mmgSolver_->iterationStep_)->setProblem(*hessianMatrix_, corr, rhs, grid_->maxLevel()+1);
+        mgStep->setProblem(*hessianMatrix_, corr, rhs, grid_->maxLevel()+1);
 
         trustRegionObstacles.back() = trustRegion.obstacles();
-        dynamic_cast<MonotoneMGStep<MatrixType, CorrectionType>*>(mmgSolver_->iterationStep_)->obstacles_ = &trustRegionObstacles;
+        mgStep->obstacles_ = &trustRegionObstacles;
 
         
-        mmgSolver_->preprocess();
+        innerSolver_->preprocess();
         
-        dynamic_cast<MultigridStep<MatrixType,CorrectionType>*>(mmgSolver_->iterationStep_)->preprocess();
+        mgStep->preprocess();
         
         
         // /////////////////////////////
         //    Solve !
         // /////////////////////////////
-        mmgSolver_->solve();
+        innerSolver_->solve();
         
-        corr = dynamic_cast<MultigridStep<MatrixType,CorrectionType>*>(mmgSolver_->iterationStep_)->getSol();
+        corr = mgStep->getSol();
         
         //std::cout << "Correction: " << std::endl << corr << std::endl;
         
diff --git a/src/riemanniantrsolver.hh b/src/riemanniantrsolver.hh
index b10e1dfef1e7e073f938f52b49aa25970e10cd1b..6ee746168be73974be08931d54f6b8daae9c8fc6 100644
--- a/src/riemanniantrsolver.hh
+++ b/src/riemanniantrsolver.hh
@@ -10,7 +10,7 @@
 
 #include <dune-solvers/boxconstraint.hh>
 #include <dune-solvers/norms/h1seminorm.hh>
-#include <dune-solvers/solvers/loopsolver.hh>
+#include <dune-solvers/solvers/solver.hh>
 
 #include "geodesicfeassembler.hh"
 
@@ -105,8 +105,8 @@ protected:
     /** \brief The assembler for the material law */
     const GeodesicFEAssembler<typename GridType::LeafGridView, TargetSpace>* assembler_;
 
-    /** \brief The multigrid solver */
-    LoopSolver<CorrectionType>* mmgSolver_;
+    /** \brief The solver for the quadratic inner problems */
+    Solver* innerSolver_;
 
     /** \brief Dummy fields containing 'true' everywhere.  The multigrid step
         expects them :-( */