From 66dfead5e2ef24e85ebaaef1f19b8977aae1be71 Mon Sep 17 00:00:00 2001
From: Dmitri Naumov <dmitri.naumov@ufz.de>
Date: Wed, 6 Nov 2019 18:03:51 +0100
Subject: [PATCH] [PL] Update the secondary vars fcts interfaces.

And only interfaces to separate the relevant code changes;
most processes do not need any changes.
---
 .../ComponentTransportFEM.h                   |  8 +--
 .../GroundwaterFlow/GroundwaterFlowFEM.h      | 12 ++--
 ProcessLib/HT/HTLocalAssemblerInterface.h     |  4 +-
 ProcessLib/HT/MonolithicHTFEM.h               |  4 +-
 ProcessLib/HT/StaggeredHTFEM-impl.h           |  9 +--
 ProcessLib/HT/StaggeredHTFEM.h                |  4 +-
 ProcessLib/HeatConduction/HeatConductionFEM.h | 36 +++++-----
 .../HydroMechanics/HydroMechanicsFEM-impl.h   |  9 +--
 ProcessLib/HydroMechanics/HydroMechanicsFEM.h | 52 +++++++-------
 .../HydroMechanics/LocalAssemblerInterface.h  | 54 +++++++-------
 .../SmallDeformationLocalAssemblerFracture.h  | 48 ++++++-------
 .../SmallDeformationLocalAssemblerInterface.h | 48 ++++++-------
 .../SmallDeformationLocalAssemblerMatrix.h    | 48 ++++++-------
 ...ormationLocalAssemblerMatrixNearFracture.h | 48 ++++++-------
 .../LiquidFlowLocalAssembler-impl.h           |  9 +--
 .../LiquidFlow/LiquidFlowLocalAssembler.h     | 12 ++--
 .../PhaseField/LocalAssemblerInterface.h      | 12 ++--
 ProcessLib/PhaseField/PhaseFieldFEM.h         |  8 +--
 .../RichardsComponentTransportFEM-impl.h      | 18 ++---
 .../RichardsComponentTransportFEM.h           | 21 +++---
 ProcessLib/RichardsFlow/RichardsFlowFEM.h     | 24 +++----
 .../LocalAssemblerInterface.h                 | 24 +++----
 .../RichardsMechanicsFEM-impl.h               | 39 +++++-----
 .../RichardsMechanics/RichardsMechanicsFEM.h  | 24 +++----
 .../LocalAssemblerInterface.h                 | 18 ++---
 .../SmallDeformation/SmallDeformationFEM.h    | 12 ++--
 .../SmallDeformationProcess.cpp               |  5 +-
 .../LocalAssemblerInterface.h                 | 36 +++++-----
 .../SmallDeformationNonlocalFEM.h             | 24 +++----
 ProcessLib/TES/TESLocalAssembler-impl.h       | 40 ++++++-----
 ProcessLib/TES/TESLocalAssembler.h            | 72 +++++++++----------
 ProcessLib/TES/TESProcess.cpp                 | 29 ++++----
 ProcessLib/TES/TESProcess.h                   | 12 ++--
 .../ThermalTwoPhaseFlowWithPPLocalAssembler.h | 24 +++----
 .../LocalAssemblerInterface.h                 | 18 ++---
 .../ThermoHydroMechanicsFEM-impl.h            | 12 ++--
 .../ThermoHydroMechanicsFEM.h                 | 12 ++--
 .../LocalAssemblerInterface.h                 | 18 ++---
 .../ThermoMechanicalPhaseFieldFEM.h           | 12 ++--
 .../ThermoMechanics/LocalAssemblerInterface.h | 13 ++--
 .../ThermoMechanics/ThermoMechanicsFEM-impl.h | 18 ++---
 .../ThermoMechanics/ThermoMechanicsFEM.h      | 12 ++--
 .../TwoPhaseFlowWithPPLocalAssembler.h        | 24 +++----
 .../TwoPhaseFlowWithPrhoLocalAssembler.h      | 24 +++----
 Tests/NumLib/TestExtrapolation.cpp            | 33 ++++-----
 45 files changed, 530 insertions(+), 513 deletions(-)

diff --git a/ProcessLib/ComponentTransport/ComponentTransportFEM.h b/ProcessLib/ComponentTransport/ComponentTransportFEM.h
index 671b8044da3..d326a0ebb93 100644
--- a/ProcessLib/ComponentTransport/ComponentTransportFEM.h
+++ b/ProcessLib/ComponentTransport/ComponentTransportFEM.h
@@ -66,8 +66,8 @@ public:
 
     virtual std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
 protected:
@@ -706,8 +706,8 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override
     {
         auto const indices = NumLib::getIndices(_element.getID(), dof_table);
diff --git a/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h b/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h
index 9b6ffce86ff..82d9b6b72bc 100644
--- a/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h
+++ b/ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h
@@ -63,10 +63,10 @@ class GroundwaterFlowLocalAssemblerInterface
 {
 public:
     virtual std::vector<double> const& getIntPtDarcyVelocity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 };
 
 template <typename ShapeFunction, typename IntegrationMethod,
@@ -175,8 +175,8 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override
     {
         auto const n_integration_points =
diff --git a/ProcessLib/HT/HTLocalAssemblerInterface.h b/ProcessLib/HT/HTLocalAssemblerInterface.h
index 54b4569465b..27057b023aa 100644
--- a/ProcessLib/HT/HTLocalAssemblerInterface.h
+++ b/ProcessLib/HT/HTLocalAssemblerInterface.h
@@ -55,8 +55,8 @@ public:
 
     virtual std::vector<double> const& getIntPtDarcyVelocity(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const = 0;
 
     Eigen::Vector3d getFlux(MathLib::Point3d const& pnt_local_coords,
diff --git a/ProcessLib/HT/MonolithicHTFEM.h b/ProcessLib/HT/MonolithicHTFEM.h
index d1b5db88c6f..8348a447d55 100644
--- a/ProcessLib/HT/MonolithicHTFEM.h
+++ b/ProcessLib/HT/MonolithicHTFEM.h
@@ -204,8 +204,8 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override
     {
         auto const indices =
diff --git a/ProcessLib/HT/StaggeredHTFEM-impl.h b/ProcessLib/HT/StaggeredHTFEM-impl.h
index 71a10acc432..e3302dc990f 100644
--- a/ProcessLib/HT/StaggeredHTFEM-impl.h
+++ b/ProcessLib/HT/StaggeredHTFEM-impl.h
@@ -293,10 +293,11 @@ template <typename ShapeFunction, typename IntegrationMethod,
           unsigned GlobalDim>
 std::vector<double> const&
 StaggeredHTFEM<ShapeFunction, IntegrationMethod, GlobalDim>::
-    getIntPtDarcyVelocity(const double t,
-                          GlobalVector const& /*current_solution*/,
-                          NumLib::LocalToGlobalIndexMap const& dof_table,
-                          std::vector<double>& cache) const
+    getIntPtDarcyVelocity(
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const
 {
     auto const indices = NumLib::getIndices(this->_element.getID(), dof_table);
     assert(!indices.empty());
diff --git a/ProcessLib/HT/StaggeredHTFEM.h b/ProcessLib/HT/StaggeredHTFEM.h
index 76a84c5a1be..8ce6c973de5 100644
--- a/ProcessLib/HT/StaggeredHTFEM.h
+++ b/ProcessLib/HT/StaggeredHTFEM.h
@@ -79,8 +79,8 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
 private:
diff --git a/ProcessLib/HeatConduction/HeatConductionFEM.h b/ProcessLib/HeatConduction/HeatConductionFEM.h
index ded19469066..ba6e9097356 100644
--- a/ProcessLib/HeatConduction/HeatConductionFEM.h
+++ b/ProcessLib/HeatConduction/HeatConductionFEM.h
@@ -34,22 +34,22 @@ class HeatConductionLocalAssemblerInterface
 {
 public:
     virtual std::vector<double> const& getIntPtHeatFluxX(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtHeatFluxY(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtHeatFluxZ(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 };
 
 template <typename ShapeFunction, typename IntegrationMethod,
@@ -171,8 +171,8 @@ public:
 
     std::vector<double> const& getIntPtHeatFluxX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(!_heat_fluxes.empty());
@@ -181,8 +181,8 @@ public:
 
     std::vector<double> const& getIntPtHeatFluxY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(_heat_fluxes.size() > 1);
@@ -191,8 +191,8 @@ public:
 
     std::vector<double> const& getIntPtHeatFluxZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(_heat_fluxes.size() > 2);
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
index 50e9740d1cd..78879632db5 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
@@ -301,10 +301,11 @@ template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
 std::vector<double> const&
 HydroMechanicsLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
                              IntegrationMethod, DisplacementDim>::
-    getIntPtDarcyVelocity(const double t,
-                          GlobalVector const& current_solution,
-                          NumLib::LocalToGlobalIndexMap const& dof_table,
-                          std::vector<double>& cache) const
+    getIntPtDarcyVelocity(
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const
 {
     auto const num_intpts = _ip_data.size();
 
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM.h
index 044b7e18781..c3cfe2fade5 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsFEM.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM.h
@@ -215,8 +215,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 0);
@@ -224,8 +224,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 1);
@@ -233,8 +233,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 2);
@@ -242,8 +242,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 3);
@@ -251,8 +251,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -261,8 +261,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -271,8 +271,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 0);
@@ -280,8 +280,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 1);
@@ -289,8 +289,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 2);
@@ -298,8 +298,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 3);
@@ -307,8 +307,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -317,8 +317,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -327,8 +327,8 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
 private:
diff --git a/ProcessLib/HydroMechanics/LocalAssemblerInterface.h b/ProcessLib/HydroMechanics/LocalAssemblerInterface.h
index 70a9db0bda2..c9f999626b5 100644
--- a/ProcessLib/HydroMechanics/LocalAssemblerInterface.h
+++ b/ProcessLib/HydroMechanics/LocalAssemblerInterface.h
@@ -31,80 +31,80 @@ struct LocalAssemblerInterface
 
     virtual std::vector<double> const& getIntPtSigmaXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtDarcyVelocity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 };
 
diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerFracture.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerFracture.h
index 680afc0d621..2ee9c736111 100644
--- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerFracture.h
+++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerFracture.h
@@ -107,8 +107,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 0);
@@ -116,8 +116,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 1);
@@ -125,8 +125,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 2);
@@ -134,8 +134,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 3);
@@ -143,8 +143,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -153,8 +153,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -163,8 +163,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 0);
@@ -172,8 +172,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 1);
@@ -181,8 +181,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 2);
@@ -190,8 +190,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 3);
@@ -199,8 +199,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -209,8 +209,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h
index 46e9c2ee88f..79a8fe37f4f 100644
--- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h
+++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h
@@ -104,74 +104,74 @@ public:
 
     virtual std::vector<double> const& getIntPtSigmaXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilonYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 
 protected:
diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h
index bced565a7e5..c517558d42b 100644
--- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h
+++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h
@@ -106,8 +106,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 0);
@@ -115,8 +115,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 1);
@@ -124,8 +124,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 2);
@@ -133,8 +133,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 3);
@@ -142,8 +142,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -152,8 +152,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -162,8 +162,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 0);
@@ -171,8 +171,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 1);
@@ -180,8 +180,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 2);
@@ -189,8 +189,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 3);
@@ -198,8 +198,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -208,8 +208,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrixNearFracture.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrixNearFracture.h
index 4e5cd530eeb..94e67497d7a 100644
--- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrixNearFracture.h
+++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrixNearFracture.h
@@ -109,8 +109,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 0);
@@ -118,8 +118,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 1);
@@ -127,8 +127,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 2);
@@ -136,8 +136,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtSigma(cache, 3);
@@ -145,8 +145,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -155,8 +155,8 @@ public:
 
     std::vector<double> const& getIntPtSigmaYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -165,8 +165,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 0);
@@ -174,8 +174,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonYY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 1);
@@ -183,8 +183,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonZZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 2);
@@ -192,8 +192,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXY(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         return getIntPtEpsilon(cache, 3);
@@ -201,8 +201,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonXZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
@@ -211,8 +211,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilonYZ(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         assert(DisplacementDim == 3);
diff --git a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h
index c1ec1f3c742..c8d810787c8 100644
--- a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h
+++ b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h
@@ -122,10 +122,11 @@ template <typename ShapeFunction, typename IntegrationMethod,
           unsigned GlobalDim>
 std::vector<double> const&
 LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
-    getIntPtDarcyVelocity(const double t,
-                          GlobalVector const& current_solution,
-                          NumLib::LocalToGlobalIndexMap const& dof_table,
-                          std::vector<double>& velocity_cache) const
+    getIntPtDarcyVelocity(
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& velocity_cache) const
 {
     auto const indices = NumLib::getIndices(_element.getID(), dof_table);
     auto const local_x = current_solution.get(indices);
diff --git a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h
index 4fe16c01425..08f6657186d 100644
--- a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h
+++ b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h
@@ -61,10 +61,10 @@ class LiquidFlowLocalAssemblerInterface
 {
 public:
     virtual std::vector<double> const& getIntPtDarcyVelocity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 };
 
 template <typename ShapeFunction, typename IntegrationMethod,
@@ -139,8 +139,8 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& velocity_cache) const override;
 
 private:
diff --git a/ProcessLib/PhaseField/LocalAssemblerInterface.h b/ProcessLib/PhaseField/LocalAssemblerInterface.h
index 6e7900c2186..eb982b0b90f 100644
--- a/ProcessLib/PhaseField/LocalAssemblerInterface.h
+++ b/ProcessLib/PhaseField/LocalAssemblerInterface.h
@@ -24,15 +24,15 @@ struct PhaseFieldLocalAssemblerInterface
       public NumLib::ExtrapolatableElement
 {
     virtual std::vector<double> const& getIntPtSigma(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilon(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual void computeCrackIntegral(
diff --git a/ProcessLib/PhaseField/PhaseFieldFEM.h b/ProcessLib/PhaseField/PhaseFieldFEM.h
index 14bd3e546c1..4e4ce9cdccb 100644
--- a/ProcessLib/PhaseField/PhaseFieldFEM.h
+++ b/ProcessLib/PhaseField/PhaseFieldFEM.h
@@ -277,8 +277,8 @@ public:
 private:
     std::vector<double> const& getIntPtSigma(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         static const int kelvin_vector_size =
@@ -303,8 +303,8 @@ private:
 
     std::vector<double> const& getIntPtEpsilon(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         auto const kelvin_vector_size =
diff --git a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM-impl.h b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM-impl.h
index 179256b9e74..4c57dd869aa 100644
--- a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM-impl.h
+++ b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM-impl.h
@@ -214,10 +214,11 @@ template <typename ShapeFunction, typename IntegrationMethod,
           unsigned GlobalDim>
 std::vector<double> const&
 LocalAssemblerData<ShapeFunction, IntegrationMethod, GlobalDim>::
-    getIntPtDarcyVelocity(const double t,
-                          GlobalVector const& current_solution,
-                          NumLib::LocalToGlobalIndexMap const& dof_table,
-                          std::vector<double>& cache) const
+    getIntPtDarcyVelocity(
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const
 {
     unsigned const n_integration_points =
         _integration_method.getNumberOfPoints();
@@ -302,10 +303,11 @@ template <typename ShapeFunction, typename IntegrationMethod,
           unsigned GlobalDim>
 std::vector<double> const&
 LocalAssemblerData<ShapeFunction, IntegrationMethod, GlobalDim>::
-    getIntPtSaturation(const double t,
-                       GlobalVector const& current_solution,
-                       NumLib::LocalToGlobalIndexMap const& dof_table,
-                       std::vector<double>& cache) const
+    getIntPtSaturation(
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const
 {
     ParameterLib::SpatialPosition pos;
     pos.setElementID(_element_id);
diff --git a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM.h b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM.h
index 040422af215..1dddf0e98b6 100644
--- a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM.h
+++ b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM.h
@@ -60,15 +60,15 @@ class RichardsComponentTransportLocalAssemblerInterface
 {
 public:
     virtual std::vector<double> const& getIntPtSaturation(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 };
 
@@ -112,8 +112,8 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
     Eigen::Map<const Eigen::RowVectorXd> getShapeMatrix(
@@ -121,9 +121,10 @@ public:
 
     std::vector<double> const& getIntPtSaturation(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
+
 private:
     unsigned const _element_id;
     RichardsComponentTransportProcessData const& _process_data;
diff --git a/ProcessLib/RichardsFlow/RichardsFlowFEM.h b/ProcessLib/RichardsFlow/RichardsFlowFEM.h
index ff890573a56..6f0e663cb89 100644
--- a/ProcessLib/RichardsFlow/RichardsFlowFEM.h
+++ b/ProcessLib/RichardsFlow/RichardsFlowFEM.h
@@ -59,16 +59,16 @@ class RichardsFlowLocalAssemblerInterface
 {
 public:
     virtual std::vector<double> const& getIntPtSaturation(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtDarcyVelocity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 };
 
 template <typename ShapeFunction, typename IntegrationMethod,
@@ -237,8 +237,8 @@ public:
 
     std::vector<double> const& getIntPtSaturation(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(!_saturation.empty());
@@ -247,8 +247,8 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override
     {
         auto const num_intpts = _shape_matrices.size();
diff --git a/ProcessLib/RichardsMechanics/LocalAssemblerInterface.h b/ProcessLib/RichardsMechanics/LocalAssemblerInterface.h
index 4a0e28ab05e..eaf661fd32f 100644
--- a/ProcessLib/RichardsMechanics/LocalAssemblerInterface.h
+++ b/ProcessLib/RichardsMechanics/LocalAssemblerInterface.h
@@ -21,27 +21,27 @@ struct LocalAssemblerInterface : public ProcessLib::LocalAssemblerInterface,
                                  public NumLib::ExtrapolatableElement
 {
     virtual std::vector<double> const& getIntPtSigma(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilon(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtDarcyVelocity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSaturation(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 };
 
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
index 3c7a31f6c25..d8a8a598abd 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
@@ -590,10 +590,11 @@ template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
 std::vector<double> const& RichardsMechanicsLocalAssembler<
     ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod,
     DisplacementDim>::
-    getIntPtSigma(const double /*t*/,
-                  GlobalVector const& /*current_solution*/,
-                  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-                  std::vector<double>& cache) const
+    getIntPtSigma(
+        const double /*t*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
+        std::vector<double>& cache) const
 {
     static const int kelvin_vector_size =
         MathLib::KelvinVector::KelvinVectorDimensions<DisplacementDim>::value;
@@ -619,10 +620,11 @@ template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
 std::vector<double> const& RichardsMechanicsLocalAssembler<
     ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod,
     DisplacementDim>::
-    getIntPtEpsilon(const double /*t*/,
-                    GlobalVector const& /*current_solution*/,
-                    NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-                    std::vector<double>& cache) const
+    getIntPtEpsilon(
+        const double /*t*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
+        std::vector<double>& cache) const
 {
     auto const kelvin_vector_size =
         MathLib::KelvinVector::KelvinVectorDimensions<DisplacementDim>::value;
@@ -646,12 +648,12 @@ template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
           typename IntegrationMethod, int DisplacementDim>
 std::vector<double> const& RichardsMechanicsLocalAssembler<
     ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod,
-    DisplacementDim>::getIntPtDarcyVelocity(const double t,
-                                            GlobalVector const&
-                                                current_solution,
-                                            NumLib::LocalToGlobalIndexMap const&
-                                                dof_table,
-                                            std::vector<double>& cache) const
+    DisplacementDim>::
+    getIntPtDarcyVelocity(
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const
 {
     auto const num_intpts = _ip_data.size();
 
@@ -709,10 +711,11 @@ template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
 std::vector<double> const& RichardsMechanicsLocalAssembler<
     ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod,
     DisplacementDim>::
-    getIntPtSaturation(const double /*t*/,
-                       GlobalVector const& /*current_solution*/,
-                       NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-                       std::vector<double>& cache) const
+    getIntPtSaturation(
+        const double /*t*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
+        std::vector<double>& cache) const
 {
     auto const num_intpts = _ip_data.size();
 
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h
index e5ade7e7472..45db4b36971 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h
@@ -152,26 +152,26 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
     std::vector<double> const& getIntPtSaturation(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const override;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const override;
 
     std::vector<double> const& getIntPtSigma(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
     std::vector<double> const& getIntPtEpsilon(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
 private:
diff --git a/ProcessLib/SmallDeformation/LocalAssemblerInterface.h b/ProcessLib/SmallDeformation/LocalAssemblerInterface.h
index 40fb04ad7e1..243180f244d 100644
--- a/ProcessLib/SmallDeformation/LocalAssemblerInterface.h
+++ b/ProcessLib/SmallDeformation/LocalAssemblerInterface.h
@@ -32,22 +32,22 @@ struct SmallDeformationLocalAssemblerInterface
         int const integration_order) = 0;
 
     virtual std::vector<double> const& getIntPtFreeEnergyDensity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> getSigma() const = 0;
     virtual std::vector<double> const& getIntPtSigma(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilon(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     // TODO move to NumLib::ExtrapolatableElement
diff --git a/ProcessLib/SmallDeformation/SmallDeformationFEM.h b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
index 040b2bcfdde..1d1356c165f 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationFEM.h
+++ b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
@@ -366,8 +366,8 @@ public:
 
     std::vector<double> const& getIntPtFreeEnergyDensity(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         cache.clear();
@@ -430,8 +430,8 @@ public:
 
     std::vector<double> const& getIntPtSigma(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         static const int kelvin_vector_size =
@@ -456,8 +456,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilon(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         auto const kelvin_vector_size =
diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp
index fbcb74cfc3b..6d58551a294 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp
+++ b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp
@@ -146,8 +146,9 @@ void SmallDeformationProcess<DisplacementDim>::initializeConcreteProcess(
             [fct, num_components](
                 LocalAssemblerInterface const& loc_asm,
                 const double /*t*/,
-                GlobalVector const& /*current_solution*/,
-                NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+                std::vector<GlobalVector*> const& /*x*/,
+                std::vector<
+                    NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
                 std::vector<double>& cache) -> std::vector<double> const& {
             const unsigned num_int_pts = loc_asm.getNumberOfIntegrationPoints();
 
diff --git a/ProcessLib/SmallDeformationNonlocal/LocalAssemblerInterface.h b/ProcessLib/SmallDeformationNonlocal/LocalAssemblerInterface.h
index 6f368a91c55..b906af5c082 100644
--- a/ProcessLib/SmallDeformationNonlocal/LocalAssemblerInterface.h
+++ b/ProcessLib/SmallDeformationNonlocal/LocalAssemblerInterface.h
@@ -40,39 +40,39 @@ struct SmallDeformationNonlocalLocalAssemblerInterface
         double& crack_volume) = 0;
 
     virtual std::vector<double> const& getIntPtEpsPV(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
     virtual std::vector<double> const& getIntPtEpsPDXX(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
     virtual std::vector<double> getKappaD() const = 0;
     virtual std::vector<double> const& getIntPtDamage(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtFreeEnergyDensity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> getSigma() const = 0;
     virtual std::vector<double> const& getIntPtSigma(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilon(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     // TODO move to NumLib::ExtrapolatableElement
diff --git a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h
index fcf5bd2a4af..950b03ad50a 100644
--- a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h
+++ b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h
@@ -607,8 +607,8 @@ public:
 
     std::vector<double> const& getIntPtFreeEnergyDensity(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         cache.clear();
@@ -624,8 +624,8 @@ public:
 
     std::vector<double> const& getIntPtEpsPV(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         cache.clear();
@@ -641,8 +641,8 @@ public:
 
     std::vector<double> const& getIntPtEpsPDXX(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         cache.clear();
@@ -658,8 +658,8 @@ public:
 
     std::vector<double> const& getIntPtSigma(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         static const int kelvin_vector_size =
@@ -684,8 +684,8 @@ public:
 
     std::vector<double> const& getIntPtEpsilon(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         auto const kelvin_vector_size =
@@ -792,8 +792,8 @@ public:
 
     std::vector<double> const& getIntPtDamage(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         cache.clear();
diff --git a/ProcessLib/TES/TESLocalAssembler-impl.h b/ProcessLib/TES/TESLocalAssembler-impl.h
index 25e75429a3a..e530a7d4fa7 100644
--- a/ProcessLib/TES/TESLocalAssembler-impl.h
+++ b/ProcessLib/TES/TESLocalAssembler-impl.h
@@ -199,10 +199,11 @@ template <typename ShapeFunction_, typename IntegrationMethod_,
           unsigned GlobalDim>
 std::vector<double> const&
 TESLocalAssembler<ShapeFunction_, IntegrationMethod_, GlobalDim>::
-    getIntPtSolidDensity(const double /*t*/,
-                         GlobalVector const& /*current_solution*/,
-                         NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-                         std::vector<double>& /*cache*/) const
+    getIntPtSolidDensity(
+        const double /*t*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
+        std::vector<double>& /*cache*/) const
 {
     return _d.getData().solid_density;
 }
@@ -211,10 +212,11 @@ template <typename ShapeFunction_, typename IntegrationMethod_,
           unsigned GlobalDim>
 std::vector<double> const&
 TESLocalAssembler<ShapeFunction_, IntegrationMethod_, GlobalDim>::
-    getIntPtLoading(const double /*t*/,
-                    GlobalVector const& /*current_solution*/,
-                    NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-                    std::vector<double>& cache) const
+    getIntPtLoading(
+        const double /*t*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
+        std::vector<double>& cache) const
 {
     auto const rho_SR = _d.getData().solid_density;
     auto const rho_SR_dry = _d.getAssemblyParameters().rho_SR_dry;
@@ -234,8 +236,8 @@ template <typename ShapeFunction_, typename IntegrationMethod_,
 std::vector<double> const&
 TESLocalAssembler<ShapeFunction_, IntegrationMethod_, GlobalDim>::
     getIntPtReactionDampingFactor(
-        const double /*t*/, GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double /*t*/, std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const
 {
     auto const fac = _d.getData().reaction_adaptor->getReactionDampingFactor();
@@ -251,10 +253,11 @@ template <typename ShapeFunction_, typename IntegrationMethod_,
           unsigned GlobalDim>
 std::vector<double> const&
 TESLocalAssembler<ShapeFunction_, IntegrationMethod_, GlobalDim>::
-    getIntPtReactionRate(const double /*t*/,
-                         GlobalVector const& /*current_solution*/,
-                         NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-                         std::vector<double>& /*cache*/) const
+    getIntPtReactionRate(
+        const double /*t*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
+        std::vector<double>& /*cache*/) const
 {
     return _d.getData().reaction_rate;
 }
@@ -263,10 +266,11 @@ template <typename ShapeFunction_, typename IntegrationMethod_,
           unsigned GlobalDim>
 std::vector<double> const&
 TESLocalAssembler<ShapeFunction_, IntegrationMethod_, GlobalDim>::
-    getIntPtDarcyVelocity(const double /*t*/,
-                          GlobalVector const& current_solution,
-                          NumLib::LocalToGlobalIndexMap const& dof_table,
-                          std::vector<double>& cache) const
+    getIntPtDarcyVelocity(
+        const double /*t*/,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const
 {
     auto const n_integration_points = _integration_method.getNumberOfPoints();
 
diff --git a/ProcessLib/TES/TESLocalAssembler.h b/ProcessLib/TES/TESLocalAssembler.h
index c34040057a7..7d92327b1c3 100644
--- a/ProcessLib/TES/TESLocalAssembler.h
+++ b/ProcessLib/TES/TESLocalAssembler.h
@@ -34,34 +34,34 @@ public:
                              std::vector<double> const& local_x_prev_ts) = 0;
 
     virtual std::vector<double> const& getIntPtSolidDensity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtLoading(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtReactionDampingFactor(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtReactionRate(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtDarcyVelocity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 };
 
 template <typename ShapeFunction_, typename IntegrationMethod_,
@@ -99,34 +99,34 @@ public:
                      std::vector<double> const& local_x_prev_ts) override;
 
     std::vector<double> const& getIntPtSolidDensity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const override;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const override;
 
     std::vector<double> const& getIntPtLoading(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
     std::vector<double> const& getIntPtReactionDampingFactor(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
     std::vector<double> const& getIntPtReactionRate(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const override;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const override;
 
     std::vector<double> const& getIntPtDarcyVelocity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const override;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const override;
 
 private:
     MeshLib::Element const& _element;
diff --git a/ProcessLib/TES/TESProcess.cpp b/ProcessLib/TES/TESProcess.cpp
index 8e29d55cd4c..9afb605ab4b 100644
--- a/ProcessLib/TES/TESProcess.cpp
+++ b/ProcessLib/TES/TESProcess.cpp
@@ -161,14 +161,15 @@ void TESProcess::initializeSecondaryVariables()
     };
 
     // creates an extrapolator
-    auto makeEx = [&](
-        unsigned const n_components,
-        std::vector<double> const& (TESLocalAssemblerInterface::*method)(
-            const double /*t*/,
-            GlobalVector const& /*current_solution*/,
-            NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-            std::vector<double>& /*cache*/)
-            const) -> SecondaryVariableFunctions {
+    auto makeEx =
+        [&](unsigned const n_components,
+            std::vector<double> const& (TESLocalAssemblerInterface::*method)(
+                const double /*t*/,
+                std::vector<GlobalVector*> const& /*x*/,
+                std::vector<
+                    NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
+                std::vector<double>& /*cache*/)
+                const) -> SecondaryVariableFunctions {
         return ProcessLib::makeExtrapolator(n_components, getExtrapolator(),
                                             _local_assemblers, method);
     };
@@ -317,8 +318,8 @@ NumLib::IterationResult TESProcess::postIterationConcreteProcess(
 
 GlobalVector const& TESProcess::computeVapourPartialPressure(
     const double /*t*/,
-    GlobalVector const& x,
-    NumLib::LocalToGlobalIndexMap const& dof_table,
+    std::vector<GlobalVector*> const& x,
+    std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
     std::unique_ptr<GlobalVector>& result_cache)
 {
     assert(&dof_table == _local_to_global_index_map.get());
@@ -349,8 +350,8 @@ GlobalVector const& TESProcess::computeVapourPartialPressure(
 
 GlobalVector const& TESProcess::computeRelativeHumidity(
     double const /*t*/,
-    GlobalVector const& x,
-    NumLib::LocalToGlobalIndexMap const& dof_table,
+    std::vector<GlobalVector*> const& xs,
+    std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
     std::unique_ptr<GlobalVector>& result_cache)
 {
     assert(&dof_table == _local_to_global_index_map.get());
@@ -386,8 +387,8 @@ GlobalVector const& TESProcess::computeRelativeHumidity(
 
 GlobalVector const& TESProcess::computeEquilibriumLoading(
     double const /*t*/,
-    GlobalVector const& x,
-    NumLib::LocalToGlobalIndexMap const& dof_table,
+    std::vector<GlobalVector*> const& xs,
+    std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
     std::unique_ptr<GlobalVector>& result_cache)
 {
     assert(&dof_table == _local_to_global_index_map.get());
diff --git a/ProcessLib/TES/TESProcess.h b/ProcessLib/TES/TESProcess.h
index bbc727937ce..4676ecac04a 100644
--- a/ProcessLib/TES/TESProcess.h
+++ b/ProcessLib/TES/TESProcess.h
@@ -73,20 +73,20 @@ private:
 
     GlobalVector const& computeVapourPartialPressure(
         const double t,
-        GlobalVector const& x,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::unique_ptr<GlobalVector>& result_cache);
 
     GlobalVector const& computeRelativeHumidity(
         const double t,
-        GlobalVector const& x,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::unique_ptr<GlobalVector>& result_cache);
 
     GlobalVector const& computeEquilibriumLoading(
         const double t,
-        GlobalVector const& x,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::unique_ptr<GlobalVector>& result_cache);
 
     std::vector<std::unique_ptr<TESLocalAssemblerInterface>> _local_assemblers;
diff --git a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPLocalAssembler.h b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPLocalAssembler.h
index 37eb7ab2a10..6a82a844458 100644
--- a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPLocalAssembler.h
+++ b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPLocalAssembler.h
@@ -62,16 +62,16 @@ class ThermalTwoPhaseFlowWithPPLocalAssemblerInterface
 {
 public:
     virtual std::vector<double> const& getIntPtSaturation(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtWettingPressure(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 };
 
 template <typename ShapeFunction, typename IntegrationMethod,
@@ -149,8 +149,8 @@ public:
 
     std::vector<double> const& getIntPtSaturation(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(!_saturation.empty());
@@ -159,8 +159,8 @@ public:
 
     std::vector<double> const& getIntPtWettingPressure(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(!_pressure_wetting.empty());
diff --git a/ProcessLib/ThermoHydroMechanics/LocalAssemblerInterface.h b/ProcessLib/ThermoHydroMechanics/LocalAssemblerInterface.h
index c528e5ed678..ef3681a5759 100644
--- a/ProcessLib/ThermoHydroMechanics/LocalAssemblerInterface.h
+++ b/ProcessLib/ThermoHydroMechanics/LocalAssemblerInterface.h
@@ -23,21 +23,21 @@ struct LocalAssemblerInterface : public ProcessLib::LocalAssemblerInterface,
     virtual std::vector<double> getSigma() const = 0;
 
     virtual std::vector<double> const& getIntPtSigma(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilon(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtDarcyVelocity(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 };
 
diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h
index c0c41764f0e..2bd21774998 100644
--- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h
+++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h
@@ -440,12 +440,12 @@ template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
           typename IntegrationMethod, int DisplacementDim>
 std::vector<double> const& ThermoHydroMechanicsLocalAssembler<
     ShapeFunctionDisplacement, ShapeFunctionPressure, IntegrationMethod,
-    DisplacementDim>::getIntPtDarcyVelocity(const double t,
-                                            GlobalVector const&
-                                                current_solution,
-                                            NumLib::LocalToGlobalIndexMap const&
-                                                dof_table,
-                                            std::vector<double>& cache) const
+    DisplacementDim>::
+    getIntPtDarcyVelocity(
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const
 {
     auto const num_intpts = _ip_data.size();
 
diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h
index 152913696a0..b8490fd83bf 100644
--- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h
+++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h
@@ -136,8 +136,8 @@ public:
 
     std::vector<double> const& getIntPtDarcyVelocity(
         const double t,
-        GlobalVector const& current_solution,
-        NumLib::LocalToGlobalIndexMap const& dof_table,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
 private:
@@ -192,8 +192,8 @@ private:
 
     std::vector<double> const& getIntPtSigma(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         static const int kelvin_vector_size =
@@ -219,8 +219,8 @@ private:
 
     virtual std::vector<double> const& getIntPtEpsilon(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         auto const kelvin_vector_size =
diff --git a/ProcessLib/ThermoMechanicalPhaseField/LocalAssemblerInterface.h b/ProcessLib/ThermoMechanicalPhaseField/LocalAssemblerInterface.h
index 6a768a8937c..82b551f145c 100644
--- a/ProcessLib/ThermoMechanicalPhaseField/LocalAssemblerInterface.h
+++ b/ProcessLib/ThermoMechanicalPhaseField/LocalAssemblerInterface.h
@@ -24,21 +24,21 @@ struct ThermoMechanicalPhaseFieldLocalAssemblerInterface
       public NumLib::ExtrapolatableElement
 {
     virtual std::vector<double> const& getIntPtSigma(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilon(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtHeatFlux(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 };
 
diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h
index e3dcf0e5bcf..0308107da7f 100644
--- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h
+++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h
@@ -275,8 +275,8 @@ public:
 private:
     std::vector<double> const& getIntPtSigma(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         static const int kelvin_vector_size =
@@ -301,8 +301,8 @@ private:
 
     std::vector<double> const& getIntPtEpsilon(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         auto const kelvin_vector_size =
@@ -327,8 +327,8 @@ private:
 
     std::vector<double> const& getIntPtHeatFlux(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         using KelvinVectorType = typename BMatricesType::KelvinVectorType;
diff --git a/ProcessLib/ThermoMechanics/LocalAssemblerInterface.h b/ProcessLib/ThermoMechanics/LocalAssemblerInterface.h
index bf1ea78dbf2..4ff9813373b 100644
--- a/ProcessLib/ThermoMechanics/LocalAssemblerInterface.h
+++ b/ProcessLib/ThermoMechanics/LocalAssemblerInterface.h
@@ -34,17 +34,16 @@ struct ThermoMechanicsLocalAssemblerInterface
     virtual std::vector<double> getEpsilonMechanical() const = 0;
 
     virtual std::vector<double> const& getIntPtSigma(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtEpsilon(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const = 0;
-
 };
 
 }  // namespace ThermoMechanics
diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h
index 7e51944ba75..91157581dd5 100644
--- a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h
+++ b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h
@@ -596,10 +596,11 @@ template <typename ShapeFunction, typename IntegrationMethod,
           int DisplacementDim>
 std::vector<double> const& ThermoMechanicsLocalAssembler<
     ShapeFunction, IntegrationMethod, DisplacementDim>::
-    getIntPtSigma(const double /*t*/,
-                  GlobalVector const& /*current_solution*/,
-                  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-                  std::vector<double>& cache) const
+    getIntPtSigma(
+        const double /*t*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
+        std::vector<double>& cache) const
 {
     static const int kelvin_vector_size =
         MathLib::KelvinVector::KelvinVectorDimensions<DisplacementDim>::value;
@@ -675,10 +676,11 @@ template <typename ShapeFunction, typename IntegrationMethod,
           int DisplacementDim>
 std::vector<double> const& ThermoMechanicsLocalAssembler<
     ShapeFunction, IntegrationMethod, DisplacementDim>::
-    getIntPtEpsilon(const double /*t*/,
-                    GlobalVector const& /*current_solution*/,
-                    NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-                    std::vector<double>& cache) const
+    getIntPtEpsilon(
+        const double /*t*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
+        std::vector<double>& cache) const
 {
     auto const kelvin_vector_size =
         MathLib::KelvinVector::KelvinVectorDimensions<DisplacementDim>::value;
diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h
index 028c1799672..144c113791a 100644
--- a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h
+++ b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h
@@ -276,9 +276,9 @@ private:
     std::vector<double> getSigma() const override;
 
     std::vector<double> const& getIntPtSigma(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
     std::size_t setEpsilon(double const* values);
@@ -286,9 +286,9 @@ private:
     std::vector<double> getEpsilon() const override;
 
     std::vector<double> const& getIntPtEpsilon(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
         std::vector<double>& cache) const override;
 
     std::size_t setEpsilonMechanical(double const* values);
diff --git a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPLocalAssembler.h b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPLocalAssembler.h
index 36c5397b053..8788660ea53 100644
--- a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPLocalAssembler.h
+++ b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPLocalAssembler.h
@@ -57,16 +57,16 @@ class TwoPhaseFlowWithPPLocalAssemblerInterface
 {
 public:
     virtual std::vector<double> const& getIntPtSaturation(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtWetPressure(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 };
 
 template <typename ShapeFunction, typename IntegrationMethod,
@@ -141,8 +141,8 @@ public:
 
     std::vector<double> const& getIntPtSaturation(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(!_saturation.empty());
@@ -151,8 +151,8 @@ public:
 
     std::vector<double> const& getIntPtWetPressure(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(!_pressure_wet.empty());
diff --git a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoLocalAssembler.h b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoLocalAssembler.h
index 0060c668218..c2fd4008d13 100644
--- a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoLocalAssembler.h
+++ b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoLocalAssembler.h
@@ -64,16 +64,16 @@ class TwoPhaseFlowWithPrhoLocalAssemblerInterface
 {
 public:
     virtual std::vector<double> const& getIntPtSaturation(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtNonWettingPressure(
-        const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
-        std::vector<double>& /*cache*/) const = 0;
+        const double t,
+        std::vector<GlobalVector*> const& x,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
+        std::vector<double>& cache) const = 0;
 };
 
 template <typename ShapeFunction, typename IntegrationMethod,
@@ -149,8 +149,8 @@ public:
 
     std::vector<double> const& getIntPtSaturation(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(!_saturation.empty());
@@ -159,8 +159,8 @@ public:
 
     std::vector<double> const& getIntPtNonWettingPressure(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         assert(!_pressure_nonwetting.empty());
diff --git a/Tests/NumLib/TestExtrapolation.cpp b/Tests/NumLib/TestExtrapolation.cpp
index de9ddc5eb7a..9bb98871d56 100644
--- a/Tests/NumLib/TestExtrapolation.cpp
+++ b/Tests/NumLib/TestExtrapolation.cpp
@@ -54,22 +54,23 @@ public:
 
     virtual std::vector<double> const& getStoredQuantity(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const = 0;
 
     virtual std::vector<double> const& getDerivedQuantity(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const = 0;
 };
 
 using IntegrationPointValuesMethod = std::vector<double> const& (
     LocalAssemblerDataInterface::*)(const double /*t*/,
-                                    GlobalVector const& /*current_solution*/,
-                                    NumLib::
-                                        LocalToGlobalIndexMap const& /*dof_table*/
+                                    std::vector<GlobalVector*> const& /*x*/,
+                                    std::vector<
+                                        NumLib::
+                                            LocalToGlobalIndexMap const*> const& /*dof_table*/
                                     ,
                                     std::vector<double>& /*cache*/) const;
 
@@ -105,8 +106,8 @@ public:
 
     std::vector<double> const& getStoredQuantity(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& /*cache*/) const override
     {
         return _int_pt_values;
@@ -114,8 +115,8 @@ public:
 
     std::vector<double> const& getDerivedQuantity(
         const double /*t*/,
-        GlobalVector const& /*current_solution*/,
-        NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+        std::vector<GlobalVector*> const& /*x*/,
+        std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
         std::vector<double>& cache) const override
     {
         cache.clear();
@@ -188,7 +189,7 @@ public:
 
     std::pair<GlobalVector const*, GlobalVector const*> extrapolate(
         IntegrationPointValuesMethod method, const double t,
-        const GlobalVector& x) const
+        std::vector<GlobalVector*> const& x) const
     {
         auto const extrapolatables =
             NumLib::makeExtrapolatable(_local_assemblers, method);
@@ -212,10 +213,10 @@ private:
     std::unique_ptr<ExtrapolatorInterface> _extrapolator;
 };
 
-void extrapolate(ExtrapolationTestProcess const& pcs,
-                 IntegrationPointValuesMethod method,
-                 GlobalVector const& expected_extrapolated_global_nodal_values,
-                 std::size_t const nnodes, std::size_t const nelements)
+void extrapolate(
+    ExtrapolationTestProcess const& pcs, IntegrationPointValuesMethod method,
+    std::vector<GlobalVector*> const& expected_extrapolated_global_nodal_values,
+    std::size_t const nnodes, std::size_t const nelements)
 {
     namespace LinAlg = MathLib::LinAlg;
 
-- 
GitLab