diff --git a/dune/tectonic/myblockproblem.hh b/dune/tectonic/myblockproblem.hh
index a6aceca491e5dfad06a3c5f043a03c8bf05ca0bd..d9b01fb0641d46627ac7107ec072d1816ab0cd52 100644
--- a/dune/tectonic/myblockproblem.hh
+++ b/dune/tectonic/myblockproblem.hh
@@ -271,7 +271,7 @@ class MyBlockProblem<MyConvexProblemTypeTEMPLATE>::IterateObject {
         if (j == m)
           localA = &(*it); // localA = A[m][m]
         else
-          Arithmetic::addProduct(localb, -1.0, *it, u[j]);
+          Arithmetic::subtractProduct(localb, *it, u[j]);
       }
       assert(localA != nullptr);
 
diff --git a/src/timestepping/eulerpair.cc b/src/timestepping/eulerpair.cc
index 85cee621ffb022705b336b6cd9338d271e505564..480f86052e75ed439305c4040b8d22e7e18956a9 100644
--- a/src/timestepping/eulerpair.cc
+++ b/src/timestepping/eulerpair.cc
@@ -28,7 +28,7 @@ void EulerPair<VectorType, MatrixType, FunctionType, dim>::setup(
 
   problem_rhs = ell;
   Arithmetic::addProduct(problem_rhs, 1.0 / tau, B, v_old);
-  Arithmetic::addProduct(problem_rhs, -1.0, A, u_old);
+  Arithmetic::subtractProduct(problem_rhs, A, u_old);
 
   // For fixed tau, we'd only really have to do this once
   Dune::MatrixIndexSet indices(A.N(), A.M());
diff --git a/src/timestepping/impliciteuler.cc b/src/timestepping/impliciteuler.cc
index f2dd8be8c7488f9804dc953977c604e6c0210feb..753689da45bb9547179174310cabdf51c14243e0 100644
--- a/src/timestepping/impliciteuler.cc
+++ b/src/timestepping/impliciteuler.cc
@@ -25,7 +25,7 @@ void ImplicitEuler<VectorType, MatrixType, FunctionType, dim>::setup(
   tau = _tau;
 
   problem_rhs = ell;
-  Arithmetic::addProduct(problem_rhs, -1.0, A, u_old);
+  Arithmetic::subtractProduct(problem_rhs, A, u_old);
 
   // For fixed tau, we'd only really have to do this once
   problem_AB = A;
diff --git a/src/timestepping/newmark.cc b/src/timestepping/newmark.cc
index a479b08f8f4df326f1fcf997398a15e147a99fa7..30ad548af61f1b9384d6230e0ce3f3ae2cd80515 100644
--- a/src/timestepping/newmark.cc
+++ b/src/timestepping/newmark.cc
@@ -30,8 +30,8 @@ void Newmark<VectorType, MatrixType, FunctionType, dim>::setup(
   problem_rhs = ell;
   Arithmetic::addProduct(problem_rhs, 2.0 / tau, B, v_old);
   Arithmetic::addProduct(problem_rhs, B, a_old);
-  Arithmetic::addProduct(problem_rhs, -1.0, A, u_old);
-  Arithmetic::addProduct(problem_rhs, -tau / 2.0, A, v_old);
+  Arithmetic::subtractProduct(problem_rhs, A, u_old);
+  Arithmetic::subtractProduct(problem_rhs, tau / 2.0, A, v_old);
 
   // For fixed tau, we'd only really have to do this once
   Dune::MatrixIndexSet indices(A.N(), A.M());
@@ -81,8 +81,8 @@ void Newmark<VectorType, MatrixType, FunctionType, dim>::postProcess(
 
   a = 0;
   Arithmetic::addProduct(a, 2.0 / tau, v);
-  Arithmetic::addProduct(a, -2.0 / tau, v_old);
-  Arithmetic::addProduct(a, -1.0, a_old);
+  Arithmetic::subtractProduct(a, 2.0 / tau, v_old);
+  Arithmetic::subtractProduct(a, 1.0, a_old);
 }
 
 template <class VectorType, class MatrixType, class FunctionType, int dim>