OGS
Property.cpp
Go to the documentation of this file.
1 
13 #include "Property.h"
14 
15 #include <string>
16 
17 #include "Component.h"
18 #include "Medium.h"
19 #include "Phase.h"
20 
21 namespace MaterialPropertyLib
22 {
23 PropertyDataType fromVector(std::vector<double> const& values)
24 {
25  switch (values.size())
26  {
27  case 1:
28  {
29  return values[0];
30  }
31  case 2:
32  {
33  return Eigen::Vector2d{values[0], values[1]};
34  }
35  case 3:
36  {
37  return Eigen::Vector3d{values[0], values[1], values[2]};
38  }
39  case 4:
40  {
41  using M = Eigen::Matrix2d;
42  return M{Eigen::Map<M const>{values.data(), 2, 2}};
43  }
44  case 6:
45  {
46  // Symmetric Tensor - xx, yy, zz, xy, xz, yz
47  using M = Eigen::Matrix<double, 6, 1>;
48  return M{Eigen::Map<M const>{values.data(), 6}};
49  }
50  case 9:
51  {
52  using M = Eigen::Matrix3d;
53  return M{Eigen::Map<M const>{values.data(), 3, 3}};
54  }
55  default:
56  {
57  OGS_FATAL(
58  "Conversion of a {:d}-vector to PropertyDataType is not "
59  "implemented.",
60  values.size());
61  }
62  }
63 }
64 
66  ParameterLib::SpatialPosition const& pos, double const t) const
67 {
68  return value(VariableArray{}, pos, t,
69  std::numeric_limits<double>::quiet_NaN());
70 }
71 
73 {
74 #ifndef NDEBUG
75  property_used = true;
76 #endif
77  return value_;
78 }
79 
81  VariableArray const& /*variable_array_prev*/,
82  ParameterLib::SpatialPosition const& /*pos*/,
83  double const /*t*/, double const /*dt*/) const
84 {
85 #ifndef NDEBUG
86  property_used = true;
87 #endif
88  return value_;
89 }
90 
93  double const t, double const dt) const
94 {
95 #ifndef NDEBUG
96  property_used = true;
97 #endif
98  return value(variable_array, VariableArray{}, pos, t, dt);
99 }
100 
102  VariableArray const& /*variable_array_prev*/,
103  Variable const /*variable*/,
104  ParameterLib::SpatialPosition const& /*pos*/,
105  double const /*t*/, double const /*dt*/) const
106 {
107 #ifndef NDEBUG
108  property_used = true;
109 #endif
110  return dvalue_;
111 }
112 
116  Variable const variable,
118  double const t, double const dt) const
119 {
120 #ifndef NDEBUG
121  property_used = true;
122 #endif
123  return dValue(variable_array, VariableArray{}, variable, pos, t, dt);
124 }
125 
128  Variable const /*variable*/,
129  Variable const /*variable*/,
130  ParameterLib::SpatialPosition const& /*pos*/,
131  double const /*t*/,
132  double const /*dt*/) const
133 {
134 #ifndef NDEBUG
135  property_used = true;
136 #endif
137  return 0.0;
138 }
139 
140 std::string Property::description() const
141 {
142  return "property '" + name_ + "' defined for " +
143  std::visit([](auto&& scale) -> std::string
144  { return scale->description(); },
145  scale_);
146 }
147 } // namespace MaterialPropertyLib
#define OGS_FATAL(...)
Definition: Error.h:26
virtual PropertyDataType d2Value(VariableArray const &variable_array, Variable const variable1, Variable const variable2, ParameterLib::SpatialPosition const &pos, double const t, double const dt) const
Default implementation: 2nd derivative of any constant property is zero.
Definition: Property.cpp:127
std::string description() const
Definition: Property.cpp:140
PropertyDataType value_
The single value of a property.
Definition: Property.h:282
PropertyDataType dvalue_
Definition: Property.h:283
virtual PropertyDataType value() const
Definition: Property.cpp:72
std::variant< Medium *, Phase *, Component * > scale_
Definition: Property.h:287
virtual PropertyDataType initialValue(ParameterLib::SpatialPosition const &pos, double const t) const
Definition: Property.cpp:65
virtual PropertyDataType dValue(VariableArray const &variable_array, VariableArray const &variable_array_prev, Variable const variable, ParameterLib::SpatialPosition const &pos, double const t, double const dt) const
Definition: Property.cpp:101
PropertyDataType fromVector(std::vector< double > const &values)
Definition: Property.cpp:23
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 > > PropertyDataType
Definition: Property.h:35
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
Definition: VariableType.h:108
void scale(PETScVector &x, double const a)
Definition: LinAlg.cpp:44