diff --git a/docs/users-guide/examples/CMakeLists.txt b/docs/users-guide/examples/CMakeLists.txt index c383b7ca1..a08a18b77 100644 --- a/docs/users-guide/examples/CMakeLists.txt +++ b/docs/users-guide/examples/CMakeLists.txt @@ -49,6 +49,9 @@ set(examples_device GuideExampleFilterDataSetWithField.cxx GuideExampleGenerateMeshConstantShape.cxx GuideExampleSimpleAlgorithm.cxx + GuideExampleSimpleHistogram.cxx + GuideExampleSumOfAngles.cxx + GuideExampleTriangleQuality.cxx GuideExampleUnknownArrayHandle.cxx GuideExampleUseWorkletMapField.cxx GuideExampleUseWorkletPointNeighborhood.cxx diff --git a/docs/users-guide/examples/GuideExampleSimpleHistogram.cxx b/docs/users-guide/examples/GuideExampleSimpleHistogram.cxx new file mode 100644 index 000000000..9d4f04180 --- /dev/null +++ b/docs/users-guide/examples/GuideExampleSimpleHistogram.cxx @@ -0,0 +1,122 @@ +//============================================================================= +// +// Copyright (c) Kitware, Inc. +// All rights reserved. +// See LICENSE.txt for details. +// +// This software is distributed WITHOUT ANY WARRANTY; without even +// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR +// PURPOSE. See the above copyright notice for more information. +// +//============================================================================= + +#include +#include +#include + +#include +#include +#include +#include + +#include + +#include + +struct SimpleHistogramStruct +{ + //// + //// BEGIN-EXAMPLE SimpleHistogram + //// + struct CountBins : vtkm::worklet::WorkletMapField + { + using ControlSignature = void(FieldIn data, AtomicArrayInOut histogramBins); + using ExecutionSignature = void(_1, _2); + using InputDomain = _1; + + vtkm::Range HistogramRange; + vtkm::Id NumberOfBins; + + VTKM_CONT + CountBins(const vtkm::Range& histogramRange, vtkm::Id& numBins) + : HistogramRange(histogramRange) + , NumberOfBins(numBins) + { + } + + template + VTKM_EXEC void operator()(T value, const AtomicArrayType& histogramBins) const + { + vtkm::Float64 interp = + (value - this->HistogramRange.Min) / this->HistogramRange.Length(); + vtkm::Id bin = static_cast(interp * this->NumberOfBins); + if (bin < 0) + { + bin = 0; + } + if (bin >= this->NumberOfBins) + { + bin = this->NumberOfBins - 1; + } + + histogramBins.Add(bin, 1); + } + }; + //// + //// END-EXAMPLE SimpleHistogram + //// + + template + VTKM_CONT static vtkm::cont::ArrayHandle Run(const InputArray& input, + vtkm::Id numberOfBins) + { + VTKM_IS_ARRAY_HANDLE(InputArray); + + // Histograms only work on scalar values + using ValueType = typename InputArray::ValueType; + VTKM_STATIC_ASSERT_MSG( + (std::is_same::HasMultipleComponents, + vtkm::VecTraitsTagSingleComponent>::value), + "Histiogram input not a scalar value."); + + vtkm::Range range = vtkm::cont::ArrayRangeCompute(input).ReadPortal().Get(0); + + // Initialize histogram to 0 + vtkm::cont::ArrayHandle histogram; + vtkm::cont::Algorithm::Copy( + vtkm::cont::ArrayHandleConstant(0, numberOfBins), histogram); + + CountBins histogramWorklet(range, numberOfBins); + + vtkm::cont::Invoker invoker; + invoker(histogramWorklet, input, histogram); + + return histogram; + } +}; + +VTKM_CONT +static inline void TrySimpleHistogram() +{ + std::cout << "Try Simple Histogram" << std::endl; + + static const vtkm::Id ARRAY_SIZE = 100; + vtkm::cont::ArrayHandle inputArray; + inputArray.Allocate(ARRAY_SIZE); + SetPortal(inputArray.WritePortal()); + + vtkm::cont::ArrayHandle histogram = + SimpleHistogramStruct::Run(inputArray, ARRAY_SIZE / 2); + + VTKM_TEST_ASSERT(histogram.GetNumberOfValues() == ARRAY_SIZE / 2, "Bad array size"); + for (vtkm::Id index = 0; index < histogram.GetNumberOfValues(); ++index) + { + vtkm::Int32 binSize = histogram.ReadPortal().Get(index); + VTKM_TEST_ASSERT(binSize == 2, "Bad bin size."); + } +} + +int GuideExampleSimpleHistogram(int argc, char* argv[]) +{ + return vtkm::cont::testing::Testing::Run(TrySimpleHistogram, argc, argv); +} diff --git a/docs/users-guide/examples/GuideExampleSumOfAngles.cxx b/docs/users-guide/examples/GuideExampleSumOfAngles.cxx new file mode 100644 index 000000000..ff17ac5a5 --- /dev/null +++ b/docs/users-guide/examples/GuideExampleSumOfAngles.cxx @@ -0,0 +1,196 @@ +//============================================================================= +// +// Copyright (c) Kitware, Inc. +// All rights reserved. +// See LICENSE.txt for details. +// +// This software is distributed WITHOUT ANY WARRANTY; without even +// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR +// PURPOSE. See the above copyright notice for more information. +// +//============================================================================= + +#include + +#include +#include +#include + +#include + +#include +#include + +#include +#include + +#include + +namespace +{ + +// This worklet computes the sum of the angles of all polygons connected +// to each point. This sum is related (but not equal to) the Gaussian +// curvature of the surface. A flat mesh will have a sum equal to 2 pi. +// A concave or convex surface will have a sum less than 2 pi. A saddle +// will have a sum greater than 2 pi. The actual Gaussian curvature is +// equal to (2 pi - angle sum)/A where A is an area of influence (which +// we are not calculating here). See +// http://computergraphics.stackexchange.com/questions/1718/what-is-the-simplest-way-to-compute-principal-curvature-for-a-mesh-triangle#1721 +// or the publication "Discrete Differential-Geometry Operators for +// Triangulated 2-Manifolds" by Mayer et al. (Equation 9). +//// +//// BEGIN-EXAMPLE SumOfAngles +//// +struct SumOfAngles : vtkm::worklet::WorkletVisitPointsWithCells +{ + using ControlSignature = void(CellSetIn inputCells, + WholeCellSetIn<>, // Same as inputCells + WholeArrayIn pointCoords, + FieldOutPoint angleSum); + using ExecutionSignature = void(CellIndices incidentCells, + InputIndex pointIndex, + _2 cellSet, + _3 pointCoordsPortal, + _4 outSum); + using InputDomain = _1; + + template + VTKM_EXEC void operator()(const IncidentCellVecType& incidentCells, + vtkm::Id pointIndex, + const CellSetType& cellSet, + const PointCoordsPortalType& pointCoordsPortal, + SumType& outSum) const + { + using CoordType = typename PointCoordsPortalType::ValueType; + + CoordType thisPoint = pointCoordsPortal.Get(pointIndex); + + outSum = 0; + for (vtkm::IdComponent incidentCellIndex = 0; + incidentCellIndex < incidentCells.GetNumberOfComponents(); + ++incidentCellIndex) + { + // Get information about incident cell. + vtkm::Id cellIndex = incidentCells[incidentCellIndex]; + typename CellSetType::CellShapeTag cellShape = cellSet.GetCellShape(cellIndex); + typename CellSetType::IndicesType cellConnections = cellSet.GetIndices(cellIndex); + vtkm::IdComponent numPointsInCell = cellSet.GetNumberOfIndices(cellIndex); + vtkm::IdComponent numEdges; + vtkm::exec::CellEdgeNumberOfEdges(numPointsInCell, cellShape, numEdges); + + // Iterate over all edges and find the first one with pointIndex. + // Use that to find the first vector. + vtkm::IdComponent edgeIndex = -1; + CoordType vec1; + while (true) + { + ++edgeIndex; + if (edgeIndex >= numEdges) + { + this->RaiseError("Bad cell. Could not find two incident edges."); + return; + } + vtkm::IdComponent2 edge; + vtkm::exec::CellEdgeLocalIndex( + numPointsInCell, 0, edgeIndex, cellShape, edge[0]); + vtkm::exec::CellEdgeLocalIndex( + numPointsInCell, 1, edgeIndex, cellShape, edge[1]); + if (cellConnections[edge[0]] == pointIndex) + { + vec1 = pointCoordsPortal.Get(cellConnections[edge[1]]) - thisPoint; + break; + } + else if (cellConnections[edge[1]] == pointIndex) + { + vec1 = pointCoordsPortal.Get(cellConnections[edge[0]]) - thisPoint; + break; + } + else + { + // Continue to next iteration of loop. + } + } + + // Continue iteration over remaining edges and find the second one with + // pointIndex. Use that to find the second vector. + CoordType vec2; + while (true) + { + ++edgeIndex; + if (edgeIndex >= numEdges) + { + this->RaiseError("Bad cell. Could not find two incident edges."); + return; + } + vtkm::IdComponent2 edge; + vtkm::exec::CellEdgeLocalIndex( + numPointsInCell, 0, edgeIndex, cellShape, edge[0]); + vtkm::exec::CellEdgeLocalIndex( + numPointsInCell, 1, edgeIndex, cellShape, edge[1]); + if (cellConnections[edge[0]] == pointIndex) + { + vec2 = pointCoordsPortal.Get(cellConnections[edge[1]]) - thisPoint; + break; + } + else if (cellConnections[edge[1]] == pointIndex) + { + vec2 = pointCoordsPortal.Get(cellConnections[edge[0]]) - thisPoint; + break; + } + else + { + // Continue to next iteration of loop. + } + } + + // The dot product of two unit vectors is equal to the cosine of the + // angle between them. + vtkm::Normalize(vec1); + vtkm::Normalize(vec2); + SumType cosine = static_cast(vtkm::Dot(vec1, vec2)); + + outSum += vtkm::ACos(cosine); + } + } +}; +//// +//// END-EXAMPLE SumOfAngles +//// + +VTKM_CONT +static void TrySumOfAngles() +{ + std::cout << "Read input data" << std::endl; + vtkm::io::VTKDataSetReader reader(vtkm::cont::testing::Testing::GetTestDataBasePath() + + "unstructured/cow.vtk"); + vtkm::cont::DataSet dataSet = reader.ReadDataSet(); + + std::cout << "Get information out of data" << std::endl; + vtkm::cont::CellSetExplicit<> cellSet; + dataSet.GetCellSet().AsCellSet(cellSet); + + auto pointCoordinates = dataSet.GetCoordinateSystem().GetData(); + + std::cout << "Run algorithm" << std::endl; + vtkm::cont::Invoker invoker; + vtkm::cont::ArrayHandle angleSums; + invoker(SumOfAngles{}, cellSet, cellSet, pointCoordinates, angleSums); + + std::cout << "Add field to data set" << std::endl; + dataSet.AddPointField("angle-sum", angleSums); + + std::cout << "Write result" << std::endl; + vtkm::io::VTKDataSetWriter writer("cow-curvature.vtk"); + writer.WriteDataSet(dataSet); +} + +} // anonymous namespace + +int GuideExampleSumOfAngles(int argc, char* argv[]) +{ + return vtkm::cont::testing::Testing::Run(TrySumOfAngles, argc, argv); +} diff --git a/docs/users-guide/examples/GuideExampleTriangleQuality.cxx b/docs/users-guide/examples/GuideExampleTriangleQuality.cxx new file mode 100644 index 000000000..cd7f2c088 --- /dev/null +++ b/docs/users-guide/examples/GuideExampleTriangleQuality.cxx @@ -0,0 +1,439 @@ +//============================================================================= +// +// Copyright (c) Kitware, Inc. +// All rights reserved. +// See LICENSE.txt for details. +// +// This software is distributed WITHOUT ANY WARRANTY; without even +// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR +// PURPOSE. See the above copyright notice for more information. +// +//============================================================================= + +#include +#include +#include + +#include + +#include +#include +#include + +//// +//// BEGIN-EXAMPLE TriangleQualityWholeArray +//// +namespace detail +{ + +static const vtkm::Id TRIANGLE_QUALITY_TABLE_DIMENSION = 8; +static const vtkm::Id TRIANGLE_QUALITY_TABLE_SIZE = + TRIANGLE_QUALITY_TABLE_DIMENSION * TRIANGLE_QUALITY_TABLE_DIMENSION; + +VTKM_CONT +vtkm::cont::ArrayHandle GetTriangleQualityTable() +{ + // Use these precomputed values for the array. A real application would + // probably use a larger array, but we are keeping it small for demonstration + // purposes. + static vtkm::Float32 triangleQualityBuffer[TRIANGLE_QUALITY_TABLE_SIZE] = { + 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0.24431f, + 0, 0, 0, 0, 0, 0, 0.43298f, 0.47059f, + 0, 0, 0, 0, 0, 0.54217f, 0.65923f, 0.66408f, + 0, 0, 0, 0, 0.57972f, 0.75425f, 0.82154f, 0.81536f, + 0, 0, 0, 0.54217f, 0.75425f, 0.87460f, 0.92567f, 0.92071f, + 0, 0, 0.43298f, 0.65923f, 0.82154f, 0.92567f, 0.97664f, 0.98100f, + 0, 0.24431f, 0.47059f, 0.66408f, 0.81536f, 0.92071f, 0.98100f, 1 + }; + + return vtkm::cont::make_ArrayHandle( + triangleQualityBuffer, TRIANGLE_QUALITY_TABLE_SIZE, vtkm::CopyFlag::Off); +} + +template +VTKM_EXEC_CONT vtkm::Vec TriangleEdgeLengths(const vtkm::Vec& point1, + const vtkm::Vec& point2, + const vtkm::Vec& point3) +{ + return vtkm::make_Vec(vtkm::Magnitude(point1 - point2), + vtkm::Magnitude(point2 - point3), + vtkm::Magnitude(point3 - point1)); +} + +VTKM_SUPPRESS_EXEC_WARNINGS +template +VTKM_EXEC_CONT static vtkm::Float32 LookupTriangleQuality( + const PortalType& triangleQualityPortal, + const vtkm::Vec& point1, + const vtkm::Vec& point2, + const vtkm::Vec& point3) +{ + vtkm::Vec edgeLengths = TriangleEdgeLengths(point1, point2, point3); + + // To reduce the size of the table, we just store the quality of triangles + // with the longest edge of size 1. The table is 2D indexed by the length + // of the other two edges. Thus, to use the table we have to identify the + // longest edge and scale appropriately. + T smallEdge1 = vtkm::Min(edgeLengths[0], edgeLengths[1]); + T tmpEdge = vtkm::Max(edgeLengths[0], edgeLengths[1]); + T smallEdge2 = vtkm::Min(edgeLengths[2], tmpEdge); + T largeEdge = vtkm::Max(edgeLengths[2], tmpEdge); + + smallEdge1 /= largeEdge; + smallEdge2 /= largeEdge; + + // Find index into array. + vtkm::Id index1 = static_cast( + vtkm::Floor(smallEdge1 * (TRIANGLE_QUALITY_TABLE_DIMENSION - 1) + 0.5)); + vtkm::Id index2 = static_cast( + vtkm::Floor(smallEdge2 * (TRIANGLE_QUALITY_TABLE_DIMENSION - 1) + 0.5)); + vtkm::Id totalIndex = index1 + index2 * TRIANGLE_QUALITY_TABLE_DIMENSION; + + return triangleQualityPortal.Get(totalIndex); +} + +} // namespace detail + +struct TriangleQualityWorklet : vtkm::worklet::WorkletVisitCellsWithPoints +{ + using ControlSignature = void(CellSetIn cells, + FieldInPoint pointCoordinates, + WholeArrayIn triangleQualityTable, + FieldOutCell triangleQuality); + using ExecutionSignature = _4(CellShape, _2, _3); + using InputDomain = _1; + + template + VTKM_EXEC vtkm::Float32 operator()( + CellShape shape, + const PointCoordinatesType& pointCoordinates, + const TriangleQualityTablePortalType& triangleQualityTable) const + { + if (shape.Id != vtkm::CELL_SHAPE_TRIANGLE) + { + this->RaiseError("Only triangles are supported for triangle quality."); + return vtkm::Nan32(); + } + else + { + return detail::LookupTriangleQuality(triangleQualityTable, + pointCoordinates[0], + pointCoordinates[1], + pointCoordinates[2]); + } + } +}; + +// +// Later in the associated Filter class... +// + +//// PAUSE-EXAMPLE +struct DemoTriangleQuality +{ + vtkm::cont::Invoker Invoke; + + template + VTKM_CONT vtkm::cont::ArrayHandle Run( + const vtkm::cont::DataSet inputDataSet, + const CoordinatesType& inputPointCoordinatesField) + { + //// RESUME-EXAMPLE + vtkm::cont::ArrayHandle triangleQualityTable = + detail::GetTriangleQualityTable(); + + vtkm::cont::ArrayHandle triangleQualities; + + this->Invoke(TriangleQualityWorklet{}, + inputDataSet.GetCellSet(), + inputPointCoordinatesField, + triangleQualityTable, + triangleQualities); + //// + //// END-EXAMPLE TriangleQualityWholeArray + //// + + return triangleQualities; + } +}; + +//// +//// BEGIN-EXAMPLE TriangleQualityExecObject +//// +template +class TriangleQualityTableExecutionObject +{ + using TableArrayType = vtkm::cont::ArrayHandle; + using TablePortalType = typename TableArrayType::ReadPortalType; + TablePortalType TablePortal; + +public: + VTKM_CONT + TriangleQualityTableExecutionObject(const TablePortalType& tablePortal) + : TablePortal(tablePortal) + { + } + + template + VTKM_EXEC vtkm::Float32 GetQuality(const vtkm::Vec& point1, + const vtkm::Vec& point2, + const vtkm::Vec& point3) const + { + return detail::LookupTriangleQuality(this->TablePortal, point1, point2, point3); + } +}; + +class TriangleQualityTable : public vtkm::cont::ExecutionObjectBase +{ +public: + template + VTKM_CONT TriangleQualityTableExecutionObject PrepareForExecution( + Device, + vtkm::cont::Token& token) const + { + return TriangleQualityTableExecutionObject( + detail::GetTriangleQualityTable().PrepareForInput(Device{}, token)); + } +}; + +struct TriangleQualityWorklet2 : vtkm::worklet::WorkletVisitCellsWithPoints +{ + using ControlSignature = void(CellSetIn cells, + FieldInPoint pointCoordinates, + ExecObject triangleQualityTable, + FieldOutCell triangleQuality); + using ExecutionSignature = _4(CellShape, _2, _3); + using InputDomain = _1; + + template + VTKM_EXEC vtkm::Float32 operator()( + CellShape shape, + const PointCoordinatesType& pointCoordinates, + const TriangleQualityTableType& triangleQualityTable) const + { + if (shape.Id != vtkm::CELL_SHAPE_TRIANGLE) + { + this->RaiseError("Only triangles are supported for triangle quality."); + return vtkm::Nan32(); + } + else + { + return triangleQualityTable.GetQuality( + pointCoordinates[0], pointCoordinates[1], pointCoordinates[2]); + } + } +}; + +// +// Later in the associated Filter class... +// + +//// PAUSE-EXAMPLE +struct DemoTriangleQuality2 +{ + vtkm::cont::Invoker Invoke; + + template + VTKM_CONT vtkm::cont::ArrayHandle Run( + const vtkm::cont::DataSet inputDataSet, + const CoordinatesType& inputPointCoordinatesField) + { + //// RESUME-EXAMPLE + TriangleQualityTable triangleQualityTable; + + vtkm::cont::ArrayHandle triangleQualities; + + this->Invoke(TriangleQualityWorklet2{}, + inputDataSet.GetCellSet(), + inputPointCoordinatesField, + triangleQualityTable, + triangleQualities); + //// + //// END-EXAMPLE TriangleQualityExecObject + //// + + return triangleQualities; + } +}; + +#include +#include + +#include + +#include + +namespace TriangleQualityNamespace +{ + +template +VTKM_EXEC T TriangleQuality(const vtkm::Vec& edgeLengths) +{ + // Heron's formula for triangle area. + T semiperimeter = (edgeLengths[0] + edgeLengths[1] + edgeLengths[2]) / 2; + T areaSquared = (semiperimeter * (semiperimeter - edgeLengths[0]) * + (semiperimeter - edgeLengths[1]) * (semiperimeter - edgeLengths[2])); + + if (areaSquared < 0) + { + // If the edge lengths do not make a valid triangle (i.e. the sum of the + // two smaller lengths is smaller than the larger length), then Heron's + // formula gives an imaginary number. If that happens, just return a + // quality of 0 for the degenerate triangle. + return 0; + } + T area = vtkm::Sqrt(areaSquared); + + // Formula for triangle quality. + return 4 * area * vtkm::Sqrt(T(3)) / vtkm::MagnitudeSquared(edgeLengths); +} + +struct ComputeTriangleQualityValues : vtkm::worklet::WorkletMapField +{ + using ControlSignature = void(FieldIn, FieldOut); + using ExecutionSignature = _2(_1); + + template + VTKM_EXEC T operator()(const vtkm::Vec& edgeLengths) const + { + return TriangleQuality(edgeLengths); + } +}; + +VTKM_CONT +vtkm::cont::ArrayHandle BuildTriangleQualityTable() +{ + // Repurpose uniform point coordinates to compute triange edge lengths. + vtkm::cont::ArrayHandleUniformPointCoordinates edgeLengths( + vtkm::Id3(detail::TRIANGLE_QUALITY_TABLE_DIMENSION, + detail::TRIANGLE_QUALITY_TABLE_DIMENSION, + 1), + vtkm::Vec3f(0, 0, 1), + vtkm::Vec3f(1.0f / (detail::TRIANGLE_QUALITY_TABLE_DIMENSION - 1), + 1.0f / (detail::TRIANGLE_QUALITY_TABLE_DIMENSION - 1), + 1.0f)); + + vtkm::cont::ArrayHandle triQualityArray; + + vtkm::cont::Invoker invoke; + invoke(ComputeTriangleQualityValues{}, edgeLengths, triQualityArray); + + return triQualityArray; +} + +template +VTKM_CONT void PrintTriangleQualityTable(const PortalType& portal) +{ + for (vtkm::Id index = 0; index < portal.GetNumberOfValues(); index++) + { + if (index % detail::TRIANGLE_QUALITY_TABLE_DIMENSION == 0) + { + std::cout << std::endl; + } + std::cout << portal.Get(index) << ", "; + } + std::cout << std::endl << std::endl; +} + +VTKM_CONT +vtkm::cont::DataSet BuildDataSet() +{ + static const vtkm::Id NUM_ROWS = 5; + + vtkm::cont::DataSetBuilderExplicitIterative dataSetBuilder; + dataSetBuilder.Begin(); + + for (vtkm::Id row = 0; row < NUM_ROWS; row++) + { + dataSetBuilder.AddPoint(0, static_cast(row * row), 0); + dataSetBuilder.AddPoint(1, static_cast(row * row), 0); + } + + for (vtkm::Id row = 0; row < NUM_ROWS - 1; row++) + { + vtkm::Id firstPoint = 2 * row; + + dataSetBuilder.AddCell(vtkm::CELL_SHAPE_TRIANGLE); + dataSetBuilder.AddCellPoint(firstPoint + 0); + dataSetBuilder.AddCellPoint(firstPoint + 1); + dataSetBuilder.AddCellPoint(firstPoint + 2); + + dataSetBuilder.AddCell(vtkm::CELL_SHAPE_TRIANGLE); + dataSetBuilder.AddCellPoint(firstPoint + 1); + dataSetBuilder.AddCellPoint(firstPoint + 3); + dataSetBuilder.AddCellPoint(firstPoint + 2); + } + + return dataSetBuilder.Create(); +} + +VTKM_CONT +void CheckQualityArray(vtkm::cont::ArrayHandle qualities) +{ + vtkm::cont::printSummary_ArrayHandle(qualities, std::cout); + std::cout << std::endl; + + vtkm::cont::ArrayHandle::ReadPortalType qualityPortal = + qualities.ReadPortal(); + + // Pairwise triangles should have the same quality. + for (vtkm::Id pairIndex = 0; pairIndex < qualities.GetNumberOfValues() / 2; + pairIndex++) + { + vtkm::Float32 q1 = qualityPortal.Get(2 * pairIndex); + vtkm::Float32 q2 = qualityPortal.Get(2 * pairIndex + 1); + VTKM_TEST_ASSERT(test_equal(q1, q2), "Isometric triangles have different quality."); + } + + // Triangle qualities should be monotonically decreasing. + vtkm::Float32 lastQuality = 1; + for (vtkm::Id triIndex = 0; triIndex < qualities.GetNumberOfValues(); triIndex++) + { + vtkm::Float32 quality = qualityPortal.Get(triIndex); + VTKM_TEST_ASSERT(test_equal(quality, lastQuality) || (quality <= lastQuality), + "Triangle quality not monotonically decreasing."); + lastQuality = quality; + } + + // The first quality should definitely be better than the last. + vtkm::Float32 firstQuality = qualityPortal.Get(0); + VTKM_TEST_ASSERT(firstQuality > lastQuality, "First quality not better than last."); +} + +VTKM_CONT +void TestTriangleQuality() +{ + std::cout << "Building triangle quality array." << std::endl; + vtkm::cont::ArrayHandle triQualityTable = BuildTriangleQualityTable(); + VTKM_TEST_ASSERT(triQualityTable.GetNumberOfValues() == + detail::TRIANGLE_QUALITY_TABLE_DIMENSION * + detail::TRIANGLE_QUALITY_TABLE_DIMENSION, + "Bad size for triangle quality array."); + PrintTriangleQualityTable(triQualityTable.ReadPortal()); + + std::cout << "Creating a data set." << std::endl; + vtkm::cont::DataSet dataSet = BuildDataSet(); + + std::cout << "Getting triangle quality using whole array argument." << std::endl; + vtkm::cont::ArrayHandle qualities = + DemoTriangleQuality().Run(dataSet, dataSet.GetCoordinateSystem().GetData()); + CheckQualityArray(qualities); + + std::cout << "Getting triangle quality using execution object argument." << std::endl; + qualities = + DemoTriangleQuality2().Run(dataSet, dataSet.GetCoordinateSystem().GetData()); + CheckQualityArray(qualities); +} + +} // namespace TriangleQualityNamespace + +int GuideExampleTriangleQuality(int argc, char* argv[]) +{ + return vtkm::cont::testing::Testing::Run( + TriangleQualityNamespace::TestTriangleQuality, argc, argv); +} diff --git a/docs/users-guide/globals.rst b/docs/users-guide/globals.rst new file mode 100644 index 000000000..a8d2315f2 --- /dev/null +++ b/docs/users-guide/globals.rst @@ -0,0 +1,167 @@ +============================== +Global Arrays and Topology +============================== + +.. index:: worklet; creating + +When writing an algorithm in |VTKm| by creating a worklet, the data each instance of the worklet has access to is intentionally limited. +This allows |VTKm| to provide safety from race conditions and other parallel programming difficulties. +However, there are times when the complexity of an algorithm requires all threads to have shared global access to a global structure. +This chapter describes worklet tags that can be used to pass data globally to all instances of a worklet. + + +------------------------------ +Whole Arrays +------------------------------ + +.. index:: + single: whole array + single: worklet; whole array + single: control signature; whole array + +A *whole array* argument to a worklet allows you to pass in a :class:`vtkm::cont::ArrayHandle`. +All instances of the worklet will have access to all the data in the :class:`vtkm::cont::ArrayHandle`. + +.. commonerrors:: + The |VTKm| worklet invoking mechanism performs many safety checks to prevent race conditions across concurrently running worklets. + Using a whole array within a worklet circumvents this guarantee of safety, so be careful when using whole arrays, especially when writing to whole arrays. + +A whole array is declared by adding a :class:`WholeArrayIn`, a :class:`WholeArrayInOut`, or a :class:`WholeArrayOut` to the \controlsignature of a worklet. +The corresponding argument to the :class:`vtkm::cont::Invoker` should be an :class:`vtkm::cont::ArrayHandle`. +The :class:`vtkm::cont::ArrayHandle` must already be allocated in all cases, including when using :class:`WholeArrayOut`. +When the data are passed to the operator of the worklet, it is passed as an array portal object. +(Array portals are discussed in :secref:`basic-array-handles:Array Portals`.) +This means that the worklet can access any entry in the array with :func:`ArrayPortal::Get` and/or :func:`ArrayPortal::Set` methods. + +We have already seen a demonstration of using a whole array in :numref:`ex:RandomArrayAccess` in :chapref:`worklet-types:Worklet Types` to perform a simple array copy. +Here we will construct a more thorough example of building functionality that requires random array access. + +Let's say we want to measure the quality of triangles in a mesh. +A common method for doing this is using the equation + +.. math:: + + q = \frac{4a\sqrt{3}}{h_1^2 + h_2^2 + h_3^2} + +where :math:`a` is the area of the triangle and :math:`h_1`, :math:`h_2`, and :math:`h_3` are the lengths of the sides. +We can easily compute this in a cell to point map, but what if we want to speed up the computations by reducing precision? +After all, we probably only care if the triangle is good, reasonable, or bad. +So instead, let's build a lookup table and then retrieve the triangle quality from that lookup table based on its sides. + +The following example demonstrates creating such a table lookup in an array and using a worklet argument tagged with :class:`WholeArrayIn` to make it accessible. + +.. load-example:: TriangleQualityWholeArray + :file: GuideExampleTriangleQuality.cxx + :caption: Using :class:`WholeArrayIn` to access a lookup table in a worklet. + + +------------------------------ +Atomic Arrays +------------------------------ + +.. index:: + single: atomic array + single: worklet; atomic array + sintle: control signature; atomic array + +One of the problems with writing to whole arrays is that it is difficult to coordinate the access to an array from multiple threads. +If multiple threads are going to write to a common index of an array, then you will probably need to use an *atomic array*. + +An atomic array allows random access into an array of data, similar to a whole array. +However, the operations on the values in the atomic array allow you to perform an operation that modifies its value that is guaranteed complete without being interrupted and potentially corrupted. + +.. commonerrors:: + Due to limitations in available atomic operations, atomic arrays can currently only contain :type:`vtkm::Int32` or :type:`vtkm::Int64` values. + + +To use an array as an atomic array, first add the :class:`AtomicArrayInOut` tag to the worklet's ``ControlSignature``. +The corresponding argument to the :class:`vtkm::cont::Invoker` should be an :class:`vtkm::cont::ArrayHandle`, which must already be allocated and initialized with values. + +When the data are passed to the operator of the worklet, it is passed in a \vtkmexec{AtomicArrayExecutionObject} structure. + +.. doxygenclass:: vtkm::exec::AtomicArrayExecutionObject + :members: + +.. commonerrors:: + Atomic arrays help resolve hazards in parallel algorithms, but they come at a cost. + Atomic operations are more costly than non-thread-safe ones, and they can slow a parallel program immensely if used incorrectly. + +.. index:: histogram + +The following example uses an atomic array to count the bins in a histogram. +It does this by making the array of histogram bins an atomic array and then using an atomic add. +Note that this is not the fastest way to create a histogram. +We gave an implementation in :secref:`worklet-types:Reduce by Key` that is generally faster (unless your histogram happens to be very sparse). +|VTKm| also comes with a histogram worklet that uses a similar approach. + +.. load-example:: SimpleHistogram + :file: GuideExampleSimpleHistogram.cxx + :caption: Using :class:`AtomicArrayInOut` to count histogram bins in a worklet. + + +------------------------------ +Whole Cell Sets +------------------------------ + +.. index:: + single: whole cell set + single: cell set; whole + single: worklet; whole cell set + single: control signature; whole cell set + +:secref:`worklet-types:Topology Map` describes how to make a topology map filter that performs an operation on cell sets. +The worklet has access to a single cell element (such as point or cell) and its immediate connections. +But there are cases when you need more general queries on a topology. +For example, you might need more detailed information than the topology map gives or you might need to trace connections from one cell to the next. +To do this |VTKm| allows you to provide a *whole cell set* argument to a worklet that provides random access to the entire topology. + +A whole cell set is declared by adding a :class:`WholeCellSetIn` to the worklet's ``ControlSignature``. +The corresponding argument to the :class:`vtkm::cont::Invoker` should be a :class:`vtkm::cont::CellSet` subclass or an :class:`vtkm::cont::UnknownCellSet` (both of which are described in :secref:`dataset:Cell Sets`). + +The :class:`WholeCellSetIn` is templated and takes two arguments: the "visit" topology type and the "incident" topology type, respectively. +These template arguments must be one of the topology element tags, but for convenience you can use :class:`Point` and :class:`Cell` in lieu of :class:`vtkm::TopologyElementTagPoint` and :class:`vtkm::TopologyElementTagCell`, respectively. +The "visit" and "incident" topology types define which topological elements can be queried (visited) and which incident elements are returned. +The semantics of the "visit" and "incident" topology is the same as that for the general topology maps described in :secref:`worklet-types:Topology Map`. +You can look up an element of the "visit" topology by index and then get all of the "incident" elements from it. + +For example, a ``WholeCellSetIn`` allows you to find all the points that are incident on each cell (as well as querying the cell shape). Likewise, a ``WholeCellSetIn`` allows you to find all the cells that are incident on each point. +The default parameters of :class:`WholeCellSetIn` are visiting cells with incident points. +That is, ``WholeCellSetIn<>`` is equivalent to ``WholeCellSetIn``. + +When the cell set is passed to the operator of the worklet, it is passed in a special connectivity object. +The actual object type depends on the cell set, but :class:`vtkm::exec::ConnectivityExplicit` and are two common examples :class:`vtkm::exec::ConnectivityStructured`. + +.. doxygenclass:: vtkm::exec::ConnectivityExplicit + :members: +.. doxygenclass:: vtkm::exec::ConnectivityStructured + :members: + +All these connectivity objects share a common interface. +In particular, the share the types ``CellShapeTag`` and ``IndicesType``. +They also share the methods ``GetNumberOfElements()``, ``GetCellShape()``, ``GetNumberOfIndices()``, and ``GetIndices()``. + +|VTKm| comes with several functions to work with the shape and index information returned from these connectivity objects. +Most of these methods are documented in :chapref:`working-with-cells:Working with Cells`. + +Let us use the whole cell set feature to help us determine the "flatness" of a polygonal mesh. +We will do this by summing up all the angles incident on each on each point. +That is, for each point, we will find each incident polygon, then find the part of that polygon using the given point, then computing the angle at that point, and then summing for all such angles. +So, for example, in the mesh fragment shown in :numref:`fig:PointIncidentAngles` one of the angles attached to the middle point is labeled :math:`\theta_{j}`. + +.. figure:: images/PointIncidentAngles.png + :width: 25% + :name: fig:PointIncidentAngles + + The angles incident around a point in a mesh. + +We want a worklet to compute :math:`\sum_{j} \theta` for all such attached angles. +This measure is related (but not the same as) the curvature of the surface. +A flat surface will have a sum of :math:`2\pi`. +Convex and concave surfaces have a value less than :math:`2\pi`, and saddle surfaces have a value greater than :math:`2\pi`. + +To do this, we create a visit points with cells worklet (:secref:`worklet-types:Visit Points with Cells`) that visits every point and gives the index of every incident cell. +The worklet then uses a whole cell set to inspect each incident cell to measure the attached angle and sum them together. + +.. load-example:: SumOfAngles + :file: GuideExampleSumOfAngles.cxx + :caption: Using :class:`WholeCellSetIn` to sum the angles around each point. diff --git a/docs/users-guide/images/PointIncidentAngles.pdf b/docs/users-guide/images/PointIncidentAngles.pdf new file mode 100644 index 000000000..16db9e44b --- /dev/null +++ b/docs/users-guide/images/PointIncidentAngles.pdf @@ -0,0 +1,3 @@ +version https://git-lfs.github.com/spec/v1 +oid sha256:54ac8dce07317525954319e2d525a1d39d521e8f073ada1590945bb92e35f18c +size 834766 diff --git a/docs/users-guide/images/PointIncidentAngles.png b/docs/users-guide/images/PointIncidentAngles.png new file mode 100644 index 000000000..1fabde773 --- /dev/null +++ b/docs/users-guide/images/PointIncidentAngles.png @@ -0,0 +1,3 @@ +version https://git-lfs.github.com/spec/v1 +oid sha256:a259ca34e02fc5133895279dfc712a5786fe3f49a899c9c17308687606daeff9 +size 19872 diff --git a/docs/users-guide/part-advanced.rst b/docs/users-guide/part-advanced.rst index 0b9439be3..16660f76e 100644 --- a/docs/users-guide/part-advanced.rst +++ b/docs/users-guide/part-advanced.rst @@ -14,3 +14,4 @@ Advanced Development working-with-cells.rst memory-layout.rst fancy-array-handles.rst + globals.rst diff --git a/docs/users-guide/worklet-types.rst b/docs/users-guide/worklet-types.rst index 4cc07b318..f773a5734 100644 --- a/docs/users-guide/worklet-types.rst +++ b/docs/users-guide/worklet-types.rst @@ -89,9 +89,9 @@ Field maps most commonly perform basic calculator arithmetic, as demonstrated in Although simple, the :class:`vtkm::worklet::WorkletMapField` worklet type can be used (and abused) as a general parallel-for/scheduling mechanism. In particular, the :class:`WorkIndex` execution signature tag can be used to get a unique index, the ``WholeArray*`` tags can be used to get random access to arrays, and the :class:`ExecObject` control signature tag can be used to pass execution objects directly to the worklet. -Whole arrays and execution objects are talked about in more detail in Chapters \ref{chap:Globals} and \ref{chap:ExecutionObjects}, respectively, in more detail, but here is a simple example that uses the random access of :class`WholeArrayOut` to make a worklet that copies an array in reverse order. +Whole arrays and execution objects are talked about in more detail in :chapref:`globals:Global Arrays and Topology` and Chapter \ref{chap:ExecutionObjects}, respectively, in more detail, but here is a simple example that uses the random access of :class:`WholeArrayOut` to make a worklet that copies an array in reverse order. -.. todo:: Fix references to globals and execution object chapters above. +.. todo:: Fix reference to execution object chapter above. .. load-example:: RandomArrayAccess :file: GuideExampleUseWorkletMapField.cxx diff --git a/vtkm/exec/AtomicArrayExecutionObject.h b/vtkm/exec/AtomicArrayExecutionObject.h index d36355447..7075b2e10 100644 --- a/vtkm/exec/AtomicArrayExecutionObject.h +++ b/vtkm/exec/AtomicArrayExecutionObject.h @@ -77,6 +77,13 @@ struct ArithType }; } +/// An object passed to a worklet when accessing an atomic array. +/// This object is created for the worklet when a `ControlSignature` argument is +/// `AtomicArrayInOut`. +/// +/// `AtomicArrayExecutionObject` behaves similar to a normal `ArrayPortal`. +/// It has similar `Get()` and `Set()` methods, but it has additional features +/// that allow atomic access as well as more methods for atomic operations. template class AtomicArrayExecutionObject { @@ -105,15 +112,16 @@ public: "GetIteratorBegin()."); } + /// @brief Retrieve the number of values in the atomic array. VTKM_SUPPRESS_EXEC_WARNINGS VTKM_EXEC vtkm::Id GetNumberOfValues() const { return this->NumberOfValues; } - /// \brief Perform an atomic load of the indexed element with acquire memory + /// @brief Perform an atomic load of the indexed element with acquire memory /// ordering. - /// \param index The index of the element to load. - /// \param order The memory ordering to use for the load operation. - /// \return The value of the atomic array at \a index. + /// @param index The index of the element to load. + /// @param order The memory ordering to use for the load operation. + /// @return The value of the atomic array at \a index. /// VTKM_SUPPRESS_EXEC_WARNINGS VTKM_EXEC @@ -127,13 +135,13 @@ public: return static_cast(vtkm::AtomicLoad(reinterpret_cast(this->Data + index), order)); } - /// \brief Peform an atomic addition with sequentially consistent memory + /// @brief Peform an atomic addition with sequentially consistent memory /// ordering. - /// \param index The index of the array element that will be added to. - /// \param value The addend of the atomic add operation. - /// \param order The memory ordering to use for the add operation. - /// \return The original value of the element at \a index (before addition). - /// \warning Overflow behavior from this operation is undefined. + /// @param index The index of the array element that will be added to. + /// @param value The addend of the atomic add operation. + /// @param order The memory ordering to use for the add operation. + /// @return The original value of the element at \a index (before addition). + /// @warning Overflow behavior from this operation is undefined. /// VTKM_SUPPRESS_EXEC_WARNINGS VTKM_EXEC @@ -153,13 +161,13 @@ public: reinterpret_cast(this->Data + index), static_cast(value), order)); } - /// \brief Peform an atomic store to memory while enforcing, at minimum, "release" + /// @brief Peform an atomic store to memory while enforcing, at minimum, "release" /// memory ordering. - /// \param index The index of the array element that will be added to. - /// \param value The value to write for the atomic store operation. - /// \param order The memory ordering to use for the store operation. - /// \warning Using something like: - /// ``` + /// @param index The index of the array element that will be added to. + /// @param value The value to write for the atomic store operation. + /// @param order The memory ordering to use for the store operation. + /// @warning Using something like: + /// ```cpp /// Set(index, Get(index)+N) /// ``` /// Should not be done as it is not thread safe, instead you should use @@ -183,16 +191,16 @@ public: reinterpret_cast(this->Data + index), static_cast(value), order); } - /// \brief Perform an atomic compare and exchange operation with sequentially consistent + /// @brief Perform an atomic compare and exchange operation with sequentially consistent /// memory ordering. - /// \param index The index of the array element that will be atomically + /// @param index The index of the array element that will be atomically /// modified. - /// \param oldValue A pointer to the expected value of the indexed element. - /// \param newValue The value to replace the indexed element with. - /// \param order The memory ordering to use for the compare and exchange operation. - /// \return If the operation is successful, \a true is returned. Otherwise, - /// \a oldValue is replaced with the current value of the indexed element, - /// the element is not modified, and \a false is returned. In either case, \a oldValue + /// @param oldValue A pointer to the expected value of the indexed element. + /// @param newValue The value to replace the indexed element with. + /// @param order The memory ordering to use for the compare and exchange operation. + /// @return If the operation is successful, `true` is returned. Otherwise, + /// `oldValue` is replaced with the current value of the indexed element, + /// the element is not modified, and `false` is returned. In either case, `oldValue` /// becomes the value that was originally in the indexed element. /// /// This operation is typically used in a loop. For example usage, @@ -209,18 +217,18 @@ public: /// } while (!atomicArray.CompareExchange(idx, ¤t, newVal)); /// ``` /// - /// The while condition here updates \a newVal what the proper multiplication + /// The `while` condition here updates `newVal` what the proper multiplication /// is given the expected current value. It then compares this to the /// value in the array. If the values match, the operation was successful and the - /// loop exits. If the values do not match, the value at \a idx was changed + /// loop exits. If the values do not match, the value at `idx` was changed /// by another thread since the initial Get, and the compare-exchange operation failed -- /// the target element was not modified by the compare-exchange call. If this happens, the - /// loop body re-executes using the new value of \a current and tries again until + /// loop body re-executes using the new value of `current` and tries again until /// it succeeds. /// /// Note that for demonstration purposes, the previous code is unnecessarily verbose. /// We can express the same atomic operation more succinctly with just two lines where - /// \a newVal is just computed in place. + /// `newVal` is just computed in place. /// /// ```cpp /// vtkm::Int32 current = atomicArray.Get(idx); // Load the current value at idx diff --git a/vtkm/exec/ConnectivityExplicit.h b/vtkm/exec/ConnectivityExplicit.h index 3cdd71895..4e88604d9 100644 --- a/vtkm/exec/ConnectivityExplicit.h +++ b/vtkm/exec/ConnectivityExplicit.h @@ -20,6 +20,11 @@ namespace vtkm namespace exec { +/// @brief A class holding information about topology connections. +/// +/// An object of `ConnectivityExplicit` is provided to a worklet when the +/// `ControlSignature` argument is `WholeCellSetIn` and the `vtkm::cont::CellSet` +/// provided is a `vtkm::cont::CellSetExplicit`. template class ConnectivityExplicit { @@ -37,22 +42,40 @@ public: { } + /// @brief Provides the number of elements in the topology. + /// + /// This number of elements is associated with the "visit" type of topology element, + /// which is the first template argument to `WholeCellSetIn`. The number of elements + /// defines the valid indices for the other methods of this class. VTKM_EXEC SchedulingRangeType GetNumberOfElements() const { return this->Shapes.GetNumberOfValues(); } + /// @brief The tag representing the cell shape of the visited elements. + /// + /// The tag type is allways `vtkm::CellShapeTagGeneric` and its id is filled with the + /// identifier for the appropriate shape. using CellShapeTag = vtkm::CellShapeTagGeneric; + /// @brief Returns a tagfor the cell shape associated with the element at the given index. + /// + /// The tag type is allways `vtkm::CellShapeTagGeneric` and its id is filled with the + /// identifier for the appropriate shape. VTKM_EXEC CellShapeTag GetCellShape(vtkm::Id index) const { return CellShapeTag(this->Shapes.Get(index)); } + /// Given the index of a visited element, returns the number of incident elements + /// touching it. VTKM_EXEC vtkm::IdComponent GetNumberOfIndices(vtkm::Id index) const { return static_cast(this->Offsets.Get(index + 1) - this->Offsets.Get(index)); } + /// @brief Type of variable that lists of incident indices will be put into. using IndicesType = vtkm::VecFromPortal; + /// Provides the indices of all elements incident to the visit element of the provided + /// index. /// Returns a Vec-like object containing the indices for the given index. /// The object returned is not an actual array, but rather an object that /// loads the indices lazily out of the connectivity array. This prevents diff --git a/vtkm/exec/ConnectivityStructured.h b/vtkm/exec/ConnectivityStructured.h index fd3fce942..2734812c2 100644 --- a/vtkm/exec/ConnectivityStructured.h +++ b/vtkm/exec/ConnectivityStructured.h @@ -21,6 +21,11 @@ namespace vtkm namespace exec { +/// @brief A class holding information about topology connections. +/// +/// An object of `ConnectivityStructured` is provided to a worklet when the +/// `ControlSignature` argument is `WholeCellSetIn` and the `vtkm::cont::CellSet` +/// provided is a `vtkm::cont::CellSetStructured`. template class ConnectivityStructured { @@ -57,43 +62,72 @@ public: ConnectivityStructured& operator=(ConnectivityStructured&& src) = default; + /// @brief Provides the number of elements in the topology. + /// + /// This number of elements is associated with the "visit" type of topology element, + /// which is the first template argument to `WholeCellSetIn`. The number of elements + /// defines the valid indices for the other methods of this class. VTKM_EXEC vtkm::Id GetNumberOfElements() const { return Helper::GetNumberOfElements(this->Internals); } + /// @brief The tag representing the cell shape of the visited elements. + /// + /// If the "visit" element is cells, then the returned tag is `vtkm::CellShapeTagHexahedron` + /// for a 3D structured grid, `vtkm::CellShapeTagQuad` for a 2D structured grid, or + /// `vtkm::CellShapeLine` for a 1D structured grid. using CellShapeTag = typename Helper::CellShapeTag; + + /// @brief Returns a tag for the cell shape associated with the element at the given index. + /// + /// If the "visit" element is cells, then the returned tag is `vtkm::CellShapeTagHexahedron` + /// for a 3D structured grid, `vtkm::CellShapeTagQuad` for a 2D structured grid, or + /// `vtkm::CellShapeLine` for a 1D structured grid. VTKM_EXEC CellShapeTag GetCellShape(vtkm::Id) const { return CellShapeTag(); } + /// Given the index of a visited element, returns the number of incident elements + /// touching it. template VTKM_EXEC vtkm::IdComponent GetNumberOfIndices(const IndexType& index) const { return Helper::GetNumberOfIndices(this->Internals, index); } + /// @brief Type of variable that lists of incident indices will be put into. using IndicesType = typename Helper::IndicesType; + /// Provides the indices of all elements incident to the visit element of the provided + /// index. template VTKM_EXEC IndicesType GetIndices(const IndexType& index) const { return Helper::GetIndices(this->Internals, index); } + /// Convenience method that converts a flat, 1D index to the visited elements to a `vtkm::Vec` + /// containing the logical indices in the grid. VTKM_EXEC_CONT SchedulingRangeType FlatToLogicalVisitIndex(vtkm::Id flatVisitIndex) const { return Helper::FlatToLogicalVisitIndex(this->Internals, flatVisitIndex); } + /// Convenience method that converts a flat, 1D index to the incident elements to a `vtkm::Vec` + /// containing the logical indices in the grid. VTKM_EXEC_CONT SchedulingRangeType FlatToLogicalIncidentIndex(vtkm::Id flatIncidentIndex) const { return Helper::FlatToLogicalIncidentIndex(this->Internals, flatIncidentIndex); } + /// Convenience method that converts logical indices in a `vtkm::Vec` of a visited element + /// to a flat, 1D index. VTKM_EXEC_CONT vtkm::Id LogicalToFlatVisitIndex( const SchedulingRangeType& logicalVisitIndex) const { return Helper::LogicalToFlatVisitIndex(this->Internals, logicalVisitIndex); } + /// Convenience method that converts logical indices in a `vtkm::Vec` of an incident element + /// to a flat, 1D index. VTKM_EXEC_CONT vtkm::Id LogicalToFlatIncidentIndex( const SchedulingRangeType& logicalIncidentIndex) const { @@ -128,12 +162,14 @@ public: return this->LogicalToFlatVisitIndex(logicalToIndex); } + /// Return the dimensions of the points in the cell set. VTKM_EXEC_CONT vtkm::Vec GetPointDimensions() const { return this->Internals.GetPointDimensions(); } + /// Return the dimensions of the points in the cell set. VTKM_EXEC_CONT vtkm::Vec GetCellDimensions() const {