89 #ifndef OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED 90 #define OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED 92 #include <openvdb/Types.h> 93 #include <openvdb/math/ConjGradient.h> 94 #include <openvdb/tree/LeafManager.h> 95 #include <openvdb/tree/Tree.h> 96 #include <openvdb/util/NullInterrupter.h> 100 #include <boost/scoped_array.hpp> 117 template<
typename TreeType>
130 inline typename TreeType::Ptr
133 template<
typename TreeType,
typename Interrupter>
134 inline typename TreeType::Ptr
140 template<
typename TreeType,
typename BoundaryOp,
typename Interrupter>
171 inline typename TreeType::Ptr
174 template<
typename PreconditionerType,
typename TreeType,
typename BoundaryOp,
typename Interrupter>
175 inline typename TreeType::Ptr
179 template<
typename PreconditionerType,
typename TreeType,
typename DomainTreeType,
typename BoundaryOp,
typename Interrupter>
180 inline typename TreeType::Ptr
192 template<
typename VIndexTreeType>
198 template<
typename TreeType>
199 inline typename TreeType::template ValueConverter<VIndex>::Type::Ptr
209 template<
typename VectorValueType,
typename SourceTreeType>
212 const SourceTreeType& source,
213 const typename SourceTreeType::template ValueConverter<VIndex>::Type& index);
223 template<
typename TreeValueType,
typename VIndexTreeType,
typename VectorValueType>
224 inline typename VIndexTreeType::template ValueConverter<TreeValueType>::Type::Ptr
227 const VIndexTreeType& index,
228 const TreeValueType& background);
235 template<
typename BoolTreeType>
238 const typename BoolTreeType::template ValueConverter<VIndex>::Type& vectorIndexTree,
239 const BoolTreeType& interiorMask);
261 template<
typename BoolTreeType,
typename BoundaryOp>
264 const typename BoolTreeType::template ValueConverter<VIndex>::Type& vectorIndexTree,
265 const BoolTreeType& interiorMask,
266 const BoundaryOp& boundaryOp,
279 template<
typename LeafType>
284 void operator()(
const LeafType& leaf,
size_t leafIdx)
const {
285 count[leafIdx] =
static_cast<VIndex
>(leaf.onVoxelCount());
292 template<
typename LeafType>
298 VIndex idx = (leafIdx == 0) ? 0 : count[leafIdx - 1];
299 for (
typename LeafType::ValueOnIter it = leaf.beginValueOn(); it; ++it) {
308 template<
typename VIndexTreeType>
312 typedef typename VIndexTreeType::LeafNodeType LeafT;
316 LeafMgrT leafManager(result);
317 const size_t leafCount = leafManager.leafCount();
320 boost::scoped_array<VIndex> perLeafCount(
new VIndex[leafCount]);
321 VIndex* perLeafCountPtr = perLeafCount.get();
326 for (
size_t i = 1; i < leafCount; ++i) {
327 perLeafCount[i] += perLeafCount[i - 1];
331 assert(
Index64(perLeafCount[leafCount-1]) == result.activeVoxelCount());
339 template<
typename TreeType>
340 inline typename TreeType::template ValueConverter<VIndex>::Type::Ptr
343 typedef typename TreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
346 const VIndex invalidIdx = -1;
347 typename VIdxTreeT::Ptr result(
351 result->voxelizeActiveTiles();
366 template<
typename VectorValueType,
typename SourceTreeType>
369 typedef typename SourceTreeType::template ValueConverter<VIndex>::Type
VIdxTreeT;
371 typedef typename SourceTreeType::LeafNodeType
LeafT;
378 CopyToVecOp(
const SourceTreeType& t, VectorT& v): tree(&t), vector(&v) {}
382 VectorT& vec = *vector;
383 if (
const LeafT* leaf = tree->probeLeaf(idxLeaf.origin())) {
386 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
387 vec[*it] = leaf->getValue(it.pos());
392 const TreeValueT& value = tree->getValue(idxLeaf.origin());
393 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
403 template<
typename VectorValueType,
typename SourceTreeType>
406 const typename SourceTreeType::template ValueConverter<VIndex>::Type& idxTree)
408 typedef typename SourceTreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
413 const size_t numVoxels = idxTree.activeVoxelCount();
414 typename VectorT::Ptr result(
new VectorT(static_cast<math::pcg::SizeType>(numVoxels)));
418 VIdxLeafMgrT leafManager(idxTree);
432 template<
typename TreeValueType,
typename VIndexTreeType,
typename VectorValueType>
435 typedef typename VIndexTreeType::template ValueConverter<TreeValueType>::Type
OutTreeT;
437 typedef typename VIndexTreeType::LeafNodeType
VIdxLeafT;
447 const VectorT& vec = *vector;
448 OutLeafT* leaf = tree->probeLeaf(idxLeaf.origin());
449 assert(leaf != NULL);
450 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
451 leaf->setValueOnly(it.pos(),
static_cast<TreeValueType
>(vec[*it]));
459 template<
typename TreeValueType,
typename VIndexTreeType,
typename VectorValueType>
460 inline typename VIndexTreeType::template ValueConverter<TreeValueType>::Type::Ptr
463 const VIndexTreeType& idxTree,
464 const TreeValueType& background)
466 typedef typename VIndexTreeType::template ValueConverter<TreeValueType>::Type OutTreeT;
470 typename OutTreeT::Ptr result(
new OutTreeT(idxTree, background,
TopologyCopy()));
471 OutTreeT& tree = *result;
475 VIdxLeafMgrT leafManager(idxTree);
489 template<
typename ValueType>
492 const Coord&,
const Coord&, ValueType&, ValueType& diag)
const { diag -= 1; }
497 template<
typename BoolTreeType,
typename BoundaryOp>
500 typedef typename BoolTreeType::template ValueConverter<VIndex>::Type
VIdxTreeT;
512 const BoolTreeType& mask,
const BoundaryOp& op, VectorT& src):
513 laplacian(&m), idxTree(&idx), interiorMask(&mask), boundaryOp(op), source(&src) {}
523 const ValueT diagonal = -6.f, offDiagonal = 1.f;
526 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
527 assert(it.getValue() > -1);
556 ValueT modifiedDiagonal = 0.f;
559 if (vectorIdx.
probeValue(ijk.offsetBy(-1, 0, 0), column)) {
561 modifiedDiagonal -= 1;
563 boundaryOp(ijk, ijk.offsetBy(-1, 0, 0), source->
at(rowNum), modifiedDiagonal);
566 if (vectorIdx.
probeValue(ijk.offsetBy(0, -1, 0), column)) {
568 modifiedDiagonal -= 1;
570 boundaryOp(ijk, ijk.offsetBy(0, -1, 0), source->
at(rowNum), modifiedDiagonal);
573 if (vectorIdx.
probeValue(ijk.offsetBy(0, 0, -1), column)) {
575 modifiedDiagonal -= 1;
577 boundaryOp(ijk, ijk.offsetBy(0, 0, -1), source->
at(rowNum), modifiedDiagonal);
580 if (vectorIdx.
probeValue(ijk.offsetBy(0, 0, 1), column)) {
582 modifiedDiagonal -= 1;
584 boundaryOp(ijk, ijk.offsetBy(0, 0, 1), source->
at(rowNum), modifiedDiagonal);
587 if (vectorIdx.
probeValue(ijk.offsetBy(0, 1, 0), column)) {
589 modifiedDiagonal -= 1;
591 boundaryOp(ijk, ijk.offsetBy(0, 1, 0), source->
at(rowNum), modifiedDiagonal);
594 if (vectorIdx.
probeValue(ijk.offsetBy(1, 0, 0), column)) {
596 modifiedDiagonal -= 1;
598 boundaryOp(ijk, ijk.offsetBy(1, 0, 0), source->
at(rowNum), modifiedDiagonal);
601 row.
setValue(rowNum, modifiedDiagonal);
610 template<
typename BoolTreeType>
613 const BoolTreeType& interiorMask)
617 static_cast<math::pcg::SizeType>(idxTree.activeVoxelCount()));
623 template<
typename BoolTreeType,
typename BoundaryOp>
626 const typename BoolTreeType::template ValueConverter<VIndex>::Type& idxTree,
627 const BoolTreeType& interiorMask,
628 const BoundaryOp& boundaryOp,
631 typedef typename BoolTreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
635 const Index64 numDoF = idxTree.activeVoxelCount();
640 LaplacianMatrix&
laplacian = *laplacianPtr;
643 VIdxLeafMgrT idxLeafManager(idxTree);
645 laplacian, idxTree, interiorMask, boundaryOp, source));
654 template<
typename TreeType>
655 inline typename TreeType::Ptr
659 return solve(inTree, state, interrupter);
663 template<
typename TreeType,
typename Interrupter>
664 inline typename TreeType::Ptr
672 template<
typename TreeType,
typename BoundaryOp,
typename Interrupter>
673 inline typename TreeType::Ptr
678 return solveWithBoundaryConditionsAndPreconditioner<DefaultPrecondT>(
679 inTree, boundaryOp, state, interrupter);
683 template<
typename PreconditionerType,
typename TreeType,
typename BoundaryOp,
typename Interrupter>
684 inline typename TreeType::Ptr
686 const BoundaryOp& boundaryOp,
math::pcg::State& state, Interrupter& interrupter)
689 return solveWithBoundaryConditionsAndPreconditioner<PreconditionerType>(inTree , inTree ,
690 boundaryOp, state, interrupter);
693 template<
typename PreconditionerType,
typename TreeType,
typename DomainTreeType,
typename BoundaryOp,
typename Interrupter>
694 inline typename TreeType::Ptr
696 const DomainTreeType& domainMask,
697 const BoundaryOp& boundaryOp,
701 typedef typename TreeType::ValueType TreeValueT;
704 typedef typename TreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
705 typedef typename TreeType::template ValueConverter<bool>::Type MaskTreeT;
711 typename VectorT::Ptr b = createVectorFromTree<VecValueT>(inTree, *idxTree);
714 typename MaskTreeT::Ptr interiorMask(
720 *idxTree, *interiorMask, boundaryOp, *b);
723 laplacian->scale(-1.0);
725 typename VectorT::Ptr x(
new VectorT(b->size(), zeroVal<VecValueT>()));
727 new PreconditionerType(*laplacian));
728 if (!precond->isValid()) {
736 return createTreeFromVector<TreeValueT>(*x, *idxTree, zeroVal<TreeValueT>());
744 #endif // OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED
boost::shared_ptr< SparseStencilMatrix > Ptr
Definition: ConjGradient.h:271
Read/write accessor to a row of this matrix.
Definition: ConjGradient.h:448
Preconditioner using incomplete Cholesky factorization.
Definition: ConjGradient.h:71
This class manages a linear array of pointers to a given tree's leaf nodes, as well as optional auxil...
Definition: LeafManager.h:114
Diagonal preconditioner.
Definition: ConjGradient.h:70
ValueType_ ValueType
Definition: ConjGradient.h:269
Index32 SizeType
Definition: ConjGradient.h:61
const ValueType & getValue(const Coord &xyz) const
Return the value of the voxel at the given coordinates.
Definition: ValueAccessor.h:256
Sparse, square matrix representing a 3D stencil operator of size STENCIL_SIZE.
Definition: ConjGradient.h:67
SizeType setValue(SizeType column, const ValueType &value)
Set the value of the entry in the specified column.
Definition: ConjGradient.h:1265
#define OPENVDB_VERSION_NAME
Definition: version.h:43
bool isValueOn(const Coord &xyz) const
Return the active state of the voxel at the given coordinates.
Definition: ValueAccessor.h:263
Dummy NOOP interrupter class defining interface.
Definition: NullInterrupter.h:52
Lightweight, variable-length vector.
Definition: ConjGradient.h:65
Information about the state of a conjugate gradient solution.
Definition: ConjGradient.h:74
RowEditor getRowEditor(SizeType row)
Return a read/write view onto the given row of this matrix.
Definition: ConjGradient.h:1109
Implementation of morphological dilation and erosion.
T & at(SizeType i)
Return the value of this vector's ith element.
Definition: ConjGradient.h:229
Definition: Exceptions.h:39
Definition: ValueAccessor.h:219
int32_t Int32
Definition: Types.h:60
boost::shared_ptr< Vector > Ptr
Definition: ConjGradient.h:170
boost::shared_ptr< Preconditioner > Ptr
Definition: ConjGradient.h:497
uint64_t Index64
Definition: Types.h:57
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:71
bool probeValue(const Coord &xyz, ValueType &value) const
Return the active state of the voxel as well as its value.
Definition: ValueAccessor.h:266