Updating temporal advection example

- renaming long directory name to something concise
- Adding Dave's Unit test for particle advection
- Fixing really trivial issues
This commit is contained in:
ayenpure 2018-04-18 22:06:40 -07:00
parent faf8c8337a
commit 4cf9942d44
8 changed files with 92 additions and 319 deletions

@ -23,21 +23,21 @@
#path so that our examples can find VTK-m
set(CMAKE_PREFIX_PATH ${VTKm_BINARY_DIR}/${VTKm_INSTALL_CONFIG_DIR})
#add_subdirectory(clipping)
#add_subdirectory(contour_tree)
#add_subdirectory(cosmotools)
#add_subdirectory(demo)
#add_subdirectory(dynamic_dispatcher)
#add_subdirectory(game_of_life)
#add_subdirectory(hello_world)
#add_subdirectory(isosurface)
#add_subdirectory(multi_backend)
#add_subdirectory(streamline)
#add_subdirectory(tetrahedra)
#add_subdirectory(particle_advection)
add_subdirectory(temporalparticleadvection)
#if(VTKm_ENABLE_RENDERING)
# add_subdirectory(rendering)
#endif()
#add_subdirectory(unified_memory)
add_subdirectory(clipping)
add_subdirectory(contour_tree)
add_subdirectory(cosmotools)
add_subdirectory(demo)
add_subdirectory(dynamic_dispatcher)
add_subdirectory(game_of_life)
add_subdirectory(hello_world)
add_subdirectory(isosurface)
add_subdirectory(multi_backend)
add_subdirectory(streamline)
add_subdirectory(tetrahedra)
add_subdirectory(particle_advection)
add_subdirectory(temporaladvection)
if(VTKm_ENABLE_RENDERING)
add_subdirectory(rendering)
endif()
add_subdirectory(unified_memory)

@ -1,47 +0,0 @@
##=============================================================================
##
## 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.
##
## Copyright 2015 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
## Copyright 2015 UT-Battelle, LLC.
## Copyright 2015 Los Alamos National Security.
##
## Under the terms of Contract DE-NA0003525 with NTESS,
## the U.S. Government retains certain rights in this software.
## Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
## Laboratory (LANL), the U.S. Government retains certain rights in
## this software.
##
##=============================================================================
#Find the VTK-m package
find_package(VTKm REQUIRED QUIET
OPTIONAL_COMPONENTS Serial CUDA TBB Rendering
)
add_executable(Particle_Advection_SERIAL ParticleAdvection.cxx)
target_include_directories(Particle_Advection_SERIAL PRIVATE ${VTKm_INCLUDE_DIRS})
target_link_libraries(Particle_Advection_SERIAL ${VTKm_LIBRARIES})
target_compile_options(Particle_Advection_SERIAL PRIVATE ${VTKm_COMPILE_OPTIONS})
if(VTKm_TBB_FOUND)
add_executable(Particle_Advection_TBB ParticleAdvectionTBB.cxx)
target_include_directories(Particle_Advection_TBB PRIVATE ${VTKm_INCLUDE_DIRS})
target_link_libraries(Particle_Advection_TBB ${VTKm_LIBRARIES})
target_compile_options(Particle_Advection_TBB PRIVATE ${VTKm_COMPILE_OPTIONS})
endif()
if(VTKm_CUDA_FOUND)
# Cuda compiles do not respect target_include_directories
cuda_include_directories(${VTKm_INCLUDE_DIRS})
cuda_add_executable(Particle_Advection_CUDA ParticleAdvection.cu)
target_include_directories(Particle_Advection_SERIAL PRIVATE ${VTKm_INCLUDE_DIRS})
target_link_libraries(Particle_Advection_CUDA PRIVATE ${VTKm_LIBRARIES})
target_compile_options(Particle_Advection_CUDA PRIVATE ${VTKm_COMPILE_OPTIONS})
endif()

@ -1,23 +0,0 @@
//============================================================================
// 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.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#define VTKM_DEVICE_ADAPTER VTKM_DEVICE_ADAPTER_CUDA
#include "ParticleAdvection.cxx"

@ -1,176 +0,0 @@
//============================================================================
// 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.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#ifndef VTKM_DEVICE_ADAPTER
#define VTKM_DEVICE_ADAPTER VTKM_DEVICE_ADAPTER_SERIAL
#endif
#include <vtkm/cont/DataSet.h>
#include <vtkm/worklet/ParticleAdvection.h>
#include <vtkm/worklet/particleadvection/Integrators.h>
#include <vtkm/worklet/particleadvection/ParticleAdvectionWorklets.h>
#include <vtkm/worklet/particleadvection/TemporalGridEvaluators.h>
#include <vtkm/cont/Timer.h>
#include <vtkm/io/reader/BOVDataSetReader.h>
#include <vtkm/io/writer/VTKDataSetWriter.h>
// For offscreen rendering
/*#include <vtkm/rendering/Actor.h>
#include <vtkm/rendering/CanvasRayTracer.h>
#include <vtkm/rendering/MapperWireframer.h>
#include <vtkm/rendering/Scene.h>
#include <vtkm/rendering/View3D.h>
*/
#include <cstdlib>
#include <vector>
int renderAndWriteDataSet(vtkm::cont::DataSet& dataset)
{
std::cout << "Trying to render the dataset" << std::endl;
/* vtkm::rendering::CanvasRayTracer canvas;
vtkm::rendering::MapperWireframer mapper;
vtkm::rendering::Actor actor(dataset.GetCellSet(),
dataset.GetCoordinateSystem(),
dataset.GetPointField(0),
vtkm::rendering::ColorTable("temperature"));
vtkm::rendering::Scene scene;
scene.AddActor(actor);
// save images.
vtkm::rendering::View3D view(scene, mapper, canvas);
view.Initialize();
view.SetBackgroundColor(vtkm::rendering::Color(1,1,1,1));
view.SetForegroundColor(vtkm::rendering::Color(0,0,0,1));
for(int i = 0; i < 16; i++)
{
std::ostringstream filename;
filename << "clipped" << i << ".ppm";
std::cout << "Writing " << filename.str() << std::endl;
view.GetCamera().Azimuth(i*45.0);
std::cout << "shifted " << std::endl;
view.GetCamera().Elevation(i*45.0);
std::cout << "lifted " << std::endl;
view.Paint();
std::cout << "painted" << std::endl;
view.SaveAs(filename.str());
}*/
vtkm::io::writer::VTKDataSetWriter writer("vtkmwritten.vtk");
writer.WriteDataSet(dataset, static_cast<vtkm::Id>(0));
return 0;
}
void RunTest(vtkm::Id numSteps, vtkm::Float32 stepSize, vtkm::Id advectType)
{
using DeviceAdapter = VTKM_DEFAULT_DEVICE_ADAPTER_TAG;
using FieldType = vtkm::Float32;
using FieldHandle = vtkm::cont::ArrayHandle<vtkm::Vec<FieldType, 3>>;
using FieldPortalConstType =
typename FieldHandle::template ExecutionTypes<DeviceAdapter>::PortalConst;
vtkm::cont::ArrayHandle<vtkm::Vec<FieldType, 3>> fieldArray1;
vtkm::cont::ArrayHandle<vtkm::Vec<FieldType, 3>> fieldArray2;
vtkm::Id numValues;
vtkm::io::reader::BOVDataSetReader reader1("slice1.bov");
vtkm::cont::DataSet ds1 = reader1.ReadDataSet();
ds1.GetField(0).GetData().CopyTo(fieldArray1);
vtkm::io::reader::BOVDataSetReader reader2("slice2.bov");
vtkm::cont::DataSet ds2 = reader2.ReadDataSet();
ds2.GetField(0).GetData().CopyTo(fieldArray2);
using GridEvaluator =
vtkm::worklet::particleadvection::TemporalGridEvaluator<FieldPortalConstType,
FieldType,
DeviceAdapter>;
using Integrator = vtkm::worklet::particleadvection::EulerIntegrator<GridEvaluator, FieldType>;
//GridEvaluator eval(ds.GetCoordinateSystem(), ds.GetCellSet(0), fieldArray);
GridEvaluator eval(ds1.GetCoordinateSystem(),
ds1.GetCellSet(0),
fieldArray1,
0,
ds2.GetCoordinateSystem(),
ds2.GetCellSet(0),
fieldArray2,
10.0);
Integrator integrator(eval, stepSize);
std::vector<vtkm::Vec<FieldType, 3>> seeds;
vtkm::Id x = 0, y = 5, z = 0;
for (int i = 0; i <= 11; i++)
{
vtkm::Vec<FieldType, 3> point;
point[0] = x;
point[1] = y;
point[2] = z++;
seeds.push_back(point);
}
vtkm::cont::ArrayHandle<vtkm::Vec<FieldType, 3>> seedArray;
seedArray = vtkm::cont::make_ArrayHandle(seeds);
if (advectType == 0)
{
vtkm::worklet::ParticleAdvection particleAdvection;
particleAdvection.Run(integrator, seedArray, numSteps, DeviceAdapter());
}
else
{
vtkm::worklet::Streamline streamline;
vtkm::worklet::StreamlineResult<FieldType> res =
streamline.Run(integrator, seedArray, numSteps, DeviceAdapter());
vtkm::cont::DataSet outData;
vtkm::cont::CoordinateSystem outputCoords("coordinates", res.positions);
outData.AddCellSet(res.polyLines);
outData.AddCoordinateSystem(outputCoords);
renderAndWriteDataSet(outData);
}
}
int main(int argc, char** argv)
{
vtkm::Id numSteps;
vtkm::Float32 stepSize;
vtkm::Id advectionType;
if (argc < 4)
{
std::cout << "Wrong number of parameters provided" << std::endl;
exit(EXIT_FAILURE);
}
numSteps = atoi(argv[1]);
stepSize = atof(argv[2]);
advectionType = atoi(argv[3]);
RunTest(numSteps, stepSize, advectionType);
return 0;
}

@ -1,26 +0,0 @@
//============================================================================
// 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.
//
// Copyright 2014 National Technology & Engineering Solutions of Sandia, LLC (NTESS).
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-NA0003525 with NTESS,
// the U.S. Government retains certain rights in this software.
//
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
// Laboratory (LANL), the U.S. Government retains certain rights in
// this software.
//============================================================================
#define VTKM_DEVICE_ADAPTER VTKM_DEVICE_ADAPTER_TBB
#define __BUILDING_TBB_VERSION__
#include "ParticleAdvection.cxx"
#include <tbb/task_scheduler_init.h>

@ -57,7 +57,7 @@ public:
}
VTKM_EXEC_CONT
bool IsWithinTemporalBoundary(const FieldType time) const { return true; }
bool IsWithinTemporalBoundary(const FieldType vtkmNotUsed(time)) const { return true; }
VTKM_EXEC_CONT
void GetSpatialBoundary(vtkm::Vec<FieldType, 3>& dir, vtkm::Vec<FieldType, 3>& boundary) const
@ -121,7 +121,7 @@ public:
}
VTKM_EXEC_CONT
bool IsWithinTemporalBoundary(const FieldType time) const { return true; }
bool IsWithinTemporalBoundary(const FieldType vtkmNotUsed(time)) const { return true; }
VTKM_EXEC_CONT
void GetSpatialBoundary(vtkm::Vec<FieldType, 3>& dir, vtkm::Vec<FieldType, 3>& boundary) const
@ -257,7 +257,7 @@ public:
}
VTKM_EXEC_CONT
bool IsWithinTemporalBoundary(const FieldType time) const { return true; }
bool IsWithinTemporalBoundary(const FieldType vtkmNotUsed(time)) const { return true; }
VTKM_EXEC_CONT
void GetSpatialBoundary(vtkm::Vec<FieldType, 3>& dir, vtkm::Vec<FieldType, 3>& boundary) const
@ -448,7 +448,7 @@ public:
}
VTKM_EXEC_CONT
bool IsWithinTemporalBoundary(const FieldType time) const { return true; }
bool IsWithinTemporalBoundary(const FieldType vtkmNotUsed(time)) const { return true; }
VTKM_EXEC_CONT
void GetSpatialBoundary(vtkm::Vec<FieldType, 3>& dir, vtkm::Vec<FieldType, 3>& boundary) const

@ -142,6 +142,8 @@ public:
time += stepLength;
return ParticleStatus::EXITED_TEMPORAL_BOUNDARY;
}
//If the control reaches here, it is an invalid case.
return ParticleStatus::STATUS_ERROR;
}
VTKM_EXEC
@ -193,14 +195,14 @@ public:
bool status3 = this->Evaluator.Evaluate(inpos + var1 * k2, var2, k3);
bool status4 = this->Evaluator.Evaluate(inpos + stepLength * k3, var3, k4);
if (status1 & status2 & status3 & status4 == ParticleStatus::STATUS_OK)
if ((status1 & status2 & status3 & status4) == ParticleStatus::STATUS_OK)
{
velocity = (k1 + 2 * k2 + 2 * k3 + k4) / 6.0f;
return ParticleStatus::STATUS_OK;
}
else
{
return ParticleStatus::EXITED_SPATIAL_BOUNDARY;
return ParticleStatus::AT_SPATIAL_BOUNDARY;
}
}
};
@ -223,7 +225,7 @@ public:
VTKM_EXEC
ParticleStatus CheckStep(const vtkm::Vec<FieldType, 3>& inpos,
FieldType stepLength,
FieldType vtkmNotUsed(stepLength),
FieldType time,
vtkm::Vec<FieldType, 3>& velocity) const
{

@ -416,7 +416,7 @@ void TestParticleWorklets()
using FieldHandle = vtkm::cont::ArrayHandle<vtkm::Vec<FieldType, 3>>;
using FieldPortalConstType = FieldHandle::template ExecutionTypes<DeviceAdapter>::PortalConst;
FieldType stepSize = 0.05f;
FieldType stepSize = 0.01f;
vtkm::cont::DataSetBuilderUniform dataSetBuilder;
@ -452,35 +452,78 @@ void TestParticleWorklets()
pts.push_back(vtkm::Vec<FieldType, 3>(2, 2, 2));
pts.push_back(vtkm::Vec<FieldType, 3>(3, 3, 3));
vtkm::Id maxSteps = 100;
vtkm::cont::ArrayHandle<vtkm::Vec<FieldType, 3>> seeds;
seeds = vtkm::cont::make_ArrayHandle(pts);
vtkm::Id maxSteps = 1000;
vtkm::Id nSeeds = static_cast<vtkm::Id>(pts.size());
std::vector<vtkm::Id> stepsTaken = { 10, 20, 600 };
if (i == 0)
{
vtkm::worklet::ParticleAdvection particleAdvection;
vtkm::worklet::ParticleAdvectionResult<FieldType> res;
res = particleAdvection.Run(rk4, seeds, maxSteps, DeviceAdapter());
VTKM_TEST_ASSERT(res.positions.GetNumberOfValues() == seeds.GetNumberOfValues(),
"Number of output particles does not match input.");
for (int j = 0; j < 2; j++)
{
vtkm::worklet::ParticleAdvection particleAdvection;
vtkm::worklet::ParticleAdvectionResult<FieldType> res;
vtkm::cont::ArrayHandle<vtkm::Vec<FieldType, 3>> seeds;
seeds.Allocate(nSeeds);
for (vtkm::Id k = 0; k < nSeeds; k++)
seeds.GetPortalControl().Set(k, pts[k]);
if (j == 0)
res = particleAdvection.Run(rk4, seeds, maxSteps, DeviceAdapter());
else if (j == 1)
{
vtkm::cont::ArrayHandle<vtkm::Id> stepsTakenArrayHandle =
vtkm::cont::make_ArrayHandle(stepsTaken);
res = particleAdvection.Run(rk4, seeds, stepsTakenArrayHandle, maxSteps, DeviceAdapter());
}
VTKM_TEST_ASSERT(res.positions.GetNumberOfValues() == seeds.GetNumberOfValues(),
"Number of output particles does not match input.");
std::cout << "Test: " << i << " " << j << std::endl;
vtkm::cont::printSummary_ArrayHandle(res.stepsTaken, std::cout);
for (vtkm::Id k = 0; k < nSeeds; k++)
VTKM_TEST_ASSERT(res.stepsTaken.GetPortalConstControl().Get(k) <= maxSteps,
"Too many steps taken in particle advection");
}
}
else if (i == 1)
{
vtkm::worklet::Streamline streamline;
vtkm::worklet::StreamlineResult<FieldType> res;
res = streamline.Run(rk4, seeds, maxSteps, DeviceAdapter());
//Make sure we have the right number of streamlines.
VTKM_TEST_ASSERT(res.polyLines.GetNumberOfCells() == seeds.GetNumberOfValues(),
"Number of output streamlines does not match input.");
//Make sure we have the right number of samples in each streamline.
vtkm::Id nSeeds = static_cast<vtkm::Id>(pts.size());
for (vtkm::Id j = 0; j < nSeeds; j++)
for (int j = 0; j < 2; j++)
{
vtkm::Id numPoints = static_cast<vtkm::Id>(res.polyLines.GetNumberOfPointsInCell(j));
vtkm::Id numSteps = res.stepsTaken.GetPortalConstControl().Get(j);
VTKM_TEST_ASSERT(numPoints == numSteps, "Invalid number of points in streamline.");
vtkm::worklet::Streamline streamline;
vtkm::worklet::StreamlineResult<FieldType> res;
vtkm::cont::ArrayHandle<vtkm::Vec<FieldType, 3>> seeds;
seeds.Allocate(nSeeds);
for (vtkm::Id k = 0; k < nSeeds; k++)
seeds.GetPortalControl().Set(k, pts[k]);
if (j == 0)
res = streamline.Run(rk4, seeds, maxSteps, DeviceAdapter());
else if (j == 1)
{
vtkm::cont::ArrayHandle<vtkm::Id> stepsTakenArrayHandle =
vtkm::cont::make_ArrayHandle(stepsTaken);
res = streamline.Run(rk4, seeds, stepsTakenArrayHandle, maxSteps, DeviceAdapter());
}
//Make sure we have the right number of streamlines.
VTKM_TEST_ASSERT(res.polyLines.GetNumberOfCells() == seeds.GetNumberOfValues(),
"Number of output streamlines does not match input.");
std::cout << "Test: " << i << " " << j << std::endl;
vtkm::cont::printSummary_ArrayHandle(res.stepsTaken, std::cout);
//Make sure we have the right number of samples in each streamline.
nSeeds = static_cast<vtkm::Id>(pts.size());
for (vtkm::Id k = 0; k < nSeeds; k++)
{
vtkm::Id numPoints = static_cast<vtkm::Id>(res.polyLines.GetNumberOfPointsInCell(k));
vtkm::Id numSteps = res.stepsTaken.GetPortalConstControl().Get(k);
VTKM_TEST_ASSERT(numPoints == numSteps, "Invalid number of points in streamline.");
VTKM_TEST_ASSERT(res.stepsTaken.GetPortalConstControl().Get(k) <= maxSteps,
"Too many steps taken in streamline");
}
}
}
}