Merge topic 'osc_refactor'

1d2c4c29 fixed single/double precision converstion
e680815d add include sstream and see if it clears the dashboard
d5973276 remove comment about DeviceAdapter.h
08ad08b1 Merge branch 'master' into osc_refactor
7add915d include DeviceAdapter.h in StorageVirtual.h
b9c2ddf9 add missing #include
5448ffb2 decouple Timer.h from device compiler
e90b186b

 do not use device compiler to compile source unit tests
...
Acked-by: Kitware Robot <kwrobot@kitware.com>
Acked-by: Kenneth Moreland <morelandkd@ornl.gov>
Merge-request: !2602

benchmarking/Benchmarker.h

@@ -10,7 +10,7 @@
 
 #ifndef vtk_m_benchmarking_Benchmarker_h
 #define vtk_m_benchmarking_Benchmarker_h
-
+#include <sstream>
 #include <vtkm/cont/Initialize.h>
 #include <vtkm/cont/Logging.h>
 #include <vtkm/cont/RuntimeDeviceTracker.h>

vtkm/cont/CellSetStructured.h

@@ -14,7 +14,6 @@
 
 #include <vtkm/TopologyElementTag.h>
 #include <vtkm/cont/CellSet.h>
-#include <vtkm/cont/DeviceAdapter.h>
 #include <vtkm/cont/ErrorBadType.h>
 #include <vtkm/exec/ConnectivityStructured.h>
 #include <vtkm/internal/ConnectivityStructuredInternals.h>

vtkm/cont/FieldRangeCompute.h

@@ -14,6 +14,8 @@
 #include <vtkm/cont/Field.h>
 #include <vtkm/cont/PartitionedDataSet.h>
 
+#include <numeric>
+
 namespace vtkm
 {
 namespace cont

vtkm/cont/StorageVirtual.cxx

@@ -10,8 +10,6 @@
 #define vtk_m_cont_StorageVirtual_cxx
 #include <vtkm/cont/StorageVirtual.h>
 
-#include <vtkm/cont/DeviceAdapter.h>
-
 namespace vtkm
 {
 namespace cont

vtkm/cont/StorageVirtual.hxx

@@ -10,6 +10,7 @@
 #ifndef vtk_m_cont_StorageVirtual_hxx
 #define vtk_m_cont_StorageVirtual_hxx
 
+#include <vtkm/cont/DeviceAdapter.h>
 #include <vtkm/cont/StorageVirtual.h>
 #include <vtkm/cont/TryExecute.h>
 

vtkm/cont/Timer.cxx

@@ -7,6 +7,8 @@
 //  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
 //  PURPOSE.  See the above copyright notice for more information.
 //============================================================================
+#include <vtkm/cont/DeviceAdapter.h>
+#include <vtkm/cont/DeviceAdapterList.h>
 #include <vtkm/cont/Logging.h>
 #include <vtkm/cont/RuntimeDeviceTracker.h>
 #include <vtkm/cont/Timer.h>

vtkm/cont/Timer.h

@@ -10,9 +10,7 @@
 #ifndef vtk_m_cont_Timer_h
 #define vtk_m_cont_Timer_h
 
-#include <vtkm/List.h>
-#include <vtkm/cont/DeviceAdapter.h>
-#include <vtkm/cont/DeviceAdapterList.h>
+#include <vtkm/cont/DeviceAdapterTag.h>
 
 #include <vtkm/cont/vtkm_cont_export.h>
 

vtkm/filter/Filter.h

@@ -12,6 +12,7 @@
 
 #include <vtkm/cont/CoordinateSystem.h>
 #include <vtkm/cont/DataSet.h>
+#include <vtkm/cont/ErrorExecution.h>
 #include <vtkm/cont/Field.h>
 #include <vtkm/cont/Invoker.h>
 #include <vtkm/cont/Logging.h>

vtkm/io/ImageUtils.cxx

@@ -7,12 +7,14 @@
 //  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
 //  PURPOSE.  See the above copyright notice for more information.
 //============================================================================
-#include <vtkm/io/ImageUtils.h>
+
+#include <vtkm/cont/ErrorExecution.h>
 
 #include <vtkm/io/FileUtils.h>
 #include <vtkm/io/ImageReaderBase.h>
 #include <vtkm/io/ImageReaderPNG.h>
 #include <vtkm/io/ImageReaderPNM.h>
+#include <vtkm/io/ImageUtils.h>
 #include <vtkm/io/ImageWriterBase.h>
 #include <vtkm/io/ImageWriterPNG.h>
 #include <vtkm/io/ImageWriterPNM.h>

vtkm/source/Oscillator.cxx

@@ -8,56 +8,188 @@
 //  PURPOSE.  See the above copyright notice for more information.
 //============================================================================
 #include <vtkm/source/Oscillator.h>
+#include <vtkm/worklet/WorkletMapField.h>
 
 namespace vtkm
 {
 namespace source
 {
+namespace internal
+{
+struct Oscillator
+{
+  vtkm::Vec3f Center;
+  vtkm::FloatDefault Radius;
+  vtkm::FloatDefault Omega;
+  vtkm::FloatDefault Zeta;
+};
+
+class OscillatorSource : public vtkm::worklet::WorkletMapField
+{
+public:
+  typedef void ControlSignature(FieldIn, FieldOut);
+  typedef _2 ExecutionSignature(_1);
+
+  VTKM_CONT
+  void AddPeriodic(vtkm::FloatDefault x,
+                   vtkm::FloatDefault y,
+                   vtkm::FloatDefault z,
+                   vtkm::FloatDefault radius,
+                   vtkm::FloatDefault omega,
+                   vtkm::FloatDefault zeta)
+  {
+    if (this->PeriodicOscillators.GetNumberOfComponents() < MAX_OSCILLATORS)
+    {
+      this->PeriodicOscillators.Append(Oscillator{ { x, y, z }, radius, omega, zeta });
+    }
+  }
+
+  VTKM_CONT
+  void AddDamped(vtkm::FloatDefault x,
+                 vtkm::FloatDefault y,
+                 vtkm::FloatDefault z,
+                 vtkm::FloatDefault radius,
+                 vtkm::FloatDefault omega,
+                 vtkm::FloatDefault zeta)
+  {
+    if (this->DampedOscillators.GetNumberOfComponents() < MAX_OSCILLATORS)
+    {
+      this->DampedOscillators.Append(Oscillator{ { x, y, z }, radius, omega, zeta });
+    }
+  }
+
+  VTKM_CONT
+  void AddDecaying(vtkm::FloatDefault x,
+                   vtkm::FloatDefault y,
+                   vtkm::FloatDefault z,
+                   vtkm::FloatDefault radius,
+                   vtkm::FloatDefault omega,
+                   vtkm::FloatDefault zeta)
+  {
+    if (this->DecayingOscillators.GetNumberOfComponents() < MAX_OSCILLATORS)
+    {
+      this->DecayingOscillators.Append(Oscillator{ { x, y, z }, radius, omega, zeta });
+    }
+  }
+
+  VTKM_CONT
+  void SetTime(vtkm::FloatDefault time) { this->Time = time; }
+
+  VTKM_EXEC
+  vtkm::FloatDefault operator()(const vtkm::Vec3f& vec) const
+  {
+    vtkm::UInt8 oIdx;
+    vtkm::FloatDefault t0, t, result = 0;
+    const internal::Oscillator* oscillator;
+
+    t0 = 0.0;
+    t = vtkm::FloatDefault(this->Time * 2 * 3.14159265358979323846);
+
+    // Compute damped
+    for (oIdx = 0; oIdx < this->DampedOscillators.GetNumberOfComponents(); oIdx++)
+    {
+      oscillator = &this->DampedOscillators[oIdx];
+
+      vtkm::Vec3f delta = oscillator->Center - vec;
+      vtkm::FloatDefault dist2 = dot(delta, delta);
+      vtkm::FloatDefault dist_damp =
+        vtkm::Exp(-dist2 / (2 * oscillator->Radius * oscillator->Radius));
+      vtkm::FloatDefault phi = vtkm::ACos(oscillator->Zeta);
+      vtkm::FloatDefault val = vtkm::FloatDefault(
+        1. -
+        vtkm::Exp(-oscillator->Zeta * oscillator->Omega * t0) *
+          (vtkm::Sin(vtkm::Sqrt(1 - oscillator->Zeta * oscillator->Zeta) * oscillator->Omega * t +
+                     phi) /
+           vtkm::Sin(phi)));
+      result += val * dist_damp;
+    }
+
+    // Compute decaying
+    for (oIdx = 0; oIdx < this->DecayingOscillators.GetNumberOfComponents(); oIdx++)
+    {
+      oscillator = &this->DecayingOscillators[oIdx];
+      t = t0 + 1 / oscillator->Omega;
+      vtkm::Vec3f delta = oscillator->Center - vec;
+      vtkm::FloatDefault dist2 = dot(delta, delta);
+      vtkm::FloatDefault dist_damp =
+        vtkm::Exp(-dist2 / (2 * oscillator->Radius * oscillator->Radius));
+      vtkm::FloatDefault val = vtkm::Sin(t / oscillator->Omega) / (oscillator->Omega * t);
+      result += val * dist_damp;
+    }
+
+    // Compute periodic
+    for (oIdx = 0; oIdx < this->PeriodicOscillators.GetNumberOfComponents(); oIdx++)
+    {
+      oscillator = &this->PeriodicOscillators[oIdx];
+      t = t0 + 1 / oscillator->Omega;
+      vtkm::Vec3f delta = oscillator->Center - vec;
+      vtkm::FloatDefault dist2 = dot(delta, delta);
+      vtkm::FloatDefault dist_damp =
+        vtkm::Exp(-dist2 / (2 * oscillator->Radius * oscillator->Radius));
+      vtkm::FloatDefault val = vtkm::Sin(t / oscillator->Omega);
+      result += val * dist_damp;
+    }
+
+    // We are done...
+    return result;
+  }
+
+private:
+  static constexpr vtkm::IdComponent MAX_OSCILLATORS = 10;
+  vtkm::VecVariable<internal::Oscillator, MAX_OSCILLATORS> PeriodicOscillators;
+  vtkm::VecVariable<internal::Oscillator, MAX_OSCILLATORS> DampedOscillators;
+  vtkm::VecVariable<internal::Oscillator, MAX_OSCILLATORS> DecayingOscillators;
+  vtkm::FloatDefault Time;
+}; // OscillatorSource
+
+} // internal
 
 //-----------------------------------------------------------------------------
 Oscillator::Oscillator(vtkm::Id3 dims)
   : Dims(dims)
-  , Worklet()
+  , Worklet(std::make_unique<internal::OscillatorSource>())
 {
 }
 
+Oscillator::~Oscillator() = default;
+
 //-----------------------------------------------------------------------------
-void Oscillator::SetTime(vtkm::Float64 time)
+void Oscillator::SetTime(vtkm::FloatDefault time)
 {
-  this->Worklet.SetTime(time);
+  this->Worklet->SetTime(time);
 }
 
 //-----------------------------------------------------------------------------
-void Oscillator::AddPeriodic(vtkm::Float64 x,
-                             vtkm::Float64 y,
-                             vtkm::Float64 z,
-                             vtkm::Float64 radius,
-                             vtkm::Float64 omega,
-                             vtkm::Float64 zeta)
+void Oscillator::AddPeriodic(vtkm::FloatDefault x,
+                             vtkm::FloatDefault y,
+                             vtkm::FloatDefault z,
+                             vtkm::FloatDefault radius,
+                             vtkm::FloatDefault omega,
+                             vtkm::FloatDefault zeta)
 {
-  this->Worklet.AddPeriodic(x, y, z, radius, omega, zeta);
+  this->Worklet->AddPeriodic(x, y, z, radius, omega, zeta);
 }
 
 //-----------------------------------------------------------------------------
-void Oscillator::AddDamped(vtkm::Float64 x,
-                           vtkm::Float64 y,
-                           vtkm::Float64 z,
-                           vtkm::Float64 radius,
-                           vtkm::Float64 omega,
-                           vtkm::Float64 zeta)
+void Oscillator::AddDamped(vtkm::FloatDefault x,
+                           vtkm::FloatDefault y,
+                           vtkm::FloatDefault z,
+                           vtkm::FloatDefault radius,
+                           vtkm::FloatDefault omega,
+                           vtkm::FloatDefault zeta)
 {
-  this->Worklet.AddDamped(x, y, z, radius, omega, zeta);
+  this->Worklet->AddDamped(x, y, z, radius, omega, zeta);
 }
 
 //-----------------------------------------------------------------------------
-void Oscillator::AddDecaying(vtkm::Float64 x,
-                             vtkm::Float64 y,
-                             vtkm::Float64 z,
-                             vtkm::Float64 radius,
-                             vtkm::Float64 omega,
-                             vtkm::Float64 zeta)
+void Oscillator::AddDecaying(vtkm::FloatDefault x,
+                             vtkm::FloatDefault y,
+                             vtkm::FloatDefault z,
+                             vtkm::FloatDefault radius,
+                             vtkm::FloatDefault omega,
+                             vtkm::FloatDefault zeta)
 {
-  this->Worklet.AddDecaying(x, y, z, radius, omega, zeta);
+  this->Worklet->AddDecaying(x, y, z, radius, omega, zeta);
 }
 
 
@@ -68,8 +200,9 @@ vtkm::cont::DataSet Oscillator::Execute() const
 
   vtkm::cont::DataSet dataSet;
 
+  const vtkm::Id3 pdims{ this->Dims + vtkm::Id3{ 1, 1, 1 } };
   vtkm::cont::CellSetStructured<3> cellSet;
-  cellSet.SetPointDimensions(this->Dims);
+  cellSet.SetPointDimensions(pdims);
   dataSet.SetCellSet(cellSet);
 
   const vtkm::Vec3f origin(0.0f, 0.0f, 0.0f);
@@ -77,15 +210,12 @@ vtkm::cont::DataSet Oscillator::Execute() const
                             1.0f / static_cast<vtkm::FloatDefault>(this->Dims[1]),
                             1.0f / static_cast<vtkm::FloatDefault>(this->Dims[2]));
 
-  const vtkm::Id3 pdims{ this->Dims + vtkm::Id3{ 1, 1, 1 } };
   vtkm::cont::ArrayHandleUniformPointCoordinates coordinates(pdims, origin, spacing);
   dataSet.AddCoordinateSystem(vtkm::cont::CoordinateSystem("coordinates", coordinates));
 
 
-  vtkm::cont::ArrayHandle<vtkm::Float64> outArray;
-  //todo, we need to use the policy to determine the valid conversions
-  //that the dispatcher should do
-  this->Invoke(this->Worklet, coordinates, outArray);
+  vtkm::cont::ArrayHandle<vtkm::FloatDefault> outArray;
+  this->Invoke(*(this->Worklet), coordinates, outArray);
   dataSet.AddField(vtkm::cont::make_FieldPoint("oscillating", outArray));
 
   return dataSet;

vtkm/source/Oscillator.h

@@ -11,12 +11,15 @@
 #define vtk_m_source_OscillatorSource_h
 
 #include <vtkm/source/Source.h>
-#include <vtkm/worklet/OscillatorSource.h>
 
 namespace vtkm
 {
 namespace source
 {
+namespace internal
+{
+class OscillatorSource;
+}
 
 /**\brief An analytical, time-varying uniform dataset with a point based array
  *
@@ -33,38 +36,44 @@ public:
   VTKM_CONT
   Oscillator(vtkm::Id3 dims);
 
+  // We can not declare default destructor here since compiler does not know how
+  // to create one for the Worklet at this point yet. However, the implementation
+  // in Oscillator.cxx does have ~Oscillator() = default;
+  VTKM_CONT
+  ~Oscillator();
+
   VTKM_CONT
-  void SetTime(vtkm::Float64 time);
+  void SetTime(vtkm::FloatDefault time);
 
   VTKM_CONT
-  void AddPeriodic(vtkm::Float64 x,
-                   vtkm::Float64 y,
-                   vtkm::Float64 z,
-                   vtkm::Float64 radius,
-                   vtkm::Float64 omega,
-                   vtkm::Float64 zeta);
+  void AddPeriodic(vtkm::FloatDefault x,
+                   vtkm::FloatDefault y,
+                   vtkm::FloatDefault z,
+                   vtkm::FloatDefault radius,
+                   vtkm::FloatDefault omega,
+                   vtkm::FloatDefault zeta);
 
   VTKM_CONT
-  void AddDamped(vtkm::Float64 x,
-                 vtkm::Float64 y,
-                 vtkm::Float64 z,
-                 vtkm::Float64 radius,
-                 vtkm::Float64 omega,
-                 vtkm::Float64 zeta);
+  void AddDamped(vtkm::FloatDefault x,
+                 vtkm::FloatDefault y,
+                 vtkm::FloatDefault z,
+                 vtkm::FloatDefault radius,
+                 vtkm::FloatDefault omega,
+                 vtkm::FloatDefault zeta);
 
   VTKM_CONT
-  void AddDecaying(vtkm::Float64 x,
-                   vtkm::Float64 y,
-                   vtkm::Float64 z,
-                   vtkm::Float64 radius,
-                   vtkm::Float64 omega,
-                   vtkm::Float64 zeta);
+  void AddDecaying(vtkm::FloatDefault x,
+                   vtkm::FloatDefault y,
+                   vtkm::FloatDefault z,
+                   vtkm::FloatDefault radius,
+                   vtkm::FloatDefault omega,
+                   vtkm::FloatDefault zeta);
 
   VTKM_CONT vtkm::cont::DataSet Execute() const;
 
 private:
   vtkm::Id3 Dims;
-  vtkm::worklet::OscillatorSource Worklet;
+  std::unique_ptr<internal::OscillatorSource> Worklet;
 };
 }
 }

vtkm/source/testing/CMakeLists.txt

@@ -9,6 +9,7 @@
 ##============================================================================
 
 set(unit_tests
+  UnitTestOscillatorSource.cxx
   UnitTestTangleSource.cxx
   UnitTestWaveletSource.cxx
  )
@@ -16,5 +17,4 @@ set(unit_tests
 vtkm_unit_tests(
   SOURCES ${unit_tests}
   LIBRARIES vtkm_source
-  ALL_BACKENDS
 )

vtkm/source/testing/UnitTestOscillatorSource.cxx

+//============================================================================
+//  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.
+//============================================================================
+
+#include <vtkm/source/Oscillator.h>
+
+#include <vtkm/cont/Timer.h>
+#include <vtkm/cont/testing/Testing.h>
+
+void OscillatorSourceTest()
+{
+  vtkm::cont::Timer timer;
+  timer.Start();
+
+  vtkm::source::Oscillator source(vtkm::Id3{ 20, 20, 20 });
+  source.SetTime(0.5);
+  source.AddDamped(0.25f, 0.25f, 0.25f, 0.5f, 0.1f, 0.2f);
+  source.AddDecaying(0.5f, 0.5f, 0.5f, 0.35f, 0.2f, 0.1f);
+  source.AddPeriodic(0.6f, 0.2f, 0.7f, 0.15f, 0.1f, 0.2f);
+
+  vtkm::cont::DataSet ds = source.Execute();
+
+  double time = timer.GetElapsedTime();
+
+  std::cout << "Default oscillator took " << time << "s.\n";
+
+  {
+    auto coords = ds.GetCoordinateSystem("coordinates");
+    auto data = coords.GetData();
+    VTKM_TEST_ASSERT(test_equal(data.GetNumberOfValues(), 9261), "Incorrect number of points.");
+  }
+
+  {
+    auto cells = ds.GetCellSet();
+    VTKM_TEST_ASSERT(test_equal(cells.GetNumberOfCells(), 8000), "Incorrect number of cells.");
+  }
+
+  // Spot check some node scalars
+  {
+    using ScalarHandleType = vtkm::cont::ArrayHandle<vtkm::FloatDefault>;
+
+    auto field = ds.GetPointField("oscillating");
+    auto dynData = field.GetData();
+    VTKM_TEST_ASSERT(dynData.IsType<ScalarHandleType>(), "Invalid scalar handle type.");
+    ScalarHandleType handle = dynData.AsArrayHandle<ScalarHandleType>();
+    auto data = handle.ReadPortal();
+
+    VTKM_TEST_ASSERT(test_equal(data.GetNumberOfValues(), 9261), "Incorrect number of scalars.");
+
+    VTKM_TEST_ASSERT(test_equal(data.Get(0), -0.0163996), "Incorrect scalar value.");
+    VTKM_TEST_ASSERT(test_equal(data.Get(16), -0.0182232), "Incorrect scalar value.");
+    VTKM_TEST_ASSERT(test_equal(data.Get(21), -0.0181952), "Incorrect scalar value.");
+    VTKM_TEST_ASSERT(test_equal(data.Get(3110), -0.0404135), "Incorrect scalar value.");
+  }
+}
+
+int UnitTestOscillatorSource(int argc, char* argv[])
+{
+  return vtkm::cont::testing::Testing::Run(OscillatorSourceTest, argc, argv);
+}

vtkm/source/testing/UnitTestWaveletSource.cxx

@@ -10,7 +10,6 @@
 
 #include <vtkm/source/Wavelet.h>
 
-#include <vtkm/cont/Timer.h>
 #include <vtkm/cont/testing/Testing.h>
 
 void WaveletSourceTest()

vtkm/worklet/CMakeLists.txt

@@ -57,7 +57,6 @@ set(headers
   OrientNormals.h
   OrientPointNormals.h
   OrientPointAndCellNormals.h
-  OscillatorSource.h
   ParticleAdvection.h
   PointAverage.h
   PointElevation.h

vtkm/worklet/OscillatorSource.h

-//============================================================================
-//  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.
-//============================================================================
-#ifndef vtk_m_worklet_OscillatorSource_h
-#define vtk_m_worklet_OscillatorSource_h
-
-#include <vtkm/Math.h>
-#include <vtkm/worklet/WorkletMapField.h>
-
-#define MAX_OSCILLATORS 10
-
-namespace vtkm
-{
-namespace worklet
-{
-namespace internal
-{
-
-struct Oscillator
-{
-  void Set(vtkm::Float64 x,
-           vtkm::Float64 y,
-           vtkm::Float64 z,
-           vtkm::Float64 radius,
-           vtkm::Float64 omega,
-           vtkm::Float64 zeta)
-  {
-    this->Center[0] = x;
-    this->Center[1] = y;
-    this->Center[2] = z;
-    this->Radius = radius;
-    this->Omega = omega;
-    this->Zeta = zeta;
-  }
-
-  vtkm::Vec3f_64 Center;
-  vtkm::Float64 Radius;
-  vtkm::Float64 Omega;
-  vtkm::Float64 Zeta;
-};
-}
-
-class OscillatorSource : public vtkm::worklet::WorkletMapField
-{
-public:
-  typedef void ControlSignature(FieldIn, FieldOut);
-  typedef _2 ExecutionSignature(_1);
-
-  VTKM_CONT
-  OscillatorSource()
-    : NumberOfPeriodics(0)
-    , NumberOfDamped(0)
-    , NumberOfDecaying(0)
-  {
-  }
-
-  VTKM_CONT
-  void AddPeriodic(vtkm::Float64 x,
-                   vtkm::Float64 y,
-                   vtkm::Float64 z,
-                   vtkm::Float64 radius,
-                   vtkm::Float64 omega,
-                   vtkm::Float64 zeta)
-  {
-    if (this->NumberOfPeriodics < MAX_OSCILLATORS)
-    {
-      this->PeriodicOscillators[this->NumberOfPeriodics].Set(x, y, z, radius, omega, zeta);
-      this->NumberOfPeriodics++;
-    }
-  }
-
-  VTKM_CONT
-  void AddDamped(vtkm::Float64 x,
-                 vtkm::Float64 y,
-                 vtkm::Float64 z,
-                 vtkm::Float64 radius,
-                 vtkm::Float64 omega,
-                 vtkm::Float64 zeta)