Merge branch 'sgde_orthogonal_adaptivity' into 'master'

Sgde orthogonal adaptivity

See merge request !25

base/src/sgpp/base/grid/generation/hashmap/HashCoarsening.hpp

@@ -15,11 +15,10 @@
 
 #include <sgpp/globaldef.hpp>
 
-#include <vector>
-#include <list>
 #include <cmath>
+#include <list>
 #include <utility>
-
+#include <vector>
 
 namespace sgpp {
 namespace base {
@@ -44,11 +43,14 @@ class HashCoarsening {
    *
    * @param storage hashmap that stores the grid points
    * @param functor a function used to determine if refinement is needed
-   * @param alpha pointer to the gridpoints' coefficients removed points must also be considered in this vector
+   * @param alpha pointer to the gridpoints' coefficients removed points must also be considered in
+   * this vector
    * @param numFirstPoints number of grid points that are regarded to be coarsened
+   * @param minIndexConsidered indices of coarsen point candidates must be higher than this
+   * parameter to be allowed to get coarsened
    */
-  void free_coarsen_NFirstOnly(GridStorage& storage, CoarseningFunctor& functor,
-                               DataVector& alpha, size_t numFirstPoints) {
+  void free_coarsen_NFirstOnly(GridStorage& storage, CoarseningFunctor& functor, DataVector& alpha,
+                               size_t numFirstPoints, size_t minIndexConsidered = 0) {
     // check if the grid has any points
     if (storage.getSize() == 0) {
       throw generation_exception("storage empty");
@@ -77,7 +79,8 @@ class HashCoarsening {
     size_t max_idx = 0;
 
     // assure that only the first numFirstPoints are checked for coarsening
-    for (size_t z = 0; z < numFirstPoints; z++) {
+    // also assure, that indices bigger than minIndexConsidered are not checked
+    for (size_t z = minIndexConsidered; z < numFirstPoints; z++) {
       GridPoint& point = storage.getPoint(z);
 
       if (point.isLeaf() && point.isInnerPoint()) {
@@ -102,10 +105,11 @@ class HashCoarsening {
     }
 
     // DEBUG : print list of removable candidates
-    // for (size_t i = 0; i < remove_num; i++)
-    // {
-    //   std::cout << "Index: " << removePoints[i].first <<
-    //   " with surplus " << removePoints[i].second << std::endl;
+    // std::cout << "list of removable candidates:\n";
+    // for (size_t i = 0; i < remove_num; i++) {
+    //   std::cout << "Index: " << removePoints[i].first << " with surplus " <<
+    //   removePoints[i].second
+    //             << std::endl;
     // }
     // std::cout << std::endl;
 
@@ -117,20 +121,19 @@ class HashCoarsening {
     // vector to save remaining points
     std::vector<size_t> remainingIndex;
 
-    // vector to stored the points that match all condition for deleting
+    // vector to store the points that match all condition for deleting
     this->deletePoints.clear();
 
     for (size_t i = 0; i < remove_num; i++) {
-      if (removePoints[i].second < initValue &&
-          removePoints[i].second <= threshold) {
+      if (removePoints[i].second < initValue && removePoints[i].second <= threshold) {
         this->deletePoints.push_back(removePoints[i].first);
       }
     }
 
     // DEBUG : print list points to delete
-    // for(std::list<size_t>::iterator iter = deletePoints.begin();
-    // iter != deletePoints.end(); iter++)
-    // {
+    // std::cout << "list of points to delete:\n";
+    // for (std::list<size_t>::iterator iter = deletePoints.begin(); iter != deletePoints.end();
+    //      iter++) {
     //   std::cout << "Index: " << *iter << std::endl;
     // }
 
@@ -162,10 +165,10 @@ class HashCoarsening {
    *
    * @param storage hashmap that stores the grid points
    * @param functor a function used to determine if refinement is needed
-   * @param alpha pointer to the gridpoints' coefficients removed points must also be considered in this vector
+   * @param alpha pointer to the gridpoints' coefficients removed points must also be considered in
+   * this vector
    */
-  void free_coarsen(GridStorage& storage, CoarseningFunctor& functor,
-                    DataVector& alpha) {
+  void free_coarsen(GridStorage& storage, CoarseningFunctor& functor, DataVector& alpha) {
     free_coarsen_NFirstOnly(storage, functor, alpha, storage.getSize());
   }
 
@@ -184,9 +187,7 @@ class HashCoarsening {
     GridPoint point;
     GridStorage::grid_map_iterator end_iter = storage.end();
 
-    for (GridStorage::grid_map_iterator iter = storage.begin();
-         iter != end_iter;
-         iter++) {
+    for (GridStorage::grid_map_iterator iter = storage.begin(); iter != end_iter; iter++) {
       point = *(iter->first);
 
       if (point.isLeaf()) {
@@ -197,9 +198,7 @@ class HashCoarsening {
     return counter;
   }
 
-  std::list<size_t> getDeletedPoints() {
-    return this->deletePoints;
-  }
+  std::list<size_t> getDeletedPoints() { return this->deletePoints; }
 
  private:
   std::list<size_t> deletePoints;

datadriven/build/pysgpp/datadriven.i

@@ -31,6 +31,7 @@
 %include "datadriven/src/sgpp/datadriven/algorithm/DBMatDensityConfiguration.hpp"
 %include "datadriven/src/sgpp/datadriven/algorithm/DBMatDMSChol.hpp"
 %include "datadriven/src/sgpp/datadriven/algorithm/DBMatDMSDenseIChol.hpp"
+%include "datadriven/src/sgpp/datadriven/algorithm/DBMatDMSOrthoAdapt.hpp"
 %ignore *::operator=;
 %include "datadriven/src/sgpp/datadriven/algorithm/DBMatOffline.hpp"
 %include "datadriven/src/sgpp/datadriven/algorithm/DBMatOfflineGE.hpp"
@@ -48,8 +49,10 @@
 %include "datadriven/src/sgpp/datadriven/algorithm/DBMatOfflineEigen.hpp"
 %ignore sgpp::datadriven::DBMatOfflineLU::DBMatOfflineLU(DBMatOfflineLU &&);
 %include "datadriven/src/sgpp/datadriven/algorithm/DBMatOfflineLU.hpp"
+%include "datadriven/src/sgpp/datadriven/algorithm/DBMatOfflineOrthoAdapt.hpp"
 %include "datadriven/src/sgpp/datadriven/algorithm/DBMatOnlineDEEigen.hpp"
 %include "datadriven/src/sgpp/datadriven/algorithm/DBMatOnlineDELU.hpp"
+%include "datadriven/src/sgpp/datadriven/algorithm/DBMatOnlineDEOrthoAdapt.hpp"
 #endif /* USE_GSL */
 
 #ifdef __AVX__

datadriven/examples/benchmark_OrthoAdapt.cpp

@@ -0,0 +1,87 @@
+// Copyright (C) 2008-today The SG++ project
+// This file is part of the SG++ project. For conditions of distribution and
+// use, please see the copyright notice provided with SG++ or at
+// sgpp.sparsegrids.org
+
+#include <sgpp/datadriven/algorithm/DBMatDMSOrthoAdapt.hpp>
+#include <sgpp/datadriven/algorithm/DBMatOfflineOrthoAdapt.hpp>
+#include <sgpp/datadriven/algorithm/DBMatOnlineDEOrthoAdapt.hpp>
+
+#include <chrono>
+#include <iostream>
+#include <vector>
+
+int main() {
+  std::cout << "ortho_adapt algorithm benchmarks: \n";
+
+  // config
+  sgpp::datadriven::DBMatDensityConfiguration config;
+  config.grid_dim_ = 4;
+  config.grid_level_ = 6;
+  config.lambda_ = 0.0001;
+  config.regularization_ = sgpp::datadriven::RegularizationType::Identity;
+  config.decomp_type_ = sgpp::datadriven::DBMatDecompostionType::OrthoAdapt;
+
+  size_t number_points_to_refine = 1;
+  size_t number_points_to_coarsen = 1;
+
+  std::cout << "dim = " << config.grid_dim_ << "\n";
+  std::cout << "lvl = " << config.grid_level_ << "\n";
+  std::cout << "lambda = " << config.lambda_ << "\n\n";
+
+  // offline phase
+  sgpp::datadriven::DBMatOfflineOrthoAdapt offline(config);
+
+  offline.buildMatrix();
+  std::cout << "initial matrix size = " << offline.getDimA();
+
+  std::cout << "\n\ndecomposition took ";
+  auto begin = std::chrono::high_resolution_clock::now();
+  offline.decomposeMatrix();
+  auto end = std::chrono::high_resolution_clock::now();
+  std::cout << std::chrono::duration_cast<std::chrono::milliseconds>(end - begin).count() << "ms"
+            << std::endl;
+
+  // online phase
+  sgpp::datadriven::DBMatOnlineDEOrthoAdapt online(offline);
+
+  size_t refine_size = offline.getDimA() + number_points_to_refine;
+  size_t coarsen_size = refine_size - number_points_to_coarsen;
+
+  // create random points to refine
+  for (size_t i = 0; i < number_points_to_refine; i++) {
+    sgpp::base::DataVector vec(refine_size);
+    for (size_t j = 0; j < refine_size; j++) {
+      double value = (static_cast<double>(rand()) / (RAND_MAX));  // values in [0, 1]
+      vec.set(j, value);
+    }
+    online.add_new_refine_point(vec);
+  }
+
+  std::cout << "\nrefining " << number_points_to_refine << " points took ";
+  begin = std::chrono::high_resolution_clock::now();
+  online.sherman_morrison_adapt(number_points_to_refine, true);
+  end = std::chrono::high_resolution_clock::now();
+  std::cout << std::chrono::duration_cast<std::chrono::milliseconds>(end - begin).count() << "ms"
+            << std::endl;
+
+  // create random indices for coarsening
+  std::vector<size_t> coarsen_indices = {};
+  for (size_t i = 0; i < number_points_to_coarsen; i++) {
+    size_t index = ((size_t)rand() % (refine_size - coarsen_size)) + offline.getDimA();
+    coarsen_indices.push_back(index);
+  }
+  // DEBUG:
+  // std::cout << "indices to coarsen are: \n";
+  // for (auto& i : coarsen_indices) {
+  //   std::cout << i << " ";
+  // }
+
+  std::cout << "\ncoarsening " << number_points_to_coarsen << " points took ";
+  online.sherman_morrison_adapt(0, false, coarsen_indices);
+  end = std::chrono::high_resolution_clock::now();
+  std::cout << std::chrono::duration_cast<std::chrono::milliseconds>(end - begin).count() << "ms"
+            << std::endl;
+
+  return 0;
+}

datadriven/examples/learnerSGDEOnOffTest.cpp

@@ -9,6 +9,7 @@
 #include <sgpp/datadriven/application/LearnerSGDEOnOff.hpp>
 #include <sgpp/datadriven/tools/ARFFTools.hpp>
 
+#include <chrono>
 #include <string>
 
 using sgpp::base::DataMatrix;
@@ -33,6 +34,7 @@ using sgpp::base::DataVector;
  */
 
 int main() {
+  auto begin = std::chrono::high_resolution_clock::now();
   /**
    * Specify the number of runs to perform.
    * If only one specific example should be executed, set
@@ -99,28 +101,35 @@ int main() {
 
       /**
        * Select the desired decomposition type for the offline step.
-       * Note: Refinement/Coarsening only possible for Cholesky decomposition.
+       * Note: Refinement/Coarsening only possible for Cholesky decomposition
+       * and OrthoAdapt
        */
       sgpp::datadriven::DBMatDecompostionType dt;
       std::string decompType;
       // choose "LU decomposition"
       // dt = DBMatDecompostionType::DBMatDecompLU;
       // decompType = "LU decomposition";
+
       // choose"Eigen decomposition"
       // dt = DBMatDecompostionType::DBMatDecompEigen;
       // decompType = "Eigen decomposition";
+
       // choose "Cholesky decomposition"
-      //      dt = sgpp::datadriven::DBMatDecompostionType::Chol;
-      //      decompType = "Cholesky decomposition";
-      //      dt = sgpp::datadriven::DBMatDecompostionType::IChol;
-      //      decompType = "Incomplete Cholesky decomposition";
-      dt = sgpp::datadriven::DBMatDecompostionType::DenseIchol;
-      decompType = "Incomplete Cholesky decomposition on Dense Matrix";
+      // dt = sgpp::datadriven::DBMatDecompostionType::Chol;
+      // decompType = "Cholesky decomposition";
+      // dt = sgpp::datadriven::DBMatDecompostionType::IChol;
+      // decompType = "Incomplete Cholesky decomposition";
+      // dt = sgpp::datadriven::DBMatDecompostionType::DenseIchol;
+
+      // choose "orthogonal Adaptivity"
+      dt = sgpp::datadriven::DBMatDecompostionType::OrthoAdapt;
+      decompType = "orthogonal Adaptivity";
+
       std::cout << "Decomposition type: " << decompType << std::endl;
 
       /**
-       * Configure adaptive refinement (if Cholesky is chosen). As refinement
-       * monitor the periodic monitor or the convergence monitor
+       * Configure adaptive refinement (if Cholesky or OrthoAdapt is chosen).
+       * As refinement monitor the periodic monitor or the convergence monitor
        * can be chosen. Possible refinement indicators are
        * surplus refinement, data-based refinement, zero-crossings-based
        * refinement.
@@ -146,19 +155,19 @@ int main() {
       std::cout << "Refinement monitor: " << refMonitor << std::endl;
       std::string refType;
       // select surplus refinement
-      // refType = "surplus";
+      refType = "surplus";
       // select data-based refinement
       // refType = "data";
       // select zero-crossings-based refinement
-      refType = "zero";
+      // refType = "zero";
       std::cout << "Refinement type: " << refType << std::endl;
       sgpp::base::AdpativityConfiguration adaptConfig;
       /**
        * Specify number of refinement steps and the max number
        * of grid points to refine each step.
        */
-      adaptConfig.numRefinements_ = 2;
-      adaptConfig.noPoints_ = 7;
+      adaptConfig.numRefinements_ = 10;
+      adaptConfig.noPoints_ = 10;
       adaptConfig.threshold_ = 0.0;  // only required for surplus refinement
 
       // initial regularization parameter lambda
@@ -200,7 +209,7 @@ int main() {
 
       // specify batch size
       // (set to 1 for processing only a single data point each iteration)
-      size_t batchSize = 1;
+      size_t batchSize = 10;
       // specify max number of passes over traininig data set
       size_t maxDataPasses = 2;
 
@@ -238,16 +247,20 @@ int main() {
     avgErrorsFolds.mult(1.0 / static_cast<double>(totalFolds));
 
     // write error evaluation to csv-file
-    /*std::ofstream output;
-    output.open("SGDEOnOff_avg_classification_error_"+std::to_string(numSets+1)+".csv");
+    /*
+    std::ofstream output;
+    output.open("SGDEOnOff_avg_classification_error_" + std::to_string(numSets + 1) + ".csv");
     if (output.fail()) {
       std::cout << "failed to create csv file!" << std::endl;
-    }
-    else {
+    } else {
       for (size_t i = 0; i < avgErrorsFolds.getSize(); i++) {
         output << avgErrorsFolds.get(i) << ";" << std::endl;
       }
       output.close();
     }*/
   }
+  auto end = std::chrono::high_resolution_clock::now();
+  std::cout << "whole learnerSGDEOnOff test took "
+            << std::chrono::duration_cast<std::chrono::milliseconds>(end - begin).count() << "ms"
+            << std::endl;
 }