From 27a5fd1b4a0b4420bde61f0e20f249a19619621b Mon Sep 17 00:00:00 2001
From: Florian Zill <florian.zill@ufz.de>
Date: Wed, 29 Jan 2020 16:56:53 +0100
Subject: [PATCH] [HM] proper calculation of compressibility
removed compressibility from process_data
---
.../CreateHydroMechanicsProcess.cpp | 11 ++-------
.../HydroMechanics/HydroMechanicsFEM-impl.h | 23 +++++++++---------
.../HydroMechanicsProcessData.h | 24 -------------------
3 files changed, 13 insertions(+), 45 deletions(-)
diff --git a/ProcessLib/HydroMechanics/CreateHydroMechanicsProcess.cpp b/ProcessLib/HydroMechanics/CreateHydroMechanicsProcess.cpp
index 8652773a330..2a5b33ecb18 100644
--- a/ProcessLib/HydroMechanics/CreateHydroMechanicsProcess.cpp
+++ b/ProcessLib/HydroMechanics/CreateHydroMechanicsProcess.cpp
@@ -169,13 +169,6 @@ std::unique_ptr<Process> createHydroMechanicsProcess(
"specific_gas_constant", std::numeric_limits<double>::quiet_NaN());
DBUG("Use 'specific_gas_constant' as specific gas constant.");
- // Fluid compressibility
- double const fluid_compressibility =
- //! \ogs_file_param{prj__processes__process__HYDRO_MECHANICS__fluid_compressibility}
- config.getConfigParameter<double>(
- "fluid_compressibility", std::numeric_limits<double>::quiet_NaN());
- DBUG("Use 'fluid_compressibility' as fluid compressibility parameter.");
-
auto const fluid_type =
FluidType::strToFluidType(
//! \ogs_file_param{prj__processes__process__HYDRO_MECHANICS__fluid_type}
@@ -193,8 +186,8 @@ std::unique_ptr<Process> createHydroMechanicsProcess(
materialIDs(mesh), std::move(media_map),
std::move(solid_constitutive_relations),
initial_stress, specific_body_force,
- fluid_compressibility, reference_temperature,
- specific_gas_constant, fluid_type};
+ reference_temperature, specific_gas_constant,
+ fluid_type};
SecondaryVariableCollection secondary_variables;
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
index 58d3263baf3..ef20b4c4919 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h
@@ -228,14 +228,14 @@ void HydroMechanicsLocalAssembler<ShapeFunctionDisplacement,
auto const K_S = solid_material.getBulkModulus(t, x_position);
auto const rho_sr = solid_phase.property(MPL::PropertyType::density)
.template value<double>(vars, x_position, t);
- // TODO (FZill) get fluid properties from GPML
- double const p_fr =
- (_process_data.fluid_type == FluidType::Fluid_Type::IDEAL_GAS)
- ? N_p.dot(p)
- : std::numeric_limits<double>::quiet_NaN();
- double const beta_p = _process_data.getFluidCompressibility(p_fr);
+
auto const rho_fr = gas_phase.property(MPL::PropertyType::density)
.template value<double>(vars, x_position, t);
+ auto const beta_p =
+ gas_phase.property(MPL::PropertyType::density)
+ .template dValue<double>(vars, MPL::Variable::phase_pressure,
+ x_position, t) / rho_fr;
+
auto const porosity = solid_phase.property(MPL::PropertyType::porosity)
.template value<double>(vars, x_position, t);
auto const& identity2 = MathLib::KelvinVector::Invariants<
@@ -450,14 +450,13 @@ void HydroMechanicsLocalAssembler<ShapeFunctionDisplacement,
auto const alpha_b =
solid_phase.property(MPL::PropertyType::biot_coefficient)
.template value<double>(vars, x_position, t);
- // TODO (FZill) get fluid properties from GPML
- double const p_fr =
- (_process_data.fluid_type == FluidType::Fluid_Type::IDEAL_GAS)
- ? N_p.dot(p)
- : std::numeric_limits<double>::quiet_NaN();
- double const beta_p = _process_data.getFluidCompressibility(p_fr);
+
auto const rho_fr = gas_phase.property(MPL::PropertyType::density)
.template value<double>(vars, x_position, t);
+ auto const beta_p =
+ gas_phase.property(MPL::PropertyType::density)
+ .template dValue<double>(vars, MPL::Variable::phase_pressure,
+ x_position, t) / rho_fr;
auto const porosity = solid_phase.property(MPL::PropertyType::porosity)
.template value<double>(vars, x_position, t);
auto const K_S = solid_material.getBulkModulus(t, x_position);
diff --git a/ProcessLib/HydroMechanics/HydroMechanicsProcessData.h b/ProcessLib/HydroMechanics/HydroMechanicsProcessData.h
index f9cd065565a..af8a6fab75d 100644
--- a/ProcessLib/HydroMechanics/HydroMechanicsProcessData.h
+++ b/ProcessLib/HydroMechanics/HydroMechanicsProcessData.h
@@ -56,11 +56,6 @@ struct HydroMechanicsProcessData
/// A vector of displacement dimension's length.
Eigen::Matrix<double, DisplacementDim, 1> const specific_body_force;
- /// Fluid's compressibility. A scalar quantity.
- /// Only used for compressible_fluid fluid_type
- double const fluid_compressibility =
- std::numeric_limits<double>::quiet_NaN();
-
/// Reference Temperature. A scalar quantity.
/// Only used for ideal_gas fluid_type
double const reference_temperature =
@@ -75,25 +70,6 @@ struct HydroMechanicsProcessData
/// incompressible_fluid, compressible_fluid, ideal_gas
FluidType::Fluid_Type const fluid_type;
-
- /// will be removed after linking with MPL
- double getFluidCompressibility(double const& p_fr)
- {
- if (fluid_type == FluidType::Fluid_Type::INCOMPRESSIBLE_FLUID)
- {
- return 0.0;
- }
- if (fluid_type == FluidType::Fluid_Type::COMPRESSIBLE_FLUID)
- {
- return fluid_compressibility;
- }
- if (fluid_type == FluidType::Fluid_Type::IDEAL_GAS)
- {
- return 1.0 / p_fr;
- }
- OGS_FATAL("unknown fluid type %d", static_cast<int> (fluid_type));
- }
-
MeshLib::PropertyVector<double>* pressure_interpolated = nullptr;
std::array<MeshLib::PropertyVector<double>*, 3> principal_stress_vector = {
nullptr, nullptr, nullptr};
--
GitLab