OGS
TransportPorosityFromMassBalance.cpp
Go to the documentation of this file.
1
11
12#include <algorithm>
13#include <cmath>
14
16
17namespace MaterialPropertyLib
18{
20{
21 if (!std::holds_alternative<Medium*>(scale_))
22 {
24 "The property 'TransportPorosityFromMassBalance' is "
25 "implemented on the 'medium' scales only.");
26 }
27}
28
30 VariableArray const& variable_array,
31 VariableArray const& variable_array_prev,
32 ParameterLib::SpatialPosition const& pos, double const t,
33 double const dt) const
34{
35 double const beta_SR = variable_array.grain_compressibility;
36 auto const alpha_b =
37 std::get<Medium*>(scale_)
39 .template value<double>(variable_array, pos, t, dt);
40
41 double const e = variable_array.volumetric_strain;
42 double const e_prev = variable_array_prev.volumetric_strain;
43 double const delta_e = e - e_prev;
44
45 double const p_eff = variable_array.effective_pore_pressure;
46 double const p_eff_prev = variable_array_prev.effective_pore_pressure;
47 double const delta_p_eff = p_eff - p_eff_prev;
48
49 double const phi = variable_array.porosity;
50
51 double const phi_tr_prev = variable_array_prev.transport_porosity;
52
53 double const w = delta_e + delta_p_eff * beta_SR;
54 return std::clamp(phi_tr_prev + (alpha_b - phi) * w, phi_min_, phi_max_);
55}
56
58 VariableArray const& /*variable_array*/,
59 ParameterLib::SpatialPosition const& /*pos*/, double const /*t*/,
60 double const /*dt*/) const
61{
63 "TransportPorosityFromMassBalance value call requires previous time "
64 "step values.");
65}
66
68 VariableArray const& /*variable_array*/, Variable const /*variable*/,
69 ParameterLib::SpatialPosition const& /*pos*/, double const /*t*/,
70 double const /*dt*/) const
71{
73 "TransportPorosityFromMassBalance derivatives are not implemented.");
74}
75
76} // namespace MaterialPropertyLib
#define OGS_FATAL(...)
Definition Error.h:26
virtual PropertyDataType value() const
Definition Property.cpp:76
std::variant< Medium *, Phase *, Component * > scale_
Definition Property.h:297
PropertyDataType dValue(VariableArray const &variable_array, Variable const variable, ParameterLib::SpatialPosition const &pos, double const t, double const dt) const override
std::variant< double, Eigen::Matrix< double, 2, 1 >, Eigen::Matrix< double, 3, 1 >, Eigen::Matrix< double, 2, 2 >, Eigen::Matrix< double, 3, 3 >, Eigen::Matrix< double, 4, 1 >, Eigen::Matrix< double, 6, 1 >, Eigen::MatrixXd > PropertyDataType
Definition Property.h:31