diff --git a/MaterialLib/MPL/Properties/IdealGasLaw.cpp b/MaterialLib/MPL/Properties/IdealGasLaw.cpp
index d362a23a1a63ccc93b756e51fb6f2f75e9fdab20..5ce307d36a2efd60c07a7c60aa280b36ce1e1c57 100644
--- a/MaterialLib/MPL/Properties/IdealGasLaw.cpp
+++ b/MaterialLib/MPL/Properties/IdealGasLaw.cpp
@@ -20,60 +20,53 @@
namespace MaterialPropertyLib
{
-double molarMass(std::variant<Medium*, Phase*, Component*> const scale,
- VariableArray const& variable_array,
- ParameterLib::SpatialPosition const& pos, double const t,
- double const dt)
-{
- return std::visit(
- [&variable_array, &pos, t, dt](auto&& s) -> double {
- return s->property(PropertyType::molar_mass)
- .template value<double>(variable_array, pos, t, dt);
- },
- scale);
-}
-
IdealGasLaw::IdealGasLaw(std::string name)
{
name_ = std::move(name);
}
-PropertyDataType IdealGasLaw::value(VariableArray const& variable_array,
- ParameterLib::SpatialPosition const& pos,
- double const t, double const dt) const
+PropertyDataType IdealGasLaw::value(
+ VariableArray const& variable_array,
+ ParameterLib::SpatialPosition const& /*pos*/, double const /*t*/,
+ double const /*dt*/) const
{
const double gas_constant = MaterialLib::PhysicalConstant::IdealGasConstant;
const double pressure = std::get<double>(
variable_array[static_cast<int>(Variable::phase_pressure)]);
const double temperature = std::get<double>(
variable_array[static_cast<int>(Variable::temperature)]);
- double molar_mass = molarMass(scale_, variable_array, pos, t, dt);
+ const double molar_mass = std::get<double>(
+ variable_array[static_cast<int>(Variable::molar_mass)]);
const double density = pressure * molar_mass / gas_constant / temperature;
return density;
}
-PropertyDataType IdealGasLaw::dValue(VariableArray const& variable_array,
- Variable const primary_variable,
- ParameterLib::SpatialPosition const& pos,
- double const t, double const dt) const
+PropertyDataType IdealGasLaw::dValue(
+ VariableArray const& variable_array, Variable const primary_variable,
+ ParameterLib::SpatialPosition const& /*pos*/, double const /*t*/,
+ double const /*dt*/) const
{
const double gas_constant = MaterialLib::PhysicalConstant::IdealGasConstant;
const double pressure = std::get<double>(
variable_array[static_cast<int>(Variable::phase_pressure)]);
const double temperature = std::get<double>(
variable_array[static_cast<int>(Variable::temperature)]);
- double molar_mass = molarMass(scale_, variable_array, pos, t, dt);
+ const double molar_mass = std::get<double>(
+ variable_array[static_cast<int>(Variable::molar_mass)]);
+ // todo: add molar mass derivatives
if (primary_variable == Variable::temperature)
{
+ // extend to take temperature-dependent molar mass into account
return -pressure * molar_mass / gas_constant / temperature /
temperature;
}
if (primary_variable == Variable::phase_pressure)
{
+ // extend to take pressure-dependent molar mass into account
return molar_mass / gas_constant / temperature;
}
@@ -84,29 +77,33 @@ PropertyDataType IdealGasLaw::dValue(VariableArray const& variable_array,
return 0.;
}
-PropertyDataType IdealGasLaw::d2Value(VariableArray const& variable_array,
- Variable const primary_variable1,
- Variable const primary_variable2,
- ParameterLib::SpatialPosition const& pos,
- double const t, double const dt) const
+PropertyDataType IdealGasLaw::d2Value(
+ VariableArray const& variable_array, Variable const primary_variable1,
+ Variable const primary_variable2,
+ ParameterLib::SpatialPosition const& /*pos*/, double const /*t*/,
+ double const /*dt*/) const
{
const double gas_constant = MaterialLib::PhysicalConstant::IdealGasConstant;
const double pressure = std::get<double>(
variable_array[static_cast<int>(Variable::phase_pressure)]);
const double temperature = std::get<double>(
variable_array[static_cast<int>(Variable::temperature)]);
- double molar_mass = molarMass(scale_, variable_array, pos, t, dt);
+ const double molar_mass = std::get<double>(
+ variable_array[static_cast<int>(Variable::molar_mass)]);
+ // todo: add molar mass derivatives
if ((primary_variable1 == Variable::phase_pressure) &&
(primary_variable2 == Variable::phase_pressure))
{
// d2rho_dp2
+ // extend to take pressure-dependent molar mass into account
return 0.;
}
if ((primary_variable1 == Variable::temperature) &&
(primary_variable2 == Variable::temperature))
{
// d2rho_dT2
+ // extend to take temperature-dependent molar mass into account
return 2. * molar_mass * pressure / gas_constant / temperature /
temperature / temperature;
}
@@ -116,6 +113,7 @@ PropertyDataType IdealGasLaw::d2Value(VariableArray const& variable_array,
(primary_variable2 == Variable::phase_pressure)))
{
// d2rho_dpdT or d2rho_dTdp
+ // extend to take pressure-temperature-dependent molar mass into account
return -molar_mass / gas_constant / temperature / temperature;
}
diff --git a/MaterialLib/MPL/Properties/IdealGasLaw.h b/MaterialLib/MPL/Properties/IdealGasLaw.h
index 17739e1c2a3fb0db8de41b6159c328ca86ae6bc4..a981b38010399655829db43073c161ee11cfac1e 100644
--- a/MaterialLib/MPL/Properties/IdealGasLaw.h
+++ b/MaterialLib/MPL/Properties/IdealGasLaw.h
@@ -57,8 +57,9 @@ public:
double const /*dt*/) const override;
PropertyDataType d2Value(VariableArray const& variable_array,
Variable const variable1, Variable const variable2,
- ParameterLib::SpatialPosition const& pos,
- double const t, double const dt) const override;
+ ParameterLib::SpatialPosition const& /*pos*/,
+ double const /*t*/,
+ double const /*dt*/) const override;
};
} // namespace MaterialPropertyLib
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
index d264fefd36df11b4dc8d0b42f9851e5395231e77..052e58df7728b2acb8c7f99bf3f13bbab7532dd0 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
@@ -240,6 +240,18 @@ void HydroMechanicsLocalAssembler<ShapeFunctionDisplacement,
auto const mu = fluid.property(MPL::PropertyType::viscosity)
.template value<double>(vars, x_position, t, dt);
+
+ // Quick workaround: If fluid density is described as ideal gas, then
+ // the molar mass must be passed to the MPL::IdealGasLaw via the
+ // varialbe_array and the fluid must have the property
+ // MPL::PropertyType::molar_mass. For other density models (e.g.
+ // Constant), it is not mandatory to specify the molar mass.
+ if (fluid.hasProperty(MPL::PropertyType::molar_mass))
+ {
+ vars[static_cast<int>(MPL::Variable::molar_mass)] =
+ fluid.property(MPL::PropertyType::molar_mass)
+ .template value<double>(vars, x_position, t, dt);
+ }
auto const rho_fr =
fluid.property(MPL::PropertyType::density)
.template value<double>(vars, x_position, t, dt);
@@ -435,6 +447,17 @@ HydroMechanicsLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
auto const mu = fluid.property(MPL::PropertyType::viscosity)
.template value<double>(vars, x_position, t, dt);
+ // Quick workaround: If fluid density is described as ideal gas, then
+ // the molar mass must be passed to the MPL::IdealGasLaw via the
+ // varialbe_array and the fluid must have the property
+ // MPL::PropertyType::molar_mass. For other density models (e.g.
+ // Constant), it is not mandatory to specify the molar mass.
+ if (fluid.hasProperty(MPL::PropertyType::molar_mass))
+ {
+ vars[static_cast<int>(MPL::Variable::molar_mass)] =
+ fluid.property(MPL::PropertyType::molar_mass)
+ .template value<double>(vars, x_position, t, dt);
+ }
auto const rho_fr =
fluid.property(MPL::PropertyType::density)
.template value<double>(vars, x_position, t, dt);
@@ -557,6 +580,17 @@ void HydroMechanicsLocalAssembler<ShapeFunctionDisplacement,
auto const mu = fluid.property(MPL::PropertyType::viscosity)
.template value<double>(vars, x_position, t, dt);
+ // Quick workaround: If fluid density is described as ideal gas, then
+ // the molar mass must be passed to the MPL::IdealGasLaw via the
+ // varialbe_array and the fluid must have the property
+ // MPL::PropertyType::molar_mass. For other density models (e.g.
+ // Constant), it is not mandatory to specify the molar mass.
+ if (fluid.hasProperty(MPL::PropertyType::molar_mass))
+ {
+ vars[static_cast<int>(MPL::Variable::molar_mass)] =
+ fluid.property(MPL::PropertyType::molar_mass)
+ .template value<double>(vars, x_position, t, dt);
+ }
auto const rho_fr =
fluid.property(MPL::PropertyType::density)
.template value<double>(vars, x_position, t, dt);
@@ -683,6 +717,17 @@ void HydroMechanicsLocalAssembler<ShapeFunctionDisplacement,
medium->property(MPL::PropertyType::porosity)
.template value<double>(vars, x_position, t, dt);
+ // Quick workaround: If fluid density is described as ideal gas, then
+ // the molar mass must be passed to the MPL::IdealGasLaw via the
+ // varialbe_array and the fluid must have the property
+ // MPL::PropertyType::molar_mass. For other density models (e.g.
+ // Constant), it is not mandatory to specify the molar mass.
+ if (fluid.hasProperty(MPL::PropertyType::molar_mass))
+ {
+ vars[static_cast<int>(MPL::Variable::molar_mass)] =
+ fluid.property(MPL::PropertyType::molar_mass)
+ .template value<double>(vars, x_position, t, dt);
+ }
auto const rho_fr =
fluid.property(MPL::PropertyType::density)
.template value<double>(vars, x_position, t, dt);
diff --git a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPLocalAssembler-impl.h b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPLocalAssembler-impl.h
index c0e0f3861a367978a7af68785bd5dae9bf8f5a06..fa778769de14a49aedc9f51c787a2b5e461e9278 100644
--- a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPLocalAssembler-impl.h
+++ b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPLocalAssembler-impl.h
@@ -122,6 +122,9 @@ void TwoPhaseFlowWithPPLocalAssembler<
pc_int_pt;
variables[static_cast<int>(MPL::Variable::temperature)] =
_process_data.temperature(t, pos)[0];
+ variables[static_cast<int>(MPL::Variable::molar_mass)] =
+ gas_phase.property(MPL::PropertyType::molar_mass)
+ .template value<double>(variables, pos, t, dt);
auto const rho_nonwet =
gas_phase.property(MPL::PropertyType::density)
diff --git a/Tests/MaterialLib/TestMPLIdealGasLaw.cpp b/Tests/MaterialLib/TestMPLIdealGasLaw.cpp
index d23cf8a61f15114215b848668d93e2188f9bc61b..36855eb94edcc6f6a06dadc0a529be31aadad5f5 100644
--- a/Tests/MaterialLib/TestMPLIdealGasLaw.cpp
+++ b/Tests/MaterialLib/TestMPLIdealGasLaw.cpp
@@ -64,6 +64,8 @@ TEST(MaterialPropertyLib, IdealGasLawOfPurePhase)
MaterialPropertyLib::Variable::phase_pressure)] = pressure_norm;
variable_array[static_cast<int>(
MaterialPropertyLib::Variable::temperature)] = temperature_norm;
+ variable_array[static_cast<int>(
+ MaterialPropertyLib::Variable::molar_mass)] = molar_mass_air;
ParameterLib::SpatialPosition const pos;
double const time = std::numeric_limits<double>::quiet_NaN();