Yuhao smp pca

benchmark/scripts/scanThreads/config_CPPSMPPCA.csv

@@ -0,0 +1,9 @@
+key,value,type
+aRow,10000,U64
+aCol,1000,U64
+bCol,500,U64
+sketchDimension,25,U64
+threads,1,U64
+ptFile,,String
+useCPP,1,U64
+cppAlgoTag,smp-pca,String

benchmark/scripts/scanThreads/config_CPPWeightedCR.csv

@@ -0,0 +1,9 @@
+key,value,type
+aRow,10000,U64
+aCol,1000,U64
+bCol,500,U64
+sketchDimension,25,U64
+threads,1,U64
+ptFile,torchscripts/WeightedCR.pt,String
+useCPP,1,U64
+cppAlgoTag,weighted-cr,String

benchmark/scripts/scanThreads/drawTogether.py

@@ -137,9 +137,9 @@ def main():
     figPath = os.path.abspath(os.path.join(os.getcwd(), "../..")) + "/figures/" + scanTag
     configTemplate = exeSpace + "config.csv"
     commonBase = resultPath + "/"
-    resultPaths = ["mm", "crs","ews","counterSketch","cofd","tow"]
-    csvTemplates = ["config_CPPMM.csv", "config_CPPCRS.csv","config_CPPEWS.csv","config_CPPCOUNTERSKETCH.csv","config_CPPCOFD.csv","config_CPPTOW.csv"]
-    methodTags = ["mm", "crs","ews","counterSketch","cofd","tow"]
+    resultPaths = ["mm", "crs","ews","counterSketch","cofd","tow","weighted-cr","smp-pca"]
+    csvTemplates = ["config_CPPMM.csv", "config_CPPCRS.csv","config_CPPEWS.csv","config_CPPCOUNTERSKETCH.csv","config_CPPCOFD.csv","config_CPPTOW.csv","config_CPPWeightedCR.csv","config_CPPSMPPCA.csv"]
+    methodTags = ["mm", "crs","ews","counterSketch","cofd","tow","weighted-cr","smp-pca"]
     valueVec = [ 2, 4, 6, 8, 10, 12]
     valueVecRun = valueVec
     print(configTemplate)

include/AMMBench.h

@@ -126,6 +126,8 @@
 #include <CPPAlgos/INT8CPPAlgo.h>
 
 #include <CPPAlgos/TugOfWarCPPAlgo.h>
+#include <CPPAlgos/WeightedCRCPPAlgo.h>
+#include <CPPAlgos/SMPPCACPPAlgo.h>
 /**
  * @}
  *

include/CPPAlgos/SMPPCACPPAlgo.h

@@ -0,0 +1,56 @@
+/*! \file SMPPCACPPAlgo.h*/
+//
+// Created by yuhao on 6/6/23.
+//
+
+#ifndef INTELLISTREAM_INCLUDE_CPPALGOS_SMPPCACppAlgo_H_
+#define INTELLISTREAM_INCLUDE_CPPALGOS_SMPPCACppAlgo_H_
+#include <CPPAlgos/AbstractCPPAlgo.h>
+namespace AMMBench {
+/**
+ * @ingroup AMMBENCH_CppAlgos The algorithms written in c++
+ * @{
+ */
+/**
+ * @class SMPPCACPPAlgo CPPAlgos/SMPPCACPPAlgo.h
+ * @brief sketch scaled JL class of c++ algos
+ *
+ */
+class SMPPCACPPAlgo : public AMMBench::AbstractCPPAlgo {
+ public:
+  SMPPCACPPAlgo() {
+
+  }
+  ~SMPPCACPPAlgo() {
+
+  }
+  /**
+   * @brief the virtual function provided for outside callers, rewrite in children classes
+   * @param A the A matrix
+   * @param B the B matrix
+   * @param sketchSize the size of sketch
+   * @return the output c matrix
+   */
+  virtual torch::Tensor amm(torch::Tensor A, torch::Tensor B, uint64_t sketchSize);
+
+};
+/**
+ * @ingroup AMMBENCH_CppAlgos
+ * @typedef AbstractMatrixCppAlgoPtr
+ * @brief The class to describe a shared pointer to @ref SMPPCACppAlgo
+
+ */
+typedef std::shared_ptr<class AMMBench::SMPPCACPPAlgo> SMPPCACPPAlgoPtr;
+/**
+ * @ingroup AMMBENCH_CppAlgos
+ * @def newSMPPCACppAlgo
+ * @brief (Macro) To creat a new @ref  SMPPCACppAlgounder shared pointer.
+ */
+#define newSMPPCACPPAlgo std::make_shared<AMMBench::SMPPCACPPAlgo>
+}
+/**
+ * @}
+ */
+
+#endif //INTELLISTREAM_INCLUDE_CPPALGOS_SMPPCACppAlgo_H_
+

include/CPPAlgos/WeightedCRCPPAlgo.h

@@ -0,0 +1,56 @@
+/*! \file WeightedCRCPPAlgo.h*/
+//
+// Created by yuhao on 27/05/23.
+//
+
+#ifndef INTELLISTREAM_INCLUDE_CPPALGOS_WeightedCRCPPAlgo_H_
+#define INTELLISTREAM_INCLUDE_CPPALGOS_WeightedCRCPPAlgo_H_
+#include <CPPAlgos/AbstractCPPAlgo.h>
+namespace AMMBench {
+/**
+ * @ingroup AMMBENCH_CPPAlgos The algorithms writtrn in c++
+ * @{
+ */
+/**
+ * @class WeightedCRCPPlgo CPPAlgos/WeightedCRCPPAlgo.h
+ * @brief The weighted cloumn row sampling class of c++ algos
+ *
+ */
+class WeightedCRCPPAlgo : public AMMBench::AbstractCPPAlgo {
+ public:
+  WeightedCRCPPAlgo() {
+
+  }
+  ~WeightedCRCPPAlgo() {
+
+  }
+  /**
+   * @brief the virtual function provided for outside callers, rewrite in children classes
+   * @param A the A matrix
+   * @param B the B matrix
+   * @param sketchSize the size of sketc or sampling
+   * @return the output c matrix
+   */
+  virtual torch::Tensor amm(torch::Tensor A, torch::Tensor B, uint64_t sketchSize);
+
+};
+/**
+ * @ingroup AMMBENCH_CPPAlgos
+ * @typedef AbstractMatrixCPPAlgoPtr
+ * @brief The class to describe a shared pointer to @ref WeightedCRCPPAlgo
+
+ */
+typedef std::shared_ptr<class AMMBench::WeightedCRCPPAlgo> WeightedCRCPPAlgoPtr;
+/**
+ * @ingroup AMMBENCH_CPPAlgos
+ * @def newWeightedCRCPPAlgo
+ * @brief (Macro) To creat a new @ref  WeightedCRCPPAlgounder shared pointer.
+ */
+#define newWeightedCRCPPAlgo std::make_shared<AMMBench::WeightedCRCPPAlgo>
+}
+/**
+ * @}
+ */
+
+#endif //INTELLISTREAM_INCLUDE_CPPALGOS_WeightedCRCPPAlgo_H_
+

src/CPPAlgos/CMakeLists.txt

@@ -1,6 +1,7 @@
 add_sources(
         AbstractCPPAlgo.cpp
         CRSCPPAlgo.cpp
+        WeightedCRCPPAlgo.cpp
         CPPAlgoTable.cpp
         CRSV2CPPAlgo.cpp
         BCRSCPPAlgo.cpp
@@ -10,5 +11,6 @@ add_sources(
         CountSketchCPPAlgo.cpp
         INT8CPPAlgo.cpp
         TugOfWarCPPAlgo.cpp
+        SMPPCACPPAlgo.cpp
 )
 

src/CPPAlgos/CPPAlgoTable.cpp

@@ -12,9 +12,9 @@
 #include <CPPAlgos/EWSCPPAlgo.h>
 #include <CPPAlgos/CoOccurringFDCPPAlgo.h>
 #include <CPPAlgos/BetaCoOFDCPPAlgo.h>
-
+#include <CPPAlgos/WeightedCRCPPAlgo.h>
+#include <CPPAlgos/SMPPCACPPAlgo.h>
 #include <CPPAlgos/INT8CPPAlgo.h>
-
 #include <CPPAlgos/TugOfWarCPPAlgo.h>
 namespace AMMBench {
 AMMBench::CPPAlgoTable::CPPAlgoTable() {
@@ -28,6 +28,8 @@ AMMBench::CPPAlgoTable::CPPAlgoTable() {
   algoMap["bcooFD"] = newBetaCoOFDCPPAlgo();
   algoMap["int8"] = newINT8CPPAlgo();
   algoMap["tugOfWar"] = newTugOfWarCPPAlgo();
+  algoMap["weighted-cr"] = newWeightedCRCPPAlgo();
+  algoMap["smp-pca"] = newSMPPCACPPAlgo();
 }
 
 } // AMMBench

src/CPPAlgos/SMPPCACPPAlgo.cpp

@@ -0,0 +1,36 @@
+//
+// Created by yuhao on 6/6/23.
+//
+
+#include <CPPAlgos/SMPPCACPPAlgo.h>
+
+namespace AMMBench {
+torch::Tensor AMMBench::SMPPCACPPAlgo::amm(torch::Tensor A, torch::Tensor B, uint64_t k) {
+   // Step 1: Input
+    A = A.t(); // d*n1
+    int64_t d = A.size(0);
+    int64_t n1 = A.size(1);
+    int64_t n2 = B.size(1);
+
+    // Step 2: Get sketched matrix
+    torch::Tensor pi = 1/std::sqrt(k) * torch::randn({k, d}); // Gaussian sketching matrix
+    torch::Tensor A_tilde = torch::matmul(pi, A); // k*n1
+    torch::Tensor B_tilde = torch::matmul(pi, B); // k*n2
+
+    torch::Tensor A_tilde_B_tilde = torch::matmul(A_tilde.t(), B_tilde);
+
+    // Compute column norms of A and B
+    torch::Tensor col_norm_A = torch::linalg::vector_norm(A, 2, {0}, false, c10::nullopt); // ||Ai|| for i in [n1]
+    torch::Tensor col_norm_B = torch::linalg::vector_norm(B, 2, {0}, false, c10::nullopt); // ||Bj|| for j in [n2]
+
+    torch::Tensor col_norm_A_tilde = torch::linalg::vector_norm(A_tilde, 2, {0}, false, c10::nullopt); // ||Ai|| for i in [n1]
+    torch::Tensor col_norm_B_tilde = torch::linalg::vector_norm(B_tilde, 2, {0}, false, c10::nullopt); // ||Bj|| for j in [n2]
+
+    // Compute M_tilde
+    torch::Tensor col_norm_A_col_norm_B = torch::matmul(col_norm_A.reshape({n1,1}), col_norm_B.reshape({1,n2}));
+    torch::Tensor col_norm_A_tilde_col_norm_B_tilde = torch::matmul(col_norm_A_tilde.reshape({n1,1}), col_norm_B_tilde.reshape({1,n2}));
+    torch::Tensor M_tilde = torch::div(torch::mul(A_tilde_B_tilde, col_norm_A_col_norm_B), col_norm_A_tilde_col_norm_B_tilde);
+
+    return M_tilde;
+}
+} // AMMBench

src/CPPAlgos/WeightedCRCPPAlgo.cpp

@@ -0,0 +1,44 @@
+//
+// Created by tony on 25/05/23.
+//
+
+#include <CPPAlgos/WeightedCRCPPAlgo.h>
+#include <ATen/ATen.h>
+
+namespace AMMBench {
+torch::Tensor AMMBench::WeightedCRCPPAlgo::amm(torch::Tensor A, torch::Tensor B, uint64_t c) {
+
+  int64_t n = A.size(1); // A: m*n, B: n*d
+
+  // Probability distribution
+  torch::Tensor probability_distribution = torch::zeros({n});
+    for (int i = 0; i < n; ++i) {
+      probability_distribution[i] = torch::norm(A.t()[i]) * torch::norm(B[i]);
+    }
+    probability_distribution /= probability_distribution.sum();
+
+  // S
+  torch::Tensor sample_indices = torch::multinomial(probability_distribution, /*num_samples*/c, /*replacement*/true);
+
+  torch::Tensor unique_indices, _, occurences;
+  std::tie(unique_indices, _, occurences) = at::_unique2(sample_indices, /*sorted*/false, /*return_inverse*/false, /*return_counts*/true);
+
+  torch::Tensor S = torch::zeros({n, unique_indices.size(0)});
+
+  for (int trial = 0; trial < unique_indices.size(0); ++trial) {
+    int index = unique_indices[trial].item<int>();
+    S[index][trial] = 1;
+  }
+
+  // D
+  torch::Tensor D = torch::diag(torch::sqrt(occurences) / torch::sqrt(torch::tensor({static_cast<long>(c)}, torch::kLong) * probability_distribution.index_select(0, unique_indices)));
+
+  // ASD(SD)^TB
+  torch::Tensor SS = torch::matmul(S, D);
+  torch::Tensor C = torch::matmul(A, SS);
+  torch::Tensor R = torch::matmul(SS.t(), B);
+  torch::Tensor weighted_CR = torch::matmul(C, R);
+
+  return weighted_CR;
+}
+} // AMMBench

test/CMakeLists.txt

@@ -23,6 +23,7 @@ add_catch_test(sketch_test SystemTest/SketchTest.cpp IntelliStream)
 add_catch_test(crs_test SystemTest/CRSTest.cpp IntelliStream)
 add_catch_test(ews_test SystemTest/EWSTest.cpp IntelliStream)
 add_catch_test(weighted_cr_test SystemTest/WeightedCRTest.cpp IntelliStream)
+add_catch_test(smp_pca_test SystemTest/SMPPCATest.cpp IntelliStream)
 add_catch_test(block_partition_test SystemTest/BlockPartitionTest.cpp IntelliStream)
 add_catch_test(tug_of_war_test SystemTest/TugOfWarTest.cpp IntelliStream)
 add_catch_test(int8_test SystemTest/INT8Test.cpp IntelliStream)

test/SystemTest/SMPPCATest.cpp

@@ -0,0 +1,21 @@
+#include <vector>
+
+#define CATCH_CONFIG_MAIN
+#include "catch.hpp"
+#include <AMMBench.h>
+using namespace std;
+using namespace INTELLI;
+using namespace torch;
+
+TEST_CASE("Test CRS in cpp", "[short]")
+{
+  torch::manual_seed(114514);
+  AMMBench::SMPPCACPPAlgo smppca;
+  auto A = torch::rand({600, 400});
+  auto B = torch::rand({400, 1000});
+  auto realC = torch::matmul(A, B);
+  auto ammC = smppca.amm(A, B, 20);
+  double froError = INTELLI::UtilityFunctions::relativeFrobeniusNorm(realC, ammC);
+  std::cout << froError << std::endl;
+  // REQUIRE(froError < 0.5);
+}

test/SystemTest/WeightedCRTest.cpp

@@ -53,4 +53,15 @@ TEST_CASE("Test the COLUMN ROW SAMPLINGS", "[short]")
   runSingleThreadTest("scripts/config_WeightedCR.csv");
   // place your test here
   REQUIRE(a == 0);
+}
+TEST_CASE("Test Weighted CR in cpp", "[short]")
+{
+    torch::manual_seed(114514);
+    AMMBench::WeightedCRCPPAlgo wcr;
+    auto A = torch::rand({400, 400});
+    auto B = torch::rand({400, 400});
+    auto realC = torch::matmul(A, B);
+    auto ammC = wcr.amm(A, B, 20);
+    double froError = INTELLI::UtilityFunctions::relativeFrobeniusNorm(realC, ammC);
+    REQUIRE(froError < 0.5);
 }

test/scripts/config_WeightedCR.csv

@@ -3,4 +3,4 @@ aRow,100,U64
 aCol,1000,U64
 bCol,500,U64
 sketchDimension,25,U64
-ptFile,torchscripts/CRS.pt,String
+ptFile,torchscripts/WeightedCR.pt,String