Module for Generalied Eigen Solver

CMakeLists.txt

@@ -49,11 +49,19 @@ foreach(dependency_i ${DEPENDENCIES})
     include(get_${dependency_i})
 endforeach()
 
+cmaize_find_or_build_dependency(
+    eigen
+    URL https://www.gitlab.com/libeigen/eigen
+    VERSION 3.4.0
+    BUILD_TARGET eigen
+    FIND_TARGET Eigen3::Eigen
+)
+
 cmaize_add_library(
     scf
     SOURCE_DIR "${CMAKE_CURRENT_LIST_DIR}/src/scf"
     INCLUDE_DIRS "${CMAKE_CURRENT_LIST_DIR}/include/scf"
-    DEPENDS "${DEPENDENCIES}"
+    DEPENDS "${DEPENDENCIES}" eigen
 )
 
 if("${BUILD_TAMM_SCF}")

src/scf/eigen_solver/eigen_generalized.cpp

@@ -0,0 +1,88 @@
+/*
+ * Copyright 2024 NWChemEx-Project
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include "eigen_solver.hpp"
+#include <Eigen/Eigen>
+#include <simde/simde.hpp>
+
+namespace scf::eigen_solver {
+
+using pt = simde::GeneralizedEigenSolve;
+
+const auto desc = R"(
+Generalized Eigen Solve via Eigen
+---------------------------------
+
+TODO: Write me!!!
+)";
+
+MODULE_CTOR(EigenGeneralized) {
+    description(desc);
+    satisfies_property_type<pt>();
+}
+
+MODULE_RUN(EigenGeneralized) {
+    auto&& [A, B] = pt::unwrap_inputs(inputs);
+
+    // Convert to Eigen buffers
+    using matrix_alloc_t = tensorwrapper::allocator::Eigen<double, 2>;
+    using vector_alloc_t = tensorwrapper::allocator::Eigen<double, 1>;
+    const auto& eigen_A  = matrix_alloc_t::rebind(A.buffer());
+    const auto& eigen_B  = matrix_alloc_t::rebind(B.buffer());
+
+    // Wrap the tensors in Eigen::Map objects to avoid copy
+    const auto* pA      = eigen_A.value().data();
+    const auto* pB      = eigen_B.value().data();
+    const auto& shape_A = eigen_A.layout().shape().as_smooth();
+    auto rows           = shape_A.extent(0);
+    auto cols           = shape_A.extent(1);
+    Eigen::Map<const Eigen::MatrixXd> A_map(pA, rows, cols);
+    Eigen::Map<const Eigen::MatrixXd> B_map(pB, rows, cols);
+
+    // Compute
+    Eigen::GeneralizedEigenSolver<Eigen::MatrixXd> ges(A_map, B_map);
+    auto eigen_values  = ges.eigenvalues();
+    auto eigen_vectors = ges.eigenvectors();
+
+    // Wrap in TensorWrapper Tensor
+    tensorwrapper::shape::Smooth vector_shape{rows};
+    tensorwrapper::shape::Smooth matrix_shape{rows, cols};
+    tensorwrapper::layout::Physical vector_layout(vector_shape);
+    tensorwrapper::layout::Physical matrix_layout(matrix_shape);
+
+    using matrix_buffer = typename matrix_alloc_t::eigen_buffer_type;
+    using vector_buffer = typename vector_alloc_t::eigen_buffer_type;
+
+    typename vector_buffer::data_type vector_tensor(rows);
+    typename matrix_buffer::data_type matrix_tensor(rows, cols);
+    for(auto i = 0; i < rows; ++i) {
+        vector_tensor(i) = eigen_values(i).real();
+        for(auto j = 0; j < cols; ++j) {
+            matrix_tensor(i, j) = eigen_vectors(i, j).real();
+        }
+    }
+
+    vector_buffer values_buffer(vector_tensor, vector_layout);
+    matrix_buffer vectors_buffer(matrix_tensor, matrix_layout);
+
+    simde::type::tensor values(vector_shape, values_buffer);
+    simde::type::tensor vectors(matrix_shape, vectors_buffer);
+
+    auto rv = results();
+    return pt::wrap_results(rv, values, vectors);
+}
+
+} // namespace scf::eigen_solver

src/scf/eigen_solver/eigen_solver.hpp

@@ -0,0 +1,28 @@
+/*
+ * Copyright 2024 NWChemEx-Project
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#pragma once
+#include <simde/simde.hpp>
+
+namespace scf::eigen_solver {
+
+DECLARE_MODULE(EigenGeneralized);
+
+inline void load_modules(pluginplay::ModuleManager& mm) {
+    mm.add_module<EigenGeneralized>("Generalized eigensolve via Eigen");
+}
+
+} // namespace scf::eigen_solver

src/scf/scf_mm.cpp

@@ -14,13 +14,15 @@
  * limitations under the License.
  */
 
+#include "eigen_solver/eigen_solver.hpp"
 #include "fock_operator/fock_operator.hpp"
 #include "scf_modules.hpp"
 #include <scf/scf_mm.hpp>
 
 namespace scf {
 
 void load_modules(pluginplay::ModuleManager& mm) {
+    eigen_solver::load_modules(mm);
     fock_operator::load_modules(mm);
     mm.add_module<CoulombsLaw>("Coulomb's Law");
 #ifdef BUILD_TAMM_SCF

tests/cxx/test_scf.hpp

@@ -17,6 +17,7 @@
 #include <catch2/catch_approx.hpp>
 #include <catch2/catch_template_test_macros.hpp>
 #include <catch2/catch_test_macros.hpp>
+#include <catch2/matchers/catch_matchers_floating_point.hpp>
 #include <simde/simde.hpp>
 
 namespace test_scf {

tests/cxx/unit_tests/eigen_solver/eigen_generalized.cpp

@@ -0,0 +1,40 @@
+/*
+ * Copyright 2024 NWChemEx-Project
+ *
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+#include "../../test_scf.hpp"
+#include <scf/scf.hpp>
+#include <simde/simde.hpp>
+
+TEST_CASE("EigenGeneralized") {
+    pluginplay::ModuleManager mm;
+    scf::load_modules(mm);
+
+    using pt = simde::GeneralizedEigenSolve;
+
+    auto& mod = mm.at("Generalized eigensolve via Eigen");
+
+    simde::type::tensor A({{1.0, 2.0}, {2.0, 3.0}});
+    simde::type::tensor B({{1.0, 0.0}, {0.0, 1.0}});
+
+    auto&& [values, vector] = mod.run_as<pt>(A, B);
+
+    using value_alloc_t      = tensorwrapper::allocator::Eigen<double, 1>;
+    const auto& eigen_values = value_alloc_t::rebind(values.buffer());
+
+    using Catch::Matchers::WithinAbs;
+    REQUIRE_THAT(eigen_values.value()(0), WithinAbs(-0.236068, 1E-6));
+    REQUIRE_THAT(eigen_values.value()(1), WithinAbs(4.236068, 1E-6));
+}