Use pybind for subsidence bindings

lib/CMakeLists.txt

@@ -325,6 +325,34 @@ set_target_properties(
   _well_info PROPERTIES LIBRARY_OUTPUT_DIRECTORY
                         ${PROJECT_SOURCE_DIR}/python/resdata/well)
 
+# -----------------------------------------------------------------
+# Pybind11 extension module: resdata.gravimetry._subsidence
+# -----------------------------------------------------------------
+pybind11_add_module(_subsidence resdata/rd_subsidence_pybind.cpp
+                    resdata/cwrap_pybind.cpp)
+target_link_libraries(_subsidence PRIVATE resdata fmt::fmt)
+target_include_directories(
+  _subsidence
+  PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/include
+          ${CMAKE_CURRENT_BINARY_DIR}/include
+          ${CMAKE_CURRENT_SOURCE_DIR}/private-include)
+if(SKBUILD)
+  set(_subsidence_install_dir "python/resdata/gravimetry")
+  if(APPLE)
+    set_target_properties(_subsidence PROPERTIES INSTALL_RPATH
+                                                 "@loader_path/../.libs")
+  else()
+    set_target_properties(_subsidence PROPERTIES INSTALL_RPATH
+                                                 "$ORIGIN/../.libs")
+  endif()
+else()
+  set(_subsidence_install_dir "${CMAKE_INSTALL_LIBDIR}")
+endif()
+install(TARGETS _subsidence LIBRARY DESTINATION ${_subsidence_install_dir})
+set_target_properties(
+  _subsidence PROPERTIES LIBRARY_OUTPUT_DIRECTORY
+                         ${PROJECT_SOURCE_DIR}/python/resdata/gravimetry)
+
 # -----------------------------------------------------------------
 # Pybind11 extension module: resdata.summary._rd_sum
 # -----------------------------------------------------------------

lib/include/resdata/rd_subsidence.hpp

@@ -1,15 +1,13 @@
-#ifndef ERT_RD_SUBSIDENCE_H
-#define ERT_RD_SUBSIDENCE_H
+#pragma once
+
+#include <optional>
+#include <string>
 
 #include <resdata/rd_file.hpp>
 #include <resdata/rd_file_view.hpp>
 #include <resdata/rd_grid.hpp>
 #include <resdata/rd_region.hpp>
 
-#ifdef __cplusplus
-extern "C" {
-#endif
-
 typedef struct rd_subsidence_struct rd_subsidence_type;
 typedef struct rd_subsidence_survey_struct rd_subsidence_survey_type;
 
@@ -18,32 +16,28 @@ rd_subsidence_type *rd_subsidence_alloc(rd_grid_type *rd_grid,
                                         const rd_file_type *init_file);
 rd_subsidence_survey_type *
 rd_subsidence_add_survey_PRESSURE(rd_subsidence_type *subsidence,
-                                  const char *name,
+                                  const std::string &name,
                                   const rd_file_view_type *restart_view);
 
 bool rd_subsidence_has_survey(const rd_subsidence_type *subsidence,
-                              const char *name);
+                              const std::string &name);
 double rd_subsidence_eval(const rd_subsidence_type *subsidence,
-                          const char *base, const char *monitor,
+                          const std::string &base,
+                          const std::optional<std::string> &monitor,
                           rd_region_type *region, double utm_x, double utm_y,
                           double depth, double compressibility,
                           double poisson_ratio);
 
 double rd_subsidence_eval_geertsma(const rd_subsidence_type *subsidence,
-                                   const char *base, const char *monitor,
+                                   const std::string &base,
+                                   const std::optional<std::string> &monitor,
                                    rd_region_type *region, double utm_x,
                                    double utm_y, double depth,
                                    double youngs_modulus, double poisson_ratio,
                                    double seabed);
 
-double rd_subsidence_eval_geertsma_rporv(const rd_subsidence_type *subsidence,
-                                         const char *base, const char *monitor,
-                                         rd_region_type *region, double utm_x,
-                                         double utm_y, double depth,
-                                         double youngs_modulus,
-                                         double poisson_ratio, double seabed);
-
-#ifdef __cplusplus
-}
-#endif
-#endif
+double rd_subsidence_eval_geertsma_rporv(
+    const rd_subsidence_type *subsidence, const std::string &base,
+    const std::optional<std::string> &monitor, rd_region_type *region,
+    double utm_x, double utm_y, double depth, double youngs_modulus,
+    double poisson_ratio, double seabed);

lib/resdata/cwrap_pybind.cpp

@@ -12,6 +12,8 @@
 #include <resdata/FortIO.hpp>
 #include <resdata/rd_file_view.hpp>
 #include <resdata/rd_type.hpp>
+#include <resdata/rd_subsidence.hpp>
+#include <resdata/rd_region.hpp>
 
 #include <ert/util/double_vector.hpp>
 #include <ert/util/stringlist.hpp>
@@ -233,3 +235,35 @@ template <> rd_file_view_type *from_cwrap<rd_file_view_type>(py::handle obj) {
                              static_cast<std::string>(py::repr(obj)));
     return cast_cwrap<rd_file_view_type>(obj);
 }
+py::object ResdataSubsidence() {
+    static py::object cls;
+    if (!cls) {
+        cls =
+            py::module_::import("resdata.gravimetry").attr("ResdataSubsidence");
+    }
+    return cls;
+}
+
+template <> rd_subsidence_type *from_cwrap<rd_subsidence_type>(py::handle obj) {
+    if (!py::isinstance(obj, ResdataSubsidence()))
+        throw py::type_error("Expected ResdataSubsidence, got " +
+                             static_cast<std::string>(py::repr(obj)));
+
+    return cast_cwrap<rd_subsidence_type>(obj);
+}
+
+py::object ResdataRegion() {
+    static py::object cls;
+    if (!cls) {
+        cls = py::module_::import("resdata.grid").attr("ResdataRegion");
+    }
+    return cls;
+}
+
+template <> rd_region_type *from_cwrap<rd_region_type>(py::handle obj) {
+    if (!py::isinstance(obj, ResdataRegion()))
+        throw py::type_error("Expected ResdataRegion, got " +
+                             static_cast<std::string>(py::repr(obj)));
+
+    return cast_cwrap<rd_region_type>(obj);
+}

lib/resdata/rd_subsidence.cpp

@@ -8,6 +8,7 @@
 #include <utility>
 #include <vector>
 #include <optional>
+#include <fmt/format.h>
 
 #include <resdata/rd_kw.hpp>
 #include <resdata/rd_util.hpp>
@@ -59,7 +60,8 @@ struct rd_subsidence_survey_struct {
     std::optional<std::vector<double>>
         dynamic_porevolume; /* Porevolume in each grid cell at survey time */
 
-    rd_subsidence_survey_struct(const rd_subsidence_type &sub, const char *name)
+    rd_subsidence_survey_struct(const rd_subsidence_type &sub,
+                                const std::string name)
         : grid_cache(sub.grid_cache.get()),
           aquifer_cell(sub.aquifer_cell.get()), name(name),
           porv(sub.grid_cache->size()), pressure(sub.grid_cache->size()) {}
@@ -68,7 +70,7 @@ struct rd_subsidence_survey_struct {
 static std::unique_ptr<rd_subsidence_survey_type>
 rd_subsidence_survey_alloc_PRESSURE(rd_subsidence_type *rd_subsidence,
                                     const rd_file_view_type *restart_view,
-                                    const char *name) {
+                                    const std::string &name) {
 
     auto survey(
         std::make_unique<rd_subsidence_survey_struct>(*rd_subsidence, name));
@@ -169,17 +171,17 @@ static double rd_subsidence_survey_eval_geertsma_rporv(
     std::vector<double> weight(size);
 
     if (!base_survey->dynamic_porevolume.has_value()) {
-        util_abort(
-            "%s: Keyword RPORV not defined in .UNRST file for %s. Please add "
+        throw std::invalid_argument(fmt::format(
+            "Keyword RPORV not defined in .UNRST file for {}. Please add "
             "RPORV keyword to output in RPTRST clause in .DATA file.\n",
-            __func__, base_survey->name.c_str());
+            base_survey->name));
     }
 
     if (monitor_survey && !monitor_survey->dynamic_porevolume.has_value()) {
-        util_abort(
-            "%s: Keyword RPORV not defined in .UNRST file for %s. Please add "
-            "RPORV keyword to output in RPTRST clause in .DATA file.\n",
-            __func__, monitor_survey->name.c_str());
+        throw std::invalid_argument(fmt::format(
+            "Keyword RPORV not defined in .UNRST file for {}. Please add "
+            "RPORV keyword to output in RPTRST clause in .DATA file.",
+            monitor_survey->name));
     }
 
     for (size_t index = 0; index < weight.size(); ++index) {
@@ -216,10 +218,8 @@ rd_subsidence_type *rd_subsidence_alloc(rd_grid_type *rd_grid,
 
 rd_subsidence_survey_type *
 rd_subsidence_add_survey_PRESSURE(rd_subsidence_type *subsidence,
-                                  const char *name,
+                                  const std::string &name,
                                   const rd_file_view_type *restart_view) {
-    if (name == nullptr)
-        throw std::invalid_argument("Name cannot be NULL");
     auto survey =
         rd_subsidence_survey_alloc_PRESSURE(subsidence, restart_view, name);
     rd_subsidence_survey_type *ret = survey.get();
@@ -228,51 +228,50 @@ rd_subsidence_add_survey_PRESSURE(rd_subsidence_type *subsidence,
 }
 
 bool rd_subsidence_has_survey(const rd_subsidence_type *subsidence,
-                              const char *name) {
-    if (name == nullptr)
-        return false;
+                              const std::string &name) {
     return subsidence->surveys.count(name) > 0;
 }
 
 double rd_subsidence_eval(const rd_subsidence_type *subsidence,
-                          const char *base, const char *monitor,
+                          const std::string &base,
+                          const std::optional<std::string> &monitor,
                           rd_region_type *region, double utm_x, double utm_y,
                           double depth, double compressibility,
                           double poisson_ratio) {
     rd_subsidence_survey_type *base_survey = subsidence->surveys.at(base).get();
     rd_subsidence_survey_type *monitor_survey = nullptr;
     if (monitor)
-        monitor_survey = subsidence->surveys.at(monitor).get();
+        monitor_survey = subsidence->surveys.at(*monitor).get();
     return rd_subsidence_survey_eval(base_survey, monitor_survey, region, utm_x,
                                      utm_y, depth, compressibility,
                                      poisson_ratio);
 }
 
 double rd_subsidence_eval_geertsma(const rd_subsidence_type *subsidence,
-                                   const char *base, const char *monitor,
+                                   const std::string &base,
+                                   const std::optional<std::string> &monitor,
                                    rd_region_type *region, double utm_x,
                                    double utm_y, double depth,
                                    double youngs_modulus, double poisson_ratio,
                                    double seabed) {
     rd_subsidence_survey_type *base_survey = subsidence->surveys.at(base).get();
     rd_subsidence_survey_type *monitor_survey = nullptr;
     if (monitor)
-        monitor_survey = subsidence->surveys.at(monitor).get();
+        monitor_survey = subsidence->surveys.at(*monitor).get();
     return rd_subsidence_survey_eval_geertsma(
         base_survey, monitor_survey, region, utm_x, utm_y, depth,
         youngs_modulus, poisson_ratio, seabed);
 }
 
-double rd_subsidence_eval_geertsma_rporv(const rd_subsidence_type *subsidence,
-                                         const char *base, const char *monitor,
-                                         rd_region_type *region, double utm_x,
-                                         double utm_y, double depth,
-                                         double youngs_modulus,
-                                         double poisson_ratio, double seabed) {
+double rd_subsidence_eval_geertsma_rporv(
+    const rd_subsidence_type *subsidence, const std::string &base,
+    const std::optional<std::string> &monitor, rd_region_type *region,
+    double utm_x, double utm_y, double depth, double youngs_modulus,
+    double poisson_ratio, double seabed) {
     rd_subsidence_survey_type *base_survey = subsidence->surveys.at(base).get();
     rd_subsidence_survey_type *monitor_survey = nullptr;
     if (monitor)
-        monitor_survey = subsidence->surveys.at(monitor).get();
+        monitor_survey = subsidence->surveys.at(*monitor).get();
     return rd_subsidence_survey_eval_geertsma_rporv(
         base_survey, monitor_survey, region, utm_x, utm_y, depth,
         youngs_modulus, poisson_ratio, seabed);

lib/resdata/rd_subsidence_pybind.cpp

@@ -0,0 +1,93 @@
+#include <cstdint>
+#include <optional>
+#include <string>
+
+#include <pybind11/pybind11.h>
+#include <pybind11/stl.h>
+
+#include <resdata/rd_file.hpp>
+#include <resdata/rd_file_view.hpp>
+#include <resdata/rd_grid.hpp>
+#include <resdata/rd_region.hpp>
+#include <resdata/rd_subsidence.hpp>
+
+#include <detail/resdata/cwrap_pybind.hpp>
+
+namespace py = pybind11;
+
+PYBIND11_MODULE(_subsidence, m) {
+    register_exceptions(m);
+    m.doc() = "pybind11 bindings for rd_subsidence.cpp";
+
+    m.def("_alloc", [](py::handle grid, py::handle init_file) -> py::object {
+        auto *subsidence =
+            rd_subsidence_alloc(from_cwrap<rd_grid_type>(grid),
+                                from_cwrap<rd_file_type>(init_file));
+
+        if (subsidence == nullptr)
+            return py::none();
+
+        return py::cast(reinterpret_cast<std::uintptr_t>(subsidence));
+    });
+
+    m.def("_free", [](py::handle self) {
+        rd_subsidence_free(from_cwrap<rd_subsidence_type>(self));
+    });
+
+    m.def("_add_survey_PRESSURE", [](py::handle self, std::string survey_name,
+                                     py::handle restart_file_view) {
+        rd_subsidence_add_survey_PRESSURE(
+            from_cwrap<rd_subsidence_type>(self), survey_name,
+            from_cwrap<rd_file_view_type>(restart_file_view));
+    });
+
+    m.def("_eval", [](py::handle self, std::string base,
+                      std::optional<std::string> monitor, py::object region,
+                      double utm_x, double utm_y, double depth,
+                      double compressibility, double poisson_ratio) {
+        rd_region_type *region_ptr = nullptr;
+        if (!region.is_none()) {
+            region_ptr = from_cwrap<rd_region_type>(region);
+        }
+
+        return rd_subsidence_eval(from_cwrap<rd_subsidence_type>(self), base,
+                                  monitor, region_ptr, utm_x, utm_y, depth,
+                                  compressibility, poisson_ratio);
+    });
+
+    m.def("_eval_geertsma", [](py::handle self, std::string base,
+                               std::optional<std::string> monitor,
+                               py::object region, double utm_x, double utm_y,
+                               double depth, double youngs_modulus,
+                               double poisson_ratio, double seabed) {
+        rd_region_type *region_ptr = nullptr;
+        if (!region.is_none()) {
+            region_ptr = from_cwrap<rd_region_type>(region);
+        }
+
+        return rd_subsidence_eval_geertsma(
+            from_cwrap<rd_subsidence_type>(self), base, monitor, region_ptr,
+            utm_x, utm_y, depth, youngs_modulus, poisson_ratio, seabed);
+    });
+
+    m.def("_eval_geertsma_rporv",
+          [](py::handle self, std::string base,
+             std::optional<std::string> monitor, py::object region,
+             double utm_x, double utm_y, double depth, double youngs_modulus,
+             double poisson_ratio, double seabed) {
+              rd_region_type *region_ptr = nullptr;
+              if (!region.is_none()) {
+                  region_ptr = from_cwrap<rd_region_type>(region);
+              }
+
+              return rd_subsidence_eval_geertsma_rporv(
+                  from_cwrap<rd_subsidence_type>(self), base, monitor,
+                  region_ptr, utm_x, utm_y, depth, youngs_modulus,
+                  poisson_ratio, seabed);
+          });
+
+    m.def("_has_survey", [](py::handle self, std::string survey_name) {
+        return rd_subsidence_has_survey(from_cwrap<rd_subsidence_type>(self),
+                                        survey_name);
+    });
+}