Skip to content

GitLab

Commit d8b826d3 authored by joergbuchwald's avatar joergbuchwald Committed by Dmitry Yu. Naumov
Browse files

move porosity, transport_porosity, permeability, biot_coefficient and storage... 

move porosity, transport_porosity, permeability, biot_coefficient and storage to medium in RM process and MPL properties
parent 883db3d5
Pipeline #3784 passed with stages
in 40 minutes and 30 seconds
Hide whitespace changes
Inline Side-by-side
MaterialLib/MPL/Properties/PermeabilityOrthotropicPowerLaw.cpp
View file @ d8b826d3
......@@ -37,20 +37,11 @@ PermeabilityOrthotropicPowerLaw<DisplacementDim>::
template <int DisplacementDim>
void PermeabilityOrthotropicPowerLaw<DisplacementDim>::checkScale() const
{
if (!std::holds_alternative<Phase*>(scale_))
if (!std::holds_alternative<Medium*>(scale_))
{
OGS_FATAL(
"The property 'PermeabilityOrthotropicPowerLaw' is "
"implemented on the 'phase' scales only.");
}
auto const phase = std::get<Phase*>(scale_);
if (phase->name != "Solid")
{
OGS_FATAL(
"The property 'PermeabilityOrthotropicPowerLaw' must be given in "
"the 'Solid' phase, not in '{:s}' phase.",
phase->name);
"implemented on the 'medium' scales only.");
}
}
template <int DisplacementDim>
......@@ -64,13 +55,13 @@ PropertyDataType PermeabilityOrthotropicPowerLaw<DisplacementDim>::value(
// TODO (naumov) The phi0 must be evaluated once upon
// creation/initialization and be stored in a local state.
// For now assume porosity's initial value does not change with time.
auto const phase = std::get<Phase*>(scale_);
auto const medium = std::get<Medium*>(scale_);
auto const phi_0 =
phase->hasProperty(PropertyType::transport_porosity)
? phase->property(PropertyType::transport_porosity)
medium->hasProperty(PropertyType::transport_porosity)
? medium->property(PropertyType::transport_porosity)
.template initialValue<double>(
pos, std::numeric_limits<double>::quiet_NaN())
: phase->property(PropertyType::porosity)
: medium->property(PropertyType::porosity)
.template initialValue<double>(
pos, std::numeric_limits<double>::quiet_NaN());
......
MaterialLib/MPL/Properties/PorosityFromMassBalance.cpp
View file @ d8b826d3
......@@ -18,19 +18,11 @@ namespace MaterialPropertyLib
{
void PorosityFromMassBalance::checkScale() const
{
if (!std::holds_alternative<Phase*>(scale_))
if (!std::holds_alternative<Medium*>(scale_))
{
OGS_FATAL(
"The property 'PorosityFromMassBalance' is "
"implemented on the 'phase' scales only.");
}
auto const phase = std::get<Phase*>(scale_);
if (phase->name != "Solid")
{
OGS_FATAL(
"The property 'PorosityFromMassBalance' must be given in the "
"'Solid' phase, not in '{:s}' phase.",
phase->name);
"implemented on the 'medium' scales only.");
}
}
......@@ -53,7 +45,7 @@ PropertyDataType PorosityFromMassBalance::value(
double const beta_SR = std::get<double>(
variable_array[static_cast<int>(Variable::grain_compressibility)]);
auto const alpha_b =
std::get<Phase*>(scale_)
std::get<Medium*>(scale_)
->property(PropertyType::biot_coefficient)
.template value<double>(variable_array, pos, t, dt);
......
MaterialLib/MPL/Properties/TransportPorosityFromMassBalance.cpp
View file @ d8b826d3
......@@ -18,19 +18,11 @@ namespace MaterialPropertyLib
{
void TransportPorosityFromMassBalance::checkScale() const
{
if (!std::holds_alternative<Phase*>(scale_))
if (!std::holds_alternative<Medium*>(scale_))
{
OGS_FATAL(
"The property 'TransportPorosityFromMassBalance' is "
"implemented on the 'phase' scales only.");
}
auto const phase = std::get<Phase*>(scale_);
if (phase->name != "Solid")
{
OGS_FATAL(
"The property 'TransportPorosityFromMassBalance' must be given in "
"the 'Solid' phase, not in '{:s}' phase.",
phase->name);
"implemented on the 'medium' scales only.");
}
}
......@@ -43,7 +35,7 @@ PropertyDataType TransportPorosityFromMassBalance::value(
double const beta_SR = std::get<double>(
variable_array[static_cast<int>(Variable::grain_compressibility)]);
auto const alpha_b =
std::get<Phase*>(scale_)
std::get<Medium*>(scale_)
->property(PropertyType::biot_coefficient)
.template value<double>(variable_array, pos, t, dt);
......
ProcessLib/RichardsMechanics/CreateRichardsMechanicsProcess.cpp
View file @ d8b826d3
......@@ -35,11 +35,12 @@ void checkMPLProperties(
MaterialPropertyLib::reference_temperature,
MaterialPropertyLib::bishops_effective_stress,
MaterialPropertyLib::relative_permeability,
MaterialPropertyLib::saturation};
MaterialPropertyLib::saturation,
MaterialPropertyLib::porosity, MaterialPropertyLib::biot_coefficient
};
std::array const required_liquid_properties = {
MaterialPropertyLib::viscosity, MaterialPropertyLib::density};
std::array const required_solid_properties = {
MaterialPropertyLib::porosity, MaterialPropertyLib::biot_coefficient,
MaterialPropertyLib::density};
for (auto const& m : media)
......
ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
View file @ d8b826d3
......@@ -66,7 +66,6 @@ RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
e.getID());
auto const& medium = _process_data.media_map->getMedium(_element.getID());
auto const& solid_phase = medium->phase("Solid");
ParameterLib::SpatialPosition x_position;
x_position.setElementID(_element.getID());
......@@ -99,17 +98,17 @@ RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
// Initial porosity. Could be read from integration point data or mesh.
ip_data.porosity =
solid_phase.property(MPL::porosity)
medium->property(MPL::porosity)
.template initialValue<double>(
x_position,
std::numeric_limits<
double>::quiet_NaN() /* t independent */);
ip_data.transport_porosity = ip_data.porosity;
if (solid_phase.hasProperty(MPL::PropertyType::transport_porosity))
if (medium->hasProperty(MPL::PropertyType::transport_porosity))
{
ip_data.transport_porosity =
solid_phase.property(MPL::transport_porosity)
medium->property(MPL::transport_porosity)
.template initialValue<double>(
x_position,
std::numeric_limits<
......@@ -370,7 +369,7 @@ void RichardsMechanicsLocalAssembler<
variables[static_cast<int>(MPL::Variable::temperature)] = temperature;
auto const alpha =
solid_phase.property(MPL::PropertyType::biot_coefficient)
medium->property(MPL::PropertyType::biot_coefficient)
.template value<double>(variables, x_position, t, dt);
auto const C_el = _ip_data[ip].computeElasticTangentStiffness(
t, x_position, dt, temperature);
......@@ -430,7 +429,7 @@ void RichardsMechanicsLocalAssembler<
{ // Porosity update
variables_prev[static_cast<int>(MPL::Variable::porosity)] =
_ip_data[ip].porosity_prev;
phi = solid_phase.property(MPL::PropertyType::porosity)
phi = medium->property(MPL::PropertyType::porosity)
.template value<double>(variables, variables_prev,
x_position, t, dt);
variables[static_cast<int>(MPL::Variable::porosity)] = phi;
......@@ -474,14 +473,14 @@ void RichardsMechanicsLocalAssembler<
identity2.transpose() * C_el.inverse() * sigma_sw_prev;
}
if (solid_phase.hasProperty(MPL::PropertyType::transport_porosity))
if (medium->hasProperty(MPL::PropertyType::transport_porosity))
{
variables_prev[static_cast<int>(
MPL::Variable::transport_porosity)] =
_ip_data[ip].transport_porosity_prev;
_ip_data[ip].transport_porosity =
solid_phase.property(MPL::PropertyType::transport_porosity)
medium->property(MPL::PropertyType::transport_porosity)
.template value<double>(variables, variables_prev,
x_position, t, dt);
variables[static_cast<int>(MPL::Variable::transport_porosity)] =
......@@ -522,7 +521,7 @@ void RichardsMechanicsLocalAssembler<
_ip_data[ip].material_state_variables->getEquivalentPlasticStrain();
auto const K_intrinsic = MPL::formEigenTensor<DisplacementDim>(
solid_phase.property(MPL::PropertyType::permeability)
medium->property(MPL::PropertyType::permeability)
.value(variables, x_position, t, dt));
GlobalDimMatrixType const rho_K_over_mu =
......@@ -738,7 +737,7 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
auto& S_L = _ip_data[ip].saturation;
auto const S_L_prev = _ip_data[ip].saturation_prev;
auto const alpha =
solid_phase.property(MPL::PropertyType::biot_coefficient)
medium->property(MPL::PropertyType::biot_coefficient)
.template value<double>(variables, x_position, t, dt);
auto const C_el = _ip_data[ip].computeElasticTangentStiffness(
......@@ -799,7 +798,7 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
variables_prev[static_cast<int>(MPL::Variable::porosity)] =
_ip_data[ip].porosity_prev;
phi = solid_phase.property(MPL::PropertyType::porosity)
phi = medium->property(MPL::PropertyType::porosity)
.template value<double>(variables, variables_prev,
x_position, t, dt);
variables[static_cast<int>(MPL::Variable::porosity)] = phi;
......@@ -894,14 +893,14 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
sigma_sw = sigma_sw_prev + delta_sigma_sw;
}
if (solid_phase.hasProperty(MPL::PropertyType::transport_porosity))
if (medium->hasProperty(MPL::PropertyType::transport_porosity))
{
variables_prev[static_cast<int>(
MPL::Variable::transport_porosity)] =
_ip_data[ip].transport_porosity_prev;
_ip_data[ip].transport_porosity =
solid_phase.property(MPL::PropertyType::transport_porosity)
medium->property(MPL::PropertyType::transport_porosity)
.template value<double>(variables, variables_prev,
x_position, t, dt);
variables[static_cast<int>(MPL::Variable::transport_porosity)] =
......@@ -934,7 +933,7 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
_ip_data[ip].material_state_variables->getEquivalentPlasticStrain();
auto const K_intrinsic = MPL::formEigenTensor<DisplacementDim>(
solid_phase.property(MPL::PropertyType::permeability)
medium->property(MPL::PropertyType::permeability)
.value(variables, x_position, t, dt));
GlobalDimMatrixType const rho_Ki_over_mu = K_intrinsic * rho_LR / mu;
......@@ -1612,7 +1611,7 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
double const chi_S_L_prev = chi(S_L_prev);
auto const alpha =
solid_phase.property(MPL::PropertyType::biot_coefficient)
medium->property(MPL::PropertyType::biot_coefficient)
.template value<double>(variables, x_position, t, dt);
auto const C_el = _ip_data[ip].computeElasticTangentStiffness(
......@@ -1640,7 +1639,7 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
{ // Porosity update
variables_prev[static_cast<int>(MPL::Variable::porosity)] =
_ip_data[ip].porosity_prev;
phi = solid_phase.property(MPL::PropertyType::porosity)
phi = medium->property(MPL::PropertyType::porosity)
.template value<double>(variables, variables_prev,
x_position, t, dt);
variables[static_cast<int>(MPL::Variable::porosity)] = phi;
......@@ -1676,14 +1675,14 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
identity2.transpose() * C_el.inverse() * sigma_sw_prev;
}
if (solid_phase.hasProperty(MPL::PropertyType::transport_porosity))
if (medium->hasProperty(MPL::PropertyType::transport_porosity))
{
variables_prev[static_cast<int>(
MPL::Variable::transport_porosity)] =
_ip_data[ip].transport_porosity_prev;
_ip_data[ip].transport_porosity =
solid_phase.property(MPL::PropertyType::transport_porosity)
medium->property(MPL::PropertyType::transport_porosity)
.template value<double>(variables, variables_prev,
x_position, t, dt);
variables[static_cast<int>(MPL::Variable::transport_porosity)] =
......@@ -1720,7 +1719,7 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
_ip_data[ip].material_state_variables->getEquivalentPlasticStrain();
auto const K_intrinsic = MPL::formEigenTensor<DisplacementDim>(
solid_phase.property(MPL::PropertyType::permeability)
medium->property(MPL::PropertyType::permeability)
.value(variables, x_position, t, dt));
double const k_rel =
......
Tests/Data/RichardsMechanics/LiakopoulosHM/liakopoulos.prj
View file @ d8b826d3
<?xml version="1.0" encoding="ISO-8859-1"?>
<?xml version='1.0' encoding='ISO-8859-1'?>
<OpenGeoSysProject>
<meshes>
<mesh axially_symmetric="true">liakopoulos.vtu</mesh>
......@@ -69,25 +69,25 @@
<type>Constant</type>
<value>1e12</value>
</property>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>1</value>
</property>
<property>
<name>porosity</name>
<type>Constant</type>
<value>0.2975</value>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.5e-13</value>
</property>
</properties>
</phase>
</phases>
<properties>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>1</value>
</property>
<property>
<name>porosity</name>
<type>Constant</type>
<value>0.2975</value>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.5e-13</value>
</property>
<property>
<name>saturation</name>
<type>SaturationLiakopoulos</type>
......
Tests/Data/RichardsMechanics/LiakopoulosHM/liakopoulos_QN.prj
View file @ d8b826d3
......@@ -89,25 +89,25 @@
<type>Constant</type>
<value>1e12</value>
</property>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>1</value>
</property>
<property>
<name>porosity</name>
<type>Constant</type>
<value>0.2975</value>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.5e-13</value>
</property>
</properties>
</phase>
</phases>
<properties>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>1</value>
</property>
<property>
<name>porosity</name>
<type>Constant</type>
<value>0.2975</value>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.5e-13</value>
</property>
<property>
<name>saturation</name>
<type>SaturationLiakopoulos</type>
......
Tests/Data/RichardsMechanics/RichardsFlow_2d_quasinewton.prj
View file @ d8b826d3
<?xml version="1.0" encoding="ISO-8859-1"?>
<?xml version='1.0' encoding='ISO-8859-1'?>
<OpenGeoSysProject>
<mesh>Richards_2d.vtu</mesh>
<geometry>Richards_2d.gml</geometry>
......@@ -57,37 +57,37 @@
<phase>
<type>Solid</type>
<properties>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>0.38000000000000006</value>
</property>
<property>
<name>density</name>
<type>Constant</type>
<value>1</value>
</property>
<property>
<name>porosity</name>
<type>PorosityFromMassBalance</type>
<initial_porosity>phi0</initial_porosity>
<minimal_porosity>0</minimal_porosity>
<maximal_porosity>1</maximal_porosity>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.46e-13</value>
</property>
<property>
<name>storage</name>
<type>Constant</type>
<value> 0.0 </value>
</property>
</properties>
</phase>
</phases>
<properties>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>0.38000000000000006</value>
</property>
<property>
<name>porosity</name>
<type>PorosityFromMassBalance</type>
<initial_porosity>phi0</initial_porosity>
<minimal_porosity>0</minimal_porosity>
<maximal_porosity>1</maximal_porosity>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.46e-13</value>
</property>
<property>
<name>storage</name>
<type>Constant</type>
<value> 0.0 </value>
</property>
<property>
<name>reference_temperature</name>
<type>Constant</type>
......
Tests/Data/RichardsMechanics/RichardsFlow_2d_richardsflow.prj
View file @ d8b826d3
<?xml version="1.0" encoding="ISO-8859-1"?>
<?xml version='1.0' encoding='ISO-8859-1'?>
<OpenGeoSysProject>
<mesh>Richards_2d_linear.vtu</mesh>
<geometry>Richards_2d.gml</geometry>
......
Tests/Data/RichardsMechanics/RichardsFlow_2d_small.prj
View file @ d8b826d3
<?xml version="1.0" encoding="ISO-8859-1"?>
<?xml version='1.0' encoding='ISO-8859-1'?>
<OpenGeoSysProject>
<mesh>RichardsFlow_2d_small.vtu</mesh>
<geometry>RichardsFlow_2d_small.gml</geometry>
......@@ -52,37 +52,37 @@
<phase>
<type>Solid</type>
<properties>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>0.38</value>
</property>
<property>
<name>density</name>
<type>Constant</type>
<value>1</value>
</property>
<property>
<name>porosity</name>
<type>PorosityFromMassBalance</type>
<initial_porosity>phi0</initial_porosity>
<minimal_porosity>0</minimal_porosity>
<maximal_porosity>1</maximal_porosity>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.46e-13</value>
</property>
<property>
<name>storage</name>
<type>Constant</type>
<value> 0.0 </value>
</property>
</properties>
</phase>
</phases>
<properties>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>0.38</value>
</property>
<property>
<name>porosity</name>
<type>PorosityFromMassBalance</type>
<initial_porosity>phi0</initial_porosity>
<minimal_porosity>0</minimal_porosity>
<maximal_porosity>1</maximal_porosity>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.46e-13</value>
</property>
<property>
<name>storage</name>
<type>Constant</type>
<value> 0.0 </value>
</property>
<property>
<name>reference_temperature</name>
<type>Constant</type>
......
Tests/Data/RichardsMechanics/RichardsFlow_2d_small_masslumping.prj
View file @ d8b826d3
<?xml version="1.0" encoding="ISO-8859-1"?>
<?xml version='1.0' encoding='ISO-8859-1'?>
<OpenGeoSysProject>
<mesh>Richards_2d.vtu</mesh>
<geometry>Richards_2d.gml</geometry>
......@@ -58,32 +58,32 @@
<type>Constant</type>
<value>1</value>
</property>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>0.38</value>
</property>
<property>
<name>porosity</name>
<type>PorosityFromMassBalance</type>
<initial_porosity>phi0</initial_porosity>
<minimal_porosity>0</minimal_porosity>
<maximal_porosity>1</maximal_porosity>
</property>
<property>
<name>storage</name>
<type>Constant</type>
<value> 0.0 </value>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.46e-13</value>
</property>
</properties>
</phase>
</phases>
<properties>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>0.38</value>
</property>
<property>
<name>porosity</name>
<type>PorosityFromMassBalance</type>
<initial_porosity>phi0</initial_porosity>
<minimal_porosity>0</minimal_porosity>
<maximal_porosity>1</maximal_porosity>
</property>
<property>
<name>storage</name>
<type>Constant</type>
<value> 0.0 </value>
</property>
<property>
<name>permeability</name>
<type>Constant</type>
<value>4.46e-13</value>
</property>
<property>
<name>reference_temperature</name>
<type>Constant</type>
......
Tests/Data/RichardsMechanics/alternative_mass_balance_anzInterval_10.prj
View file @ d8b826d3
<?xml version="1.0" encoding="ISO-8859-1"?>
<?xml version='1.0' encoding='ISO-8859-1'?>
<OpenGeoSysProject>
<mesh>square_1x1_quad_1e0.vtu</mesh>
<geometry>square_1x1.gml</geometry>
......@@ -57,27 +57,27 @@
<type>Constant</type>
<value>2.6e3</value>
</property>
<property>
<name>biot_coefficient</name>
<type>Constant</type>
<value>0.72222222222222222222</value>
</property>
<property>
<name>porosity</name>
<type>PorosityFromMassBalance</type>
<initial_porosity>phi0</initial_porosity>
<minimal_porosity>0</minimal_porosity>