diff --git a/dune/gfe/riemanniantrsolver.cc b/dune/gfe/riemanniantrsolver.cc
index ffd77f1e6cf38589a12951c689ecd875288ee035..4315233003d2fc346e5b31dca508748542d49478 100644
--- a/dune/gfe/riemanniantrsolver.cc
+++ b/dune/gfe/riemanniantrsolver.cc
@@ -26,10 +26,10 @@
 #include <dune/gfe/parallel/matrixcommunicator.hh>
 #include <dune/gfe/parallel/vectorcommunicator.hh>
 
-template <class GridType, class TargetSpace>
-void RiemannianTrustRegionSolver<GridType,TargetSpace>::
+template <class Basis, class TargetSpace>
+void RiemannianTrustRegionSolver<Basis,TargetSpace>::
 setup(const GridType& grid,
-      const GeodesicFEAssembler<BasisType, TargetSpace>* assembler,
+      const GeodesicFEAssembler<Basis, TargetSpace>* assembler,
          const SolutionType& x,
          const Dune::BitSetVector<blocksize>& dirichletNodes,
          double tolerance,
@@ -112,10 +112,10 @@ setup(const GridType& grid,
     //   Assemble a Laplace matrix to create a norm that's equivalent to the H1-norm
     // //////////////////////////////////////////////////////////////////////////////////////
 
-    BasisType basis(grid.leafGridView());
-    OperatorAssembler<BasisType,BasisType> operatorAssembler(basis, basis);
+    Basis basis(grid.leafGridView());
+    OperatorAssembler<Basis,Basis> operatorAssembler(basis, basis);
 
-    LaplaceAssembler<GridType, typename BasisType::LocalFiniteElement, typename BasisType::LocalFiniteElement> laplaceStiffness;
+    LaplaceAssembler<GridType, typename Basis::LocalFiniteElement, typename Basis::LocalFiniteElement> laplaceStiffness;
     typedef Dune::BCRSMatrix<Dune::FieldMatrix<double,1,1> > ScalarMatrixType;
     ScalarMatrixType localA;
 
@@ -148,7 +148,7 @@ setup(const GridType& grid,
     //   This will be used to monitor the gradient
     // //////////////////////////////////////////////////////////////////////////////////////
 
-    MassAssembler<GridType, typename BasisType::LocalFiniteElement, typename BasisType::LocalFiniteElement> massStiffness;
+    MassAssembler<GridType, typename Basis::LocalFiniteElement, typename Basis::LocalFiniteElement> massStiffness;
     ScalarMatrixType localMassMatrix;
 
     operatorAssembler.assemble(massStiffness, localMassMatrix);
@@ -191,7 +191,7 @@ setup(const GridType& grid,
         TransferOperatorType pkToP1TransferMatrix;
         assembleBasisInterpolationMatrix<TransferOperatorType,
                                          P1NodalBasis<typename GridType::LeafGridView,double>,
-                                         BasisType>(pkToP1TransferMatrix,p1Basis,basis);
+                                         Basis>(pkToP1TransferMatrix,p1Basis,basis);
 
         // If we are on more than 1 processors, join all local transfer matrices on rank 0,
         // and construct a single global transfer operator there.
@@ -261,8 +261,8 @@ setup(const GridType& grid,
 }
 
 
-template <class GridType, class TargetSpace>
-void RiemannianTrustRegionSolver<GridType,TargetSpace>::solve()
+template <class Basis, class TargetSpace>
+void RiemannianTrustRegionSolver<Basis,TargetSpace>::solve()
 {
     int argc = 0;
     char** argv;
diff --git a/dune/gfe/riemanniantrsolver.hh b/dune/gfe/riemanniantrsolver.hh
index 1880f8a4a30fb61517ea8f173d1d8bcf09aa2b0d..118707d6f1ee8d39e5919cf98d7073bac2e10271 100644
--- a/dune/gfe/riemanniantrsolver.hh
+++ b/dune/gfe/riemanniantrsolver.hh
@@ -51,11 +51,13 @@ struct MapperFactory<P3NodalBasis<GridView,field_type> >
 };
 
 /** \brief Riemannian trust-region solver for geodesic finite-element problems */
-template <class GridType, class TargetSpace>
+template <class Basis, class TargetSpace>
 class RiemannianTrustRegionSolver
     : public IterativeSolver<std::vector<TargetSpace>,
                              Dune::BitSetVector<TargetSpace::TangentVector::dimension> >
 {
+    typedef typename Basis::GridView::Grid GridType;
+
     const static int blocksize = TargetSpace::TangentVector::dimension;
 
     const static int gridDim = GridType::dimension;
@@ -68,15 +70,8 @@ class RiemannianTrustRegionSolver
     typedef Dune::BlockVector<Dune::FieldVector<field_type, blocksize> >           CorrectionType;
     typedef std::vector<TargetSpace>                                               SolutionType;
 
-#ifdef THIRD_ORDER
-    typedef P3NodalBasis<typename GridType::LeafGridView,double> BasisType;
-#elif defined SECOND_ORDER
-    typedef P2NodalBasis<typename GridType::LeafGridView,double> BasisType;
-#else
-    typedef P1NodalBasis<typename GridType::LeafGridView,double> BasisType;
-#endif
-    typedef typename MapperFactory<BasisType>::GlobalMapper GlobalMapper;
-    typedef typename MapperFactory<BasisType>::LocalMapper LocalMapper;
+    typedef typename MapperFactory<Basis>::GlobalMapper GlobalMapper;
+    typedef typename MapperFactory<Basis>::LocalMapper LocalMapper;
 
 
 public:
@@ -90,7 +85,7 @@ public:
 
     /** \brief Set up the solver using a monotone multigrid method as the inner solver */
     void setup(const GridType& grid,
-               const GeodesicFEAssembler<BasisType, TargetSpace>* assembler,
+               const GeodesicFEAssembler<Basis, TargetSpace>* assembler,
                const SolutionType& x,
                const Dune::BitSetVector<blocksize>& dirichletNodes,
                double tolerance,
@@ -109,7 +104,7 @@ public:
     {
       scaling_(scaling);
     }
-    
+
     void setIgnoreNodes(const Dune::BitSetVector<blocksize>& ignoreNodes)
     {
         ignoreNodes_ = &ignoreNodes;
@@ -159,7 +154,7 @@ protected:
     std::auto_ptr<MatrixType> hessianMatrix_;
 
     /** \brief The assembler for the material law */
-    const GeodesicFEAssembler<BasisType, TargetSpace>* assembler_;
+    const GeodesicFEAssembler<Basis, TargetSpace>* assembler_;
 
     /** \brief The solver for the quadratic inner problems */
     std::shared_ptr<Solver> innerSolver_;
diff --git a/src/cosserat-continuum.cc b/src/cosserat-continuum.cc
index 1a4bf2d2f98c2a0f86b97e0b1b64e53107ff785d..f58b44964a1c18e659adec1709a0f04c14a2db5f 100644
--- a/src/cosserat-continuum.cc
+++ b/src/cosserat-continuum.cc
@@ -271,7 +271,7 @@ int main (int argc, char *argv[]) try
     //   Create a Riemannian trust-region solver
     // /////////////////////////////////////////////////
 
-    RiemannianTrustRegionSolver<GridType,TargetSpace> solver;
+    RiemannianTrustRegionSolver<FEBasis,TargetSpace> solver;
     solver.setup(*grid,
                  &assembler,
                  x,
diff --git a/src/harmonicmaps.cc b/src/harmonicmaps.cc
index 52bd3b8e5688da0e68424a26a545b889ecf024a1..df0abdfc6d5033e035e2d6806364b13d881c2880 100644
--- a/src/harmonicmaps.cc
+++ b/src/harmonicmaps.cc
@@ -185,7 +185,7 @@ int main (int argc, char *argv[]) try
     //   Create a Riemannian trust-region solver
     // /////////////////////////////////////////////////
 
-    RiemannianTrustRegionSolver<GridType,TargetSpace> solver;
+    RiemannianTrustRegionSolver<FEBasis,TargetSpace> solver;
     solver.setup(*grid,
                  &assembler,
                  x,