Shuenhoy/vtk-m2
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merge topic 'add_vorticity_and_qcriterion'

Browse Source 
f71f3ed0 Allow the gradient filter to compute Vorticity and QCriterion.
efb119ea Sort the worklet headers, as they should be in alphabetical order.

Acked-by: Kitware Robot <kwrobot@kitware.com>
Merge-request: !638
This commit is contained in:
Robert Maynard 2016-12-19 15:31:33 -05:00 committed by Kitware Robot
commit d5c80802c4
7 changed files with 229 additions and 34 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

32
vtkm/filter/Gradient.h
View File

@ -22,7 +22,6 @@
#define vtk_m_filter_Gradient_h
#include <vtkm/filter/FilterCell.h>
#include <vtkm/worklet/Gradient.h>
namespace vtkm {
namespace filter {
@ -38,7 +37,7 @@ namespace filter {
class Gradient : public vtkm::filter::FilterCell<Gradient>
{
public:
Gradient(): ComputePointGradient(false) {}
Gradient();
/// When this flag is on (default is off), the gradient filter will provide a
/// point based gradients, which are significantly more costly since for each
@ -46,6 +45,30 @@ public:
void SetComputePointGradient(bool enable) { ComputePointGradient=enable; }
bool GetComputePointGradient() const { return ComputePointGradient; }
/// Add voriticity/curl field to the output data. The name of the array
/// will be Vorticity and will be a cell field unless \c ComputePointGradient
/// is enabled. The input array must have 3 components in order to
/// compute this. The default is off.
void SetComputeVorticity(bool enable) { ComputeVorticity=enable; }
bool GetComputeVorticity() const { return ComputeVorticity; }
/// Add Q-criterion field to the output data. The name of the array
/// will be QCriterion and will be a cell field unless \c ComputePointGradient
/// is enabled. The input array must have 3 components in order to
/// compute this. The default is off.
void SetComputeQCriterion(bool enable) { ComputeQCriterion=enable; }
bool GetComputeQCriterion() const { return ComputeQCriterion; }
void SetVorticityName( const std::string& name )
{ this->VorticityName = name; }
const std::string& GetVorticityName() const
{ return this->VorticityName; }
void SetQCriterionName( const std::string& name )
{ this->QCriterionName = name; }
const std::string& GetQCriterionName() const
{ return this->QCriterionName; }
template<typename T, typename StorageType, typename DerivedPolicy, typename DeviceAdapter>
vtkm::filter::ResultField DoExecute(const vtkm::cont::DataSet &input,
@ -56,6 +79,11 @@ public:
private:
bool ComputePointGradient;
bool ComputeVorticity;
bool ComputeQCriterion;
std::string VorticityName;
std::string QCriterionName;
};
template<>

86
vtkm/filter/Gradient.hxx
View File

@ -22,8 +22,12 @@
#include <vtkm/worklet/DispatcherMapTopology.h>
#include <vtkm/worklet/Gradient.h>
#include <vtkm/worklet/Vorticity.h>
namespace {
//-----------------------------------------------------------------------------
template<typename DerivedPolicy, typename T, typename S>
struct PointGrad
{
@ -53,12 +57,83 @@ struct PointGrad
vtkm::cont::ArrayHandle< vtkm::Vec<T,3> >* Result;
};
//-----------------------------------------------------------------------------
template<typename HandleType>
inline
void add_field( vtkm::filter::ResultField& result,
const HandleType& handle,
const std::string name)
{
const vtkm::cont::Field::AssociationEnum assoc = result.GetField().GetAssociation();
if ((assoc == vtkm::cont::Field::ASSOC_WHOLE_MESH) ||
(assoc == vtkm::cont::Field::ASSOC_POINTS))
{
vtkm::cont::Field field(name, assoc, handle);
result.GetDataSet().AddField(field);
}
else
{
vtkm::cont::Field field(name,
assoc,
result.GetField().GetAssocCellSet(),
handle);
result.GetDataSet().AddField(field);
}
}
//-----------------------------------------------------------------------------
template<typename T, typename S, typename DeviceAdapter>
inline
void add_extra_vec_fields( const vtkm::cont::ArrayHandle< vtkm::Vec< vtkm::Vec<T,3>, 3>, S> &inField,
const vtkm::filter::Gradient* const filter,
vtkm::filter::ResultField& result,
const DeviceAdapter&)
{
if(filter->GetComputeVorticity())
{
vtkm::cont::ArrayHandle< vtkm::Vec<T,3> > vorticity;
vtkm::worklet::DispatcherMapField< vtkm::worklet::Vorticity, DeviceAdapter > dispatcher;
dispatcher.Invoke(inField, vorticity);
add_field(result, vorticity, filter->GetVorticityName());
}
if(filter->GetComputeQCriterion())
{
vtkm::cont::ArrayHandle< T > qc;
vtkm::worklet::DispatcherMapField< vtkm::worklet::QCriterion, DeviceAdapter > dispatcher;
dispatcher.Invoke(inField, qc);
add_field(result, qc, filter->GetQCriterionName());
}
}
template<typename T, typename S, typename DeviceAdapter>
inline
void add_extra_vec_fields( const vtkm::cont::ArrayHandle< T, S> &,
const vtkm::filter::Gradient* const,
vtkm::filter::ResultField&,
const DeviceAdapter&)
{
//not a vector array handle so add nothing
}
} //namespace
namespace vtkm {
namespace filter {
//-----------------------------------------------------------------------------
Gradient::Gradient():
ComputePointGradient(false),
ComputeVorticity(false),
ComputeQCriterion(false),
VorticityName("Vorticity"),
QCriterionName("QCriterion")
{
}
//-----------------------------------------------------------------------------
template<typename T,
typename StorageType,
@ -70,7 +145,7 @@ vtkm::filter::ResultField Gradient::DoExecute(
const vtkm::cont::ArrayHandle<T, StorageType> &inField,
const vtkm::filter::FieldMetadata &fieldMetadata,
const vtkm::filter::PolicyBase<DerivedPolicy>& policy,
const DeviceAdapter&)
const DeviceAdapter& adapter)
{
if(!fieldMetadata.IsPointField())
{
@ -115,12 +190,17 @@ vtkm::filter::ResultField Gradient::DoExecute(
fieldAssociation = vtkm::cont::Field::ASSOC_CELL_SET;
}
return vtkm::filter::ResultField(input,
vtkm::filter::ResultField result(input,
outArray,
outputName,
fieldAssociation,
cells.GetName());
//Add the vorticity and qcriterion fields if they are enabled to the result
add_extra_vec_fields(outArray, this, result, adapter);
return result;
}

12
vtkm/filter/ResultBase.h
View File

@ -50,6 +50,10 @@ public:
VTKM_CONT
const vtkm::cont::DataSet &GetDataSet() const { return this->Data; }
/// Returns the results of the filter in terms of a writable \c DataSet.
VTKM_CONT
vtkm::cont::DataSet &GetDataSet() { return this->Data; }
protected:
VTKM_CONT
ResultBase(): Valid(false) { }
@ -71,14 +75,6 @@ protected:
this->SetValid(true);
}
/// Returns a writable reference to the data set.
///
VTKM_CONT
vtkm::cont::DataSet &GetDataSetReference()
{
return this->Data;
}
private:
bool Valid;
vtkm::cont::DataSet Data;

12
vtkm/filter/ResultDataSet.h
View File

@ -45,18 +45,6 @@ public:
VTKM_CONT
ResultDataSet(const vtkm::cont::DataSet &dataSet)
: ResultBase(dataSet) { }
VTKM_CONT
const vtkm::cont::DataSet &GetDataSet() const
{
return this->ResultBase::GetDataSet();
}
VTKM_CONT
vtkm::cont::DataSet &GetDataSet()
{
return this->ResultBase::GetDataSetReference();
}
};
}

24
vtkm/filter/testing/UnitTestGradient.cxx
View File

@ -36,16 +36,26 @@ void TestCellGradientUniform3D()
vtkm::filter::ResultField result;
vtkm::filter::Gradient gradient;
gradient.SetOutputFieldName("gradient");
gradient.SetOutputFieldName("Gradient");
gradient.SetComputeVorticity(true); //this wont work as we have a scalar field
gradient.SetComputeQCriterion(true); //this wont work as we have a scalar field
result = gradient.Execute( dataSet, dataSet.GetField("pointvar"));
VTKM_TEST_ASSERT(result.GetField().GetName() == "gradient",
VTKM_TEST_ASSERT(result.GetField().GetName() == "Gradient",
"Field was given the wrong name.");
VTKM_TEST_ASSERT(result.GetField().GetAssociation() ==
vtkm::cont::Field::ASSOC_CELL_SET,
"Field was given the wrong association.");
//verify that the vorticity and qcriterion fields don't exist
const vtkm::cont::DataSet& outputDS = result.GetDataSet();
VTKM_TEST_ASSERT(outputDS.HasField("Vorticity") == false,
"scalar gradients can't generate vorticity");
VTKM_TEST_ASSERT(outputDS.HasField("QCriterion") == false,
"scalar gradients can't generate qcriterion");
vtkm::cont::ArrayHandle< vtkm::Vec<vtkm::Float32,3> > resultArrayHandle;
const bool valid = result.FieldAs(resultArrayHandle);
VTKM_TEST_ASSERT( valid, "result of gradient is not expected type");
@ -89,6 +99,8 @@ void TestCellGradientUniform3DWithVectorField()
vtkm::filter::ResultField result;
vtkm::filter::Gradient gradient;
gradient.SetOutputFieldName("vec_gradient");
gradient.SetComputeVorticity(true);
gradient.SetComputeQCriterion(true);
result = gradient.Execute( dataSet, dataSet.GetField("vec_pointvar"));
@ -99,6 +111,14 @@ void TestCellGradientUniform3DWithVectorField()
vtkm::cont::Field::ASSOC_CELL_SET,
"Field was given the wrong association.");
//verify that the vorticity and qcriterion fields DO exist
const vtkm::cont::DataSet& outputDS = result.GetDataSet();
VTKM_TEST_ASSERT(outputDS.HasField("Vorticity") == true,
"vec gradients should generate vorticity");
VTKM_TEST_ASSERT(outputDS.HasField("QCriterion") == true,
"vec gradients should generate qcriterion");
vtkm::cont::ArrayHandle< vtkm::Vec< vtkm::Vec<vtkm::Float64,3>, 3> > resultArrayHandle;
const bool valid = result.FieldAs(resultArrayHandle);
VTKM_TEST_ASSERT( valid, "result of gradient is not expected type");

15
vtkm/worklet/CMakeLists.txt
View File

@ -19,15 +19,13 @@
##============================================================================
set(headers
DispatcherMapField.h
DispatcherStreamingMapField.h
WorkletMapField.h
DispatcherMapTopology.h
WorkletMapTopology.h
AverageByKey.h
CellAverage.h
CellDeepCopy.h
Clip.h
DispatcherMapField.h
DispatcherMapTopology.h
DispatcherStreamingMapField.h
ExternalFaces.h
FieldHistogram.h
FieldStatistics.h
@ -36,8 +34,8 @@ set(headers
Magnitude.h
MarchingCubes.h
MarchingCubesDataTables.h
PointElevation.h
PointAverage.h
PointElevation.h
RemoveUnusedPoints.h
ScatterCounting.h
ScatterIdentity.h
@ -45,11 +43,14 @@ set(headers
StreamLineUniformGrid.h
TetrahedralizeExplicitGrid.h
TetrahedralizeUniformGrid.h
Threshold.h
TriangulateExplicitGrid.h
TriangulateUniformGrid.h
Threshold.h
VertexClustering.h
Vorticity.h
WaveletCompressor.h
WorkletMapField.h
WorkletMapTopology.h
)
#-----------------------------------------------------------------------------

82
vtkm/worklet/Vorticity.h Normal file
View File

@ -0,0 +1,82 @@
//============================================================================
// 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 Sandia Corporation.
// Copyright 2014 UT-Battelle, LLC.
// Copyright 2014 Los Alamos National Security.
//
// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
// 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 vtk_m_worklet_Vorticity_h
#define vtk_m_worklet_Vorticity_h
#include <vtkm/worklet/WorkletMapField.h>
namespace vtkm {
namespace worklet {
struct VorticityInputTypes
: vtkm::ListTagBase<
vtkm::Vec< vtkm::Vec<vtkm::Float32,3>, 3>,
vtkm::Vec< vtkm::Vec<vtkm::Float64,3>, 3>
> { };
struct Vorticity : public vtkm::worklet::WorkletMapField
{
typedef void ControlSignature(FieldIn<VorticityInputTypes> input,
FieldOut<Vec3> output);
typedef void ExecutionSignature(_1,_2);
typedef _1 InputDomain;
VTKM_EXEC
template<typename InputType, typename OutputType>
void operator()(const InputType &input, OutputType &vorticity) const
{
vorticity[0] = input[2][1] - input[1][2];
vorticity[1] = input[0][2] - input[2][0];
vorticity[2] = input[1][0] - input[0][1];
}
};
struct QCriterion : public vtkm::worklet::WorkletMapField
{
typedef void ControlSignature(FieldIn<VorticityInputTypes> input,
FieldOut<Scalar> output);
typedef void ExecutionSignature(_1,_2);
typedef _1 InputDomain;
VTKM_EXEC
template<typename InputType, typename OutputType>
void operator()(const InputType &input, OutputType &qcriterion) const
{
OutputType t1 =
((input[2][1] - input[1][2]) * (input[2][1] - input[1][2]) +
(input[1][0] - input[0][1]) * (input[1][0] - input[0][1]) +
(input[0][2] - input[2][0]) * (input[0][2] - input[2][0])) / 2.0f;
OutputType t2 =
input[0][0] * input[0][0] + input[1][1] * input[1][1] +
input[2][2] * input[2][2] +
((input[1][0] + input[0][1]) * (input[1][0] + input[0][1]) +
(input[2][0] + input[0][2]) * (input[2][0] + input[0][2]) +
(input[2][1] + input[1][2]) * (input[2][1] + input[1][2])) / 2.0f;
qcriterion = (t1 - t2) / 2.0f;
}
};
}
} // namespace vtkm::worklet
#endif