From 9f72c2eeb9b144d150a6ffbbba1da28ddd5be302 Mon Sep 17 00:00:00 2001
From: Wenqing Wang <wenqing.wang@ufz.de>
Date: Thu, 10 Oct 2024 13:59:43 +0200
Subject: [PATCH] [LF] Set phase_variable in LiquidFlowData

---
 .../LiquidFlow/CreateLiquidFlowProcess.cpp    | 10 +++++++++-
 ProcessLib/LiquidFlow/LiquidFlowData.h        |  3 +++
 .../LiquidFlowLocalAssembler-impl.h           | 20 ++++++++-----------
 3 files changed, 20 insertions(+), 13 deletions(-)

diff --git a/ProcessLib/LiquidFlow/CreateLiquidFlowProcess.cpp b/ProcessLib/LiquidFlow/CreateLiquidFlowProcess.cpp
index 42ada7a720c..31042261674 100644
--- a/ProcessLib/LiquidFlow/CreateLiquidFlowProcess.cpp
+++ b/ProcessLib/LiquidFlow/CreateLiquidFlowProcess.cpp
@@ -188,6 +188,13 @@ std::unique_ptr<Process> createLiquidFlowProcess(
             *aperture_config, "parameter", parameters, 1);
     }
 
+    // The uniqueness of phase has already been checked in
+    // `checkMPLProperties`.
+    MaterialPropertyLib::Variable const phase_variable =
+        (*ranges::begin(media_map.media()))->hasPhase("Gas")
+            ? MaterialPropertyLib::Variable::gas_phase_pressure
+            : MaterialPropertyLib::Variable::liquid_phase_pressure;
+
     LiquidFlowData process_data{
         covertEquationBalanceTypeFromString(equation_balance_type_str),
         std::move(media_map),
@@ -197,7 +204,8 @@ std::unique_ptr<Process> createLiquidFlowProcess(
         std::move(specific_body_force),
         has_gravity,
         *aperture_size_parameter,
-        NumLib::ShapeMatrixCache{integration_order}};
+        NumLib::ShapeMatrixCache{integration_order},
+        phase_variable};
 
     return std::make_unique<LiquidFlowProcess>(
         std::move(name), mesh, std::move(jacobian_assembler), parameters,
diff --git a/ProcessLib/LiquidFlow/LiquidFlowData.h b/ProcessLib/LiquidFlow/LiquidFlowData.h
index c0f5b8640c0..9d63e2faef1 100644
--- a/ProcessLib/LiquidFlow/LiquidFlowData.h
+++ b/ProcessLib/LiquidFlow/LiquidFlowData.h
@@ -13,6 +13,7 @@
 #include <Eigen/Core>
 
 #include "MaterialLib/MPL/MaterialSpatialDistributionMap.h"
+#include "MaterialLib/MPL/VariableType.h"
 #include "NumLib/Fem/ShapeMatrixCache.h"
 #include "ParameterLib/Parameter.h"
 
@@ -54,6 +55,8 @@ struct LiquidFlowData final
 
     /// caches for each mesh element type the shape matrix
     NumLib::ShapeMatrixCache shape_matrix_cache;
+
+    MaterialPropertyLib::Variable const phase_variable;
 };
 
 }  // namespace LiquidFlow
diff --git a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h
index 6774e9ba19c..0bbf1b0aeca 100644
--- a/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h
+++ b/ProcessLib/LiquidFlow/LiquidFlowLocalAssembler-impl.h
@@ -134,14 +134,10 @@ void LiquidFlowLocalAssembler<ShapeFunction, GlobalDim>::
 
     auto const& medium = *_process_data.media_map.getMedium(_element.getID());
     auto const& fluid_phase = fluidPhase(medium);
-    bool is_gas_phase = medium.hasPhase("Gas");
-    auto const fluid_pressure_variable =
-        is_gas_phase ? MaterialPropertyLib::Variable::gas_phase_pressure
-                     : MaterialPropertyLib::Variable::liquid_phase_pressure;
 
     MaterialPropertyLib::VariableArray vars;
-    auto& phase_pressue =
-        is_gas_phase ? vars.gas_phase_pressure : vars.liquid_phase_pressure;
+    auto& phase_pressure = medium.hasPhase("Gas") ? vars.gas_phase_pressure
+                                                  : vars.liquid_phase_pressure;
 
     vars.temperature =
         medium[MaterialPropertyLib::PropertyType::reference_temperature]
@@ -160,7 +156,7 @@ void LiquidFlowLocalAssembler<ShapeFunction, GlobalDim>::
         auto const& ip_data = _ip_data[ip];
         auto const& N = Ns[ip];
 
-        phase_pressue = N.dot(local_p_vec);
+        phase_pressure = N.dot(local_p_vec);
 
         auto const [fluid_density, viscosity] =
             getFluidDensityAndViscosity(t, dt, pos, fluid_phase, vars);
@@ -175,8 +171,8 @@ void LiquidFlowLocalAssembler<ShapeFunction, GlobalDim>::
         auto const get_drho_dp = [&]()
         {
             return fluid_phase[MaterialPropertyLib::PropertyType::density]
-                .template dValue<double>(vars, fluid_pressure_variable, pos, t,
-                                         dt);
+                .template dValue<double>(vars, _process_data.phase_variable,
+                                         pos, t, dt);
         };
         auto const storage =
             _process_data.equation_balance_type == EquationBalanceType::volume
@@ -295,8 +291,8 @@ void LiquidFlowLocalAssembler<ShapeFunction, GlobalDim>::
     auto const& fluid_phase = fluidPhase(medium);
 
     MaterialPropertyLib::VariableArray vars;
-    auto& phase_pressue = medium.hasPhase("Gas") ? vars.gas_phase_pressure
-                                                 : vars.liquid_phase_pressure;
+    auto& phase_pressure = medium.hasPhase("Gas") ? vars.gas_phase_pressure
+                                                  : vars.liquid_phase_pressure;
 
     vars.temperature =
         medium[MaterialPropertyLib::PropertyType::reference_temperature]
@@ -315,7 +311,7 @@ void LiquidFlowLocalAssembler<ShapeFunction, GlobalDim>::
         auto const& ip_data = _ip_data[ip];
         auto const& N = Ns[ip];
 
-        phase_pressue = N.dot(local_p_vec);
+        phase_pressure = N.dot(local_p_vec);
 
         auto const [fluid_density, viscosity] =
             getFluidDensityAndViscosity(t, dt, pos, fluid_phase, vars);