From 1f9b6325243e4d3f09be5fceb538e69c93c61eba Mon Sep 17 00:00:00 2001
From: renchao_lu <renchao.lu@gmail.com>
Date: Tue, 22 Dec 2020 16:52:24 +0100
Subject: [PATCH] [PL] Move into solveReactionEquation.

---
 .../ComponentTransportProcess.cpp             | 25 +++++++++++++++++++
 .../ComponentTransportProcess.h               |  3 +++
 ProcessLib/Process.h                          |  5 ++++
 ProcessLib/TimeLoop.cpp                       | 19 +-------------
 4 files changed, 34 insertions(+), 18 deletions(-)

diff --git a/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp b/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp
index 6736b76e350..d29078e8d83 100644
--- a/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp
+++ b/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp
@@ -197,6 +197,31 @@ void ComponentTransportProcess::
         _local_assemblers, pv.getActiveElementIDs(), _coupled_solutions);
 }
 
+void ComponentTransportProcess::solveReactionEquation(
+    std::vector<GlobalVector*>& x, double const t, double const dt)
+{
+    if (!_chemical_solver_interface)
+    {
+        return;
+    }
+
+    // Sequential non-iterative approach applied here to perform water
+    // chemistry calculation followed by resolving component transport
+    // process.
+    // TODO: move into a global loop to consider both mass balance over
+    // space and localized chemical equilibrium between solutes.
+    BaseLib::RunTime time_phreeqc;
+    time_phreeqc.start();
+
+    _chemical_solver_interface->doWaterChemistryCalculation(
+        interpolateNodalValuesToIntegrationPoints(x), dt);
+
+    extrapolateIntegrationPointValuesToNodes(
+        t, _chemical_solver_interface->getIntPtProcessSolutions(), x);
+
+    INFO("[time] Phreeqc took {:g} s.", time_phreeqc.elapsed());
+}
+
 std::vector<GlobalVector>
 ComponentTransportProcess::interpolateNodalValuesToIntegrationPoints(
     std::vector<GlobalVector*> const& nodal_values_vectors) const
diff --git a/ProcessLib/ComponentTransport/ComponentTransportProcess.h b/ProcessLib/ComponentTransport/ComponentTransportProcess.h
index 873eed1e112..3f1a89caeeb 100644
--- a/ProcessLib/ComponentTransport/ComponentTransportProcess.h
+++ b/ProcessLib/ComponentTransport/ComponentTransportProcess.h
@@ -124,6 +124,9 @@ public:
     void setCoupledTermForTheStaggeredSchemeToLocalAssemblers(
         int const process_id) override;
 
+    void solveReactionEquation(std::vector<GlobalVector*>& x, double const t,
+                               double const dt) override;
+
     std::vector<GlobalVector> interpolateNodalValuesToIntegrationPoints(
         std::vector<GlobalVector*> const& nodal_values_vectors) const override;
 
diff --git a/ProcessLib/Process.h b/ProcessLib/Process.h
index c01cf00daa5..ce4db97ebb0 100644
--- a/ProcessLib/Process.h
+++ b/ProcessLib/Process.h
@@ -179,6 +179,11 @@ public:
         return Eigen::Vector3d{};
     }
 
+    virtual void solveReactionEquation(std::vector<GlobalVector*>& /*x*/,
+                                       double const /*t*/, double const /*dt*/)
+    {
+    }
+
 protected:
     NumLib::Extrapolator& getExtrapolator() const
     {
diff --git a/ProcessLib/TimeLoop.cpp b/ProcessLib/TimeLoop.cpp
index 3e232971801..b8114d76bcc 100644
--- a/ProcessLib/TimeLoop.cpp
+++ b/ProcessLib/TimeLoop.cpp
@@ -823,26 +823,9 @@ TimeLoop::solveCoupledEquationSystemsByStaggeredScheme(
             timestep_id, t);
     }
 
-    if (_chemical_solver_interface)
     {
-        // Sequential non-iterative approach applied here to perform water
-        // chemistry calculation followed by resolving component transport
-        // process.
-        // TODO: move into a global loop to consider both mass balance over
-        // space and localized chemical equilibrium between solutes.
-        BaseLib::RunTime time_phreeqc;
-        time_phreeqc.start();
-
         auto& pcs = _per_process_data[0]->process;
-        _chemical_solver_interface->doWaterChemistryCalculation(
-            pcs.interpolateNodalValuesToIntegrationPoints(_process_solutions),
-            dt);
-
-        pcs.extrapolateIntegrationPointValuesToNodes(
-            t, _chemical_solver_interface->getIntPtProcessSolutions(),
-            _process_solutions);
-
-        INFO("[time] Phreeqc took {:g} s.", time_phreeqc.elapsed());
+        pcs.solveReactionEquation(_process_solutions, t, dt);
     }
 
     return nonlinear_solver_status;
-- 
GitLab