From abe3174b498393d025995a3a3f1baf80391ade49 Mon Sep 17 00:00:00 2001
From: Wenqing Wang <wenqing.wang@ufz.de>
Date: Mon, 10 Oct 2016 16:02:44 +0200
Subject: [PATCH] [PCS] Corrected the computation of velocity for
 HeatConduction

---
 ProcessLib/HeatConduction/HeatConductionFEM.h | 29 +++++++++++++++++--
 .../HeatConduction/HeatConductionProcess.cpp  |  9 ++++++
 Tests/Data                                    |  2 +-
 3 files changed, 36 insertions(+), 4 deletions(-)

diff --git a/ProcessLib/HeatConduction/HeatConductionFEM.h b/ProcessLib/HeatConduction/HeatConductionFEM.h
index c9de972da31..e2fad4f8534 100644
--- a/ProcessLib/HeatConduction/HeatConductionFEM.h
+++ b/ProcessLib/HeatConduction/HeatConductionFEM.h
@@ -116,10 +116,33 @@ public:
             local_M.noalias() += sm.N.transpose() * density * heat_capacity *
                                  sm.N * sm.detJ * wp.getWeight() *
                                  sm.integralMeasure;
+        }
+    }
+
+    void computeSecondaryVariableConcrete(
+                                    const double t,
+                                    std::vector<double> const& local_x) override
+    {
+        auto const local_matrix_size = local_x.size();
+        // This assertion is valid only if all nodal d.o.f. use the same shape
+        // matrices.
+        assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
+
+        unsigned const n_integration_points =
+            _integration_method.getNumberOfPoints();
+
+        SpatialPosition pos;
+        pos.setElementID(_element.getID());
+        const auto local_x_vec =
+        MathLib::toVector<NodalVectorType>(local_x, local_matrix_size);
+
+        for (unsigned ip = 0; ip < n_integration_points; ip++)
+        {
+            pos.setIntegrationPoint(ip);
+            auto const& sm = _shape_matrices[ip];
+            auto const k = _process_data.thermal_conductivity(t, pos)[0];
             // heat flux only computed for output.
-            GlobalDimVectorType const heat_flux =
-                -k * sm.dNdx * Eigen::Map<const NodalVectorType>(
-                                   local_x.data(), ShapeFunction::NPOINTS);
+            GlobalDimVectorType const heat_flux = -k * sm.dNdx * local_x_vec;
 
             for (unsigned d = 0; d < GlobalDim; ++d)
             {
diff --git a/ProcessLib/HeatConduction/HeatConductionProcess.cpp b/ProcessLib/HeatConduction/HeatConductionProcess.cpp
index 617b90ff0ca..04cbe6d9471 100644
--- a/ProcessLib/HeatConduction/HeatConductionProcess.cpp
+++ b/ProcessLib/HeatConduction/HeatConductionProcess.cpp
@@ -96,5 +96,14 @@ void HeatConductionProcess::assembleWithJacobianConcreteProcess(
         dx_dx, M, K, b, Jac);
 }
 
+void HeatConductionProcess::computeSecondaryVariableConcrete(const double t,
+                                                         GlobalVector const& x)
+{
+    DBUG("Compute the velocity for LiquidFlowProcess.");
+    GlobalExecutor::executeMemberOnDereferenced(
+            &HeatConductionLocalAssemblerInterface::computeSecondaryVariable,
+            _local_assemblers, *_local_to_global_index_map, t, x);
+}
+
 }  // namespace HeatConduction
 }  // namespace ProcessLib
diff --git a/Tests/Data b/Tests/Data
index 782b7f47fbf..8806cf31aa3 160000
--- a/Tests/Data
+++ b/Tests/Data
@@ -1 +1 @@
-Subproject commit 782b7f47fbf89f7ed72dbe034f372966a680c953
+Subproject commit 8806cf31aa300ee9c631c73674e218e9a60297cc
-- 
GitLab