mirror of
https://gitlab.kitware.com/vtk/vtk-m
synced 2024-10-05 01:49:02 +00:00
Particle advection tests with file.
This commit is contained in:
parent
da57431d9c
commit
755137a822
3
data/data/rectilinear/fishtank.vtk
Normal file
3
data/data/rectilinear/fishtank.vtk
Normal file
@ -0,0 +1,3 @@
|
||||
version https://git-lfs.github.com/spec/v1
|
||||
oid sha256:ef3dfd79f0c8d18780d0749014d71c0226134041283d33de0bcd994e343dd421
|
||||
size 2001070
|
3
data/data/rectilinear/fishtank_double_ascii.vtk
Normal file
3
data/data/rectilinear/fishtank_double_ascii.vtk
Normal file
@ -0,0 +1,3 @@
|
||||
version https://git-lfs.github.com/spec/v1
|
||||
oid sha256:2bb3d36ea5ecef5e7ef1057d0dddebbc590424915083091ead3dac2928000524
|
||||
size 2904465
|
3
data/data/rectilinear/fishtank_double_big_endian.vtk
Normal file
3
data/data/rectilinear/fishtank_double_big_endian.vtk
Normal file
@ -0,0 +1,3 @@
|
||||
version https://git-lfs.github.com/spec/v1
|
||||
oid sha256:bffad7dae3dd6ef018ad7a9e109464ced0f3b9bc15cf1fb5d555f6d0d00b621f
|
||||
size 3001624
|
3
data/data/rectilinear/fusion.vtk
Normal file
3
data/data/rectilinear/fusion.vtk
Normal file
@ -0,0 +1,3 @@
|
||||
version https://git-lfs.github.com/spec/v1
|
||||
oid sha256:2cbdf56fd5445ddc5b6bc05507b8825fb8d74fe1ccce894bde03e5ff2ecf5fb6
|
||||
size 525141
|
@ -64,7 +64,6 @@
|
||||
#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>
|
||||
@ -459,17 +458,17 @@ int main(int argc, char* argv[])
|
||||
// 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::worklet::contourtree_augmented::GetRowsColsSlices temp;
|
||||
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()));
|
||||
values.resize(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));
|
||||
values[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");
|
||||
|
@ -12,6 +12,8 @@
|
||||
#include <vtkm/cont/testing/Testing.h>
|
||||
#include <vtkm/filter/Pathline.h>
|
||||
#include <vtkm/filter/Streamline.h>
|
||||
#include <vtkm/io/VTKDataSetReader.h>
|
||||
#include <vtkm/io/VTKDataSetWriter.h>
|
||||
|
||||
namespace
|
||||
{
|
||||
@ -102,10 +104,99 @@ void TestPathline()
|
||||
VTKM_TEST_ASSERT(dcells.GetNumberOfCells() == 3, "Wrong number of cells");
|
||||
}
|
||||
|
||||
void TestStreamlineFile(const std::string& fname,
|
||||
const std::vector<vtkm::Vec3f>& pts,
|
||||
vtkm::FloatDefault stepSize,
|
||||
vtkm::Id maxSteps,
|
||||
const std::vector<vtkm::Vec3f>& endPts)
|
||||
{
|
||||
vtkm::io::VTKDataSetReader reader(fname);
|
||||
vtkm::cont::DataSet ds;
|
||||
try
|
||||
{
|
||||
ds = reader.ReadDataSet();
|
||||
}
|
||||
catch (vtkm::io::ErrorIO& e)
|
||||
{
|
||||
std::string message("Error reading: ");
|
||||
message += fname;
|
||||
message += ", ";
|
||||
message += e.GetMessage();
|
||||
|
||||
VTKM_TEST_FAIL(message.c_str());
|
||||
}
|
||||
vtkm::Id numPoints = static_cast<vtkm::Id>(pts.size());
|
||||
|
||||
std::vector<vtkm::Particle> seeds;
|
||||
for (vtkm::Id i = 0; i < numPoints; i++)
|
||||
seeds.push_back(vtkm::Particle(pts[static_cast<std::size_t>(i)], i));
|
||||
auto seedArray = vtkm::cont::make_ArrayHandle(seeds);
|
||||
|
||||
vtkm::filter::Streamline streamline;
|
||||
streamline.SetStepSize(stepSize);
|
||||
streamline.SetNumberOfSteps(maxSteps);
|
||||
streamline.SetSeeds(seedArray);
|
||||
|
||||
streamline.SetActiveField("vec");
|
||||
auto output = streamline.Execute(ds);
|
||||
|
||||
auto coords = output.GetCoordinateSystem().GetData();
|
||||
vtkm::cont::DynamicCellSet dcells = output.GetCellSet();
|
||||
VTKM_TEST_ASSERT(dcells.GetNumberOfCells() == numPoints, "Wrong number of cells");
|
||||
VTKM_TEST_ASSERT(dcells.IsType<vtkm::cont::CellSetExplicit<>>(), "Wrong cell type");
|
||||
|
||||
auto cells = dcells.Cast<vtkm::cont::CellSetExplicit<>>();
|
||||
auto cPortal = coords.ReadPortal();
|
||||
const vtkm::FloatDefault eps = 1e-3;
|
||||
|
||||
for (vtkm::Id i = 0; i < numPoints; i++)
|
||||
{
|
||||
vtkm::Id numPts = cells.GetNumberOfPointsInCell(i);
|
||||
std::vector<vtkm::Id> ids(static_cast<std::size_t>(numPts));
|
||||
cells.GetCellPointIds(i, ids.data());
|
||||
|
||||
vtkm::Vec3f e = endPts[static_cast<std::size_t>(i)];
|
||||
vtkm::Vec3f pt = cPortal.Get(ids[ids.size() - 1]);
|
||||
if (vtkm::Magnitude(pt - e) > eps)
|
||||
{
|
||||
VTKM_LOG_S(vtkm::cont::LogLevel::Error,
|
||||
"Expected magnitude of error <=" << eps << ", but got error "
|
||||
<< vtkm::Magnitude(pt - e));
|
||||
VTKM_TEST_ASSERT(false, "Streamline end point is wrong.");
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void TestStreamlineFilters()
|
||||
{
|
||||
TestStreamline();
|
||||
TestPathline();
|
||||
|
||||
std::string basePath = vtkm::cont::testing::Testing::GetTestDataBasePath();
|
||||
|
||||
//Fusion test.
|
||||
std::vector<vtkm::Vec3f> fusionPts, fusionEndPts;
|
||||
fusionPts.push_back(vtkm::Vec3f(0.8f, 0.6f, 0.6f));
|
||||
fusionPts.push_back(vtkm::Vec3f(0.8f, 0.8f, 0.6f));
|
||||
fusionPts.push_back(vtkm::Vec3f(0.8f, 0.8f, 0.3f));
|
||||
fusionEndPts.push_back(vtkm::Vec3f(0.5335789918f, 0.87112802267f, 0.6723330020f));
|
||||
fusionEndPts.push_back(vtkm::Vec3f(0.5601879954f, 0.91389900446f, 0.43989110522f));
|
||||
fusionEndPts.push_back(vtkm::Vec3f(0.7004770041f, 0.63193398714f, 0.64524400234f));
|
||||
vtkm::FloatDefault fusionStep = 0.005f;
|
||||
std::string fusionFile = basePath + "/rectilinear/fusion.vtk";
|
||||
TestStreamlineFile(fusionFile, fusionPts, fusionStep, 1000, fusionEndPts);
|
||||
|
||||
//Fishtank test.
|
||||
std::vector<vtkm::Vec3f> fishPts, fishEndPts;
|
||||
fishPts.push_back(vtkm::Vec3f(0.75f, 0.5f, 0.01f));
|
||||
fishPts.push_back(vtkm::Vec3f(0.4f, 0.2f, 0.7f));
|
||||
fishPts.push_back(vtkm::Vec3f(0.5f, 0.3f, 0.8f));
|
||||
fishEndPts.push_back(vtkm::Vec3f(0.7734669447f, 0.4870159328f, 0.8979591727f));
|
||||
fishEndPts.push_back(vtkm::Vec3f(0.7257543206f, 0.1277695596f, 0.7468645573f));
|
||||
fishEndPts.push_back(vtkm::Vec3f(0.8347796798f, 0.1276152730f, 0.4985143244f));
|
||||
vtkm::FloatDefault fishStep = 0.001f;
|
||||
std::string fishFile = basePath + "/rectilinear/fishtank.vtk";
|
||||
TestStreamlineFile(fishFile, fishPts, fishStep, 100, fishEndPts);
|
||||
}
|
||||
}
|
||||
|
||||
|
@ -46,27 +46,34 @@ void VTKRectilinearGridReader::Read()
|
||||
this->DataFile->Stream >> dim[0] >> dim[1] >> dim[2] >> std::ws;
|
||||
|
||||
//Read the points.
|
||||
std::string dataType;
|
||||
std::string fileStorageDataType;
|
||||
std::size_t numPoints[3];
|
||||
vtkm::cont::VariantArrayHandle X, Y, Z;
|
||||
|
||||
// Always read coordinates as vtkm::FloatDefault
|
||||
std::string readDataType = vtkm::io::internal::DataTypeName<vtkm::FloatDefault>::Name();
|
||||
|
||||
this->DataFile->Stream >> tag >> numPoints[0] >> dataType >> std::ws;
|
||||
this->DataFile->Stream >> tag >> numPoints[0] >> fileStorageDataType >> std::ws;
|
||||
if (tag != "X_COORDINATES")
|
||||
throw vtkm::io::ErrorIO("X_COORDINATES tag not found");
|
||||
X = this->DoReadArrayVariant(vtkm::cont::Field::Association::ANY, readDataType, numPoints[0], 1);
|
||||
|
||||
this->DataFile->Stream >> tag >> numPoints[1] >> dataType >> std::ws;
|
||||
// In binary mode, we must read the data as they are stored in the file.
|
||||
// In text mode we can parse as FloatDefault no matter the precision of the storage.
|
||||
X = this->DoReadArrayVariant(
|
||||
vtkm::cont::Field::Association::ANY, fileStorageDataType, numPoints[0], 1);
|
||||
|
||||
|
||||
this->DataFile->Stream >> tag >> numPoints[1] >> fileStorageDataType >> std::ws;
|
||||
if (tag != "Y_COORDINATES")
|
||||
throw vtkm::io::ErrorIO("Y_COORDINATES tag not found");
|
||||
Y = this->DoReadArrayVariant(vtkm::cont::Field::Association::ANY, readDataType, numPoints[1], 1);
|
||||
|
||||
this->DataFile->Stream >> tag >> numPoints[2] >> dataType >> std::ws;
|
||||
Y = this->DoReadArrayVariant(
|
||||
vtkm::cont::Field::Association::ANY, fileStorageDataType, numPoints[1], 1);
|
||||
|
||||
this->DataFile->Stream >> tag >> numPoints[2] >> fileStorageDataType >> std::ws;
|
||||
if (tag != "Z_COORDINATES")
|
||||
throw vtkm::io::ErrorIO("Z_COORDINATES tag not found");
|
||||
Z = this->DoReadArrayVariant(vtkm::cont::Field::Association::ANY, readDataType, numPoints[2], 1);
|
||||
|
||||
Z = this->DoReadArrayVariant(
|
||||
vtkm::cont::Field::Association::ANY, fileStorageDataType, numPoints[2], 1);
|
||||
|
||||
|
||||
if (dim != vtkm::Id3(static_cast<vtkm::Id>(numPoints[0]),
|
||||
static_cast<vtkm::Id>(numPoints[1]),
|
||||
@ -78,14 +85,46 @@ void VTKRectilinearGridReader::Read()
|
||||
vtkm::cont::ArrayHandle<vtkm::FloatDefault>>
|
||||
coords;
|
||||
|
||||
// We need to store all coordinate arrays as FloatDefault.
|
||||
vtkm::cont::ArrayHandle<vtkm::FloatDefault> Xc, Yc, Zc;
|
||||
X.CopyTo(Xc);
|
||||
Y.CopyTo(Yc);
|
||||
Z.CopyTo(Zc);
|
||||
// But the VariantArrayHandle has type fileStorageDataType.
|
||||
// If the fileStorageDataType is the same as the storage type, no problem:
|
||||
if (fileStorageDataType == vtkm::io::internal::DataTypeName<vtkm::FloatDefault>::Name())
|
||||
{
|
||||
X.CopyTo(Xc);
|
||||
Y.CopyTo(Yc);
|
||||
Z.CopyTo(Zc);
|
||||
}
|
||||
else
|
||||
{
|
||||
// Two cases if the data in the file differs from FloatDefault:
|
||||
if (fileStorageDataType == "float")
|
||||
{
|
||||
vtkm::cont::ArrayHandle<float> Xcf, Ycf, Zcf;
|
||||
X.CopyTo(Xcf);
|
||||
Y.CopyTo(Ycf);
|
||||
Z.CopyTo(Zcf);
|
||||
vtkm::cont::ArrayCopy(Xcf, Xc);
|
||||
vtkm::cont::ArrayCopy(Ycf, Yc);
|
||||
vtkm::cont::ArrayCopy(Zcf, Zc);
|
||||
}
|
||||
else
|
||||
{
|
||||
vtkm::cont::ArrayHandle<double> Xcd, Ycd, Zcd;
|
||||
X.CopyTo(Xcd);
|
||||
Y.CopyTo(Ycd);
|
||||
Z.CopyTo(Zcd);
|
||||
vtkm::cont::ArrayCopy(Xcd, Xc);
|
||||
vtkm::cont::ArrayCopy(Ycd, Yc);
|
||||
vtkm::cont::ArrayCopy(Zcd, Zc);
|
||||
}
|
||||
}
|
||||
// As a postscript to this somewhat branchy code, I thought that X.CopyTo(Xc) should be able to cast to the value_type of Xc.
|
||||
// But that would break the ability to make the cast lazy, and hence if you change it you induce performance bugs.
|
||||
|
||||
coords = vtkm::cont::make_ArrayHandleCartesianProduct(Xc, Yc, Zc);
|
||||
vtkm::cont::CoordinateSystem coordSys("coordinates", coords);
|
||||
this->DataSet.AddCoordinateSystem(coordSys);
|
||||
|
||||
this->DataSet.SetCellSet(internal::CreateCellSetStructured(dim));
|
||||
|
||||
// Read points and cell attributes
|
||||
|
@ -578,6 +578,254 @@ void TestReadingStructuredGridBin()
|
||||
"Incorrect cellset type");
|
||||
}
|
||||
|
||||
void TestReadingFishTank()
|
||||
{
|
||||
std::string fishtank =
|
||||
vtkm::cont::testing::Testing::GetTestDataBasePath() + "/rectilinear/fishtank.vtk";
|
||||
vtkm::cont::DataSet ds = readVTKDataSet(fishtank.c_str());
|
||||
|
||||
// This is information you can glean by running 'strings' on fishtank.vtk:
|
||||
VTKM_TEST_ASSERT(ds.GetCellSet().IsType<vtkm::cont::CellSetStructured<3>>(),
|
||||
"Incorrect cellset type");
|
||||
VTKM_TEST_ASSERT(ds.GetNumberOfPoints() == 50 * 50 * 50, "Incorrect number of points");
|
||||
VTKM_TEST_ASSERT(ds.GetCellSet().GetNumberOfPoints() == 50 * 50 * 50,
|
||||
"Incorrect number of points (from cell set)");
|
||||
VTKM_TEST_ASSERT(ds.GetNumberOfFields() == 2, "Incorrect number of fields");
|
||||
VTKM_TEST_ASSERT(ds.HasField("vec"), "The vtk file has a field 'vec', but the dataset does not.");
|
||||
VTKM_TEST_ASSERT(ds.HasField("vec_magnitude"),
|
||||
"The vtk file has a field 'vec_magnitude', but the dataset does not.");
|
||||
|
||||
// I believe the coordinate system is implicitly given by the first element of X_COORDINATES:
|
||||
VTKM_TEST_ASSERT(ds.GetNumberOfCoordinateSystems() == 1,
|
||||
"Need one and only one coordinate system.");
|
||||
// In order to get the data from the coordinate system, I used the following workflow:
|
||||
// First, I deleted all ascii header lines just past 'X_COORDINATES 50 float'.
|
||||
// Once this is done, I can get the binary data from
|
||||
// $ od -tfF --endian=big fishtank_copy.vtk
|
||||
// The result is:
|
||||
// 0 0.020408163 ... 0.9591837 0.97959185 1
|
||||
// So monotone increasing, bound [0,1].
|
||||
const vtkm::cont::CoordinateSystem& coordinateSystem = ds.GetCoordinateSystem();
|
||||
vtkm::Vec<vtkm::Range, 3> ranges = coordinateSystem.GetRange();
|
||||
vtkm::Range xRange = ranges[0];
|
||||
VTKM_TEST_ASSERT(xRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(xRange.Max == 1);
|
||||
// Do the same past 'Y_COORDINATES 50 float'.
|
||||
// You get exactly the same as the x data.
|
||||
vtkm::Range yRange = ranges[1];
|
||||
VTKM_TEST_ASSERT(yRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(yRange.Max == 1);
|
||||
// And finally, do it past 'Z_COORDINATES 50 float':
|
||||
vtkm::Range zRange = ranges[2];
|
||||
VTKM_TEST_ASSERT(zRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(zRange.Max == 1);
|
||||
|
||||
// Now delete the text up to LOOKUP TABLE default.
|
||||
// I see:
|
||||
// 0 0 0 0 3.5267966 . . .
|
||||
// This is a vector magnitude, so all values must be >= 0.
|
||||
// A cursory glance shows that 124.95 is a large value, so we can sanity check the data with the bounds
|
||||
// [0, ~130].
|
||||
// And if we open the file in Paraview, we can observe the bounds [0, 156.905].
|
||||
const vtkm::cont::Field& vec_magnitude = ds.GetField("vec_magnitude");
|
||||
VTKM_TEST_ASSERT(vec_magnitude.GetName() == "vec_magnitude");
|
||||
VTKM_TEST_ASSERT(vec_magnitude.IsFieldPoint());
|
||||
|
||||
vtkm::Range mag_range;
|
||||
vec_magnitude.GetRange(&mag_range);
|
||||
VTKM_TEST_ASSERT(mag_range.Min == 0);
|
||||
VTKM_TEST_ASSERT(mag_range.Max <= 156.906);
|
||||
|
||||
// This info was gleaned from the Paraview Information panel:
|
||||
const vtkm::cont::Field& vec = ds.GetField("vec");
|
||||
VTKM_TEST_ASSERT(vec.GetName() == "vec");
|
||||
VTKM_TEST_ASSERT(vec.IsFieldPoint());
|
||||
// Bounds from Information panel:
|
||||
// [-65.3147, 86.267], [-88.0325, 78.7217], [-67.0969, 156.867]
|
||||
const vtkm::cont::ArrayHandle<vtkm::Range>& vecRanges = vec.GetRange();
|
||||
VTKM_TEST_ASSERT(vecRanges.GetNumberOfValues() == 3);
|
||||
auto vecRangesReadPortal = vecRanges.ReadPortal();
|
||||
|
||||
auto xVecRange = vecRangesReadPortal.Get(0);
|
||||
VTKM_TEST_ASSERT(xVecRange.Min >= -65.3148 && xVecRange.Min <= -65.3146);
|
||||
VTKM_TEST_ASSERT(xVecRange.Max >= 86.26 && xVecRange.Min <= 86.268);
|
||||
|
||||
auto yVecRange = vecRangesReadPortal.Get(1);
|
||||
VTKM_TEST_ASSERT(yVecRange.Min >= -88.0326 && yVecRange.Min <= -88.0324);
|
||||
VTKM_TEST_ASSERT(yVecRange.Max >= 78.721);
|
||||
VTKM_TEST_ASSERT(yVecRange.Max <= 78.7218);
|
||||
|
||||
auto zVecRange = vecRangesReadPortal.Get(2);
|
||||
VTKM_TEST_ASSERT(zVecRange.Min >= -67.097 && zVecRange.Min <= -67.096);
|
||||
VTKM_TEST_ASSERT(zVecRange.Max >= 156.866 && zVecRange.Max <= 156.868);
|
||||
}
|
||||
|
||||
void TestReadingDoublePrecisionFishTank()
|
||||
{
|
||||
std::string fishtank = vtkm::cont::testing::Testing::GetTestDataBasePath() +
|
||||
"/rectilinear/fishtank_double_big_endian.vtk";
|
||||
vtkm::cont::DataSet ds = readVTKDataSet(fishtank.c_str());
|
||||
|
||||
// This is information you can glean by running 'strings' on fishtank.vtk:
|
||||
VTKM_TEST_ASSERT(ds.GetCellSet().IsType<vtkm::cont::CellSetStructured<3>>(),
|
||||
"Incorrect cellset type");
|
||||
VTKM_TEST_ASSERT(ds.GetNumberOfPoints() == 50 * 50 * 50, "Incorrect number of points");
|
||||
VTKM_TEST_ASSERT(ds.GetCellSet().GetNumberOfPoints() == 50 * 50 * 50,
|
||||
"Incorrect number of points (from cell set)");
|
||||
|
||||
VTKM_TEST_ASSERT(ds.HasField("vec"), "The vtk file has a field 'vec', but the dataset does not.");
|
||||
|
||||
|
||||
VTKM_TEST_ASSERT(ds.GetNumberOfCoordinateSystems() == 1,
|
||||
"fishtank has one and only one coordinate system.");
|
||||
// See the single precision version for info:
|
||||
const vtkm::cont::CoordinateSystem& coordinateSystem = ds.GetCoordinateSystem();
|
||||
vtkm::Vec<vtkm::Range, 3> ranges = coordinateSystem.GetRange();
|
||||
vtkm::Range xRange = ranges[0];
|
||||
VTKM_TEST_ASSERT(xRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(xRange.Max == 1);
|
||||
vtkm::Range yRange = ranges[1];
|
||||
VTKM_TEST_ASSERT(yRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(yRange.Max == 1);
|
||||
vtkm::Range zRange = ranges[2];
|
||||
VTKM_TEST_ASSERT(zRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(zRange.Max == 1);
|
||||
|
||||
// This info was gleaned from the Paraview Information panel:
|
||||
const vtkm::cont::Field& vec = ds.GetField("vec");
|
||||
VTKM_TEST_ASSERT(vec.GetName() == "vec");
|
||||
VTKM_TEST_ASSERT(vec.IsFieldPoint());
|
||||
// Bounds from Information panel:
|
||||
// [-65.3147, 86.267], [-88.0325, 78.7217], [-67.0969, 156.867]
|
||||
const vtkm::cont::ArrayHandle<vtkm::Range>& vecRanges = vec.GetRange();
|
||||
VTKM_TEST_ASSERT(vecRanges.GetNumberOfValues() == 3);
|
||||
auto vecRangesReadPortal = vecRanges.ReadPortal();
|
||||
|
||||
auto xVecRange = vecRangesReadPortal.Get(0);
|
||||
VTKM_TEST_ASSERT(xVecRange.Min >= -65.3148 && xVecRange.Min <= -65.3146);
|
||||
VTKM_TEST_ASSERT(xVecRange.Max >= 86.26 && xVecRange.Min <= 86.268);
|
||||
|
||||
auto yVecRange = vecRangesReadPortal.Get(1);
|
||||
VTKM_TEST_ASSERT(yVecRange.Min >= -88.0326 && yVecRange.Min <= -88.0324);
|
||||
VTKM_TEST_ASSERT(yVecRange.Max >= 78.721);
|
||||
VTKM_TEST_ASSERT(yVecRange.Max <= 78.7218);
|
||||
|
||||
auto zVecRange = vecRangesReadPortal.Get(2);
|
||||
VTKM_TEST_ASSERT(zVecRange.Min >= -67.097 && zVecRange.Min <= -67.096);
|
||||
VTKM_TEST_ASSERT(zVecRange.Max >= 156.866 && zVecRange.Max <= 156.868);
|
||||
}
|
||||
|
||||
void TestReadingASCIIFishTank()
|
||||
{
|
||||
std::string fishtank =
|
||||
vtkm::cont::testing::Testing::GetTestDataBasePath() + "/rectilinear/fishtank_double_ascii.vtk";
|
||||
vtkm::cont::DataSet ds = readVTKDataSet(fishtank.c_str());
|
||||
VTKM_TEST_ASSERT(ds.GetCellSet().IsType<vtkm::cont::CellSetStructured<3>>(),
|
||||
"Incorrect cellset type");
|
||||
VTKM_TEST_ASSERT(ds.GetNumberOfPoints() == 50 * 50 * 50, "Incorrect number of points");
|
||||
VTKM_TEST_ASSERT(ds.GetCellSet().GetNumberOfPoints() == 50 * 50 * 50,
|
||||
"Incorrect number of points (from cell set)");
|
||||
VTKM_TEST_ASSERT(ds.HasField("vec"), "The vtk file has a field 'vec', but the dataset does not.");
|
||||
VTKM_TEST_ASSERT(ds.GetNumberOfCoordinateSystems() == 1,
|
||||
"fishtank has one and only one coordinate system.");
|
||||
const vtkm::cont::CoordinateSystem& coordinateSystem = ds.GetCoordinateSystem();
|
||||
vtkm::Vec<vtkm::Range, 3> ranges = coordinateSystem.GetRange();
|
||||
vtkm::Range xRange = ranges[0];
|
||||
VTKM_TEST_ASSERT(xRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(xRange.Max == 1);
|
||||
vtkm::Range yRange = ranges[1];
|
||||
VTKM_TEST_ASSERT(yRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(yRange.Max == 1);
|
||||
vtkm::Range zRange = ranges[2];
|
||||
VTKM_TEST_ASSERT(zRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(zRange.Max == 1);
|
||||
|
||||
const vtkm::cont::Field& vec = ds.GetField("vec");
|
||||
VTKM_TEST_ASSERT(vec.GetName() == "vec");
|
||||
VTKM_TEST_ASSERT(vec.IsFieldPoint());
|
||||
// Bounds from Paraview information panel:
|
||||
// [-65.3147, 86.267], [-88.0325, 78.7217], [-67.0969, 156.867]
|
||||
const vtkm::cont::ArrayHandle<vtkm::Range>& vecRanges = vec.GetRange();
|
||||
VTKM_TEST_ASSERT(vecRanges.GetNumberOfValues() == 3);
|
||||
auto vecRangesReadPortal = vecRanges.ReadPortal();
|
||||
auto xVecRange = vecRangesReadPortal.Get(0);
|
||||
VTKM_TEST_ASSERT(xVecRange.Min >= -65.3148 && xVecRange.Min <= -65.3146);
|
||||
VTKM_TEST_ASSERT(xVecRange.Max >= 86.26 && xVecRange.Min <= 86.268);
|
||||
|
||||
auto yVecRange = vecRangesReadPortal.Get(1);
|
||||
VTKM_TEST_ASSERT(yVecRange.Min >= -88.0326 && yVecRange.Min <= -88.0324);
|
||||
VTKM_TEST_ASSERT(yVecRange.Max >= 78.721);
|
||||
VTKM_TEST_ASSERT(yVecRange.Max <= 78.7218);
|
||||
|
||||
auto zVecRange = vecRangesReadPortal.Get(2);
|
||||
VTKM_TEST_ASSERT(zVecRange.Min >= -67.097 && zVecRange.Min <= -67.096);
|
||||
VTKM_TEST_ASSERT(zVecRange.Max >= 156.866 && zVecRange.Max <= 156.868);
|
||||
}
|
||||
|
||||
void TestReadingFusion()
|
||||
{
|
||||
std::string fusion =
|
||||
vtkm::cont::testing::Testing::GetTestDataBasePath() + "/rectilinear/fusion.vtk";
|
||||
vtkm::cont::DataSet ds = readVTKDataSet(fusion.c_str());
|
||||
|
||||
VTKM_TEST_ASSERT(ds.GetCellSet().IsType<vtkm::cont::CellSetStructured<3>>(),
|
||||
"Incorrect cellset type");
|
||||
VTKM_TEST_ASSERT(ds.GetNumberOfPoints() == 32 * 32 * 32, "Incorrect number of points");
|
||||
VTKM_TEST_ASSERT(ds.GetCellSet().GetNumberOfPoints() == 32 * 32 * 32,
|
||||
"Incorrect number of points (from cell set)");
|
||||
VTKM_TEST_ASSERT(ds.HasField("vec_magnitude"),
|
||||
"The vtk file has a field 'vec_magnitude', but the dataset does not.");
|
||||
VTKM_TEST_ASSERT(ds.HasField("vec"), "The vtk file has a field 'vec', but the dataset does not.");
|
||||
VTKM_TEST_ASSERT(ds.GetNumberOfCoordinateSystems() == 1,
|
||||
"The vtk file has a field 'vec', but the dataset does not.");
|
||||
|
||||
// Taken from Paraview + clicking Data Axes Grid:
|
||||
const vtkm::cont::CoordinateSystem& coordinateSystem = ds.GetCoordinateSystem();
|
||||
vtkm::Vec<vtkm::Range, 3> ranges = coordinateSystem.GetRange();
|
||||
vtkm::Range xRange = ranges[0];
|
||||
VTKM_TEST_ASSERT(xRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(xRange.Max == 1);
|
||||
vtkm::Range yRange = ranges[1];
|
||||
VTKM_TEST_ASSERT(yRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(yRange.Max == 1);
|
||||
vtkm::Range zRange = ranges[2];
|
||||
VTKM_TEST_ASSERT(zRange.Min == 0);
|
||||
VTKM_TEST_ASSERT(zRange.Max == 1);
|
||||
|
||||
// Paraview Information Panel of this file:
|
||||
// vec_magnitude [0, 3.73778]
|
||||
vtkm::cont::Field vec_magnitude = ds.GetField("vec_magnitude");
|
||||
VTKM_TEST_ASSERT(vec_magnitude.GetName() == "vec_magnitude");
|
||||
VTKM_TEST_ASSERT(vec_magnitude.IsFieldPoint());
|
||||
|
||||
vtkm::Range mag_range;
|
||||
vec_magnitude.GetRange(&mag_range);
|
||||
VTKM_TEST_ASSERT(mag_range.Min == 0);
|
||||
VTKM_TEST_ASSERT(mag_range.Max <= 3.73779);
|
||||
VTKM_TEST_ASSERT(mag_range.Max >= 3.73777);
|
||||
|
||||
vtkm::cont::Field vec = ds.GetField("vec");
|
||||
VTKM_TEST_ASSERT(vec.GetName() == "vec");
|
||||
VTKM_TEST_ASSERT(vec.IsFieldPoint());
|
||||
const vtkm::cont::ArrayHandle<vtkm::Range>& vecRanges = vec.GetRange();
|
||||
VTKM_TEST_ASSERT(vecRanges.GetNumberOfValues() == 3);
|
||||
auto vecRangesReadPortal = vecRanges.ReadPortal();
|
||||
|
||||
// vec float [-3.41054, 3.40824], [-3.41018, 3.41036], [-0.689022, 0.480726]
|
||||
auto xVecRange = vecRangesReadPortal.Get(0);
|
||||
VTKM_TEST_ASSERT(test_equal(xVecRange.Min, -3.41054));
|
||||
VTKM_TEST_ASSERT(test_equal(xVecRange.Max, 3.40824));
|
||||
|
||||
auto yVecRange = vecRangesReadPortal.Get(1);
|
||||
|
||||
VTKM_TEST_ASSERT(test_equal(yVecRange.Min, -3.41018));
|
||||
VTKM_TEST_ASSERT(test_equal(yVecRange.Max, 3.41036));
|
||||
|
||||
auto zVecRange = vecRangesReadPortal.Get(2);
|
||||
VTKM_TEST_ASSERT(test_equal(zVecRange.Min, -0.689022));
|
||||
VTKM_TEST_ASSERT(test_equal(zVecRange.Max, 0.480726));
|
||||
}
|
||||
|
||||
void TestReadingVTKDataSet()
|
||||
{
|
||||
std::cout << "Test reading VTK Polydata file in ASCII" << std::endl;
|
||||
@ -607,6 +855,14 @@ void TestReadingVTKDataSet()
|
||||
TestReadingStructuredGridASCII();
|
||||
std::cout << "Test reading VTK StructuredGrid file in BINARY" << std::endl;
|
||||
TestReadingStructuredGridBin();
|
||||
std::cout << "Test reading float precision fishtank" << std::endl;
|
||||
TestReadingFishTank();
|
||||
std::cout << "Test reading double precision fishtank" << std::endl;
|
||||
TestReadingDoublePrecisionFishTank();
|
||||
std::cout << "Test ASCII fishtank" << std::endl;
|
||||
TestReadingASCIIFishTank();
|
||||
std::cout << "Test reading fusion" << std::endl;
|
||||
TestReadingFusion();
|
||||
}
|
||||
|
||||
int UnitTestVTKDataSetReader(int argc, char* argv[])
|
||||
|
@ -89,7 +89,7 @@ set(unit_tests
|
||||
|
||||
vtkm_unit_tests(
|
||||
SOURCES ${unit_tests}
|
||||
LIBRARIES vtkm_source vtkm_worklet vtkm_filter
|
||||
LIBRARIES vtkm_source vtkm_worklet vtkm_filter vtkm_io
|
||||
ALL_BACKENDS
|
||||
USE_VTKM_JOB_POOL
|
||||
)
|
||||
|
@ -16,6 +16,7 @@
|
||||
#include <vtkm/cont/DataSetBuilderRectilinear.h>
|
||||
#include <vtkm/cont/DataSetBuilderUniform.h>
|
||||
#include <vtkm/cont/testing/Testing.h>
|
||||
#include <vtkm/io/VTKDataSetReader.h>
|
||||
#include <vtkm/worklet/ParticleAdvection.h>
|
||||
#include <vtkm/worklet/particleadvection/GridEvaluators.h>
|
||||
#include <vtkm/worklet/particleadvection/Integrators.h>
|
||||
@ -847,6 +848,111 @@ void TestWorkletsBasic()
|
||||
}
|
||||
}
|
||||
|
||||
template <class ResultType>
|
||||
void ValidateResult(const ResultType& res,
|
||||
vtkm::Id maxSteps,
|
||||
const std::vector<vtkm::Vec3f>& endPts)
|
||||
{
|
||||
const vtkm::FloatDefault eps = 1e-3;
|
||||
vtkm::Id numPts = static_cast<vtkm::Id>(endPts.size());
|
||||
|
||||
VTKM_TEST_ASSERT(res.Particles.GetNumberOfValues() == numPts,
|
||||
"Wrong number of points in particle advection result.");
|
||||
|
||||
auto portal = res.Particles.ReadPortal();
|
||||
bool fail = false;
|
||||
for (vtkm::Id i = 0; i < 3; i++)
|
||||
{
|
||||
vtkm::Vec3f p = portal.Get(i).Pos;
|
||||
vtkm::Vec3f e = endPts[static_cast<std::size_t>(i)];
|
||||
if (vtkm::Magnitude(p - e) > eps)
|
||||
{
|
||||
std::cout << std::setprecision(15) << "P_" << i << " p " << p << " e " << e << " diff "
|
||||
<< vtkm::Magnitude(p - e) << " eps= " << eps << std::endl;
|
||||
fail = true;
|
||||
}
|
||||
|
||||
//VTKM_TEST_ASSERT(vtkm::Magnitude(p - e) <= eps, "Particle advection point is wrong");
|
||||
VTKM_TEST_ASSERT(portal.Get(i).NumSteps == maxSteps, "Particle advection NumSteps is wrong");
|
||||
VTKM_TEST_ASSERT(portal.Get(i).Status.CheckOk(), "Particle advection Status is wrong");
|
||||
VTKM_TEST_ASSERT(portal.Get(i).Status.CheckTerminate(),
|
||||
"Particle advection particle did not terminate");
|
||||
}
|
||||
VTKM_TEST_ASSERT(fail == false, "Particle advection point is wrong");
|
||||
}
|
||||
|
||||
|
||||
void TestParticleAdvectionFile(const std::string& fname,
|
||||
const std::vector<vtkm::Vec3f>& pts,
|
||||
vtkm::FloatDefault stepSize,
|
||||
vtkm::Id maxSteps,
|
||||
const std::vector<vtkm::Vec3f>& endPts)
|
||||
{
|
||||
VTKM_LOG_S(vtkm::cont::LogLevel::Info, "Testing particle advection on file " << fname);
|
||||
vtkm::io::VTKDataSetReader reader(fname);
|
||||
vtkm::cont::DataSet ds;
|
||||
try
|
||||
{
|
||||
ds = reader.ReadDataSet();
|
||||
}
|
||||
catch (vtkm::io::ErrorIO& e)
|
||||
{
|
||||
std::string message("Error reading: ");
|
||||
message += fname;
|
||||
message += ", ";
|
||||
message += e.GetMessage();
|
||||
|
||||
VTKM_TEST_FAIL(message.c_str());
|
||||
}
|
||||
|
||||
using FieldHandle = vtkm::cont::ArrayHandle<vtkm::Vec3f_32>;
|
||||
using GridEvalType = vtkm::worklet::particleadvection::GridEvaluator<FieldHandle>;
|
||||
using RK4Type = vtkm::worklet::particleadvection::RK4Integrator<GridEvalType>;
|
||||
|
||||
VTKM_TEST_ASSERT(ds.HasField("vec"));
|
||||
vtkm::cont::Field& field = ds.GetField("vec");
|
||||
auto fieldData = field.GetData();
|
||||
|
||||
if (!fieldData.IsType<FieldHandle>())
|
||||
{
|
||||
VTKM_LOG_S(vtkm::cont::LogLevel::Error,
|
||||
"The field data is of type "
|
||||
<< vtkm::cont::TypeToString<decltype(fieldData)>()
|
||||
<< ", but we expect type vtkm::cont::ArrayHandle<vtkm::Vec3f>");
|
||||
VTKM_TEST_FAIL("No field with correct type found.");
|
||||
}
|
||||
|
||||
|
||||
FieldHandle fieldArray = fieldData.Cast<FieldHandle>();
|
||||
GridEvalType eval(ds.GetCoordinateSystem(), ds.GetCellSet(), fieldArray);
|
||||
RK4Type rk4(eval, stepSize);
|
||||
|
||||
for (int i = 0; i < 2; i++)
|
||||
{
|
||||
std::vector<vtkm::Particle> seeds;
|
||||
for (size_t j = 0; j < pts.size(); j++)
|
||||
seeds.push_back(vtkm::Particle(pts[j], static_cast<vtkm::Id>(j)));
|
||||
auto seedArray = vtkm::cont::make_ArrayHandle(seeds);
|
||||
|
||||
if (i == 0)
|
||||
{
|
||||
vtkm::worklet::ParticleAdvection pa;
|
||||
vtkm::worklet::ParticleAdvectionResult res;
|
||||
|
||||
res = pa.Run(rk4, seedArray, maxSteps);
|
||||
ValidateResult(res, maxSteps, endPts);
|
||||
}
|
||||
else if (i == 1)
|
||||
{
|
||||
vtkm::worklet::Streamline s;
|
||||
vtkm::worklet::StreamlineResult res;
|
||||
|
||||
res = s.Run(rk4, seedArray, maxSteps);
|
||||
ValidateResult(res, maxSteps, endPts);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
void TestParticleAdvection()
|
||||
{
|
||||
TestIntegrators();
|
||||
@ -854,6 +960,33 @@ void TestParticleAdvection()
|
||||
TestParticleStatus();
|
||||
TestWorkletsBasic();
|
||||
TestParticleWorkletsWithDataSetTypes();
|
||||
|
||||
std::string basePath = vtkm::cont::testing::Testing::GetTestDataBasePath();
|
||||
|
||||
//Fusion test.
|
||||
std::vector<vtkm::Vec3f> fusionPts, fusionEndPts;
|
||||
fusionPts.push_back(vtkm::Vec3f(0.8f, 0.6f, 0.6f));
|
||||
fusionPts.push_back(vtkm::Vec3f(0.8f, 0.8f, 0.6f));
|
||||
fusionPts.push_back(vtkm::Vec3f(0.8f, 0.8f, 0.3f));
|
||||
fusionEndPts.push_back(vtkm::Vec3f(0.5335789918f, 0.87112802267f, 0.6723330020f));
|
||||
fusionEndPts.push_back(vtkm::Vec3f(0.5601879954f, 0.91389900446f, 0.43989110522f));
|
||||
fusionEndPts.push_back(vtkm::Vec3f(0.7004770041f, 0.63193398714f, 0.64524400234f));
|
||||
vtkm::FloatDefault fusionStep = 0.005f;
|
||||
std::string fusionFile = basePath + "/rectilinear/fusion.vtk";
|
||||
TestParticleAdvectionFile(fusionFile, fusionPts, fusionStep, 1000, fusionEndPts);
|
||||
|
||||
//Fishtank test.
|
||||
std::vector<vtkm::Vec3f> fishPts, fishEndPts;
|
||||
fishPts.push_back(vtkm::Vec3f(0.75f, 0.5f, 0.01f));
|
||||
fishPts.push_back(vtkm::Vec3f(0.4f, 0.2f, 0.7f));
|
||||
fishPts.push_back(vtkm::Vec3f(0.5f, 0.3f, 0.8f));
|
||||
fishEndPts.push_back(vtkm::Vec3f(0.7734669447f, 0.4870159328f, 0.8979591727f));
|
||||
fishEndPts.push_back(vtkm::Vec3f(0.7257543206f, 0.1277695596f, 0.7468645573f));
|
||||
fishEndPts.push_back(vtkm::Vec3f(0.8347796798f, 0.1276152730f, 0.4985143244f));
|
||||
|
||||
vtkm::FloatDefault fishStep = 0.001f;
|
||||
std::string fishFile = basePath + "/rectilinear/fishtank.vtk";
|
||||
TestParticleAdvectionFile(fishFile, fishPts, fishStep, 100, fishEndPts);
|
||||
}
|
||||
|
||||
int UnitTestParticleAdvection(int argc, char* argv[])
|
||||
|
Loading…
Reference in New Issue
Block a user