diff --git a/ProcessLib/ComponentTransport/ComponentTransportFEM.h b/ProcessLib/ComponentTransport/ComponentTransportFEM.h
index d0a90d31864427a2d73d1fbdcd7666d81375b86b..974ca4be42e03cc2cd8d96309fef83c5f8a47632 100644
--- a/ProcessLib/ComponentTransport/ComponentTransportFEM.h
+++ b/ProcessLib/ComponentTransport/ComponentTransportFEM.h
@@ -1285,8 +1285,6 @@ public:
     void assembleWithJacobianForStaggeredScheme(
         double const t, double const dt, Eigen::VectorXd const& local_x,
         Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& /*local_M_data*/,
-        std::vector<double>& /*local_K_data*/,
         std::vector<double>& local_b_data,
         std::vector<double>& local_Jac_data) override
     {
diff --git a/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp b/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp
index 55c2cdf83c11681452cce6849e49c9683739d59f..b635727bfa35b3e7623bf71833e9c0af8ab14386 100644
--- a/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp
+++ b/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp
@@ -276,7 +276,7 @@ void ComponentTransportProcess::assembleConcreteProcess(
 void ComponentTransportProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian ComponentTransportProcess.");
 
@@ -296,7 +296,7 @@ void ComponentTransportProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     // b is the negated residumm used in the Newton's method.
     // Here negating b is to recover the primitive residuum.
diff --git a/ProcessLib/ComponentTransport/ComponentTransportProcess.h b/ProcessLib/ComponentTransport/ComponentTransportProcess.h
index 224933c61b98d74893c7f56d82b985c03a75ddff..6ed3b112afdc66c7d5ac9e6bb267c63acc730b37 100644
--- a/ProcessLib/ComponentTransport/ComponentTransportProcess.h
+++ b/ProcessLib/ComponentTransport/ComponentTransportProcess.h
@@ -166,8 +166,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preOutputConcreteProcess(const double t, double const dt,
                                   std::vector<GlobalVector*> const& x,
diff --git a/ProcessLib/HT/HTProcess.cpp b/ProcessLib/HT/HTProcess.cpp
index 63d6ddfcc24f164074ee915fed1f0783b525c43a..2f41777e6cea82a0a160f2ad8a984c7fc7ca41b7 100644
--- a/ProcessLib/HT/HTProcess.cpp
+++ b/ProcessLib/HT/HTProcess.cpp
@@ -112,7 +112,7 @@ void HTProcess::assembleConcreteProcess(
 void HTProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian HTProcess.");
 
@@ -131,7 +131,7 @@ void HTProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 
 std::tuple<NumLib::LocalToGlobalIndexMap*, bool>
diff --git a/ProcessLib/HT/HTProcess.h b/ProcessLib/HT/HTProcess.h
index 950793e04118ff890aaa88f5a11d913d39f00c68..2e7c65b9fa38770512258c0f2d52120ff4a4396f 100644
--- a/ProcessLib/HT/HTProcess.h
+++ b/ProcessLib/HT/HTProcess.h
@@ -98,8 +98,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     /**
      * @copydoc ProcessLib::Process::getDOFTableForExtrapolatorData()
diff --git a/ProcessLib/HeatConduction/HeatConductionFEM.h b/ProcessLib/HeatConduction/HeatConductionFEM.h
index bb93cc4520983ae710bc04980ecc20fe1735994d..4096af46a59809857c7e5a34c42e706b886a2368 100644
--- a/ProcessLib/HeatConduction/HeatConductionFEM.h
+++ b/ProcessLib/HeatConduction/HeatConductionFEM.h
@@ -152,8 +152,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override
     {
diff --git a/ProcessLib/HeatConduction/HeatConductionProcess.cpp b/ProcessLib/HeatConduction/HeatConductionProcess.cpp
index a12e8a11fa0a364749a8564a078fea51b4f44759..727736df4ee2778ac231cf365e060c0fe59c63cc 100644
--- a/ProcessLib/HeatConduction/HeatConductionProcess.cpp
+++ b/ProcessLib/HeatConduction/HeatConductionProcess.cpp
@@ -111,7 +111,7 @@ void HeatConductionProcess::assembleConcreteProcess(
 void HeatConductionProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian HeatConductionProcess.");
 
@@ -121,7 +121,7 @@ void HeatConductionProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     transformVariableFromGlobalVector(b, 0 /*variable id*/,
                                       *_local_to_global_index_map, *_heat_flux,
diff --git a/ProcessLib/HeatConduction/HeatConductionProcess.h b/ProcessLib/HeatConduction/HeatConductionProcess.h
index c0a461d4377e160fb85a317c07baf9d250ca57d2..dc60006530c4daa5649232339b705867f905351f 100644
--- a/ProcessLib/HeatConduction/HeatConductionProcess.h
+++ b/ProcessLib/HeatConduction/HeatConductionProcess.h
@@ -70,8 +70,7 @@ private:
         const double t, double const /*dt*/,
         std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preOutputConcreteProcess(const double t, double const dt,
                                   std::vector<GlobalVector*> const& x,
diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp
index 03878612a980c2b6458edb4107d37387b0508c30..5b2a54b1110f5fd3d1fa6c6ca23693c812efce4a 100644
--- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp
+++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp
@@ -177,7 +177,7 @@ void HeatTransportBHEProcess::assembleConcreteProcess(
 void HeatTransportBHEProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian HeatTransportBHE process.");
 
@@ -188,7 +188,7 @@ void HeatTransportBHEProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 
 void HeatTransportBHEProcess::computeSecondaryVariableConcrete(
diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h
index 328e27a5995cb7ca00ab1aa243ea2bddda50234f..e59d3af65fcc6ecc1cc3cc5d45273c5ef0ea6865 100644
--- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h
+++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h
@@ -63,8 +63,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void createBHEBoundaryConditionTopBottom(
         std::vector<std::vector<MeshLib::Node*>> const& all_bhe_nodes);
diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h
index f26a5e5d950e6a6e2dc2c2e271f834a37480b10e..e4fb2d485fe66eb8490b213f7d3355f195ea4337 100644
--- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h
+++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h
@@ -219,8 +219,6 @@ void HeatTransportBHELocalAssemblerBHE<ShapeFunction, BHEType>::
     assembleWithJacobian(double const t, double const dt,
                          std::vector<double> const& local_x,
                          std::vector<double> const& local_x_prev,
-                         std::vector<double>& local_M_data,
-                         std::vector<double>& local_K_data,
                          std::vector<double>& local_rhs_data,
                          std::vector<double>& local_Jac_data)
 {
@@ -236,6 +234,8 @@ void HeatTransportBHELocalAssemblerBHE<ShapeFunction, BHEType>::
     auto local_rhs = MathLib::createZeroedVector<BheLocalVectorType>(
         local_rhs_data, local_matrix_size);
 
+    std::vector<double> local_M_data(local_Jac_data.size());
+    std::vector<double> local_K_data(local_Jac_data.size());
     assemble(t, dt, local_x, local_x_prev, local_M_data, local_K_data,
              local_rhs_data /*not going to be used*/);
 
diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE.h
index bef12205b7e08942bd714df7c78a390a6f9790ab..2664222a5f696cdae260ddc9d8efcdb71575954e 100644
--- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE.h
+++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE.h
@@ -73,8 +73,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h
index 56351127972e70d75d6d0b4b01468fec850d0d8a..fc1443e76f99470bc426e1ed73dbb11e93b0359c 100644
--- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h
+++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h
@@ -195,9 +195,8 @@ void HeatTransportBHELocalAssemblerSoil<ShapeFunction>::assemble(
 template <typename ShapeFunction>
 void HeatTransportBHELocalAssemblerSoil<ShapeFunction>::assembleWithJacobian(
     double const t, double const dt, std::vector<double> const& local_x,
-    std::vector<double> const& local_x_prev, std::vector<double>& local_M_data,
-    std::vector<double>& local_K_data, std::vector<double>& local_rhs_data,
-    std::vector<double>& local_Jac_data)
+    std::vector<double> const& local_x_prev,
+    std::vector<double>& local_rhs_data, std::vector<double>& local_Jac_data)
 {
     assert(local_x.size() == ShapeFunction::NPOINTS);
     auto const local_matrix_size = local_x.size();
@@ -212,6 +211,8 @@ void HeatTransportBHELocalAssemblerSoil<ShapeFunction>::assembleWithJacobian(
     auto local_rhs = MathLib::createZeroedVector<NodalVectorType>(
         local_rhs_data, local_matrix_size);
 
+    std::vector<double> local_M_data(local_Jac_data.size());
+    std::vector<double> local_K_data(local_Jac_data.size());
     assemble(t, dt, local_x, local_x_prev, local_M_data, local_K_data,
              local_rhs_data /*not going to be used*/);
 
diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h
index 82488e738cdc0d22f2333ef6e4d88219c2960d60..8eed921cf2ed7b6f3bfbdb258c2aad972ac258c6 100644
--- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h
+++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h
@@ -62,8 +62,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
index 07094dc53be39cf1b8325ba54537c86ada07406c..867359cb049095e6f49fe40e8ed028d2c36315c2 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
@@ -103,8 +103,6 @@ void HydroMechanicsLocalAssembler<ShapeFunctionDisplacement,
     assembleWithJacobian(double const t, double const dt,
                          std::vector<double> const& local_x,
                          std::vector<double> const& local_x_prev,
-                         std::vector<double>& /*local_M_data*/,
-                         std::vector<double>& /*local_K_data*/,
                          std::vector<double>& local_rhs_data,
                          std::vector<double>& local_Jac_data)
 {
@@ -791,12 +789,12 @@ template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
           int DisplacementDim>
 void HydroMechanicsLocalAssembler<ShapeFunctionDisplacement,
                                   ShapeFunctionPressure, DisplacementDim>::
-    assembleWithJacobianForStaggeredScheme(
-        const double t, double const dt, Eigen::VectorXd const& local_x,
-        Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& /*local_M_data*/,
-        std::vector<double>& /*local_K_data*/,
-        std::vector<double>& local_b_data, std::vector<double>& local_Jac_data)
+    assembleWithJacobianForStaggeredScheme(const double t, double const dt,
+                                           Eigen::VectorXd const& local_x,
+                                           Eigen::VectorXd const& local_x_prev,
+                                           int const process_id,
+                                           std::vector<double>& local_b_data,
+                                           std::vector<double>& local_Jac_data)
 {
     // For the equations with pressure
     if (process_id == _process_data.hydraulic_process_id)
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM.h
index 45dcdc49d94aaabf990bb8ee22e77e5b847233d6..fdb3010d97a06f564c65adaa8569c7117ec757e7 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsFEM.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM.h
@@ -215,15 +215,12 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override;
 
     void assembleWithJacobianForStaggeredScheme(
         const double t, double const dt, Eigen::VectorXd const& local_x,
         Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& local_M_data, std::vector<double>& local_K_data,
         std::vector<double>& local_b_data,
         std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp b/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp
index 2addf268dd63e359bd9776b69e35bb122ebf1b6a..d5a123a71abe403aaf7f3fe70f58f45e304454b4 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp
+++ b/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp
@@ -308,7 +308,7 @@ void HydroMechanicsProcess<DisplacementDim>::
     assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     // For the monolithic scheme
     bool const use_monolithic_scheme = _process_data.isMonolithicSchemeUsed();
@@ -339,7 +339,7 @@ void HydroMechanicsProcess<DisplacementDim>::
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     auto copyRhs = [&](int const variable_id, auto& output_vector)
     {
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsProcess.h b/ProcessLib/HydroMechanics/HydroMechanicsProcess.h
index 70df7d70b7ec1d568412da189a38faed2fb9c880..775fcd7f49ecdc3ffa59b3964e61e98feedada22 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsProcess.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsProcess.h
@@ -85,8 +85,7 @@ private:
         const double t, double const /*dt*/,
         std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                     double const t, double const dt,
diff --git a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp
index 976676c913d36406a9ead8bd6577d23af75b74c6..dda3ef2d93f3bc3c95adbcaf466c28e32c821d67 100644
--- a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp
+++ b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp
@@ -461,7 +461,7 @@ template <int GlobalDim>
 void HydroMechanicsProcess<GlobalDim>::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian HydroMechanicsProcess.");
 
@@ -471,7 +471,7 @@ void HydroMechanicsProcess<GlobalDim>::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     auto copyRhs = [&](int const variable_id, auto& output_vector)
     {
diff --git a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.h b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.h
index fe8b621d6eebbe8fe088aeef6561a297513cd57d..884f224705aed0bb08a5ca6ca6eaa875c5331745 100644
--- a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.h
+++ b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.h
@@ -73,8 +73,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                     double const t, double const dt,
                                     int const process_id) override;
diff --git a/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerInterface.h b/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerInterface.h
index 34f37086cec87a5c3baccc17142f2b7be43643be..7bb3281a40ae32f50e3c548c4559ea08a37e1777 100644
--- a/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerInterface.h
+++ b/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerInterface.h
@@ -67,8 +67,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x_,
                               std::vector<double> const& local_x_prev_,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_b_data,
                               std::vector<double>& local_Jac_data) override
     {
diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h
index a9c7999bcbc3ec660a87d831b23cf4c44215989d..4c1c98d4539eab80abb452b8641f5173725bc15e 100644
--- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h
+++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h
@@ -41,8 +41,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x_,
                               std::vector<double> const& /*local_x_prev*/,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_b_data,
                               std::vector<double>& local_Jac_data) override
     {
diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix-impl.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix-impl.h
index fc19e7c026c0c95e23d9f7ab3d93f67b5eea3d82..750193bc5ce0e6ecc9c299d2e648614902b76286 100644
--- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix-impl.h
+++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix-impl.h
@@ -94,8 +94,6 @@ void SmallDeformationLocalAssemblerMatrix<ShapeFunction, DisplacementDim>::
     assembleWithJacobian(double const t, double const dt,
                          std::vector<double> const& local_x,
                          std::vector<double> const& /*local_x_prev*/,
-                         std::vector<double>& /*local_M_data*/,
-                         std::vector<double>& /*local_K_data*/,
                          std::vector<double>& local_b_data,
                          std::vector<double>& local_Jac_data)
 {
diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h
index b624622d5ea32c2ddb5e997c481e1b13ddfb44ef..9ea94d6e3153d2aca058138a81299a073f29132e 100644
--- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h
+++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h
@@ -74,8 +74,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& /*local_x_prev*/,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_b_data,
                               std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp
index 6cf8ae896524869f1064143e4ef02fe09ae98526..ce7ae11c8e73b36d6f746837cb46d3f60a407683 100644
--- a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp
+++ b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp
@@ -429,7 +429,7 @@ void SmallDeformationProcess<DisplacementDim>::
     assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian SmallDeformationProcess.");
 
@@ -439,7 +439,7 @@ void SmallDeformationProcess<DisplacementDim>::
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 template <int DisplacementDim>
 void SmallDeformationProcess<DisplacementDim>::preTimestepConcreteProcess(
diff --git a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.h b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.h
index 61f76db9a4d5b893e25f2f9c33c17e7565ba9137..5083459097081ead5b6316f8384c17701f74ef79 100644
--- a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.h
+++ b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.h
@@ -70,8 +70,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                     double const t, double const dt,
diff --git a/ProcessLib/LargeDeformation/LargeDeformationFEM.h b/ProcessLib/LargeDeformation/LargeDeformationFEM.h
index 31c20f977d01a3a6738e7695c0693cd604e579fe..bde4331a96d787a532f57b0927814bdf12dc99cf 100644
--- a/ProcessLib/LargeDeformation/LargeDeformationFEM.h
+++ b/ProcessLib/LargeDeformation/LargeDeformationFEM.h
@@ -242,8 +242,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_b_data,
                               std::vector<double>& local_Jac_data) override
     {
diff --git a/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp b/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp
index 169ad6ed5e88618b2fab3af2e05945422f967fe7..bf5a13101ff98b31c240f041c689c7a87ebb72a4 100644
--- a/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp
+++ b/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp
@@ -145,7 +145,7 @@ void LargeDeformationProcess<DisplacementDim>::assembleConcreteProcess(
         getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K,
         &b);
 
-    _global_output(t, process_id, M, K, b);
+    _global_output(t, process_id, &M, &K, b);
 }
 
 template <int DisplacementDim>
@@ -153,7 +153,7 @@ void LargeDeformationProcess<DisplacementDim>::
     assembleWithJacobianConcreteProcess(
         double const t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian LargeDeformationProcess.");
 
@@ -163,12 +163,12 @@ void LargeDeformationProcess<DisplacementDim>::
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     transformVariableFromGlobalVector(b, 0, *_local_to_global_index_map,
                                       *_nodal_forces, std::negate<double>());
 
-    _global_output(t, process_id, M, K, b, &Jac);
+    _global_output(t, process_id, nullptr, nullptr, b, &Jac);
 }
 
 template <int DisplacementDim>
diff --git a/ProcessLib/LargeDeformation/LargeDeformationProcess.h b/ProcessLib/LargeDeformation/LargeDeformationProcess.h
index c30f829d583554086a97e4fc5982301cbd299101..27e41b7e1a846de86bbf7c91e71ef43a18084b3e 100644
--- a/ProcessLib/LargeDeformation/LargeDeformationProcess.h
+++ b/ProcessLib/LargeDeformation/LargeDeformationProcess.h
@@ -58,8 +58,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void postTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                      std::vector<GlobalVector*> const& x_prev,
diff --git a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp
index af379eab6fcef297bd1415749729ba0b5f0b5e48..8a7ed8a18b84273c6a9078b0e4ece3d603c992ce 100644
--- a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp
+++ b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp
@@ -97,7 +97,7 @@ void LiquidFlowProcess::assembleConcreteProcess(
 void LiquidFlowProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian LiquidFlowProcess.");
 
@@ -108,7 +108,7 @@ void LiquidFlowProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 
 void LiquidFlowProcess::computeSecondaryVariableConcrete(
diff --git a/ProcessLib/LiquidFlow/LiquidFlowProcess.h b/ProcessLib/LiquidFlow/LiquidFlowProcess.h
index 8ac5910005eeb4ec6a5c8ea37724280b5ff5cbfa..b61365e5b6bee6a270174c69522f4f43735780c3 100644
--- a/ProcessLib/LiquidFlow/LiquidFlowProcess.h
+++ b/ProcessLib/LiquidFlow/LiquidFlowProcess.h
@@ -100,8 +100,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     LiquidFlowData _process_data;
 
diff --git a/ProcessLib/PhaseField/PhaseFieldFEM-impl.h b/ProcessLib/PhaseField/PhaseFieldFEM-impl.h
index cccd13e0d994f72d7b7388f5830049aa1632eb0f..56b62e26e9688b6e1f59fed60fd4af319dd7a11b 100644
--- a/ProcessLib/PhaseField/PhaseFieldFEM-impl.h
+++ b/ProcessLib/PhaseField/PhaseFieldFEM-impl.h
@@ -22,8 +22,6 @@ void PhaseFieldLocalAssembler<ShapeFunction, DisplacementDim>::
     assembleWithJacobianForStaggeredScheme(
         double const t, double const dt, Eigen::VectorXd const& local_x,
         Eigen::VectorXd const& /*local_x_prev*/, int const process_id,
-        std::vector<double>& /*local_M_data*/,
-        std::vector<double>& /*local_K_data*/,
         std::vector<double>& local_b_data, std::vector<double>& local_Jac_data)
 {
     // For the equations with phase field.
diff --git a/ProcessLib/PhaseField/PhaseFieldFEM.h b/ProcessLib/PhaseField/PhaseFieldFEM.h
index a3d1e3170c7a7146b86772ad2df2babbf2250178..17bfe2b9f1e8764c671c5c3227814b97890a2185 100644
--- a/ProcessLib/PhaseField/PhaseFieldFEM.h
+++ b/ProcessLib/PhaseField/PhaseFieldFEM.h
@@ -277,7 +277,6 @@ public:
     void assembleWithJacobianForStaggeredScheme(
         double const t, double const dt, Eigen::VectorXd const& local_x,
         Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& local_M_data, std::vector<double>& local_K_data,
         std::vector<double>& local_b_data,
         std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/PhaseField/PhaseFieldProcess.cpp b/ProcessLib/PhaseField/PhaseFieldProcess.cpp
index 1a35046a8fd5642594f668726eb373b27dfc461c..38233cf3db2ae52ae7e4165409d4562329d80eb2 100644
--- a/ProcessLib/PhaseField/PhaseFieldProcess.cpp
+++ b/ProcessLib/PhaseField/PhaseFieldProcess.cpp
@@ -218,7 +218,7 @@ template <int DisplacementDim>
 void PhaseFieldProcess<DisplacementDim>::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
 
@@ -242,7 +242,7 @@ void PhaseFieldProcess<DisplacementDim>::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     if (process_id == 0)
     {
diff --git a/ProcessLib/PhaseField/PhaseFieldProcess.h b/ProcessLib/PhaseField/PhaseFieldProcess.h
index 1d3548728711cde0e60819e7315a7d5b71384aee..9f2a2db958e439bdcf1a7735c60137c63f5727ef 100644
--- a/ProcessLib/PhaseField/PhaseFieldProcess.h
+++ b/ProcessLib/PhaseField/PhaseFieldProcess.h
@@ -70,8 +70,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                     double const t, double const dt,
diff --git a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp
index f034a1c1360a48ec3327643ed722f5535bb38282..69c43e6d5b12ce473c679cefd9a22a1572fe026b 100644
--- a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp
+++ b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp
@@ -84,7 +84,7 @@ void RichardsComponentTransportProcess::assembleConcreteProcess(
 void RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian RichardsComponentTransportProcess.");
 
@@ -95,7 +95,7 @@ void RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 
 }  // namespace RichardsComponentTransport
diff --git a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.h b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.h
index d84d50f6318221c52d0dcfb4d0b1f871d912d5cd..9e9143fa126a78a476d78a79f50cac3caab6362c 100644
--- a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.h
+++ b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.h
@@ -137,8 +137,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     RichardsComponentTransportProcessData _process_data;
 
diff --git a/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp b/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp
index e00e2ac3da6301b48f1ee5c2a412cd46bd6b7623..6b41d1f8394bda1cdcc911ba9e26ab461081b4e4 100644
--- a/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp
+++ b/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp
@@ -78,7 +78,7 @@ void RichardsFlowProcess::assembleConcreteProcess(
 void RichardsFlowProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian RichardsFlowProcess.");
 
@@ -89,7 +89,7 @@ void RichardsFlowProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 
 }  // namespace RichardsFlow
diff --git a/ProcessLib/RichardsFlow/RichardsFlowProcess.h b/ProcessLib/RichardsFlow/RichardsFlowProcess.h
index 70bee188df1e69cdb9fcc3e695f663f5828898a0..dae1ee4bdb119c7069984412cc1c90b32cdedd72 100644
--- a/ProcessLib/RichardsFlow/RichardsFlowProcess.h
+++ b/ProcessLib/RichardsFlow/RichardsFlowProcess.h
@@ -56,8 +56,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     RichardsFlowProcessData _process_data;
 
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
index b225b114f766bfff34116877690666df0b1d13bd..84be48989be244f3402b24a0fc6a967262893bd8 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
@@ -1036,8 +1036,6 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
     assembleWithJacobian(double const t, double const dt,
                          std::vector<double> const& local_x,
                          std::vector<double> const& local_x_prev,
-                         std::vector<double>& /*local_M_data*/,
-                         std::vector<double>& /*local_K_data*/,
                          std::vector<double>& local_rhs_data,
                          std::vector<double>& local_Jac_data)
 {
@@ -1451,8 +1449,6 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
         const double /*t*/, double const /*dt*/,
         Eigen::VectorXd const& /*local_x*/,
         Eigen::VectorXd const& /*local_x_prev*/,
-        std::vector<double>& /*local_M_data*/,
-        std::vector<double>& /*local_K_data*/,
         std::vector<double>& /*local_b_data*/,
         std::vector<double>& /*local_Jac_data*/)
 {
@@ -1467,8 +1463,6 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
         const double /*t*/, double const /*dt*/,
         Eigen::VectorXd const& /*local_x*/,
         Eigen::VectorXd const& /*local_x_prev*/,
-        std::vector<double>& /*local_M_data*/,
-        std::vector<double>& /*local_K_data*/,
         std::vector<double>& /*local_b_data*/,
         std::vector<double>& /*local_Jac_data*/)
 {
@@ -1479,24 +1473,23 @@ template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
           int DisplacementDim>
 void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
                                      ShapeFunctionPressure, DisplacementDim>::
-    assembleWithJacobianForStaggeredScheme(
-        double const t, double const dt, Eigen::VectorXd const& local_x,
-        Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& local_M_data, std::vector<double>& local_K_data,
-        std::vector<double>& local_b_data, std::vector<double>& local_Jac_data)
+    assembleWithJacobianForStaggeredScheme(double const t, double const dt,
+                                           Eigen::VectorXd const& local_x,
+                                           Eigen::VectorXd const& local_x_prev,
+                                           int const process_id,
+                                           std::vector<double>& local_b_data,
+                                           std::vector<double>& local_Jac_data)
 {
     // For the equations with pressure
     if (process_id == 0)
     {
         assembleWithJacobianForPressureEquations(t, dt, local_x, local_x_prev,
-                                                 local_M_data, local_K_data,
                                                  local_b_data, local_Jac_data);
         return;
     }
 
     // For the equations with deformation
     assembleWithJacobianForDeformationEquations(t, dt, local_x, local_x_prev,
-                                                local_M_data, local_K_data,
                                                 local_b_data, local_Jac_data);
 }
 
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h
index f11da1d504d2338eb261a5e7b19fdd3042c0e6f3..4cec7a14c8e6836cc0867b6d09d585cd295c3056 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h
@@ -104,15 +104,12 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override;
 
     void assembleWithJacobianForStaggeredScheme(
         double const t, double const dt, Eigen::VectorXd const& local_x,
         Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& local_M_data, std::vector<double>& local_K_data,
         std::vector<double>& local_b_data,
         std::vector<double>& local_Jac_data) override;
 
@@ -186,19 +183,13 @@ private:
      * @param dt              Time increment
      * @param local_x         Nodal values of \f$x\f$ of an element.
      * @param local_x_prev    Nodal values of \f$x_{prev}\f$ of an element.
-     * @param local_M_data    Mass matrix of an element, which takes the form of
-     *                        \f$ \int N^T N\mathrm{d}\Omega\f$. Not used.
-     * @param local_K_data    Laplacian matrix of an element, which takes the
-     *         form of \f$ \int (\nabla N)^T K \nabla N\mathrm{d}\Omega\f$.
-     *                        Not used.
      * @param local_b_data    Right hand side vector of an element.
      * @param local_Jac_data  Element Jacobian matrix for the Newton-Raphson
      *                        method.
      */
     void assembleWithJacobianForDeformationEquations(
         double const t, double const dt, Eigen::VectorXd const& local_x,
-        Eigen::VectorXd const& local_x_prev, std::vector<double>& local_M_data,
-        std::vector<double>& local_K_data, std::vector<double>& local_b_data,
+        Eigen::VectorXd const& local_x_prev, std::vector<double>& local_b_data,
         std::vector<double>& local_Jac_data);
 
     /**
@@ -218,19 +209,13 @@ private:
      * @param dt              Time increment
      * @param local_x         Nodal values of \f$x\f$ of an element.
      * @param local_x_prev    Nodal values of \f$x_{prev}\f$ of an element.
-     * @param local_M_data    Mass matrix of an element, which takes the form of
-     *                        \f$ \int N^T N\mathrm{d}\Omega\f$. Not used.
-     * @param local_K_data    Laplacian matrix of an element, which takes the
-     *         form of \f$ \int (\nabla N)^T K \nabla N\mathrm{d}\Omega\f$.
-     *                        Not used.
      * @param local_b_data    Right hand side vector of an element.
      * @param local_Jac_data  Element Jacobian matrix for the Newton-Raphson
      *                        method.
      */
     void assembleWithJacobianForPressureEquations(
         double const t, double const dt, Eigen::VectorXd const& local_x,
-        Eigen::VectorXd const& local_x_prev, std::vector<double>& local_M_data,
-        std::vector<double>& local_K_data, std::vector<double>& local_b_data,
+        Eigen::VectorXd const& local_x_prev, std::vector<double>& local_b_data,
         std::vector<double>& local_Jac_data);
 
 private:
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp
index d04fe1f174a4f9f841e327a5111f5ca3cdcbbc23..7d7877c301514d26a757a0d0c3b05d8cf9e296d9 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp
@@ -292,7 +292,7 @@ void RichardsMechanicsProcess<DisplacementDim>::
     assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     // For the monolithic scheme
     if (_use_monolithic_scheme)
@@ -322,7 +322,7 @@ void RichardsMechanicsProcess<DisplacementDim>::
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     auto copyRhs = [&](int const variable_id, auto& output_vector)
     {
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.h b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.h
index 4bde5cffffda9ac37416b338edaca8b912249472..11212f11c619bee0394e0d7f8f54f55be335cd03 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.h
@@ -88,8 +88,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void postTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                      std::vector<GlobalVector*> const& x_prev,
diff --git a/ProcessLib/SmallDeformation/SmallDeformationFEM.h b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
index 4a436dfd6dafd834b2459ad953b11cf1ddf7087c..646d18e319edc5f638dc5a44d86b96d70cacfa33 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationFEM.h
+++ b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
@@ -229,8 +229,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_b_data,
                               std::vector<double>& local_Jac_data) override
     {
diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp
index dc768172fa7bb7de0a8d8dc19ffabced9985aede..b588aef7a42083ef6a8f6ad4dd4d371b011d12cf 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp
+++ b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp
@@ -147,7 +147,7 @@ void SmallDeformationProcess<DisplacementDim>::assembleConcreteProcess(
         getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K,
         &b);
 
-    _global_output(t, process_id, M, K, b);
+    _global_output(t, process_id, &M, &K, b);
 }
 
 template <int DisplacementDim>
@@ -155,7 +155,7 @@ void SmallDeformationProcess<DisplacementDim>::
     assembleWithJacobianConcreteProcess(
         double const t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian SmallDeformationProcess.");
 
@@ -165,12 +165,12 @@ void SmallDeformationProcess<DisplacementDim>::
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     transformVariableFromGlobalVector(b, 0, *_local_to_global_index_map,
                                       *_nodal_forces, std::negate<double>());
 
-    _global_output(t, process_id, M, K, b, &Jac);
+    _global_output(t, process_id, nullptr, nullptr, b, &Jac);
 }
 
 template <int DisplacementDim>
diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcess.h b/ProcessLib/SmallDeformation/SmallDeformationProcess.h
index ec5a51ca5a10bf636299f47120e639ec5366acb5..3693283872a2488db11cfe1773af11d9d73c0da0 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationProcess.h
+++ b/ProcessLib/SmallDeformation/SmallDeformationProcess.h
@@ -58,8 +58,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void postTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                      std::vector<GlobalVector*> const& x_prev,
diff --git a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h
index 5c0b788c637039e6f51855f9919e3de6a5f16108..6661b54d76d67046cb862af609232b829668bd52 100644
--- a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h
+++ b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h
@@ -440,8 +440,6 @@ public:
     void assembleWithJacobian(double const t, double const /*dt*/,
                               std::vector<double> const& local_x,
                               std::vector<double> const& /*local_x_prev*/,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_b_data,
                               std::vector<double>& local_Jac_data) override
     {
diff --git a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp
index cbf799fa30070c81db0148b2081e62517e024eef..83751f61293c3597e39a6d5be897c9e65b959bcd 100644
--- a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp
+++ b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp
@@ -213,7 +213,7 @@ void SmallDeformationNonlocalProcess<DisplacementDim>::
     assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian SmallDeformationNonlocalProcess.");
 
@@ -224,7 +224,7 @@ void SmallDeformationNonlocalProcess<DisplacementDim>::
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     b.copyValues(*_nodal_forces);
     std::transform(_nodal_forces->begin(), _nodal_forces->end(),
diff --git a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.h b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.h
index f505faffa9e5176bfa38dfe3e44c026e3d972034..545feeab1015cd4144e543def680c87269d696e2 100644
--- a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.h
+++ b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.h
@@ -62,8 +62,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void postTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                      std::vector<GlobalVector*> const& x_prev,
diff --git a/ProcessLib/StokesFlow/StokesFlowProcess.cpp b/ProcessLib/StokesFlow/StokesFlowProcess.cpp
index 3fc37e549b22851b1e750f7bc4b4ab189833acdf..b34d904f47c582ea9cabb55167d8573df3a57ecb 100644
--- a/ProcessLib/StokesFlow/StokesFlowProcess.cpp
+++ b/ProcessLib/StokesFlow/StokesFlowProcess.cpp
@@ -152,8 +152,7 @@ void StokesFlowProcess<GlobalDim>::assembleWithJacobianConcreteProcess(
     const double /*t*/, double const /*dt*/,
     std::vector<GlobalVector*> const& /*x*/,
     std::vector<GlobalVector*> const& /*x_prev*/, int const /*process_id*/,
-    GlobalMatrix& /*M*/, GlobalMatrix& /*K*/, GlobalVector& /*b*/,
-    GlobalMatrix& /*Jac*/)
+    GlobalVector& /*b*/, GlobalMatrix& /*Jac*/)
 {
     OGS_FATAL(
         "Assembly of Jacobian matrix has not yet been implemented for "
diff --git a/ProcessLib/StokesFlow/StokesFlowProcess.h b/ProcessLib/StokesFlow/StokesFlowProcess.h
index e6a3126c32d1890f7b5a62d17eb49c7b45ec56f1..fae43b5f7f9e819703f0fac3d26d081f506cb2b3 100644
--- a/ProcessLib/StokesFlow/StokesFlowProcess.h
+++ b/ProcessLib/StokesFlow/StokesFlowProcess.h
@@ -80,8 +80,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     std::vector<MeshLib::Node*> _base_nodes;
     std::unique_ptr<MeshLib::MeshSubset const> _mesh_subset_base_nodes;
diff --git a/ProcessLib/TES/TESProcess.cpp b/ProcessLib/TES/TESProcess.cpp
index 466026afe1554d57fd8c2a5b3e0817a2604e51cf..c2f59694243167588355e201274df44c8d3a0ec6 100644
--- a/ProcessLib/TES/TESProcess.cpp
+++ b/ProcessLib/TES/TESProcess.cpp
@@ -223,7 +223,7 @@ void TESProcess::assembleConcreteProcess(
 void TESProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     std::vector<NumLib::LocalToGlobalIndexMap const*> dof_table = {
         _local_to_global_index_map.get()};
@@ -232,7 +232,7 @@ void TESProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 
 void TESProcess::preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
diff --git a/ProcessLib/TES/TESProcess.h b/ProcessLib/TES/TESProcess.h
index 249021ab363edfac0bb84c077575910d9f1944eb..ab05be02beb4aefab5cb232bb696deb3dcd50825 100644
--- a/ProcessLib/TES/TESProcess.h
+++ b/ProcessLib/TES/TESProcess.h
@@ -68,8 +68,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     GlobalVector const& computeVapourPartialPressure(
         const double t,
diff --git a/ProcessLib/TH2M/TH2MFEM-impl.h b/ProcessLib/TH2M/TH2MFEM-impl.h
index 42d127276f9a65aaa8cf626ac7393123d3c2988d..a39ae5a5af62fa000e698952a8094b6aef37d6da 100644
--- a/ProcessLib/TH2M/TH2MFEM-impl.h
+++ b/ProcessLib/TH2M/TH2MFEM-impl.h
@@ -1210,8 +1210,6 @@ void TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
     assembleWithJacobian(double const t, double const dt,
                          std::vector<double> const& local_x,
                          std::vector<double> const& local_x_prev,
-                         std::vector<double>& /*local_M_data*/,
-                         std::vector<double>& /*local_K_data*/,
                          std::vector<double>& local_rhs_data,
                          std::vector<double>& local_Jac_data)
 {
diff --git a/ProcessLib/TH2M/TH2MFEM.h b/ProcessLib/TH2M/TH2MFEM.h
index bad9ef32fd5fc973b2f19b53b8b47deff84e9d0c..ac63333df313315de7dc6ed23b2f602fb34bb5d9 100644
--- a/ProcessLib/TH2M/TH2MFEM.h
+++ b/ProcessLib/TH2M/TH2MFEM.h
@@ -109,8 +109,6 @@ private:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/TH2M/TH2MProcess.cpp b/ProcessLib/TH2M/TH2MProcess.cpp
index 3162a1876a4a49b6cb207d3697d1ba905f1fdad5..c7f66f32522cb681da24cf5ce46d044c069c9b1a 100644
--- a/ProcessLib/TH2M/TH2MProcess.cpp
+++ b/ProcessLib/TH2M/TH2MProcess.cpp
@@ -252,7 +252,7 @@ template <int DisplacementDim>
 void TH2MProcess<DisplacementDim>::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     if (!_use_monolithic_scheme)
     {
@@ -260,7 +260,7 @@ void TH2MProcess<DisplacementDim>::assembleWithJacobianConcreteProcess(
     }
 
     AssemblyMixin<TH2MProcess<DisplacementDim>>::assembleWithJacobian(
-        t, dt, x, x_prev, process_id, M, K, b, Jac);
+        t, dt, x, x_prev, process_id, b, Jac);
 }
 
 template <int DisplacementDim>
diff --git a/ProcessLib/TH2M/TH2MProcess.h b/ProcessLib/TH2M/TH2MProcess.h
index c1831c7056440fc33644c61ad7e2e44b2ea2efd7..6b8942efd4411dcd75074fe61452df056788ba4a 100644
--- a/ProcessLib/TH2M/TH2MProcess.h
+++ b/ProcessLib/TH2M/TH2MProcess.h
@@ -91,8 +91,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                     double const t, double const dt,
diff --git a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp
index 474a0426a91f2c76bf409f1164dee669300e330f..20165029fdda7c147ce11078f14ef6a307fdab55 100644
--- a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp
+++ b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp
@@ -96,13 +96,13 @@ void ThermalTwoPhaseFlowWithPPProcess::assembleConcreteProcess(
         getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K,
         &b);
 
-    _global_output(t, process_id, M, K, b);
+    _global_output(t, process_id, &M, &K, b);
 }
 
 void ThermalTwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian ThermalTwoPhaseFlowWithPPProcess.");
 
@@ -113,9 +113,9 @@ void ThermalTwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
-    _global_output(t, process_id, M, K, b, &Jac);
+    _global_output(t, process_id, nullptr, nullptr, b, &Jac);
 }
 
 void ThermalTwoPhaseFlowWithPPProcess::preTimestepConcreteProcess(
diff --git a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.h b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.h
index d58e6fc1405bbe62abb3eaaef6c90cf61d719484..9b41114ee3a6c612b83667fe31404cc6881e2229 100644
--- a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.h
+++ b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.h
@@ -71,8 +71,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                     const double t, const double delta_t,
diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h
index ae3e8597c753f64e65bb441d01dfa5781c5b29c0..1475f9c878f6e2c09ed5c121bae0546e4bda26f4 100644
--- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h
+++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h
@@ -549,14 +549,13 @@ ConstitutiveRelationsValues<DisplacementDim> ThermoHydroMechanicsLocalAssembler<
 template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
           int DisplacementDim>
 void ThermoHydroMechanicsLocalAssembler<
-    ShapeFunctionDisplacement, ShapeFunctionPressure, DisplacementDim>::
-    assembleWithJacobian(double const t, double const dt,
-                         std::vector<double> const& local_x,
-                         std::vector<double> const& local_x_prev,
-                         std::vector<double>& /*local_M_data*/,
-                         std::vector<double>& /*local_K_data*/,
-                         std::vector<double>& local_rhs_data,
-                         std::vector<double>& local_Jac_data)
+    ShapeFunctionDisplacement, ShapeFunctionPressure,
+    DisplacementDim>::assembleWithJacobian(double const t, double const dt,
+                                           std::vector<double> const& local_x,
+                                           std::vector<double> const&
+                                               local_x_prev,
+                                           std::vector<double>& local_rhs_data,
+                                           std::vector<double>& local_Jac_data)
 {
     assert(local_x.size() ==
            pressure_size + displacement_size + temperature_size);
diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h
index 2fe17a9a636d6a9c69dbee8731a66861dd4173eb..335757843ea5eda6ecd82139359acfef2b777d7b 100644
--- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h
+++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h
@@ -109,8 +109,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp
index bd2f596277ada1c5e09584edcb98964a1fabd264..2fe1b4fd4cb295445baf7d78c5db42600a59cf8a 100644
--- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp
+++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp
@@ -355,7 +355,7 @@ void ThermoHydroMechanicsProcess<DisplacementDim>::
     assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     // For the monolithic scheme
     if (_use_monolithic_scheme)
@@ -392,7 +392,7 @@ void ThermoHydroMechanicsProcess<DisplacementDim>::
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     auto copyRhs = [&](int const variable_id, auto& output_vector)
     {
diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.h b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.h
index afc726448f27490212bcaf4cea4602a77ef14583..c71b18e4ff8b399fc7642f8fd536e4e8e9d642c8 100644
--- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.h
+++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.h
@@ -87,8 +87,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                     double const t, double const dt,
diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM-impl.h b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM-impl.h
index d2f503df1cf5a7ec517c31b4b63fbe2118c0223b..4f676c7bb7368f987c10e1c2524e40101544440d 100644
--- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM-impl.h
+++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM-impl.h
@@ -19,12 +19,12 @@ namespace ThermoMechanicalPhaseField
 {
 template <typename ShapeFunction, int DisplacementDim>
 void ThermoMechanicalPhaseFieldLocalAssembler<ShapeFunction, DisplacementDim>::
-    assembleWithJacobianForStaggeredScheme(
-        double const t, double const dt, Eigen::VectorXd const& local_x,
-        Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& /*local_M_data*/,
-        std::vector<double>& /*local_K_data*/,
-        std::vector<double>& local_b_data, std::vector<double>& local_Jac_data)
+    assembleWithJacobianForStaggeredScheme(double const t, double const dt,
+                                           Eigen::VectorXd const& local_x,
+                                           Eigen::VectorXd const& local_x_prev,
+                                           int const process_id,
+                                           std::vector<double>& local_b_data,
+                                           std::vector<double>& local_Jac_data)
 {
     if (process_id == _phase_field_process_id)
     {
diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h
index ec0c5fe6fcca55fb2dd4cd72aa4421f734a04401..d40b691dee4258abaf39771b2a0225af77adf409 100644
--- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h
+++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h
@@ -237,7 +237,6 @@ public:
     void assembleWithJacobianForStaggeredScheme(
         double const t, double const dt, Eigen::VectorXd const& local_x,
         Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& local_M_data, std::vector<double>& local_K_data,
         std::vector<double>& local_b_data,
         std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp
index 72a414084ac4b1af7dbefa20e854210095150b2c..fcd0a3d5f9427ffb2fda1b6c2e63fa6775ea3f4d 100644
--- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp
+++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp
@@ -213,7 +213,7 @@ void ThermoMechanicalPhaseFieldProcess<DisplacementDim>::
     assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     // For the staggered scheme
     if (process_id == _mechanics_related_process_id)
@@ -247,7 +247,7 @@ void ThermoMechanicalPhaseFieldProcess<DisplacementDim>::
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 
 template <int DisplacementDim>
diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.h b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.h
index 62f54a89c1d1a5ae4a1a34797f49b67c15eb5894..ec88268d9b10ca476c542c106dae71683bbced00 100644
--- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.h
+++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.h
@@ -101,8 +101,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                     double const t, double const dt,
diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h
index f9e4f0673630710f42eb1e4e2d17c71cb5e347cc..125ac03d23cb4c0b3258338b303697d3dac80e17 100644
--- a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h
+++ b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h
@@ -116,8 +116,6 @@ void ThermoMechanicsLocalAssembler<ShapeFunction, DisplacementDim>::
     assembleWithJacobian(double const t, double const dt,
                          std::vector<double> const& local_x,
                          std::vector<double> const& local_x_prev,
-                         std::vector<double>& /*local_M_data*/,
-                         std::vector<double>& /*local_K_data*/,
                          std::vector<double>& local_rhs_data,
                          std::vector<double>& local_Jac_data)
 {
@@ -329,12 +327,12 @@ void ThermoMechanicsLocalAssembler<ShapeFunction, DisplacementDim>::
 
 template <typename ShapeFunction, int DisplacementDim>
 void ThermoMechanicsLocalAssembler<ShapeFunction, DisplacementDim>::
-    assembleWithJacobianForStaggeredScheme(
-        const double t, double const dt, Eigen::VectorXd const& local_x,
-        Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& /*local_M_data*/,
-        std::vector<double>& /*local_K_data*/,
-        std::vector<double>& local_b_data, std::vector<double>& local_Jac_data)
+    assembleWithJacobianForStaggeredScheme(const double t, double const dt,
+                                           Eigen::VectorXd const& local_x,
+                                           Eigen::VectorXd const& local_x_prev,
+                                           int const process_id,
+                                           std::vector<double>& local_b_data,
+                                           std::vector<double>& local_Jac_data)
 {
     // For the equations with pressure
     if (process_id == _process_data.heat_conduction_process_id)
diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h
index a8a45231d3aeb27ffc70a3aeb21c3af5890e18a6..681b6043ec5c3e987b2763a9d9f13a69dbfa747d 100644
--- a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h
+++ b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h
@@ -182,15 +182,12 @@ public:
     void assembleWithJacobianForStaggeredScheme(
         double const t, double const dt, Eigen::VectorXd const& local_x,
         Eigen::VectorXd const& local_x_prev, int const process_id,
-        std::vector<double>& local_M_data, std::vector<double>& local_K_data,
         std::vector<double>& local_b_data,
         std::vector<double>& local_Jac_data) override;
 
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& local_M_data,
-                              std::vector<double>& local_K_data,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp
index 21fa72010ca777648da7aede5b4a731a9ad7c889..3160fca6646aa153e853fc7010e942907354d041 100644
--- a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp
+++ b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp
@@ -225,7 +225,7 @@ void ThermoMechanicsProcess<DisplacementDim>::
     assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleJacobian ThermoMechanicsProcess.");
 
@@ -273,7 +273,7 @@ void ThermoMechanicsProcess<DisplacementDim>::
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 
     // TODO (naumov): Refactor the copy rhs part. This is copy from HM.
     auto copyRhs = [&](int const variable_id, auto& output_vector)
diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h
index ff35d33c61a2dc874dea53fcb8c460318d4157ce..bd5930b7ac71e990d07c06c7036bf8f10a8eb984 100644
--- a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h
+++ b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h
@@ -77,8 +77,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                     double const t, double const dt,
diff --git a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM-impl.h b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM-impl.h
index 58c8b2724fa54b6312c582450d468d5f6d1a7bd8..3ccc832d9db04f01fcf508671f94a9d5adbf4c3c 100644
--- a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM-impl.h
+++ b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM-impl.h
@@ -155,8 +155,6 @@ void ThermoRichardsFlowLocalAssembler<ShapeFunction, GlobalDim>::
     assembleWithJacobian(double const t, double const dt,
                          std::vector<double> const& local_x,
                          std::vector<double> const& local_x_prev,
-                         std::vector<double>& /*local_M_data*/,
-                         std::vector<double>& /*local_K_data*/,
                          std::vector<double>& local_rhs_data,
                          std::vector<double>& local_Jac_data)
 {
diff --git a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM.h b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM.h
index 8adfb0e0cb83f5e9fe672030b69f4c8514fd0ae3..9bc94a65e9d8acf44cf3e431980c80887b2b15f7 100644
--- a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM.h
+++ b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM.h
@@ -70,8 +70,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp
index 8c88b5ddc93eba166a5471a3304a6889b050877c..762fae976350b3aeb43941726234e9c7350092b0 100644
--- a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp
+++ b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp
@@ -144,7 +144,7 @@ void ThermoRichardsFlowProcess::assembleConcreteProcess(
 void ThermoRichardsFlowProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
 
@@ -154,8 +154,8 @@ void ThermoRichardsFlowProcess::assembleWithJacobianConcreteProcess(
     dof_tables.emplace_back(_local_to_global_index_map.get());
 
     _pvma.assembleWithJacobian(_local_assemblers, getActiveElementIDs(),
-                               dof_tables, t, dt, x, x_prev, process_id, M, K,
-                               b, Jac);
+                               dof_tables, t, dt, x, x_prev, process_id, b,
+                               Jac);
 
     auto copyRhs = [&](int const variable_id, auto& output_vector)
     {
diff --git a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.h b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.h
index 12cbffd2d99955c37767aeea75d74918b9fb0f5d..c4c3cb94266877c58993f03f1902e002377c6059 100644
--- a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.h
+++ b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.h
@@ -122,8 +122,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void postTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
                                      std::vector<GlobalVector*> const& x_prev,
diff --git a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM-impl.h b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM-impl.h
index 3a9c9776320d3184d3037b359871fe0f2ee182c9..421998d8795b94167409a673f0dc9a3a6ed2c2bd 100644
--- a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM-impl.h
+++ b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM-impl.h
@@ -202,8 +202,6 @@ void ThermoRichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
     assembleWithJacobian(double const t, double const dt,
                          std::vector<double> const& local_x,
                          std::vector<double> const& local_x_prev,
-                         std::vector<double>& /*local_M_data*/,
-                         std::vector<double>& /*local_K_data*/,
                          std::vector<double>& local_rhs_data,
                          std::vector<double>& local_Jac_data)
 {
diff --git a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM.h b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM.h
index 410cad16ad4db31795188cb1a3652c78475ffeb8..0bb2215d1ff9188c11d4337d8aefb1d85cac7344 100644
--- a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM.h
+++ b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM.h
@@ -244,8 +244,6 @@ public:
     void assembleWithJacobian(double const t, double const dt,
                               std::vector<double> const& local_x,
                               std::vector<double> const& local_x_prev,
-                              std::vector<double>& /*local_M_data*/,
-                              std::vector<double>& /*local_K_data*/,
                               std::vector<double>& local_rhs_data,
                               std::vector<double>& local_Jac_data) override;
 
diff --git a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp
index 3232236c0772e3406aa4ec7be32e91edd94e01f9..3813235dd1bfc89445f697f364d6ab01ef2ed13a 100644
--- a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp
+++ b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp
@@ -218,14 +218,13 @@ void ThermoRichardsMechanicsProcess<DisplacementDim, ConstitutiveTraits>::
     assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+        GlobalVector& b, GlobalMatrix& Jac)
 {
     AssemblyMixin<ThermoRichardsMechanicsProcess<
         DisplacementDim, ConstitutiveTraits>>::assembleWithJacobian(t, dt, x,
                                                                     x_prev,
                                                                     process_id,
-                                                                    M, K, b,
-                                                                    Jac);
+                                                                    b, Jac);
 }
 
 template <int DisplacementDim, typename ConstitutiveTraits>
diff --git a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.h b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.h
index 45a99603e9c7eb5869641046fffee21c8e543837..a13c2b0d33896f9cb9a390e6d55d541bc53e9925 100644
--- a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.h
+++ b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.h
@@ -184,8 +184,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     void preTimestepConcreteProcess(std::vector<GlobalVector*> const& /*x*/,
                                     const double /*t*/,
diff --git a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp
index 2ea49f2db016aa588a51578e39ec3f434653308c..268193e4644720377e33d844f224edc24cbddacd 100644
--- a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp
+++ b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp
@@ -78,7 +78,7 @@ void TwoPhaseFlowWithPPProcess::assembleConcreteProcess(
 void TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian TwoPhaseFlowWithPPProcess.");
 
@@ -89,7 +89,7 @@ void TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 
 }  // namespace TwoPhaseFlowWithPP
diff --git a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.h b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.h
index bcbf185c829e25560f59caf95052181622107c87..d3e85ca1f42a20cf30797f7693af3808382e137c 100644
--- a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.h
+++ b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.h
@@ -67,8 +67,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     TwoPhaseFlowWithPPProcessData _process_data;
 
diff --git a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp
index 1bbd135aafcfb0da466130541b06865b2872165a..6026da1846e1a5509da11f4fa3112c376573c799 100644
--- a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp
+++ b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp
@@ -78,7 +78,7 @@ void TwoPhaseFlowWithPrhoProcess::assembleConcreteProcess(
 void TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess(
     const double t, double const dt, std::vector<GlobalVector*> const& x,
     std::vector<GlobalVector*> const& x_prev, int const process_id,
-    GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac)
+    GlobalVector& b, GlobalMatrix& Jac)
 {
     DBUG("AssembleWithJacobian TwoPhaseFlowWithPrhoProcess.");
 
@@ -89,7 +89,7 @@ void TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess(
     GlobalExecutor::executeSelectedMemberDereferenced(
         _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
         _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
-        process_id, &M, &K, &b, &Jac);
+        process_id, &b, &Jac);
 }
 
 }  // namespace TwoPhaseFlowWithPrho
diff --git a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.h b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.h
index 14b0b8d8f709e2727e9d05f6911abf5792eaed2f..111a0f5c8a7365f4256e613e8bde6ad36bc75c72 100644
--- a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.h
+++ b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.h
@@ -65,8 +65,7 @@ private:
     void assembleWithJacobianConcreteProcess(
         const double t, double const dt, std::vector<GlobalVector*> const& x,
         std::vector<GlobalVector*> const& x_prev, int const process_id,
-        GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b,
-        GlobalMatrix& Jac) override;
+        GlobalVector& b, GlobalMatrix& Jac) override;
 
     TwoPhaseFlowWithPrhoProcessData _process_data;