perf: avoid unnecessary copies, fix: small coding style issues

gpu/benchmark/Benchmark.hpp

@@ -6,16 +6,36 @@
 #include <psDomain.hpp>
 #include <rayParticle.hpp>
 
-#define MAKE_GEO Trench
+#define MAKE_GEO Hole
 #define DEFAULT_GRID_DELTA 0.1
 #define DEFAULT_STICKING 0.1
-#define DIM 2
+#define DIM 3
 #define FIXED_RAYS false
 
-constexpr int particleType = 1;
 using TranslatorType = std::unordered_map<unsigned long, unsigned long>;
 using namespace viennaps;
 
+template <class NumericType, int N>
+consteval std::array<NumericType, N> linspace(NumericType start,
+                                              NumericType end) {
+  std::array<NumericType, N> arr{};
+  NumericType step = (end - start) / static_cast<NumericType>(N - 1);
+  for (int i = 0; i < N; ++i) {
+    arr[i] = start + i * step;
+  }
+  return arr;
+}
+
+constexpr int particleType = 0;
+constexpr bool runDisk = false;
+constexpr bool runLine = false;
+constexpr bool runTriangle = true;
+
+constexpr auto gridDeltaValues = linspace<float, 8>(0.08f, 0.4f);
+constexpr int numRuns = 10;
+constexpr int raysPerPoint = 1000;
+constexpr int numRays = int(1.4e8);
+
 template <class NumericType>
 auto Trench(NumericType gridDelta = DEFAULT_GRID_DELTA) {
   NumericType xExtent = 20.;
@@ -42,16 +62,6 @@ auto Hole(NumericType gridDelta = DEFAULT_GRID_DELTA) {
   return domain;
 }
 
-template <class NumericType, int N>
-std::array<NumericType, N> linspace(NumericType start, NumericType end) {
-  std::array<NumericType, N> arr{};
-  NumericType step = (end - start) / static_cast<NumericType>(N - 1);
-  for (int i = 0; i < N; ++i) {
-    arr[i] = start + i * step;
-  }
-  return arr;
-}
-
 template <typename NumericType> auto getIBEParameters() {
   viennaps::IBEParameters<NumericType> params;
   params.redepositionRate = 1.0;
@@ -199,8 +209,9 @@ void setupTriangleGeometry(
     elementKdTree_->setPoints(triangleCenters);
     elementKdTree_->build();
   } else {
-    triangleMesh = CreateTriangleMesh(
-        static_cast<float>(domain->getGridDelta()), surfaceMesh_);
+    CopyTriangleMesh(static_cast<float>(domain->getGridDelta()), surfaceMesh_,
+                     triangleMesh);
+    // preserves surfaceMesh_
   }
 
   rayTracer_.setGeometry(triangleMesh);

gpu/benchmark/CPU_Benchmark.cpp

@@ -7,27 +7,18 @@
 
 #include "Benchmark.hpp"
 
-using namespace viennaps;
-
 int main() {
   omp_set_num_threads(16);
   using NumericType = float;
   constexpr int D = DIM;
 
-  // const std::array<NumericType, 4> gridDeltaValues = {0.05f, 0.1f, 0.2f,
-  // 0.4f};
-  auto gridDeltaValues = linspace<NumericType, 8>(0.01f, 0.4f);
-  const int numRuns = 20;
-  const int raysPerPoint = 1000;
-  const int numRays = int(1.4e8);
-
   auto particle = makeCPUParticle<NumericType, D>();
   std::vector<std::unique_ptr<viennaray::AbstractParticle<NumericType>>>
       particles;
   particles.push_back(particle->clone());
   std::string fluxLabel = particleType == 0 ? "flux" : "ionFlux";
 
-  { // Disk
+  if constexpr (runDisk) { // Disk
     std::ofstream file("CPU_Benchmark_Disk.txt");
     file << "Meshing;Tracing;Postprocessing;GridDelta\n";
 
@@ -112,7 +103,7 @@ int main() {
     file.close();
   }
 
-  { // Triangle
+  if constexpr (runTriangle) { // Triangle
     std::ofstream file("CPU_Benchmark_Triangle.txt");
     file << "Meshing;Tracing;Postprocessing;GridDelta\n";
 
@@ -125,7 +116,7 @@ int main() {
     tracer.setUseRandomSeeds(false);
     tracer.setParticleType(particle);
 
-    std::vector<std::string> dataLabels = particle->getLocalDataLabels();
+    const auto &dataLabels = particle->getLocalDataLabels();
 
     for (int i = 0; i < gridDeltaValues.size(); i++) {
       std::cout << "  Grid Delta: " << gridDeltaValues[i] << "\n";

gpu/benchmark/GPU_Benchmark.cpp

@@ -13,21 +13,10 @@
 
 #include "Benchmark.hpp"
 
-using namespace viennaps;
-
 int main() {
   omp_set_num_threads(16);
   using NumericType = float;
   constexpr int D = DIM;
-
-  const NumericType cSticking = DEFAULT_STICKING;
-  // const std::array<NumericType, 4> gridDeltaValues = {0.05f, 0.1f, 0.2f,
-  // 0.4f};
-  auto gridDeltaValues = linspace<NumericType, 8>(0.01f, 0.4f);
-  const int numRuns = 20;
-  const int raysPerPoint = 1000;
-  const int numRays = int(1.4e8);
-
   auto context = DeviceContext::createContext();
 
   CudaBuffer deviceParamsBuffer;
@@ -36,7 +25,7 @@ int main() {
     deviceParamsBuffer.allocUploadSingle(deviceParams);
   }
 
-  { // Triangle
+  if constexpr (runTriangle) { // Triangle
     std::ofstream file("GPU_Benchmark_Triangle.txt");
     file << "Meshing;Tracing;Postprocessing;GridDelta\n";
 
@@ -55,8 +44,7 @@ int main() {
     }
     tracer.prepareParticlePrograms();
 
-    std::vector<std::string> dataLabels =
-        std::get<0>(particleConfig).dataLabels;
+    const auto &dataLabels = std::get<0>(particleConfig).dataLabels;
 
     std::cout << "Starting Triangle Benchmark\n";
 
@@ -138,7 +126,7 @@ int main() {
     file.close();
   }
 
-  { // Disk
+  if constexpr (runDisk) { // Disk
     std::ofstream file("GPU_Benchmark_Disk.txt");
     file << "Meshing;Tracing;Postprocessing;GridDelta\n";
 
@@ -186,7 +174,6 @@ int main() {
       for (int j = 0; j < numRuns; j++) {
         std::cout << "    Process Step: " << j + 1 << "\n";
         advectionKernel.prepareLS();
-        file << cSticking << ";";
 
         Timer timer;
 
@@ -235,7 +222,7 @@ int main() {
     file.close();
   }
 
-  if constexpr (D == 2) { // Line
+  if constexpr (D == 2 && runLine) { // Line
     std::ofstream file("GPU_Benchmark_Line.txt");
     file << "Meshing;Tracing;Postprocessing;GridDelta\n";
 

gpu/benchmark/holeCompare.cpp

@@ -23,7 +23,7 @@ void runBenchmark(SmartPointer<Domain<NumericType, D>> geometry,
   process.apply();
   timer.finish();
 
-  std::cout << "Flux Engine: " << to_string(fluxEngine)
+  std::cout << "Flux Engine: " << util::toString(fluxEngine)
             << ", Time taken: " << timer.currentDuration / 1e9 << " s"
             << std::endl;
   copy->saveSurfaceMesh(outputFilename);
@@ -61,7 +61,7 @@ int main(int argc, char *argv[]) {
 
     for (const auto &fluxEngine : fluxEngines) {
       runBenchmark<NumericType, D>(geometry, model, fluxEngine,
-                                   "result_PE_" + to_string(fluxEngine) +
+                                   "result_PE_" + util::toString(fluxEngine) +
                                        ".vtp");
     }
   }
@@ -72,7 +72,7 @@ int main(int argc, char *argv[]) {
 
     for (const auto &fluxEngine : fluxEngines) {
       runBenchmark<NumericType, D>(geometry, model, fluxEngine,
-                                   "result_SP_" + to_string(fluxEngine) +
+                                   "result_SP_" + util::toString(fluxEngine) +
                                        ".vtp");
     }
   }
@@ -90,7 +90,7 @@ int main(int argc, char *argv[]) {
 
     for (const auto &fluxEngine : fluxEngines) {
       runBenchmark<NumericType, D>(geometry, model, fluxEngine,
-                                   "result_IBE_" + to_string(fluxEngine) +
+                                   "result_IBE_" + util::toString(fluxEngine) +
                                        ".vtp");
     }
   }
@@ -114,7 +114,7 @@ int main(int argc, char *argv[]) {
 
     for (const auto &fluxEngine : fluxEngines) {
       runBenchmark<NumericType, D>(geometry, model, fluxEngine,
-                                   "result_MP_" + to_string(fluxEngine) +
+                                   "result_MP_" + util::toString(fluxEngine) +
                                        ".vtp");
     }
   }

gpu/benchmark/holeCompareGrid.cpp

@@ -28,7 +28,7 @@ void runBenchmark(SmartPointer<Domain<NumericType, D>> geometry,
   process.apply();
   timer.finish();
 
-  std::cout << "Flux Engine: " << to_string(fluxEngine)
+  std::cout << "Flux Engine: " << util::toString(fluxEngine)
             << ", Time taken: " << timer.currentDuration / 1e9 << " s"
             << std::endl;
   copy->saveSurfaceMesh(outputFilename);
@@ -65,7 +65,7 @@ int main(int argc, char *argv[]) {
 
       for (const auto &fluxEngine : fluxEngines) {
         runBenchmark<NumericType, D>(geometry, model, fluxEngine,
-                                     "result_PE_" + to_string(fluxEngine) +
+                                     "result_PE_" + util::toString(fluxEngine) +
                                          "_" + std::to_string(i) + ".vtp",
                                      0.2);
       }
@@ -92,8 +92,9 @@ int main(int argc, char *argv[]) {
 
       for (const auto &fluxEngine : fluxEngines) {
         runBenchmark<NumericType, D>(geometry, model, fluxEngine,
-                                     "result_IBE_" + to_string(fluxEngine) +
-                                         "_" + std::to_string(i) + ".vtp");
+                                     "result_IBE_" +
+                                         util::toString(fluxEngine) + "_" +
+                                         std::to_string(i) + ".vtp");
       }
     }
   }

include/viennaps/process/psCPUTriangleEngine.hpp

@@ -114,8 +114,8 @@ class CPUTriangleEngine final : public FluxEngine<NumericType, D> {
       elementKdTree_->setPoints(triangleCenters);
       elementKdTree_->build();
     } else {
-      triangleMesh = CreateTriangleMesh(
-          static_cast<float>(context.domain->getGridDelta()), surfaceMesh_);
+      CopyTriangleMesh(static_cast<float>(context.domain->getGridDelta()),
+                       surfaceMesh_, triangleMesh);
     }
 
     rayTracer_.setGeometry(triangleMesh);

include/viennaps/process/psGPUTriangleEngine.hpp

@@ -136,8 +136,9 @@ class GPUTriangleEngine final : public FluxEngine<NumericType, D> {
       elementKdTree_->setPoints(triangleCenters);
       elementKdTree_->build();
     } else {
-      triangleMesh = CreateTriangleMesh(
-          static_cast<float>(context.domain->getGridDelta()), surfaceMesh_);
+      CopyTriangleMesh(static_cast<float>(context.domain->getGridDelta()),
+                       surfaceMesh_, triangleMesh);
+      // preserves surfaceMesh_
     }
 
     rayTracer_.setGeometry(triangleMesh);
@@ -228,7 +229,8 @@ class GPUTriangleEngine final : public FluxEngine<NumericType, D> {
         downloadCoverages(d_coverages, surfaceMesh_->getCellData(), coverages,
                           surfaceMesh_->getElements<3>().size());
       }
-      saveResultsToPointData(surfaceMesh_->getCellData());
+      saveResultsToPointData(
+          surfaceMesh_->getCellData()); // save fluxes to elements
       viennals::VTKWriter<float>(
           surfaceMesh_, context.intermediateOutputPath +
                             context.getProcessName() + "_flux_" +
@@ -269,22 +271,20 @@ class GPUTriangleEngine final : public FluxEngine<NumericType, D> {
     assert(numPoints == surfaceMesh_->getElements<3>().size());
     auto const &results = rayTracer_.getResults();
     auto particles = rayTracer_.getParticles();
+    const auto &dataLabels = model_->getParticleDataLabels();
 
-    int offset = 0;
-    for (int pIdx = 0; pIdx < particles.size(); pIdx++) {
-      for (int dIdx = 0; dIdx < particles[pIdx].dataLabels.size(); dIdx++) {
-        auto name = particles[pIdx].dataLabels[dIdx];
-        assert(offset + dIdx < results.size());
-        const auto &data = results[offset + dIdx];
+    assert(dataLabels.size() == results.size());
 
-        std::vector<float> values(numPoints);
-        for (unsigned i = 0; i < numPoints; ++i) {
-          values[i] = static_cast<float>(data[i]);
-        }
+    for (int dIdx = 0; dIdx < dataLabels.size(); dIdx++) {
+      const auto &name = dataLabels[dIdx];
+      const auto &data = results[dIdx];
 
-        pointData.insertReplaceScalarData(std::move(values), name);
+      std::vector<float> values(numPoints);
+      for (unsigned i = 0; i < numPoints; ++i) {
+        values[i] = static_cast<float>(data[i]);
       }
-      offset += particles[pIdx].dataLabels.size();
+
+      pointData.insertReplaceScalarData(std::move(values), name);
     }
   }
 

include/viennaps/process/psProcess.hpp

@@ -99,7 +99,7 @@ template <typename NumericType, int D> class Process {
     VIENNACORE_LOG_DEBUG("Using strategy: " + std::string(strategy->name()));
 
     if (strategy->requiresFluxEngine()) {
-      VIENNACORE_LOG_DEBUG("Setting up " + to_string(fluxEngineType_) +
+      VIENNACORE_LOG_DEBUG("Setting up " + util::toString(fluxEngineType_) +
                            " flux engine for strategy.");
       strategy->setFluxEngine(createFluxEngine());
     }
@@ -191,7 +191,7 @@ template <typename NumericType, int D> class Process {
     assert(fluxEngineType_ != FluxEngineType::AUTO &&
            "Flux engine type must be specified before creation.");
     VIENNACORE_LOG_DEBUG("Creating flux engine of type: " +
-                         to_string(fluxEngineType_));
+                         util::toString(fluxEngineType_));
     // Create CPU engine
     if (fluxEngineType_ == FluxEngineType::CPU_DISK) {
       return std::make_unique<CPUDiskEngine<NumericType, D>>();
@@ -295,7 +295,7 @@ template <typename NumericType, int D> class Process {
         fluxEngineType_ = FluxEngineType::CPU_DISK;
       }
       VIENNACORE_LOG_DEBUG("Auto-selected flux engine type: " +
-                           to_string(fluxEngineType_));
+                           util::toString(fluxEngineType_));
     }
 
 #ifdef VIENNACORE_COMPILE_GPU

include/viennaps/process/psProcessModel.hpp

@@ -174,8 +174,6 @@ class ProcessModelCPU : public ProcessModelBase<NumericType, D> {
 #ifdef VIENNACORE_COMPILE_GPU
 namespace viennaps::gpu {
 
-using namespace viennacore;
-
 template <class NumericType, int D>
 class ProcessModelGPU : public ProcessModelBase<NumericType, D> {
 private:
@@ -187,7 +185,7 @@ class ProcessModelGPU : public ProcessModelBase<NumericType, D> {
   std::vector<viennaray::gpu::CallableConfig> callableMap_;
 
 public:
-  CudaBuffer processData;
+  viennacore::CudaBuffer processData;
   auto &getParticleTypes() { return particles; }
   auto getProcessDataDPtr() const { return processData.dPointer(); }
   bool useMaterialIds() const { return materialIds; }

include/viennaps/psCreateSurfaceMesh.hpp

@@ -24,7 +24,18 @@ CreateTriangleMesh(const float gridDelta,
   triangleMesh.maximumExtent = mesh->maximumExtent;
   triangleMesh.normals = *mesh->getCellData().getVectorData("Normals");
 
-  return triangleMesh;
+  return std::move(triangleMesh);
+}
+
+inline void CopyTriangleMesh(const float gridDelta,
+                             const SmartPointer<viennals::Mesh<float>> &mesh,
+                             viennaray::TriangleMesh &triangleMesh) {
+  triangleMesh.gridDelta = gridDelta;
+  triangleMesh.triangles = mesh->triangles;
+  triangleMesh.nodes = mesh->nodes;
+  triangleMesh.minimumExtent = mesh->minimumExtent;
+  triangleMesh.maximumExtent = mesh->maximumExtent;
+  triangleMesh.normals = *mesh->getCellData().getVectorData("Normals");
 }
 
 template <class LsNT, class MeshNT = LsNT, int D = 3> class CreateSurfaceMesh {