88 #ifndef OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED 89 #define OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED 97 #include <boost/scoped_array.hpp> 114 template<
typename TreeType>
127 inline typename TreeType::Ptr
130 template<
typename TreeType,
typename Interrupter>
131 inline typename TreeType::Ptr
137 template<
typename TreeType,
typename BoundaryOp,
typename Interrupter>
168 inline typename TreeType::Ptr
174 bool staggered =
false);
177 typename PreconditionerType,
180 typename Interrupter>
181 inline typename TreeType::Ptr
187 bool staggered =
false);
190 typename PreconditionerType,
192 typename DomainTreeType,
194 typename Interrupter>
195 inline typename TreeType::Ptr
198 const DomainTreeType&,
202 bool staggered =
false);
212 template<
typename VIndexTreeType>
218 template<
typename TreeType>
219 inline typename TreeType::template ValueConverter<VIndex>::Type::Ptr
229 template<
typename VectorValueType,
typename SourceTreeType>
232 const SourceTreeType& source,
233 const typename SourceTreeType::template ValueConverter<VIndex>::Type& index);
243 template<
typename TreeValueType,
typename VIndexTreeType,
typename VectorValueType>
244 inline typename VIndexTreeType::template ValueConverter<TreeValueType>::Type::Ptr
247 const VIndexTreeType& index,
248 const TreeValueType& background);
255 template<
typename BoolTreeType>
258 const typename BoolTreeType::template ValueConverter<VIndex>::Type& vectorIndexTree,
260 bool staggered =
false);
282 template<
typename BoolTreeType,
typename BoundaryOp>
285 const typename BoolTreeType::template ValueConverter<VIndex>::Type& vectorIndexTree,
287 const BoundaryOp& boundaryOp,
289 bool staggered =
false);
295 template<
typename ValueType>
314 template<
typename LeafType>
319 void operator()(
const LeafType& leaf,
size_t leafIdx)
const {
320 count[leafIdx] =
static_cast<VIndex
>(leaf.onVoxelCount());
327 template<
typename LeafType>
333 VIndex idx = (leafIdx == 0) ? 0 : count[leafIdx - 1];
334 for (
typename LeafType::ValueOnIter it = leaf.beginValueOn(); it; ++it) {
343 template<
typename VIndexTreeType>
347 typedef typename VIndexTreeType::LeafNodeType LeafT;
351 LeafMgrT leafManager(result);
352 const size_t leafCount = leafManager.leafCount();
354 if (leafCount == 0)
return;
357 boost::scoped_array<VIndex> perLeafCount(
new VIndex[leafCount]);
358 VIndex* perLeafCountPtr = perLeafCount.get();
363 for (
size_t i = 1; i < leafCount; ++i) {
364 perLeafCount[i] += perLeafCount[i - 1];
368 assert(
Index64(perLeafCount[leafCount-1]) == result.activeVoxelCount());
376 template<
typename TreeType>
377 inline typename TreeType::template ValueConverter<VIndex>::Type::Ptr
380 typedef typename TreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
383 const VIndex invalidIdx = -1;
384 typename VIdxTreeT::Ptr result(
388 result->voxelizeActiveTiles();
403 template<
typename VectorValueType,
typename SourceTreeType>
406 typedef typename SourceTreeType::template ValueConverter<VIndex>::Type
VIdxTreeT;
408 typedef typename SourceTreeType::LeafNodeType
LeafT;
415 CopyToVecOp(
const SourceTreeType& t, VectorT& v): tree(&t), vector(&v) {}
419 VectorT& vec = *vector;
420 if (
const LeafT* leaf = tree->probeLeaf(idxLeaf.origin())) {
423 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
424 vec[*it] = leaf->getValue(it.pos());
429 const TreeValueT& value = tree->getValue(idxLeaf.origin());
430 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
440 template<
typename VectorValueType,
typename SourceTreeType>
443 const typename SourceTreeType::template ValueConverter<VIndex>::Type& idxTree)
445 typedef typename SourceTreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
450 const size_t numVoxels = idxTree.activeVoxelCount();
451 typename VectorT::Ptr result(
new VectorT(static_cast<math::pcg::SizeType>(numVoxels)));
455 VIdxLeafMgrT leafManager(idxTree);
469 template<
typename TreeValueType,
typename VIndexTreeType,
typename VectorValueType>
472 typedef typename VIndexTreeType::template ValueConverter<TreeValueType>::Type
OutTreeT;
474 typedef typename VIndexTreeType::LeafNodeType
VIdxLeafT;
484 const VectorT& vec = *vector;
485 OutLeafT* leaf = tree->probeLeaf(idxLeaf.origin());
486 assert(leaf != NULL);
487 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
488 leaf->setValueOnly(it.pos(),
static_cast<TreeValueType
>(vec[*it]));
496 template<
typename TreeValueType,
typename VIndexTreeType,
typename VectorValueType>
497 inline typename VIndexTreeType::template ValueConverter<TreeValueType>::Type::Ptr
500 const VIndexTreeType& idxTree,
501 const TreeValueType& background)
503 typedef typename VIndexTreeType::template ValueConverter<TreeValueType>::Type OutTreeT;
507 typename OutTreeT::Ptr result(
new OutTreeT(idxTree, background,
TopologyCopy()));
508 OutTreeT& tree = *result;
512 VIdxLeafMgrT leafManager(idxTree);
526 template<
typename BoolTreeType,
typename BoundaryOp>
529 typedef typename BoolTreeType::template ValueConverter<VIndex>::Type
VIdxTreeT;
541 const BoolTreeType& mask,
const BoundaryOp& op, VectorT& src):
542 laplacian(&m), idxTree(&idx), interiorMask(&mask), boundaryOp(op), source(&src) {}
552 const ValueT diagonal = -6.f, offDiagonal = 1.f;
555 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
556 assert(it.getValue() > -1);
585 ValueT modifiedDiagonal = 0.f;
590 modifiedDiagonal -= 1;
592 boundaryOp(ijk, ijk.
offsetBy(-1, 0, 0), source->
at(rowNum), modifiedDiagonal);
597 modifiedDiagonal -= 1;
599 boundaryOp(ijk, ijk.
offsetBy(0, -1, 0), source->
at(rowNum), modifiedDiagonal);
604 modifiedDiagonal -= 1;
606 boundaryOp(ijk, ijk.
offsetBy(0, 0, -1), source->
at(rowNum), modifiedDiagonal);
611 modifiedDiagonal -= 1;
613 boundaryOp(ijk, ijk.
offsetBy(0, 0, 1), source->
at(rowNum), modifiedDiagonal);
618 modifiedDiagonal -= 1;
620 boundaryOp(ijk, ijk.
offsetBy(0, 1, 0), source->
at(rowNum), modifiedDiagonal);
625 modifiedDiagonal -= 1;
627 boundaryOp(ijk, ijk.
offsetBy(1, 0, 0), source->
at(rowNum), modifiedDiagonal);
630 row.
setValue(rowNum, modifiedDiagonal);
640 #define OPENVDB_TOOLS_POISSON_LAPLACIAN_STENCIL 2 643 template<
typename VIdxTreeT,
typename BoundaryOp>
655 ISLaplacianOp(LaplacianMatrix& m,
const VIdxTreeT& idx,
const BoundaryOp& op, VectorT& src):
656 laplacian(&m), idxTree(&idx), boundaryOp(op), source(&src) {}
662 const int kNumOffsets = 6;
663 const Coord ijkOffset[kNumOffsets] = {
664 #if OPENVDB_TOOLS_POISSON_LAPLACIAN_STENCIL == 1 665 Coord(-1,0,0),
Coord(1,0,0),
Coord(0,-1,0),
Coord(0,1,0),
Coord(0,0,-1),
Coord(0,0,1)
667 Coord(-2,0,0),
Coord(2,0,0),
Coord(0,-2,0),
Coord(0,2,0),
Coord(0,0,-2),
Coord(0,0,2)
672 for (
typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
673 assert(it.getValue() > -1);
675 const Coord ijk = it.getCoord();
680 ValueT modifiedDiagonal = 0.f;
683 for (
int dir = 0; dir < kNumOffsets; ++dir) {
684 const Coord neighbor = ijk + ijkOffset[dir];
689 #if OPENVDB_TOOLS_POISSON_LAPLACIAN_STENCIL == 1 690 const bool ijkIsInterior = (vectorIdx.
probeValue(neighbor + ijkOffset[dir], column)
693 const bool ijkIsInterior = vectorIdx.
probeValue(neighbor, column);
699 modifiedDiagonal -= 1.f;
703 boundaryOp(ijk, neighbor, source->
at(rowNum), modifiedDiagonal);
707 row.
setValue(rowNum, modifiedDiagonal);
715 template<
typename BoolTreeType>
722 static_cast<math::pcg::SizeType>(idxTree.activeVoxelCount()));
728 template<
typename BoolTreeType,
typename BoundaryOp>
731 const typename BoolTreeType::template ValueConverter<VIndex>::Type& idxTree,
733 const BoundaryOp& boundaryOp,
737 typedef typename BoolTreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
741 const Index64 numDoF = idxTree.activeVoxelCount();
746 LaplacianMatrix&
laplacian = *laplacianPtr;
749 VIdxLeafMgrT idxLeafManager(idxTree);
752 laplacian, idxTree, interiorMask, boundaryOp, source));
755 laplacian, idxTree, boundaryOp, source));
765 template<
typename TreeType>
766 inline typename TreeType::Ptr
770 return solve(inTree, state, interrupter, staggered);
774 template<
typename TreeType,
typename Interrupter>
775 inline typename TreeType::Ptr
783 template<
typename TreeType,
typename BoundaryOp,
typename Interrupter>
784 inline typename TreeType::Ptr
789 return solveWithBoundaryConditionsAndPreconditioner<DefaultPrecondT>(
790 inTree, boundaryOp, state, interrupter, staggered);
795 typename PreconditionerType,
798 typename Interrupter>
799 inline typename TreeType::Ptr
801 const TreeType& inTree,
802 const BoundaryOp& boundaryOp,
804 Interrupter& interrupter,
807 return solveWithBoundaryConditionsAndPreconditioner<PreconditionerType>(
808 inTree, inTree, boundaryOp, state, interrupter, staggered);
812 typename PreconditionerType,
814 typename DomainTreeType,
816 typename Interrupter>
817 inline typename TreeType::Ptr
819 const TreeType& inTree,
820 const DomainTreeType& domainMask,
821 const BoundaryOp& boundaryOp,
823 Interrupter& interrupter,
826 typedef typename TreeType::ValueType TreeValueT;
829 typedef typename TreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
830 typedef typename TreeType::template ValueConverter<bool>::Type MaskTreeT;
836 typename VectorT::Ptr b = createVectorFromTree<VecValueT>(inTree, *idxTree);
846 *idxTree, *interiorMask, boundaryOp, *b, staggered);
849 laplacian->scale(-1.0);
851 typename VectorT::Ptr x(
new VectorT(b->size(), zeroVal<VecValueT>()));
853 new PreconditionerType(*laplacian));
854 if (!precond->isValid()) {
862 return createTreeFromVector<TreeValueT>(*x, *idxTree, zeroVal<TreeValueT>());
870 #endif // OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED Dummy NOOP interrupter class defining interface.
Definition: NullInterrupter.h:52
uint64_t Index64
Definition: Types.h:59
bool isValueOn(const Coord &xyz) const
Return the active state of the voxel at the given coordinates.
Definition: ValueAccessor.h:263
bool probeValue(const Coord &xyz, ValueType &value) const
Return the active state of the voxel as well as its value.
Definition: ValueAccessor.h:266
Preconditioned conjugate gradient solver (solves Ax = b using the conjugate gradient method with one ...
SharedPtr< SparseStencilMatrix > Ptr
Definition: ConjGradient.h:269
Definition: ValueAccessor.h:219
Signed (x, y, z) 32-bit integer coordinates.
Definition: Coord.h:51
SharedPtr< Vector > Ptr
Definition: ConjGradient.h:169
Preconditioner using incomplete Cholesky factorization.
Definition: ConjGradient.h:70
RowEditor getRowEditor(SizeType row)
Return a read/write view onto the given row of this matrix.
Definition: ConjGradient.h:1096
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h:136
SharedPtr< Preconditioner > Ptr
Definition: ConjGradient.h:495
Diagonal preconditioner.
Definition: ConjGradient.h:69
Implementation of morphological dilation and erosion.
const ValueType & getValue(const Coord &xyz) const
Return the value of the voxel at the given coordinates.
Definition: ValueAccessor.h:256
Definition: Exceptions.h:39
Lightweight, variable-length vector.
Definition: ConjGradient.h:64
This class manages a linear array of pointers to a given tree's leaf nodes, as well as optional auxil...
Definition: LeafManager.h:110
Coord offsetBy(Int32 dx, Int32 dy, Int32 dz) const
Definition: Coord.h:118
Index32 SizeType
Definition: ConjGradient.h:60
Information about the state of a conjugate gradient solution.
Definition: ConjGradient.h:73
Tag dispatch class that distinguishes topology copy constructors from deep copy constructors.
Definition: Types.h:508
T & at(SizeType i)
Return the value of this vector's ith element.
Definition: ConjGradient.h:228
Read/write accessor to a row of this matrix.
Definition: ConjGradient.h:446
SizeType setValue(SizeType column, const ValueType &value)
Set the value of the entry in the specified column.
Definition: ConjGradient.h:1252
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:188
A LeafManager manages a linear array of pointers to a given tree's leaf nodes, as well as optional au...
Sparse, square matrix representing a 3D stencil operator of size STENCIL_SIZE.
Definition: ConjGradient.h:66
ValueType_ ValueType
Definition: ConjGradient.h:267
int32_t Int32
Definition: Types.h:62