diff --git a/dune/tectonic/localnonlinearity.hh b/dune/tectonic/localnonlinearity.hh
index 97f61a1136ad45e0ff2344970bd683fb7ca228ba..a6ce8e3b27aa091c1a0b560ee7b12d6364ec30de 100644
--- a/dune/tectonic/localnonlinearity.hh
+++ b/dune/tectonic/localnonlinearity.hh
@@ -81,6 +81,7 @@ template <int dimension> class LocalNonlinearity {
       A[k][k] += (h2 - h1ox) * x[k] * x[k] / normX2 + h1ox;
   }
 
+  // TODO: do not evaluate at zero
   void addGradient(VectorType const &x, VectorType &y) const {
     if (x == VectorType(0))
       return;
@@ -90,11 +91,13 @@ template <int dimension> class LocalNonlinearity {
     y += tmp;
   }
 
+  // TODO: do not evaluate at zero
   void upperGradient(VectorType const &x, VectorType &ret) const {
     ret = x;
     ret *= func_->rightDifferential(x.two_norm()) / x.two_norm();
   }
 
+  // TODO: do not evaluate at zero
   void lowerGradient(VectorType const &x, VectorType &ret) const {
     ret = x;
     ret *= func_->leftDifferential(x.two_norm()) / x.two_norm();
diff --git a/dune/tectonic/myblockproblem.hh b/dune/tectonic/myblockproblem.hh
index ecd7b7b9b18906627d3e72f956a4bcd19809448f..f79ee7b212d681c64e14bf2963e09c06710fec59 100644
--- a/dune/tectonic/myblockproblem.hh
+++ b/dune/tectonic/myblockproblem.hh
@@ -86,6 +86,7 @@ template <class MyConvexProblemTypeTEMPLATE> class MyBlockProblem {
     Interval<double> ab;
     for (int i = 0; i < u.size(); ++i) {
       for (int j = 0; j < block_size; ++j) {
+        // TODO: handle smooth domain
         // problem.phi.smoothDomain(i, u[i][j], regularity_tol, ab, j);
 
         // // if (phi is not regular at u) or (u is not in smooth domain) or
@@ -102,6 +103,7 @@ template <class MyConvexProblemTypeTEMPLATE> class MyBlockProblem {
     // construct sparsity pattern for linearization
     Dune::MatrixIndexSet indices(problem.A.N(), problem.A.M());
     indices.import(problem.A);
+    // TODO: handle smooth domain
     // problem.phi.addHessianIndices(indices);
 
     // construct matrix from pattern and initialize it
diff --git a/src/TODO.org b/src/TODO.org
new file mode 100644
index 0000000000000000000000000000000000000000..012f29b2d59743b278015f0eb6abd48c1baea911
--- /dev/null
+++ b/src/TODO.org
@@ -0,0 +1,4 @@
+* compute damping parameter
+* investigate why laursen fails
+* handle smooth domain properly
+* fix up octave bindings (low-pri)
diff --git a/src/one-body-sample.cc b/src/one-body-sample.cc
index 478d22586f566a952168e645626706e505ced5f9..564b7f5d0744435331c3dc5bd14b49c08a187fe6 100644
--- a/src/one-body-sample.cc
+++ b/src/one-body-sample.cc
@@ -244,13 +244,13 @@ int main(int argc, char *argv[]) {
     nonlinearGSStep.ignoreNodes_ = &ignoreNodes;
 
     VectorType u1(grid.size(grid.maxLevel(), dim));
-    u1 = 0;
+    u1 = 0.0; // Has to be zero!
     VectorType u2 = u1;
     VectorType u3 = u1;
     VectorType u4 = u1;
 
     VectorType u1_diff_new(grid.size(grid.maxLevel(), dim));
-    u1_diff_new = 0;
+    u1_diff_new = 0.0; // Has to be zero!
     VectorType u2_diff_new = u1_diff_new;
     VectorType u3_diff_new = u1_diff_new;
     VectorType u4_diff_new = u1_diff_new;
diff --git a/src/one-body-sample.parset b/src/one-body-sample.parset
index d49f051d8a2c6f3e46e59aef31a109ea2285973d..1a1c33701a0054f95d41aee3d2f413b04b64b708 100644
--- a/src/one-body-sample.parset
+++ b/src/one-body-sample.parset
@@ -39,7 +39,7 @@ normalstress = 0.1 # laursen depends a lot more on this
 #    http://earthquake.usgs.gov/research/physics/lab/prediction.pdf
 mu = 0.5
 eta = 1
-model = Ruina
+model = Ruina # TODO: Laursen doesn't work with TNNMG
 
 [boundary.friction.ruina]
 # "For rocks, typical values of A and B range from 0.005 to 0.015"