diff --git a/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h b/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h
index 072090486b64fafd0179298130c4239d888999d1..64b7e28402c88dd592ef2cfb1bb1da85216d28ef 100644
--- a/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h
+++ b/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h
@@ -166,7 +166,7 @@ public:
                 local_x, ShapeFunction::NPOINTS);
 
         cache.clear();
-        auto cache_vec = MathLib::createZeroedMatrix<
+        auto cache_mat = MathLib::createZeroedMatrix<
             Eigen::Matrix<double, GlobalDim, Eigen::Dynamic, Eigen::RowMajor>>(
             cache, GlobalDim, num_intpts);
 
@@ -178,7 +178,7 @@ public:
             pos.setIntegrationPoint(i);
             auto const k = _process_data.hydraulic_conductivity(t, pos)[0];
             // dimensions: (d x 1) = (d x n) * (n x 1)
-            cache_vec.col(i).noalias() =
+            cache_mat.col(i).noalias() =
                 -k * _shape_matrices[i].dNdx * local_x_vec;
         }
 
diff --git a/ProcessLib/HT/HTFEM.h b/ProcessLib/HT/HTFEM.h
index c39ef860601c8935cc5a701ed2ab3fbbdb4e8f0b..bb0a60d26d3d079bb7bd5b6bb13191ba83cb033b 100644
--- a/ProcessLib/HT/HTFEM.h
+++ b/ProcessLib/HT/HTFEM.h
@@ -280,7 +280,7 @@ public:
         auto const local_x = current_solution.get(indices);
 
         cache.clear();
-        auto cache_vec = MathLib::createZeroedMatrix<
+        auto cache_mat = MathLib::createZeroedMatrix<
             Eigen::Matrix<double, GlobalDim, Eigen::Dynamic, Eigen::RowMajor>>(
             cache, GlobalDim, num_intpts);
 
@@ -317,7 +317,7 @@ public:
                 MaterialLib::Fluid::FluidPropertyType::Viscosity, vars);
             GlobalDimMatrixType const K_over_mu = K / mu;
 
-            cache_vec.col(ip).noalias() = -K_over_mu * dNdx * p_nodal_values;
+            cache_mat.col(ip).noalias() = -K_over_mu * dNdx * p_nodal_values;
 
             if (_process_data.has_gravity)
             {
@@ -325,7 +325,7 @@ public:
                     MaterialLib::Fluid::FluidPropertyType::Density, vars);
                 auto const b = _process_data.specific_body_force;
                 // here it is assumed that the vector b is directed 'downwards'
-                cache_vec.col(ip).noalias() += K_over_mu * rho_w * b;
+                cache_mat.col(ip).noalias() += K_over_mu * rho_w * b;
             }
         }
 
diff --git a/ProcessLib/TES/TESLocalAssembler-impl.h b/ProcessLib/TES/TESLocalAssembler-impl.h
index 746938045eb0ed4153a2e7ae3413b02e6fbd750a..1f562bd092e25b5a8a21aba50217c4fb243f55d2 100644
--- a/ProcessLib/TES/TESLocalAssembler-impl.h
+++ b/ProcessLib/TES/TESLocalAssembler-impl.h
@@ -264,14 +264,12 @@ TESLocalAssembler<ShapeFunction_, IntegrationMethod_, GlobalDim>::
     assert(!indices.empty());
     auto const local_x = current_solution.get(indices);
     // local_x is ordered by component, local_x_mat is row major
-    // TODO fixed size matrix ShapeFunction_::NPOINTS columns. Conflicts with
-    // RowMajor for (presumably) 0D element.
     auto const local_x_mat = MathLib::toMatrix<
         Eigen::Matrix<double, NODAL_DOF, Eigen::Dynamic, Eigen::RowMajor>>(
         local_x, NODAL_DOF, ShapeFunction_::NPOINTS);
 
     cache.clear();
-    auto cache_vec = MathLib::createZeroedMatrix<
+    auto cache_mat = MathLib::createZeroedMatrix<
         Eigen::Matrix<double, GlobalDim, Eigen::Dynamic, Eigen::RowMajor>>(
         cache, GlobalDim, num_intpts);
 
@@ -283,7 +281,7 @@ TESLocalAssembler<ShapeFunction_, IntegrationMethod_, GlobalDim>::
         const double eta_GR = fluid_viscosity(p, T, x);
 
         auto const& k = _d.getAssemblyParameters().solid_perm_tensor;
-        cache_vec.col(i).noalias() =
+        cache_mat.col(i).noalias() =
             k * (_shape_matrices[i].dNdx *
                  local_x_mat.row(COMPONENT_ID_PRESSURE).transpose()) /
             -eta_GR;