BenchmarkTopologyAlgorithms.cxx 14.5 KB
Newer Older
1 2 3 4 5 6 7 8
//============================================================================
//  Copyright (c) Kitware, Inc.
//  All rights reserved.
//  See LICENSE.txt for details.
//  This software is distributed WITHOUT ANY WARRANTY; without even
//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
//  PURPOSE.  See the above copyright notice for more information.
//
Kenneth Moreland's avatar
Kenneth Moreland committed
9
//  Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
10 11 12
//  Copyright 2014 UT-Battelle, LLC.
//  Copyright 2014 Los Alamos National Security.
//
Kenneth Moreland's avatar
Kenneth Moreland committed
13
//  Under the terms of Contract DE-NA0003525 with NTESS,
14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
//  the U.S. Government retains certain rights in this software.
//
//  Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
//  Laboratory (LANL), the U.S. Government retains certain rights in
//  this software.
//============================================================================
#include <vtkm/Math.h>
#include <vtkm/VectorAnalysis.h>

#include <vtkm/cont/ArrayHandle.h>
#include <vtkm/cont/CellSetStructured.h>
#include <vtkm/cont/Timer.h>

#include <vtkm/worklet/DispatcherMapField.h>
#include <vtkm/worklet/DispatcherMapTopology.h>
29 30
#include <vtkm/worklet/WorkletMapField.h>
#include <vtkm/worklet/WorkletMapTopology.h>
31

32
#include "Benchmarker.h"
33
#include <vtkm/cont/testing/Testing.h>
34

35
#include <cctype>
36
#include <random>
37 38
#include <string>

39 40 41 42
namespace vtkm
{
namespace benchmarking
{
43 44

#define CUBE_SIZE 256
45
static const std::string DIVIDER(40, '-');
46

47 48
enum BenchmarkName
{
49 50 51
  CELL_TO_POINT = 1 << 1,
  POINT_TO_CELL = 1 << 2,
  MC_CLASSIFY = 1 << 3,
52
  ALL = CELL_TO_POINT | POINT_TO_CELL | MC_CLASSIFY
53 54 55 56 57
};

class AveragePointToCell : public vtkm::worklet::WorkletMapPointToCell
{
public:
58
  using ControlSignature = void(FieldInPoint<> inPoints,
59
                                CellSetIn cellset,
60
                                FieldOutCell<> outCells);
61
  using ExecutionSignature = void(_1, PointCount, _3);
62
  using InputDomain = _2;
63

64 65
  template <typename PointValueVecType, typename OutType>
  VTKM_EXEC void operator()(const PointValueVecType& pointValues,
66 67
                            const vtkm::IdComponent& numPoints,
                            OutType& average) const
68 69 70
  {
    OutType sum = static_cast<OutType>(pointValues[0]);
    for (vtkm::IdComponent pointIndex = 1; pointIndex < numPoints; ++pointIndex)
71
    {
72
      sum = sum + static_cast<OutType>(pointValues[pointIndex]);
73
    }
74 75 76 77 78 79 80 81

    average = sum / static_cast<OutType>(numPoints);
  }
};

class AverageCellToPoint : public vtkm::worklet::WorkletMapCellToPoint
{
public:
82 83
  using ControlSignature = void(FieldInCell<> inCells, CellSetIn topology, FieldOut<> outPoints);
  using ExecutionSignature = void(_1, _3, CellCount);
84
  using InputDomain = _2;
85

86
  template <typename CellVecType, typename OutType>
87 88
  VTKM_EXEC void operator()(const CellVecType& cellValues,
                            OutType& avgVal,
89
                            const vtkm::IdComponent& numCellIDs) const
90 91
  {
    //simple functor that returns the average cell Value.
Sujin Philip's avatar
Sujin Philip committed
92 93
    avgVal = vtkm::TypeTraits<OutType>::ZeroInitialization();
    if (numCellIDs != 0)
94
    {
Sujin Philip's avatar
Sujin Philip committed
95 96 97 98 99
      for (vtkm::IdComponent cellIndex = 0; cellIndex < numCellIDs; ++cellIndex)
      {
        avgVal += static_cast<OutType>(cellValues[cellIndex]);
      }
      avgVal = avgVal / static_cast<OutType>(numCellIDs);
100 101 102 103 104
    }
  }
};

// -----------------------------------------------------------------------------
105
template <typename T>
106 107 108
class Classification : public vtkm::worklet::WorkletMapPointToCell
{
public:
109
  using ControlSignature = void(FieldInPoint<> inNodes,
110
                                CellSetIn cellset,
111
                                FieldOutCell<IdComponentType> outCaseId);
112
  using ExecutionSignature = void(_1, _3);
113
  using InputDomain = _2;
114 115 116

  T IsoValue;

117
  VTKM_CONT
118 119
  Classification(T isovalue)
    : IsoValue(isovalue)
120 121 122
  {
  }

123 124
  template <typename FieldInType>
  VTKM_EXEC void operator()(const FieldInType& fieldIn, vtkm::IdComponent& caseNumber) const
125
  {
126
    using FieldType = typename vtkm::VecTraits<FieldInType>::ComponentType;
127 128
    const FieldType iso = static_cast<FieldType>(this->IsoValue);

129 130 131
    caseNumber = ((fieldIn[0] > iso) | (fieldIn[1] > iso) << 1 | (fieldIn[2] > iso) << 2 |
                  (fieldIn[3] > iso) << 3 | (fieldIn[4] > iso) << 4 | (fieldIn[5] > iso) << 5 |
                  (fieldIn[6] > iso) << 6 | (fieldIn[7] > iso) << 7);
132 133 134
  }
};

135 136 137 138
struct ValueTypes
  : vtkm::ListTagBase<vtkm::UInt32, vtkm::Int32, vtkm::Int64, vtkm::Float32, vtkm::Float64>
{
};
139
using StorageListTag = ::vtkm::cont::StorageListTagBasic;
140 141 142

/// This class runs a series of micro-benchmarks to measure
/// performance of different field operations
143 144 145
template <class DeviceAdapterTag>
class BenchmarkTopologyAlgorithms
{
146
  using StorageTag = vtkm::cont::StorageTagBasic;
147

148
  using Algorithm = vtkm::cont::DeviceAdapterAlgorithm<DeviceAdapterTag>;
149

150
  using Timer = vtkm::cont::Timer<DeviceAdapterTag>;
151

152
  using ValueDynamicHandle = vtkm::cont::DynamicArrayHandleBase<ValueTypes, StorageListTag>;
153

154
private:
155 156 157 158
  template <typename T, typename Enable = void>
  struct NumberGenerator
  {
  };
159

160 161
  template <typename T>
  struct NumberGenerator<T, typename std::enable_if<std::is_floating_point<T>::value>::type>
162 163 164
  {
    std::mt19937 rng;
    std::uniform_real_distribution<T> distribution;
165 166 167 168 169
    NumberGenerator(T low, T high)
      : rng()
      , distribution(low, high)
    {
    }
170 171 172
    T next() { return distribution(rng); }
  };

173 174
  template <typename T>
  struct NumberGenerator<T, typename std::enable_if<!std::is_floating_point<T>::value>::type>
175 176 177 178
  {
    std::mt19937 rng;
    std::uniform_int_distribution<T> distribution;

179 180 181 182 183
    NumberGenerator(T low, T high)
      : rng()
      , distribution(low, high)
    {
    }
184 185 186
    T next() { return distribution(rng); }
  };

187 188 189 190 191
  template <typename Value>
  struct BenchCellToPointAvg
  {
    std::vector<Value> input;
    vtkm::cont::ArrayHandle<Value, StorageTag> InputHandle;
192
    std::size_t DomainSize;
193

194
    VTKM_CONT
195 196
    BenchCellToPointAvg()
    {
197
      NumberGenerator<Value> generator(static_cast<Value>(1.0), static_cast<Value>(100.0));
198
      //cube size is points in each dim
199 200 201
      this->DomainSize = (CUBE_SIZE - 1) * (CUBE_SIZE - 1) * (CUBE_SIZE - 1);
      this->input.resize(DomainSize);
      for (std::size_t i = 0; i < DomainSize; ++i)
202
      {
203
        this->input[i] = generator.next();
204 205 206 207
      }
      this->InputHandle = vtkm::cont::make_ArrayHandle(this->input);
    }

208
    VTKM_CONT
209 210 211
    vtkm::Float64 operator()()
    {
      vtkm::cont::CellSetStructured<3> cellSet;
212 213
      cellSet.SetPointDimensions(vtkm::Id3(CUBE_SIZE, CUBE_SIZE, CUBE_SIZE));
      vtkm::cont::ArrayHandle<Value, StorageTag> result;
214 215 216

      Timer timer;

217 218
      vtkm::worklet::DispatcherMapTopology<AverageCellToPoint> dispatcher;
      dispatcher.Invoke(this->InputHandle, cellSet, result);
219
      //result.SyncControlArray();
220 221 222 223

      return timer.GetElapsedTime();
    }

224 225
    virtual std::string Type() const { return std::string("Static"); }

226
    VTKM_CONT
227 228
    std::string Description() const
    {
229

230
      std::stringstream description;
231 232 233
      description << "Computing Cell To Point Average "
                  << "[" << this->Type() << "] "
                  << "with a domain size of: " << this->DomainSize;
234 235 236 237
      return description.str();
    }
  };

238 239 240
  template <typename Value>
  struct BenchCellToPointAvgDynamic : public BenchCellToPointAvg<Value>
  {
241

242
    VTKM_CONT
243 244 245
    vtkm::Float64 operator()()
    {
      vtkm::cont::CellSetStructured<3> cellSet;
246
      cellSet.SetPointDimensions(vtkm::Id3(CUBE_SIZE, CUBE_SIZE, CUBE_SIZE));
247 248

      ValueDynamicHandle dinput(this->InputHandle);
249
      vtkm::cont::ArrayHandle<Value, StorageTag> result;
250 251 252

      Timer timer;

253
      vtkm::worklet::DispatcherMapTopology<AverageCellToPoint> dispatcher;
254

255
      dispatcher.Invoke(dinput, cellSet, result);
256
      //result.SyncControlArray();
257 258 259 260 261 262 263

      return timer.GetElapsedTime();
    }

    virtual std::string Type() const { return std::string("Dynamic"); }
  };

264
  VTKM_MAKE_BENCHMARK(CellToPointAvg, BenchCellToPointAvg);
265
  VTKM_MAKE_BENCHMARK(CellToPointAvgDynamic, BenchCellToPointAvgDynamic);
266

267 268 269 270 271
  template <typename Value>
  struct BenchPointToCellAvg
  {
    std::vector<Value> input;
    vtkm::cont::ArrayHandle<Value, StorageTag> InputHandle;
272
    std::size_t DomainSize;
273

274
    VTKM_CONT
275 276
    BenchPointToCellAvg()
    {
277
      NumberGenerator<Value> generator(static_cast<Value>(1.0), static_cast<Value>(100.0));
278

279 280 281
      this->DomainSize = (CUBE_SIZE) * (CUBE_SIZE) * (CUBE_SIZE);
      this->input.resize(DomainSize);
      for (std::size_t i = 0; i < DomainSize; ++i)
282
      {
283
        this->input[i] = generator.next();
284 285 286 287
      }
      this->InputHandle = vtkm::cont::make_ArrayHandle(this->input);
    }

288
    VTKM_CONT
289 290 291
    vtkm::Float64 operator()()
    {
      vtkm::cont::CellSetStructured<3> cellSet;
292 293
      cellSet.SetPointDimensions(vtkm::Id3(CUBE_SIZE, CUBE_SIZE, CUBE_SIZE));
      vtkm::cont::ArrayHandle<Value, StorageTag> result;
294 295 296

      Timer timer;

297 298
      vtkm::worklet::DispatcherMapTopology<AveragePointToCell> dispatcher;
      dispatcher.Invoke(this->InputHandle, cellSet, result);
299
      //result.SyncControlArray();
300 301 302 303

      return timer.GetElapsedTime();
    }

304 305
    virtual std::string Type() const { return std::string("Static"); }

306
    VTKM_CONT
307 308
    std::string Description() const
    {
309

310
      std::stringstream description;
311 312 313
      description << "Computing Point To Cell Average "
                  << "[" << this->Type() << "] "
                  << "with a domain size of: " << this->DomainSize;
314 315 316 317
      return description.str();
    }
  };

318 319 320
  template <typename Value>
  struct BenchPointToCellAvgDynamic : public BenchPointToCellAvg<Value>
  {
321

322
    VTKM_CONT
323 324 325
    vtkm::Float64 operator()()
    {
      vtkm::cont::CellSetStructured<3> cellSet;
326
      cellSet.SetPointDimensions(vtkm::Id3(CUBE_SIZE, CUBE_SIZE, CUBE_SIZE));
327 328

      ValueDynamicHandle dinput(this->InputHandle);
329
      vtkm::cont::ArrayHandle<Value, StorageTag> result;
330 331 332

      Timer timer;

333 334
      vtkm::worklet::DispatcherMapTopology<AveragePointToCell> dispatcher;
      dispatcher.Invoke(dinput, cellSet, result);
335
      //result.SyncControlArray();
336 337 338 339 340 341 342

      return timer.GetElapsedTime();
    }

    virtual std::string Type() const { return std::string("Dynamic"); }
  };

343
  VTKM_MAKE_BENCHMARK(PointToCellAvg, BenchPointToCellAvg);
344
  VTKM_MAKE_BENCHMARK(PointToCellAvgDynamic, BenchPointToCellAvgDynamic);
345

346 347 348 349 350
  template <typename Value>
  struct BenchClassification
  {
    std::vector<Value> input;
    vtkm::cont::ArrayHandle<Value, StorageTag> InputHandle;
351
    Value IsoValue;
352
    size_t DomainSize;
353

354
    VTKM_CONT
355 356
    BenchClassification()
    {
357
      NumberGenerator<Value> generator(static_cast<Value>(1.0), static_cast<Value>(100.0));
358

359 360 361
      this->DomainSize = (CUBE_SIZE) * (CUBE_SIZE) * (CUBE_SIZE);
      this->input.resize(DomainSize);
      for (std::size_t i = 0; i < DomainSize; ++i)
362
      {
363
        this->input[i] = generator.next();
364 365
      }
      this->InputHandle = vtkm::cont::make_ArrayHandle(this->input);
366
      this->IsoValue = generator.next();
367 368
    }

369
    VTKM_CONT
370 371 372
    vtkm::Float64 operator()()
    {
      vtkm::cont::CellSetStructured<3> cellSet;
373 374
      cellSet.SetPointDimensions(vtkm::Id3(CUBE_SIZE, CUBE_SIZE, CUBE_SIZE));
      vtkm::cont::ArrayHandle<vtkm::IdComponent, StorageTag> result;
375

376 377
      ValueDynamicHandle dinput(this->InputHandle);

378 379 380
      Timer timer;

      Classification<Value> worklet(this->IsoValue);
381 382
      vtkm::worklet::DispatcherMapTopology<Classification<Value>> dispatcher(worklet);
      dispatcher.Invoke(dinput, cellSet, result);
383
      //result.SyncControlArray();
384 385 386 387

      return timer.GetElapsedTime();
    }

388 389
    virtual std::string Type() const { return std::string("Static"); }

390
    VTKM_CONT
391 392
    std::string Description() const
    {
393

394
      std::stringstream description;
395 396 397
      description << "Computing Marching Cubes Classification "
                  << "[" << this->Type() << "] "
                  << "with a domain size of: " << this->DomainSize;
398 399 400 401
      return description.str();
    }
  };

402 403 404
  template <typename Value>
  struct BenchClassificationDynamic : public BenchClassification<Value>
  {
405
    VTKM_CONT
406 407 408
    vtkm::Float64 operator()()
    {
      vtkm::cont::CellSetStructured<3> cellSet;
409 410
      cellSet.SetPointDimensions(vtkm::Id3(CUBE_SIZE, CUBE_SIZE, CUBE_SIZE));
      vtkm::cont::ArrayHandle<vtkm::IdComponent, StorageTag> result;
411 412 413 414

      Timer timer;

      Classification<Value> worklet(this->IsoValue);
415 416
      vtkm::worklet::DispatcherMapTopology<Classification<Value>> dispatcher(worklet);
      dispatcher.Invoke(this->InputHandle, cellSet, result);
417
      //result.SyncControlArray();
418 419 420 421 422 423 424

      return timer.GetElapsedTime();
    }

    virtual std::string Type() const { return std::string("Dynamic"); }
  };

425
  VTKM_MAKE_BENCHMARK(Classification, BenchClassification);
426
  VTKM_MAKE_BENCHMARK(ClassificationDynamic, BenchClassificationDynamic);
427 428

public:
429 430
  static VTKM_CONT int Run(int benchmarks)
  {
431 432
    std::cout << DIVIDER << "\nRunning Topology Algorithm benchmarks\n";

433 434
    if (benchmarks & CELL_TO_POINT)
    {
435 436
      std::cout << DIVIDER << "\nBenchmarking Cell To Point Average\n";
      VTKM_RUN_BENCHMARK(CellToPointAvg, ValueTypes());
437
      VTKM_RUN_BENCHMARK(CellToPointAvgDynamic, ValueTypes());
438 439
    }

440 441
    if (benchmarks & POINT_TO_CELL)
    {
442 443
      std::cout << DIVIDER << "\nBenchmarking Point to Cell Average\n";
      VTKM_RUN_BENCHMARK(PointToCellAvg, ValueTypes());
444
      VTKM_RUN_BENCHMARK(PointToCellAvgDynamic, ValueTypes());
445 446
    }

447 448
    if (benchmarks & MC_CLASSIFY)
    {
449 450
      std::cout << DIVIDER << "\nBenchmarking Hex/Voxel MC Classification\n";
      VTKM_RUN_BENCHMARK(Classification, ValueTypes());
451
      VTKM_RUN_BENCHMARK(ClassificationDynamic, ValueTypes());
452 453 454 455 456 457 458 459 460 461
    }

    return 0;
  }
};

#undef ARRAY_SIZE
}
} // namespace vtkm::benchmarking

462
int main(int argc, char* argv[])
463 464
{
  int benchmarks = 0;
465 466
  if (argc < 2)
  {
467 468
    benchmarks = vtkm::benchmarking::ALL;
  }
469 470 471 472
  else
  {
    for (int i = 1; i < argc; ++i)
    {
473
      std::string arg = argv[i];
474 475 476
      std::transform(arg.begin(), arg.end(), arg.begin(), [](char c) {
        return static_cast<char>(std::tolower(static_cast<unsigned char>(c)));
      });
477 478 479 480 481 482 483 484 485 486 487 488
      if (arg == "celltopoint")
      {
        benchmarks |= vtkm::benchmarking::CELL_TO_POINT;
      }
      else if (arg == "pointtocell")
      {
        benchmarks |= vtkm::benchmarking::POINT_TO_CELL;
      }
      else if (arg == "classify")
      {
        benchmarks |= vtkm::benchmarking::MC_CLASSIFY;
      }
489 490
      else
      {
491 492 493 494 495 496 497
        std::cout << "Unrecognized benchmark: " << argv[i] << std::endl;
        return 1;
      }
    }
  }

  //now actually execute the benchmarks
498 499
  return vtkm::benchmarking::BenchmarkTopologyAlgorithms<VTKM_DEFAULT_DEVICE_ADAPTER_TAG>::Run(
    benchmarks);
500
}