diff --git a/dune/tectonic/globalnonlinearity.hh b/dune/tectonic/globalnonlinearity.hh
index 5a2f582666f3216a053b7e4e5f5983a0559b72fe..ffcd617ec301cb98d00c81ef2a69e0e3f1187a44 100644
--- a/dune/tectonic/globalnonlinearity.hh
+++ b/dune/tectonic/globalnonlinearity.hh
@@ -18,10 +18,23 @@ template <int dim, class OuterFunctionType> class GlobalNonlinearity {
         normalStress(normalStress),
         nodalIntegrals(nodalIntegrals) {}
 
+  /*
+    Return a restriction of the outer function to the i'th node. If
+    mu and sigma_n denote the coefficient of friction and the normal
+    stress, respectively, at the i'th node, this function is given
+    by
+
+    sigma_n [(1/eta \bar Gamma)* + mu id]
+
+    TODO: We chose Gamma = id, so that (\bar Gamma)* = \Gamma^{-1}
+    = id^{-1} = id. The factor 1/eta cancels in this special case, leaving us
+    with
+
+    sigma_n [id + mu id] = sigma_n (1 + mu) id
+  */
   void restriction(int i, OuterFunctionType &f) const {
-    double coefficient = normalStress[i] * nodalIntegrals[i];
-    // FIXME: Assume Gamma = id and h_{n+1} = 1 for now;
-    // We then only have to evaluate (1 + F_xi)(|x|)
+    double coefficient = nodalIntegrals[i][0];
+    coefficient *= normalStress[i];
     coefficient *= 1 + coefficientOfFriction[i];
     f = OuterFunctionType(coefficient);
   }