10#include <pybind11/eigen.h>
11#include <pybind11/pybind11.h>
12#include <pybind11/stl.h>
14#include <range/v3/range/conversion.hpp>
15#include <range/v3/view/transform.hpp>
25namespace py = pybind11;
36 py::object
const& coordinates)
38 if (coordinates.is_none())
42 return std::make_optional<MathLib::Point3d>(
43 py::cast<std::array<double, 3>>(coordinates));
54 return py::array_t<double>(3, opt_coords->data());
61 py::class_<SpatialPosition>(m,
"SpatialPosition",
63 Describes a spatial position within a mesh or domain.
65 A SpatialPosition may refer to a node (via node ID), an element (via element ID),
66 or a physical point in space (via coordinates). It is typically used when evaluating
67 material properties that depend on the spatial context within a simulation.
69 The position may be partially specified. For example:
70 - Only coordinates for geometric evaluation
71 - Only node ID for nodal properties
72 - Only element ID for element-wise values
73 - Or any combination of the above
75 The class is used as a parameter to property evaluations in the OGS material properties library.
78 [](std::optional<std::size_t> node_id,
79 std::optional<std::size_t>
81 py::object const& coordinates)
87 py::arg(
"node_id") = std::nullopt,
88 py::arg(
"element_id") = std::nullopt,
89 py::arg(
"coords") = py::none(),
91 SpatialPosition(node_id=None, element_id=None, coords=None)
94 node_id (int, optional): Node ID
95 element_id (int, optional): Element ID
96 coords (array-like of 3 floats, optional): Coordinates
105 Node ID of the spatial position.
107 This property can be read and set from Python. Setting to None is not supported.
116 Element ID of the spatial position.
118 This property can be read and set from Python. Setting to None is not supported.
130 Coordinates of the spatial position as a 3-element array.
132 This property can be read and set from Python. Use a 3-element list, tuple,
133 or NumPy array. Setting to None is not supported.
142 auto const element_id =
147 return "<SpatialPosition(" + node_id +
", " + element_id +
148 ", " + py::str(coords).cast<std::string>() +
")>";
151 Return string representation of the SpatialPosition.
157 py::class_<VariableArray>(m,
"VariableArray")
161 .def_readwrite(
"deformation_gradient",
164 .def_readwrite(
"effective_pore_pressure",
167 .def_readwrite(
"enthalpy_of_evaporation",
169 .def_readwrite(
"equivalent_plastic_strain",
172 .def_readwrite(
"grain_compressibility",
174 .def_readwrite(
"liquid_phase_pressure",
179 .def_readwrite(
"molar_mass_derivative",
184 .def_readwrite(
"solid_grain_pressure",
197 py::enum_<Variable>(m,
"Variable")
198 .value(
"capillary_pressure", Variable::capillary_pressure)
199 .value(
"concentration", Variable::concentration)
200 .value(
"deformation_gradient", Variable::deformation_gradient)
201 .value(
"density", Variable::density)
202 .value(
"effective_pore_pressure", Variable::effective_pore_pressure)
203 .value(
"enthalpy", Variable::enthalpy)
204 .value(
"enthalpy_of_evaporation", Variable::enthalpy_of_evaporation)
205 .value(
"equivalent_plastic_strain", Variable::equivalent_plastic_strain)
206 .value(
"fracture_aperture", Variable::fracture_aperture)
207 .value(
"grain_compressibility", Variable::grain_compressibility)
208 .value(
"liquid_phase_pressure", Variable::liquid_phase_pressure)
209 .value(
"liquid_saturation", Variable::liquid_saturation)
210 .value(
"mechanical_strain", Variable::mechanical_strain)
211 .value(
"molar_mass", Variable::molar_mass)
212 .value(
"molar_mass_derivative", Variable::molar_mass_derivative)
213 .value(
"molar_fraction", Variable::molar_fraction)
214 .value(
"gas_phase_pressure", Variable::gas_phase_pressure)
215 .value(
"porosity", Variable::porosity)
216 .value(
"solid_grain_pressure", Variable::solid_grain_pressure)
217 .value(
"stress", Variable::stress)
218 .value(
"temperature", Variable::temperature)
219 .value(
"total_strain", Variable::total_strain)
220 .value(
"total_stress", Variable::total_stress)
221 .value(
"transport_porosity", Variable::transport_porosity)
222 .value(
"vapour_pressure", Variable::vapour_pressure)
223 .value(
"volumetric_strain", Variable::volumetric_strain)
229 py::class_<Property>(m,
"Property", R
"pbdoc(
230 Base class for material properties.
236 {
return p.
value(va, pos, t, dt); },
237 py::arg(
"variable_array"),
242 Evaluate the property value.
245 variable_array: Current time step values.
246 pos: Spatial position.
255 double dt) {
return p.
value(va, vap, pos, t, dt); },
256 py::arg(
"variable_array"),
257 py::arg(
"variable_array_prev"),
262 Evaluate the property value with previous time step data.
265 variable_array: Current time step values.
266 variable_array_prev: Previous time step values.
267 pos: Spatial position.
277 {
return p.
dValue(va, vap, var, pos, t, dt); },
278 py::arg(
"variable_array"),
279 py::arg(
"variable_array_prev"),
285 Evaluate the derivative of the property with respect to a variable,
286 using previous time step data.
289 variable_array: Current time step values.
290 variable_array_prev: Previous time step values.
291 variable: Variable to differentiate with respect to.
292 pos: Spatial position.
301 {
return p.
dValue(va, var, pos, t, dt); },
302 py::arg(
"variable_array"),
308 Evaluate the derivative of the property with respect to a variable.
311 variable_array: Current time step values.
312 variable: Variable to differentiate with respect to.
313 pos: Spatial position.
321 py::class_<Constant, Property>(m,
"Constant", R
"pbdoc(
322 Constant property with a fixed value.
324 .def(py::init<std::string, PropertyDataType>(),
328 Construct a Constant property.
331 name: Name of the property.
332 value: Constant value.
338 py::class_<Linear, Property>(m,
"Linear", R
"pbdoc(
339 Linearly interpolated property dependent on independent variables.
348 {
return std::make_from_tuple<IndependentVariable>(iv); };
350 auto independent_variables =
351 ivs | ranges::views::transform(makeIV) |
352 ranges::to<std::vector>();
353 return std::make_unique<Linear>(
name, reference_value,
354 independent_variables);
357 py::arg(
"reference_value"),
358 py::arg(
"independent_variables"),
360 Construct a Linear property.
363 name: Name of the property.
364 reference_value: Base property value.
365 independent_variables: List of tuples (variable, reference condition, slope), each defining a linear dependency on a variable.
378 m.attr(
"__name__") =
"ogs.mpl";
379 m.doc() =
"pybind11 ogs MPL bindings";
Definition of the Point3d class.
virtual PropertyDataType value() const
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
KelvinVector mechanical_strain
DeformationGradient deformation_gradient
KelvinVector total_stress
double solid_grain_pressure
double transport_porosity
double grain_compressibility
double molar_mass_derivative
KelvinVector total_strain
double gas_phase_pressure
double effective_pore_pressure
double equivalent_plastic_strain
double capillary_pressure
double enthalpy_of_evaporation
double liquid_phase_pressure
std::optional< std::size_t > getNodeID() const
void setNodeID(std::size_t node_id)
std::optional< std::size_t > getElementID() const
void setCoordinates(MathLib::Point3d const &coordinates)
void setElementID(std::size_t element_id)
std::optional< MathLib::Point3d > const getCoordinates() const
void initOGSLogger(std::string const &log_level)
std::variant< std::monostate, double, Eigen::Vector< double, 4 >, Eigen::Vector< double, 5 >, Eigen::Vector< double, 6 >, Eigen::Vector< double, 9 > > VariableType
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
void bindSpatialPosition(py::module_ &m)
void bindConstant(py::module_ &m)
void bindLinear(py::module_ &m)
void bindVariableArray(py::module_ &m)
void bindVariableEnum(py::module_ &m)
std::optional< MathLib::Point3d > pythonToSpatialPositionCoords(py::object const &coordinates)
py::object spatialPositionCoordsToPython(SpatialPosition const &pos)
void bindProperty(py::module_ &m)