Merge topic 'add/hyperstructstats'

3273f1c48 Update to AddPointField of dataset
eef30c232 Merge branch 'master' into add/hyperstructstats
1014dbf04 Fix BOV reader support in ContourTreeApp
e908826c1 Temporary disabled BOV reader in ContourTreeApp pending #517
8a27ff510 Merge branch 'master' into add/hyperstructstats
43e470de7 Fix g++ conversion warning in ContourTreeApp
a9a4736c3 Add BOV reader support and fix TBB settings in contour tree app
38da8dc43 Fix TBB check bug in ContourTreeApp
...

Acked-by: Kitware Robot <kwrobot@kitware.com>
Acked-by: Kenneth Moreland <kmorel@sandia.gov>
Merge-request: !2002
This commit is contained in:
Kenneth Moreland 2020-05-28 17:24:58 +00:00 committed by Kitware Robot
commit da6fc5e378
3 changed files with 226 additions and 99 deletions

@ -59,7 +59,7 @@ find_package(VTKm REQUIRED QUIET)
# Serial
####################################
add_executable(ContourTree_Augmented ContourTreeApp.cxx)
target_link_libraries(ContourTree_Augmented vtkm_filter)
target_link_libraries(ContourTree_Augmented vtkm_filter vtkm_io)
vtkm_add_target_information(ContourTree_Augmented
DROP_UNUSED_SYMBOLS MODIFY_CUDA_FLAGS
DEVICE_SOURCES ContourTreeApp.cxx)
@ -78,7 +78,7 @@ endif()
####################################
if (VTKm_ENABLE_MPI)
add_executable(ContourTree_Augmented_MPI ContourTreeApp.cxx)
target_link_libraries(ContourTree_Augmented_MPI vtkm_filter)
target_link_libraries(ContourTree_Augmented_MPI vtkm_filter vtkm_io)
vtkm_add_target_information(ContourTree_Augmented_MPI
MODIFY_CUDA_FLAGS
DEVICE_SOURCES ContourTreeApp.cxx)

@ -64,9 +64,13 @@
#include <vtkm/cont/ArrayHandle.h>
#include <vtkm/cont/DataSet.h>
#include <vtkm/cont/DataSetBuilderUniform.h>
#include <vtkm/cont/DataSetFieldAdd.h>
#include <vtkm/cont/DeviceAdapterTag.h>
#include <vtkm/cont/Initialize.h>
#include <vtkm/cont/RuntimeDeviceTracker.h>
#include <vtkm/cont/Timer.h>
#include <vtkm/io/BOVDataSetReader.h>
#include <vtkm/filter/ContourTreeUniformAugmented.h>
#include <vtkm/worklet/contourtree_augmented/PrintVectors.h>
#include <vtkm/worklet/contourtree_augmented/ProcessContourTree.h>
@ -228,18 +232,21 @@ int main(int argc, char* argv[])
int numThreads = tbb::task_scheduler_init::default_num_threads();
if (parser.hasOption("--numThreads"))
{
// Print warning about mismatch between the --numThreads and -d/--device opton
VTKM_LOG_IF_S(vtkm::cont::LogLevel::Warn,
device != vtkm::cont::DeviceAdapterTagTBB{},
"WARNING: Mismatch between --numThreads and -d/--device option."
"numThreads option requires the use of TBB as device. "
"Ignoring the numThread option.");
bool deviceIsTBB = (device.GetName() == "TBB");
// Set the number of threads to be used for TBB
if (device == vtkm::cont::DeviceAdapterTagTBB{})
if (deviceIsTBB)
{
numThreads = std::stoi(parser.getOption("--numThreads"));
tbb::task_scheduler_init schedulerInit(numThreads);
}
// Print warning about mismatch between the --numThreads and -d/--device option
else
{
VTKM_LOG_S(vtkm::cont::LogLevel::Warn,
"WARNING: Mismatch between --numThreads and -d/--device option."
"numThreads option requires the use of TBB as device. "
"Ignoring the numThread option.");
}
}
#endif
@ -432,84 +439,99 @@ int main(int argc, char* argv[])
///////////////////////////////////////////////
// Read the input data
///////////////////////////////////////////////
std::ifstream inFile(filename);
if (inFile.bad())
return 0;
// Read the dimensions of the mesh, i.e,. number of elementes in x, y, and z
std::vector<std::size_t> dims;
std::string line;
getline(inFile, line);
std::istringstream linestream(line);
std::size_t dimVertices;
while (linestream >> dimVertices)
{
dims.push_back(dimVertices);
}
// Compute the number of vertices, i.e., xdim * ydim * zdim
unsigned short nDims = static_cast<unsigned short>(dims.size());
std::size_t numVertices = static_cast<std::size_t>(
std::accumulate(dims.begin(), dims.end(), std::size_t(1), std::multiplies<std::size_t>()));
// Print the mesh metadata
if (rank == 0)
{
VTKM_LOG_S(vtkm::cont::LogLevel::Info,
std::endl
<< " ---------------- Input Mesh Properties --------------"
<< std::endl
<< " Number of dimensions: "
<< nDims
<< std::endl
<< " Number of mesh vertices: "
<< numVertices
<< std::endl);
}
// Check for fatal input errors
// Check the the number of dimensiosn is either 2D or 3D
bool invalidNumDimensions = (nDims < 2 || nDims > 3);
// Check if marching cubes is enabled for non 3D data
bool invalidMCOption = (useMarchingCubes && nDims != 3);
// Log any errors if found on rank 0
VTKM_LOG_IF_S(vtkm::cont::LogLevel::Error,
invalidNumDimensions && (rank == 0),
"The input mesh is " << nDims << "D. "
"The input data must be either 2D or 3D.");
VTKM_LOG_IF_S(
vtkm::cont::LogLevel::Error,
invalidMCOption && (rank == 0),
"The input mesh is " << nDims << "D. "
<< "Contour tree using marching cubes is only supported for 3D data.");
// If we found any errors in the setttings than finalize MPI and exit the execution
if (invalidNumDimensions || invalidMCOption)
vtkm::Float64 dataReadTime = 0;
vtkm::Float64 buildDatasetTime = 0;
std::vector<vtkm::Float32>::size_type nDims = 0;
vtkm::cont::DataSet inDataSet;
std::vector<ValueType> values;
std::vector<vtkm::Id> dims;
if (filename.compare(filename.length() - 3, 3, "bov") == 0)
{
std::cout << "Reading BOV file" << std::endl;
vtkm::io::BOVDataSetReader reader(filename);
inDataSet = reader.ReadDataSet();
nDims = 3;
currTime = totalTime.GetElapsedTime();
dataReadTime = currTime - prevTime;
prevTime = currTime;
#ifdef WITH_MPI
// Copy the data into the values array so we can construct a multiblock dataset
// TODO All we should need to do to implement BOV support is to copy the values
// in the values vector and copy the dimensions in the dims vector
vtkm::Id nRows, nCols, nSlices;
vtkm::filter::GetRowsColsSlices temp;
temp(inDataSet.GetCellSet(), nRows, nCols, nSlices);
dims[0] = nRows;
dims[1] = nCols;
dims[2] = nSlices;
auto tempField = inDataSet.GetField("values").GetData();
values.resize(static_cast<std::size_t>(tempField.GetNumberOfValues()));
auto tempFieldHandle = tempField.AsVirtual<ValueType>().ReadPortal();
for (vtkm::Id i = 0; i < tempField.GetNumberOfValues(); i++)
{
values[static_cast<std::size_t>(i)] = static_cast<ValueType>(tempFieldHandle.Get(i));
}
VTKM_LOG_S(vtkm::cont::LogLevel::Error,
"BOV reader not yet support in MPI mode by this example");
MPI_Finalize();
#endif
return EXIT_SUCCESS;
#endif
}
// Read data
std::vector<ValueType> values(numVertices);
for (std::size_t vertex = 0; vertex < numVertices; ++vertex)
else // Read ASCII data input
{
inFile >> values[vertex];
}
std::cout << "Reading ASCII file" << std::endl;
std::ifstream inFile(filename);
if (inFile.bad())
return 0;
// finish reading the data
inFile.close();
// Read the dimensions of the mesh, i.e,. number of elementes in x, y, and z
std::string line;
getline(inFile, line);
std::istringstream linestream(line);
vtkm::Id dimVertices;
while (linestream >> dimVertices)
{
dims.push_back(dimVertices);
}
currTime = totalTime.GetElapsedTime();
vtkm::Float64 dataReadTime = currTime - prevTime;
prevTime = currTime;
// Compute the number of vertices, i.e., xdim * ydim * zdim
nDims = static_cast<unsigned short>(dims.size());
std::size_t numVertices = static_cast<std::size_t>(
std::accumulate(dims.begin(), dims.end(), std::size_t(1), std::multiplies<std::size_t>()));
vtkm::cont::DataSetBuilderUniform dsb;
#ifndef WITH_MPI // construct regular, single-block VTK-M input dataset
vtkm::cont::DataSet inDataSet; // Single block dataset
{
// Check for fatal input errors
// Check the the number of dimensiosn is either 2D or 3D
bool invalidNumDimensions = (nDims < 2 || nDims > 3);
// Log any errors if found on rank 0
VTKM_LOG_IF_S(vtkm::cont::LogLevel::Error,
invalidNumDimensions && (rank == 0),
"The input mesh is " << nDims << "D. The input data must be either 2D or 3D.");
// If we found any errors in the setttings than finalize MPI and exit the execution
if (invalidNumDimensions)
{
#ifdef WITH_MPI
MPI_Finalize();
#endif
return EXIT_SUCCESS;
}
// Read data
values.resize(numVertices);
for (std::size_t vertex = 0; vertex < numVertices; ++vertex)
{
inFile >> values[vertex];
}
// finish reading the data
inFile.close();
currTime = totalTime.GetElapsedTime();
dataReadTime = currTime - prevTime;
prevTime = currTime;
#ifndef WITH_MPI // We only need the inDataSet if are not using MPI otherwise we'll constructe a multi-block dataset
// build the input dataset
vtkm::cont::DataSetBuilderUniform dsb;
// 2D data
if (nDims == 2)
{
@ -528,9 +550,41 @@ int main(int argc, char* argv[])
inDataSet = dsb.Create(vdims);
}
inDataSet.AddPointField("values", values);
#endif
} // END ASCII Read
// Print the mesh metadata
if (rank == 0)
{
VTKM_LOG_S(vtkm::cont::LogLevel::Info,
std::endl
<< " ---------------- Input Mesh Properties --------------"
<< std::endl
<< " Number of dimensions: "
<< nDims);
}
#else // Create a multi-block dataset for multi-block DIY-paralle processing
vtkm::cont::PartitionedDataSet inDataSet; // Partitioned variant of the input dataset
// Check if marching cubes is enabled for non 3D data
bool invalidMCOption = (useMarchingCubes && nDims != 3);
VTKM_LOG_IF_S(
vtkm::cont::LogLevel::Error,
invalidMCOption && (rank == 0),
"The input mesh is " << nDims << "D. "
<< "Contour tree using marching cubes is only supported for 3D data.");
// If we found any errors in the setttings than finalize MPI and exit the execution
if (invalidMCOption)
{
#ifdef WITH_MPI
MPI_Finalize();
#endif
return EXIT_SUCCESS;
}
#ifndef WITH_MPI // construct regular, single-block VTK-M input dataset
vtkm::cont::DataSet useDataSet = inDataSet; // Single block dataset
#else // Create a multi-block dataset for multi-block DIY-paralle processing
vtkm::cont::PartitionedDataSet useDataSet; // Partitioned variant of the input dataset
vtkm::Id3 blocksPerDim =
nDims == 3 ? vtkm::Id3(1, 1, numBlocks) : vtkm::Id3(1, numBlocks, 1); // Decompose the data into
vtkm::Id3 globalSize = nDims == 3 ? vtkm::Id3(static_cast<vtkm::Id>(dims[0]),
@ -584,6 +638,7 @@ int main(int argc, char* argv[])
}
vtkm::Id currBlockSize = (vtkm::Id)((blockEnd - blockStart) / blockSliceSize);
vtkm::cont::DataSetBuilderUniform dsb;
vtkm::cont::DataSet ds;
// 2D data
@ -625,13 +680,13 @@ int main(int argc, char* argv[])
(values.begin() + blockEnd));
ds.AddPointField("values", subValues);
inDataSet.AppendPartition(ds);
useDataSet.AppendPartition(ds);
}
}
#endif // WITH_MPI construct input dataset
#endif // WITH_MPI construct input dataset
currTime = totalTime.GetElapsedTime();
vtkm::Float64 buildDatasetTime = currTime - prevTime;
buildDatasetTime = currTime - prevTime;
prevTime = currTime;
// Convert the mesh of values into contour tree, pairs of vertex ids
@ -645,7 +700,7 @@ int main(int argc, char* argv[])
// Execute the contour tree analysis. NOTE: If MPI is used the result will be
// a vtkm::cont::PartitionedDataSet instead of a vtkm::cont::DataSet
auto result = filter.Execute(inDataSet);
auto result = filter.Execute(useDataSet);
currTime = totalTime.GetElapsedTime();
vtkm::Float64 computeContourTreeTime = currTime - prevTime;
@ -710,7 +765,7 @@ int main(int argc, char* argv[])
bool dataFieldIsSorted = true;
#else
vtkm::cont::ArrayHandle<ValueType> dataField;
inDataSet.GetField(0).GetData().CopyTo(dataField);
useDataSet.GetField(0).GetData().CopyTo(dataField);
bool dataFieldIsSorted = false;
#endif
@ -916,6 +971,12 @@ int main(int argc, char* argv[])
<< ": "
<< ct.Hyperarcs.GetNumberOfValues()
<< std::endl);
// Print hyperstructure statistics
VTKM_LOG_S(vtkm::cont::LogLevel::Info,
std::endl
<< ct.PrintHyperStructureStatistics(false)
<< std::endl);
// Flush ouput streams just to make sure everything has been logged (in particular when using MPI)
std::cout << std::flush;
std::cerr << std::flush;

@ -57,6 +57,8 @@
#include <algorithm>
#include <iomanip>
#include <iostream>
#include <sstream>
#include <string>
// local includes
#include <vtkm/worklet/contourtree_augmented/PrintVectors.h>
@ -161,15 +163,14 @@ public:
inline void Init(vtkm::Id dataSize);
// debug routine
inline void DebugPrint(const char* message, const char* fileName, long lineNum);
inline void DebugPrint(const char* message, const char* fileName, long lineNum) const;
// print contents
inline void PrintContent() const;
// print routines
//void PrintDotHyperStructure();
inline void PrintDotSuperStructure();
//void PrintDotRegularStructure();
inline void PrintDotSuperStructure() const;
inline std::string PrintHyperStructureStatistics(bool print = true) const;
}; // class ContourTree
@ -217,7 +218,7 @@ inline void ContourTree::PrintContent() const
PrintIndices("Augmentarcs", this->Augmentarcs);
}
void ContourTree::DebugPrint(const char* message, const char* fileName, long lineNum)
void ContourTree::DebugPrint(const char* message, const char* fileName, long lineNum) const
{ // DebugPrint()
#ifdef DEBUG_PRINT
std::cout << "---------------------------" << std::endl;
@ -237,13 +238,7 @@ void ContourTree::DebugPrint(const char* message, const char* fileName, long lin
#endif
} // DebugPrint()
// print routines
// void ContourTree::PrintDotHyperStructure()
// { // PrintDotHyperStructure()
// } // PrintDotHyperStructure()
void ContourTree::PrintDotSuperStructure()
void ContourTree::PrintDotSuperStructure() const
{ // PrintDotSuperStructure()
// print the header information
printf("digraph G\n\t{\n");
@ -330,9 +325,80 @@ void ContourTree::PrintDotSuperStructure()
printf("\t}\n");
} // PrintDotSuperStructure()
// void ContourTree::PrintDotRegularStructure()
// { // PrintDotRegularStructure()
// } // PrintDotRegularStructure()
std::string ContourTree::PrintHyperStructureStatistics(bool print) const
{ // PrintHyperStructureStatistics()
// arrays for collecting statistics
std::vector<vtkm::Id> minPath;
std::vector<vtkm::Id> maxPath;
std::vector<vtkm::Id> supernodeCount;
std::vector<vtkm::Id> hypernodeCount;
auto whenTransferredPortal = this->WhenTransferred.ReadPortal();
auto hypernodesPortal = this->Hypernodes.ReadPortal();
// set an initial iteration number to negative to get it started
long whichIteration = -1;
// loop through the hypernodes
for (vtkm::Id hypernode = 0; hypernode < this->Hypernodes.GetNumberOfValues(); hypernode++)
{ // per hypernode
// retrieve corresponding supernode ID
vtkm::Id supernodeID = hypernodesPortal.Get(hypernode);
// and the iteration of transfer
vtkm::Id iterationNo = MaskedIndex(whenTransferredPortal.Get(supernodeID));
// if it doesn't match, we've hit a boundary
if (whichIteration != iterationNo)
{ // new iteration
// initialise the next iteration
// this one is larger than the maximum possible to force minimum
minPath.push_back(static_cast<vtkm::Id>(this->Supernodes.GetNumberOfValues() + 1));
maxPath.push_back(0);
supernodeCount.push_back(0);
hypernodeCount.push_back(0);
// and increment the iteration ID
whichIteration++;
} // new iteration
// now compute the new path length - default to off the end
vtkm::Id pathLength = static_cast<vtkm::Id>(this->Supernodes.GetNumberOfValues() - supernodeID);
// for all except the last, take the next one
if (hypernode != this->Hypernodes.GetNumberOfValues() - 1)
{
pathLength = hypernodesPortal.Get(hypernode + 1) - supernodeID;
}
// update the statistics
if (pathLength < minPath[static_cast<std::size_t>(whichIteration)])
{
minPath[static_cast<std::size_t>(whichIteration)] = pathLength;
}
if (pathLength > maxPath[static_cast<std::size_t>(whichIteration)])
{
maxPath[static_cast<std::size_t>(whichIteration)] = pathLength;
}
supernodeCount[static_cast<std::size_t>(whichIteration)] += pathLength;
hypernodeCount[static_cast<std::size_t>(whichIteration)]++;
} // per hypernode
// now print out the statistics
std::stringstream resultString;
for (std::size_t iteration = 0; iteration < minPath.size(); iteration++)
{ // per iteration
double averagePath = static_cast<double>(supernodeCount[iteration]) /
static_cast<double>(hypernodeCount[iteration]);
resultString << "Iteration: " << iteration << " Hyper: " << hypernodeCount[iteration]
<< " Super: " << supernodeCount[iteration] << " Min: " << minPath[iteration]
<< " Avg: " << averagePath << " Max: " << maxPath[iteration] << std::endl;
} // per iteration
resultString << "Total Hypernodes: " << this->Hypernodes.GetNumberOfValues()
<< " Supernodes: " << this->Supernodes.GetNumberOfValues() << std::endl;
if (print)
{
std::cout << resultString.str() << std::endl;
}
return resultString.str();
} // PrintHyperStructureStatistics()
} // namespace contourtree_augmented
} // worklet