test: Added boilerplate code for different Yang-Mills correlation fun…

…ctions and used these to rewrite the benchmarks to make a more realistic benchmark. Signed-off-by: Franz R. Sattler <[email protected]>
satfra · Jan 10, 2025 · 12db594 · 12db594
1 parent f4f439f
commit 12db594
Show file tree

Hide file tree

Showing 15 changed files with 1,335 additions and 446 deletions.
diff --git a/DiFfRG/tests/physics/boilerplate/1_angle/AA/ZA.hh b/DiFfRG/tests/physics/boilerplate/1_angle/AA/ZA.hh
@@ -0,0 +1,74 @@
+#pragma once
+
+#include "../def.hh"
+#include "ZA.kernel"
+
+#include <future>
+#include <memory>
+
+namespace DiFfRG
+{
+  namespace Flows
+  {
+    template <template <int, typename, typename> class INT> class ZA_integrator
+    {
+    public:
+      ZA_integrator(QuadratureProvider &quadrature_provider, std::array<uint, 2> grid_sizes, const double x_extent,
+                    const JSONValue &json)
+          : quadrature_provider(quadrature_provider), grid_sizes(grid_sizes), x_extent(x_extent),
+            jacobian_quadrature_factor(json.get_double("/integration/jacobian_quadrature_factor")), json(json)
+      {
+        integrator =
+            std::make_unique<INT<4, double, ZA_kernel<__REGULATOR__>>>(quadrature_provider, grid_sizes, x_extent, json);
+      }
+      ~ZA_integrator() = default;
+
+      template <typename NT, typename... T> std::future<NT> request(T &&...t)
+      {
+        static_assert(std::is_same_v<NT, double>, "Unknown type requested of ZA_integrator::request");
+        if constexpr (std::is_same_v<NT, double>) return request_CT(std::forward<T>(t)...);
+      }
+
+      template <typename NT, typename... T> NT get(T &&...t)
+      {
+        static_assert(std::is_same_v<NT, double>, "Unknown type requested of ZA_integrator::request");
+        if constexpr (std::is_same_v<NT, double>) return get_CT(std::forward<T>(t)...);
+      }
+
+    private:
+      std::future<double> request_CT(const double k, const double p,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA3,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZAcbc,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA4,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &dtZc,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &Zc,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &dtZA,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA,
+                                     const double m2A)
+      {
+        return integrator->request(k, p, ZA3, ZAcbc, ZA4, dtZc, Zc, dtZA, ZA, m2A);
+      }
+
+      double get_CT(const double k, const double p,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA3,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZAcbc,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA4,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &dtZc,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &Zc,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &dtZA,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA, const double m2A)
+      {
+        return integrator->get(k, p, ZA3, ZAcbc, ZA4, dtZc, Zc, dtZA, ZA, m2A);
+      }
+
+      QuadratureProvider &quadrature_provider;
+      const std::array<uint, 2> grid_sizes;
+      std::array<uint, 2> jac_grid_sizes;
+      const double x_extent;
+      const double jacobian_quadrature_factor;
+      const JSONValue json;
+
+      std::unique_ptr<INT<4, double, ZA_kernel<__REGULATOR__>>> integrator;
+    };
+  } // namespace Flows
+} // namespace DiFfRG
diff --git a/DiFfRG/tests/physics/boilerplate/1_angle/AA/ZA.kernel b/DiFfRG/tests/physics/boilerplate/1_angle/AA/ZA.kernel
@@ -0,0 +1,57 @@
+template <typename REG> class ZA_kernel
+{
+public:
+  template <typename T1, typename T2, typename T3, typename T4, typename T5, typename T6, typename T7, typename T8, typename T9>
+  static __forceinline__ __host__ __device__ auto kernel(const double q, const double cos1, const double k, const T1 p, const T2 &ZA3, const T3 &ZAcbc, const T4 &ZA4,
+                                                         const T5 &dtZc, const T6 &Zc, const T7 &dtZA, const T8 &ZA, const T9 m2A)
+  {
+    // optimization definitions
+    const double DiFfRGREPL1 = dtZA(pow(1. + powr<6>(k), 0.16666666666666666));
+    const double DiFfRGREPL2 = dtZc(k);
+    const double DiFfRGREPL3 = RB(powr<2>(k), powr<2>(q));
+    const double DiFfRGREPL4 = RB(powr<2>(k), -2. * cos1 * p * q + powr<2>(p) + powr<2>(q));
+    const double DiFfRGREPL5 = RBdot(powr<2>(k), powr<2>(q));
+    const double DiFfRGREPL6 = ZA(pow(1. + powr<6>(k), 0.16666666666666666));
+    const double DiFfRGREPL7 = ZA(1.02 * pow(1. + powr<6>(k), 0.16666666666666666));
+    const double DiFfRGREPL8 = ZA(q);
+    const double DiFfRGREPL9 = Zc(k);
+    const double DiFfRGREPL10 = Zc(1.02 * k);
+    const double DiFfRGREPL11 = Zc(q);
+
+    // flow
+    return powr<-2>(p) *
+           (-4. * (-1. + powr<2>(cos1)) *
+                (-6. * cos1 * q * powr<3>(p) + 3. * powr<4>(p) + (8. + powr<2>(cos1)) * powr<2>(p) * powr<2>(q) - 6. * cos1 * p * powr<3>(q) + 3. * powr<4>(q)) *
+                (DiFfRGREPL5 * DiFfRGREPL6 + DiFfRGREPL3 * (DiFfRGREPL1 + 50. * (-1. * DiFfRGREPL6 + DiFfRGREPL7) * powr<6>(k) * powr<-1>(1. + powr<6>(k)))) *
+                powr<-1>(-2. * cos1 * p * q + powr<2>(p) + powr<2>(q)) * powr<-2>(DiFfRGREPL3 * DiFfRGREPL6 + DiFfRGREPL8 * powr<2>(q)) *
+                powr<-1>(DiFfRGREPL4 * DiFfRGREPL6 + (-2. * cos1 * p * q + powr<2>(p) + powr<2>(q)) * ZA((sqrt(-2. * cos1 * p * q + powr<2>(p) + powr<2>(q))))) *
+                powr<2>(ZA3(0.816496580927726 * (sqrt(-1. * cos1 * p * q + powr<2>(p) + powr<2>(q))))) +
+            (DiFfRGREPL3 * (50. * DiFfRGREPL10 + DiFfRGREPL2 - 50. * DiFfRGREPL9) + DiFfRGREPL5 * DiFfRGREPL9) * (-1. + powr<2>(cos1)) * powr<2>(q) *
+                powr<-2>(DiFfRGREPL3 * DiFfRGREPL9 + DiFfRGREPL11 * powr<2>(q)) * powr<2>(ZAcbc(0.816496580927726 * (sqrt(-1. * cos1 * p * q + powr<2>(p) + powr<2>(q))))) *
+                powr<-1>(DiFfRGREPL4 * DiFfRGREPL9 + (-2. * cos1 * p * q + powr<2>(p) + powr<2>(q)) * Zc((sqrt(-2. * cos1 * p * q + powr<2>(p) + powr<2>(q))))) +
+            (DiFfRGREPL3 * (50. * DiFfRGREPL10 + DiFfRGREPL2 - 50. * DiFfRGREPL9) + DiFfRGREPL5 * DiFfRGREPL9) * (-1. + powr<2>(cos1)) * powr<2>(q) *
+                powr<-2>(DiFfRGREPL3 * DiFfRGREPL9 + DiFfRGREPL11 * powr<2>(q)) * powr<2>(ZAcbc(0.816496580927726 * (sqrt(cos1 * p * q + powr<2>(p) + powr<2>(q))))) *
+                powr<-1>(DiFfRGREPL9 * RB(powr<2>(k), 2. * cos1 * p * q + powr<2>(p) + powr<2>(q)) +
+                         (2. * cos1 * p * q + powr<2>(p) + powr<2>(q)) * Zc((sqrt(2. * cos1 * p * q + powr<2>(p) + powr<2>(q))))) +
+            (-7. + powr<2>(cos1)) * (DiFfRGREPL5 * DiFfRGREPL6 + DiFfRGREPL3 * (DiFfRGREPL1 + 50. * (-1. * DiFfRGREPL6 + DiFfRGREPL7) * powr<6>(k) * powr<-1>(1. + powr<6>(k)))) *
+                powr<-2>(DiFfRGREPL3 * DiFfRGREPL6 + DiFfRGREPL8 * powr<2>(q)) * ZA4(0.7071067811865475 * (sqrt(powr<2>(p) + powr<2>(q)))));
+  }
+
+  template <typename T1, typename T2, typename T3, typename T4, typename T5, typename T6, typename T7, typename T8, typename T9>
+  static __forceinline__ __host__ __device__ auto constant(const double k, const T1 p, const T2 &ZA3, const T3 &ZAcbc, const T4 &ZA4, const T5 &dtZc, const T6 &Zc, const T7 &dtZA,
+                                                           const T8 &ZA, const T9 m2A)
+  {
+    // flow
+    return 0.;
+  }
+
+private:
+  static __forceinline__ __device__ __host__ auto RB(const auto k2, const auto p2) { return REG::RB(k2, p2); }
+  static __forceinline__ __device__ __host__ auto RF(const auto k2, const auto p2) { return REG::RF(k2, p2); }
+
+  static __forceinline__ __device__ __host__ auto RBdot(const auto k2, const auto p2) { return REG::RBdot(k2, p2); }
+  static __forceinline__ __device__ __host__ auto RFdot(const auto k2, const auto p2) { return REG::RFdot(k2, p2); }
+
+  static __forceinline__ __device__ __host__ auto dq2RB(const auto k2, const auto p2) { return REG::dq2RB(k2, p2); }
+  static __forceinline__ __device__ __host__ auto dq2RF(const auto k2, const auto p2) { return REG::dq2RF(k2, p2); }
+};
diff --git a/DiFfRG/tests/physics/boilerplate/1_angle/def.hh b/DiFfRG/tests/physics/boilerplate/1_angle/def.hh
@@ -0,0 +1,11 @@
+#pragma once
+
+#include <DiFfRG/common/utils.hh>
+#include <DiFfRG/physics/integration.hh>
+#include <DiFfRG/physics/interpolation.hh>
+#include <DiFfRG/physics/regulators.hh>
+#include <DiFfRG/physics/thermodynamics.hh>
+
+using namespace ::DiFfRG;
+
+#define __REGULATOR__ ::DiFfRG::PolynomialExpRegulator<>
diff --git a/DiFfRG/tests/physics/boilerplate/1_angle/flows.hh b/DiFfRG/tests/physics/boilerplate/1_angle/flows.hh
@@ -0,0 +1,30 @@
+#pragma once
+
+#include "AA/ZA.hh"
+
+#include "def.hh"
+#include <DiFfRG/physics/flow_equations.hh>
+
+template <template <int, typename, typename> class INT> class YangMillsFlowEquations : public FlowEquations
+{
+public:
+  YangMillsFlowEquations(const JSONValue &json)
+      : FlowEquations(json,
+                      [](double x) { return powr<-1>(x + __REGULATOR__::RB(1., x)) * __REGULATOR__::RBdot(1., x); }),
+        grid_size_int{{x_quadrature_order}}, grid_sizes_angle_int{{x_quadrature_order, 2 * angle_quadrature_order}},
+        grid_sizes_3D_int{{x_quadrature_order, angle_quadrature_order, angle_quadrature_order}},
+        grid_sizes_4D_int{{x_quadrature_order, angle_quadrature_order, angle_quadrature_order, angle_quadrature_order}},
+        ZA_integrator(quadrature_provider, grid_sizes_angle_int, x_extent, json)
+  {
+  }
+
+private:
+  const std::array<uint, 1> grid_size_int;
+  const std::array<uint, 2> grid_sizes_angle_int;
+  const std::array<uint, 3> grid_sizes_3D_int;
+  const std::array<uint, 4> grid_sizes_4D_int;
+
+public:
+  ::DiFfRG::QuadratureProvider quadrature_provider;
+  ::DiFfRG::Flows::ZA_integrator<INT> ZA_integrator;
+};
diff --git a/DiFfRG/tests/physics/boilerplate/2_angle/AAA/ZA3.hh b/DiFfRG/tests/physics/boilerplate/2_angle/AAA/ZA3.hh
@@ -0,0 +1,76 @@
+#pragma once
+
+#include "../def.hh"
+#include "ZA3.kernel"
+
+#include <future>
+#include <memory>
+
+namespace DiFfRG
+{
+  namespace Flows
+  {
+    template <template <typename, typename> class INT> class ZA3_integrator
+    {
+    public:
+      ZA3_integrator(QuadratureProvider &quadrature_provider, std::array<uint, 3> grid_sizes, const double x_extent,
+                     const JSONValue &json)
+          : quadrature_provider(quadrature_provider), grid_sizes(grid_sizes), x_extent(x_extent),
+            jacobian_quadrature_factor(json.get_double("/integration/jacobian_quadrature_factor")), json(json)
+      {
+        integrator =
+            std::make_unique<INT<double, ZA3_kernel<__REGULATOR__>>>(quadrature_provider, grid_sizes, x_extent, json);
+      }
+
+      ~ZA3_integrator() = default;
+
+      template <typename NT, typename... T> std::future<NT> request(T &&...t)
+      {
+        static_assert(std::is_same_v<NT, double>, "Unknown type requested of ZA3_integrator::request");
+        if constexpr (std::is_same_v<NT, double>) return request_CT(std::forward<T>(t)...);
+      }
+
+      template <typename NT, typename... T> NT get(T &&...t)
+      {
+        static_assert(std::is_same_v<NT, double>, "Unknown type requested of ZA3_integrator::request");
+        if constexpr (std::is_same_v<NT, double>) return get_CT(std::forward<T>(t)...);
+      }
+
+    private:
+      std::future<double> request_CT(const double k, const double p,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA3,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZAcbc,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA4,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &dtZc,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &Zc,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &dtZA,
+                                     const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA,
+                                     const double m2A)
+      {
+        return integrator->request(k, p, ZA3, ZAcbc, ZA4, dtZc, Zc, dtZA, ZA, m2A);
+      }
+
+      double get_CT(const double k, const double p,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA3,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZAcbc,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA4,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &dtZc,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &Zc,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &dtZA,
+                    const TexLinearInterpolator1D<double, LogarithmicCoordinates1D<float>> &ZA, const double m2A)
+      {
+        return integrator->get(k, p, ZA3, ZAcbc, ZA4, dtZc, Zc, dtZA, ZA, m2A);
+      }
+
+      QuadratureProvider &quadrature_provider;
+      const std::array<uint, 3> grid_sizes;
+      std::array<uint, 3> jac_grid_sizes;
+      const double x_extent;
+      const double jacobian_quadrature_factor;
+      const JSONValue json;
+
+      std::unique_ptr<INT<double, ZA3_kernel<__REGULATOR__>>> integrator;
+    };
+
+  } // namespace Flows
+} // namespace DiFfRG