diff --git a/cosserat-continuum.cc b/cosserat-continuum.cc
index 6ceed2e1d4191de150d94b7dd3cf6890ec772927..f02c99d5bac9fd3eab66505ef0f7ff8d54a64deb 100644
--- a/cosserat-continuum.cc
+++ b/cosserat-continuum.cc
@@ -228,10 +228,15 @@ int main (int argc, char *argv[]) try
     //   Create an assembler for the energy functional
     // ////////////////////////////////////////////////////////////
 
+    typedef P1NodalBasis<typename GridType::LeafGridView,double> P1Basis;
+    P1Basis p1Basis(grid.leafView());
+    
     const ParameterTree& materialParameters = parameterSet.sub("materialParameters");
-    CosseratEnergyLocalStiffness<GridType::LeafGridView,3> cosseratEnergyLocalStiffness(materialParameters);
+    CosseratEnergyLocalStiffness<GridType::LeafGridView,
+                                 typename P1Basis::LocalFiniteElement,
+                                 3> cosseratEnergyLocalStiffness(materialParameters);
 
-    GeodesicFEAssembler<GridType::LeafGridView,TargetSpace> assembler(grid.leafView(),
+    GeodesicFEAssembler<P1Basis,TargetSpace> assembler(grid.leafView(),
                                                                       &cosseratEnergyLocalStiffness);
 
     // /////////////////////////////////////////////////
diff --git a/dune/gfe/cosseratenergystiffness.hh b/dune/gfe/cosseratenergystiffness.hh
index 707297f55ff91807b77623f8fa19e9f54ebf78a1..6e2fb19f32550a872a1c361c3ddba7fadfc063d3 100644
--- a/dune/gfe/cosseratenergystiffness.hh
+++ b/dune/gfe/cosseratenergystiffness.hh
@@ -10,9 +10,9 @@
 #include <dune/gfe/rigidbodymotion.hh>
 
 
-template<class GridView, int dim>
+template<class GridView, class LocalFiniteElement, int dim>
 class CosseratEnergyLocalStiffness 
-    : public LocalGeodesicFEStiffness<GridView,RigidBodyMotion<dim> >
+    : public LocalGeodesicFEStiffness<GridView,LocalFiniteElement,RigidBodyMotion<dim> >
 {
     // grid types
     typedef typename GridView::Grid::ctype DT;
@@ -120,7 +120,8 @@ public:
 
     /** \brief Assemble the energy for a single element */
     RT energy (const Entity& e,
-               const Dune::array<TargetSpace, gridDim+1>& localSolution) const;
+               const LocalFiniteElement& localFiniteElement,
+               const std::vector<TargetSpace>& localSolution) const;
                
     RT quadraticMembraneEnergy(const Dune::FieldMatrix<double,3,3>& U) const
     {
@@ -182,16 +183,17 @@ public:
     
 };
 
-template <class GridView, int dim>
-typename CosseratEnergyLocalStiffness<GridView, dim>::RT CosseratEnergyLocalStiffness<GridView, dim>::
+template <class GridView, class LocalFiniteElement, int dim>
+typename CosseratEnergyLocalStiffness<GridView,LocalFiniteElement,dim>::RT
+CosseratEnergyLocalStiffness<GridView,LocalFiniteElement,dim>::
 energy(const Entity& element,
-       const Dune::array<RigidBodyMotion<dim>, gridDim+1>& localSolution) const
+       const LocalFiniteElement& localFiniteElement,
+       const std::vector<RigidBodyMotion<dim> >& localSolution) const
 {
     RT energy = 0;
 
-    assert(element.type().isSimplex());
-    
-    LocalGeodesicFEFunction<gridDim, double, TargetSpace> localGeodesicFEFunction(localSolution);
+    LocalGeodesicFEFunction<gridDim, double, LocalFiniteElement, TargetSpace> localGeodesicFEFunction(localFiniteElement,
+                                                                                                      localSolution);
 
     int quadOrder = 1;//gridDim;
 
diff --git a/dune/gfe/localgeodesicfefunction.hh b/dune/gfe/localgeodesicfefunction.hh
index 5f4f729996b10db302598437fe5fc7cad09fbe39..9efccb2583cbe8e87278f07b89691bddff35b71c 100644
--- a/dune/gfe/localgeodesicfefunction.hh
+++ b/dune/gfe/localgeodesicfefunction.hh
@@ -31,16 +31,6 @@ class LocalGeodesicFEFunction
 
 public:
 
-    /** \brief Constructor */
-    LocalGeodesicFEFunction(const LocalFiniteElement& localFiniteElement,
-                            const Dune::array<TargetSpace,dim+1>& coefficients) DUNE_DEPRECATED
-        : localFiniteElement_(localFiniteElement),
-        coefficients_(coefficients.size())
-    {
-        for (size_t i=0; i<coefficients.size(); i++)
-            coefficients_[i] = coefficients[i];
-    }
-
     /** \brief Constructor 
      \param type Type of the reference element
      */
@@ -567,16 +557,19 @@ class LocalGeodesicFEFunction<dim,ctype,LocalFiniteElement,RigidBodyMotion<3> >
 public:
 
     /** \brief Constructor */
-    LocalGeodesicFEFunction(const Dune::array<TargetSpace,dim+1>& coefficients)
+    LocalGeodesicFEFunction(const LocalFiniteElement& localFiniteElement,
+                            const std::vector<TargetSpace>& coefficients)
     {
-        for (int i=0; i<dim+1; i++)
+        assert(localFiniteElement.localBasis().size() == coefficients.size());
+        
+        for (int i=0; i<coefficients.size(); i++)
             translationCoefficients_[i] = coefficients[i].r;
         
-        Dune::array<Rotation<3,ctype>,dim+1> orientationCoefficients;
-        for (int i=0; i<dim+1; i++)
+        std::vector<Rotation<3,ctype> > orientationCoefficients(coefficients.size());
+        for (int i=0; i<coefficients.size(); i++)
             orientationCoefficients[i] = coefficients[i].q;
         
-        orientationFEFunction_ = std::auto_ptr<LocalGeodesicFEFunction<dim,ctype,LocalFiniteElement,Rotation<3,double> > > (new LocalGeodesicFEFunction<dim,ctype,LocalFiniteElement,Rotation<3,double> >(orientationCoefficients));
+        orientationFEFunction_ = std::auto_ptr<LocalGeodesicFEFunction<dim,ctype,LocalFiniteElement,Rotation<3,double> > > (new LocalGeodesicFEFunction<dim,ctype,LocalFiniteElement,Rotation<3,double> >(localFiniteElement,orientationCoefficients));
         
     }
 
diff --git a/dune/gfe/riemanniantrsolver.hh b/dune/gfe/riemanniantrsolver.hh
index 9912e07dd2d7403aaf84ea3f25e14cf79780860b..2a1d4642262c4419ea5e4703f2c540e168c813d1 100644
--- a/dune/gfe/riemanniantrsolver.hh
+++ b/dune/gfe/riemanniantrsolver.hh
@@ -13,6 +13,8 @@
 #include <dune/solvers/solvers/iterativesolver.hh>
 #include <dune/solvers/solvers/loopsolver.hh>
 
+#include <dune/fufem/functionspacebases/p1nodalbasis.hh>
+
 #include "geodesicfeassembler.hh"
 
 /** \brief Riemannian trust-region solver for geodesic finite-element problems */