diff --git a/src/one-body-sample.cc b/src/one-body-sample.cc
index 0c01069ee8e2cdfa4afce641b63bc9800441edf0..60d96e82394bc0aac2635cf162708e6973a6cf13 100644
--- a/src/one-body-sample.cc
+++ b/src/one-body-sample.cc
@@ -325,6 +325,15 @@ int main(int argc, char *argv[]) {
     auto myGlobalNonlinearity = assemble_nonlinearity<MatrixType, VectorType>(
         frictionalNodes, *nodalIntegrals, frictionData);
 
+    // Problem formulation: right-hand side
+    auto const createRHS = [&](double _time, VectorType &_ell) {
+      assemble_neumann(leafView, p1Assembler, neumannNodes, _ell,
+                       neumannFunction, _time);
+      _ell += gravityFunctional;
+    };
+    VectorType ell(finestSize);
+
+    // {{{ Initial conditions
     VectorType u_initial(finestSize);
     u_initial = 0.0;
     VectorType v_initial(finestSize);
@@ -371,11 +380,6 @@ int main(int argc, char *argv[]) {
 
     auto const timesteps = parset.get<size_t>("timeSteps");
     auto const tau = parset.get<double>("endOfTime") / timesteps;
-    auto const createRHS = [&](double _time, VectorType &_ell) {
-      assemble_neumann(leafView, p1Assembler, neumannNodes, _ell,
-                       neumannFunction, _time);
-      _ell += gravityFunctional;
-    };
 
     VectorType v = v_initial;
     SingletonVectorType alpha = alpha_initial;
@@ -389,8 +393,6 @@ int main(int argc, char *argv[]) {
       timeSteppingScheme->nextTimeStep();
 
       double const time = tau * run;
-
-      VectorType ell(finestSize);
       createRHS(time, ell);
 
       MatrixType problem_AB;