Shuenhoy/vtk-m
1
0
Fork 0
mirror of https://gitlab.kitware.com/vtk/vtk-m synced 2024-09-08 13:23:51 +00:00
Code Issues 2 Actions Packages Projects Releases Wiki Activity

Adding support for evaluating 2D meshes.

Browse Source 
- Adding support for evaluating 2D meshes using Rectilinear Grids
- Adding support for particle advection for 2D meshes
This commit is contained in:
Abhishek Yenpure 2019-07-11 10:45:30 -07:00
parent b25599cb2e
commit c2fa0467c1
6 changed files with 99 additions and 41 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

65
vtkm/cont/CellLocatorRectilinearGrid.cxx
View File

@ -24,7 +24,8 @@ CellLocatorRectilinearGrid::CellLocatorRectilinearGrid() = default;
CellLocatorRectilinearGrid::~CellLocatorRectilinearGrid() = default;
using StructuredType = vtkm::cont::CellSetStructured<3>;
using Structured2DType = vtkm::cont::CellSetStructured<2>;
using Structured3DType = vtkm::cont::CellSetStructured<3>;
using AxisHandle = vtkm::cont::ArrayHandle<vtkm::FloatDefault>;
using RectilinearType = vtkm::cont::ArrayHandleCartesianProduct<AxisHandle, AxisHandle, AxisHandle>;
@ -35,18 +36,33 @@ void CellLocatorRectilinearGrid::Build()
if (!coords.GetData().IsType<RectilinearType>())
throw vtkm::cont::ErrorInternal("Coordinates are not rectilinear.");
if (!cellSet.IsSameType(StructuredType()))
throw vtkm::cont::ErrorInternal("Cells are not 3D structured.");
vtkm::Vec<vtkm::Id, 3> celldims =
cellSet.Cast<StructuredType>().GetSchedulingRange(vtkm::TopologyElementTagCell());
this->PlaneSize = celldims[0] * celldims[1];
this->RowSize = celldims[0];
if (cellSet.IsSameType(Structured2DType()))
{
this->Is3D = false;
vtkm::Vec<vtkm::Id, 2> celldims =
cellSet.Cast<Structured2DType>().GetSchedulingRange(vtkm::TopologyElementTagCell());
this->PlaneSize = celldims[0] * celldims[1];
this->RowSize = celldims[0];
}
else if (cellSet.IsSameType(Structured3DType()))
{
this->Is3D = true;
vtkm::Vec<vtkm::Id, 3> celldims =
cellSet.Cast<Structured3DType>().GetSchedulingRange(vtkm::TopologyElementTagCell());
this->PlaneSize = celldims[0] * celldims[1];
this->RowSize = celldims[0];
}
else
{
throw vtkm::cont::ErrorInternal("Cells are not structured.");
}
}
namespace
{
template <vtkm::IdComponent dimensions>
struct CellLocatorRectilinearGridPrepareForExecutionFunctor
{
template <typename DeviceAdapter, typename... Args>
@ -54,7 +70,7 @@ struct CellLocatorRectilinearGridPrepareForExecutionFunctor
vtkm::cont::VirtualObjectHandle<vtkm::exec::CellLocator>& execLocator,
Args&&... args) const
{
using ExecutionType = vtkm::exec::CellLocatorRectilinearGrid<DeviceAdapter>;
using ExecutionType = vtkm::exec::CellLocatorRectilinearGrid<DeviceAdapter, dimensions>;
ExecutionType* execObject = new ExecutionType(std::forward<Args>(args)..., DeviceAdapter());
execLocator.Reset(execObject);
return true;
@ -65,14 +81,29 @@ struct CellLocatorRectilinearGridPrepareForExecutionFunctor
const vtkm::exec::CellLocator* CellLocatorRectilinearGrid::PrepareForExecution(
vtkm::cont::DeviceAdapterId device) const
{
const bool success = vtkm::cont::TryExecuteOnDevice(
device,
CellLocatorRectilinearGridPrepareForExecutionFunctor(),
this->ExecutionObjectHandle,
this->PlaneSize,
this->RowSize,
this->GetCellSet().template Cast<StructuredType>(),
this->GetCoordinates().GetData().template Cast<RectilinearType>());
bool success = false;
if (this->Is3D)
{
success = vtkm::cont::TryExecuteOnDevice(
device,
CellLocatorRectilinearGridPrepareForExecutionFunctor<3>(),
this->ExecutionObjectHandle,
this->PlaneSize,
this->RowSize,
this->GetCellSet().template Cast<Structured3DType>(),
this->GetCoordinates().GetData().template Cast<RectilinearType>());
}
else
{
success = vtkm::cont::TryExecuteOnDevice(
device,
CellLocatorRectilinearGridPrepareForExecutionFunctor<2>(),
this->ExecutionObjectHandle,
this->PlaneSize,
this->RowSize,
this->GetCellSet().template Cast<Structured2DType>(),
this->GetCoordinates().GetData().template Cast<RectilinearType>());
}
if (!success)
{
throwFailedRuntimeDeviceTransfer("CellLocatorRectilinearGrid", device);

1
vtkm/cont/CellLocatorRectilinearGrid.h
View File

@ -35,6 +35,7 @@ private:
vtkm::Bounds Bounds;
vtkm::Id PlaneSize;
vtkm::Id RowSize;
bool Is3D = true;
mutable vtkm::cont::VirtualObjectHandle<vtkm::exec::CellLocator> ExecutionObjectHandle;
};

26
vtkm/exec/CellLocatorRectilinearGrid.h
View File

@ -28,7 +28,7 @@ namespace vtkm
namespace exec
{
template <typename DeviceAdapter>
template <typename DeviceAdapter, vtkm::IdComponent dimensions>
class VTKM_ALWAYS_EXPORT CellLocatorRectilinearGrid final : public vtkm::exec::CellLocator
{
private:
@ -36,7 +36,7 @@ private:
using ToType = vtkm::TopologyElementTagCell;
using CellSetPortal = vtkm::exec::ConnectivityStructured<vtkm::TopologyElementTagPoint,
vtkm::TopologyElementTagCell,
3>;
dimensions>;
using AxisHandle = vtkm::cont::ArrayHandle<vtkm::FloatDefault>;
using RectilinearType =
vtkm::cont::ArrayHandleCartesianProduct<AxisHandle, AxisHandle, AxisHandle>;
@ -48,7 +48,7 @@ public:
VTKM_CONT
CellLocatorRectilinearGrid(const vtkm::Id planeSize,
const vtkm::Id rowSize,
const vtkm::cont::CellSetStructured<3>& cellSet,
const vtkm::cont::CellSetStructured<dimensions>& cellSet,
const RectilinearType& coords,
DeviceAdapter)
: PlaneSize(planeSize)
@ -58,16 +58,18 @@ public:
, PointDimensions(cellSet.GetPointDimensions())
{
this->AxisPortals[0] = Coords.GetFirstPortal();
this->AxisPortals[1] = Coords.GetSecondPortal();
this->AxisPortals[2] = Coords.GetThirdPortal();
this->MinPoint[0] = coords.GetPortalConstControl().GetFirstPortal().Get(0);
this->MinPoint[1] = coords.GetPortalConstControl().GetSecondPortal().Get(0);
this->MinPoint[2] = coords.GetPortalConstControl().GetThirdPortal().Get(0);
this->MaxPoint[0] = coords.GetPortalConstControl().GetFirstPortal().Get(PointDimensions[0] - 1);
this->AxisPortals[1] = Coords.GetSecondPortal();
this->MinPoint[1] = coords.GetPortalConstControl().GetSecondPortal().Get(0);
this->MaxPoint[1] =
coords.GetPortalConstControl().GetSecondPortal().Get(PointDimensions[1] - 1);
if
constexpr(dimensions == 2) return;
this->AxisPortals[2] = Coords.GetThirdPortal();
this->MinPoint[2] = coords.GetPortalConstControl().GetThirdPortal().Get(0);
this->MaxPoint[2] = coords.GetPortalConstControl().GetThirdPortal().Get(PointDimensions[2] - 1);
}
@ -85,6 +87,8 @@ public:
inside = false;
if (point[1] < this->MinPoint[1] || point[1] > this->MaxPoint[1])
inside = false;
if
constexpr(dimensions == 2) return inside;
if (point[2] < this->MinPoint[2] || point[2] > this->MaxPoint[2])
inside = false;
return inside;
@ -104,7 +108,7 @@ public:
// Get the Cell Id from the point.
vtkm::Vec<vtkm::Id, 3> logicalCell(0, 0, 0);
for (vtkm::Int32 dim = 0; dim < 3; ++dim)
for (vtkm::Int32 dim = 0; dim < dimensions; ++dim)
{
//
// When searching for points, we consider the max value of the cell
@ -167,7 +171,7 @@ private:
CellSetPortal CellSet;
RectilinearPortalType Coords;
AxisPortalType AxisPortals[3];
vtkm::Id3 PointDimensions;
vtkm::Vec<vtkm::Id, dimensions> PointDimensions;
vtkm::Vec<vtkm::FloatDefault, 3> MinPoint;
vtkm::Vec<vtkm::FloatDefault, 3> MaxPoint;
};

39
vtkm/worklet/particleadvection/CellInterpolationHelper.h
View File

@ -49,9 +49,10 @@ public:
StructuredCellInterpolationHelper() = default;
VTKM_CONT
StructuredCellInterpolationHelper(vtkm::Id3 cellDims, vtkm::Id3 pointDims)
StructuredCellInterpolationHelper(vtkm::Id3 cellDims, vtkm::Id3 pointDims, bool is3D)
: CellDims(cellDims)
, PointDims(pointDims)
, Is3D(is3D)
{
}
@ -71,11 +72,16 @@ public:
indices.Append(indices[0] + 1);
indices.Append(indices[1] + PointDims[0]);
indices.Append(indices[2] - 1);
if (!this->Is3D)
{
cellShape = static_cast<vtkm::UInt8>(vtkm::CELL_SHAPE_QUAD);
numVerts = 4;
return;
}
indices.Append(indices[0] + PointDims[0] * PointDims[1]);
indices.Append(indices[4] + 1);
indices.Append(indices[5] + PointDims[0]);
indices.Append(indices[6] - 1);
cellShape = static_cast<vtkm::UInt8>(vtkm::CELL_SHAPE_HEXAHEDRON);
numVerts = 8;
}
@ -83,6 +89,7 @@ public:
private:
vtkm::Id3 CellDims;
vtkm::Id3 PointDims;
bool Is3D = true;
};
template <typename DeviceAdapter>
@ -196,21 +203,34 @@ public:
class StructuredCellInterpolationHelper : public vtkm::cont::CellInterpolationHelper
{
public:
using StructuredType = vtkm::cont::CellSetStructured<3>;
using Structured2DType = vtkm::cont::CellSetStructured<2>;
using Structured3DType = vtkm::cont::CellSetStructured<3>;
StructuredCellInterpolationHelper() = default;
VTKM_CONT
StructuredCellInterpolationHelper(const vtkm::cont::DynamicCellSet& cellSet)
{
if (cellSet.IsSameType(StructuredType()))
if (cellSet.IsSameType(Structured2DType()))
{
CellDims = cellSet.Cast<StructuredType>().GetSchedulingRange(vtkm::TopologyElementTagCell());
PointDims =
cellSet.Cast<StructuredType>().GetSchedulingRange(vtkm::TopologyElementTagPoint());
this->Is3D = false;
vtkm::Id2 cellDims =
cellSet.Cast<Structured2DType>().GetSchedulingRange(vtkm::TopologyElementTagCell());
vtkm::Id2 pointDims =
cellSet.Cast<Structured2DType>().GetSchedulingRange(vtkm::TopologyElementTagPoint());
this->CellDims = vtkm::Id3(cellDims[0], cellDims[1], 0);
this->PointDims = vtkm::Id3(pointDims[0], pointDims[1], 1);
}
else if (cellSet.IsSameType(Structured3DType()))
{
this->Is3D = true;
this->CellDims =
cellSet.Cast<Structured3DType>().GetSchedulingRange(vtkm::TopologyElementTagCell());
this->PointDims =
cellSet.Cast<Structured3DType>().GetSchedulingRange(vtkm::TopologyElementTagPoint());
}
else
throw vtkm::cont::ErrorBadType("Cell set is not 3D structured type");
throw vtkm::cont::ErrorBadType("Cell set is not structured type");
}
VTKM_CONT
@ -225,7 +245,7 @@ public:
}
using ExecutionType = vtkm::exec::StructuredCellInterpolationHelper;
ExecutionType* execObject = new ExecutionType(this->CellDims, this->PointDims);
ExecutionType* execObject = new ExecutionType(this->CellDims, this->PointDims, this->Is3D);
this->ExecHandle.Reset(execObject);
return this->ExecHandle.PrepareForExecution(deviceId);
@ -234,6 +254,7 @@ public:
private:
vtkm::Id3 CellDims;
vtkm::Id3 PointDims;
bool Is3D = true;
mutable HandleType ExecHandle;
};

8
vtkm/worklet/particleadvection/GridEvaluators.h
View File

@ -127,7 +127,8 @@ public:
using AxisHandle = vtkm::cont::ArrayHandle<vtkm::FloatDefault>;
using RectilinearType =
vtkm::cont::ArrayHandleCartesianProduct<AxisHandle, AxisHandle, AxisHandle>;
using StructuredType = vtkm::cont::CellSetStructured<3>;
using Structured2DType = vtkm::cont::CellSetStructured<2>;
using Structured3DType = vtkm::cont::CellSetStructured<3>;
VTKM_CONT
GridEvaluator() = default;
@ -139,7 +140,7 @@ public:
: Vectors(field)
, Bounds(coordinates.GetBounds())
{
if (cellset.IsSameType(StructuredType()))
if (cellset.IsSameType(Structured2DType()) || cellset.IsSameType(Structured3DType()))
{
if (coordinates.GetData().IsType<UniformType>())
{
@ -149,8 +150,7 @@ public:
locator.Update();
this->Locator = std::make_shared<vtkm::cont::CellLocatorUniformGrid>(locator);
}
else if (coordinates.GetData().IsType<RectilinearType>() &&
cellset.IsSameType(StructuredType()))
else if (coordinates.GetData().IsType<RectilinearType>())
{
vtkm::cont::CellLocatorRectilinearGrid locator;
locator.SetCoordinates(coordinates);

1
vtkm/worklet/particleadvection/ParticleAdvectionWorklets.h
View File

@ -43,6 +43,7 @@ public:
IntegralCurveType& integralCurve) const
{
vtkm::Vec<ScalarType, 3> inpos = integralCurve.GetPos(idx);
std::cout << "Advecting index " << idx << " point " << inpos << std::endl;
vtkm::Vec<ScalarType, 3> outpos;
ScalarType time = integralCurve.GetTime(idx);
ParticleStatus status;