mirror of
https://gitlab.kitware.com/vtk/vtk-m
synced 2024-09-19 18:45:43 +00:00
Parallelised half of the hypersweep.
This commit is contained in:
Petar Hristov
parent
40b3c7c2b8
commit
791843214c
@ -1405,104 +1405,16 @@ public:
|
||||
//{
|
||||
//printf("The min value in the subtree defined by the edge (%lld, %lld) is %lld.\n", i, maskedIndex(minParents.GetPortalConstControl().Get(i)), minValues.GetPortalConstControl().Get(i));
|
||||
//}
|
||||
|
||||
//
|
||||
// Make array for each arc (i, j, subtreeHeight, subtreeMin), then for with that first and second criteria
|
||||
//
|
||||
class EdgeData
|
||||
{
|
||||
public:
|
||||
Id i;
|
||||
Id j;
|
||||
Id subtreeMin;
|
||||
Id subtreeMax;
|
||||
bool upEdge;
|
||||
Float64 subtreeHeight;
|
||||
|
||||
bool operator<(const EdgeData& b) const
|
||||
{
|
||||
if (this->i == b.i)
|
||||
{
|
||||
if (this->upEdge == b.upEdge)
|
||||
{
|
||||
if (this->subtreeHeight == b.subtreeHeight)
|
||||
{
|
||||
return this->subtreeMin < b.subtreeMin;
|
||||
}
|
||||
else
|
||||
{
|
||||
return this->subtreeHeight > b.subtreeHeight;
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
return this->upEdge < b.upEdge;
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
return this->i < b.i;
|
||||
}
|
||||
}
|
||||
};
|
||||
|
||||
vtkm::cont::ArrayHandle<EdgeData> arcs;
|
||||
vtkm::cont::ArrayHandle<vtkm::worklet::contourtree_augmented::EdgeData> arcs;
|
||||
arcs.Allocate(contourTree.Superarcs.GetNumberOfValues() * 2 - 2);
|
||||
|
||||
timer.Reset();
|
||||
timer.Start();
|
||||
|
||||
//
|
||||
// Set arc endpoints. Parallelisable. No HS.
|
||||
//
|
||||
int edges = 0;
|
||||
for (vtkm::Id i = 0; i < contourTree.Superarcs.GetNumberOfValues(); i++)
|
||||
{
|
||||
Id parent = MaskedIndex(contourTree.Superarcs.ReadPortal().Get(i));
|
||||
|
||||
if (parent == 0)
|
||||
continue;
|
||||
|
||||
EdgeData edge;
|
||||
edge.i = i;
|
||||
edge.j = parent;
|
||||
edge.upEdge = IsAscending((contourTree.Superarcs.ReadPortal().Get(i)));
|
||||
|
||||
EdgeData oppositeEdge;
|
||||
oppositeEdge.i = parent;
|
||||
oppositeEdge.j = i;
|
||||
oppositeEdge.upEdge = !edge.upEdge;
|
||||
|
||||
|
||||
// Is it in the direction of the minRootedTree?
|
||||
if (MaskedIndex(minParents.ReadPortal().Get(edge.j)) == edge.i)
|
||||
{
|
||||
edge.subtreeMin = minValues.ReadPortal().Get(edge.j);
|
||||
oppositeEdge.subtreeMin = 0;
|
||||
}
|
||||
else
|
||||
{
|
||||
oppositeEdge.subtreeMin = minValues.ReadPortal().Get(oppositeEdge.j);
|
||||
edge.subtreeMin = 0;
|
||||
}
|
||||
|
||||
// Is it in the direction of the maxRootedTree?
|
||||
if (MaskedIndex(maxParents.ReadPortal().Get(edge.j)) == edge.i)
|
||||
{
|
||||
edge.subtreeMax = maxValues.ReadPortal().Get(edge.j);
|
||||
oppositeEdge.subtreeMax = contourTree.Arcs.GetNumberOfValues() - 1;
|
||||
}
|
||||
else
|
||||
{
|
||||
oppositeEdge.subtreeMax = maxValues.ReadPortal().Get(oppositeEdge.j);
|
||||
edge.subtreeMax = contourTree.Arcs.GetNumberOfValues() - 1;
|
||||
}
|
||||
|
||||
// We cannot use i here because one of the vertices is skipped (the root one and we don't know where it is)
|
||||
arcs.WritePortal().Set(edges++, edge);
|
||||
arcs.WritePortal().Set(edges++, oppositeEdge);
|
||||
}
|
||||
vtkm::worklet::contourtree_augmented::process_contourtree_inc::InitialiseArcs initArcs(
|
||||
0, contourTree.Arcs.GetNumberOfValues() - 1, minPath[minPath.size() - 1]);
|
||||
|
||||
Invoke(initArcs, minParents, maxParents, minValues, maxValues, contourTree.Superarcs, arcs);
|
||||
|
||||
//
|
||||
// Set whether an arc is up or down arcs. Parallelisable. No HS.
|
||||
@ -1515,31 +1427,15 @@ public:
|
||||
<< " seconds." << std::endl;
|
||||
}
|
||||
|
||||
//
|
||||
// Compute the height of all subtrees using the min and max
|
||||
//
|
||||
timer.Reset();
|
||||
timer.Start();
|
||||
|
||||
//
|
||||
// Set subtree height based on the best min & max. Parallelisable.
|
||||
//
|
||||
for (vtkm::Id i = 0; i < arcs.GetNumberOfValues(); i++)
|
||||
{
|
||||
EdgeData edge = arcs.ReadPortal().Get(i);
|
||||
|
||||
Float64 minIsoval =
|
||||
fieldValues.ReadPortal().Get(ctSortOrder.ReadPortal().Get(edge.subtreeMin));
|
||||
Float64 maxIsoval =
|
||||
fieldValues.ReadPortal().Get(ctSortOrder.ReadPortal().Get(edge.subtreeMax));
|
||||
Float64 vertexIsoval = fieldValues.ReadPortal().Get(
|
||||
ctSortOrder.ReadPortal().Get(contourTree.Supernodes.ReadPortal().Get(edge.i)));
|
||||
|
||||
// We need to incorporate the value of the vertex into the height of the tree (otherwise leafs edges have 0 persistence)
|
||||
minIsoval = minOperator(minIsoval, vertexIsoval);
|
||||
maxIsoval = maxOperator(maxIsoval, vertexIsoval);
|
||||
|
||||
edge.subtreeHeight = maxIsoval - minIsoval;
|
||||
|
||||
arcs.WritePortal().Set(i, edge);
|
||||
}
|
||||
vtkm::worklet::contourtree_augmented::process_contourtree_inc::ComputeSubtreeHeight
|
||||
computeSubtreeHeight;
|
||||
Invoke(computeSubtreeHeight, fieldValues, ctSortOrder, contourTree.Supernodes, arcs);
|
||||
|
||||
timer.Stop();
|
||||
if (true == printTime)
|
||||
@ -1548,6 +1444,9 @@ public:
|
||||
<< " seconds." << std::endl;
|
||||
}
|
||||
|
||||
//
|
||||
// Sort to find the bestUp and Down
|
||||
//
|
||||
timer.Reset();
|
||||
timer.Start();
|
||||
|
||||
@ -1560,41 +1459,15 @@ public:
|
||||
<< std::endl;
|
||||
}
|
||||
|
||||
timer.Reset();
|
||||
timer.Start();
|
||||
|
||||
//
|
||||
// Set bestUp and bestDown. Parallelisable
|
||||
//
|
||||
for (vtkm::Id i = 0; i < arcs.GetNumberOfValues(); i++)
|
||||
{
|
||||
if (i == 0)
|
||||
{
|
||||
if (arcs.ReadPortal().Get(0).upEdge == 0)
|
||||
{
|
||||
bestDownward.WritePortal().Set(arcs.ReadPortal().Get(0).i, arcs.ReadPortal().Get(0).j);
|
||||
}
|
||||
else
|
||||
{
|
||||
bestUpward.WritePortal().Set(arcs.ReadPortal().Get(0).i, arcs.ReadPortal().Get(0).j);
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
if (arcs.ReadPortal().Get(i).upEdge == 0 &&
|
||||
arcs.ReadPortal().Get(i).i != arcs.ReadPortal().Get(i - 1).i)
|
||||
{
|
||||
bestDownward.WritePortal().Set(arcs.ReadPortal().Get(i).i, arcs.ReadPortal().Get(i).j);
|
||||
}
|
||||
timer.Reset();
|
||||
timer.Start();
|
||||
|
||||
if (arcs.ReadPortal().Get(i).upEdge == 1 &&
|
||||
(arcs.ReadPortal().Get(i).i != arcs.ReadPortal().Get(i - 1).i ||
|
||||
arcs.ReadPortal().Get(i - 1).upEdge == 0))
|
||||
{
|
||||
bestUpward.WritePortal().Set(arcs.ReadPortal().Get(i).i, arcs.ReadPortal().Get(i).j);
|
||||
}
|
||||
}
|
||||
}
|
||||
vtkm::worklet::contourtree_augmented::process_contourtree_inc::SetBestUpDown setBestUpDown;
|
||||
Invoke(setBestUpDown, bestUpward, bestDownward, arcs);
|
||||
|
||||
timer.Stop();
|
||||
if (true == printTime)
|
||||
@ -1603,12 +1476,6 @@ public:
|
||||
<< " seconds." << std::endl;
|
||||
}
|
||||
|
||||
//for (int i = 0 ; i < bestUpward.GetNumberOfValues() ; i++)
|
||||
//{
|
||||
//printf("The bestUP of %2lld is %2lld.\n", i, bestUpward.GetPortalConstControl().Get(i));
|
||||
//printf("The bestDown of %2lld is %2lld.\n\n", i, bestDownward.GetPortalConstControl().Get(i));
|
||||
//}
|
||||
|
||||
timer.Reset();
|
||||
timer.Start();
|
||||
|
||||
|
||||
@ -146,6 +146,42 @@ inline std::string FlagString(vtkm::Id flaggedIndex)
|
||||
return fString;
|
||||
} // FlagString()
|
||||
|
||||
class EdgeData
|
||||
{
|
||||
public:
|
||||
Id i;
|
||||
Id j;
|
||||
Id subtreeMin;
|
||||
Id subtreeMax;
|
||||
bool upEdge;
|
||||
Float64 subtreeHeight;
|
||||
|
||||
bool operator<(const EdgeData& b) const
|
||||
{
|
||||
if (this->i == b.i)
|
||||
{
|
||||
if (this->upEdge == b.upEdge)
|
||||
{
|
||||
if (this->subtreeHeight == b.subtreeHeight)
|
||||
{
|
||||
return this->subtreeMin < b.subtreeMin;
|
||||
}
|
||||
else
|
||||
{
|
||||
return this->subtreeHeight > b.subtreeHeight;
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
return this->upEdge < b.upEdge;
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
return this->i < b.i;
|
||||
}
|
||||
}
|
||||
};
|
||||
|
||||
|
||||
} // namespace contourtree_augmented
|
||||
|
||||
@ -70,90 +70,173 @@ namespace process_contourtree_inc
|
||||
{
|
||||
|
||||
|
||||
|
||||
|
||||
template <typename Operator>
|
||||
class IncorporateEdge : public vtkm::worklet::WorkletMapField
|
||||
class InitialiseArcs : public vtkm::worklet::WorkletMapField
|
||||
{
|
||||
public:
|
||||
typedef void ControlSignature(WholeArrayIn minParents,
|
||||
WholeArrayIn supernodes,
|
||||
WholeArrayOut minMaxValues);
|
||||
typedef void ExecutionSignature(InputIndex, _1, _2, _3);
|
||||
typedef void ControlSignature(WholeArrayIn,
|
||||
WholeArrayIn,
|
||||
WholeArrayIn,
|
||||
WholeArrayIn,
|
||||
WholeArrayIn,
|
||||
WholeArrayInOut);
|
||||
typedef void ExecutionSignature(InputIndex, _1, _2, _3, _4, _5, _6);
|
||||
using InputDomain = _1;
|
||||
|
||||
Operator op;
|
||||
VTKM_EXEC_CONT IncorporateEdge(Operator _op)
|
||||
: op(_op)
|
||||
vtkm::Id globalMinSortedIndex, globalMaxSortedIndex, rootSupernodeId;
|
||||
|
||||
VTKM_EXEC_CONT InitialiseArcs(vtkm::Id _globalMinSortedIndex,
|
||||
vtkm::Id _globalMaxSortedIndex,
|
||||
vtkm::Id _rootSupernodeId)
|
||||
: globalMinSortedIndex(_globalMinSortedIndex)
|
||||
, globalMaxSortedIndex(_globalMaxSortedIndex)
|
||||
, rootSupernodeId(_rootSupernodeId)
|
||||
{
|
||||
}
|
||||
|
||||
template <typename IdWholeArrayInPortalType, typename IdWholeArrayOutPortalType>
|
||||
template <typename IdWholeArrayInPortalType, typename EdgeWholeArrayInOutPortal>
|
||||
VTKM_EXEC void operator()(const vtkm::Id currentId,
|
||||
const IdWholeArrayInPortalType& parentsPortal,
|
||||
const IdWholeArrayInPortalType& supernodesPortal,
|
||||
const IdWholeArrayOutPortalType& minMaxValuesPortal) const
|
||||
const IdWholeArrayInPortalType& minParentsPortal,
|
||||
const IdWholeArrayInPortalType& maxParentsPortal,
|
||||
const IdWholeArrayInPortalType& minValuesPortal,
|
||||
const IdWholeArrayInPortalType& maxValuesPortal,
|
||||
const IdWholeArrayInPortalType& superarcsPortal,
|
||||
const EdgeWholeArrayInOutPortal& arcsPortal) const
|
||||
{
|
||||
Id parent = MaskedIndex(parentsPortal.Get(currentId));
|
||||
Id subtreeValue = minMaxValuesPortal.Get(currentId);
|
||||
Id parentValue = MaskedIndex(supernodesPortal.Get(parent));
|
||||
minMaxValuesPortal.Set(currentId, op(parentValue, subtreeValue));
|
||||
Id i = currentId;
|
||||
Id parent = MaskedIndex(superarcsPortal.Get(i));
|
||||
if (parent == 0)
|
||||
return;
|
||||
|
||||
EdgeData edge;
|
||||
edge.i = i;
|
||||
edge.j = parent;
|
||||
edge.upEdge = IsAscending((superarcsPortal.Get(i)));
|
||||
|
||||
EdgeData oppositeEdge;
|
||||
oppositeEdge.i = parent;
|
||||
oppositeEdge.j = i;
|
||||
oppositeEdge.upEdge = !edge.upEdge;
|
||||
|
||||
|
||||
// Is it in the direction of the minRootedTree?
|
||||
if (MaskedIndex(minParentsPortal.Get(edge.j)) == edge.i)
|
||||
{
|
||||
edge.subtreeMin = minValuesPortal.Get(edge.j);
|
||||
oppositeEdge.subtreeMin = globalMinSortedIndex;
|
||||
}
|
||||
else
|
||||
{
|
||||
oppositeEdge.subtreeMin = minValuesPortal.Get(oppositeEdge.j);
|
||||
edge.subtreeMin = globalMinSortedIndex;
|
||||
}
|
||||
|
||||
// Is it in the direction of the maxRootedTree?
|
||||
if (MaskedIndex(maxParentsPortal.Get(edge.j)) == edge.i)
|
||||
{
|
||||
edge.subtreeMax = maxValuesPortal.Get(edge.j);
|
||||
oppositeEdge.subtreeMax = globalMaxSortedIndex;
|
||||
}
|
||||
else
|
||||
{
|
||||
oppositeEdge.subtreeMax = maxValuesPortal.Get(oppositeEdge.j);
|
||||
edge.subtreeMax = globalMaxSortedIndex;
|
||||
}
|
||||
|
||||
// Compensate for the missing edge where the root is
|
||||
if (i > rootSupernodeId)
|
||||
{
|
||||
i--;
|
||||
}
|
||||
|
||||
// We cannot use i here because one of the vertices is skipped (the root one and we don't know where it is)
|
||||
arcsPortal.Set(i * 2, edge);
|
||||
arcsPortal.Set(i * 2 + 1, oppositeEdge);
|
||||
}
|
||||
}; // ComputeMinMaxValues
|
||||
|
||||
|
||||
|
||||
|
||||
class GetOppositeValue : public vtkm::worklet::WorkletMapField
|
||||
class ComputeSubtreeHeight : public vtkm::worklet::WorkletMapField
|
||||
{
|
||||
public:
|
||||
typedef void ControlSignature(WholeArrayIn minParents,
|
||||
WholeArrayIn maxParents,
|
||||
WholeArrayIn minValues,
|
||||
WholeArrayIn maxValues,
|
||||
WholeArrayInOut arcs);
|
||||
typedef void ExecutionSignature(InputIndex, _1, _2, _3, _4, _5);
|
||||
using InputDomain = _1;
|
||||
typedef void ControlSignature(WholeArrayIn, WholeArrayIn, WholeArrayIn, WholeArrayInOut);
|
||||
typedef void ExecutionSignature(InputIndex, _1, _2, _3, _4);
|
||||
using InputDomain = _4;
|
||||
|
||||
vtkm::Id globalMinSortedIndex, globalMaxSortedIndex;
|
||||
VTKM_EXEC_CONT GetOppositeValue(vtkm::Id _globalMinSortedIndex, vtkm::Id _globalMaxSortedIndex)
|
||||
: globalMinSortedIndex(_globalMinSortedIndex)
|
||||
, globalMaxSortedIndex(_globalMaxSortedIndex)
|
||||
VTKM_EXEC_CONT ComputeSubtreeHeight() {}
|
||||
|
||||
template <typename Float64WholeArrayInPortalType,
|
||||
typename IdWholeArrayInPortalType,
|
||||
typename EdgeWholeArrayInOutPortal>
|
||||
VTKM_EXEC void operator()(const vtkm::Id currentId,
|
||||
const Float64WholeArrayInPortalType& fieldValuesPortal,
|
||||
const IdWholeArrayInPortalType& ctSortOrderPortal,
|
||||
const IdWholeArrayInPortalType& supernodesPortal,
|
||||
const EdgeWholeArrayInOutPortal& arcsPortal) const
|
||||
{
|
||||
Id i = currentId;
|
||||
EdgeData edge = arcsPortal.Get(i);
|
||||
|
||||
Float64 minIsoval = fieldValuesPortal.Get(ctSortOrderPortal.Get(edge.subtreeMin));
|
||||
Float64 maxIsoval = fieldValuesPortal.Get(ctSortOrderPortal.Get(edge.subtreeMax));
|
||||
Float64 vertexIsoval =
|
||||
fieldValuesPortal.Get(ctSortOrderPortal.Get(supernodesPortal.Get(edge.i)));
|
||||
|
||||
// We need to incorporate the value of the vertex into the height of the tree (otherwise leafs edges have 0 persistence)
|
||||
minIsoval = vtkm::Minimum()(minIsoval, vertexIsoval);
|
||||
maxIsoval = vtkm::Maximum()(maxIsoval, vertexIsoval);
|
||||
|
||||
edge.subtreeHeight = maxIsoval - minIsoval;
|
||||
|
||||
arcsPortal.Set(i, edge);
|
||||
}
|
||||
}; // ComputeMinMaxValues
|
||||
|
||||
|
||||
|
||||
|
||||
class SetBestUpDown : public vtkm::worklet::WorkletMapField
|
||||
{
|
||||
public:
|
||||
typedef void ControlSignature(WholeArrayInOut, WholeArrayInOut, WholeArrayIn);
|
||||
typedef void ExecutionSignature(InputIndex, _1, _2, _3);
|
||||
using InputDomain = _3;
|
||||
|
||||
VTKM_EXEC_CONT SetBestUpDown() {}
|
||||
|
||||
template <typename IdWholeArrayInPortalType, typename EdgeWholeArrayInOutPortal>
|
||||
VTKM_EXEC void operator()(const vtkm::Id currentId,
|
||||
const IdWholeArrayInPortalType& minParents,
|
||||
const IdWholeArrayInPortalType& maxParents,
|
||||
const IdWholeArrayInPortalType& minValues,
|
||||
const IdWholeArrayInPortalType& maxValues,
|
||||
const EdgeWholeArrayInOutPortal& arcs) const
|
||||
const IdWholeArrayInPortalType& bestUpwardPortal,
|
||||
const IdWholeArrayInPortalType& bestDownwardPortal,
|
||||
const EdgeWholeArrayInOutPortal& arcsPortal) const
|
||||
{
|
||||
auto i = currentId;
|
||||
auto edge = arcs.Get(i);
|
||||
vtkm::Id i = currentId;
|
||||
|
||||
// Is it in the direction of the minRootedTree?
|
||||
if (MaskedIndex(minParents.ReadPortal().Get(edge.j)) == edge.i)
|
||||
if (i == 0)
|
||||
{
|
||||
edge.subtreeMin = minValues.ReadPortal().Get(edge.j);
|
||||
if (arcsPortal.Get(0).upEdge == 0)
|
||||
{
|
||||
bestDownwardPortal.Set(arcsPortal.Get(0).i, arcsPortal.Get(0).j);
|
||||
}
|
||||
else
|
||||
{
|
||||
bestUpwardPortal.Set(arcsPortal.Get(0).i, arcsPortal.Get(0).j);
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
edge.subtreeMin = globalMinSortedIndex;
|
||||
}
|
||||
if (arcsPortal.Get(i).upEdge == 0 && arcsPortal.Get(i).i != arcsPortal.Get(i - 1).i)
|
||||
{
|
||||
bestDownwardPortal.Set(arcsPortal.Get(i).i, arcsPortal.Get(i).j);
|
||||
}
|
||||
|
||||
// Is it in the direction of the maxRootedTree?
|
||||
if (MaskedIndex(maxParents.ReadPortal().Get(edge.j)) == edge.i)
|
||||
{
|
||||
edge.subtreeMax = maxValues.ReadPortal().Get(edge.j);
|
||||
if (arcsPortal.Get(i).upEdge == 1 &&
|
||||
(arcsPortal.Get(i).i != arcsPortal.Get(i - 1).i || arcsPortal.Get(i - 1).upEdge == 0))
|
||||
{
|
||||
bestUpwardPortal.Set(arcsPortal.Get(i).i, arcsPortal.Get(i).j);
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
edge.subtreeMax = globalMinSortedIndex;
|
||||
}
|
||||
|
||||
arcs.Set(i, edge);
|
||||
}
|
||||
}; // ComputeMinMaxValues
|
||||
|
||||
|
||||
Loading…
Reference in New Issue
Block a user