consider the liquid phase and change the calculation of the volumtric heat capacity accordingly

ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h

@@ -14,6 +14,7 @@
 #include <vector>
 
 #include "MaterialLib/MPL/Medium.h"
+#include "MaterialLib/MPL/Utils/FormEffectiveThermalConductivity.h"
 #include "MaterialLib/MPL/Utils/FormEigenTensor.h"
 #include "MathLib/LinAlg/Eigen/EigenMapTools.h"
 #include "NumLib/Fem/ShapeMatrixPolicy.h"
@@ -113,37 +114,45 @@ void HeatTransportBHELocalAssemblerSoil<ShapeFunction, IntegrationMethod>::
                     MaterialPropertyLib::PropertyType::thermal_conductivity)
                 .template value<double>(vars, pos, t);
 
+        auto const density_s =
+            solid_phase.property(MaterialPropertyLib::PropertyType::density)
+                .template value<double>(vars, pos, t);
+
         auto const heat_capacity_s =
             solid_phase
                 .property(
                     MaterialPropertyLib::PropertyType::specific_heat_capacity)
                 .template value<double>(vars, pos, t);
 
-        auto const heat_capacity_f =
+        auto const k_f =
             liquid_phase
                 .property(
-                    MaterialPropertyLib::PropertyType::specific_heat_capacity)
-                .template value<double>(vars, pos, t);
-
-        auto const density_s =
-            solid_phase.property(MaterialPropertyLib::PropertyType::density)
+                    MaterialPropertyLib::PropertyType::thermal_conductivity)
                 .template value<double>(vars, pos, t);
 
         auto const density_f =
             liquid_phase.property(MaterialPropertyLib::PropertyType::density)
                 .template value<double>(vars, pos, t);
-        // for now only using the solid phase parameters
 
-        /*auto const velocity =
+        auto const heat_capacity_f =
             liquid_phase
-                .property(MaterialPropertyLib::PropertyType::phase_velocity)
-                .template value<GlobalDimVectorType>(vars, pos, t);*/
+                .property(
+                    MaterialPropertyLib::PropertyType::specific_heat_capacity)
+                .template value<double>(vars, pos, t);
+
+        auto const porosity =
+            medium.property(MaterialPropertyLib::PropertyType::porosity)
+                .template value<double>(vars, pos, t);
+
+        // for now only using the solid and liquid phase parameters
+
         auto const velocity_arr =
             liquid_phase
                 .property(MaterialPropertyLib::PropertyType::phase_velocity)
                 .value(vars, pos, t);
         auto velocity = Eigen::Map<Eigen::Matrix<double, 1, 3> const>(
             std::get<MaterialPropertyLib::Vector>(velocity_arr).data(), 1, 3);
+
         // assemble Conductance matrix
         local_K.noalias() +=
             (dNdx.transpose() * dNdx * k_s +
@@ -151,8 +160,9 @@ void HeatTransportBHELocalAssemblerSoil<ShapeFunction, IntegrationMethod>::
             w;
 
         // assemble Mass matrix
-        local_M.noalias() +=
-            N.transpose() * N * w * density_s * heat_capacity_s;
+        local_M.noalias() += N.transpose() * N * w *
+                             (density_s * heat_capacity_s * (1 - porosity) +
+                              density_f * heat_capacity_f * porosity);
     }
 
     // debugging