adding support and tests for cell-to-point topology worklet.

This commit is contained in:
Jeremy Meredith 2015-08-28 13:51:42 -04:00
parent f3e3ecb68e
commit 2987142fa6
6 changed files with 250 additions and 12 deletions

@ -26,6 +26,9 @@
#include <vtkm/cont/internal/ConnectivityExplicitInternals.h>
#include <vtkm/exec/ConnectivityExplicit.h>
#include <map>
#include <utility>
namespace vtkm {
namespace cont {
@ -59,17 +62,20 @@ public:
typedef vtkm::Id SchedulingRangeType;
VTKM_CONT_EXPORT
CellSetExplicit(const std::string &name = std::string(),
CellSetExplicit(vtkm::Id numpoints = 0,
const std::string &name = std::string(),
vtkm::IdComponent dimensionality = 3)
: CellSet(name, dimensionality),
: NumberOfPoints(numpoints),
CellSet(name, dimensionality),
ConnectivityLength(-1),
NumberOfCells(-1)
{
}
VTKM_CONT_EXPORT
CellSetExplicit(int dimensionality)
: CellSet(std::string(), dimensionality),
CellSetExplicit(vtkm::Id numpoints, int dimensionality)
: NumberOfPoints(numpoints),
CellSet(std::string(), dimensionality),
ConnectivityLength(-1),
NumberOfCells(-1)
{
@ -80,12 +86,23 @@ public:
return this->PointToCell.GetNumberOfElements();
}
virtual vtkm::Id GetNumberOfPoints() const
{
return this->NumberOfPoints;
}
VTKM_CONT_EXPORT
vtkm::Id GetSchedulingRange(vtkm::TopologyElementTagCell) const
{
return this->GetNumberOfCells();
}
VTKM_CONT_EXPORT
vtkm::Id GetSchedulingRange(vtkm::TopologyElementTagPoint) const
{
return this->GetNumberOfPoints();
}
VTKM_CONT_EXPORT
vtkm::Id GetNumberOfPointsInCell(vtkm::Id cellIndex) const
{
@ -236,6 +253,76 @@ public:
connectivity.IndexOffsets.PrepareForInput(Device()));
}
void CreateCellToPointConnectivity()
{
std::multimap<int,int> cells_of_nodes;
vtkm::Id maxNodeID = 0;
vtkm::Id numCells = GetNumberOfCells();
vtkm::Id numPoints = GetNumberOfPoints();
for (vtkm::Id cell = 0, cindex = 0; cell < numCells; ++cell)
{
vtkm::Id npts = this->PointToCell.NumIndices.GetPortalControl().Get(cell);
for (int pt=0; pt<npts; ++pt)
{
vtkm::Id index = this->PointToCell.Connectivity.GetPortalControl().Get(cindex++);
if (index > maxNodeID)
maxNodeID = index;
cells_of_nodes.insert(std::pair<int,int>(index,cell));
}
}
std::vector<vtkm::Id> shapes;
std::vector<vtkm::Id> numindices;
std::vector<vtkm::Id> conn;
vtkm::Id filled_array_to_node = 0;
for (std::multimap<int,int>::iterator iter = cells_of_nodes.begin();
iter != cells_of_nodes.end(); iter++)
{
int node = iter->first;
while (filled_array_to_node <= node)
{
// add empty spots to skip nodes not referenced by our cells
// also add a spot for the current one
++filled_array_to_node;
shapes.push_back(VTKM_VERTEX);
numindices.push_back(0);
}
vtkm::Id cell = iter->second;
conn.push_back(cell);
++numindices[numindices.size()-1];
}
while (filled_array_to_node < numPoints)
{
// add empty spots for tail nodes not referenced by our cells
++filled_array_to_node;
shapes.push_back(VTKM_VERTEX);
numindices.push_back(0);
}
///\todo: THIS IS A COPY OF "FillViaCopy", because that method
/// is specific to PointToCell. Should make it non-specific and call
/// it instead.
this->CellToPoint.Shapes.Allocate( static_cast<vtkm::Id>(shapes.size()) );
std::copy(shapes.begin(), shapes.end(),
vtkm::cont::ArrayPortalToIteratorBegin(
this->CellToPoint.Shapes.GetPortalControl()));
this->CellToPoint.NumIndices.Allocate( static_cast<vtkm::Id>(numindices.size()) );
std::copy(numindices.begin(), numindices.end(),
vtkm::cont::ArrayPortalToIteratorBegin(
this->CellToPoint.NumIndices.GetPortalControl()));
this->CellToPoint.Connectivity.Allocate( static_cast<vtkm::Id>(conn.size()) );
std::copy(conn.begin(), conn.end(),
vtkm::cont::ArrayPortalToIteratorBegin(
this->CellToPoint.Connectivity.GetPortalControl()));
this->CellToPoint.ElementsValid = true;
this->CellToPoint.IndexOffsetsValid = false;
}
virtual void PrintSummary(std::ostream &out) const
{
out << " ExplicitCellSet: " << this->Name
@ -306,6 +393,7 @@ private:
// cells.
vtkm::Id ConnectivityLength;
vtkm::Id NumberOfCells;
vtkm::Id NumberOfPoints;
};
namespace detail {

@ -150,8 +150,9 @@ MakeTestDataSet::Make3DExplicitDataSet0()
conn.push_back(3);
conn.push_back(4);
vtkm::cont::CellSetExplicit<> cellSet("cells", 2);
vtkm::cont::CellSetExplicit<> cellSet(nVerts, "cells", 2);
cellSet.FillViaCopy(shapes, numindices, conn);
cellSet.CreateCellToPointConnectivity();
dataSet.AddCellSet(cellSet);
@ -184,14 +185,15 @@ MakeTestDataSet::Make3DExplicitDataSet1()
vtkm::Float32 cellvar[2] = {100.1f, 100.2f};
dataSet.AddField(Field("cellvar", 1, vtkm::cont::Field::ASSOC_CELL_SET, "cells", cellvar, 2));
vtkm::cont::CellSetExplicit<> cellSet("cells", 2);
vtkm::cont::CellSetExplicit<> cellSet(nVerts, "cells", 2);
cellSet.PrepareToAddCells(2, 7);
cellSet.AddCell(vtkm::CELL_SHAPE_TRIANGLE, 3, make_Vec<vtkm::Id>(0,1,2));
cellSet.AddCell(vtkm::CELL_SHAPE_QUAD, 4, make_Vec<vtkm::Id>(2,1,3,4));
cellSet.CompleteAddingCells();
//todo this need to be a reference/shared_ptr style class
cellSet.CreateCellToPointConnectivity();
dataSet.AddCellSet(cellSet);
return dataSet;
@ -251,7 +253,7 @@ MakeTestDataSet::Make3DExplicitDataSetCowNose(double *pBounds)
dataSet.AddCoordinateSystem(
vtkm::cont::CoordinateSystem("coordinates", 1, coordinates, nVerts));
vtkm::cont::CellSetExplicit<> cellSet("cells", 2);
vtkm::cont::CellSetExplicit<> cellSet(nVerts, "cells", 2);
cellSet.PrepareToAddCells(nPointIds/3, nPointIds);
for (vtkm::Id i=0; i<nPointIds/3; i++)
@ -264,7 +266,8 @@ MakeTestDataSet::Make3DExplicitDataSetCowNose(double *pBounds)
}
cellSet.CompleteAddingCells();
//todo this need to be a reference/shared_ptr style class
cellSet.CreateCellToPointConnectivity();
dataSet.AddCellSet(cellSet);
// copy bounds

@ -53,7 +53,7 @@ vtkm::cont::DataSet MakePointElevationTestDataSet()
dataSet.AddCoordinateSystem(
vtkm::cont::CoordinateSystem("coordinates", 1, coordinates));
vtkm::cont::CellSetExplicit<> cellSet("cells", 3);
vtkm::cont::CellSetExplicit<> cellSet(numVerts, "cells", 3);
cellSet.PrepareToAddCells(numCells, numCells * 4);
for (vtkm::Id j = 0; j < dim - 1; ++j)
{

@ -88,7 +88,7 @@ vtkm::cont::DataSet RunVertexClustering(vtkm::cont::DataSet &dataSet,
//typedef typename vtkm::cont::ArrayHandleConstant<vtkm::Id>::StorageTag ConstantStorage;
//typedef typename vtkm::cont::ArrayHandleImplicit<vtkm::Id, CounterOfThree>::StorageTag CountingStorage;
vtkm::cont::CellSetExplicit<> newCellSet("cells", 0);
vtkm::cont::CellSetExplicit<> newCellSet(pointArray.GetNumberOfValues(), "cells", 0);
newCellSet.Fill(
copyFromImplicit(vtkm::cont::make_ArrayHandleConstant<vtkm::Id>(vtkm::CELL_SHAPE_TRIANGLE, cells)),

@ -96,12 +96,43 @@ public:
}
};
class AverageCellToPointValue :
public vtkm::worklet::WorkletMapTopology<vtkm::TopologyElementTagCell,
vtkm::TopologyElementTagPoint>
{
public:
typedef void ControlSignature(FieldInFrom<Scalar> inCells,
TopologyIn topology,
FieldOut<Scalar> outPoints);
typedef void ExecutionSignature(_1, _3, FromCount);
typedef _2 InputDomain;
VTKM_CONT_EXPORT
AverageCellToPointValue() { }
template<typename CellVecType, typename OutType>
VTKM_EXEC_EXPORT
void operator()(const CellVecType &cellValues,
OutType &avgVal,
const vtkm::IdComponent &numCellIDs) const
{
//simple functor that returns the average cell Value.
avgVal = static_cast<OutType>(cellValues[0]);
for (vtkm::IdComponent cellIndex = 1; cellIndex < numCellIDs; ++cellIndex)
{
avgVal += static_cast<OutType>(cellValues[cellIndex]);
}
avgVal = avgVal / static_cast<OutType>(numCellIDs);
}
};
}
namespace {
static void TestMaxPointOrCell();
static void TestAvgPointToCell();
static void TestAvgCellToPoint();
void TestWorkletMapTopologyExplicit()
{
@ -112,6 +143,7 @@ void TestWorkletMapTopologyExplicit()
TestMaxPointOrCell();
TestAvgPointToCell();
TestAvgCellToPoint();
}
@ -195,6 +227,46 @@ TestAvgPointToCell()
"Wrong result for PointToCellAverage worklet");
}
static void
TestAvgCellToPoint()
{
std::cout<<"Testing AvgCellToPoint worklet"<<std::endl;
vtkm::cont::testing::MakeTestDataSet testDataSet;
vtkm::cont::DataSet dataSet = testDataSet.Make3DExplicitDataSet1();
//Run a worklet to populate a point centered field.
//Here, we're filling it with test values.
vtkm::cont::Field f("outpointvar",
1,
vtkm::cont::Field::ASSOC_POINTS,
vtkm::Float32());
dataSet.AddField(f);
VTKM_TEST_ASSERT(dataSet.GetNumberOfCellSets() == 1,
"Incorrect number of cell sets");
VTKM_TEST_ASSERT(dataSet.GetNumberOfFields() == 6,
"Incorrect number of fields");
vtkm::worklet::DispatcherMapTopology< ::test_explicit::AverageCellToPointValue > dispatcher;
dispatcher.Invoke(dataSet.GetField("cellvar").GetData(),
dataSet.GetCellSet(),
dataSet.GetField("outpointvar").GetData());
//make sure we got the right answer.
vtkm::cont::ArrayHandle<vtkm::Float32> res;
res = dataSet.GetField("outpointvar").GetData().CastToArrayHandle(vtkm::Float32(),
VTKM_DEFAULT_STORAGE_TAG());
VTKM_TEST_ASSERT(test_equal(res.GetPortalConstControl().Get(0), 100.1f),
"Wrong result for CellToPointAverage worklet");
VTKM_TEST_ASSERT(test_equal(res.GetPortalConstControl().Get(1), 100.15f),
"Wrong result for CellToPointAverage worklet");
}
} // anonymous namespace
int UnitTestWorkletMapTopologyExplicit(int, char *[])

@ -96,12 +96,42 @@ public:
}
};
class AverageCellToPointValue :
public vtkm::worklet::WorkletMapTopology<vtkm::TopologyElementTagCell,
vtkm::TopologyElementTagPoint>
{
public:
typedef void ControlSignature(FieldInFrom<Scalar> inCells,
TopologyIn topology,
FieldOut<Scalar> outPoints);
typedef void ExecutionSignature(_1, _3, FromCount);
typedef _2 InputDomain;
VTKM_CONT_EXPORT
AverageCellToPointValue() { }
template<typename CellVecType, typename OutType>
VTKM_EXEC_EXPORT
void operator()(const CellVecType &cellValues,
OutType &avgVal,
const vtkm::IdComponent &numCellIDs) const
{
//simple functor that returns the average cell Value.
avgVal = static_cast<OutType>(cellValues[0]);
for (vtkm::IdComponent cellIndex = 1; cellIndex < numCellIDs; ++cellIndex)
{
avgVal += static_cast<OutType>(cellValues[cellIndex]);
}
avgVal = avgVal / static_cast<OutType>(numCellIDs);
}
};
}
namespace {
static void TestMaxPointOrCell();
static void TestAvgPointToCell();
static void TestAvgCellToPoint();
void TestWorkletMapTopologyRegular()
{
@ -112,6 +142,7 @@ void TestWorkletMapTopologyRegular()
TestMaxPointOrCell();
TestAvgPointToCell();
TestAvgCellToPoint();
}
static void
@ -168,7 +199,7 @@ TestAvgPointToCell()
//Run a worklet to populate a cell centered field.
//Here, we're filling it with test values.
vtkm::cont::Field f("outcellvar",
1,
0,
vtkm::cont::Field::ASSOC_CELL_SET,
std::string("cells"),
vtkm::Float32());
@ -201,6 +232,50 @@ TestAvgPointToCell()
"Wrong result for PointToCellAverage worklet");
}
static void
TestAvgCellToPoint()
{
std::cout<<"Testing AvgCellToPoint worklet"<<std::endl;
vtkm::cont::testing::MakeTestDataSet testDataSet;
vtkm::cont::DataSet dataSet = testDataSet.Make2DRegularDataSet0();
//Run a worklet to populate a point centered field.
//Here, we're filling it with test values.
vtkm::cont::Field f("outpointvar",
1,
vtkm::cont::Field::ASSOC_POINTS,
vtkm::Float32());
dataSet.AddField(f);
VTKM_TEST_ASSERT(test_equal(dataSet.GetNumberOfCellSets(), 1),
"Incorrect number of cell sets");
VTKM_TEST_ASSERT(test_equal(dataSet.GetNumberOfFields(), 5),
"Incorrect number of fields");
vtkm::worklet::DispatcherMapTopology< ::test_regular::AverageCellToPointValue > dispatcher;
dispatcher.Invoke(dataSet.GetField("cellvar").GetData(),
// We know that the cell set is a structured 2D grid and
// The worklet does not work with general types because
// of the way we get cell indices. We need to make that
// part more flexible.
dataSet.GetCellSet(0).ResetCellSetList(
vtkm::cont::CellSetListTagStructured2D()),
dataSet.GetField("outpointvar").GetData());
//make sure we got the right answer.
vtkm::cont::ArrayHandle<vtkm::Float32> res;
res = dataSet.GetField(4).GetData().CastToArrayHandle(vtkm::Float32(),
VTKM_DEFAULT_STORAGE_TAG());
VTKM_TEST_ASSERT(test_equal(res.GetPortalConstControl().Get(0), 100.1f),
"Wrong result for CellToPointAverage worklet");
VTKM_TEST_ASSERT(test_equal(res.GetPortalConstControl().Get(1), 150.1f),
"Wrong result for CellToPointAverage worklet");
}
} // anonymous namespace
int UnitTestWorkletMapTopologyRegular(int, char *[])