OGS
ogs_mpl.cpp
Go to the documentation of this file.
1
10#include <pybind11/eigen.h>
11#include <pybind11/pybind11.h>
12#include <pybind11/stl.h>
13
14#include <range/v3/range/conversion.hpp>
15#include <range/v3/view/transform.hpp>
16
17#include "BaseLib/Logging.h"
18#include "BaseLib/MPI.h"
22#include "MathLib/Point3d.h"
24
25namespace py = pybind11;
26
27using namespace MaterialPropertyLib;
28using namespace ParameterLib;
29
30namespace
31{
32// Convert NumPy array to Point3d if not None. This function is there to avoid
33// pulling the MeshLib::Point3d object into python, but use std::array (which
34// has all automatic conversions in pybind11).
35std::optional<MathLib::Point3d> pythonToSpatialPositionCoords(
36 py::object const& coordinates)
37{
38 if (coordinates.is_none())
39 {
40 return std::nullopt;
41 }
42 return std::make_optional<MathLib::Point3d>(
43 py::cast<std::array<double, 3>>(coordinates));
44}
45
46// Convert SpatialPositon to python array if coordinates are set. This function
47// is there to avoid pulling the MeshLib::Point3d object into python, but use
48// python arrays.
50{
51 auto const& opt_coords = pos.getCoordinates();
52 if (opt_coords)
53 {
54 return py::array_t<double>(3, opt_coords->data());
55 }
56 return py::none();
57}
58
59void bindSpatialPosition(py::module_& m)
60{
61 py::class_<SpatialPosition>(m, "SpatialPosition",
62 R"pbdoc(
63 Describes a spatial position within a mesh or domain.
64
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.
68
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
74
75 The class is used as a parameter to property evaluations in the OGS material properties library.
76)pbdoc")
77 .def(py::init(
78 [](std::optional<std::size_t> node_id,
79 std::optional<std::size_t>
80 element_id,
81 py::object const& coordinates)
82 {
83 return SpatialPosition(
84 node_id, element_id,
86 }),
87 py::arg("node_id") = std::nullopt,
88 py::arg("element_id") = std::nullopt,
89 py::arg("coords") = py::none(),
90 R"pbdoc(
91 SpatialPosition(node_id=None, element_id=None, coords=None)
92
93 Parameters:
94 node_id (int, optional): Node ID
95 element_id (int, optional): Element ID
96 coords (array-like of 3 floats, optional): Coordinates
97 )pbdoc")
98
99 .def_property(
100 "node_id",
101 [](SpatialPosition const& pos) { return pos.getNodeID(); },
102 [](SpatialPosition& pos, std::size_t const id)
103 { pos.setNodeID(id); },
104 R"pbdoc(
105 Node ID of the spatial position.
106
107 This property can be read and set from Python. Setting to None is not supported.
108 )pbdoc")
109
110 .def_property(
111 "element_id",
112 [](SpatialPosition const& pos) { return pos.getElementID(); },
113 [](SpatialPosition& pos, std::size_t const id)
114 { pos.setElementID(id); },
115 R"pbdoc(
116 Element ID of the spatial position.
117
118 This property can be read and set from Python. Setting to None is not supported.
119 )pbdoc")
120
121 .def_property(
122 "coordinates",
124 [](SpatialPosition& pos, py::object const& coordinates)
125 {
126 pos.setCoordinates(
127 pythonToSpatialPositionCoords(coordinates).value());
128 },
129 R"pbdoc(
130 Coordinates of the spatial position as a 3-element array.
131
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.
134 )pbdoc")
135
136 .def(
137 "__repr__",
138 [](SpatialPosition const& pos)
139 {
140 auto const node_id =
141 pos.getNodeID() ? std::to_string(*pos.getNodeID()) : "None";
142 auto const element_id =
143 pos.getElementID() ? std::to_string(*pos.getElementID())
144 : "None";
145 auto const coords = spatialPositionCoordsToPython(pos);
146
147 return "<SpatialPosition(" + node_id + ", " + element_id +
148 ", " + py::str(coords).cast<std::string>() + ")>";
149 },
150 R"pbdoc(
151 Return string representation of the SpatialPosition.
152 )pbdoc");
153}
154
155void bindVariableArray(py::module_& m)
156{
157 py::class_<VariableArray>(m, "VariableArray")
158 .def(py::init<>())
159 .def_readwrite("capillary_pressure", &VariableArray::capillary_pressure)
160 .def_readwrite("concentration", &VariableArray::concentration)
161 .def_readwrite("deformation_gradient",
163 .def_readwrite("density", &VariableArray::density)
164 .def_readwrite("effective_pore_pressure",
166 .def_readwrite("enthalpy", &VariableArray::enthalpy)
167 .def_readwrite("enthalpy_of_evaporation",
169 .def_readwrite("equivalent_plastic_strain",
171 .def_readwrite("fracture_aperture", &VariableArray::fracture_aperture)
172 .def_readwrite("grain_compressibility",
174 .def_readwrite("liquid_phase_pressure",
176 .def_readwrite("liquid_saturation", &VariableArray::liquid_saturation)
177 .def_readwrite("mechanical_strain", &VariableArray::mechanical_strain)
178 .def_readwrite("molar_mass", &VariableArray::molar_mass)
179 .def_readwrite("molar_mass_derivative",
181 .def_readwrite("molar_fraction", &VariableArray::molar_fraction)
182 .def_readwrite("gas_phase_pressure", &VariableArray::gas_phase_pressure)
183 .def_readwrite("porosity", &VariableArray::porosity)
184 .def_readwrite("solid_grain_pressure",
186 .def_readwrite("stress", &VariableArray::stress)
187 .def_readwrite("temperature", &VariableArray::temperature)
188 .def_readwrite("total_strain", &VariableArray::total_strain)
189 .def_readwrite("total_stress", &VariableArray::total_stress)
190 .def_readwrite("transport_porosity", &VariableArray::transport_porosity)
191 .def_readwrite("vapour_pressure", &VariableArray::vapour_pressure)
192 .def_readwrite("volumetric_strain", &VariableArray::volumetric_strain);
193}
194
195void bindVariableEnum(py::module_& m)
196{
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)
224 .export_values();
225}
226
227void bindProperty(py::module_& m)
228{
229 py::class_<Property>(m, "Property", R"pbdoc(
230 Base class for material properties.
231 )pbdoc")
232 .def(
233 "value",
234 [](const Property& p, const VariableArray& va,
235 const SpatialPosition& pos, double t, double dt)
236 { return p.value(va, pos, t, dt); },
237 py::arg("variable_array"),
238 py::arg("pos"),
239 py::arg("t"),
240 py::arg("dt"),
241 R"pbdoc(
242 Evaluate the property value.
243
244 Parameters:
245 variable_array: Current time step values.
246 pos: Spatial position.
247 t: Current time.
248 dt: Time step size.
249 )pbdoc")
250
251 .def(
252 "value",
253 [](const Property& p, const VariableArray& va,
254 const VariableArray& vap, const SpatialPosition& pos, double t,
255 double dt) { return p.value(va, vap, pos, t, dt); },
256 py::arg("variable_array"),
257 py::arg("variable_array_prev"),
258 py::arg("pos"),
259 py::arg("t"),
260 py::arg("dt"),
261 R"pbdoc(
262 Evaluate the property value with previous time step data.
263
264 Parameters:
265 variable_array: Current time step values.
266 variable_array_prev: Previous time step values.
267 pos: Spatial position.
268 t: Current time.
269 dt: Time step size.
270 )pbdoc")
271
272 .def(
273 "dValue",
274 [](const Property& p, const VariableArray& va,
275 const VariableArray& vap, Variable var,
276 const SpatialPosition& pos, double t, double dt)
277 { return p.dValue(va, vap, var, pos, t, dt); },
278 py::arg("variable_array"),
279 py::arg("variable_array_prev"),
280 py::arg("variable"),
281 py::arg("pos"),
282 py::arg("t"),
283 py::arg("dt"),
284 R"pbdoc(
285 Evaluate the derivative of the property with respect to a variable,
286 using previous time step data.
287
288 Parameters:
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.
293 t: Current time.
294 dt: Time step size.
295 )pbdoc")
296
297 .def(
298 "dValue",
299 [](const Property& p, const VariableArray& va, Variable var,
300 const SpatialPosition& pos, double t, double dt)
301 { return p.dValue(va, var, pos, t, dt); },
302 py::arg("variable_array"),
303 py::arg("variable"),
304 py::arg("pos"),
305 py::arg("t"),
306 py::arg("dt"),
307 R"pbdoc(
308 Evaluate the derivative of the property with respect to a variable.
309
310 Parameters:
311 variable_array: Current time step values.
312 variable: Variable to differentiate with respect to.
313 pos: Spatial position.
314 t: Current time.
315 dt: Time step size.
316 )pbdoc");
317}
318
319void bindConstant(py::module_& m)
320{
321 py::class_<Constant, Property>(m, "Constant", R"pbdoc(
322 Constant property with a fixed value.
323 )pbdoc")
324 .def(py::init<std::string, PropertyDataType>(),
325 py::arg("name"),
326 py::arg("value"),
327 R"pbdoc(
328 Construct a Constant property.
329
330 Parameters:
331 name: Name of the property.
332 value: Constant value.
333 )pbdoc");
334}
335
336void bindLinear(py::module_& m)
337{
338 py::class_<Linear, Property>(m, "Linear", R"pbdoc(
339 Linearly interpolated property dependent on independent variables.
340 )pbdoc")
341 .def(py::init(
342 [](std::string name, PropertyDataType reference_value,
343 std::vector<std::tuple<Variable, VariableType,
344 VariableType>> const& ivs)
345 {
346 auto makeIV = [](std::tuple<Variable, VariableType,
347 VariableType> const& iv)
348 { return std::make_from_tuple<IndependentVariable>(iv); };
349
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);
355 }),
356 py::arg("name"),
357 py::arg("reference_value"),
358 py::arg("independent_variables"),
359 R"pbdoc(
360 Construct a Linear property.
361
362 Parameters:
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.
366 )pbdoc");
367}
368} // namespace
369
371{
372#ifndef NDEBUG
374#else // NDEBUG
376#endif // NDEBUG
377
378 m.attr("__name__") = "ogs.mpl";
379 m.doc() = "pybind11 ogs MPL bindings";
380
384 bindProperty(m);
385 bindConstant(m);
386 bindLinear(m);
387}
Definition of the Point3d class.
virtual PropertyDataType value() const
Definition Property.cpp:76
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:105
DeformationGradient deformation_gradient
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)
Definition Logging.cpp:64
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
Definition Property.h:31
void bindSpatialPosition(py::module_ &m)
Definition ogs_mpl.cpp:59
void bindConstant(py::module_ &m)
Definition ogs_mpl.cpp:246
void bindLinear(py::module_ &m)
Definition ogs_mpl.cpp:257
void bindVariableArray(py::module_ &m)
Definition ogs_mpl.cpp:119
void bindVariableEnum(py::module_ &m)
Definition ogs_mpl.cpp:159
std::optional< MathLib::Point3d > pythonToSpatialPositionCoords(py::object const &coordinates)
Definition ogs_mpl.cpp:35
py::object spatialPositionCoordsToPython(SpatialPosition const &pos)
Definition ogs_mpl.cpp:49
void bindProperty(py::module_ &m)
Definition ogs_mpl.cpp:191
PYBIND11_MODULE(mpl, m)
Definition ogs_mpl.cpp:370