diff --git a/dirneucoupling.cc b/dirneucoupling.cc
index 31bffb14fd5f763d8c06aa85d09062014da5cd0d..55dfb41f10993e9c9ff20c85aa9923fa9458eea6 100644
--- a/dirneucoupling.cc
+++ b/dirneucoupling.cc
@@ -392,7 +392,7 @@ int main (int argc, char *argv[]) try
             // Using that index 0 is always the left boundary for a uniformly refined OneDGrid
             computeAveragePressure<GridType::LevelGridView>(resultantForce, resultantTorque, 
                                               interfaceBoundary[grid.maxLevel()], 
-                                              rodX[0],
+                                              rodX[0].r,
                                               neumannValues);
 
             rhs3d = 0;
@@ -520,7 +520,7 @@ int main (int argc, char *argv[]) try
                 // Using that index 0 is always the left boundary for a uniformly refined OneDGrid
                 computeAveragePressure<GridType::LevelGridView>(resultantForce, resultantTorque, 
                                               interfaceBoundary[grid.maxLevel()], 
-                                              rodX[0],
+                                              rodX[0].r,
                                               neumannValues);
 
                 rhs3d = 0;
diff --git a/dune/gfe/averageinterface.hh b/dune/gfe/averageinterface.hh
index 944633a34f5ad5c6b79a89688a40add6c7c4691e..ef479816d852486a0959cc8c513129f318b0972c 100644
--- a/dune/gfe/averageinterface.hh
+++ b/dune/gfe/averageinterface.hh
@@ -472,7 +472,7 @@ template <class GridView>
 void computeAveragePressure(const Dune::FieldVector<double,GridView::dimension>& resultantForce,
                             const Dune::FieldVector<double,GridView::dimension>& resultantTorque,
                             const BoundaryPatchBase<GridView>& interface,
-                            const RigidBodyMotion<3>& crossSection,
+                            const Dune::FieldVector<double,GridView::dimension>& centerOfTorque,
                             Dune::BlockVector<Dune::FieldVector<double, GridView::dimension> >& pressure)
 {
     const GridView& gridView                    = interface.gridView();
@@ -567,7 +567,7 @@ void computeAveragePressure(const Dune::FieldVector<double,GridView::dimension>&
                         phi_i[j] = shapeFunctionValues[i];
                         
                         mu_tilde[i][j].axpy(quad[qp].weight() * integrationElement,
-                                            crossProduct(Dune::FieldVector<double,dim>(worldPos-crossSection.r), phi_i));
+                                            crossProduct(Dune::FieldVector<double,dim>(worldPos-centerOfTorque), phi_i));
 
                     }
                     
@@ -653,7 +653,7 @@ void computeAveragePressure(const Dune::FieldVector<double,GridView::dimension>&
 
     Dune::FieldVector<double,3> outputForce(0), outputTorque(0);
 
-    computeTotalForceAndTorque(interface, pressure, crossSection.r, outputForce, outputTorque);
+    computeTotalForceAndTorque(interface, pressure, centerOfTorque, outputForce, outputTorque);
 
     outputForce  -= resultantForce;
     outputTorque -= resultantTorque;