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>
24namespace py = pybind11;
35 py::object
const& coordinates)
37 if (coordinates.is_none())
41 return std::make_optional<MathLib::Point3d>(
42 py::cast<std::array<double, 3>>(coordinates));
53 return py::array_t<double>(3, opt_coords->data());
60 py::class_<SpatialPosition>(m,
"SpatialPosition",
62 Describes a spatial position within a mesh or domain.
64 A SpatialPosition may refer to a node (via node ID), an element (via element ID),
65 or a physical point in space (via coordinates). It is typically used when evaluating
66 material properties that depend on the spatial context within a simulation.
68 The position may be partially specified. For example:
69 - Only coordinates for geometric evaluation
70 - Only node ID for nodal properties
71 - Only element ID for element-wise values
72 - Or any combination of the above
74 The class is used as a parameter to property evaluations in the OGS material properties library.
77 [](std::optional<std::size_t> node_id,
78 std::optional<std::size_t>
80 py::object const& coordinates)
86 py::arg(
"node_id") = std::nullopt,
87 py::arg(
"element_id") = std::nullopt,
88 py::arg(
"coords") = py::none(),
90 SpatialPosition(node_id=None, element_id=None, coords=None)
93 node_id (int, optional): Node ID
94 element_id (int, optional): Element ID
95 coords (array-like of 3 floats, optional): Coordinates
104 Node ID of the spatial position.
106 This property can be read and set from Python. Setting to None is not supported.
115 Element ID of the spatial position.
117 This property can be read and set from Python. Setting to None is not supported.
129 Coordinates of the spatial position as a 3-element array.
131 This property can be read and set from Python. Use a 3-element list, tuple,
132 or NumPy array. Setting to None is not supported.
141 auto const element_id =
146 return "<SpatialPosition(" + node_id +
", " + element_id +
147 ", " + py::str(coords).cast<std::string>() +
")>";
150 Return string representation of the SpatialPosition.
156 py::class_<VariableArray>(m,
"VariableArray")
160 .def_readwrite(
"deformation_gradient",
163 .def_readwrite(
"effective_pore_pressure",
166 .def_readwrite(
"enthalpy_of_evaporation",
168 .def_readwrite(
"equivalent_plastic_strain",
171 .def_readwrite(
"grain_compressibility",
173 .def_readwrite(
"liquid_phase_pressure",
178 .def_readwrite(
"molar_mass_derivative",
183 .def_readwrite(
"solid_grain_pressure",
196 py::enum_<Variable>(m,
"Variable")
197 .value(
"capillary_pressure", Variable::capillary_pressure)
198 .value(
"concentration", Variable::concentration)
199 .value(
"deformation_gradient", Variable::deformation_gradient)
200 .value(
"density", Variable::density)
201 .value(
"effective_pore_pressure", Variable::effective_pore_pressure)
202 .value(
"enthalpy", Variable::enthalpy)
203 .value(
"enthalpy_of_evaporation", Variable::enthalpy_of_evaporation)
204 .value(
"equivalent_plastic_strain", Variable::equivalent_plastic_strain)
205 .value(
"fracture_aperture", Variable::fracture_aperture)
206 .value(
"grain_compressibility", Variable::grain_compressibility)
207 .value(
"liquid_phase_pressure", Variable::liquid_phase_pressure)
208 .value(
"liquid_saturation", Variable::liquid_saturation)
209 .value(
"mechanical_strain", Variable::mechanical_strain)
210 .value(
"molar_mass", Variable::molar_mass)
211 .value(
"molar_mass_derivative", Variable::molar_mass_derivative)
212 .value(
"molar_fraction", Variable::molar_fraction)
213 .value(
"gas_phase_pressure", Variable::gas_phase_pressure)
214 .value(
"porosity", Variable::porosity)
215 .value(
"solid_grain_pressure", Variable::solid_grain_pressure)
216 .value(
"stress", Variable::stress)
217 .value(
"temperature", Variable::temperature)
218 .value(
"total_strain", Variable::total_strain)
219 .value(
"total_stress", Variable::total_stress)
220 .value(
"transport_porosity", Variable::transport_porosity)
221 .value(
"vapour_pressure", Variable::vapour_pressure)
222 .value(
"volumetric_strain", Variable::volumetric_strain)
228 py::class_<Property>(m,
"Property", R
"pbdoc(
229 Base class for material properties.
235 {
return p.value(va, pos, t, dt); },
236 py::arg(
"variable_array"),
241 Evaluate the property value.
244 variable_array: Current time step values.
245 pos: Spatial position.
254 double dt) {
return p.value(va, vap, pos, t, dt); },
245 pos: Spatial position. {
…}
255 py::arg(
"variable_array"),
256 py::arg(
"variable_array_prev"),
261 Evaluate the property value with previous time step data.
264 variable_array: Current time step values.
265 variable_array_prev: Previous time step values.
266 pos: Spatial position.
276 {
return p.dValue(va, vap, var, pos, t, dt); },
277 py::arg(
"variable_array"),
278 py::arg(
"variable_array_prev"),
256 py::arg(
"variable_array_prev"), {
…}
284 Evaluate the derivative of the property with respect to a variable,
285 using previous time step data.
288 variable_array: Current time step values.
289 variable_array_prev: Previous time step values.
290 variable: Variable to differentiate with respect to.
291 pos: Spatial position.
300 {
return p.dValue(va, var, pos, t, dt); },
301 py::arg(
"variable_array"),
307 Evaluate the derivative of the property with respect to a variable.
310 variable_array: Current time step values.
311 variable: Variable to differentiate with respect to.
312 pos: Spatial position.
320 py::class_<Constant, Property>(m,
"Constant", R
"pbdoc(
321 Constant property with a fixed value.
323 .def(py::init<std::string, PropertyDataType>(),
327 Construct a Constant property.
330 name: Name of the property.
331 value: Constant value.
337 py::class_<Linear, Property>(m,
"Linear", R
"pbdoc(
338 Linearly interpolated property dependent on independent variables.
347 {
return std::make_from_tuple<IndependentVariable>(iv); };
349 auto independent_variables =
350 ivs | ranges::views::transform(makeIV) |
351 ranges::to<std::vector>();
352 return std::make_unique<Linear>(
name, reference_value,
353 independent_variables);
356 py::arg(
"reference_value"),
357 py::arg(
"independent_variables"),
359 Construct a Linear property.
362 name: Name of the property.
363 reference_value: Base property value.
364 independent_variables: List of tuples (variable, reference condition, slope), each defining a linear dependency on a variable.
371 m.attr(
"__name__") =
"ogs.mpl";
372 m.doc() =
"pybind11 ogs MPL bindings";
Definition of the Point3d class.
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
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)