Slater function D=1,2,3 implemented and tested

src/functions/Slater.cpp

@@ -45,7 +45,8 @@ Slater<D>::Slater(double a, double c, const Coord<D> &r)
 
 template<> double Slater<1>::calcSquareNorm() const {
     double c2 = this->coef * this->coef;
-    return c2 / this->alpha;
+    double square_norm = c2 / this->alpha;
+    return square_norm;
 }
 
 template<> double Slater<2>::calcSquareNorm() const {

src/functions/Slater.h

@@ -77,10 +77,10 @@ template <int D> class Slater : public RepresentableFunction<D> {
     double alpha;                /**< exponent  */
     Coord<D> pos;                /**< center  */
 
-    bool isVisibleAtScale(int scale, int nQuadPts) const;
-    bool isZeroOnInterval(const double *a, const double *b) const;
+    bool isVisibleAtScale(int scale, int nQuadPts) const {return true;}
+    bool isZeroOnInterval(const double *a, const double *b) const {return false;}
 
-    virtual std::ostream &print(std::ostream &o) const = 0;
+    //    virtual std::ostream &print(std::ostream &o) const = 0;
 };
 
 } // namespace mrcpp

tests/functions/CMakeLists.txt

@@ -2,12 +2,14 @@ target_sources(mrcpp-tests PRIVATE
     ${CMAKE_CURRENT_SOURCE_DIR}/legendre_poly.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/polynomial.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/gaussians.cpp
+    ${CMAKE_CURRENT_SOURCE_DIR}/slater.cpp
     ${CMAKE_CURRENT_SOURCE_DIR}/periodify_gaussians.cpp
     )
 
-add_Catch_test(NAME legendre_poly           LABELS legendre_poly)
-add_Catch_test(NAME polynomials             LABELS polynomials)
-add_Catch_test(NAME gaussians               LABELS gaussians)
-add_Catch_test(NAME periodic_narrow_gaussian       LABELS periodic_narrow_gaussian)
-add_Catch_test(NAME periodic_wide_gaussian       LABELS periodic_wide_gaussian)
-add_Catch_test(NAME periodic_gaussExp       LABELS periodic_gausExp)
+add_Catch_test(NAME legendre_poly            LABELS legendre_poly)
+add_Catch_test(NAME polynomials              LABELS polynomials)
+add_Catch_test(NAME gaussians                LABELS gaussians)
+add_Catch_test(NAME slater                   LABELS slater)
+add_Catch_test(NAME periodic_narrow_gaussian LABELS periodic_narrow_gaussian)
+add_Catch_test(NAME periodic_wide_gaussian   LABELS periodic_wide_gaussian)
+add_Catch_test(NAME periodic_gaussExp        LABELS periodic_gausExp)

tests/functions/slater.cpp

@@ -0,0 +1,126 @@
+/*
+ * MRCPP, a numerical library based on multiresolution analysis and
+ * the multiwavelet basis which provide low-scaling algorithms as well as
+ * rigorous error control in numerical computations.
+ * Copyright (C) 2021 Stig Rune Jensen, Jonas Juselius, Luca Frediani and contributors.
+ *
+ * This file is part of MRCPP.
+ *
+ * MRCPP is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * MRCPP is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public License
+ * along with MRCPP.  If not, see <https://www.gnu.org/licenses/>.
+ *
+ * For information on the complete list of contributors to MRCPP, see:
+ * <https://mrcpp.readthedocs.io/>
+ */
+
+#include "catch2/catch_all.hpp"
+
+#include "factory_functions.h"
+
+#include "functions/Slater.h"
+#include "treebuilders/grid.h"
+#include "treebuilders/project.h"
+#include "utils/math_utils.h"
+
+using namespace mrcpp;
+
+namespace gaussians {
+
+template <int D> void testSlaterFunction();
+
+SCENARIO("Slater", "[slater]") {
+    GIVEN("Slater Function in 1D") { testSlaterFunction<1>(); }
+    GIVEN("Slater Function in 2D") { testSlaterFunction<2>(); }
+    GIVEN("Slater Function in 3D") { testSlaterFunction<3>(); }
+}
+
+template <int D> void testSlaterFunction() {
+    const auto prec = 1.0e-4;
+    // Making normalized gaussian centered at the origin
+    const auto alpha = 2.;
+    const auto coeff = 3.;   
+    double pos_data[3] = {0.1, 0.2, 0.3};
+    double ref_data[3] = {-0.2, 0.5, 1.0};
+    auto pos = mrcpp::details::convert_to_std_array<double, D>(pos_data);
+    auto ref = mrcpp::details::convert_to_std_array<double, D>(ref_data);
+
+    auto slater = Slater<D>(alpha, coeff, pos);
+
+    // Making ref value
+    auto dist = math_utils::calc_distance<D>(slater.getPos(), ref);
+    const auto ref_val = slater.getCoef() * std::exp(-slater.getAlpha() * std::abs(dist));
+
+    // Making MRA
+    const auto min_scale = -4;
+    auto corner = std::array<int, D>{};
+    auto boxes = std::array<int, D>{};
+    corner.fill(-1);
+    boxes.fill(2);
+    auto world = BoundingBox<D>(min_scale, corner, boxes);
+    const auto basis = InterpolatingBasis(7);
+    auto MRA = MultiResolutionAnalysis<D>(world, basis, 25);
+
+    // Making a function tree and projects the slater fcn onto it
+    FunctionTree<D> f_tree(MRA);
+    build_grid<D>(f_tree, slater);
+    project(prec, f_tree, slater);
+    f_tree.calcSquareNorm();
+
+    WHEN("Slater function is projected") {
+        THEN("The Slater function can be evaluated") { REQUIRE(slater.evalf(ref) == Catch::Approx(ref_val)); }
+        THEN("The projected Slater function can be evaluated") { REQUIRE(f_tree.evalf(ref) == Catch::Approx(ref_val)); }
+        THEN("the square norm matches the analytical value") { REQUIRE(f_tree.getSquareNorm() == Catch::Approx(slater.calcSquareNorm())); }
+    }
+}
+
+    /*
+SCENARIO("Slater_2D", "[slater_2d]") {
+    const auto D = 2;
+    const auto prec = 1.0e-4;
+    // Making normalized gaussian centered at the origin
+    const auto pos = Coord<D>{0.0, 0.0};
+    const auto alpha = 2.2;
+    const auto coeff = 3.3;   
+    auto slater = Slater<D>(alpha, coeff, pos);
+
+    // Making ref value
+    const auto r_ref = std::array<double, D>{.1, 0.2};
+    double dist = math_utils::calc_distance<D>(r_ref, slater.getPos());
+    const auto ref_val = slater.getCoef() * std::exp(-slater.getAlpha() * std::abs(dist));
+
+    // Making MRA
+    const auto min_scale = -4;
+    const auto corner = std::array<int, D>{-1,-1};
+    const auto boxes = std::array<int, D>{2,2};
+    auto world = BoundingBox<D>(min_scale, corner, boxes);
+    const auto basis = InterpolatingBasis(7);
+    auto MRA = MultiResolutionAnalysis<D>(world, basis, 25);
+
+    // Making a function tree and projects the slater fcn onto it
+    FunctionTree<D> f_tree(MRA);
+    build_grid<D>(f_tree, slater);
+    project(prec, f_tree, slater);
+    f_tree.calcSquareNorm();
+
+    WHEN("Slater function is projected") {
+        THEN("The Slater function can be evaluated") { REQUIRE(slater.evalf(r_ref) == Catch::Approx(ref_val)); }
+
+        THEN("The projected Slater function can be evaluated") { REQUIRE(f_tree.evalf(r_ref) == Catch::Approx(ref_val)); }
+
+        THEN("the square norm matches") { REQUIRE(f_tree.getSquareNorm() == Catch::Approx(slater.calcSquareNorm())); }
+
+    }
+}
+    */
+
+} // namespace gaussians