[LF] Added an optional input of reference temperature for LiquidFlow

Documentation/ProjectFile/prj/processes/process/LIQUID_FLOW/t_reference_temperature.md

+An optional input of a reference temperature for liquid flow process when the
+ flow process is not coupled to heat transport process, and a value other than
+ the default one (room temperature of 273.15+18 K) is selected.

ProcessLib/LiquidFlow/CreateLiquidFlowProcess.cpp

@@ -15,6 +15,7 @@
 
 #include "MeshLib/MeshGenerators/MeshGenerator.h"
 
+#include "MaterialLib/PhysicalConstant.h"
 #include "ProcessLib/Utils/ParseSecondaryVariables.h"
 #include "ProcessLib/Utils/ProcessUtils.h"
 
@@ -67,6 +68,13 @@ std::unique_ptr<Process> createLiquidFlowProcess(
     assert(g >= 0.);
     const int gravity_axis_id = (g == 0.) ? -1 : gravity_axis_id_input;
 
+    auto const& refT =
+        //! \ogs_file_param{prj__processes__process__LIQUID_FLOW__reference_temperature}
+        config.getConfigParameterOptional<double>("reference_temperature");
+    const double reference_temperature =
+        refT ? *refT
+             : MaterialLib::PhysicalConstant::CelsiusZeroInKelvin + 18.0;
+
     //! \ogs_file_param{prj__processes__process__LIQUID_FLOW__material_property}
     auto const& mat_config = config.getConfigSubtree("material_property");
 
@@ -80,7 +88,7 @@ std::unique_ptr<Process> createLiquidFlowProcess(
             mesh, std::move(jacobian_assembler), parameters, integration_order,
             std::move(process_variables), std::move(secondary_variables),
             std::move(named_function_caller), *mat_ids, has_material_ids,
-            gravity_axis_id, g, mat_config}};
+            gravity_axis_id, g, reference_temperature, mat_config}};
     }
     else
     {
@@ -101,7 +109,8 @@ std::unique_ptr<Process> createLiquidFlowProcess(
             mesh, std::move(jacobian_assembler), parameters, integration_order,
             std::move(process_variables), std::move(secondary_variables),
             std::move(named_function_caller), *dummy_property_vector,
-            has_material_ids, gravity_axis_id, g, mat_config}};
+            has_material_ids, gravity_axis_id, g, reference_temperature,
+            mat_config}};
     }
 }
 

ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h

@@ -14,7 +14,6 @@
 
 #include "LiquidFlowLocalAssembler.h"
 
-#include "MaterialLib/PhysicalConstant.h"
 #include "NumLib/Function/Interpolation.h"
 
 #include "ProcessLib/HeatConduction/HeatConductionProcess.h"
@@ -133,10 +132,6 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
     //       the integration loop for non-constant porosity and storage models.
     double porosity_variable = 0.;
     double storage_variable = 0.;
-    // Reference temperature for fluid property models because of no coupled
-    // heat transport process. TODO: Add an optional input for this.
-    const double temperature =
-        MaterialLib::PhysicalConstant::CelsiusZeroInKelvin + 18.0;
     for (unsigned ip = 0; ip < n_integration_points; ip++)
     {
         auto const& sm = _shape_matrices[ip];
@@ -151,15 +146,16 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
         // Assemble mass matrix, M
         local_M.noalias() += _material_properties.getMassCoefficient(
                                  material_id, t, pos, porosity_variable,
-                                 storage_variable, p, temperature) *
+                                 storage_variable, p, _reference_temperature) *
                              sm.N.transpose() * sm.N * integration_factor;
 
         // Compute density:
         const double rho_g =
-            _material_properties.getLiquidDensity(p, temperature) *
+            _material_properties.getLiquidDensity(p, _reference_temperature) *
             _gravitational_acceleration;
         // Compute viscosity:
-        const double mu = _material_properties.getViscosity(p, temperature);
+        const double mu =
+            _material_properties.getViscosity(p, _reference_temperature);
 
         // Assemble Laplacian, K, and RHS by the gravitational term
         LaplacianGravityVelocityCalculator::calculateLaplacianAndGravityTerm(
@@ -286,8 +282,6 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
     unsigned const n_integration_points =
         _integration_method.getNumberOfPoints();
 
-    const double temperature =
-        MaterialLib::PhysicalConstant::CelsiusZeroInKelvin + 18.0;
     for (unsigned ip = 0; ip < n_integration_points; ip++)
     {
         auto const& sm = _shape_matrices[ip];
@@ -295,10 +289,11 @@ void LiquidFlowLocalAssembler<ShapeFunction, IntegrationMethod, GlobalDim>::
         NumLib::shapeFunctionInterpolate(local_x, sm.N, p);
 
         const double rho_g =
-            _material_properties.getLiquidDensity(p, temperature) *
+            _material_properties.getLiquidDensity(p, _reference_temperature) *
             _gravitational_acceleration;
         // Compute viscosity:
-        const double mu = _material_properties.getViscosity(p, temperature);
+        const double mu =
+            _material_properties.getViscosity(p, _reference_temperature);
 
         LaplacianGravityVelocityCalculator::calculateVelocity(
             _darcy_velocities, local_p_vec, sm, permeability, ip, mu, rho_g,

ProcessLib/LiquidFlow/LiquidFlowLocalAssembler.h

@@ -71,6 +71,7 @@ public:
         unsigned const integration_order,
         int const gravitational_axis_id,
         double const gravitational_acceleration,
+        double const reference_temperature,
         LiquidFlowMaterialProperties const& material_propertries)
         : _element(element),
           _integration_method(integration_order),
@@ -79,6 +80,7 @@ public:
               element, is_axially_symmetric, _integration_method)),
           _gravitational_axis_id(gravitational_axis_id),
           _gravitational_acceleration(gravitational_acceleration),
+          _reference_temperature(reference_temperature),
           _material_properties(material_propertries)
     {
     }
@@ -219,6 +221,7 @@ private:
 
     const int _gravitational_axis_id;
     const double _gravitational_acceleration;
+    const double _reference_temperature;
     const LiquidFlowMaterialProperties& _material_properties;
 };
 

ProcessLib/LiquidFlow/LiquidFlowProcess.cpp

@@ -41,12 +41,14 @@ LiquidFlowProcess::LiquidFlowProcess(
     bool const has_material_ids,
     int const gravitational_axis_id,
     double const gravitational_acceleration,
+    double const reference_temperature,
     BaseLib::ConfigTree const& config)
     : Process(mesh, std::move(jacobian_assembler), parameters,
               integration_order, std::move(process_variables),
               std::move(secondary_variables), std::move(named_function_caller)),
       _gravitational_axis_id(gravitational_axis_id),
       _gravitational_acceleration(gravitational_acceleration),
+      _reference_temperature(reference_temperature),
       _material_properties(createLiquidFlowMaterialProperties(
           config, parameters, has_material_ids, material_ids))
 {
@@ -63,7 +65,8 @@ void LiquidFlowProcess::initializeConcreteProcess(
         mesh.getDimension(), mesh.getElements(), dof_table,
         pv.getShapeFunctionOrder(), _local_assemblers,
         mesh.isAxiallySymmetric(), integration_order, _gravitational_axis_id,
-        _gravitational_acceleration, *_material_properties);
+        _gravitational_acceleration, _reference_temperature,
+        *_material_properties);
 
     _secondary_variables.addSecondaryVariable(
         "darcy_velocity_x", 1,

ProcessLib/LiquidFlow/LiquidFlowProcess.h

@@ -70,6 +70,7 @@ public:
         bool const has_material_ids,
         int const gravitational_axis_id,
         double const gravitational_acceleration,
+        double const reference_temperature,
         BaseLib::ConfigTree const& config);
 
     void computeSecondaryVariableConcrete(
@@ -106,6 +107,7 @@ private:
 
     const int _gravitational_axis_id;
     const double _gravitational_acceleration;
+    const double _reference_temperature;
     const std::unique_ptr<LiquidFlowMaterialProperties> _material_properties;
 
     std::vector<std::unique_ptr<LiquidFlowLocalAssemblerInterface>>