diff --git a/dune/gfe/coupling/rodcontinuumcomplex.hh b/dune/gfe/coupling/rodcontinuumcomplex.hh
index 8c500de3a51af3d63f03282f9b570a1f03626661..d70bc3da07bf6f9c5bebb8f0e57aba7046dbaf38 100644
--- a/dune/gfe/coupling/rodcontinuumcomplex.hh
+++ b/dune/gfe/coupling/rodcontinuumcomplex.hh
@@ -27,9 +27,9 @@ class RodContinuumComplex
     /** \brief Holds all data for a rod/continuum coupling */
     struct Coupling
     {
-        LeafBoundaryPatch<RodGrid> rodInterfaceBoundary_;
+        BoundaryPatch<typename RodGrid::LeafGridView> rodInterfaceBoundary_;
 
-        LeafBoundaryPatch<ContinuumGrid> continuumInterfaceBoundary_;
+        BoundaryPatch<typename ContinuumGrid::LeafGridView> continuumInterfaceBoundary_;
         
         /** \brief The orientation of the interface in the reference configuration */
         RigidBodyMotion<dim> referenceInterface_;
@@ -40,7 +40,7 @@ class RodContinuumComplex
     {
         Dune::shared_ptr<RodGrid> grid_;
         
-        LeafBoundaryPatch<RodGrid> dirichletBoundary_;
+        BoundaryPatch<typename RodGrid::LeafGridView> dirichletBoundary_;
         
         RodConfiguration dirichletValues_;
     };
@@ -50,7 +50,7 @@ class RodContinuumComplex
     {
         Dune::shared_ptr<ContinuumGrid> grid_;
         
-        LeafBoundaryPatch<ContinuumGrid> dirichletBoundary_;
+        BoundaryPatch<typename ContinuumGrid::LeafGridView> dirichletBoundary_;
         
         ContinuumConfiguration dirichletValues_;
     };
@@ -125,4 +125,4 @@ public:
     
 };
 
-#endif    // ROD_CONTINUUM_COMPLEX_HH
\ No newline at end of file
+#endif    // ROD_CONTINUUM_COMPLEX_HH
diff --git a/dune/gfe/coupling/rodcontinuumddstep.hh b/dune/gfe/coupling/rodcontinuumddstep.hh
index cc7de07f657f8c853e6c210767ac6cd9061cd339..0579c48f2267f3ffab1d572aa4f759b2e3232b74 100644
--- a/dune/gfe/coupling/rodcontinuumddstep.hh
+++ b/dune/gfe/coupling/rodcontinuumddstep.hh
@@ -177,7 +177,7 @@ mergeRodDirichletAndCouplingBoundaries()
         dirichletAndCouplingNodes.resize(complex_.rodGrid(name)->size(1));
         
         // first copy the true Dirichlet boundary
-        const LeafBoundaryPatch<RodGridType>& dirichletBoundary = complex_.rods_.find(name)->second.dirichletBoundary_;
+        const BoundaryPatch<typename RodGridType::LeafGridView>& dirichletBoundary = complex_.rods_.find(name)->second.dirichletBoundary_;
 
         for (int i=0; i<dirichletAndCouplingNodes.size(); i++)
             dirichletAndCouplingNodes[i] = dirichletBoundary.containsVertex(i);
@@ -189,7 +189,7 @@ mergeRodDirichletAndCouplingBoundaries()
              cIt != continuumNames.end();
              ++cIt) {
 
-            const LeafBoundaryPatch<RodGridType>& rodInterfaceBoundary 
+            const BoundaryPatch<typename RodGridType::LeafGridView>& rodInterfaceBoundary 
                     = complex_.coupling(std::make_pair(name,*cIt)).rodInterfaceBoundary_;
 
             /** \todo Use the BoundaryPatch iterator here, for increased efficiency */
@@ -228,7 +228,7 @@ mergeContinuumDirichletAndCouplingBoundaries()
         dirichletAndCouplingNodes.resize(complex_.continuumGrid(name)->size(dim));
         
         // first copy the true Dirichlet boundary
-        const LeafBoundaryPatch<ContinuumGridType>& dirichletBoundary = complex_.continua_.find(name)->second.dirichletBoundary_;
+        const BoundaryPatch<typename ContinuumGridType::LeafGridView>& dirichletBoundary = complex_.continua_.find(name)->second.dirichletBoundary_;
 
         for (int i=0; i<dirichletAndCouplingNodes.size(); i++)
             dirichletAndCouplingNodes[i] = dirichletBoundary.containsVertex(i);
@@ -240,7 +240,7 @@ mergeContinuumDirichletAndCouplingBoundaries()
              rIt != rodNames.end();
              ++rIt) {
 
-            const LeafBoundaryPatch<ContinuumGridType>& continuumInterfaceBoundary 
+            const BoundaryPatch<typename ContinuumGridType::LeafGridView>& continuumInterfaceBoundary 
                     = complex_.coupling(std::make_pair(*rIt,name)).continuumInterfaceBoundary_;
 
             /** \todo Use the BoundaryPatch iterator here, for increased efficiency */
diff --git a/dune/gfe/coupling/rodcontinuumfixedpointstep.hh b/dune/gfe/coupling/rodcontinuumfixedpointstep.hh
index 1db6365f74f0ccbdab69e868a8ce37af78f6bf30..9c3375f035b10c6a298b99b0b6dfef2762bdcbfe 100644
--- a/dune/gfe/coupling/rodcontinuumfixedpointstep.hh
+++ b/dune/gfe/coupling/rodcontinuumfixedpointstep.hh
@@ -39,6 +39,9 @@ class RodContinuumFixedPointStep
     typedef Dune::BCRSMatrix<Dune::FieldMatrix<double,3,3> > MatrixType;
                 
     typedef typename P1NodalBasis<typename ContinuumGridType::LeafGridView,double>::LocalFiniteElement ContinuumLocalFiniteElement;
+
+    typedef BoundaryPatch<typename RodGridType::LeafGridView> RodLeafBoundaryPatch;
+    typedef BoundaryPatch<typename ContinuumGridType::LeafGridView> ContinuumLeafBoundaryPatch;
     
 public:
     
@@ -117,7 +120,12 @@ protected:
     std::map<std::pair<std::string,std::string>,RigidBodyMotion<3>::TangentVector> 
         rodDirichletToNeumannMap(const std::string& rodName, 
                                  const std::map<std::pair<std::string,std::string>,RigidBodyMotion<3> >& lambda) const;
-    
+
+    std::map<std::pair<std::string,std::string>,RigidBodyMotion<3>::TangentVector> 
+        rodDirichletToNeumannMap(const std::string& rodName, 
+                                 const std::map<std::pair<std::string,std::string>,RigidBodyMotion<3> >& lambda,
+                                 const RodConfigurationType& initialRodX) const; 
+
     std::map<std::pair<std::string,std::string>,RigidBodyMotion<3> >
     continuumNeumannToDirichletMap(const std::string& continuumName,
                                const std::map<std::pair<std::string,std::string>,RigidBodyMotion<3>::TangentVector>& forceTorque,
@@ -215,7 +223,7 @@ mergeRodDirichletAndCouplingBoundaries()
         dirichletAndCouplingNodes.resize(this->complex_.rodGrid(name)->size(1));
         
         // first copy the true Dirichlet boundary
-        const LeafBoundaryPatch<RodGridType>& dirichletBoundary = this->complex_.rods_.find(name)->second.dirichletBoundary_;
+        const RodLeafBoundaryPatch& dirichletBoundary = this->complex_.rods_.find(name)->second.dirichletBoundary_;
 
         for (int i=0; i<dirichletAndCouplingNodes.size(); i++)
             dirichletAndCouplingNodes[i] = dirichletBoundary.containsVertex(i);
@@ -227,7 +235,7 @@ mergeRodDirichletAndCouplingBoundaries()
              cIt != continuumNames.end();
              ++cIt) {
 
-            const LeafBoundaryPatch<RodGridType>& rodInterfaceBoundary 
+            const RodLeafBoundaryPatch& rodInterfaceBoundary 
                     = this->complex_.coupling(std::make_pair(name,*cIt)).rodInterfaceBoundary_;
 
             /** \todo Use the BoundaryPatch iterator here, for increased efficiency */
@@ -266,7 +274,7 @@ mergeContinuumDirichletAndCouplingBoundaries()
         dirichletAndCouplingNodes.resize(this->complex_.continuumGrid(name)->size(dim));
         
         // first copy the true Dirichlet boundary
-        const LeafBoundaryPatch<ContinuumGridType>& dirichletBoundary = this->complex_.continua_.find(name)->second.dirichletBoundary_;
+        const ContinuumLeafBoundaryPatch& dirichletBoundary = this->complex_.continua_.find(name)->second.dirichletBoundary_;
 
         for (int i=0; i<dirichletAndCouplingNodes.size(); i++)
             dirichletAndCouplingNodes[i] = dirichletBoundary.containsVertex(i);
@@ -278,7 +286,7 @@ mergeContinuumDirichletAndCouplingBoundaries()
              rIt != rodNames.end();
              ++rIt) {
 
-            const LeafBoundaryPatch<ContinuumGridType>& continuumInterfaceBoundary 
+            const ContinuumLeafBoundaryPatch& continuumInterfaceBoundary 
                     = this->complex_.coupling(std::make_pair(*rIt,name)).continuumInterfaceBoundary_;
 
             /** \todo Use the BoundaryPatch iterator here, for increased efficiency */
@@ -307,7 +315,7 @@ rodDirichletToNeumannMap(const std::string& rodName,
     ///////////////////////////////////////////////////////////
     //  Set the complete set of Dirichlet values
     ///////////////////////////////////////////////////////////
-    const LeafBoundaryPatch<RodGridType>& dirichletBoundary = this->complex_.rod(rodName).dirichletBoundary_;
+    const RodLeafBoundaryPatch& dirichletBoundary = this->complex_.rod(rodName).dirichletBoundary_;
     const RodConfigurationType& dirichletValues = this->complex_.rod(rodName).dirichletValues_;
     
     for (size_t i=0; i<rodX.size(); i++)
@@ -323,7 +331,7 @@ rodDirichletToNeumannMap(const std::string& rodName,
             continue;
     
         // Use \lambda as a Dirichlet value for the rod
-        const LeafBoundaryPatch<RodGridType>& interfaceBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
+        const RodLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
 
         /** \todo Use the BoundaryPatch iterator, which will be a lot faster
          * once we use EntitySeed for its implementation
@@ -361,7 +369,84 @@ rodDirichletToNeumannMap(const std::string& rodName,
         if (couplingName.first != rodName)
             continue;
         
-        const LeafBoundaryPatch<RodGridType>& couplingBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
+        const RodLeafBoundaryPatch& couplingBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
+
+        result[couplingName] = rod(rodName).assembler_->getResultantForce(couplingBoundary, rodX);
+    }
+    
+    return result;
+}
+
+
+
+
+template <class RodGridType, class ContinuumGridType>
+std::map<std::pair<std::string,std::string>,RigidBodyMotion<3>::TangentVector> RodContinuumFixedPointStep<RodGridType,ContinuumGridType>::
+rodDirichletToNeumannMap(const std::string& rodName, 
+                         const std::map<std::pair<std::string,std::string>,RigidBodyMotion<3> >& lambda,
+                         const RodConfigurationType& initialRodX ) const
+{
+    // container for the subdomain solution
+    RodConfigurationType& rodX = this->rodSubdomainSolutions_[rodName];
+    rodX.resize(this->complex_.rodGrid(rodName)->size(1));
+    
+    rodX = initialRodX;
+    
+    ///////////////////////////////////////////////////////////
+    //  Set the complete set of Dirichlet values
+    ///////////////////////////////////////////////////////////
+    const RodLeafBoundaryPatch& dirichletBoundary = this->complex_.rod(rodName).dirichletBoundary_;
+    const RodConfigurationType& dirichletValues = this->complex_.rod(rodName).dirichletValues_;
+    
+    for (size_t i=0; i<rodX.size(); i++)
+        if (dirichletBoundary.containsVertex(i))
+            rodX[i] = dirichletValues[i];
+    
+    typename std::map<std::pair<std::string,std::string>,RigidBodyMotion<3> >::const_iterator it = lambda.begin();
+    for (; it!=lambda.end(); ++it) {
+        
+        const std::pair<std::string,std::string>& couplingName = it->first;
+        
+        if (couplingName.first != rodName)
+            continue;
+    
+        // Use \lambda as a Dirichlet value for the rod
+        const RodLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
+
+        /** \todo Use the BoundaryPatch iterator, which will be a lot faster
+         * once we use EntitySeed for its implementation
+         */
+        for (size_t i=0; i<rodX.size(); i++)
+            if (interfaceBoundary.containsVertex(i))
+                rodX[i] = it->second;
+    }
+    
+    // Set the correct Dirichlet nodes
+    rod(rodName).solver_->setIgnoreNodes(rod(rodName).dirichletAndCouplingNodes_);
+
+    ////////////////////////////////////////////////////////////////////////////////
+    //  Solve the Dirichlet problem
+    ////////////////////////////////////////////////////////////////////////////////
+           
+
+    rod(rodName).solver_->setInitialSolution(rodX);
+    
+    // Solve the Dirichlet problem
+    rod(rodName).solver_->solve();
+
+    rodX = rod(rodName).solver_->getSol();
+
+    //   Extract Neumann values
+    std::map<std::pair<std::string,std::string>, RigidBodyMotion<3>::TangentVector > result;
+
+    for (it = lambda.begin(); it!=lambda.end(); ++it) {
+        
+        const std::pair<std::string,std::string>& couplingName = it->first;
+        
+        if (couplingName.first != rodName)
+            continue;
+        
+        const RodLeafBoundaryPatch& couplingBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
 
         result[couplingName] = rod(rodName).assembler_->getResultantForce(couplingBoundary, rodX);
     }
@@ -381,7 +466,7 @@ continuumNeumannToDirichletMap(const std::string& continuumName,
     ////////////////////////////////////////////////////
 
     typedef P1NodalBasis<typename ContinuumGridType::LeafGridView,double> P1Basis;
-    const LeafBoundaryPatch<ContinuumGridType>& dirichletBoundary = this->complex_.continuum(continuumName).dirichletBoundary_;
+    const ContinuumLeafBoundaryPatch& dirichletBoundary = this->complex_.continuum(continuumName).dirichletBoundary_;
     P1Basis basis(dirichletBoundary.gridView());
 
     const MatrixType& stiffnessMatrix = *continuum(continuumName).stiffnessMatrix_;
@@ -405,7 +490,7 @@ continuumNeumannToDirichletMap(const std::string& continuumName,
             continue;
         
         // Use 'forceTorque' as a Neumann value for the rod
-        const LeafBoundaryPatch<ContinuumGridType>& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
+        const ContinuumLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
         
         // 
         VectorType localNeumannValues;
@@ -455,7 +540,7 @@ continuumNeumannToDirichletMap(const std::string& continuumName,
             continue;
         
         // Use 'forceTorque' as a Neumann value for the rod
-        const LeafBoundaryPatch<ContinuumGridType>& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
+        const ContinuumLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
         
         computeAverageInterface(interfaceBoundary, x, averageInterface[couplingName]);
         
diff --git a/dune/gfe/coupling/rodcontinuumsteklovpoincarestep.hh b/dune/gfe/coupling/rodcontinuumsteklovpoincarestep.hh
index 627a6dd459ae094743afda0bf07a6b3937e9c5ab..00c2da917a2c9c3c98d7f77468678e96ba195ce0 100644
--- a/dune/gfe/coupling/rodcontinuumsteklovpoincarestep.hh
+++ b/dune/gfe/coupling/rodcontinuumsteklovpoincarestep.hh
@@ -40,7 +40,10 @@ class RodContinuumSteklovPoincareStep
     typedef Dune::BCRSMatrix<Dune::FieldMatrix<double,3,3> > MatrixType;
                 
     typedef typename P1NodalBasis<typename ContinuumGridType::LeafGridView,double>::LocalFiniteElement ContinuumLocalFiniteElement;
-    
+
+    typedef BoundaryPatch<typename RodGridType::LeafGridView> RodLeafBoundaryPatch;
+    typedef BoundaryPatch<typename ContinuumGridType::LeafGridView> ContinuumLeafBoundaryPatch;
+   
 public:
     
     /** \brief Constructor for a complex with one rod and one continuum */
@@ -53,7 +56,7 @@ public:
                                     RiemannianTrustRegionSolver<RodGridType,RigidBodyMotion<3> >* rodSolver,
                                     const MatrixType* stiffnessMatrix3d,
                                     const Dune::shared_ptr< ::LoopSolver<VectorType> > solver,
-                                    LinearLocalAssembler<ContinuumGridType, 
+                                    LocalOperatorAssembler<ContinuumGridType, 
                                                          ContinuumLocalFiniteElement, 
                                                          ContinuumLocalFiniteElement,
                                                          Dune::FieldMatrix<double,dim,dim> >* localAssembler)
@@ -85,7 +88,7 @@ public:
                                     const std::map<std::string,RiemannianTrustRegionSolver<RodGridType,RigidBodyMotion<3> >*>& rodSolver,
                                     const std::map<std::string,const MatrixType*>& stiffnessMatrix3d,
                                     const std::map<std::string, const Dune::shared_ptr< ::LoopSolver<VectorType> > >& solver,
-                                    const std::map<std::string,LinearLocalAssembler<ContinuumGridType, 
+                                    const std::map<std::string,LocalOperatorAssembler<ContinuumGridType, 
                                                                                     ContinuumLocalFiniteElement, 
                                                                                     ContinuumLocalFiniteElement,
                                                                                     Dune::FieldMatrix<double,dim,dim> >*>& localAssembler)
@@ -282,7 +285,7 @@ mergeRodDirichletAndCouplingBoundaries()
         dirichletAndCouplingNodes.resize(this->complex_.rodGrid(name)->size(1));
         
         // first copy the true Dirichlet boundary
-        const LeafBoundaryPatch<RodGridType>& dirichletBoundary = this->complex_.rods_.find(name)->second.dirichletBoundary_;
+        const RodLeafBoundaryPatch& dirichletBoundary = this->complex_.rods_.find(name)->second.dirichletBoundary_;
 
         for (int i=0; i<dirichletAndCouplingNodes.size(); i++)
             dirichletAndCouplingNodes[i] = dirichletBoundary.containsVertex(i);
@@ -294,7 +297,7 @@ mergeRodDirichletAndCouplingBoundaries()
              cIt != continuumNames.end();
              ++cIt) {
 
-            const LeafBoundaryPatch<RodGridType>& rodInterfaceBoundary 
+            const RodLeafBoundaryPatch& rodInterfaceBoundary 
                     = this->complex_.coupling(std::make_pair(name,*cIt)).rodInterfaceBoundary_;
 
             /** \todo Use the BoundaryPatch iterator here, for increased efficiency */
@@ -333,7 +336,7 @@ mergeContinuumDirichletAndCouplingBoundaries()
         dirichletAndCouplingNodes.resize(this->complex_.continuumGrid(name)->size(dim));
         
         // first copy the true Dirichlet boundary
-        const LeafBoundaryPatch<ContinuumGridType>& dirichletBoundary = this->complex_.continua_.find(name)->second.dirichletBoundary_;
+        const ContinuumLeafBoundaryPatch& dirichletBoundary = this->complex_.continua_.find(name)->second.dirichletBoundary_;
 
         for (int i=0; i<dirichletAndCouplingNodes.size(); i++)
             dirichletAndCouplingNodes[i] = dirichletBoundary.containsVertex(i);
@@ -345,7 +348,7 @@ mergeContinuumDirichletAndCouplingBoundaries()
              rIt != rodNames.end();
              ++rIt) {
 
-            const LeafBoundaryPatch<ContinuumGridType>& continuumInterfaceBoundary 
+            const ContinuumLeafBoundaryPatch& continuumInterfaceBoundary 
                     = this->complex_.coupling(std::make_pair(*rIt,name)).continuumInterfaceBoundary_;
 
             /** \todo Use the BoundaryPatch iterator here, for increased efficiency */
@@ -374,7 +377,7 @@ rodDirichletToNeumannMap(const std::string& rodName,
     ///////////////////////////////////////////////////////////
     //  Set the complete set of Dirichlet values
     ///////////////////////////////////////////////////////////
-    const LeafBoundaryPatch<RodGridType>& dirichletBoundary = this->complex_.rod(rodName).dirichletBoundary_;
+    const RodLeafBoundaryPatch& dirichletBoundary = this->complex_.rod(rodName).dirichletBoundary_;
     const RodConfigurationType& dirichletValues = this->complex_.rod(rodName).dirichletValues_;
     
     for (size_t i=0; i<rodX.size(); i++)
@@ -390,7 +393,7 @@ rodDirichletToNeumannMap(const std::string& rodName,
             continue;
     
         // Use \lambda as a Dirichlet value for the rod
-        const LeafBoundaryPatch<RodGridType>& interfaceBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
+        const RodLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
 
         /** \todo Use the BoundaryPatch iterator, which will be a lot faster
          * once we use EntitySeed for its implementation
@@ -428,7 +431,7 @@ rodDirichletToNeumannMap(const std::string& rodName,
         if (couplingName.first != rodName)
             continue;
         
-        const LeafBoundaryPatch<RodGridType>& couplingBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
+        const RodLeafBoundaryPatch& couplingBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
 
         result[couplingName] = rod(rodName).assembler_->getResultantForce(couplingBoundary, rodX);
     }
@@ -448,7 +451,7 @@ continuumDirichletToNeumannMap(const std::string& continuumName,
     x3d = 0;
 
     // Copy the true Dirichlet values into it
-    const LeafBoundaryPatch<ContinuumGridType>& dirichletBoundary = this->complex_.continuum(continuumName).dirichletBoundary_;
+    const ContinuumLeafBoundaryPatch& dirichletBoundary = this->complex_.continuum(continuumName).dirichletBoundary_;
     const VectorType& dirichletValues = this->complex_.continuum(continuumName).dirichletValues_;
     
     for (size_t i=0; i<x3d.size(); i++)
@@ -464,7 +467,7 @@ continuumDirichletToNeumannMap(const std::string& continuumName,
             continue;
     
         // Turn \lambda \in TSE(3) into a Dirichlet value for the continuum
-        const LeafBoundaryPatch<ContinuumGridType>& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
+        const ContinuumLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
         const RigidBodyMotion<dim>& referenceInterface = this->complex_.coupling(couplingName).referenceInterface_;
         
         setRotation(interfaceBoundary, x3d, referenceInterface, it->second);
@@ -509,7 +512,7 @@ continuumDirichletToNeumannMap(const std::string& continuumName,
         if (couplingName.second != continuumName)
             continue;
         
-        const LeafBoundaryPatch<ContinuumGridType>& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
+        const ContinuumLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
 
         VectorType neumannForces(residual.size());
         neumannForces = 0;
@@ -553,7 +556,7 @@ linearizedRodNeumannToDirichletMap(const std::string& rodName,
     //  Assemble the linearized rod problem
     ////////////////////////////////////////////////////
 
-    const LeafBoundaryPatch<RodGridType>& dirichletBoundary = this->complex_.rod(rodName).dirichletBoundary_;
+    const RodLeafBoundaryPatch& dirichletBoundary = this->complex_.rod(rodName).dirichletBoundary_;
 
     typedef Dune::BCRSMatrix<Dune::FieldMatrix<double,6,6> > MatrixType;
     GeodesicFEAssembler<typename RodGridType::LeafGridView, RigidBodyMotion<3> > assembler(dirichletBoundary.gridView(),
@@ -613,11 +616,11 @@ linearizedRodNeumannToDirichletMap(const std::string& rodName,
         }
         
         // Use 'forceTorque' as a Neumann value for the rod
-        const LeafBoundaryPatch<RodGridType>& interfaceBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
+        const RodLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
 
         const typename RodGridType::LeafGridView::IndexSet& indexSet = interfaceBoundary.gridView().indexSet();
         
-        for (typename LeafBoundaryPatch<RodGridType>::iterator bIt = interfaceBoundary.begin();
+        for (typename RodLeafBoundaryPatch::iterator bIt = interfaceBoundary.begin();
              bIt != interfaceBoundary.end();
              ++bIt) {
             
@@ -684,11 +687,11 @@ linearizedRodNeumannToDirichletMap(const std::string& rodName,
         if (couplingName.first != rodName)
             continue;
         
-        const LeafBoundaryPatch<RodGridType>& interfaceBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
+        const RodLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).rodInterfaceBoundary_;
         
         const typename RodGridType::LeafGridView::IndexSet& indexSet = interfaceBoundary.gridView().indexSet();
         
-        for (typename LeafBoundaryPatch<RodGridType>::iterator bIt = interfaceBoundary.begin();
+        for (typename RodLeafBoundaryPatch::iterator bIt = interfaceBoundary.begin();
              bIt != interfaceBoundary.end();
              ++bIt) {
             
@@ -721,7 +724,7 @@ linearizedContinuumNeumannToDirichletMap(const std::string& continuumName,
     * i.e. the linearization at zero
     */
     typedef P1NodalBasis<typename ContinuumGridType::LeafGridView,double> P1Basis;
-    const LeafBoundaryPatch<ContinuumGridType>& dirichletBoundary = this->complex_.continuum(continuumName).dirichletBoundary_;
+    const ContinuumLeafBoundaryPatch& dirichletBoundary = this->complex_.continuum(continuumName).dirichletBoundary_;
     P1Basis basis(dirichletBoundary.gridView());
     OperatorAssembler<P1Basis,P1Basis> assembler(basis, basis);
 
@@ -748,7 +751,7 @@ linearizedContinuumNeumannToDirichletMap(const std::string& continuumName,
             continue;
         
         // Use 'forceTorque' as a Neumann value for the rod
-        const LeafBoundaryPatch<ContinuumGridType>& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
+        const ContinuumLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
         
         // 
         VectorType localNeumannValues;
@@ -798,7 +801,7 @@ linearizedContinuumNeumannToDirichletMap(const std::string& continuumName,
             continue;
         
         // Use 'forceTorque' as a Neumann value for the rod
-        const LeafBoundaryPatch<ContinuumGridType>& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
+        const ContinuumLeafBoundaryPatch& interfaceBoundary = this->complex_.coupling(couplingName).continuumInterfaceBoundary_;
         
         RigidBodyMotion<3> averageInterface;
         computeAverageInterface(interfaceBoundary, x, averageInterface);