88 #ifndef OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED 89 #define OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED 113 template<
typename TreeType>
126 inline typename TreeType::Ptr
129 template<
typename TreeType,
typename Interrupter>
130 inline typename TreeType::Ptr
136 template<
typename TreeType,
typename BoundaryOp,
typename Interrupter>
167 inline typename TreeType::Ptr
173 bool staggered =
false);
176 typename PreconditionerType,
179 typename Interrupter>
180 inline typename TreeType::Ptr
186 bool staggered =
false);
189 typename PreconditionerType,
191 typename DomainTreeType,
193 typename Interrupter>
194 inline typename TreeType::Ptr
197 const DomainTreeType&,
201 bool staggered =
false);
211 template<
typename VIndexTreeType>
217 template<
typename TreeType>
218 inline typename TreeType::template ValueConverter<VIndex>::Type::Ptr
228 template<
typename VectorValueType,
typename SourceTreeType>
231 const SourceTreeType& source,
232 const typename SourceTreeType::template ValueConverter<VIndex>::Type& index);
242 template<
typename TreeValueType,
typename VIndexTreeType,
typename VectorValueType>
243 inline typename VIndexTreeType::template ValueConverter<TreeValueType>::Type::Ptr
246 const VIndexTreeType& index,
247 const TreeValueType& background);
254 template<
typename BoolTreeType>
257 const typename BoolTreeType::template ValueConverter<VIndex>::Type& vectorIndexTree,
259 bool staggered =
false);
281 template<
typename BoolTreeType,
typename BoundaryOp>
284 const typename BoolTreeType::template ValueConverter<VIndex>::Type& vectorIndexTree,
286 const BoundaryOp& boundaryOp,
288 bool staggered =
false);
294 template<
typename ValueType>
313 template<
typename LeafType>
318 void operator()(
const LeafType& leaf,
size_t leafIdx)
const {
319 count[leafIdx] = static_cast<VIndex>(leaf.onVoxelCount());
326 template<
typename LeafType>
332 VIndex idx = (leafIdx == 0) ? 0 : count[leafIdx - 1];
333 for (
typename LeafType::ValueOnIter it = leaf.beginValueOn(); it; ++it) {
342 template<
typename VIndexTreeType>
346 using LeafT =
typename VIndexTreeType::LeafNodeType;
350 LeafMgrT leafManager(result);
351 const size_t leafCount = leafManager.leafCount();
353 if (leafCount == 0)
return;
356 std::unique_ptr<VIndex[]> perLeafCount(
new VIndex[leafCount]);
357 VIndex* perLeafCountPtr = perLeafCount.get();
362 for (
size_t i = 1; i < leafCount; ++i) {
363 perLeafCount[i] += perLeafCount[i - 1];
367 assert(
Index64(perLeafCount[leafCount-1]) == result.activeVoxelCount());
375 template<
typename TreeType>
376 inline typename TreeType::template ValueConverter<VIndex>::Type::Ptr
379 using VIdxTreeT =
typename TreeType::template ValueConverter<VIndex>::Type;
382 const VIndex invalidIdx = -1;
383 typename VIdxTreeT::Ptr result(
387 result->voxelizeActiveTiles();
402 template<
typename VectorValueType,
typename SourceTreeType>
405 using VIdxTreeT =
typename SourceTreeType::template ValueConverter<VIndex>::Type;
407 using LeafT =
typename SourceTreeType::LeafNodeType;
419 if (
const LeafT* leaf = tree->probeLeaf(idxLeaf.origin())) {
422 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
423 vec[*it] = leaf->getValue(it.pos());
428 const TreeValueT& value = tree->getValue(idxLeaf.origin());
429 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
439 template<
typename VectorValueType,
typename SourceTreeType>
442 const typename SourceTreeType::template ValueConverter<VIndex>::Type& idxTree)
444 using VIdxTreeT =
typename SourceTreeType::template ValueConverter<VIndex>::Type;
449 const size_t numVoxels = idxTree.activeVoxelCount();
450 typename VectorT::Ptr result(
new VectorT(static_cast<math::pcg::SizeType>(numVoxels)));
454 VIdxLeafMgrT leafManager(idxTree);
468 template<
typename TreeValueType,
typename VIndexTreeType,
typename VectorValueType>
471 using OutTreeT =
typename VIndexTreeType::template ValueConverter<TreeValueType>::Type;
473 using VIdxLeafT =
typename VIndexTreeType::LeafNodeType;
484 OutLeafT* leaf = tree->probeLeaf(idxLeaf.origin());
485 assert(leaf !=
nullptr);
486 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
487 leaf->setValueOnly(it.pos(), static_cast<TreeValueType>(vec[*it]));
495 template<
typename TreeValueType,
typename VIndexTreeType,
typename VectorValueType>
496 inline typename VIndexTreeType::template ValueConverter<TreeValueType>::Type::Ptr
499 const VIndexTreeType& idxTree,
500 const TreeValueType& background)
502 using OutTreeT =
typename VIndexTreeType::template ValueConverter<TreeValueType>::Type;
506 typename OutTreeT::Ptr result(
new OutTreeT(idxTree, background,
TopologyCopy()));
507 OutTreeT& tree = *result;
511 VIdxLeafMgrT leafManager(idxTree);
525 template<
typename BoolTreeType,
typename BoundaryOp>
528 using VIdxTreeT =
typename BoolTreeType::template ValueConverter<VIndex>::Type;
540 const BoolTreeType& mask,
const BoundaryOp& op,
VectorT& src):
551 const ValueT diagonal = -6.f, offDiagonal = 1.f;
554 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
555 assert(it.getValue() > -1);
572 row.setValue(rowNum, diagonal);
584 ValueT modifiedDiagonal = 0.f;
588 row.setValue(column, offDiagonal);
589 modifiedDiagonal -= 1;
591 boundaryOp(ijk, ijk.
offsetBy(-1, 0, 0), source->at(rowNum), modifiedDiagonal);
595 row.setValue(column, offDiagonal);
596 modifiedDiagonal -= 1;
598 boundaryOp(ijk, ijk.
offsetBy(0, -1, 0), source->at(rowNum), modifiedDiagonal);
602 row.setValue(column, offDiagonal);
603 modifiedDiagonal -= 1;
605 boundaryOp(ijk, ijk.
offsetBy(0, 0, -1), source->at(rowNum), modifiedDiagonal);
609 row.setValue(column, offDiagonal);
610 modifiedDiagonal -= 1;
612 boundaryOp(ijk, ijk.
offsetBy(0, 0, 1), source->at(rowNum), modifiedDiagonal);
616 row.setValue(column, offDiagonal);
617 modifiedDiagonal -= 1;
619 boundaryOp(ijk, ijk.
offsetBy(0, 1, 0), source->at(rowNum), modifiedDiagonal);
623 row.setValue(column, offDiagonal);
624 modifiedDiagonal -= 1;
626 boundaryOp(ijk, ijk.
offsetBy(1, 0, 0), source->at(rowNum), modifiedDiagonal);
629 row.setValue(rowNum, modifiedDiagonal);
639 #define OPENVDB_TOOLS_POISSON_LAPLACIAN_STENCIL 2 642 template<
typename VIdxTreeT,
typename BoundaryOp>
655 laplacian(&m), idxTree(&idx), boundaryOp(op), source(&src) {}
661 const int kNumOffsets = 6;
662 const Coord ijkOffset[kNumOffsets] = {
663 #if OPENVDB_TOOLS_POISSON_LAPLACIAN_STENCIL == 1 664 Coord(-1,0,0),
Coord(1,0,0),
Coord(0,-1,0),
Coord(0,1,0),
Coord(0,0,-1),
Coord(0,0,1)
666 Coord(-2,0,0),
Coord(2,0,0),
Coord(0,-2,0),
Coord(0,2,0),
Coord(0,0,-2),
Coord(0,0,2)
671 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
672 assert(it.getValue() > -1);
674 const Coord ijk = it.getCoord();
679 ValueT modifiedDiagonal = 0.f;
682 for (
int dir = 0; dir < kNumOffsets; ++dir) {
683 const Coord neighbor = ijk + ijkOffset[dir];
688 #if OPENVDB_TOOLS_POISSON_LAPLACIAN_STENCIL == 1 689 const bool ijkIsInterior = (vectorIdx.
probeValue(neighbor + ijkOffset[dir], column)
692 const bool ijkIsInterior = vectorIdx.
probeValue(neighbor, column);
697 row.setValue(column, 1.f);
698 modifiedDiagonal -= 1.f;
702 boundaryOp(ijk, neighbor, source->at(rowNum), modifiedDiagonal);
706 row.setValue(rowNum, modifiedDiagonal);
714 template<
typename BoolTreeType>
715 inline LaplacianMatrix::Ptr
721 static_cast<math::pcg::SizeType>(idxTree.activeVoxelCount()));
727 template<
typename BoolTreeType,
typename BoundaryOp>
728 inline LaplacianMatrix::Ptr
730 const typename BoolTreeType::template ValueConverter<VIndex>::Type& idxTree,
732 const BoundaryOp& boundaryOp,
736 using VIdxTreeT =
typename BoolTreeType::template ValueConverter<VIndex>::Type;
740 const Index64 numDoF = idxTree.activeVoxelCount();
748 VIdxLeafMgrT idxLeafManager(idxTree);
754 laplacian, idxTree, boundaryOp, source));
764 template<
typename TreeType>
765 inline typename TreeType::Ptr
769 return solve(inTree, state, interrupter, staggered);
773 template<
typename TreeType,
typename Interrupter>
774 inline typename TreeType::Ptr
782 template<
typename TreeType,
typename BoundaryOp,
typename Interrupter>
783 inline typename TreeType::Ptr
788 return solveWithBoundaryConditionsAndPreconditioner<DefaultPrecondT>(
789 inTree, boundaryOp, state, interrupter, staggered);
794 typename PreconditionerType,
797 typename Interrupter>
798 inline typename TreeType::Ptr
800 const TreeType& inTree,
801 const BoundaryOp& boundaryOp,
803 Interrupter& interrupter,
806 return solveWithBoundaryConditionsAndPreconditioner<PreconditionerType>(
807 inTree, inTree, boundaryOp, state, interrupter, staggered);
811 typename PreconditionerType,
813 typename DomainTreeType,
815 typename Interrupter>
816 inline typename TreeType::Ptr
818 const TreeType& inTree,
819 const DomainTreeType& domainMask,
820 const BoundaryOp& boundaryOp,
822 Interrupter& interrupter,
825 using TreeValueT =
typename TreeType::ValueType;
828 using VIdxTreeT =
typename TreeType::template ValueConverter<VIndex>::Type;
829 using MaskTreeT =
typename TreeType::template ValueConverter<bool>::Type;
835 typename VectorT::Ptr b = createVectorFromTree<VecValueT>(inTree, *idxTree);
850 typename VectorT::Ptr x(
new VectorT(b->size(), zeroVal<VecValueT>()));
861 return createTreeFromVector<TreeValueT>(*x, *idxTree, zeroVal<TreeValueT>());
869 #endif // OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED Preconditioner using incomplete Cholesky factorization.
Definition: ConjGradient.h:70
uint64_t Index64
Definition: Types.h:60
bool isValueOn(const Coord &xyz) const
Return the active state of the voxel at the given coordinates.
Definition: ValueAccessor.h:264
bool probeValue(const Coord &xyz, ValueType &value) const
Return the active state of the voxel as well as its value.
Definition: ValueAccessor.h:267
Tag dispatch class that distinguishes topology copy constructors from deep copy constructors.
Definition: Types.h:518
Preconditioned conjugate gradient solver (solves Ax = b using the conjugate gradient method with one ...
Index32 SizeType
Definition: ConjGradient.h:60
SharedPtr< Preconditioner > Ptr
Definition: ConjGradient.h:495
Signed (x, y, z) 32-bit integer coordinates.
Definition: Coord.h:51
int32_t Int32
Definition: Types.h:63
Diagonal preconditioner.
Definition: ConjGradient.h:69
Lightweight, variable-length vector.
Definition: ConjGradient.h:64
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h:136
const ValueType & getValue(const Coord &xyz) const
Return the value of the voxel at the given coordinates.
Definition: ValueAccessor.h:257
Implementation of morphological dilation and erosion.
Definition: Exceptions.h:40
Information about the state of a conjugate gradient solution.
Definition: ConjGradient.h:73
virtual bool isValid() const
Definition: ConjGradient.h:500
This class manages a linear array of pointers to a given tree's leaf nodes, as well as optional auxil...
Definition: LeafManager.h:110
Dummy NOOP interrupter class defining interface.
Definition: NullInterrupter.h:52
Definition: ValueAccessor.h:220
ValueType_ ValueType
Definition: ConjGradient.h:267
Sparse, square matrix representing a 3D stencil operator of size STENCIL_SIZE.
Definition: ConjGradient.h:66
Read/write accessor to a row of this matrix.
Definition: ConjGradient.h:446
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:188
SharedPtr< SparseStencilMatrix > Ptr
Definition: ConjGradient.h:269
A LeafManager manages a linear array of pointers to a given tree's leaf nodes, as well as optional au...
SharedPtr< Vector > Ptr
Definition: ConjGradient.h:169
Coord offsetBy(Int32 dx, Int32 dy, Int32 dz) const
Definition: Coord.h:118