Update bindings for new libROM version

CMakeLists.txt

@@ -5,6 +5,8 @@ project(_pylibROM)
 set(CMAKE_BUILD_TYPE Debug)
 set(PYBIND11_FINDPYTHON ON)
 
+option (USE_MFEM "Build pylibROM with MFEM" OFF)
+
 #=================== ScaLAPACK (optional) ==================
 option(BUILD_SCALAPACK "Build static ScaLAPACK for libROM" OFF)
 
@@ -50,14 +52,18 @@ if (BUILD_LIBROM)
 	# )
 	# add_custom_target(RUN_LIBROM_BUILD ALL DEPENDS LIBROM_BUILD)
 
-	ExternalProject_Add(
-        	libROM
-        	SOURCE_DIR ${LIBROM_SCRIPTS_DIR}
-        	CONFIGURE_COMMAND ""
-        	BINARY_DIR ${LIBROM_DIR}
-        	BUILD_COMMAND ${LIBROM_SCRIPTS_DIR}/compile.sh -m -g -t ${LIBROM_DIR}/cmake/toolchains/simple.cmake
-        	INSTALL_COMMAND ""
-	)	
+    set(LIBROM_BUILD_CMD "${LIBROM_SCRIPTS_DIR}/compile.sh -t ${LIBROM_DIR}/cmake/toolchains/simple.cmake" CACHE STRING  "Command used to build libROM and dependencies")
+    if (USE_MFEM)
+        set(LIBROM_BUILD_CMD "${LIBROM_BUILD_CMD} -m -g")
+    endif()
+	# ExternalProject_Add(
+    #     	libROM
+    #     	SOURCE_DIR ${LIBROM_SCRIPTS_DIR}
+    #     	CONFIGURE_COMMAND ""
+    #     	BINARY_DIR ${LIBROM_DIR}
+    #         BUILD_COMMAND ${LIBROM_BUILD_CMD}
+    #         INSTALL_COMMAND ""
+    # )
 	message("Building libROM dependency...")
 endif(BUILD_LIBROM)
 
@@ -72,63 +78,58 @@ execute_process(COMMAND python3 -c "import mpi4py; print(mpi4py.get_include())"
 # # TODO(kevin): We do not bind mfem-related functions until we figure out how to type-cast SWIG Object.
 # #              Until then, mfem-related functions need to be re-implemented on python-end, using PyMFEM.
 
-find_library(MFEM mfem
-	"$ENV{MFEM_DIR}/lib"
-	"$ENV{MFEM_DIR}"
-	"${LIBROM_DIR}/dependencies/mfem")
-find_library(HYPRE HYPRE
-	"$ENV{HYPRE_DIR}/lib"
-	"${LIBROM_DIR}/dependencies/hypre/src/hypre/lib")
-find_library(PARMETIS parmetis
-	"$ENV{PARMETIS_DIR}/lib"
-	"$ENV{PARMETIS_DIR}/build/lib/libparmetis"
-	"${LIBROM_DIR}/dependencies/parmetis-4.0.3/build/lib/libparmetis")
-find_library(METIS metis
-	"$ENV{METIS_DIR}/lib"
-	"$ENV{PARMETIS_DIR}/build/lib/libmetis"
-	"${LIBROM_DIR}/dependencies/parmetis-4.0.3/build/lib/libmetis")
-find_path(MFEM_INCLUDES mfem.hpp
-	"$ENV{MFEM_DIR}/include"
-	"$ENV{MFEM_DIR}"
-	"${LIBROM_DIR}/dependencies/mfem")
-find_path(HYPRE_INCLUDES HYPRE.h
-	"$ENV{HYPRE_DIR}/include"
-	"${LIBROM_DIR}/dependencies/hypre/src/hypre/include")
-find_path(PARMETIS_INCLUDES metis.h
-	"$ENV{PARMETIS_DIR}/metis/include"
-	"${LIBROM_DIR}/dependencies/parmetis-4.0.3/metis/include")
-
+if (USE_MFEM)
+    find_library(MFEM mfem
+        "$ENV{MFEM_DIR}/lib"
+        "$ENV{MFEM_DIR}"
+        "${LIBROM_DIR}/dependencies/mfem")
+    find_library(HYPRE HYPRE
+        "$ENV{HYPRE_DIR}/lib"
+        "${LIBROM_DIR}/dependencies/hypre/src/hypre/lib")
+    find_library(PARMETIS parmetis
+        "$ENV{PARMETIS_DIR}/lib"
+        "$ENV{PARMETIS_DIR}/build/lib/libparmetis"
+        "${LIBROM_DIR}/dependencies/parmetis-4.0.3/build/lib/libparmetis")
+    find_library(METIS metis
+        "$ENV{METIS_DIR}/lib"
+        "$ENV{PARMETIS_DIR}/build/lib/libmetis"
+        "${LIBROM_DIR}/dependencies/parmetis-4.0.3/build/lib/libmetis")
+    find_path(MFEM_INCLUDES mfem.hpp
+        "$ENV{MFEM_DIR}/include"
+        "$ENV{MFEM_DIR}"
+        "${LIBROM_DIR}/dependencies/mfem")
+    find_path(HYPRE_INCLUDES HYPRE.h
+        "$ENV{HYPRE_DIR}/include"
+        "${LIBROM_DIR}/dependencies/hypre/src/hypre/include")
+    find_path(PARMETIS_INCLUDES metis.h
+        "$ENV{PARMETIS_DIR}/metis/include"
+        "${LIBROM_DIR}/dependencies/parmetis-4.0.3/metis/include")
+    set(PYLIBROM_HAS_MFEM 1)
+endif()
 #===================== pylibROM =============================
 
 
 set(CMAKE_CXX_STANDARD 14)
 
 find_package(MPI REQUIRED)
 
-set(SOURCE_DIR "bindings/pylibROM") 
+configure_file(${CMAKE_CURRENT_SOURCE_DIR}/bindings/pylibROM/pylibROM_config.h.in
+               ${CMAKE_CURRENT_SOURCE_DIR}/bindings/pylibROM/pylibROM_config.h)
+
+set(SOURCE_DIR "bindings/pylibROM")
 include_directories(
 	${SOURCE_DIR} 
 	${LIBROM_INCLUDE_DIR} 
 	${MPI_INCLUDE_PATH}
 	${MPI4PY}
 	${HDF5_C_INCLUDE_DIRS}
-	${MFEM_INCLUDES}
-	${HYPRE_INCLUDES}
-	${PARMETIS_INCLUDES}
-	${MFEM_C_INCLUDE_DIRS}
-)
-link_libraries(
-	${HDF5_LIBRARIES}
-	${MFEM}
-	${HYPRE}
-	${PARMETIS}
-	${METIS}
 )
+link_libraries(${HDF5_LIBRARIES})
 
 add_subdirectory("extern/pybind11")
 
-pybind11_add_module(_pylibROM
-    bindings/pylibROM/pylibROM.cpp
+set(PYLIBROM_SOURCES
+	bindings/pylibROM/pylibROM.cpp
 
 	bindings/pylibROM/linalg/pyMatrix.cpp
 	bindings/pylibROM/linalg/pyVector.cpp
@@ -163,14 +164,38 @@ pybind11_add_module(_pylibROM
 	bindings/pylibROM/utils/pyDatabase.hpp
 	bindings/pylibROM/utils/pyDatabase.cpp
 	bindings/pylibROM/utils/pyHDFDatabase.cpp
+	bindings/pylibROM/utils/pyHDFDatabaseMPIO.cpp
 	bindings/pylibROM/utils/pyCSVDatabase.cpp
 
-	bindings/pylibROM/mfem/pyUtilities.cpp
-	bindings/pylibROM/mfem/pyPointwiseSnapshot.cpp
-	bindings/pylibROM/mfem/pySampleMesh.cpp
-
 	bindings/pylibROM/python_utils/cpp_utils.hpp
 )
+
+if (USE_MFEM)
+    set(PYLIBROM_SOURCES ${PYLIBROM_SOURCES}
+        bindings/pylibROM/mfem/pyUtilities.cpp
+        bindings/pylibROM/mfem/pyPointwiseSnapshot.cpp
+        bindings/pylibROM/mfem/pySampleMesh.cpp)
+endif()
+
+pybind11_add_module(_pylibROM ${PYLIBROM_SOURCES})
 message("building pylibROM...")
 
+if (USE_MFEM)
+    target_include_directories(
+        _pylibROM
+        PUBLIC
+        ${MFEM_INCLUDES}
+        ${HYPRE_INCLUDES}
+        ${PARMETIS_INCLUDES}
+        ${MFEM_C_INCLUDE_DIRS})
+
+    target_link_libraries(
+        _pylibROM
+        PUBLIC
+        ${MFEM}
+        ${HYPRE}
+        ${PARMETIS}
+        ${METIS})
+endif()
+
 target_link_libraries(_pylibROM PRIVATE ROM)

bindings/pylibROM/__init__.py

@@ -6,4 +6,14 @@
 # either define/import the python routine in this file.
 # This will combine both c++ bindings/pure python routines into this module.
 
-from _pylibROM import *
+from _pylibROM.algo import *
+from _pylibROM.hyperreduction import *
+from _pylibROM.linalg import *
+
+try:
+    import _pylibROM.mfem
+    from _pylibROM.mfem import *
+except:
+    pass
+
+from _pylibROM.utils import *

bindings/pylibROM/linalg/pyBasisGenerator.cpp

@@ -20,9 +20,9 @@ void init_BasisGenerator(pybind11::module_ &m) {
         )
         .def("isNextSample", (bool (BasisGenerator::*)(double)) &BasisGenerator::isNextSample)
         .def("updateRightSV", (bool (BasisGenerator::*)()) &BasisGenerator::updateRightSV)
-        .def("takeSample", [](BasisGenerator& self, py::array_t<double> &u_in, double time, double dt, bool add_without_increase = false) {
-            return self.takeSample(getVectorPointer(u_in), time, dt, add_without_increase);
-        }, py::arg("u_in"), py::arg("time"), py::arg("dt"), py::arg("add_without_increase") = false)
+        .def("takeSample", [](BasisGenerator& self, py::array_t<double> &u_in, bool add_without_increase = false) {
+            return self.takeSample(getVectorPointer(u_in), add_without_increase);
+        }, py::arg("u_in"), py::arg("add_without_increase") = false)
         .def("endSamples", &BasisGenerator::endSamples, py::arg("kind") = "basis")
         .def("writeSnapshot", (void (BasisGenerator::*)()) &BasisGenerator::writeSnapshot)
         .def("loadSamples", (void (BasisGenerator::*)(const std::string&, const std::string&, int, Database::formats)) &BasisGenerator::loadSamples,
@@ -39,8 +39,6 @@ void init_BasisGenerator(pybind11::module_ &m) {
         .def("getTemporalBasis", (const Matrix* (BasisGenerator::*)()) &BasisGenerator::getTemporalBasis,py::return_value_policy::reference)
         .def("getSingularValues", (const Vector* (BasisGenerator::*)()) &BasisGenerator::getSingularValues,py::return_value_policy::reference)
         .def("getSnapshotMatrix", (const Matrix* (BasisGenerator::*)()) &BasisGenerator::getSnapshotMatrix,py::return_value_policy::reference)
-        .def("getNumBasisTimeIntervals", (int (BasisGenerator::*)() const) &BasisGenerator::getNumBasisTimeIntervals)
-        .def("getBasisIntervalStartTime", (double (BasisGenerator::*)(int) const) &BasisGenerator::getBasisIntervalStartTime, py::arg("which_interval"))
         .def("getNumSamples",(int (BasisGenerator::*)() const) &BasisGenerator::getNumSamples)
         .def("__del__", [](BasisGenerator& self) { self.~BasisGenerator(); }); // Destructor
 

bindings/pylibROM/linalg/pyBasisReader.cpp

@@ -12,54 +12,38 @@ using namespace CAROM;
 
 void init_BasisReader(pybind11::module_ &m) {
      py::class_<BasisReader>(m, "BasisReader")
-        .def(py::init<const std::string&, Database::formats>(),
+        .def(py::init<const std::string&, Database::formats, const int>(),
             py::arg("base_file_name"),
-            py::arg("file_format") = Database::formats::HDF5
+            py::arg("file_format") = Database::formats::HDF5,
+            py::arg("dim") = -1
         )
-        .def("isNewBasis",(bool (BasisReader::*)(double)) &BasisReader::isNewBasis,
-            py::arg("time"))
-        .def("getSpatialBasis",(Matrix* (BasisReader::*)(double)) &BasisReader::getSpatialBasis,
-            py::arg("time"))
-        .def("getSpatialBasis",(Matrix* (BasisReader::*)(double,int)) &BasisReader::getSpatialBasis,
-            py::arg("time"),
+        .def("getSpatialBasis",(Matrix* (BasisReader::*)()) &BasisReader::getSpatialBasis)
+        .def("getSpatialBasis",(Matrix* (BasisReader::*)(int)) &BasisReader::getSpatialBasis,
             py::arg("n"))
-        .def("getSpatialBasis",(Matrix* (BasisReader::*)(double,int,int)) &BasisReader::getSpatialBasis,
-            py::arg("time"),
+        .def("getSpatialBasis",(Matrix* (BasisReader::*)(int,int)) &BasisReader::getSpatialBasis,
             py::arg("start_col"),
             py::arg("end_col"))
-        .def("getSpatialBasis",(Matrix* (BasisReader::*)(double,double)) &BasisReader::getSpatialBasis,
-            py::arg("time"),
+        .def("getSpatialBasis",(Matrix* (BasisReader::*)(double)) &BasisReader::getSpatialBasis,
             py::arg("ef").noconvert())
-        .def("getTemporalBasis",(Matrix* (BasisReader::*)(double)) &BasisReader::getTemporalBasis,
-            py::arg("time"))
-        .def("getTemporalBasis",(Matrix* (BasisReader::*)(double,int)) &BasisReader::getTemporalBasis,
-            py::arg("time"),
+        .def("getTemporalBasis",(Matrix* (BasisReader::*)()) &BasisReader::getTemporalBasis)
+        .def("getTemporalBasis",(Matrix* (BasisReader::*)(int)) &BasisReader::getTemporalBasis,
             py::arg("n"))
-        .def("getTemporalBasis",(Matrix* (BasisReader::*)(double,int,int)) &BasisReader::getTemporalBasis,
-            py::arg("time"),
+        .def("getTemporalBasis",(Matrix* (BasisReader::*)(int,int)) &BasisReader::getTemporalBasis,
             py::arg("start_col"),
             py::arg("end_col"))
-        .def("getTemporalBasis",(Matrix* (BasisReader::*)(double,double)) &BasisReader::getTemporalBasis,
-            py::arg("time"),
+        .def("getTemporalBasis",(Matrix* (BasisReader::*)(double)) &BasisReader::getTemporalBasis,
             py::arg("ef").noconvert())
+        .def("getSingularValues",(Vector* (BasisReader::*)()) &BasisReader::getSingularValues)
         .def("getSingularValues",(Vector* (BasisReader::*)(double)) &BasisReader::getSingularValues,
-            py::arg("time"))
-        .def("getSingularValues",(Vector* (BasisReader::*)(double,double)) &BasisReader::getSingularValues,
-            py::arg("time"),
             py::arg("ef"))
-        .def("getDim", (int (BasisReader::*)(const std::string,double)) &BasisReader::getDim,
-            py::arg("kind"),
-            py::arg("time"))
-        .def("getNumSamples", (int (BasisReader::*)(const std::string,double)) &BasisReader::getNumSamples,
-            py::arg("kind"),
-            py::arg("time"))
-        .def("getSnapshotMatrix",(Matrix* (BasisReader::*)(double)) &BasisReader::getSnapshotMatrix,
-            py::arg("time"))
-        .def("getSnapshotMatrix",(Matrix* (BasisReader::*)(double,int)) &BasisReader::getSnapshotMatrix,
-            py::arg("time"),
+        .def("getDim", (int (BasisReader::*)(const std::string)) &BasisReader::getDim,
+            py::arg("kind"))
+        .def("getNumSamples", (int (BasisReader::*)(const std::string)) &BasisReader::getNumSamples,
+            py::arg("kind"))
+        .def("getSnapshotMatrix",(Matrix* (BasisReader::*)()) &BasisReader::getSnapshotMatrix)
+        .def("getSnapshotMatrix",(Matrix* (BasisReader::*)(int)) &BasisReader::getSnapshotMatrix,
             py::arg("n"))
-        .def("getSnapshotMatrix",(Matrix* (BasisReader::*)(double,int,int)) &BasisReader::getSnapshotMatrix,
-            py::arg("time"),
+        .def("getSnapshotMatrix",(Matrix* (BasisReader::*)(int,int)) &BasisReader::getSnapshotMatrix,
             py::arg("start_col"),
             py::arg("end_col"))
         .def("__del__", [](BasisReader& self) { self.~BasisReader(); }); // Destructor

bindings/pylibROM/linalg/pyMatrix.cpp

@@ -177,7 +177,11 @@ void init_matrix(pybind11::module_ &m) {
 
         .def("transposePseudoinverse",(void (Matrix::*)()) &Matrix::transposePseudoinverse)
 
-        .def("qr_factorize",(Matrix* (Matrix::*)() const) &Matrix::qr_factorize,py::return_value_policy::take_ownership)
+        .def("qr_factorize", [](const Matrix& self) -> std::vector<std::unique_ptr<Matrix>> {
+            std::vector<std::unique_ptr<Matrix>> qr;
+            self.qr_factorize(qr);
+            return qr;
+        })
 
         // TODO (kevin): due to the difference between python and c++, technically we should not take
         //               row_pivot and row_pivot_owner as input parameters, just returning them in the end as outputs.
@@ -189,7 +193,8 @@ void init_matrix(pybind11::module_ &m) {
             return std::make_tuple(row_pivot, row_pivot_owner);
         })
 
-        .def("orthogonalize", (void (Matrix::*)()) &Matrix::orthogonalize)
+        .def("orthogonalize", (void (Matrix::*)(bool, double)) &Matrix::orthogonalize, py::arg("double_pass") = false, py::arg("zero_tol") = 1.0e-15)
+        .def("orthogonalize_last", (void (Matrix::*)(int, bool, double)) &Matrix::orthogonalize_last, py::arg("ncols") = -1, py::arg("double_pass") = false, py::arg("zero_tol") = 1.0e-15)
 
         .def("item", (const double& (Matrix::*)(int, int) const) &Matrix::item)
         .def("__getitem__", [](Matrix& self, int row, int col) { 

bindings/pylibROM/linalg/pyOptions.cpp

@@ -10,10 +10,9 @@ using namespace CAROM;
 void init_Options(pybind11::module_ &m) {
 
     py::class_<Options>(m, "Options")
-        .def(py::init<int, int, int, bool, bool>(), py::arg("dim_"), py::arg("samples_per_time_interval_"),py::arg("max_time_intervals_") = -1,py::arg("update_right_SV_") = false, py::arg("write_snapshots_") = false)
+        .def(py::init<int, int, bool, bool>(), py::arg("dim_"), py::arg("max_num_samples_"),py::arg("update_right_SV_") = false, py::arg("write_snapshots_") = false)
         .def_readwrite("dim", &Options::dim)
-        .def_readwrite("samples_per_time_interval", &Options::samples_per_time_interval)
-        .def_readwrite("max_time_intervals", &Options::max_time_intervals)
+        .def_readwrite("max_num_samples", &Options::max_num_samples)
         .def_readwrite("update_right_SV", &Options::update_right_SV)
         .def_readwrite("write_snapshots", &Options::write_snapshots)
         .def_readwrite("max_basis_dimension", &Options::max_basis_dimension)
@@ -33,6 +32,7 @@ void init_Options(pybind11::module_ &m) {
         .def_readwrite("min_sampling_time_step_scale", &Options::min_sampling_time_step_scale)
         .def_readwrite("sampling_time_step_scale", &Options::sampling_time_step_scale)
         .def_readwrite("max_sampling_time_step_scale", &Options::max_sampling_time_step_scale)
+        .def_readwrite("static_svd_preserve_snapshot", &Options::static_svd_preserve_snapshot)
         .def("setMaxBasisDimension", &Options::setMaxBasisDimension, py::arg("max_basis_dimension_"))
         .def("setSingularValueTol", &Options::setSingularValueTol, py::arg("singular_value_tol_"))
         .def("setDebugMode", &Options::setDebugMode,  py::arg("debug_algorithm_"))

bindings/pylibROM/linalg/svd/pyIncrementalSVD.cpp

@@ -31,8 +31,8 @@ class PyIncrementalSVD : public IncrementalSVD {
         PYBIND11_OVERRIDE(const Matrix*, IncrementalSVD, getSnapshotMatrix );
     }
 
-    bool takeSample(double* u_in, double time, bool add_without_increase) override {
-        PYBIND11_OVERLOAD(bool, IncrementalSVD, takeSample, u_in, time, add_without_increase);
+    bool takeSample(double* u_in, bool add_without_increase) override {
+        PYBIND11_OVERLOAD(bool, IncrementalSVD, takeSample, u_in, add_without_increase);
     }
 
     ~PyIncrementalSVD() override {
@@ -41,8 +41,8 @@ class PyIncrementalSVD : public IncrementalSVD {
 
 protected:
 
-    void buildInitialSVD(double* u, double time) override {
-        PYBIND11_OVERRIDE_PURE(void, IncrementalSVD, buildInitialSVD, u, time);
+    void buildInitialSVD(double* u) override {
+        PYBIND11_OVERRIDE_PURE(void, IncrementalSVD, buildInitialSVD, u);
     }
 
     void computeBasis() override {
@@ -62,18 +62,15 @@ class PyIncrementalSVD : public IncrementalSVD {
 void init_IncrementalSVD(pybind11::module_ &m) {
     py::class_<IncrementalSVD,PyIncrementalSVD>(m, "IncrementalSVD")
         .def(py::init<Options, const std::string&>())
-        .def("takeSample", [](IncrementalSVD& self, py::array_t<double> &u_in, double time,bool add_without_increase = false) {
-            bool result = self.takeSample(getVectorPointer(u_in), time, add_without_increase);
+        .def("takeSample", [](IncrementalSVD& self, py::array_t<double> &u_in, bool add_without_increase = false) {
+            bool result = self.takeSample(getVectorPointer(u_in), add_without_increase);
             return result;
-        }, py::arg("u_in"), py::arg("time"),py::arg("add_without_increase") = false)
+        }, py::arg("u_in"), py::arg("add_without_increase") = false)
         .def("getSpatialBasis", (const Matrix* (IncrementalSVD::*)()) &IncrementalSVD::getSpatialBasis)
         .def("getTemporalBasis", (const Matrix* (IncrementalSVD::*)()) &IncrementalSVD::getTemporalBasis)
         .def("getSingularValues", (const Vector* (IncrementalSVD::*)()) &IncrementalSVD::getSingularValues)
         .def("getSnapshotMatrix", (const Matrix* (IncrementalSVD::*)()) &IncrementalSVD::getSnapshotMatrix)
         .def("getDim", (int (IncrementalSVD::*)() const) &IncrementalSVD::getDim)
-        .def("getNumBasisTimeIntervals", (int (IncrementalSVD::*)() const) &IncrementalSVD::getNumBasisTimeIntervals)
-        .def("getBasisIntervalStartTime", (double (IncrementalSVD::*)(int) const) &IncrementalSVD::getBasisIntervalStartTime)
-        .def("isNewTimeInterval", (bool (IncrementalSVD::*)() const) &IncrementalSVD::isNewTimeInterval)
-        .def("increaseTimeInterval", (void (IncrementalSVD::*)()) &IncrementalSVD::increaseTimeInterval)
+        .def("getMaxNumSamples", (int (IncrementalSVD::*)() const) &IncrementalSVD::getMaxNumSamples)
         .def("getNumSamples", (int (IncrementalSVD::*)() const) &IncrementalSVD::getNumSamples);
 }