Merge topic 'temporal_particle_advection'

f878ba8da Move the check for 0 inputs. With mpi, 0 input is ok for a rank.
4b2dbfbad Remove unused header.
b08082472 Clean up the code a bit.
c4341adfe update cmakelists.txt. The merge somehow missed this..
90bed8d0c refactor particle adv-based filters and make a temporal particle adv filter.

Acked-by: Kitware Robot <kwrobot@kitware.com>
Merge-request: !2515
This commit is contained in:
Dave Pugmire 2021-06-24 11:26:16 +00:00 committed by Kitware Robot
commit f7826d1efd
19 changed files with 721 additions and 426 deletions

@ -21,6 +21,8 @@ set(common_headers
FilterDataSetWithField.h FilterDataSetWithField.h
FilterField.h FilterField.h
Filter.h Filter.h
FilterParticleAdvection.h
FilterTemporalParticleAdvection.h
FilterTraits.h FilterTraits.h
MapFieldMergeAverage.h MapFieldMergeAverage.h
MapFieldPermutation.h MapFieldPermutation.h
@ -41,6 +43,8 @@ set(common_header_template_sources
FilterDataSetWithField.hxx FilterDataSetWithField.hxx
FilterField.hxx FilterField.hxx
Filter.hxx Filter.hxx
FilterParticleAdvection.hxx
FilterTemporalParticleAdvection.hxx
PointAverage.hxx PointAverage.hxx
Threshold.hxx Threshold.hxx
ThresholdPoints.hxx ThresholdPoints.hxx
@ -88,6 +92,7 @@ set(extra_headers
ParticleDensityNearestGridPoint.h ParticleDensityNearestGridPoint.h
ParticleAdvection.h ParticleAdvection.h
Pathline.h Pathline.h
PathParticle.h
PointAverage.h PointAverage.h
PointElevation.h PointElevation.h
PointTransform.h PointTransform.h
@ -141,6 +146,7 @@ set(extra_header_template_sources
ParticleDensityNearestGridPoint.hxx ParticleDensityNearestGridPoint.hxx
ParticleAdvection.hxx ParticleAdvection.hxx
Pathline.hxx Pathline.hxx
PathParticle.hxx
PointElevation.hxx PointElevation.hxx
PointTransform.hxx PointTransform.hxx
Probe.hxx Probe.hxx

@ -0,0 +1,93 @@
//============================================================================
// 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_filter_FilterParticleAdvection_h
#define vtk_m_filter_FilterParticleAdvection_h
#include <vtkm/Particle.h>
#include <vtkm/filter/FilterDataSetWithField.h>
#include <vtkm/filter/particleadvection/BoundsMap.h>
#include <vtkm/filter/particleadvection/DataSetIntegrator.h>
namespace vtkm
{
namespace filter
{
/// \brief base class for advecting particles in a vector field.
/// Takes as input a vector field and seed locations and advects the seeds
/// through the flow field.
template <class Derived>
class FilterParticleAdvection : public vtkm::filter::FilterDataSetWithField<Derived>
{
public:
using SupportedTypes = vtkm::TypeListFieldVec3;
VTKM_CONT
FilterParticleAdvection();
VTKM_CONT
void SetStepSize(vtkm::FloatDefault s) { this->StepSize = s; }
VTKM_CONT
void SetNumberOfSteps(vtkm::Id n) { this->NumberOfSteps = n; }
VTKM_CONT
void SetSeeds(const std::vector<vtkm::Particle>& seeds,
vtkm::CopyFlag copyFlag = vtkm::CopyFlag::On)
{
this->Seeds = vtkm::cont::make_ArrayHandle(seeds, copyFlag);
}
VTKM_CONT
void SetSeeds(vtkm::cont::ArrayHandle<vtkm::Particle>& seeds) { this->Seeds = seeds; }
VTKM_CONT
bool GetUseThreadedAlgorithm() { return this->UseThreadedAlgorithm; }
VTKM_CONT
void SetUseThreadedAlgorithm(bool val) { this->UseThreadedAlgorithm = val; }
template <typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet PrepareForExecution(const vtkm::cont::DataSet& input,
vtkm::filter::PolicyBase<DerivedPolicy> policy);
template <typename DerivedPolicy>
VTKM_CONT bool MapFieldOntoOutput(vtkm::cont::DataSet&,
const vtkm::cont::Field&,
vtkm::filter::PolicyBase<DerivedPolicy>)
{
return false;
}
protected:
VTKM_CONT virtual void ValidateOptions() const;
using DSIType = vtkm::filter::particleadvection::DataSetIntegrator;
VTKM_CONT std::vector<DSIType> CreateDataSetIntegrators(
const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::particleadvection::BoundsMap& boundsMap) const;
vtkm::Id NumberOfSteps;
vtkm::FloatDefault StepSize;
vtkm::cont::ArrayHandle<vtkm::Particle> Seeds;
bool UseThreadedAlgorithm;
private:
};
}
} // namespace vtkm::filter
#ifndef vtk_m_filter_FilterParticleAdvection_hxx
#include <vtkm/filter/FilterParticleAdvection.hxx>
#endif
#endif // vtk_m_filter_FilterParticleAdvection_h

@ -0,0 +1,80 @@
//============================================================================
// 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_filter_FilterParticleAdvection_hxx
#define vtk_m_filter_FilterParticleAdvection_hxx
#include <vtkm/filter/FilterParticleAdvection.h>
namespace vtkm
{
namespace filter
{
//-----------------------------------------------------------------------------
template <typename Derived>
inline VTKM_CONT FilterParticleAdvection<Derived>::FilterParticleAdvection()
: vtkm::filter::FilterDataSetWithField<Derived>()
, NumberOfSteps(0)
, StepSize(0)
, UseThreadedAlgorithm(false)
{
}
template <typename Derived>
void FilterParticleAdvection<Derived>::ValidateOptions() const
{
if (this->GetUseCoordinateSystemAsField())
throw vtkm::cont::ErrorFilterExecution("Coordinate system as field not supported");
if (this->Seeds.GetNumberOfValues() == 0)
throw vtkm::cont::ErrorFilterExecution("No seeds provided.");
if (this->NumberOfSteps == 0)
throw vtkm::cont::ErrorFilterExecution("Number of steps not specified.");
if (this->StepSize == 0)
throw vtkm::cont::ErrorFilterExecution("Step size not specified.");
}
template <typename Derived>
std::vector<vtkm::filter::particleadvection::DataSetIntegrator>
FilterParticleAdvection<Derived>::CreateDataSetIntegrators(
const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::particleadvection::BoundsMap& boundsMap) const
{
std::vector<vtkm::filter::particleadvection::DataSetIntegrator> dsi;
if (boundsMap.GetTotalNumBlocks() == 0)
throw vtkm::cont::ErrorFilterExecution("No input datasets.");
std::string activeField = this->GetActiveFieldName();
for (vtkm::Id i = 0; i < input.GetNumberOfPartitions(); i++)
{
vtkm::Id blockId = boundsMap.GetLocalBlockId(i);
auto ds = input.GetPartition(i);
if (!ds.HasPointField(activeField))
throw vtkm::cont::ErrorFilterExecution("Unsupported field assocation");
dsi.push_back(DSIType(ds, blockId, activeField));
}
return dsi;
}
//-----------------------------------------------------------------------------
template <typename Derived>
template <typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet FilterParticleAdvection<Derived>::PrepareForExecution(
const vtkm::cont::DataSet& input,
vtkm::filter::PolicyBase<DerivedPolicy> policy)
{
return (static_cast<Derived*>(this))->DoExecute(input, policy);
}
}
} // namespace vtkm::filter
#endif

@ -0,0 +1,74 @@
//============================================================================
// 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_filter_TemporalFilterParticleAdvection_h
#define vtk_m_filter_TemporalFilterParticleAdvection_h
#include <vtkm/Particle.h>
#include <vtkm/filter/FilterParticleAdvection.h>
#include <vtkm/filter/particleadvection/BoundsMap.h>
#include <vtkm/filter/particleadvection/DataSetIntegrator.h>
namespace vtkm
{
namespace filter
{
/// \brief base class for advecting particles in a vector field.
/// Takes as input a vector field and seed locations and advects the seeds
/// through the flow field.
template <class Derived>
class FilterTemporalParticleAdvection : public vtkm::filter::FilterParticleAdvection<Derived>
{
public:
VTKM_CONT
FilterTemporalParticleAdvection();
VTKM_CONT
void SetPreviousTime(vtkm::FloatDefault t) { this->PreviousTime = t; }
VTKM_CONT
void SetNextTime(vtkm::FloatDefault t) { this->NextTime = t; }
VTKM_CONT
void SetNextDataSet(const vtkm::cont::DataSet& ds)
{
this->NextDataSet = vtkm::cont::PartitionedDataSet(ds);
}
VTKM_CONT
void SetNextDataSet(const vtkm::cont::PartitionedDataSet& pds) { this->NextDataSet = pds; }
template <typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet PrepareForExecution(const vtkm::cont::DataSet& input,
vtkm::filter::PolicyBase<DerivedPolicy> policy);
protected:
VTKM_CONT void ValidateOptions(const vtkm::cont::PartitionedDataSet& input) const;
using DSIType = vtkm::filter::particleadvection::TemporalDataSetIntegrator;
VTKM_CONT std::vector<DSIType> CreateDataSetIntegrators(
const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::particleadvection::BoundsMap& boundsMap) const;
vtkm::FloatDefault PreviousTime;
vtkm::FloatDefault NextTime;
vtkm::cont::PartitionedDataSet NextDataSet;
private:
};
}
} // namespace vtkm::filter
#ifndef vtk_m_filter_FilterTemporalParticleAdvection_hxx
#include <vtkm/filter/FilterTemporalParticleAdvection.hxx>
#endif
#endif // vtk_m_filter_FilterTemporalParticleAdvection_h

@ -0,0 +1,79 @@
//============================================================================
// 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_filter_FilterTemporalParticleAdvection_hxx
#define vtk_m_filter_FilterTemporalParticleAdvection_hxx
#include <vtkm/filter/FilterTemporalParticleAdvection.h>
namespace vtkm
{
namespace filter
{
//-----------------------------------------------------------------------------
template <typename Derived>
inline VTKM_CONT FilterTemporalParticleAdvection<Derived>::FilterTemporalParticleAdvection()
: vtkm::filter::FilterParticleAdvection<Derived>()
, PreviousTime(0)
, NextTime(0)
{
}
template <typename Derived>
void FilterTemporalParticleAdvection<Derived>::ValidateOptions(
const vtkm::cont::PartitionedDataSet& input) const
{
this->vtkm::filter::FilterParticleAdvection<Derived>::ValidateOptions();
if (this->NextDataSet.GetNumberOfPartitions() != input.GetNumberOfPartitions())
throw vtkm::cont::ErrorFilterExecution("Number of partitions do not match");
if (!(this->PreviousTime < this->NextTime))
throw vtkm::cont::ErrorFilterExecution("Previous time must be less than Next time.");
}
template <typename Derived>
std::vector<vtkm::filter::particleadvection::TemporalDataSetIntegrator>
FilterTemporalParticleAdvection<Derived>::CreateDataSetIntegrators(
const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::particleadvection::BoundsMap& boundsMap) const
{
std::vector<DSIType> dsi;
std::string activeField = this->GetActiveFieldName();
if (boundsMap.GetTotalNumBlocks() == 0)
throw vtkm::cont::ErrorFilterExecution("No input datasets.");
for (vtkm::Id i = 0; i < input.GetNumberOfPartitions(); i++)
{
vtkm::Id blockId = boundsMap.GetLocalBlockId(i);
auto dsPrev = input.GetPartition(i);
auto dsNext = this->NextDataSet.GetPartition(i);
if (!dsPrev.HasPointField(activeField) || !dsNext.HasPointField(activeField))
throw vtkm::cont::ErrorFilterExecution("Unsupported field assocation");
dsi.push_back(
DSIType(dsPrev, this->PreviousTime, dsNext, this->NextTime, blockId, activeField));
}
return dsi;
}
//-----------------------------------------------------------------------------
template <typename Derived>
template <typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::DataSet FilterTemporalParticleAdvection<Derived>::PrepareForExecution(
const vtkm::cont::DataSet& input,
vtkm::filter::PolicyBase<DerivedPolicy> policy)
{
return (static_cast<Derived*>(this))->DoExecute(input, policy);
}
}
} // namespace vtkm::filter
#endif

@ -12,54 +12,26 @@
#define vtk_m_filter_ParticleAdvection_h #define vtk_m_filter_ParticleAdvection_h
#include <vtkm/Particle.h> #include <vtkm/Particle.h>
#include <vtkm/filter/FilterDataSetWithField.h> #include <vtkm/filter/FilterParticleAdvection.h>
namespace vtkm namespace vtkm
{ {
namespace filter namespace filter
{ {
/// \brief advect particles in a vector field. /// \brief Advect particles in a vector field.
/// Takes as input a vector field and seed locations and generates the /// Takes as input a vector field and seed locations and generates the
/// end points for each seed through the vector field. /// end points for each seed through the vector field.
class ParticleAdvection : public vtkm::filter::FilterDataSetWithField<ParticleAdvection>
class ParticleAdvection : public vtkm::filter::FilterParticleAdvection<ParticleAdvection>
{ {
public: public:
using SupportedTypes = vtkm::TypeListFieldVec3; VTKM_CONT ParticleAdvection();
VTKM_CONT
ParticleAdvection();
VTKM_CONT
void SetStepSize(vtkm::FloatDefault s) { this->StepSize = s; }
VTKM_CONT
void SetNumberOfSteps(vtkm::Id n) { this->NumberOfSteps = n; }
VTKM_CONT
void SetSeeds(vtkm::cont::ArrayHandle<vtkm::Particle>& seeds);
VTKM_CONT
bool GetUseThreadedAlgorithm() { return this->UseThreadedAlgorithm; }
VTKM_CONT
void SetUseThreadedAlgorithm(bool val) { this->UseThreadedAlgorithm = val; }
template <typename DerivedPolicy> template <typename DerivedPolicy>
vtkm::cont::PartitionedDataSet PrepareForExecution( vtkm::cont::PartitionedDataSet PrepareForExecution(
const vtkm::cont::PartitionedDataSet& input, const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy); const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename DerivedPolicy>
VTKM_CONT bool MapFieldOntoOutput(vtkm::cont::DataSet& result,
const vtkm::cont::Field& field,
vtkm::filter::PolicyBase<DerivedPolicy> policy);
private:
vtkm::Id NumberOfSteps;
vtkm::FloatDefault StepSize;
vtkm::cont::ArrayHandle<vtkm::Particle> Seeds;
bool UseThreadedAlgorithm;
}; };
} }
} // namespace vtkm::filter } // namespace vtkm::filter

@ -10,13 +10,8 @@
#ifndef vtk_m_filter_ParticleAdvection_hxx #ifndef vtk_m_filter_ParticleAdvection_hxx
#define vtk_m_filter_ParticleAdvection_hxx #define vtk_m_filter_ParticleAdvection_hxx
#include <vtkm/filter/ParticleAdvection.h>
#include <vtkm/cont/ArrayCopy.h>
#include <vtkm/cont/ArrayHandleIndex.h>
#include <vtkm/cont/CellSetSingleType.h>
#include <vtkm/cont/ErrorFilterExecution.h> #include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/cont/ParticleArrayCopy.h> #include <vtkm/filter/ParticleAdvection.h>
#include <vtkm/filter/particleadvection/BoundsMap.h> #include <vtkm/filter/particleadvection/BoundsMap.h>
#include <vtkm/filter/particleadvection/DataSetIntegrator.h> #include <vtkm/filter/particleadvection/DataSetIntegrator.h>
#include <vtkm/filter/particleadvection/ParticleAdvectionAlgorithm.h> #include <vtkm/filter/particleadvection/ParticleAdvectionAlgorithm.h>
@ -28,46 +23,23 @@ namespace filter
//----------------------------------------------------------------------------- //-----------------------------------------------------------------------------
inline VTKM_CONT ParticleAdvection::ParticleAdvection() inline VTKM_CONT ParticleAdvection::ParticleAdvection()
: vtkm::filter::FilterDataSetWithField<ParticleAdvection>() : vtkm::filter::FilterParticleAdvection<ParticleAdvection>()
, UseThreadedAlgorithm(false)
{ {
} }
//-----------------------------------------------------------------------------
inline VTKM_CONT void ParticleAdvection::SetSeeds(vtkm::cont::ArrayHandle<vtkm::Particle>& seeds)
{
this->Seeds = seeds;
}
//----------------------------------------------------------------------------- //-----------------------------------------------------------------------------
template <typename DerivedPolicy> template <typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::PartitionedDataSet ParticleAdvection::PrepareForExecution( inline VTKM_CONT vtkm::cont::PartitionedDataSet ParticleAdvection::PrepareForExecution(
const vtkm::cont::PartitionedDataSet& input, const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>&) const vtkm::filter::PolicyBase<DerivedPolicy>&)
{ {
if (this->GetUseCoordinateSystemAsField())
throw vtkm::cont::ErrorFilterExecution("Coordinate system as field not supported");
if (this->Seeds.GetNumberOfValues() == 0)
throw vtkm::cont::ErrorFilterExecution("No seeds provided.");
std::string activeField = this->GetActiveFieldName();
vtkm::filter::particleadvection::BoundsMap boundsMap(input);
using DSIType = vtkm::filter::particleadvection::DataSetIntegrator;
std::vector<DSIType> dsi;
for (vtkm::Id i = 0; i < input.GetNumberOfPartitions(); i++)
{
vtkm::Id blockId = boundsMap.GetLocalBlockId(i);
auto ds = input.GetPartition(i);
if (!ds.HasPointField(activeField))
throw vtkm::cont::ErrorFilterExecution("Unsupported field assocation");
dsi.push_back(DSIType(ds, blockId, activeField));
}
using AlgorithmType = vtkm::filter::particleadvection::ParticleAdvectionAlgorithm; using AlgorithmType = vtkm::filter::particleadvection::ParticleAdvectionAlgorithm;
using ThreadedAlgorithmType = vtkm::filter::particleadvection::ParticleAdvectionThreadedAlgorithm; using ThreadedAlgorithmType = vtkm::filter::particleadvection::ParticleAdvectionThreadedAlgorithm;
this->ValidateOptions();
vtkm::filter::particleadvection::BoundsMap boundsMap(input);
auto dsi = this->CreateDataSetIntegrators(input, boundsMap);
if (this->GetUseThreadedAlgorithm()) if (this->GetUseThreadedAlgorithm())
return vtkm::filter::particleadvection::RunAlgo<DSIType, ThreadedAlgorithmType>( return vtkm::filter::particleadvection::RunAlgo<DSIType, ThreadedAlgorithmType>(
boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds); boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds);
@ -76,14 +48,6 @@ inline VTKM_CONT vtkm::cont::PartitionedDataSet ParticleAdvection::PrepareForExe
boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds); boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds);
} }
//-----------------------------------------------------------------------------
template <typename DerivedPolicy>
inline VTKM_CONT bool ParticleAdvection::MapFieldOntoOutput(vtkm::cont::DataSet&,
const vtkm::cont::Field&,
vtkm::filter::PolicyBase<DerivedPolicy>)
{
return false;
}
} }
} // namespace vtkm::filter } // namespace vtkm::filter
#endif #endif

@ -0,0 +1,46 @@
//============================================================================
// 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_filter_PathParticle_h
#define vtk_m_filter_PathParticle_h
#include <vtkm/filter/FilterTemporalParticleAdvection.h>
namespace vtkm
{
namespace filter
{
/// \brief Advect particles in a time varying vector field.
/// Takes as input a vector field and seed locations and generates the
/// paths taken by the seeds through the vector field.
class PathParticle : public vtkm::filter::FilterTemporalParticleAdvection<PathParticle>
{
public:
VTKM_CONT
PathParticle();
template <typename DerivedPolicy>
vtkm::cont::PartitionedDataSet PrepareForExecution(
const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
protected:
private:
};
}
} // namespace vtkm::filter
#ifndef vtk_m_filter_PathParticle_hxx
#include <vtkm/filter/PathParticle.hxx>
#endif
#endif // vtk_m_filter_PathParticle_h

@ -0,0 +1,55 @@
//============================================================================
// 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_filter_PathParticle_hxx
#define vtk_m_filter_PathParticle_hxx
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/filter/PathParticle.h>
#include <vtkm/filter/particleadvection/BoundsMap.h>
#include <vtkm/filter/particleadvection/DataSetIntegrator.h>
#include <vtkm/filter/particleadvection/ParticleAdvectionAlgorithm.h>
namespace vtkm
{
namespace filter
{
//-----------------------------------------------------------------------------
inline VTKM_CONT PathParticle::PathParticle()
: vtkm::filter::FilterTemporalParticleAdvection<PathParticle>()
{
}
//-----------------------------------------------------------------------------
template <typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::PartitionedDataSet PathParticle::PrepareForExecution(
const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
using AlgorithmType = vtkm::filter::particleadvection::PathParticleAlgorithm;
using ThreadedAlgorithmType = vtkm::filter::particleadvection::PathParticleThreadedAlgorithm;
this->ValidateOptions(input);
vtkm::filter::particleadvection::BoundsMap boundsMap(input);
auto dsi = this->CreateDataSetIntegrators(input, boundsMap);
if (this->GetUseThreadedAlgorithm())
return vtkm::filter::particleadvection::RunAlgo<DSIType, ThreadedAlgorithmType>(
boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds);
else
return vtkm::filter::particleadvection::RunAlgo<DSIType, AlgorithmType>(
boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds);
}
}
} // namespace vtkm::filter
#endif

@ -11,79 +11,36 @@
#ifndef vtk_m_filter_Pathline_h #ifndef vtk_m_filter_Pathline_h
#define vtk_m_filter_Pathline_h #define vtk_m_filter_Pathline_h
#include <vtkm/filter/FilterDataSetWithField.h> #include <vtkm/filter/FilterTemporalParticleAdvection.h>
#include <vtkm/worklet/ParticleAdvection.h>
#include <vtkm/worklet/particleadvection/GridEvaluators.h>
#include <vtkm/worklet/particleadvection/Stepper.h>
namespace vtkm namespace vtkm
{ {
namespace filter namespace filter
{ {
/// \brief generate streamlines from a vector field. /// \brief generate pathlines from a time sequence of vector fields.
/// Takes as input a vector field and seed locations and generates the /// Takes as input a vector field and seed locations and generates the
/// paths taken by the seeds through the vector field. /// paths taken by the seeds through the vector field.
class Pathline : public vtkm::filter::FilterDataSetWithField<Pathline> class Pathline : public vtkm::filter::FilterTemporalParticleAdvection<Pathline>
{ {
public: public:
using SupportedTypes = vtkm::TypeListFieldVec3;
VTKM_CONT VTKM_CONT
Pathline(); Pathline();
VTKM_CONT
void SetPreviousTime(vtkm::FloatDefault t) { this->PreviousTime = t; }
VTKM_CONT
void SetNextTime(vtkm::FloatDefault t) { this->NextTime = t; }
VTKM_CONT
void SetNextDataSet(const vtkm::cont::DataSet& ds)
{
this->NextDataSet = vtkm::cont::PartitionedDataSet(ds);
}
VTKM_CONT
void SetNextDataSet(const vtkm::cont::PartitionedDataSet& pds) { this->NextDataSet = pds; }
VTKM_CONT
void SetStepSize(vtkm::FloatDefault s) { this->StepSize = s; }
VTKM_CONT
void SetNumberOfSteps(vtkm::Id n) { this->NumberOfSteps = n; }
VTKM_CONT
void SetSeeds(vtkm::cont::ArrayHandle<vtkm::Particle>& seeds);
VTKM_CONT
bool GetUseThreadedAlgorithm() { return this->UseThreadedAlgorithm; }
VTKM_CONT
void SetUseThreadedAlgorithm(bool val) { this->UseThreadedAlgorithm = val; }
template <typename DerivedPolicy> template <typename DerivedPolicy>
vtkm::cont::PartitionedDataSet PrepareForExecution( vtkm::cont::PartitionedDataSet PrepareForExecution(
const vtkm::cont::PartitionedDataSet& input, const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy); const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename DerivedPolicy> protected:
VTKM_CONT bool MapFieldOntoOutput(vtkm::cont::DataSet& result,
const vtkm::cont::Field& field,
vtkm::filter::PolicyBase<DerivedPolicy> policy);
private: private:
vtkm::FloatDefault StepSize;
vtkm::FloatDefault PreviousTime;
vtkm::FloatDefault NextTime;
vtkm::cont::PartitionedDataSet NextDataSet;
vtkm::Id NumberOfSteps;
vtkm::cont::ArrayHandle<vtkm::Particle> Seeds;
bool UseThreadedAlgorithm;
}; };
} }
} // namespace vtkm::filter } // namespace vtkm::filter
#ifndef vtk_m_filter_Pathline_hxx
#include <vtkm/filter/Pathline.hxx> #include <vtkm/filter/Pathline.hxx>
#endif
#endif // vtk_m_filter_Pathline_h #endif // vtk_m_filter_Pathline_h

@ -11,14 +11,11 @@
#ifndef vtk_m_filter_Pathline_hxx #ifndef vtk_m_filter_Pathline_hxx
#define vtk_m_filter_Pathline_hxx #define vtk_m_filter_Pathline_hxx
#include <vtkm/filter/Pathline.h>
#include <vtkm/cont/ArrayCopy.h>
#include <vtkm/cont/ArrayHandleIndex.h>
#include <vtkm/cont/ErrorFilterExecution.h> #include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/filter/Pathline.h>
#include <vtkm/filter/particleadvection/BoundsMap.h> #include <vtkm/filter/particleadvection/BoundsMap.h>
#include <vtkm/filter/particleadvection/DataSetIntegrator.h> #include <vtkm/filter/particleadvection/DataSetIntegrator.h>
#include <vtkm/filter/particleadvection/StreamlineAlgorithm.h> #include <vtkm/filter/particleadvection/ParticleAdvectionAlgorithm.h>
namespace vtkm namespace vtkm
{ {
@ -27,51 +24,24 @@ namespace filter
//----------------------------------------------------------------------------- //-----------------------------------------------------------------------------
inline VTKM_CONT Pathline::Pathline() inline VTKM_CONT Pathline::Pathline()
: vtkm::filter::FilterDataSetWithField<Pathline>() : vtkm::filter::FilterTemporalParticleAdvection<Pathline>()
, UseThreadedAlgorithm(false)
{ {
} }
//-----------------------------------------------------------------------------
inline VTKM_CONT void Pathline::SetSeeds(vtkm::cont::ArrayHandle<vtkm::Particle>& seeds)
{
this->Seeds = seeds;
}
//----------------------------------------------------------------------------- //-----------------------------------------------------------------------------
template <typename DerivedPolicy> template <typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::PartitionedDataSet Pathline::PrepareForExecution( inline VTKM_CONT vtkm::cont::PartitionedDataSet Pathline::PrepareForExecution(
const vtkm::cont::PartitionedDataSet& input, const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>&) const vtkm::filter::PolicyBase<DerivedPolicy>&)
{ {
if (this->GetUseCoordinateSystemAsField())
throw vtkm::cont::ErrorFilterExecution("Coordinate system as field not supported");
if (this->Seeds.GetNumberOfValues() == 0)
throw vtkm::cont::ErrorFilterExecution("No seeds provided.");
if (this->NextDataSet.GetNumberOfPartitions() != input.GetNumberOfPartitions())
throw vtkm::cont::ErrorFilterExecution("Number of partitions do not match");
if (!(this->PreviousTime < this->NextTime))
throw vtkm::cont::ErrorFilterExecution("Previous time must be less than Next time.");
std::string activeField = this->GetActiveFieldName();
vtkm::filter::particleadvection::BoundsMap boundsMap(input);
using DSIType = vtkm::filter::particleadvection::TemporalDataSetIntegrator;
std::vector<DSIType> dsi;
for (vtkm::Id i = 0; i < input.GetNumberOfPartitions(); i++)
{
vtkm::Id blockId = boundsMap.GetLocalBlockId(i);
auto dsPrev = input.GetPartition(i);
auto dsNext = this->NextDataSet.GetPartition(i);
if (!dsPrev.HasPointField(activeField) || !dsNext.HasPointField(activeField))
throw vtkm::cont::ErrorFilterExecution("Unsupported field assocation");
dsi.push_back(
DSIType(dsPrev, this->PreviousTime, dsNext, this->NextTime, blockId, activeField));
}
using AlgorithmType = vtkm::filter::particleadvection::PathlineAlgorithm; using AlgorithmType = vtkm::filter::particleadvection::PathlineAlgorithm;
using ThreadedAlgorithmType = vtkm::filter::particleadvection::PathlineThreadedAlgorithm; using ThreadedAlgorithmType = vtkm::filter::particleadvection::PathlineThreadedAlgorithm;
this->ValidateOptions(input);
vtkm::filter::particleadvection::BoundsMap boundsMap(input);
auto dsi = this->CreateDataSetIntegrators(input, boundsMap);
if (this->GetUseThreadedAlgorithm()) if (this->GetUseThreadedAlgorithm())
return vtkm::filter::particleadvection::RunAlgo<DSIType, ThreadedAlgorithmType>( return vtkm::filter::particleadvection::RunAlgo<DSIType, ThreadedAlgorithmType>(
boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds); boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds);
@ -80,14 +50,6 @@ inline VTKM_CONT vtkm::cont::PartitionedDataSet Pathline::PrepareForExecution(
boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds); boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds);
} }
//-----------------------------------------------------------------------------
template <typename DerivedPolicy>
inline VTKM_CONT bool Pathline::MapFieldOntoOutput(vtkm::cont::DataSet&,
const vtkm::cont::Field&,
vtkm::filter::PolicyBase<DerivedPolicy>)
{
return false;
}
} }
} // namespace vtkm::filter } // namespace vtkm::filter
#endif #endif

@ -11,55 +11,28 @@
#ifndef vtk_m_filter_Streamline_h #ifndef vtk_m_filter_Streamline_h
#define vtk_m_filter_Streamline_h #define vtk_m_filter_Streamline_h
#include <vtkm/Particle.h> #include <vtkm/filter/FilterParticleAdvection.h>
#include <vtkm/filter/FilterDataSetWithField.h>
namespace vtkm namespace vtkm
{ {
namespace filter namespace filter
{ {
/// \brief generate streamlines from a vector field. /// \brief Generate streamlines from a vector field.
/// Takes as input a vector field and seed locations and generates the /// Takes as input a vector field and seed locations and generates the
/// paths taken by the seeds through the vector field. /// paths taken by the seeds through the vector field.
class Streamline : public vtkm::filter::FilterDataSetWithField<Streamline> class Streamline : public vtkm::filter::FilterParticleAdvection<Streamline>
{ {
public: public:
using SupportedTypes = vtkm::TypeListFieldVec3;
VTKM_CONT VTKM_CONT
Streamline(); Streamline();
VTKM_CONT
void SetStepSize(vtkm::FloatDefault s) { this->StepSize = s; }
VTKM_CONT
void SetNumberOfSteps(vtkm::Id n) { this->NumberOfSteps = n; }
VTKM_CONT
void SetSeeds(vtkm::cont::ArrayHandle<vtkm::Particle>& seeds);
VTKM_CONT
bool GetUseThreadedAlgorithm() { return this->UseThreadedAlgorithm; }
VTKM_CONT
void SetUseThreadedAlgorithm(bool val) { this->UseThreadedAlgorithm = val; }
template <typename DerivedPolicy> template <typename DerivedPolicy>
vtkm::cont::PartitionedDataSet PrepareForExecution( vtkm::cont::PartitionedDataSet PrepareForExecution(
const vtkm::cont::PartitionedDataSet& input, const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy); const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename DerivedPolicy>
VTKM_CONT bool MapFieldOntoOutput(vtkm::cont::DataSet& result,
const vtkm::cont::Field& field,
vtkm::filter::PolicyBase<DerivedPolicy> policy);
private: private:
vtkm::Id NumberOfSteps;
vtkm::FloatDefault StepSize;
vtkm::cont::ArrayHandle<vtkm::Particle> Seeds;
bool UseThreadedAlgorithm;
}; };
} }
} // namespace vtkm::filter } // namespace vtkm::filter

@ -10,16 +10,11 @@
#ifndef vtk_m_filter_Streamline_hxx #ifndef vtk_m_filter_Streamline_hxx
#define vtk_m_filter_Streamline_hxx #define vtk_m_filter_Streamline_hxx
#include <vtkm/filter/Streamline.h>
#include <vtkm/cont/ArrayCopy.h>
#include <vtkm/cont/ArrayHandleIndex.h>
#include <vtkm/cont/ErrorFilterExecution.h> #include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/cont/ParticleArrayCopy.h> #include <vtkm/filter/Streamline.h>
#include <vtkm/filter/particleadvection/BoundsMap.h> #include <vtkm/filter/particleadvection/BoundsMap.h>
#include <vtkm/filter/particleadvection/DataSetIntegrator.h> #include <vtkm/filter/particleadvection/DataSetIntegrator.h>
#include <vtkm/filter/particleadvection/ParticleAdvectionAlgorithm.h>
#include <vtkm/filter/particleadvection/StreamlineAlgorithm.h>
namespace vtkm namespace vtkm
{ {
@ -28,45 +23,23 @@ namespace filter
//----------------------------------------------------------------------------- //-----------------------------------------------------------------------------
inline VTKM_CONT Streamline::Streamline() inline VTKM_CONT Streamline::Streamline()
: vtkm::filter::FilterDataSetWithField<Streamline>() : vtkm::filter::FilterParticleAdvection<Streamline>()
, UseThreadedAlgorithm(false)
{ {
} }
//-----------------------------------------------------------------------------
inline VTKM_CONT void Streamline::SetSeeds(vtkm::cont::ArrayHandle<vtkm::Particle>& seeds)
{
this->Seeds = seeds;
}
//----------------------------------------------------------------------------- //-----------------------------------------------------------------------------
template <typename DerivedPolicy> template <typename DerivedPolicy>
inline VTKM_CONT vtkm::cont::PartitionedDataSet Streamline::PrepareForExecution( inline VTKM_CONT vtkm::cont::PartitionedDataSet Streamline::PrepareForExecution(
const vtkm::cont::PartitionedDataSet& input, const vtkm::cont::PartitionedDataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>&) const vtkm::filter::PolicyBase<DerivedPolicy>&)
{ {
if (this->GetUseCoordinateSystemAsField())
throw vtkm::cont::ErrorFilterExecution("Coordinate system as field not supported");
if (this->Seeds.GetNumberOfValues() == 0)
throw vtkm::cont::ErrorFilterExecution("No seeds provided.");
std::string activeField = this->GetActiveFieldName();
vtkm::filter::particleadvection::BoundsMap boundsMap(input);
using DSIType = vtkm::filter::particleadvection::DataSetIntegrator;
std::vector<DSIType> dsi;
for (vtkm::Id i = 0; i < input.GetNumberOfPartitions(); i++)
{
vtkm::Id blockId = boundsMap.GetLocalBlockId(i);
auto ds = input.GetPartition(i);
if (!ds.HasPointField(activeField))
throw vtkm::cont::ErrorFilterExecution("Unsupported field assocation");
dsi.push_back(DSIType(ds, blockId, activeField));
}
using AlgorithmType = vtkm::filter::particleadvection::StreamlineAlgorithm; using AlgorithmType = vtkm::filter::particleadvection::StreamlineAlgorithm;
using ThreadedAlgorithmType = vtkm::filter::particleadvection::StreamlineThreadedAlgorithm; using ThreadedAlgorithmType = vtkm::filter::particleadvection::StreamlineThreadedAlgorithm;
this->ValidateOptions();
vtkm::filter::particleadvection::BoundsMap boundsMap(input);
auto dsi = this->CreateDataSetIntegrators(input, boundsMap);
if (this->GetUseThreadedAlgorithm()) if (this->GetUseThreadedAlgorithm())
return vtkm::filter::particleadvection::RunAlgo<DSIType, ThreadedAlgorithmType>( return vtkm::filter::particleadvection::RunAlgo<DSIType, ThreadedAlgorithmType>(
boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds); boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds);
@ -75,14 +48,6 @@ inline VTKM_CONT vtkm::cont::PartitionedDataSet Streamline::PrepareForExecution(
boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds); boundsMap, dsi, this->NumberOfSteps, this->StepSize, this->Seeds);
} }
//-----------------------------------------------------------------------------
template <typename DerivedPolicy>
inline VTKM_CONT bool Streamline::MapFieldOntoOutput(vtkm::cont::DataSet&,
const vtkm::cont::Field&,
vtkm::filter::PolicyBase<DerivedPolicy>)
{
return false;
}
} }
} // namespace vtkm::filter } // namespace vtkm::filter
#endif #endif

@ -72,8 +72,8 @@ public:
int FindRank(vtkm::Id blockId) const int FindRank(vtkm::Id blockId) const
{ {
auto it = BlockToRankMap.find(blockId); auto it = this->BlockToRankMap.find(blockId);
if (it == BlockToRankMap.end()) if (it == this->BlockToRankMap.end())
return -1; return -1;
return it->second; return it->second;
} }
@ -104,6 +104,9 @@ public:
return blockIDs; return blockIDs;
} }
vtkm::Id GetTotalNumBlocks() const { return this->TotalNumBlocks; }
vtkm::Id GetLocalNumBlocks() const { return this->LocalNumBlocks; }
private: private:
void Init(const std::vector<vtkm::cont::DataSet>& dataSets) void Init(const std::vector<vtkm::cont::DataSet>& dataSets)
{ {

@ -18,7 +18,6 @@ set(headers
Messenger.h Messenger.h
ParticleAdvectionAlgorithm.h ParticleAdvectionAlgorithm.h
ParticleMessenger.h ParticleMessenger.h
StreamlineAlgorithm.h
) )
#----------------------------------------------------------------------------- #-----------------------------------------------------------------------------

@ -34,12 +34,12 @@ template <>
template <> template <>
inline void DataSetIntegratorBase<GridEvalType>::DoAdvect( inline void DataSetIntegratorBase<GridEvalType>::DoAdvect(
vtkm::cont::ArrayHandle<vtkm::Particle>& seeds, vtkm::cont::ArrayHandle<vtkm::Particle>& seeds,
const Stepper& rk4, const Stepper& stepper,
vtkm::Id maxSteps, vtkm::Id maxSteps,
vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>& result) const vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>& result) const
{ {
vtkm::worklet::ParticleAdvection Worklet; vtkm::worklet::ParticleAdvection Worklet;
result = Worklet.Run(rk4, seeds, maxSteps); result = Worklet.Run(stepper, seeds, maxSteps);
} }
//----- //-----
@ -48,12 +48,26 @@ template <>
template <> template <>
inline void DataSetIntegratorBase<GridEvalType>::DoAdvect( inline void DataSetIntegratorBase<GridEvalType>::DoAdvect(
vtkm::cont::ArrayHandle<vtkm::Particle>& seeds, vtkm::cont::ArrayHandle<vtkm::Particle>& seeds,
const Stepper& rk4, const Stepper& stepper,
vtkm::Id maxSteps, vtkm::Id maxSteps,
vtkm::worklet::StreamlineResult<vtkm::Particle>& result) const vtkm::worklet::StreamlineResult<vtkm::Particle>& result) const
{ {
vtkm::worklet::Streamline Worklet; vtkm::worklet::Streamline Worklet;
result = Worklet.Run(rk4, seeds, maxSteps); result = Worklet.Run(stepper, seeds, maxSteps);
}
//-----
// Specialization for PathParticle worklet
template <>
template <>
inline void DataSetIntegratorBase<TemporalGridEvalType>::DoAdvect(
vtkm::cont::ArrayHandle<vtkm::Particle>& seeds,
const TemporalStepper& stepper,
vtkm::Id maxSteps,
vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>& result) const
{
vtkm::worklet::ParticleAdvection Worklet;
result = Worklet.Run(stepper, seeds, maxSteps);
} }
//----- //-----
@ -62,12 +76,12 @@ template <>
template <> template <>
inline void DataSetIntegratorBase<TemporalGridEvalType>::DoAdvect( inline void DataSetIntegratorBase<TemporalGridEvalType>::DoAdvect(
vtkm::cont::ArrayHandle<vtkm::Particle>& seeds, vtkm::cont::ArrayHandle<vtkm::Particle>& seeds,
const TemporalStepper& rk4, const TemporalStepper& stepper,
vtkm::Id maxSteps, vtkm::Id maxSteps,
vtkm::worklet::StreamlineResult<vtkm::Particle>& result) const vtkm::worklet::StreamlineResult<vtkm::Particle>& result) const
{ {
vtkm::worklet::Streamline Worklet; vtkm::worklet::Streamline Worklet;
result = Worklet.Run(rk4, seeds, maxSteps); result = Worklet.Run(stepper, seeds, maxSteps);
} }
} }

@ -23,7 +23,12 @@ namespace particleadvection
{ {
using DSIType = vtkm::filter::particleadvection::DataSetIntegrator; using DSIType = vtkm::filter::particleadvection::DataSetIntegrator;
using TDSIType = vtkm::filter::particleadvection::TemporalDataSetIntegrator;
//-------------------------------------------------------------------------------------------
//Steady state advection algorithms
//Particle advection
class VTKM_ALWAYS_EXPORT ParticleAdvectionAlgorithm class VTKM_ALWAYS_EXPORT ParticleAdvectionAlgorithm
: public AdvectorBaseAlgorithm<DSIType, vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>> : public AdvectorBaseAlgorithm<DSIType, vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>>
{ {
@ -36,6 +41,7 @@ public:
} }
}; };
//Threaded particle advection
class VTKM_ALWAYS_EXPORT ParticleAdvectionThreadedAlgorithm class VTKM_ALWAYS_EXPORT ParticleAdvectionThreadedAlgorithm
: public AdvectorBaseThreadedAlgorithm<DSIType, : public AdvectorBaseThreadedAlgorithm<DSIType,
vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>> vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>>
@ -50,6 +56,94 @@ public:
} }
}; };
//Streamline
class VTKM_ALWAYS_EXPORT StreamlineAlgorithm
: public AdvectorBaseAlgorithm<DSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>
{
public:
StreamlineAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<DSIType>& blocks)
: AdvectorBaseAlgorithm<DSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>(bm, blocks)
{
}
};
//Threaded streamline
class VTKM_ALWAYS_EXPORT StreamlineThreadedAlgorithm
: public AdvectorBaseThreadedAlgorithm<DSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>
{
public:
StreamlineThreadedAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<DSIType>& blocks)
: AdvectorBaseThreadedAlgorithm<DSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>(
bm,
blocks)
{
}
};
//-------------------------------------------------------------------------------------------
//Unsteady state advection algorithms
//PathParticle
class VTKM_ALWAYS_EXPORT PathParticleAlgorithm
: public AdvectorBaseAlgorithm<TDSIType, vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>>
{
public:
PathParticleAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<TDSIType>& blocks)
: AdvectorBaseAlgorithm<TDSIType, vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>>(
bm,
blocks)
{
}
};
//Threaded PathParticle
class VTKM_ALWAYS_EXPORT PathParticleThreadedAlgorithm
: public AdvectorBaseThreadedAlgorithm<TDSIType,
vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>>
{
public:
PathParticleThreadedAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<TDSIType>& blocks)
: AdvectorBaseThreadedAlgorithm<TDSIType,
vtkm::worklet::ParticleAdvectionResult<vtkm::Particle>>(bm,
blocks)
{
}
};
//Pathline
class VTKM_ALWAYS_EXPORT PathlineAlgorithm
: public AdvectorBaseAlgorithm<TDSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>
{
public:
PathlineAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<TDSIType>& blocks)
: AdvectorBaseAlgorithm<TDSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>(bm, blocks)
{
}
};
//Threaded pathline
class VTKM_ALWAYS_EXPORT PathlineThreadedAlgorithm
: public AdvectorBaseThreadedAlgorithm<TDSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>
{
public:
PathlineThreadedAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<TDSIType>& blocks)
: AdvectorBaseThreadedAlgorithm<TDSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>(
bm,
blocks)
{
}
};
} }
} }
} // namespace vtkm::filter::particleadvection } // namespace vtkm::filter::particleadvection

@ -1,82 +0,0 @@
//============================================================================
// 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_filter_particleadvection_StreamlineAlgorithm_h
#define vtk_m_filter_particleadvection_StreamlineAlgorithm_h
#include <vtkm/filter/particleadvection/AdvectorBaseAlgorithm.h>
#include <vtkm/filter/particleadvection/AdvectorBaseThreadedAlgorithm.h>
#include <vtkm/filter/particleadvection/DataSetIntegrator.h>
namespace vtkm
{
namespace filter
{
namespace particleadvection
{
using DSIType = vtkm::filter::particleadvection::DataSetIntegrator;
using TDSIType = vtkm::filter::particleadvection::TemporalDataSetIntegrator;
class VTKM_ALWAYS_EXPORT StreamlineAlgorithm
: public AdvectorBaseAlgorithm<DSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>
{
public:
StreamlineAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<DSIType>& blocks)
: AdvectorBaseAlgorithm<DSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>(bm, blocks)
{
}
};
class VTKM_ALWAYS_EXPORT StreamlineThreadedAlgorithm
: public AdvectorBaseThreadedAlgorithm<DSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>
{
public:
StreamlineThreadedAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<DSIType>& blocks)
: AdvectorBaseThreadedAlgorithm<DSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>(
bm,
blocks)
{
}
};
//pathline
class VTKM_ALWAYS_EXPORT PathlineAlgorithm
: public AdvectorBaseAlgorithm<TDSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>
{
public:
PathlineAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<TDSIType>& blocks)
: AdvectorBaseAlgorithm<TDSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>(bm, blocks)
{
}
};
class VTKM_ALWAYS_EXPORT PathlineThreadedAlgorithm
: public AdvectorBaseThreadedAlgorithm<TDSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>
{
public:
PathlineThreadedAlgorithm(const vtkm::filter::particleadvection::BoundsMap& bm,
const std::vector<TDSIType>& blocks)
: AdvectorBaseThreadedAlgorithm<TDSIType, vtkm::worklet::StreamlineResult<vtkm::Particle>>(
bm,
blocks)
{
}
};
}
}
} // namespace vtkm::filter::particleadvection
#endif //vtk_m_filter_particleadvection_StreamlineAlgorithm_h

@ -11,6 +11,7 @@
#include <vtkm/cont/DataSetBuilderUniform.h> #include <vtkm/cont/DataSetBuilderUniform.h>
#include <vtkm/cont/testing/Testing.h> #include <vtkm/cont/testing/Testing.h>
#include <vtkm/filter/ParticleAdvection.h> #include <vtkm/filter/ParticleAdvection.h>
#include <vtkm/filter/PathParticle.h>
#include <vtkm/filter/Pathline.h> #include <vtkm/filter/Pathline.h>
#include <vtkm/filter/Streamline.h> #include <vtkm/filter/Streamline.h>
#include <vtkm/io/VTKDataSetReader.h> #include <vtkm/io/VTKDataSetReader.h>
@ -19,6 +20,14 @@
namespace namespace
{ {
enum FilterType
{
PARTICLE_ADVECTION,
STREAMLINE,
PATHLINE,
PATH_PARTICLE
};
vtkm::cont::ArrayHandle<vtkm::Vec3f> CreateConstantVectorField(vtkm::Id num, const vtkm::Vec3f& vec) vtkm::cont::ArrayHandle<vtkm::Vec3f> CreateConstantVectorField(vtkm::Id num, const vtkm::Vec3f& vec)
{ {
vtkm::cont::ArrayHandleConstant<vtkm::Vec3f> vecConst; vtkm::cont::ArrayHandleConstant<vtkm::Vec3f> vecConst;
@ -75,52 +84,17 @@ void TestStreamline()
} }
} }
void TestPathlineSimple()
{
vtkm::cont::DataSetBuilderUniform dataSetBuilder;
vtkm::cont::DataSet inData1 = dataSetBuilder.Create(vtkm::Id3(5, 5, 5));
vtkm::cont::DataSet inData2 = dataSetBuilder.Create(vtkm::Id3(5, 5, 5));
vtkm::Id numPoints = inData1.GetCellSet().GetNumberOfPoints();
vtkm::cont::ArrayHandle<vtkm::Vec3f> vectorField1;
vtkm::cont::ArrayCopy(vtkm::cont::make_ArrayHandleConstant(vtkm::Vec3f(1, 0, 0), numPoints),
vectorField1);
inData1.AddPointField("vectorvar", vectorField1);
vtkm::cont::ArrayHandle<vtkm::Vec3f> vectorField2;
vtkm::cont::ArrayCopy(vtkm::cont::make_ArrayHandleConstant(vtkm::Vec3f(0, 1, 0), numPoints),
vectorField2);
inData2.AddPointField("vectorvar", vectorField2);
vtkm::filter::Pathline pathlines;
// Specify the seeds.
vtkm::cont::ArrayHandle<vtkm::Particle> seedArray;
seedArray.Allocate(2);
seedArray.WritePortal().Set(0, vtkm::Particle({ 0, 0, 0 }, 0));
seedArray.WritePortal().Set(1, vtkm::Particle({ 1, 1, 1 }, 1));
pathlines.SetActiveField("vectorvar");
pathlines.SetStepSize(0.1f);
pathlines.SetNumberOfSteps(100);
pathlines.SetSeeds(seedArray);
pathlines.SetPreviousTime(0.0f);
pathlines.SetNextTime(1.0f);
pathlines.SetNextDataSet(inData2);
auto output = pathlines.Execute(inData1);
//Validate the result is correct.
vtkm::cont::CoordinateSystem coords = output.GetCoordinateSystem();
VTKM_TEST_ASSERT(coords.GetNumberOfPoints() == 77, "Wrong number of coordinates");
vtkm::cont::DynamicCellSet dcells = output.GetCellSet();
VTKM_TEST_ASSERT(dcells.GetNumberOfCells() == 2, "Wrong number of cells");
}
void TestPathline() void TestPathline()
{ {
const vtkm::Id3 dims(5, 5, 5); const vtkm::Id3 dims(5, 5, 5);
const vtkm::Vec3f vecX(1, 0, 0); const vtkm::Vec3f vecX(1, 0, 0);
const vtkm::Vec3f vecY(0, 1, 0); const vtkm::Vec3f vecY(0, 1, 0);
const vtkm::Bounds bounds(0, 4, 0, 4, 0, 4); const vtkm::Bounds bounds(0, 4, 0, 4, 0, 4);
std::string fieldName = "vec"; std::string var = "vec";
//test pathline and pathparticle filters.
for (int fType = 0; fType < 2; fType++)
{
auto dataSets1 = vtkm::worklet::testing::CreateAllDataSets(bounds, dims, false); auto dataSets1 = vtkm::worklet::testing::CreateAllDataSets(bounds, dims, false);
auto dataSets2 = vtkm::worklet::testing::CreateAllDataSets(bounds, dims, false); auto dataSets2 = vtkm::worklet::testing::CreateAllDataSets(bounds, dims, false);
@ -132,34 +106,57 @@ void TestPathline()
auto vecField1 = CreateConstantVectorField(ds1.GetNumberOfPoints(), vecX); auto vecField1 = CreateConstantVectorField(ds1.GetNumberOfPoints(), vecX);
auto vecField2 = CreateConstantVectorField(ds1.GetNumberOfPoints(), vecY); auto vecField2 = CreateConstantVectorField(ds1.GetNumberOfPoints(), vecY);
ds1.AddPointField(fieldName, vecField1); ds1.AddPointField(var, vecField1);
ds2.AddPointField(fieldName, vecField2); ds2.AddPointField(var, vecField2);
vtkm::cont::ArrayHandle<vtkm::Particle> seedArray = vtkm::cont::ArrayHandle<vtkm::Particle> seedArray =
vtkm::cont::make_ArrayHandle({ vtkm::Particle(vtkm::Vec3f(.2f, 1.0f, .2f), 0), vtkm::cont::make_ArrayHandle({ vtkm::Particle(vtkm::Vec3f(.2f, 1.0f, .2f), 0),
vtkm::Particle(vtkm::Vec3f(.2f, 2.0f, .2f), 1), vtkm::Particle(vtkm::Vec3f(.2f, 2.0f, .2f), 1),
vtkm::Particle(vtkm::Vec3f(.2f, 3.0f, .2f), 2) }); vtkm::Particle(vtkm::Vec3f(.2f, 3.0f, .2f), 2) });
vtkm::filter::Pathline pathline; const vtkm::FloatDefault stepSize = 0.1;
const vtkm::FloatDefault t0 = 0, t1 = 1;
const vtkm::Id numSteps = 20;
pathline.SetPreviousTime(0.0f); vtkm::cont::DataSet output;
pathline.SetNextTime(1.0f); vtkm::Id numExpectedPoints;
pathline.SetNextDataSet(ds2); if (fType == 0)
pathline.SetStepSize(static_cast<vtkm::FloatDefault>(0.05f)); {
pathline.SetNumberOfSteps(20); vtkm::filter::Pathline filt;
pathline.SetSeeds(seedArray); filt.SetActiveField(var);
filt.SetStepSize(stepSize);
pathline.SetActiveField(fieldName); filt.SetNumberOfSteps(numSteps);
auto output = pathline.Execute(ds1); filt.SetSeeds(seedArray);
filt.SetPreviousTime(t0);
filt.SetNextTime(t1);
filt.SetNextDataSet(ds2);
output = filt.Execute(ds1);
numExpectedPoints = 63;
}
else
{
vtkm::filter::PathParticle filt;
filt.SetActiveField(var);
filt.SetStepSize(stepSize);
filt.SetNumberOfSteps(numSteps);
filt.SetSeeds(seedArray);
filt.SetPreviousTime(t0);
filt.SetNextTime(t1);
filt.SetNextDataSet(ds2);
output = filt.Execute(ds1);
numExpectedPoints = 3;
}
//Validate the result is correct. //Validate the result is correct.
vtkm::cont::CoordinateSystem coords = output.GetCoordinateSystem(); vtkm::cont::CoordinateSystem coords = output.GetCoordinateSystem();
VTKM_TEST_ASSERT(coords.GetNumberOfPoints() == 63, "Wrong number of coordinates"); VTKM_TEST_ASSERT(coords.GetNumberOfPoints() == numExpectedPoints,
"Wrong number of coordinates");
vtkm::cont::DynamicCellSet dcells = output.GetCellSet(); vtkm::cont::DynamicCellSet dcells = output.GetCellSet();
VTKM_TEST_ASSERT(dcells.GetNumberOfCells() == 3, "Wrong number of cells"); VTKM_TEST_ASSERT(dcells.GetNumberOfCells() == 3, "Wrong number of cells");
} }
} }
}
void TestAMRStreamline(bool useSL) void TestAMRStreamline(bool useSL)
{ {
@ -225,7 +222,7 @@ void TestAMRStreamline(bool useSL)
{ {
vtkm::filter::Streamline filter; vtkm::filter::Streamline filter;
filter.SetStepSize(0.1f); filter.SetStepSize(0.1f);
filter.SetNumberOfSteps(1000); filter.SetNumberOfSteps(100000);
filter.SetSeeds(seedArray); filter.SetSeeds(seedArray);
filter.SetActiveField(fieldName); filter.SetActiveField(fieldName);
@ -297,7 +294,7 @@ void TestAMRStreamline(bool useSL)
{ {
vtkm::filter::ParticleAdvection filter; vtkm::filter::ParticleAdvection filter;
filter.SetStepSize(0.1f); filter.SetStepSize(0.1f);
filter.SetNumberOfSteps(1000); filter.SetNumberOfSteps(100000);
filter.SetSeeds(seedArray); filter.SetSeeds(seedArray);
filter.SetActiveField(fieldName); filter.SetActiveField(fieldName);
@ -325,7 +322,7 @@ void TestAMRStreamline(bool useSL)
} }
} }
void TestPartitionedDataSet(vtkm::Id num, bool useGhost, bool useSL) void TestPartitionedDataSet(vtkm::Id num, bool useGhost, FilterType fType)
{ {
vtkm::Id numDims = 5; vtkm::Id numDims = 5;
vtkm::FloatDefault x0 = 0; vtkm::FloatDefault x0 = 0;
@ -354,14 +351,17 @@ void TestPartitionedDataSet(vtkm::Id num, bool useGhost, bool useSL)
x1 += dx; x1 += dx;
} }
std::vector<vtkm::cont::PartitionedDataSet> allPDs; std::vector<vtkm::cont::PartitionedDataSet> allPDs, allPDs2;
const vtkm::Id3 dims(numDims, numDims, numDims); const vtkm::Id3 dims(numDims, numDims, numDims);
allPDs = vtkm::worklet::testing::CreateAllDataSets(bounds, dims, useGhost); allPDs = vtkm::worklet::testing::CreateAllDataSets(bounds, dims, useGhost);
if (fType == FilterType::PATHLINE || fType == FilterType::PATH_PARTICLE)
allPDs2 = vtkm::worklet::testing::CreateAllDataSets(bounds, dims, useGhost);
vtkm::Vec3f vecX(1, 0, 0); vtkm::Vec3f vecX(1, 0, 0);
std::string fieldName = "vec"; std::string fieldName = "vec";
for (auto& pds : allPDs) for (std::size_t idx = 0; idx < allPDs.size(); idx++)
{ {
auto pds = allPDs[idx];
AddVectorFields(pds, fieldName, vecX); AddVectorFields(pds, fieldName, vecX);
vtkm::cont::ArrayHandle<vtkm::Particle> seedArray; vtkm::cont::ArrayHandle<vtkm::Particle> seedArray;
@ -369,20 +369,38 @@ void TestPartitionedDataSet(vtkm::Id num, bool useGhost, bool useSL)
vtkm::Particle(vtkm::Vec3f(.2f, 2.0f, .2f), 1) }); vtkm::Particle(vtkm::Vec3f(.2f, 2.0f, .2f), 1) });
vtkm::Id numSeeds = seedArray.GetNumberOfValues(); vtkm::Id numSeeds = seedArray.GetNumberOfValues();
if (useSL) if (fType == FilterType::STREAMLINE || fType == FilterType::PATHLINE)
{
vtkm::cont::PartitionedDataSet out;
if (fType == FilterType::STREAMLINE)
{ {
vtkm::filter::Streamline streamline; vtkm::filter::Streamline streamline;
streamline.SetStepSize(0.1f); streamline.SetStepSize(0.1f);
streamline.SetNumberOfSteps(1000); streamline.SetNumberOfSteps(100000);
streamline.SetSeeds(seedArray); streamline.SetSeeds(seedArray);
streamline.SetActiveField(fieldName); streamline.SetActiveField(fieldName);
auto out = streamline.Execute(pds); out = streamline.Execute(pds);
}
else
{
auto pds2 = allPDs2[idx];
AddVectorFields(pds2, fieldName, vecX);
vtkm::filter::Pathline pathline;
pathline.SetPreviousTime(0);
pathline.SetNextTime(1000);
pathline.SetNextDataSet(pds2);
pathline.SetStepSize(0.1f);
pathline.SetNumberOfSteps(100000);
pathline.SetSeeds(seedArray);
pathline.SetActiveField(fieldName);
out = pathline.Execute(pds);
}
for (vtkm::Id i = 0; i < num; i++) for (vtkm::Id i = 0; i < num; i++)
{ {
auto inputDS = pds.GetPartition(i);
auto outputDS = out.GetPartition(i); auto outputDS = out.GetPartition(i);
VTKM_TEST_ASSERT(outputDS.GetNumberOfCoordinateSystems() == 1, VTKM_TEST_ASSERT(outputDS.GetNumberOfCoordinateSystems() == 1,
"Wrong number of coordinate systems in the output dataset"); "Wrong number of coordinate systems in the output dataset");
@ -415,16 +433,36 @@ void TestPartitionedDataSet(vtkm::Id num, bool useGhost, bool useSL)
} }
} }
} }
else else if (fType == FilterType::PARTICLE_ADVECTION || fType == FilterType::PATH_PARTICLE)
{
vtkm::cont::PartitionedDataSet out;
if (fType == FilterType::PARTICLE_ADVECTION)
{ {
vtkm::filter::ParticleAdvection particleAdvection; vtkm::filter::ParticleAdvection particleAdvection;
particleAdvection.SetStepSize(0.1f); particleAdvection.SetStepSize(0.1f);
particleAdvection.SetNumberOfSteps(1000); particleAdvection.SetNumberOfSteps(100000);
particleAdvection.SetSeeds(seedArray); particleAdvection.SetSeeds(seedArray);
particleAdvection.SetActiveField(fieldName); particleAdvection.SetActiveField(fieldName);
auto out = particleAdvection.Execute(pds); out = particleAdvection.Execute(pds);
}
else
{
auto pds2 = allPDs2[idx];
AddVectorFields(pds2, fieldName, vecX);
vtkm::filter::PathParticle pathParticle;
pathParticle.SetPreviousTime(0);
pathParticle.SetNextTime(1000);
pathParticle.SetNextDataSet(pds2);
pathParticle.SetStepSize(0.1f);
pathParticle.SetNumberOfSteps(100000);
pathParticle.SetSeeds(seedArray);
pathParticle.SetActiveField(fieldName);
out = pathParticle.Execute(pds);
}
VTKM_TEST_ASSERT(out.GetNumberOfPartitions() == 1, "Wrong number of partitions in output"); VTKM_TEST_ASSERT(out.GetNumberOfPartitions() == 1, "Wrong number of partitions in output");
auto ds = out.GetPartition(0); auto ds = out.GetPartition(0);
@ -541,18 +579,21 @@ void TestStreamlineFile(const std::string& fname,
void TestStreamlineFilters() void TestStreamlineFilters()
{ {
std::vector<bool> flags = { true, false }; std::vector<bool> flags = { true, false };
std::vector<FilterType> fTypes = { FilterType::PARTICLE_ADVECTION,
FilterType::STREAMLINE,
FilterType::PATHLINE,
FilterType::PATH_PARTICLE };
for (int n = 1; n < 3; n++) for (int n = 1; n < 3; n++)
{ {
for (auto useGhost : flags) for (auto useGhost : flags)
for (auto useSL : flags) for (auto ft : fTypes)
{ {
useSL = false; TestPartitionedDataSet(n, useGhost, ft);
TestPartitionedDataSet(n, useGhost, useSL);
} }
} }
TestStreamline(); TestStreamline();
TestPathlineSimple();
TestPathline(); TestPathline();
for (auto useSL : flags) for (auto useSL : flags)