Shuenhoy/vtk-m
1
0
Fork 0
mirror of https://gitlab.kitware.com/vtk/vtk-m synced 2024-09-16 17:22:55 +00:00
Code Issues 2 Actions Packages Projects Releases Wiki Activity

Merge topic 'filter-enhancements'

Browse Source 
d37c1fee adding mpi-only filter example.
4610e8b1 increate timeout for CopyrightStatement test.
89041e4f Compute correct histogram for distributed multiblock cases.
0af0271f Filter API enhancements.
b31af29a New constructor for AssignerMultiBlock.
dd3709c3 make MultiBlock constructors explicit.
7fbd0204 FieldHistogram: support using precomputed range
ddd14f25 mpi: avoid MPI calls when not using MPI.
...

Acked-by: Kitware Robot <kwrobot@kitware.com>
Acked-by: Kenneth Moreland <kmorel@sandia.gov>
Acked-by: Matt Larsen <mlarsen@cs.uoregon.edu>
Merge-request: !1138
This commit is contained in:
Utkarsh Ayachit 2018-04-13 20:31:03 +00:00 committed by Kitware Robot
commit 94db252c62
24 changed files with 1347 additions and 106 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
CMakeLists.txt
View File

@ -169,6 +169,9 @@ if (VTKm_ENABLE_TESTING)
add_test(NAME CopyrightStatement
COMMAND ${CMAKE_COMMAND} "-DVTKm_SOURCE_DIR=${VTKm_SOURCE_DIR}" -P "${VTKm_SOURCE_DIR}/CMake/VTKmCheckCopyright.cmake"
)
# increase timeout since on some machines CopyrightStatement test takes a long time.
set_tests_properties(CopyrightStatement PROPERTIES TIMEOUT 300)
add_test(NAME SourceInBuild
COMMAND ${CMAKE_COMMAND} "-DVTKm_SOURCE_DIR=${VTKm_SOURCE_DIR}" -P "${VTKm_SOURCE_DIR}/CMake/VTKmCheckSourceInBuild.cmake"
)

1
examples/CMakeLists.txt
View File

@ -30,6 +30,7 @@ add_subdirectory(demo)
add_subdirectory(dynamic_dispatcher)
add_subdirectory(game_of_life)
add_subdirectory(hello_world)
add_subdirectory(histogram)
add_subdirectory(isosurface)
add_subdirectory(multi_backend)
add_subdirectory(particle_advection)

48
examples/histogram/CMakeLists.txt Normal file
View File

@ -0,0 +1,48 @@
##=============================================================================
##
## 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.
##
## Copyright 2015 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
## Copyright 2015 UT-Battelle, LLC.
## Copyright 2015 Los Alamos National Security.
##
## Under the terms of Contract DE-NA0003525 with NTESS,
## 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.
##
##=============================================================================
cmake_minimum_required(VERSION 3.3 FATAL_ERROR)
project(Histogram CXX)
#Find the VTK-m package
find_package(VTKm REQUIRED QUIET)
if (VTKm_ENABLE_MPI)
add_executable(Histogram_SERIAL Histogram.cxx HistogramMPI.h HistogramMPI.hxx)
target_compile_definitions(Histogram_SERIAL PRIVATE
"VTKM_DEVICE_ADAPTER=VTKM_DEVICE_ADAPTER_SERIAL")
target_link_libraries(Histogram_SERIAL PRIVATE vtkm_cont)
if(TARGET vtkm::tbb)
add_executable(Histogram_TBB Histogram.cxx HistogramMPI.h HistogramMPI.hxx)
target_compile_definitions(Histogram_TBB PRIVATE
"VTKM_DEVICE_ADAPTER=VTKM_DEVICE_ADAPTER_TBB")
target_link_libraries(Histogram_TBB PRIVATE vtkm_cont)
endif()
if(TARGET vtkm::cuda)
vtkm_compile_as_cuda(cudaSource Histogram.cxx HistogramMPI.h HistogramMPI.hxx)
add_executable(Histogram_CUDA ${cudaSource})
target_compile_definitions(Histogram_CUDA PRIVATE
"VTKM_DEVICE_ADAPTER=VTKM_DEVICE_ADAPTER_CUDA")
target_link_libraries(Histogram_CUDA PRIVATE vtkm_cont)
endif()
endif()

124
examples/histogram/Histogram.cxx Normal file
View File

@ -0,0 +1,124 @@
//============================================================================
// 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.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// 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.
//============================================================================
/*
* This example demonstrates how one can write a filter that uses MPI
* for hybrid-parallelism. The `vtkm::filter::Histogram` is another approach for
* implementing the same that uses DIY. This example doesn't use DIY, instead
* uses MPI calls directly.
*/
#include "HistogramMPI.h"
#include <vtkm/cont/ArrayPortalToIterators.h>
#include <vtkm/cont/DataSetFieldAdd.h>
#include <vtkm/cont/EnvironmentTracker.h>
// clang-format off
VTKM_THIRDPARTY_PRE_INCLUDE
#include VTKM_DIY(diy/mpi.hpp)
VTKM_THIRDPARTY_POST_INCLUDE
// clang-format on
#include <mpi.h>
#include <algorithm>
#include <numeric>
#include <random>
#include <utility>
#include <vector>
namespace
{
template <typename T>
VTKM_CONT vtkm::cont::ArrayHandle<T> CreateArray(T min, T max, vtkm::Id numVals)
{
std::mt19937 gen;
std::uniform_real_distribution<double> dis(static_cast<double>(min), static_cast<double>(max));
vtkm::cont::ArrayHandle<T> handle;
handle.Allocate(numVals);
std::generate(vtkm::cont::ArrayPortalToIteratorBegin(handle.GetPortalControl()),
vtkm::cont::ArrayPortalToIteratorEnd(handle.GetPortalControl()),
[&]() { return static_cast<T>(dis(gen)); });
return handle;
}
}
int main(int argc, char* argv[])
{
// setup MPI environment.
MPI_Init(&argc, &argv);
// tell VTK-m the communicator to use.
vtkm::cont::EnvironmentTracker::SetCommunicator(diy::mpi::communicator(MPI_COMM_WORLD));
int rank, size;
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
MPI_Comm_size(MPI_COMM_WORLD, &size);
if (argc != 2)
{
if (rank == 0)
{
std::cout << "Usage: " << std::endl << "$ " << argv[0] << " <num-bins>" << std::endl;
}
MPI_Finalize();
return EXIT_FAILURE;
}
const vtkm::Id num_bins = static_cast<vtkm::Id>(std::atoi(argv[1]));
const vtkm::Id numVals = 1024;
vtkm::cont::MultiBlock mb;
vtkm::cont::DataSet ds;
vtkm::cont::DataSetFieldAdd::AddPointField(ds, "pointvar", CreateArray(-1024, 1024, numVals));
mb.AddBlock(ds);
example::HistogramMPI histogram;
histogram.SetActiveField("pointvar");
histogram.SetNumberOfBins(std::max<vtkm::Id>(1, num_bins));
vtkm::cont::MultiBlock result = histogram.Execute(mb);
vtkm::cont::ArrayHandle<vtkm::Id> bins;
result.GetBlock(0).GetField("histogram").GetData().CopyTo(bins);
auto binPortal = bins.GetPortalConstControl();
if (rank == 0)
{
// print histogram.
std::cout << "Histogram (" << num_bins << ")" << std::endl;
vtkm::Id count = 0;
for (vtkm::Id cc = 0; cc < num_bins; ++cc)
{
std::cout << " bin[" << cc << "] = " << binPortal.Get(cc) << std::endl;
count += binPortal.Get(cc);
}
if (count != numVals * size)
{
std::cout << "ERROR: bins mismatched!" << std::endl;
MPI_Finalize();
return EXIT_FAILURE;
}
}
MPI_Finalize();
return EXIT_SUCCESS;
}

115
examples/histogram/HistogramMPI.h Normal file
View File

@ -0,0 +1,115 @@
//============================================================================
// 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.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// 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.
//============================================================================
#ifndef vtk_m_examples_histogram_HistogramMPI_h
#define vtk_m_examples_histogram_HistogramMPI_h
#include <vtkm/filter/FilterField.h>
#include <vtkm/filter/FilterTraits.h>
namespace example
{
/// \brief Construct the HistogramMPI of a given Field
///
/// Construct a HistogramMPI with a default of 10 bins.
///
class HistogramMPI : public vtkm::filter::FilterField<HistogramMPI>
{
public:
//Construct a HistogramMPI with a default of 10 bins
VTKM_CONT
HistogramMPI();
VTKM_CONT
void SetNumberOfBins(vtkm::Id count) { this->NumberOfBins = count; }
VTKM_CONT
vtkm::Id GetNumberOfBins() const { return this->NumberOfBins; }
//@{
/// Get/Set the range to use to generate the HistogramMPI. If range is set to
/// empty, the field's global range (computed using `vtkm::cont::FieldRangeGlobalCompute`)
/// will be used.
VTKM_CONT
void SetRange(const vtkm::Range& range) { this->Range = range; }
VTKM_CONT
const vtkm::Range& GetRange() const { return this->Range; }
//@}
/// Returns the bin delta of the last computed field.
VTKM_CONT
vtkm::Float64 GetBinDelta() const { return this->BinDelta; }
/// Returns the range used for most recent execute. If `SetRange` is used to
/// specify and non-empty range, then this will be same as the range after
/// the `Execute` call.
VTKM_CONT
vtkm::Range GetComputedRange() const { return this->ComputedRange; }
template <typename T, typename StorageType, typename DerivedPolicy, typename DeviceAdapter>
VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet& input,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMeta,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy,
const DeviceAdapter& tag);
//@{
/// when operating on vtkm::cont::MultiBlock, we
/// want to do processing across ranks as well. Just adding pre/post handles
/// for the same does the trick.
template <typename DerivedPolicy>
VTKM_CONT void PreExecute(const vtkm::cont::MultiBlock& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename DerivedPolicy>
VTKM_CONT void PostExecute(const vtkm::cont::MultiBlock& input,
vtkm::cont::MultiBlock& output,
const vtkm::filter::PolicyBase<DerivedPolicy>&);
//@}
private:
vtkm::Id NumberOfBins;
vtkm::Float64 BinDelta;
vtkm::Range ComputedRange;
vtkm::Range Range;
};
} // namespace example
namespace vtkm
{
namespace filter
{
template <>
class FilterTraits<example::HistogramMPI>
{ //currently the HistogramMPI filter only works on scalar data.
//this mainly has to do with getting the ranges for each bin
//would require returning a more complex value type
public:
using InputFieldTypeList = vtkm::TypeListTagScalarAll;
};
}
} // namespace vtkm::filter
#include "HistogramMPI.hxx"
#endif // vtk_m_filter_Histogram_h

204
examples/histogram/HistogramMPI.hxx Normal file
View File

@ -0,0 +1,204 @@
//============================================================================
// 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.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// 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/worklet/DispatcherMapField.h>
#include <vtkm/worklet/FieldHistogram.h>
#include <vtkm/cont/Algorithm.h>
#include <vtkm/cont/ArrayPortalToIterators.h>
#include <vtkm/cont/AssignerMultiBlock.h>
#include <vtkm/cont/EnvironmentTracker.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/cont/FieldRangeGlobalCompute.h>
// clang-format off
VTKM_THIRDPARTY_PRE_INCLUDE
#include VTKM_DIY(diy/mpi.hpp)
VTKM_THIRDPARTY_POST_INCLUDE
// clang-format on
namespace example
{
namespace detail
{
class DistributedHistogram
{
std::vector<vtkm::cont::ArrayHandle<vtkm::Id>> LocalBlocks;
public:
DistributedHistogram(vtkm::Id numLocalBlocks)
: LocalBlocks(static_cast<size_t>(numLocalBlocks))
{
}
void SetLocalHistogram(vtkm::Id index, const vtkm::cont::ArrayHandle<vtkm::Id>& bins)
{
this->LocalBlocks[static_cast<size_t>(index)] = bins;
}
void SetLocalHistogram(vtkm::Id index, const vtkm::cont::Field& field)
{
this->SetLocalHistogram(index, field.GetData().Cast<vtkm::cont::ArrayHandle<vtkm::Id>>());
}
vtkm::cont::ArrayHandle<vtkm::Id> ReduceAll(const vtkm::Id numBins) const
{
const vtkm::Id numLocalBlocks = static_cast<vtkm::Id>(this->LocalBlocks.size());
auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
if (comm.size() == 1 && numLocalBlocks <= 1)
{
// no reduction necessary.
return numLocalBlocks == 0 ? vtkm::cont::ArrayHandle<vtkm::Id>() : this->LocalBlocks[0];
}
// reduce local bins first.
vtkm::cont::ArrayHandle<vtkm::Id> local;
local.Allocate(numBins);
std::fill(vtkm::cont::ArrayPortalToIteratorBegin(local.GetPortalControl()),
vtkm::cont::ArrayPortalToIteratorEnd(local.GetPortalControl()),
static_cast<vtkm::Id>(0));
for (const auto& lbins : this->LocalBlocks)
{
vtkm::cont::Algorithm::Transform(local, lbins, local, vtkm::Add());
}
// now reduce across ranks using MPI.
// converting to std::vector
std::vector<vtkm::Id> send_buf(numBins);
std::copy(vtkm::cont::ArrayPortalToIteratorBegin(local.GetPortalConstControl()),
vtkm::cont::ArrayPortalToIteratorEnd(local.GetPortalConstControl()),
send_buf.begin());
std::vector<vtkm::Id> recv_buf(numBins);
MPI_Reduce(&send_buf[0],
&recv_buf[0],
static_cast<int>(numBins),
sizeof(vtkm::Id) == 4 ? MPI_INT : MPI_LONG,
MPI_SUM,
0,
comm);
if (comm.rank() == 0)
{
local.Allocate(numBins);
std::copy(recv_buf.begin(),
recv_buf.end(),
vtkm::cont::ArrayPortalToIteratorBegin(local.GetPortalControl()));
return local;
}
return vtkm::cont::ArrayHandle<vtkm::Id>();
}
};
} // namespace detail
//-----------------------------------------------------------------------------
inline VTKM_CONT HistogramMPI::HistogramMPI()
: NumberOfBins(10)
, BinDelta(0)
, ComputedRange()
, Range()
{
this->SetOutputFieldName("histogram");
}
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy, typename DeviceAdapter>
inline VTKM_CONT vtkm::cont::DataSet HistogramMPI::DoExecute(
const vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&,
const DeviceAdapter& device)
{
vtkm::cont::ArrayHandle<vtkm::Id> binArray;
T delta;
vtkm::worklet::FieldHistogram worklet;
if (this->ComputedRange.IsNonEmpty())
{
worklet.Run(field,
this->NumberOfBins,
static_cast<T>(this->ComputedRange.Min),
static_cast<T>(this->ComputedRange.Max),
delta,
binArray,
device);
}
else
{
worklet.Run(field, this->NumberOfBins, this->ComputedRange, delta, binArray, device);
}
this->BinDelta = static_cast<vtkm::Float64>(delta);
vtkm::cont::DataSet output;
vtkm::cont::Field rfield(
this->GetOutputFieldName(), vtkm::cont::Field::ASSOC_WHOLE_MESH, binArray);
output.AddField(rfield);
return output;
}
//-----------------------------------------------------------------------------
template <typename DerivedPolicy>
inline VTKM_CONT void HistogramMPI::PreExecute(const vtkm::cont::MultiBlock& input,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
if (this->Range.IsNonEmpty())
{
this->ComputedRange = this->Range;
}
else
{
auto handle = vtkm::cont::FieldRangeGlobalCompute(
input, this->GetActiveFieldName(), this->GetActiveFieldAssociation());
if (handle.GetNumberOfValues() != 1)
{
throw vtkm::cont::ErrorFilterExecution("expecting scalar field.");
}
this->ComputedRange = handle.GetPortalConstControl().Get(0);
}
}
//-----------------------------------------------------------------------------
template <typename DerivedPolicy>
inline VTKM_CONT void HistogramMPI::PostExecute(const vtkm::cont::MultiBlock&,
vtkm::cont::MultiBlock& result,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
// iterate and compute HistogramMPI for each local block.
detail::DistributedHistogram helper(result.GetNumberOfBlocks());
for (vtkm::Id cc = 0; cc < result.GetNumberOfBlocks(); ++cc)
{
auto& ablock = result.GetBlock(cc);
helper.SetLocalHistogram(cc, ablock.GetField(this->GetOutputFieldName()));
}
vtkm::cont::DataSet output;
vtkm::cont::Field rfield(this->GetOutputFieldName(),
vtkm::cont::Field::ASSOC_WHOLE_MESH,
helper.ReduceAll(this->NumberOfBins));
output.AddField(rfield);
result = vtkm::cont::MultiBlock(output);
}
} // namespace example

21
vtkm/cont/AssignerMultiBlock.cxx
View File

@ -38,15 +38,26 @@ namespace cont
VTKM_CONT
AssignerMultiBlock::AssignerMultiBlock(const vtkm::cont::MultiBlock& mb)
: AssignerMultiBlock(mb.GetNumberOfBlocks())
{
}
VTKM_CONT
AssignerMultiBlock::AssignerMultiBlock(vtkm::Id num_blocks)
: diy::Assigner(vtkm::cont::EnvironmentTracker::GetCommunicator().size(), 1)
, IScanBlockCounts()
{
auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
const auto num_blocks = mb.GetNumberOfBlocks();
vtkm::Id iscan;
diy::mpi::scan(comm, num_blocks, iscan, std::plus<vtkm::Id>());
diy::mpi::all_gather(comm, iscan, this->IScanBlockCounts);
if (comm.size() > 1)
{
vtkm::Id iscan;
diy::mpi::scan(comm, num_blocks, iscan, std::plus<vtkm::Id>());
diy::mpi::all_gather(comm, iscan, this->IScanBlockCounts);
}
else
{
this->IScanBlockCounts.push_back(num_blocks);
}
this->set_nblocks(static_cast<int>(this->IScanBlockCounts.back()));
}

3
vtkm/cont/AssignerMultiBlock.h
View File

@ -68,6 +68,9 @@ public:
VTKM_CONT
AssignerMultiBlock(const vtkm::cont::MultiBlock& mb);
VTKM_CONT
AssignerMultiBlock(vtkm::Id num_blocks);
///@{
/// diy::Assigner API implementation.
VTKM_CONT

11
vtkm/cont/EnvironmentTracker.cxx
View File

@ -41,6 +41,17 @@ void EnvironmentTracker::SetCommunicator(const diy::mpi::communicator& comm)
const diy::mpi::communicator& EnvironmentTracker::GetCommunicator()
{
#ifndef DIY_NO_MPI
int flag;
MPI_Initialized(&flag);
if (!flag)
{
int argc = 0;
char** argv = nullptr;
MPI_Init(&argc, &argv);
internal::GlobalCommuncator = diy::mpi::communicator(MPI_COMM_WORLD);
}
#endif
return vtkm::cont::internal::GlobalCommuncator;
}
} // namespace vtkm::cont

4
vtkm/cont/MultiBlock.h
View File

@ -44,10 +44,10 @@ public:
MultiBlock(const vtkm::cont::MultiBlock& src);
/// create a new MultiBlock with a DataSet vector "mblocks"
VTKM_CONT
MultiBlock(const std::vector<vtkm::cont::DataSet>& mblocks);
explicit MultiBlock(const std::vector<vtkm::cont::DataSet>& mblocks);
/// create a new MultiBlock with the capacity set to be "size"
VTKM_CONT
MultiBlock(vtkm::Id size);
explicit MultiBlock(vtkm::Id size);
VTKM_CONT
MultiBlock();

160
vtkm/filter/Filter.h
View File

@ -27,13 +27,153 @@
#include <vtkm/cont/RuntimeDeviceTracker.h>
#include <vtkm/filter/FieldSelection.h>
#include <vtkm/filter/FilterTraits.h>
#include <vtkm/filter/PolicyBase.h>
namespace vtkm
{
namespace filter
{
/// \brief base class for all filters.
///
/// This is the base class for all filters. To add a new filter, one can
/// subclass this (or any of the existing subclasses e.g. FilterField,
/// FilterDataSet, FilterDataSetWithField, etc. and implement relevant methods.
///
/// \section FilterUsage Usage
///
/// To execute a filter, one typically calls the `auto result =
/// filter.Execute(input)`. Typical usage is as follows:
///
/// \code{cpp}
///
/// // create the concrete subclass (e.g. MarchingCubes).
/// vtkm::filter::MarchingCubes marchingCubes;
///
/// // select fieds to map to the output, if different from default which is to
/// // map all input fields.
/// marchingCubes.SetFieldToPass({"var1", "var2"});
///
/// // execute the filter on vtkm::cont::DataSet.
/// vtkm::cont::DataSet dsInput = ...
/// auto outputDS = filter.Execute(dsInput);
///
/// // or, execute on a vtkm::cont::MultiBlock
/// vtkm::cont::MultiBlock mbInput = ...
/// auto outputMB = filter.Execute(mbInput);
/// \endcode
///
/// `Execute` methods take in the input dataset or multiblock to process and
/// return the result. The type of the result is same as the input type, thus
/// `Execute(DataSet&)` returns a DataSet while `Execute(MultiBlock&)` returns a
/// MultiBlock.
///
/// The implementation for `Execute(DataSet&)` is merely provided for
/// convenience. Internally, it creates MultiBlock with a single block for the
/// input and then forwards the call to `Execute(MultiBlock&)`. The method
/// returns the first block, if any, from the MultiBlock returned by the
/// forwarded call. If the MultiBlock returned has more than 1 block, then
/// `vtkm::cont::ErrorFilterExecution` will be thrown.
///
/// \section FilterSubclassing Subclassing
///
/// Typically, one subclasses one of the immediate subclasses of this class such as
/// FilterField, FilterDataSet, FilterDataSetWithField, etc. Those may impose
/// additional constraints on the methods to implement in the subclasses.
/// Here, we describes the things to consider when directly subclassing
/// vtkm::filter::Filter.
///
/// \subsection FilterPreExecutePostExecute PreExecute and PostExecute
///
/// Subclasses may provide implementations for either or both of the following
/// methods.
///
/// \code{cpp}
///
/// template <typename DerivedPolicy>
/// void PreExecute(const vtkm::cont::MultiBlock& input,
/// const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
///
/// template <typename DerivedPolicy>
/// void PostExecute(const vtkm::cont::MultiBlock& input, vtkm::cont::MultiBlock& output
/// const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
///
/// \endcode
///
/// As the name suggests, these are called and the beginning and before the end
/// of an `Filter::Execute` call. Most filters that don't need to handle
/// mutliblock datasets specially, e.g. clip, cut, iso-contour, need not worry
/// about these methods or provide any implementation. If, however, your filter
/// needs do to some initialization e.g. allocation buffers to accumulate
/// results, or finalization e.g. reduce results across all blocks, then these
/// methods provide convenient hooks for the same.
///
/// \subsection FilterPrepareForExecution PrepareForExecution
///
/// A concrete subclass of Filter must provide `PrepareForExecution`
/// implementation that provides the meat for the filter i.e. the implementation
/// for the filter's data processing logic. There are two signatures
/// available; which one to implement depends on the nature of the filter.
///
/// Let's consider simple filters that do not need to do anything special to
/// handle multiblock datasets e.g. clip, contour, etc. These are the filters
/// where executing the filter on a MultiBlock simply means executing the filter
/// on one block at a time and packing the output for each iteration info the
/// result MultiBlock. For such filters, one must implement the following
/// signature.
///
/// \code{cpp}
///
/// template <typename DerivedPolicy>
/// vtkm::cont::DataSet PrepareForExecution(
/// const vtkm::cont::DataSet& input,
/// const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
///
/// \endcode
///
/// The role of this method is to execute on the input dataset and generate the
/// result and return it. If there are any errors, the subclass must throw an
/// exception (e.g. `vtkm::cont::ErrorFilterExecution`).
///
/// In this case, the Filter superclass handles iterating over multiple blocks
/// in the input MultiBlock and calling `PrepareForExecution` iteratively.
///
/// The aforementioned approach is also suitable for filters that need special
/// handling for multiblock datasets which can be modelled as PreExecute and
/// PostExecute steps (e.g. `vtkm::filter::Histogram`).
///
/// For more complex filters, like streamlines, particle tracking, where the
/// processing of multiblock datasets cannot be modelled as a reduction of the
/// results, one can implement the following signature.
///
/// \code{cpp}
/// template <typename DerivedPolicy>
/// vtkm::cont::MultiBlock PrepareForExecution(
/// const vtkm::cont::MultiBlock& input,
/// const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
/// \endcode
///
/// The responsibility of this method is the same, except now the subclass is
/// given full control over the execution, including any mapping of fields to
/// output (described in next sub-section).
///
/// \subsection FilterMapFieldOntoOutput MapFieldOntoOutput
///
/// Subclasses may provide `MapFieldOntoOutput` method with the following
/// signature:
///
/// \code{cpp}
///
/// template <typename DerivedPolicy>
/// VTKM_CONT bool MapFieldOntoOutput(vtkm::cont::DataSet& result,
/// const vtkm::cont::Field& field,
/// const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
///
/// \endcode
///
/// When present, this method will be called after each block execution to map
/// an input field from the corresponding input block to the output block.
///
template <typename Derived>
class Filter
{
@ -102,8 +242,7 @@ public:
/// Executes the filter on the input and produces a result dataset.
///
/// On success, this the dataset produced. On error, vtkm::cont::ErrorExecution will be thrown.
VTKM_CONT
vtkm::cont::DataSet Execute(const vtkm::cont::DataSet& input);
VTKM_CONT vtkm::cont::DataSet Execute(const vtkm::cont::DataSet& input);
template <typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet Execute(const vtkm::cont::DataSet& input,
@ -111,15 +250,24 @@ public:
//@}
//@{
/// MultiBlock variants of execute.
VTKM_CONT
vtkm::cont::MultiBlock Execute(const vtkm::cont::MultiBlock& input);
/// Executes the filter on the input MultiBlock and produces a result MultiBlock.
///
/// On success, this the dataset produced. On error, vtkm::cont::ErrorExecution will be thrown.
VTKM_CONT vtkm::cont::MultiBlock Execute(const vtkm::cont::MultiBlock& input);
template <typename DerivedPolicy>
VTKM_CONT vtkm::cont::MultiBlock Execute(const vtkm::cont::MultiBlock& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
//@}
/// Map fields from input dataset to output.
/// This is not intended for external use. Subclasses of Filter, however, may
/// use this method to map fields.
template <typename DerivedPolicy>
VTKM_CONT void MapFieldsToPass(const vtkm::cont::DataSet& input,
vtkm::cont::DataSet& output,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
vtkm::cont::RuntimeDeviceTracker Tracker;
vtkm::filter::FieldSelection FieldsToPass;

273
vtkm/filter/Filter.hxx
View File

@ -19,15 +19,12 @@
//============================================================================
#include <vtkm/filter/FieldMetadata.h>
#include <vtkm/filter/FilterTraits.h>
#include <vtkm/filter/PolicyDefault.h>
#include <vtkm/filter/internal/ResolveFieldTypeAndExecute.h>
#include <vtkm/filter/internal/ResolveFieldTypeAndMap.h>
#include <vtkm/cont/Error.h>
#include <vtkm/cont/ErrorBadAllocation.h>
#include <vtkm/cont/ErrorExecution.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/cont/Field.h>
#include <vtkm/cont/cuda/DeviceAdapterCuda.h>
@ -38,6 +35,205 @@ namespace vtkm
namespace filter
{
namespace internal
{
template <typename T, typename InputType, typename DerivedPolicy>
struct SupportsPreExecute
{
template <typename U,
typename S = decltype(std::declval<U>().PreExecute(
std::declval<InputType>(),
std::declval<vtkm::filter::PolicyBase<DerivedPolicy>>()))>
static std::true_type has(int);
template <typename U>
static std::false_type has(...);
using type = decltype(has<T>(0));
};
template <typename T, typename InputType, typename DerivedPolicy>
struct SupportsPostExecute
{
template <typename U,
typename S = decltype(std::declval<U>().PostExecute(
std::declval<InputType>(),
std::declval<InputType&>(),
std::declval<vtkm::filter::PolicyBase<DerivedPolicy>>()))>
static std::true_type has(int);
template <typename U>
static std::false_type has(...);
using type = decltype(has<T>(0));
};
template <typename T, typename InputType, typename DerivedPolicy>
struct SupportsPrepareForExecution
{
template <typename U,
typename S = decltype(std::declval<U>().PrepareForExecution(
std::declval<InputType>(),
std::declval<vtkm::filter::PolicyBase<DerivedPolicy>>()))>
static std::true_type has(int);
template <typename U>
static std::false_type has(...);
using type = decltype(has<T>(0));
};
template <typename T, typename DerivedPolicy>
struct SupportsMapFieldOntoOutput
{
template <typename U,
typename S = decltype(std::declval<U>().MapFieldOntoOutput(
std::declval<vtkm::cont::DataSet&>(),
std::declval<vtkm::cont::Field>(),
std::declval<vtkm::filter::PolicyBase<DerivedPolicy>>()))>
static std::true_type has(int);
template <typename U>
static std::false_type has(...);
using type = decltype(has<T>(0));
};
//--------------------------------------------------------------------------------
template <typename Derived, typename... Args>
void CallPreExecuteInternal(std::true_type, Derived* self, Args&&... args)
{
return self->PreExecute(std::forward<Args>(args)...);
}
//--------------------------------------------------------------------------------
template <typename Derived, typename... Args>
void CallPreExecuteInternal(std::false_type, Derived*, Args&&...)
{
}
//--------------------------------------------------------------------------------
template <typename Derived, typename InputType, typename DerivedPolicy>
void CallPreExecute(Derived* self,
const InputType& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
using call_supported_t = typename SupportsPreExecute<Derived, InputType, DerivedPolicy>::type;
CallPreExecuteInternal(call_supported_t(), self, input, policy);
}
//--------------------------------------------------------------------------------
template <typename Derived, typename DerivedPolicy>
void CallMapFieldOntoOutputInternal(std::true_type,
Derived* self,
const vtkm::cont::DataSet& input,
vtkm::cont::DataSet& output,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
for (vtkm::IdComponent cc = 0; cc < input.GetNumberOfFields(); ++cc)
{
auto field = input.GetField(cc);
if (self->GetFieldsToPass().IsFieldSelected(field))
{
self->MapFieldOntoOutput(output, field, policy);
}
}
}
//--------------------------------------------------------------------------------
template <typename Derived, typename DerivedPolicy>
void CallMapFieldOntoOutputInternal(std::false_type,
Derived*,
const vtkm::cont::DataSet&,
vtkm::cont::DataSet&,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
// do nothing.
}
//--------------------------------------------------------------------------------
template <typename Derived, typename DerivedPolicy>
void CallMapFieldOntoOutput(Derived* self,
const vtkm::cont::DataSet& input,
vtkm::cont::DataSet& output,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
using call_supported_t = typename SupportsMapFieldOntoOutput<Derived, DerivedPolicy>::type;
CallMapFieldOntoOutputInternal(call_supported_t(), self, input, output, policy);
}
//--------------------------------------------------------------------------------
// forward declare.
template <typename Derived, typename InputType, typename DerivedPolicy>
InputType CallPrepareForExecution(Derived* self,
const InputType& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
//--------------------------------------------------------------------------------
template <typename Derived, typename InputType, typename DerivedPolicy>
InputType CallPrepareForExecutionInternal(std::true_type,
Derived* self,
const InputType& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
return self->PrepareForExecution(input, policy);
}
//--------------------------------------------------------------------------------
// specialization for MultiBlock input when `PrepareForExecution` is not provided
// by the subclass. we iterate over blocks and execute for each block
// individually.
template <typename Derived, typename DerivedPolicy>
vtkm::cont::MultiBlock CallPrepareForExecutionInternal(
std::false_type,
Derived* self,
const vtkm::cont::MultiBlock& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
vtkm::cont::MultiBlock output;
for (const auto& inBlock : input)
{
vtkm::cont::DataSet outBlock = CallPrepareForExecution(self, inBlock, policy);
CallMapFieldOntoOutput(self, inBlock, outBlock, policy);
output.AddBlock(outBlock);
}
return output;
}
//--------------------------------------------------------------------------------
template <typename Derived, typename InputType, typename DerivedPolicy>
InputType CallPrepareForExecution(Derived* self,
const InputType& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
using call_supported_t =
typename SupportsPrepareForExecution<Derived, InputType, DerivedPolicy>::type;
return CallPrepareForExecutionInternal(call_supported_t(), self, input, policy);
}
//--------------------------------------------------------------------------------
template <typename Derived, typename InputType, typename DerivedPolicy>
void CallPostExecuteInternal(std::true_type,
Derived* self,
const InputType& input,
InputType& output,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
self->PostExecute(input, output, policy);
}
//--------------------------------------------------------------------------------
template <typename Derived, typename... Args>
void CallPostExecuteInternal(std::false_type, Derived*, Args&&...)
{
}
//--------------------------------------------------------------------------------
template <typename Derived, typename InputType, typename DerivedPolicy>
void CallPostExecute(Derived* self,
const InputType& input,
InputType& output,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
using call_supported_t = typename SupportsPostExecute<Derived, InputType, DerivedPolicy>::type;
CallPostExecuteInternal(call_supported_t(), self, input, output, policy);
}
}
//----------------------------------------------------------------------------
template <typename Derived>
inline VTKM_CONT Filter<Derived>::Filter()
@ -59,6 +255,15 @@ inline VTKM_CONT vtkm::cont::DataSet Filter<Derived>::Execute(const vtkm::cont::
return this->Execute(input, vtkm::filter::PolicyDefault());
}
//----------------------------------------------------------------------------
template <typename Derived>
inline VTKM_CONT vtkm::cont::MultiBlock Filter<Derived>::Execute(
const vtkm::cont::MultiBlock& input)
{
return this->Execute(input, vtkm::filter::PolicyDefault());
}
//----------------------------------------------------------------------------
template <typename Derived>
template <typename DerivedPolicy>
@ -67,29 +272,12 @@ inline VTKM_CONT vtkm::cont::DataSet Filter<Derived>::Execute(
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
Derived* self = static_cast<Derived*>(this);
vtkm::cont::DataSet output = self->PrepareForExecution(input, policy);
// the above will throw `vtkm::cont::ErrorExecution` if the execution fails
// so we don't need to do anything special to validate the output.
for (vtkm::IdComponent cc = 0; cc < input.GetNumberOfFields(); ++cc)
vtkm::cont::MultiBlock output = self->Execute(vtkm::cont::MultiBlock(input), policy);
if (output.GetNumberOfBlocks() > 1)
{
auto field = input.GetField(cc);
if (this->GetFieldsToPass().IsFieldSelected(field))
{
// todo: should we thrown ErrorExecution if mapping failed?
self->MapFieldOntoOutput(output, field, policy);
}
throw vtkm::cont::ErrorFilterExecution("Expecting at most 1 block.");
}
return output;
}
//----------------------------------------------------------------------------
template <typename Derived>
inline VTKM_CONT vtkm::cont::MultiBlock Filter<Derived>::Execute(
const vtkm::cont::MultiBlock& input)
{
return this->Execute(input, vtkm::filter::PolicyDefault());
return output.GetNumberOfBlocks() == 1 ? output.GetBlock(0) : vtkm::cont::DataSet();
}
//----------------------------------------------------------------------------
@ -99,13 +287,36 @@ inline VTKM_CONT vtkm::cont::MultiBlock Filter<Derived>::Execute(
const vtkm::cont::MultiBlock& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
vtkm::cont::MultiBlock output;
for (auto& inDataSet : input)
{
vtkm::cont::DataSet outDataSet = this->Execute(inDataSet, policy);
output.AddBlock(outDataSet);
}
Derived* self = static_cast<Derived*>(this);
// Call `void Derived::PreExecute<DerivedPolicy>(input, policy)`, if defined.
internal::CallPreExecute(self, input, policy);
// Call `PrepareForExecution` (which should probably be renamed at some point)
vtkm::cont::MultiBlock output = internal::CallPrepareForExecution(self, input, policy);
// Call `Derived::PostExecute<DerivedPolicy>(input, output, policy)` if defined.
internal::CallPostExecute(self, input, output, policy);
return output;
}
//----------------------------------------------------------------------------
template <typename Derived>
template <typename DerivedPolicy>
inline VTKM_CONT void Filter<Derived>::MapFieldsToPass(
const vtkm::cont::DataSet& input,
vtkm::cont::DataSet& output,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy)
{
Derived* self = static_cast<Derived*>(this);
for (vtkm::IdComponent cc = 0; cc < input.GetNumberOfFields(); ++cc)
{
auto field = input.GetField(cc);
if (this->GetFieldsToPass().IsFieldSelected(field))
{
internal::CallMapFieldOntoOutput(self, output, field, policy);
}
}
}
}
}

3
vtkm/filter/FilterDataSet.h
View File

@ -58,7 +58,6 @@ public:
VTKM_CONT
vtkm::Id GetActiveCoordinateSystemIndex() const { return this->CoordinateSystemIndex; }
private:
/// These are provided to satisfy the Filter API requirements.
//From the field we can extract the association component
@ -72,12 +71,12 @@ private:
const vtkm::cont::Field& field,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet PrepareForExecution(
const vtkm::cont::DataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
vtkm::Id CellSetIndex;
vtkm::Id CoordinateSystemIndex;

2
vtkm/filter/FilterDataSetWithField.h
View File

@ -85,7 +85,6 @@ public:
bool GetUseCoordinateSystemAsField() const { return this->UseCoordinateSystemAsField; }
//@}
private:
//From the field we can extract the association component
// ASSOC_ANY -> unable to map
// ASSOC_WHOLE_MESH -> (I think this is points)
@ -102,6 +101,7 @@ private:
const vtkm::cont::DataSet& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
template <typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet PrepareForExecution(
const vtkm::cont::DataSet& input,

11
vtkm/filter/FilterField.h
View File

@ -89,16 +89,7 @@ public:
vtkm::Id GetActiveCoordinateSystemIndex() const { return this->CoordinateSystemIndex; }
//@}
private:
/// These are provided to satisfy the Filter API requirements.
/// They simply pass the field to the result dataset as there's no mapping
/// needed for FilterField.
template <typename DerivedPolicy>
VTKM_CONT bool MapFieldOntoOutput(vtkm::cont::DataSet& result,
const vtkm::cont::Field& field,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename DerivedPolicy>
VTKM_CONT vtkm::cont::DataSet PrepareForExecution(
const vtkm::cont::DataSet& input,
@ -116,7 +107,7 @@ private:
const vtkm::cont::CoordinateSystem& field,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
private:
std::string OutputFieldName;
vtkm::Id CoordinateSystemIndex;
std::string ActiveFieldName;

12
vtkm/filter/FilterField.hxx
View File

@ -123,17 +123,5 @@ inline VTKM_CONT vtkm::cont::DataSet FilterField<Derived>::PrepareForExecution(
return result;
}
//-----------------------------------------------------------------------------
template <typename Derived>
template <typename DerivedPolicy>
inline VTKM_CONT bool FilterField<Derived>::MapFieldOntoOutput(
vtkm::cont::DataSet& result,
const vtkm::cont::Field& field,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
result.AddField(field);
return true;
}
}
}

1
vtkm/filter/FilterTraits.h
View File

@ -27,6 +27,7 @@ namespace vtkm
{
namespace filter
{
template <typename Filter>
class FilterTraits
{

40
vtkm/filter/Histogram.h
View File

@ -27,9 +27,10 @@ namespace vtkm
{
namespace filter
{
/// \brief Construct the histogram of a given Field
///
/// Construct a histogram with a default of 10 bins
/// Construct a histogram with a default of 10 bins.
///
class Histogram : public vtkm::filter::FilterField<Histogram>
{
@ -44,14 +45,26 @@ public:
VTKM_CONT
vtkm::Id GetNumberOfBins() const { return this->NumberOfBins; }
//Returns the bin delta of the last computed field, be it from DoExecute
//or from MapField
//@{
/// Get/Set the range to use to generate the histogram. If range is set to
/// empty, the field's global range (computed using `vtkm::cont::FieldRangeGlobalCompute`)
/// will be used.
VTKM_CONT
void SetRange(const vtkm::Range& range) { this->Range = range; }
VTKM_CONT
const vtkm::Range& GetRange() const { return this->Range; }
//@}
/// Returns the bin delta of the last computed field.
VTKM_CONT
vtkm::Float64 GetBinDelta() const { return this->BinDelta; }
//Returns the the min and max values for that last computed field
/// Returns the range used for most recent execute. If `SetRange` is used to
/// specify and non-empty range, then this will be same as the range after
/// the `Execute` call.
VTKM_CONT
vtkm::Range GetDataRange() const { return this->DataRange; }
vtkm::Range GetComputedRange() const { return this->ComputedRange; }
template <typename T, typename StorageType, typename DerivedPolicy, typename DeviceAdapter>
VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet& input,
@ -60,10 +73,25 @@ public:
const vtkm::filter::PolicyBase<DerivedPolicy>& policy,
const DeviceAdapter& tag);
//@{
/// when operating on vtkm::cont::MultiBlock, we
/// want to do processing across ranks as well. Just adding pre/post handles
/// for the same does the trick.
template <typename DerivedPolicy>
VTKM_CONT void PreExecute(const vtkm::cont::MultiBlock& input,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy);
template <typename DerivedPolicy>
VTKM_CONT void PostExecute(const vtkm::cont::MultiBlock& input,
vtkm::cont::MultiBlock& output,
const vtkm::filter::PolicyBase<DerivedPolicy>&);
//@}
private:
vtkm::Id NumberOfBins;
vtkm::Float64 BinDelta;
vtkm::Range DataRange;
vtkm::Range ComputedRange;
vtkm::Range Range;
};
template <>

231
vtkm/filter/Histogram.hxx
View File

@ -18,20 +18,175 @@
// this software.
//============================================================================
#include <vtkm/filter/internal/CreateResult.h>
#include <vtkm/worklet/DispatcherMapField.h>
#include <vtkm/worklet/FieldHistogram.h>
#include <vtkm/cont/Algorithm.h>
#include <vtkm/cont/ArrayCopy.h>
#include <vtkm/cont/AssignerMultiBlock.h>
#include <vtkm/cont/EnvironmentTracker.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/cont/FieldRangeGlobalCompute.h>
#include <vtkm/cont/diy/Serialization.h>
#include <vtkm/filter/internal/CreateResult.h>
// clang-format off
VTKM_THIRDPARTY_PRE_INCLUDE
#include <vtkm/thirdparty/diy/Configure.h>
#include VTKM_DIY(diy/decomposition.hpp)
#include VTKM_DIY(diy/master.hpp)
#include VTKM_DIY(diy/partners/broadcast.hpp)
#include VTKM_DIY(diy/partners/swap.hpp)
#include VTKM_DIY(diy/reduce.hpp)
VTKM_THIRDPARTY_POST_INCLUDE
// clang-format on
namespace vtkm
{
namespace filter
{
namespace detail
{
class DistributedHistogram
{
class Reducer
{
public:
void operator()(vtkm::cont::ArrayHandle<vtkm::Id>* result,
const diy::ReduceProxy& srp,
const diy::RegularMergePartners&) const
{
const auto selfid = srp.gid();
// 1. dequeue.
std::vector<int> incoming;
srp.incoming(incoming);
for (const int gid : incoming)
{
if (gid != selfid)
{
vtkm::cont::ArrayHandle<vtkm::Id> in;
srp.dequeue(gid, in);
if (result->GetNumberOfValues() == 0)
{
*result = in;
}
else
{
vtkm::cont::Algorithm::Transform(*result, in, *result, vtkm::Add());
}
}
}
// 2. enqueue
for (int cc = 0; cc < srp.out_link().size(); ++cc)
{
auto target = srp.out_link().target(cc);
if (target.gid != selfid)
{
srp.enqueue(target, *result);
}
}
}
};
std::vector<vtkm::cont::ArrayHandle<vtkm::Id>> LocalBlocks;
public:
DistributedHistogram(vtkm::Id numLocalBlocks)
: LocalBlocks(static_cast<size_t>(numLocalBlocks))
{
}
void SetLocalHistogram(vtkm::Id index, const vtkm::cont::ArrayHandle<vtkm::Id>& bins)
{
this->LocalBlocks[static_cast<size_t>(index)] = bins;
}
void SetLocalHistogram(vtkm::Id index, const vtkm::cont::Field& field)
{
this->SetLocalHistogram(index, field.GetData().Cast<vtkm::cont::ArrayHandle<vtkm::Id>>());
}
vtkm::cont::ArrayHandle<vtkm::Id> ReduceAll() const
{
using ArrayType = vtkm::cont::ArrayHandle<vtkm::Id>;
const vtkm::Id numLocalBlocks = static_cast<vtkm::Id>(this->LocalBlocks.size());
auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
if (comm.size() == 1 && numLocalBlocks <= 1)
{
// no reduction necessary.
return numLocalBlocks == 0 ? ArrayType() : this->LocalBlocks[0];
}
diy::Master master(
comm,
/*threads*/ 1,
/*limit*/ -1,
[]() -> void* { return new vtkm::cont::ArrayHandle<vtkm::Id>(); },
[](void* ptr) { delete static_cast<vtkm::cont::ArrayHandle<vtkm::Id>*>(ptr); });
vtkm::cont::AssignerMultiBlock assigner(numLocalBlocks);
diy::RegularDecomposer<diy::DiscreteBounds> decomposer(
/*dims*/ 1, diy::interval(0, assigner.nblocks() - 1), assigner.nblocks());
decomposer.decompose(comm.rank(), assigner, master);
assert(static_cast<vtkm::Id>(master.size()) == numLocalBlocks);
for (vtkm::Id cc = 0; cc < numLocalBlocks; ++cc)
{
*master.block<ArrayType>(static_cast<int>(cc)) = this->LocalBlocks[static_cast<size_t>(cc)];
}
diy::RegularMergePartners partners(decomposer, /*k=*/2);
// reduce to block-0.
diy::reduce(master, assigner, partners, Reducer());
ArrayType result;
if (master.local(0))
{
result = *master.block<ArrayType>(master.lid(0));
}
this->Broadcast(result);
return result;
}
private:
void Broadcast(vtkm::cont::ArrayHandle<vtkm::Id>& data) const
{
// broadcast to all ranks (and not blocks).
auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
if (comm.size() > 1)
{
using ArrayType = vtkm::cont::ArrayHandle<vtkm::Id>;
diy::Master master(
comm,
/*threads*/ 1,
/*limit*/ -1,
[]() -> void* { return new vtkm::cont::ArrayHandle<vtkm::Id>(); },
[](void* ptr) { delete static_cast<vtkm::cont::ArrayHandle<vtkm::Id>*>(ptr); });
diy::ContiguousAssigner assigner(comm.size(), comm.size());
diy::RegularDecomposer<diy::DiscreteBounds> decomposer(
1, diy::interval(0, comm.size() - 1), comm.size());
decomposer.decompose(comm.rank(), assigner, master);
assert(master.size() == 1); // number of local blocks should be 1 per rank.
*master.block<ArrayType>(0) = data;
diy::RegularBroadcastPartners partners(decomposer, /*k=*/2);
diy::reduce(master, assigner, partners, Reducer());
data = *master.block<ArrayType>(0);
}
}
};
} // namespace detail
//-----------------------------------------------------------------------------
inline VTKM_CONT Histogram::Histogram()
: NumberOfBins(10)
, BinDelta(0)
, DataRange()
, ComputedRange()
, Range()
{
this->SetOutputFieldName("histogram");
}
@ -39,9 +194,9 @@ inline VTKM_CONT Histogram::Histogram()
//-----------------------------------------------------------------------------
template <typename T, typename StorageType, typename DerivedPolicy, typename DeviceAdapter>
inline VTKM_CONT vtkm::cont::DataSet Histogram::DoExecute(
const vtkm::cont::DataSet& inDataSet,
const vtkm::cont::DataSet&,
const vtkm::cont::ArrayHandle<T, StorageType>& field,
const vtkm::filter::FieldMetadata& fieldMetadata,
const vtkm::filter::FieldMetadata&,
const vtkm::filter::PolicyBase<DerivedPolicy>&,
const DeviceAdapter& device)
{
@ -49,14 +204,70 @@ inline VTKM_CONT vtkm::cont::DataSet Histogram::DoExecute(
T delta;
vtkm::worklet::FieldHistogram worklet;
worklet.Run(field, this->NumberOfBins, this->DataRange, delta, binArray, device);
if (this->ComputedRange.IsNonEmpty())
{
worklet.Run(field,
this->NumberOfBins,
static_cast<T>(this->ComputedRange.Min),
static_cast<T>(this->ComputedRange.Max),
delta,
binArray,
device);
}
else
{
worklet.Run(field, this->NumberOfBins, this->ComputedRange, delta, binArray, device);
}
this->BinDelta = static_cast<vtkm::Float64>(delta);
return internal::CreateResult(inDataSet,
binArray,
this->GetOutputFieldName(),
fieldMetadata.GetAssociation(),
fieldMetadata.GetCellSetName());
vtkm::cont::DataSet output;
vtkm::cont::Field rfield(
this->GetOutputFieldName(), vtkm::cont::Field::ASSOC_WHOLE_MESH, binArray);
output.AddField(rfield);
return output;
}
//-----------------------------------------------------------------------------
template <typename DerivedPolicy>
inline VTKM_CONT void Histogram::PreExecute(const vtkm::cont::MultiBlock& input,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
if (this->Range.IsNonEmpty())
{
this->ComputedRange = this->Range;
}
else
{
auto handle = vtkm::cont::FieldRangeGlobalCompute(
input, this->GetActiveFieldName(), this->GetActiveFieldAssociation());
if (handle.GetNumberOfValues() != 1)
{
throw vtkm::cont::ErrorFilterExecution("expecting scalar field.");
}
this->ComputedRange = handle.GetPortalConstControl().Get(0);
}
}
//-----------------------------------------------------------------------------
template <typename DerivedPolicy>
inline VTKM_CONT void Histogram::PostExecute(const vtkm::cont::MultiBlock&,
vtkm::cont::MultiBlock& result,
const vtkm::filter::PolicyBase<DerivedPolicy>&)
{
// iterate and compute histogram for each local block.
detail::DistributedHistogram helper(result.GetNumberOfBlocks());
for (vtkm::Id cc = 0; cc < result.GetNumberOfBlocks(); ++cc)
{
auto& ablock = result.GetBlock(cc);
helper.SetLocalHistogram(cc, ablock.GetField(this->GetOutputFieldName()));
}
vtkm::cont::DataSet output;
vtkm::cont::Field rfield(
this->GetOutputFieldName(), vtkm::cont::Field::ASSOC_WHOLE_MESH, helper.ReduceAll());
output.AddField(rfield);
result = vtkm::cont::MultiBlock(output);
}
}
} // namespace vtkm::filter

1
vtkm/filter/testing/CMakeLists.txt
View File

@ -40,6 +40,7 @@ set(unit_tests
UnitTestMaskFilter.cxx
UnitTestMaskPointsFilter.cxx
UnitTestMultiBlockFilters.cxx
UnitTestMultiBlockHistogramFilter.cxx
UnitTestNDEntropyFilter.cxx
UnitTestNDHistogramFilter.cxx
UnitTestPointAverageFilter.cxx

8
vtkm/filter/testing/UnitTestHistogramFilter.cxx
View File

@ -311,28 +311,28 @@ void TestHistogram()
histogram.SetActiveField("p_poisson");
auto result = histogram.Execute(ds);
delta = histogram.GetBinDelta();
range = histogram.GetDataRange();
range = histogram.GetComputedRange();
VerifyHistogram(result, histogram.GetNumberOfBins(), range, delta);
histogram.SetNumberOfBins(100);
histogram.SetActiveField("p_normal");
result = histogram.Execute(ds);
delta = histogram.GetBinDelta();
range = histogram.GetDataRange();
range = histogram.GetComputedRange();
VerifyHistogram(result, histogram.GetNumberOfBins(), range, delta, false);
histogram.SetNumberOfBins(1);
histogram.SetActiveField("p_chiSquare");
result = histogram.Execute(ds);
delta = histogram.GetBinDelta();
range = histogram.GetDataRange();
range = histogram.GetComputedRange();
VerifyHistogram(result, histogram.GetNumberOfBins(), range, delta);
histogram.SetNumberOfBins(1000000);
histogram.SetActiveField("p_uniform");
result = histogram.Execute(ds);
delta = histogram.GetBinDelta();
range = histogram.GetDataRange();
range = histogram.GetComputedRange();
VerifyHistogram(result, histogram.GetNumberOfBins(), range, delta, false);
} // TestFieldHistogram

9
vtkm/filter/testing/UnitTestMultiBlockFilters.cxx
View File

@ -33,7 +33,6 @@
#include <vtkm/cont/testing/MakeTestDataSet.h>
#include <vtkm/cont/testing/Testing.h>
#include <vtkm/filter/CellAverage.h>
#include <vtkm/filter/Histogram.h>
template <typename T>
@ -85,7 +84,7 @@ vtkm::cont::MultiBlock MultiBlockBuilder(std::size_t BlockNum, std::string Field
}
template <typename D>
void Result_Verify(const vtkm::cont::MultiBlock& Result,
D Filter,
D& Filter,
const vtkm::cont::MultiBlock& Blocks,
std::string FieldName)
{
@ -124,12 +123,6 @@ void TestMultiBlockFilters()
vtkm::cont::MultiBlock result;
vtkm::cont::MultiBlock Blocks;
Blocks = MultiBlockBuilder<vtkm::Float64>(BlockNum, "cellvar");
vtkm::filter::Histogram histogram;
histogram.SetActiveField("cellvar");
result = histogram.Execute(Blocks);
Result_Verify(result, histogram, Blocks, std::string("cellvar"));
Blocks = MultiBlockBuilder<vtkm::Id>(BlockNum, "pointvar");
vtkm::filter::CellAverage cellAverage;
cellAverage.SetOutputFieldName("average");

136
vtkm/filter/testing/UnitTestMultiBlockHistogramFilter.cxx Normal file
View File

@ -0,0 +1,136 @@
//============================================================================
// 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.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// 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/filter/Histogram.h>
#include <vtkm/cont/DataSet.h>
#include <vtkm/cont/MultiBlock.h>
#include <vtkm/cont/testing/Testing.h>
#include <algorithm>
#include <numeric>
#include <random>
#include <utility>
namespace
{
static unsigned int uid = 1;
template <typename T>
vtkm::cont::ArrayHandle<T> CreateArrayHandle(T min, T max, vtkm::Id numVals)
{
std::mt19937 gen(uid++);
std::uniform_real_distribution<double> dis(static_cast<double>(min), static_cast<double>(max));
vtkm::cont::ArrayHandle<T> handle;
handle.Allocate(numVals);
std::generate(vtkm::cont::ArrayPortalToIteratorBegin(handle.GetPortalControl()),
vtkm::cont::ArrayPortalToIteratorEnd(handle.GetPortalControl()),
[&]() { return static_cast<T>(dis(gen)); });
return handle;
}
template <typename T, int size>
vtkm::cont::ArrayHandle<vtkm::Vec<T, size>> CreateArrayHandle(const vtkm::Vec<T, size>& min,
const vtkm::Vec<T, size>& max,
vtkm::Id numVals)
{
std::mt19937 gen(uid++);
std::uniform_real_distribution<double> dis[size];
for (int cc = 0; cc < size; ++cc)
{
dis[cc] = std::uniform_real_distribution<double>(static_cast<double>(min[cc]),
static_cast<double>(max[cc]));
}
vtkm::cont::ArrayHandle<T> handle;
handle.Allocate(numVals);
std::generate(vtkm::cont::ArrayPortalToIteratorBegin(handle.GetPortalControl()),
vtkm::cont::ArrayPortalToIteratorEnd(handle.GetPortalControl()),
[&]() {
vtkm::Vec<T, size> val;
for (int cc = 0; cc < size; ++cc)
{
val[cc] = static_cast<T>(dis[cc](gen));
}
return val;
});
return handle;
}
template <typename T>
void AddField(vtkm::cont::DataSet& dataset,
const T& min,
const T& max,
vtkm::Id numVals,
const std::string& name,
vtkm::cont::Field::AssociationEnum assoc = vtkm::cont::Field::ASSOC_POINTS)
{
auto ah = CreateArrayHandle(min, max, numVals);
dataset.AddField(vtkm::cont::Field(name, assoc, ah));
}
}
static void TestMultiBlockHistogram()
{
// init random seed.
std::srand(100);
vtkm::cont::MultiBlock mb;
vtkm::cont::DataSet block0;
AddField<double>(block0, 0.0, 100.0, 1024, "double");
mb.AddBlock(block0);
vtkm::cont::DataSet block1;
AddField<int>(block1, 100, 1000, 1024, "double");
mb.AddBlock(block1);
vtkm::cont::DataSet block2;
AddField<double>(block2, 100.0, 500.0, 1024, "double");
mb.AddBlock(block2);
vtkm::filter::Histogram histogram;
histogram.SetActiveField("double");
auto result = histogram.Execute(mb);
VTKM_TEST_ASSERT(result.GetNumberOfBlocks() == 1, "Expecting 1 block.");
auto bins =
result.GetBlock(0).GetField("histogram").GetData().Cast<vtkm::cont::ArrayHandle<vtkm::Id>>();
VTKM_TEST_ASSERT(bins.GetNumberOfValues() == 10, "Expecting 10 bins.");
auto count = std::accumulate(vtkm::cont::ArrayPortalToIteratorBegin(bins.GetPortalConstControl()),
vtkm::cont::ArrayPortalToIteratorEnd(bins.GetPortalConstControl()),
vtkm::Id(0),
vtkm::Add());
VTKM_TEST_ASSERT(count == 1024 * 3, "Expecting 3072 values");
std::cout << "Values [" << count << "] =";
for (int cc = 0; cc < 10; ++cc)
{
std::cout << " " << bins.GetPortalConstControl().Get(cc);
}
std::cout << std::endl;
};
int UnitTestMultiBlockHistogramFilter(int, char* [])
{
return vtkm::cont::testing::Testing::Run(TestMultiBlockHistogram);
}

32
vtkm/worklet/FieldHistogram.h
View File

@ -127,22 +127,37 @@ public:
vtkm::Range& rangeOfValues,
FieldType& binDelta,
vtkm::cont::ArrayHandle<vtkm::Id>& binArray,
DeviceAdapter vtkmNotUsed(device))
DeviceAdapter device)
{
using DeviceAlgorithms = typename vtkm::cont::DeviceAdapterAlgorithm<DeviceAdapter>;
//todo: need to have a signature that can use an input range so we can
//leverage fields that have already computed there range
const vtkm::Id numberOfValues = fieldArray.GetNumberOfValues();
const vtkm::Vec<FieldType, 2> initValue(fieldArray.GetPortalConstControl().Get(0));
vtkm::Vec<FieldType, 2> result =
DeviceAlgorithms::Reduce(fieldArray, initValue, vtkm::MinAndMax<FieldType>());
const FieldType& fieldMinValue = result[0];
const FieldType& fieldMaxValue = result[1];
this->Run(fieldArray, numberOfBins, result[0], result[1], binDelta, binArray, device);
//update the users data
rangeOfValues = vtkm::Range(result[0], result[1]);
}
// Execute the histogram binning filter given data and number of bins, min,
// max values.
// Returns:
// number of values in each bin
template <typename FieldType, typename Storage, typename DeviceAdapter>
void Run(vtkm::cont::ArrayHandle<FieldType, Storage> fieldArray,
vtkm::Id numberOfBins,
FieldType fieldMinValue,
FieldType fieldMaxValue,
FieldType& binDelta,
vtkm::cont::ArrayHandle<vtkm::Id>& binArray,
DeviceAdapter vtkmNotUsed(device))
{
using DeviceAlgorithms = typename vtkm::cont::DeviceAdapterAlgorithm<DeviceAdapter>;
const vtkm::Id numberOfValues = fieldArray.GetNumberOfValues();
const FieldType fieldDelta = compute_delta(fieldMinValue, fieldMaxValue, numberOfBins);
// Worklet fills in the bin belonging to each value
@ -168,7 +183,6 @@ public:
dispatcher.Invoke(binCounter, totalCount, binArray);
//update the users data
rangeOfValues = vtkm::Range(fieldMinValue, fieldMaxValue);
binDelta = fieldDelta;
}
};