[GenVectorX] Add test Invariant Masses

math/experimental/genvectorx/test/CMakeLists.txt

@@ -8,14 +8,22 @@
 
 
     ROOT_EXECUTABLE(testGenvectorSYCL testGenVectorSYCL.cxx) # LIBRARIES GenVectorSYCL)
-    if (oneapi)
-        target_link_options(testGenvectorSYCL PUBLIC -E)
-    endif()
     target_link_libraries(testGenvectorSYCL PUBLIC GenVectorSYCL )
     add_sycl_to_root_target(TARGET testGenvectorSYCL 
     SOURCES testGenVectorSYCL.cxx 
     COMPILE_DEFINITIONS ROOT_MATH_SYCL
     DEPENDENCIES GenVectorSYCL
     )
+
+
+    ROOT_EXECUTABLE(testGenvectorSYCLInvMasses testGenVectorSYCLInvMasses.cxx) # LIBRARIES GenVectorSYCL)
+    target_link_libraries(testGenvectorSYCLInvMasses PUBLIC GenVectorSYCL )
+    add_sycl_to_root_target(TARGET testGenvectorSYCLInvMasses 
+    SOURCES testGenvectorSYCLInvMasses.cxx 
+    COMPILE_DEFINITIONS ROOT_MATH_SYCL
+    DEPENDENCIES GenVectorSYCL
+    )
+
     ROOT_ADD_TEST(test-genvector-genvectorsycl COMMAND testGenvectorSYCL)
+    ROOT_ADD_TEST(test-genvector-genvectorsyclinvmasses COMMAND testGenvectorSYCLInvMasses)
 

math/experimental/genvectorx/test/testGenVectorSYCL.cxx

@@ -105,12 +105,15 @@ int testVector3D()
 {
    int iret_host = 0;
 
-   std::cout << "testing Vector3D   \t:\t";
+   std::cout << "testing Vector3D   \t:\n";
 
    sycl::buffer<int, 1> iret_buf(&iret_host, sycl::range<1>(1));
    sycl::default_selector device_selector;
    sycl::queue queue(device_selector);
 
+   std::cout << "sycl::queue check - selected device:\n"
+             << queue.get_device().get_info<sycl::info::device::name>() << std::endl;
+
    {
       queue.submit([&](sycl::handler &cgh) {
          auto iret = iret_buf.get_access<sycl::access::mode::read_write>(cgh);
@@ -167,12 +170,15 @@ int testPoint3D()
 
    int iret_host = 0;
 
-   std::cout << "testing Point3D    \t:\t";
+   std::cout << "testing Point3D    \t:\n";
 
    sycl::buffer<int, 1> iret_buf(&iret_host, sycl::range<1>(1));
    sycl::default_selector device_selector;
    sycl::queue queue(device_selector);
 
+   std::cout << "sycl::queue check - selected device:\n"
+             << queue.get_device().get_info<sycl::info::device::name>() << std::endl;
+
    {
       queue.submit([&](sycl::handler &cgh) {
          auto iret = iret_buf.get_access<sycl::access::mode::read_write>(cgh);
@@ -240,12 +246,15 @@ int testVector2D()
 
    int iret_host = 0;
 
-   std::cout << "testing Vector2D   \t:\t";
+   std::cout << "testing Vector2D   \t:\n";
 
    sycl::buffer<int, 1> iret_buf(&iret_host, sycl::range<1>(1));
    sycl::default_selector device_selector;
    sycl::queue queue(device_selector);
 
+   std::cout << "sycl::queue check - selected device:\n"
+             << queue.get_device().get_info<sycl::info::device::name>() << std::endl;
+
    {
       queue.submit([&](sycl::handler &cgh) {
          auto iret = iret_buf.get_access<sycl::access::mode::read_write>(cgh);
@@ -316,11 +325,14 @@ int testPoint2D()
 
    int iret_host = 0;
 
-   std::cout << "testing Point2D    \t:\t";
+   std::cout << "testing Point2D    \t:\n";
 
    sycl::buffer<int, 1> iret_buf(&iret_host, sycl::range<1>(1));
    sycl::queue queue{sycl::default_selector_v};
 
+   std::cout << "sycl::queue check - selected device:\n"
+             << queue.get_device().get_info<sycl::info::device::name>() << std::endl;
+
    {
       queue.submit([&](sycl::handler &cgh) {
          auto iret = iret_buf.get_access<sycl::access::mode::read_write>(cgh);
@@ -391,12 +403,15 @@ int testRotations3D()
 {
 
    int iret_host = 0;
-   std::cout << "testing 3D Rotations\t:\t";
+   std::cout << "testing 3D Rotations\t:\n";
 
    sycl::buffer<int, 1> iret_buf(&iret_host, sycl::range<1>(1));
    sycl::default_selector device_selector;
    sycl::queue queue(device_selector);
 
+   std::cout << "sycl::queue check - selected device:\n"
+             << queue.get_device().get_info<sycl::info::device::name>() << std::endl;
+
    {
       queue.submit([&](sycl::handler &cgh) {
          auto iret = iret_buf.get_access<sycl::access::mode::read_write>(cgh);
@@ -529,12 +544,15 @@ int testTransform3D()
 {
 
    int iret_host = 0;
-   std::cout << "testing 3D Transform\t:\t";
+   std::cout << "testing 3D Transform\t:\n";
 
    sycl::buffer<int, 1> iret_buf(&iret_host, sycl::range<1>(1));
    sycl::default_selector device_selector;
    sycl::queue queue(device_selector);
 
+   std::cout << "sycl::queue check - selected device:\n"
+             << queue.get_device().get_info<sycl::info::device::name>() << std::endl;
+
    {
       queue.submit([&](sycl::handler &cgh) {
          auto iret = iret_buf.get_access<sycl::access::mode::read_write>(cgh);
@@ -757,12 +775,15 @@ int testTransform3D()
 int testVectorUtil()
 {
    int iret_host = 0;
-   std::cout << "testing VectorUtil  \t:\t";
+   std::cout << "testing VectorUtil  \t:\n";
 
    sycl::buffer<int, 1> iret_buf(&iret_host, sycl::range<1>(1));
    sycl::default_selector device_selector;
    sycl::queue queue(device_selector);
 
+   std::cout << "sycl::queue check - selected device:\n"
+             << queue.get_device().get_info<sycl::info::device::name>() << std::endl;
+
    {
       queue.submit([&](sycl::handler &cgh) {
          auto iret = iret_buf.get_access<sycl::access::mode::read_write>(cgh);
@@ -834,20 +855,36 @@ int testVectorUtil()
 
 int testLorentzVector()
 {
-   std::cout << "testing LorentzVector  \t:\t";
-   int iret = 0;
-   LorentzVector<PtEtaPhiM4D<float>> v1(1, 2, 3, 4);
-   LorentzVector<PtEtaPhiM4D<float>> v2(5, 6, 7, 8);
-   iret |= compare(v1.DeltaR(v2), 4.60575f);
-   // Result cross-validated using:
-   // TLorentzVector t1, t2;
-   // t1.SetPtEtaPhiE(1,2,3,4); t2.SetPtEtaPhiE(5,6,7,8);
-   // t1.DeltaR(t2)
-   if (iret == 0)
-      std::cout << "\t\t\tOK\n";
+   int iret_host = 0;
+   std::cout << "testing LorentzVector  \t:\n";
+
+   sycl::buffer<int, 1> iret_buf(&iret_host, sycl::range<1>(1));
+   sycl::default_selector device_selector;
+   sycl::queue queue(device_selector);
+
+   std::cout << "sycl::queue check - selected device:\n"
+             << queue.get_device().get_info<sycl::info::device::name>() << std::endl;
+
+   {
+      queue.submit([&](sycl::handler &cgh) {
+         auto iret = iret_buf.get_access<sycl::access::mode::read_write>(cgh);
+         cgh.single_task<class testRotations3D>([=]() {
+            LorentzVector<PtEtaPhiM4D<float>> v1(1, 2, 3, 4);
+            LorentzVector<PtEtaPhiM4D<float>> v2(5, 6, 7, 8);
+            iret[0] |= compare(v1.DeltaR(v2), 4.60575f);
+
+            LorentzVector<PtEtaPhiM4D<float>> v = v1 + v2;
+            iret[0] |= compare(v.M(), 62.03058f);
+         });
+      });
+   }
+
+   if (iret_host == 0)
+      std::cout << "\tOK\n";
    else
-      std::cout << "\t\t\t\t\t\tFAILED\n";
-   return iret;
+      std::cout << "\t FAILED\n";
+
+   return iret_host;
 }
 
 int testGenVector()

math/experimental/genvectorx/test/testGenVectorSYCLInvMasses.cxx

@@ -0,0 +1,102 @@
+#include "MathX/LorentzRotation.h"
+#include "MathX/PtEtaPhiM4D.h"
+#include "MathX/LorentzVector.h"
+
+#include "Math/LorentzRotation.h"
+#include "Math/PtEtaPhiM4D.h"
+#include "Math/LorentzVector.h"
+
+#include "MathX/GenVectorX/AccHeaders.h"
+
+#include <sycl/sycl.hpp>
+
+#include <vector>
+
+using namespace ROOT::ROOT_MATH_ARCH;
+
+typedef LorentzVector<PtEtaPhiM4D<double>> vec4dSYCL;
+typedef ROOT::Math::LorentzVector<ROOT::Math::PtEtaPhiM4D<double>> vec4d;
+
+template <class vec4d>
+vec4d *GenVectors(int n)
+{
+   vec4d *vectors = new vec4d[n];
+
+   // generate n -4 momentum quantities
+   for (int i = 0; i < n; ++i) {
+      // fill vectors
+      vectors[i] = {1., 1., 1., 1.};
+   }
+
+   return vectors;
+}
+
+double *InvariantMassesSYCL(int N)
+{
+   sycl::default_selector device_selector;
+   sycl::queue queue(device_selector);
+
+   auto v1 = GenVectors<LorentzVector<PtEtaPhiM4D<double>>>(N);
+   auto v2 = GenVectors<LorentzVector<PtEtaPhiM4D<double>>>(N);
+   double *invMasses = new double[N];
+
+   std::cout << "sycl::queue check - selected device:\n"
+             << queue.get_device().get_info<sycl::info::device::name>() << std::endl;
+
+   {
+
+      vec4dSYCL *d_v1 = sycl::malloc_device<vec4dSYCL>(N, queue);
+      vec4dSYCL *d_v2 = sycl::malloc_device<vec4dSYCL>(N, queue);
+      double *d_invMasses = sycl::malloc_device<double>(N, queue);
+
+      queue.memcpy(d_v1, v1, N * sizeof(vec4dSYCL));
+      queue.memcpy(d_v2, v2, N * sizeof(vec4dSYCL));
+      queue.wait();
+
+      queue.submit([&](sycl::handler &cgh) {
+         cgh.parallel_for(sycl::range<1>(N), [=](sycl::id<1> indx) {
+            vec4dSYCL v = d_v2[indx] + d_v2[indx];
+            d_invMasses[indx] = v.M();
+         });
+      });
+
+      queue.wait();
+      queue.memcpy(invMasses, d_invMasses, N * sizeof(double));
+      queue.wait();
+   }
+
+   return invMasses;
+}
+
+double *InvariantMasses(int N)
+{
+   auto v1 = GenVectors<ROOT::Math::LorentzVector<ROOT::Math::PtEtaPhiM4D<double>>>(N);
+   auto v2 = GenVectors<ROOT::Math::LorentzVector<ROOT::Math::PtEtaPhiM4D<double>>>(N);
+   double *invMasses = new double[N];
+
+   for (int i = 0; i < N; i++) {
+      invMasses[i] = (v1[i] + v2[i]).M();
+   }
+
+   return invMasses;
+}
+
+int main()
+{
+   int n = 128;
+   std::cout << "testing Invariant Masses Computation \t:\n";
+   double *msycl = InvariantMassesSYCL(n);
+   double *m = InvariantMasses(n);
+
+   int iret = 0;
+
+   for (int i = 0; i < n; i++) {
+      iret += (std::abs(msycl[i] - m[i]) > 1e-5);
+   }
+
+   if (iret)
+      std::cerr << "\t\t\t\t\t FAILED\n " << std::endl;
+   else
+      std::cout << "\t\t\t\t\t OK\n " << std::endl;
+   return iret;
+}

math/genvector/test/testGenVector.cxx

@@ -814,9 +814,9 @@ int testLorentzVector()
    iret |= compare(v1.DeltaR(v2), 4.60575f);
    // Result cross-validated using:
    // TLorentzVector t1, t2;
-   // t1.SetPtEtaPhiE(1,2,3,4); t2.SetPtEtaPhiE(5,6,7,8);
+   // t1.SetPtEtaPhiM(1,2,3,4); t2.SetPtEtaPhiM(5,6,7,8);
    // t1.DeltaR(t2)
-   
+
    // test structured bindings
    auto const [x, y, z, E] = v1;
    iret |= (x != v1.X());