diff --git a/src/integrals/ao_integrals/ao_integrals.hpp b/src/integrals/ao_integrals/ao_integrals.hpp index f7b1a8a5..029cda17 100644 --- a/src/integrals/ao_integrals/ao_integrals.hpp +++ b/src/integrals/ao_integrals/ao_integrals.hpp @@ -22,17 +22,32 @@ namespace integrals::ao_integrals { DECLARE_MODULE(AOIntegralsDriver); DECLARE_MODULE(JFourCenter); DECLARE_MODULE(KFourCenter); +DECLARE_MODULE(JDensityFitted); +DECLARE_MODULE(KDensityFitted); +DECLARE_MODULE(DFIntegral); +DECLARE_MODULE(CoulombMetric); inline void set_defaults(pluginplay::ModuleManager& mm) { - const auto ao_driver = "AO integral driver"; - mm.change_submod(ao_driver, "Coulomb matrix", "Four center J builder"); - mm.change_submod(ao_driver, "Exchange matrix", "Four center K builder"); + mm.change_submod("AO integral driver", "Coulomb matrix", + "Four center J builder"); + mm.change_submod("AO integral driver", "Exchange matrix", + "Four center K builder"); + mm.change_submod("Density Fitted J builder", "DF ERI", + "Density Fitting Integral"); + mm.change_submod("Density Fitted K builder", "DF ERI", + "Density Fitting Integral"); + mm.change_submod("Density Fitting Integral", "Coulomb Metric", + "Coulomb Metric"); } inline void load_modules(pluginplay::ModuleManager& mm) { mm.add_module("AO integral driver"); mm.add_module("Four center J builder"); mm.add_module("Four center K builder"); + mm.add_module("Density Fitted J builder"); + mm.add_module("Density Fitted K builder"); + mm.add_module("Density Fitting Integral"); + mm.add_module("Coulomb Metric"); set_defaults(mm); } diff --git a/src/integrals/ao_integrals/coulomb_metric.cpp b/src/integrals/ao_integrals/coulomb_metric.cpp new file mode 100644 index 00000000..f565fcea --- /dev/null +++ b/src/integrals/ao_integrals/coulomb_metric.cpp @@ -0,0 +1,74 @@ +/* + * Copyright 2025 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 "ao_integrals.hpp" +#include + +namespace integrals::ao_integrals { + +using pt = simde::ERI2; + +namespace { + +auto desc = R"( +Inverse Coulomb Metric +--------------------- +)"; + +} + +MODULE_CTOR(CoulombMetric) { + description(desc); + satisfies_property_type(); + add_submodule("Two-center ERI"); +} + +MODULE_RUN(CoulombMetric) { + const auto& [braket] = pt::unwrap_inputs(inputs); + auto& eri2_mod = submods.at("Two-center ERI"); + + // Get ERI2 + const auto& I = eri2_mod.run_as(braket); + + // Cholesky Decomp + tensorwrapper::allocator::Eigen allocator(get_runtime()); + const auto& eigen_I = allocator.rebind(I.buffer()); + const auto* pI = eigen_I.data(); + const auto& shape_I = eigen_I.layout().shape().as_smooth(); + auto rows = shape_I.extent(0); + auto cols = shape_I.extent(1); + + using emat_t = Eigen::MatrixXd; + Eigen::Map I_map(pI, rows, cols); + Eigen::LLT lltOfI(I_map); + emat_t U = lltOfI.matrixU(); + emat_t Linv = U.inverse().transpose(); + + // Wrap result + tensorwrapper::shape::Smooth matrix_shape{rows, cols}; + tensorwrapper::layout::Physical matrix_layout(matrix_shape); + auto pM_buffer = allocator.allocate(matrix_layout); + + for(auto i = 0; i < rows; ++i) { + for(auto j = 0; j < cols; ++j) { pM_buffer->at(i, j) = Linv(i, j); } + } + simde::type::tensor M(matrix_shape, std::move(pM_buffer)); + + auto rv = results(); + return pt::wrap_results(rv, M); +} + +} // namespace integrals::ao_integrals \ No newline at end of file diff --git a/src/integrals/ao_integrals/df_integral.cpp b/src/integrals/ao_integrals/df_integral.cpp new file mode 100644 index 00000000..27315596 --- /dev/null +++ b/src/integrals/ao_integrals/df_integral.cpp @@ -0,0 +1,59 @@ +/* + * Copyright 2025 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 "ao_integrals.hpp" + +namespace integrals::ao_integrals { + +using pt = simde::ERI3; +using pt_2c = simde::ERI2; + +namespace { + +auto desc = R"( +Three-index ERI with Coulomb metric transformation +--------------------- +)"; + +} + +MODULE_CTOR(DFIntegral) { + description(desc); + satisfies_property_type(); + add_submodule("Three-center ERI"); + add_submodule("Coulomb Metric"); +} + +MODULE_RUN(DFIntegral) { + const auto& [braket] = pt::unwrap_inputs(inputs); + auto bra = braket.bra(); + auto ket = braket.ket(); + auto& op = braket.op(); + auto& eri2_mod = submods.at("Coulomb Metric"); + auto& eri3_mod = submods.at("Three-center ERI"); + + chemist::braket::BraKet aux_v_aux(bra, op, bra); + const auto& M = eri2_mod.run_as(aux_v_aux); + const auto& I = eri3_mod.run_as(braket); + + simde::type::tensor L; + L("i,k,l") = M("i,j") * I("j,k,l"); + + auto rv = results(); + return pt::wrap_results(rv, std::move(L)); +} + +} // namespace integrals::ao_integrals \ No newline at end of file diff --git a/src/integrals/ao_integrals/j_density_fitted.cpp b/src/integrals/ao_integrals/j_density_fitted.cpp new file mode 100644 index 00000000..9d4fe6a8 --- /dev/null +++ b/src/integrals/ao_integrals/j_density_fitted.cpp @@ -0,0 +1,68 @@ +/* + * 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 "ao_integrals.hpp" + +namespace integrals::ao_integrals { + +using pt = simde::aos_j_e_aos; +using pt_3c = simde::ERI3; +using aos_t = simde::type::aos; + +namespace { + +auto desc = R"( +Density Fitted J Builder +--------------------- +)"; + +} +MODULE_CTOR(JDensityFitted) { + description(desc); + satisfies_property_type(); + add_submodule("DF ERI"); + add_input("Auxiliary Basis Set"); +} + +MODULE_RUN(JDensityFitted) { + const auto&& [braket] = pt::unwrap_inputs(inputs); + + auto bra = braket.bra(); + auto ket = braket.ket(); + const auto& j_hat = braket.op(); + const auto& rho = j_hat.rhs_particle(); + const auto& p = rho.value(); + auto rho_aos = rho.basis_set(); + auto aux = inputs.at("Auxiliary Basis Set").value(); + auto& eri_mod = submods.at("DF ERI"); + + simde::type::v_ee_type v_ee; + simde::type::aos_squared ij_pair(bra, ket); + simde::type::aos_squared kl_pair(rho_aos, rho_aos); + chemist::braket::BraKet aux_v_ij(aux, v_ee, ij_pair); + chemist::braket::BraKet aux_v_kl(aux, v_ee, kl_pair); + const auto& I_akl = eri_mod.run_as(std::move(aux_v_kl)); + const auto& I_aij = eri_mod.run_as(std::move(aux_v_ij)); + + simde::type::tensor j; + j("a") = p("k,l") * I_akl("a,k,l"); + j("i,j") = j("a") * I_aij("a,i,j"); + + auto rv = results(); + return pt::wrap_results(rv, std::move(j)); +} + +} // namespace integrals::ao_integrals \ No newline at end of file diff --git a/src/integrals/ao_integrals/k_density_fitted.cpp b/src/integrals/ao_integrals/k_density_fitted.cpp new file mode 100644 index 00000000..f401b28c --- /dev/null +++ b/src/integrals/ao_integrals/k_density_fitted.cpp @@ -0,0 +1,69 @@ +/* + * 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 "ao_integrals.hpp" + +namespace integrals::ao_integrals { + +using pt = simde::aos_k_e_aos; +using pt_3c = simde::ERI3; +using aos_t = simde::type::aos; + +namespace { + +auto desc = R"( +Density Fitted K Builder +--------------------- +)"; + +} + +MODULE_CTOR(KDensityFitted) { + description(desc); + satisfies_property_type(); + add_submodule("DF ERI"); + add_input("Auxiliary Basis Set"); +} + +MODULE_RUN(KDensityFitted) { + const auto&& [braket] = pt::unwrap_inputs(inputs); + + auto bra = braket.bra(); + auto ket = braket.ket(); + const auto& j_hat = braket.op(); + const auto& rho = j_hat.rhs_particle(); + const auto& p = rho.value(); + auto rho_aos = rho.basis_set(); + auto aux = inputs.at("Auxiliary Basis Set").value(); + auto& eri_mod = submods.at("DF ERI"); + + simde::type::v_ee_type v_ee; + simde::type::aos_squared ik_pair(bra, rho_aos); + simde::type::aos_squared jl_pair(ket, rho_aos); + chemist::braket::BraKet aux_v_ik(aux, v_ee, ik_pair); + chemist::braket::BraKet aux_v_jl(aux, v_ee, jl_pair); + const auto& I_aik = eri_mod.run_as(std::move(aux_v_ik)); + const auto& I_ajl = eri_mod.run_as(std::move(aux_v_jl)); + + simde::type::tensor k; + k("a,i,l") = p("k,l") * I_aik("a,i,k"); + k("i,j") = k("a,i,l") * I_ajl("a,j,l"); + + auto rv = results(); + return pt::wrap_results(rv, std::move(k)); +} + +} // namespace integrals::ao_integrals \ No newline at end of file diff --git a/src/integrals/integrals_mm.cpp b/src/integrals/integrals_mm.cpp index 4b3a1ab6..cd1471b2 100644 --- a/src/integrals/integrals_mm.cpp +++ b/src/integrals/integrals_mm.cpp @@ -27,11 +27,13 @@ namespace integrals { * @throw none No throw guarantee */ void set_defaults(pluginplay::ModuleManager& mm) { - const auto ao_driver = "AO integral driver"; - mm.change_submod(ao_driver, "Kinetic", "Kinetic"); - mm.change_submod(ao_driver, "Electron-Nuclear attraction", "Nuclear"); + mm.change_submod("AO integral driver", "Kinetic", "Kinetic"); + mm.change_submod("AO integral driver", "Electron-Nuclear attraction", + "Nuclear"); mm.change_submod("Four center J builder", "Four-center ERI", "ERI4"); mm.change_submod("Four center K builder", "Four-center ERI", "ERI4"); + mm.change_submod("Density Fitting Integral", "Three-center ERI", "ERI3"); + mm.change_submod("Coulomb Metric", "Two-center ERI", "ERI2"); } void load_modules(pluginplay::ModuleManager& mm) { diff --git a/tests/cxx/unit/integrals/ao_integrals/coulomb_metric.cpp b/tests/cxx/unit/integrals/ao_integrals/coulomb_metric.cpp new file mode 100644 index 00000000..6ceea87d --- /dev/null +++ b/tests/cxx/unit/integrals/ao_integrals/coulomb_metric.cpp @@ -0,0 +1,47 @@ +/* + * Copyright 2022 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 "../testing.hpp" + +TEST_CASE("Coulomb Metric") { + using test_pt = simde::ERI2; + + pluginplay::ModuleManager mm; + integrals::load_modules(mm); + REQUIRE(mm.count("Coulomb Metric")); + + // Get basis set + auto mol = test::h2_molecule(); + auto aobs = test::h2_sto3g_basis_set(); + + // Make AOS object + simde::type::aos aos(aobs); + + // Make Operator + simde::type::v_ee_type op{}; + + // Make BraKet Input + chemist::braket::BraKet braket(aos, op, aos); + + // Call module + auto T = mm.at("Coulomb Metric").run_as(braket); + + auto t = test::eigen_tensor<2>(T.buffer()); + REQUIRE(t(0, 0) == Catch::Approx(0.15824726).margin(1E-6)); + REQUIRE(t(0, 1) == Catch::Approx(0.0).margin(1E-6)); + REQUIRE(t(1, 0) == Catch::Approx(-0.23097095).margin(1E-6)); + REQUIRE(t(1, 1) == Catch::Approx(0.27998174).margin(1E-6)); +} diff --git a/tests/cxx/unit/integrals/ao_integrals/df_integral.cpp b/tests/cxx/unit/integrals/ao_integrals/df_integral.cpp new file mode 100644 index 00000000..9123e0c2 --- /dev/null +++ b/tests/cxx/unit/integrals/ao_integrals/df_integral.cpp @@ -0,0 +1,52 @@ +/* + * Copyright 2025 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 "../testing.hpp" + +TEST_CASE("Density Fitting Integral") { + using test_pt = simde::ERI3; + + pluginplay::ModuleManager mm; + integrals::load_modules(mm); + REQUIRE(mm.count("Density Fitting Integral")); + + // Get basis set + auto mol = test::h2_molecule(); + auto aobs = test::h2_sto3g_basis_set(); + + // Make AOS object + simde::type::aos aos(aobs); + simde::type::aos_squared aos_squared(aos, aos); + + // Make Operator + simde::type::v_ee_type op{}; + + // Make BraKet Input + chemist::braket::BraKet braket(aos, op, aos_squared); + + // Call module + auto T = mm.at("Density Fitting Integral").run_as(braket); + + auto t = test::eigen_tensor<3>(T.buffer()); + REQUIRE(t(0, 0, 0) == Catch::Approx(0.81362039).margin(1E-6)); + REQUIRE(t(0, 0, 1) == Catch::Approx(0.31266336).margin(1E-6)); + REQUIRE(t(0, 1, 0) == Catch::Approx(0.31266336).margin(1E-6)); + REQUIRE(t(0, 1, 1) == Catch::Approx(0.60009419).margin(1E-6)); + REQUIRE(t(1, 0, 0) == Catch::Approx(-0.12580195).margin(1E-6)); + REQUIRE(t(1, 0, 1) == Catch::Approx(0.09683503).margin(1E-6)); + REQUIRE(t(1, 1, 0) == Catch::Approx(0.09683503).margin(1E-6)); + REQUIRE(t(1, 1, 1) == Catch::Approx(0.56364026).margin(1E-6)); +} \ No newline at end of file diff --git a/tests/cxx/unit/integrals/ao_integrals/j_density_fitted.cpp b/tests/cxx/unit/integrals/ao_integrals/j_density_fitted.cpp new file mode 100644 index 00000000..ecd43771 --- /dev/null +++ b/tests/cxx/unit/integrals/ao_integrals/j_density_fitted.cpp @@ -0,0 +1,48 @@ +/* + * 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 "../testing.hpp" + +TEST_CASE("Density Fitted J builder") { + using pt = simde::aos_j_e_aos; + + pluginplay::ModuleManager mm; + integrals::load_modules(mm); + REQUIRE(mm.count("Density Fitted J builder")); + + // Get basis set + auto mol = test::h2_molecule(); + auto aobs = test::h2_sto3g_basis_set(); + + // Make AOS object + simde::type::aos aos(aobs); + + // Make Operator + simde::type::j_e_type op(simde::type::electron{}, test::h2_density()); + + // Make BraKet Input + chemist::braket::BraKet braket(aos, op, aos); + + // Call module + mm.change_input("Density Fitted J builder", "Auxiliary Basis Set", aos); + const auto& T = mm.at("Density Fitted J builder").run_as(braket); + + auto t = test::eigen_tensor<2>(T.buffer()); + REQUIRE(t(0, 0) == Catch::Approx(0.50515668).margin(1E-6)); + REQUIRE(t(0, 1) == Catch::Approx(0.22344536).margin(1E-6)); + REQUIRE(t(1, 0) == Catch::Approx(0.22344536).margin(1E-6)); + REQUIRE(t(1, 1) == Catch::Approx(0.50515668).margin(1E-6)); +} \ No newline at end of file diff --git a/tests/cxx/unit/integrals/ao_integrals/k_density_fitted.cpp b/tests/cxx/unit/integrals/ao_integrals/k_density_fitted.cpp new file mode 100644 index 00000000..04e64e9d --- /dev/null +++ b/tests/cxx/unit/integrals/ao_integrals/k_density_fitted.cpp @@ -0,0 +1,48 @@ +/* + * 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 "../testing.hpp" + +TEST_CASE("Density Fitted K builder") { + using pt = simde::aos_k_e_aos; + + pluginplay::ModuleManager mm; + integrals::load_modules(mm); + REQUIRE(mm.count("Density Fitted K builder")); + + // Get basis set + auto mol = test::h2_molecule(); + auto aobs = test::h2_sto3g_basis_set(); + + // Make AOS object + simde::type::aos aos(aobs); + + // Make Operator + simde::type::k_e_type op(simde::type::electron{}, test::h2_density()); + + // Make BraKet Input + chemist::braket::BraKet braket(aos, op, aos); + + // Call module + mm.change_input("Density Fitted K builder", "Auxiliary Basis Set", aos); + const auto& T = mm.at("Density Fitted K builder").run_as(braket); + + auto t = test::eigen_tensor<2>(T.buffer()); + REQUIRE(t(0, 0) == Catch::Approx(0.40594955).margin(1E-6)); + REQUIRE(t(0, 1) == Catch::Approx(0.32265250).margin(1E-6)); + REQUIRE(t(1, 0) == Catch::Approx(0.32265250).margin(1E-6)); + REQUIRE(t(1, 1) == Catch::Approx(0.40594955).margin(1E-6)); +} \ No newline at end of file