OpenVDB  3.1.0
PoissonSolver.h
Go to the documentation of this file.
1 //
3 // Copyright (c) 2012-2015 DreamWorks Animation LLC
4 //
5 // All rights reserved. This software is distributed under the
6 // Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ )
7 //
8 // Redistributions of source code must retain the above copyright
9 // and license notice and the following restrictions and disclaimer.
10 //
11 // * Neither the name of DreamWorks Animation nor the names of
12 // its contributors may be used to endorse or promote products derived
13 // from this software without specific prior written permission.
14 //
15 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
16 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
17 // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
18 // A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
19 // OWNER OR CONTRIBUTORS BE LIABLE FOR ANY INDIRECT, INCIDENTAL,
20 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
21 // LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
22 // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
23 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
24 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
25 // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26 // IN NO EVENT SHALL THE COPYRIGHT HOLDERS' AND CONTRIBUTORS' AGGREGATE
27 // LIABILITY FOR ALL CLAIMS REGARDLESS OF THEIR BASIS EXCEED US$250.00.
28 //
30 //
88 
89 #ifndef OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED
90 #define OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED
91 
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>
97 
98 #include "Morphology.h" // for erodeVoxels
99 
100 #include <boost/scoped_array.hpp>
101 
102 
103 namespace openvdb {
105 namespace OPENVDB_VERSION_NAME {
106 namespace tools {
107 namespace poisson {
108 
109 // This type should be at least as wide as math::pcg::SizeType.
110 typedef Int32 VIndex;
111 
114 
115 
117 template<typename TreeType>
130 inline typename TreeType::Ptr
131 solve(const TreeType&, math::pcg::State&);
132 
133 template<typename TreeType, typename Interrupter>
134 inline typename TreeType::Ptr
135 solve(const TreeType&, math::pcg::State&, Interrupter&);
137 
138 
140 template<typename TreeType, typename BoundaryOp, typename Interrupter>
171 inline typename TreeType::Ptr
172 solveWithBoundaryConditions(const TreeType&, const BoundaryOp&, math::pcg::State&, Interrupter&);
173 
174 template<typename PreconditionerType, typename TreeType, typename BoundaryOp, typename Interrupter>
175 inline typename TreeType::Ptr
176 solveWithBoundaryConditionsAndPreconditioner(const TreeType&, const BoundaryOp&,
177  math::pcg::State&, Interrupter&);
178 
179 template<typename PreconditionerType, typename TreeType, typename DomainTreeType, typename BoundaryOp, typename Interrupter>
180 inline typename TreeType::Ptr
181 solveWithBoundaryConditionsAndPreconditioner(const TreeType&, const DomainTreeType&, const BoundaryOp&,
182  math::pcg::State&, Interrupter&);
184 
185 
187 
188 // The following are low-level routines that can be used to assemble custom solvers.
189 
192 template<typename VIndexTreeType>
193 inline void populateIndexTree(VIndexTreeType&);
194 
198 template<typename TreeType>
199 inline typename TreeType::template ValueConverter<VIndex>::Type::Ptr
200 createIndexTree(const TreeType&);
201 
202 
209 template<typename VectorValueType, typename SourceTreeType>
212  const SourceTreeType& source,
213  const typename SourceTreeType::template ValueConverter<VIndex>::Type& index);
214 
215 
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);
229 
230 
235 template<typename BoolTreeType>
238  const typename BoolTreeType::template ValueConverter<VIndex>::Type& vectorIndexTree,
239  const BoolTreeType& interiorMask);
240 
241 
261 template<typename BoolTreeType, typename BoundaryOp>
264  const typename BoolTreeType::template ValueConverter<VIndex>::Type& vectorIndexTree,
265  const BoolTreeType& interiorMask,
266  const BoundaryOp& boundaryOp,
268 
270 
271 
273 
274 
275 namespace internal {
276 
279 template<typename LeafType>
281 {
282  VIndex* count;
283  LeafCountOp(VIndex* count_): count(count_) {}
284  void operator()(const LeafType& leaf, size_t leafIdx) const {
285  count[leafIdx] = static_cast<VIndex>(leaf.onVoxelCount());
286  }
287 };
288 
289 
292 template<typename LeafType>
294 {
295  const VIndex* count;
296  LeafIndexOp(const VIndex* count_): count(count_) {}
297  void operator()(LeafType& leaf, size_t leafIdx) const {
298  VIndex idx = (leafIdx == 0) ? 0 : count[leafIdx - 1];
299  for (typename LeafType::ValueOnIter it = leaf.beginValueOn(); it; ++it) {
300  it.setValue(idx++);
301  }
302  }
303 };
304 
305 } // namespace internal
306 
307 
308 template<typename VIndexTreeType>
309 inline void
310 populateIndexTree(VIndexTreeType& result)
311 {
312  typedef typename VIndexTreeType::LeafNodeType LeafT;
313  typedef typename tree::LeafManager<VIndexTreeType> LeafMgrT;
314 
315  // Linearize the tree.
316  LeafMgrT leafManager(result);
317  const size_t leafCount = leafManager.leafCount();
318 
319  // Count the number of active voxels in each leaf node.
320  boost::scoped_array<VIndex> perLeafCount(new VIndex[leafCount]);
321  VIndex* perLeafCountPtr = perLeafCount.get();
322  leafManager.foreach(internal::LeafCountOp<LeafT>(perLeafCountPtr));
323 
324  // The starting index for each leaf node is the total number
325  // of active voxels in all preceding leaf nodes.
326  for (size_t i = 1; i < leafCount; ++i) {
327  perLeafCount[i] += perLeafCount[i - 1];
328  }
329 
330  // The last accumulated value should be the total of all active voxels.
331  assert(Index64(perLeafCount[leafCount-1]) == result.activeVoxelCount());
332 
333  // Parallelize over the leaf nodes of the tree, storing a unique index
334  // in each active voxel.
335  leafManager.foreach(internal::LeafIndexOp<LeafT>(perLeafCountPtr));
336 }
337 
338 
339 template<typename TreeType>
340 inline typename TreeType::template ValueConverter<VIndex>::Type::Ptr
341 createIndexTree(const TreeType& tree)
342 {
343  typedef typename TreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
344 
345  // Construct an output tree with the same active voxel topology as the input tree.
346  const VIndex invalidIdx = -1;
347  typename VIdxTreeT::Ptr result(
348  new VIdxTreeT(tree, /*background=*/invalidIdx, TopologyCopy()));
349 
350  // All active voxels are degrees of freedom, including voxels contained in active tiles.
351  result->voxelizeActiveTiles();
352 
353  populateIndexTree(*result);
354 
355  return result;
356 }
357 
358 
360 
361 
362 namespace internal {
363 
366 template<typename VectorValueType, typename SourceTreeType>
368 {
369  typedef typename SourceTreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
370  typedef typename VIdxTreeT::LeafNodeType VIdxLeafT;
371  typedef typename SourceTreeType::LeafNodeType LeafT;
372  typedef typename SourceTreeType::ValueType TreeValueT;
374 
375  const SourceTreeType* tree;
376  VectorT* vector;
377 
378  CopyToVecOp(const SourceTreeType& t, VectorT& v): tree(&t), vector(&v) {}
379 
380  void operator()(const VIdxLeafT& idxLeaf, size_t /*leafIdx*/) const
381  {
382  VectorT& vec = *vector;
383  if (const LeafT* leaf = tree->probeLeaf(idxLeaf.origin())) {
384  // If a corresponding leaf node exists in the source tree,
385  // copy voxel values from the source node to the output vector.
386  for (typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
387  vec[*it] = leaf->getValue(it.pos());
388  }
389  } else {
390  // If no corresponding leaf exists in the source tree,
391  // fill the vector with a uniform value.
392  const TreeValueT& value = tree->getValue(idxLeaf.origin());
393  for (typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
394  vec[*it] = value;
395  }
396  }
397  }
398 };
399 
400 } // namespace internal
401 
402 
403 template<typename VectorValueType, typename SourceTreeType>
405 createVectorFromTree(const SourceTreeType& tree,
406  const typename SourceTreeType::template ValueConverter<VIndex>::Type& idxTree)
407 {
408  typedef typename SourceTreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
409  typedef tree::LeafManager<const VIdxTreeT> VIdxLeafMgrT;
410  typedef typename math::pcg::Vector<VectorValueType> VectorT;
411 
412  // Allocate the vector.
413  const size_t numVoxels = idxTree.activeVoxelCount();
414  typename VectorT::Ptr result(new VectorT(static_cast<math::pcg::SizeType>(numVoxels)));
415 
416  // Parallelize over the leaf nodes of the index tree, filling the output vector
417  // with values from corresponding voxels of the source tree.
418  VIdxLeafMgrT leafManager(idxTree);
419  leafManager.foreach(internal::CopyToVecOp<VectorValueType, SourceTreeType>(tree, *result));
420 
421  return result;
422 }
423 
424 
426 
427 
428 namespace internal {
429 
432 template<typename TreeValueType, typename VIndexTreeType, typename VectorValueType>
434 {
435  typedef typename VIndexTreeType::template ValueConverter<TreeValueType>::Type OutTreeT;
436  typedef typename OutTreeT::LeafNodeType OutLeafT;
437  typedef typename VIndexTreeType::LeafNodeType VIdxLeafT;
439 
440  const VectorT* vector;
441  OutTreeT* tree;
442 
443  CopyFromVecOp(const VectorT& v, OutTreeT& t): vector(&v), tree(&t) {}
444 
445  void operator()(const VIdxLeafT& idxLeaf, size_t /*leafIdx*/) const
446  {
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]));
452  }
453  }
454 };
455 
456 } // namespace internal
457 
458 
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)
465 {
466  typedef typename VIndexTreeType::template ValueConverter<TreeValueType>::Type OutTreeT;
467  typedef typename tree::LeafManager<const VIndexTreeType> VIdxLeafMgrT;
468 
469  // Construct an output tree with the same active voxel topology as the index tree.
470  typename OutTreeT::Ptr result(new OutTreeT(idxTree, background, TopologyCopy()));
471  OutTreeT& tree = *result;
472 
473  // Parallelize over the leaf nodes of the index tree, populating voxels
474  // of the output tree with values from the input vector.
475  VIdxLeafMgrT leafManager(idxTree);
476  leafManager.foreach(
478 
479  return result;
480 }
481 
482 
484 
485 
486 namespace internal {
487 
489 template<typename ValueType>
490 struct DirichletOp {
491  inline void operator()(
492  const Coord&, const Coord&, ValueType&, ValueType& diag) const { diag -= 1; }
493 };
494 
495 
497 template<typename BoolTreeType, typename BoundaryOp>
499 {
500  typedef typename BoolTreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
501  typedef typename VIdxTreeT::LeafNodeType VIdxLeafT;
504 
505  LaplacianMatrix* laplacian;
506  const VIdxTreeT* idxTree;
507  const BoolTreeType* interiorMask;
508  const BoundaryOp boundaryOp;
509  VectorT* source;
510 
511  ISLaplacianOp(LaplacianMatrix& m, const VIdxTreeT& idx,
512  const BoolTreeType& mask, const BoundaryOp& op, VectorT& src):
513  laplacian(&m), idxTree(&idx), interiorMask(&mask), boundaryOp(op), source(&src) {}
514 
515  void operator()(const VIdxLeafT& idxLeaf, size_t /*leafIdx*/) const
516  {
517  // Local accessors
518  typename tree::ValueAccessor<const BoolTreeType> interior(*interiorMask);
519  typename tree::ValueAccessor<const VIdxTreeT> vectorIdx(*idxTree);
520 
521  Coord ijk;
522  VIndex column;
523  const ValueT diagonal = -6.f, offDiagonal = 1.f;
524 
525  // Loop over active voxels in this leaf.
526  for (typename VIdxLeafT::ValueOnCIter it = idxLeaf.cbeginValueOn(); it; ++it) {
527  assert(it.getValue() > -1);
528  const math::pcg::SizeType rowNum = static_cast<math::pcg::SizeType>(it.getValue());
529 
530  LaplacianMatrix::RowEditor row = laplacian->getRowEditor(rowNum);
531 
532  ijk = it.getCoord();
533  if (interior.isValueOn(ijk)) {
534  // The current voxel is an interior voxel.
535  // All of its neighbors are in the solution domain.
536 
537  // -x direction
538  row.setValue(vectorIdx.getValue(ijk.offsetBy(-1, 0, 0)), offDiagonal);
539  // -y direction
540  row.setValue(vectorIdx.getValue(ijk.offsetBy(0, -1, 0)), offDiagonal);
541  // -z direction
542  row.setValue(vectorIdx.getValue(ijk.offsetBy(0, 0, -1)), offDiagonal);
543  // diagonal
544  row.setValue(rowNum, diagonal);
545  // +z direction
546  row.setValue(vectorIdx.getValue(ijk.offsetBy(0, 0, 1)), offDiagonal);
547  // +y direction
548  row.setValue(vectorIdx.getValue(ijk.offsetBy(0, 1, 0)), offDiagonal);
549  // +x direction
550  row.setValue(vectorIdx.getValue(ijk.offsetBy(1, 0, 0)), offDiagonal);
551 
552  } else {
553  // The current voxel is a boundary voxel.
554  // At least one of its neighbors is outside the solution domain.
555 
556  ValueT modifiedDiagonal = 0.f;
557 
558  // -x direction
559  if (vectorIdx.probeValue(ijk.offsetBy(-1, 0, 0), column)) {
560  row.setValue(column, offDiagonal);
561  modifiedDiagonal -= 1;
562  } else {
563  boundaryOp(ijk, ijk.offsetBy(-1, 0, 0), source->at(rowNum), modifiedDiagonal);
564  }
565  // -y direction
566  if (vectorIdx.probeValue(ijk.offsetBy(0, -1, 0), column)) {
567  row.setValue(column, offDiagonal);
568  modifiedDiagonal -= 1;
569  } else {
570  boundaryOp(ijk, ijk.offsetBy(0, -1, 0), source->at(rowNum), modifiedDiagonal);
571  }
572  // -z direction
573  if (vectorIdx.probeValue(ijk.offsetBy(0, 0, -1), column)) {
574  row.setValue(column, offDiagonal);
575  modifiedDiagonal -= 1;
576  } else {
577  boundaryOp(ijk, ijk.offsetBy(0, 0, -1), source->at(rowNum), modifiedDiagonal);
578  }
579  // +z direction
580  if (vectorIdx.probeValue(ijk.offsetBy(0, 0, 1), column)) {
581  row.setValue(column, offDiagonal);
582  modifiedDiagonal -= 1;
583  } else {
584  boundaryOp(ijk, ijk.offsetBy(0, 0, 1), source->at(rowNum), modifiedDiagonal);
585  }
586  // +y direction
587  if (vectorIdx.probeValue(ijk.offsetBy(0, 1, 0), column)) {
588  row.setValue(column, offDiagonal);
589  modifiedDiagonal -= 1;
590  } else {
591  boundaryOp(ijk, ijk.offsetBy(0, 1, 0), source->at(rowNum), modifiedDiagonal);
592  }
593  // +x direction
594  if (vectorIdx.probeValue(ijk.offsetBy(1, 0, 0), column)) {
595  row.setValue(column, offDiagonal);
596  modifiedDiagonal -= 1;
597  } else {
598  boundaryOp(ijk, ijk.offsetBy(1, 0, 0), source->at(rowNum), modifiedDiagonal);
599  }
600  // diagonal
601  row.setValue(rowNum, modifiedDiagonal);
602  }
603  } // end loop over voxels
604  }
605 };
606 
607 } // namespace internal
608 
609 
610 template<typename BoolTreeType>
611 inline LaplacianMatrix::Ptr
612 createISLaplacian(const typename BoolTreeType::template ValueConverter<VIndex>::Type& idxTree,
613  const BoolTreeType& interiorMask)
614 {
615  typedef LaplacianMatrix::ValueType ValueT;
617  static_cast<math::pcg::SizeType>(idxTree.activeVoxelCount()));
619  return createISLaplacianWithBoundaryConditions(idxTree, interiorMask, op, unused);
620 }
621 
622 
623 template<typename BoolTreeType, typename BoundaryOp>
624 inline LaplacianMatrix::Ptr
626  const typename BoolTreeType::template ValueConverter<VIndex>::Type& idxTree,
627  const BoolTreeType& interiorMask,
628  const BoundaryOp& boundaryOp,
630 {
631  typedef typename BoolTreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
632  typedef typename tree::LeafManager<const VIdxTreeT> VIdxLeafMgrT;
633 
634  // The number of active voxels is the number of degrees of freedom.
635  const Index64 numDoF = idxTree.activeVoxelCount();
636 
637  // Construct the matrix.
638  LaplacianMatrix::Ptr laplacianPtr(
639  new LaplacianMatrix(static_cast<math::pcg::SizeType>(numDoF)));
640  LaplacianMatrix& laplacian = *laplacianPtr;
641 
642  // Populate the matrix using a second-order, 7-point CD stencil.
643  VIdxLeafMgrT idxLeafManager(idxTree);
645  laplacian, idxTree, interiorMask, boundaryOp, source));
646 
647  return laplacianPtr;
648 }
649 
650 
652 
653 
654 template<typename TreeType>
655 inline typename TreeType::Ptr
656 solve(const TreeType& inTree, math::pcg::State& state)
657 {
658  util::NullInterrupter interrupter;
659  return solve(inTree, state, interrupter);
660 }
661 
662 
663 template<typename TreeType, typename Interrupter>
664 inline typename TreeType::Ptr
665 solve(const TreeType& inTree, math::pcg::State& state, Interrupter& interrupter)
666 {
668  return solveWithBoundaryConditions(inTree, boundaryOp, state, interrupter);
669 }
670 
671 
672 template<typename TreeType, typename BoundaryOp, typename Interrupter>
673 inline typename TreeType::Ptr
674 solveWithBoundaryConditions(const TreeType& inTree, const BoundaryOp& boundaryOp,
675  math::pcg::State& state, Interrupter& interrupter)
676 {
678  return solveWithBoundaryConditionsAndPreconditioner<DefaultPrecondT>(
679  inTree, boundaryOp, state, interrupter);
680 }
681 
682 
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)
687 {
688 
689  return solveWithBoundaryConditionsAndPreconditioner<PreconditionerType>(inTree /*source*/, inTree /*domain mask*/,
690  boundaryOp, state, interrupter);
691 }
692 
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,
698  math::pcg::State& state, Interrupter& interrupter)
699 {
700 
701  typedef typename TreeType::ValueType TreeValueT;
702  typedef LaplacianMatrix::ValueType VecValueT;
703  typedef typename math::pcg::Vector<VecValueT> VectorT;
704  typedef typename TreeType::template ValueConverter<VIndex>::Type VIdxTreeT;
705  typedef typename TreeType::template ValueConverter<bool>::Type MaskTreeT;
706 
707  // 1. Create a mapping from active voxels of the input tree to elements of a vector.
708  typename VIdxTreeT::ConstPtr idxTree = createIndexTree(domainMask);
709 
710  // 2. Populate a vector with values from the input tree.
711  typename VectorT::Ptr b = createVectorFromTree<VecValueT>(inTree, *idxTree);
712 
713  // 3. Create a mask of the interior voxels of the input tree (from the densified index tree).
714  typename MaskTreeT::Ptr interiorMask(
715  new MaskTreeT(*idxTree, /*background=*/false, TopologyCopy()));
716  tools::erodeVoxels(*interiorMask, /*iterations=*/1, tools::NN_FACE);
717 
718  // 4. Create the Laplacian matrix.
720  *idxTree, *interiorMask, boundaryOp, *b);
721 
722  // 5. Solve the Poisson equation.
723  laplacian->scale(-1.0); // matrix is negative-definite; solve -M x = -b
724  b->scale(-1.0);
725  typename VectorT::Ptr x(new VectorT(b->size(), zeroVal<VecValueT>()));
727  new PreconditionerType(*laplacian));
728  if (!precond->isValid()) {
729  precond.reset(new math::pcg::JacobiPreconditioner<LaplacianMatrix>(*laplacian));
730  }
731 
732  state = math::pcg::solve(*laplacian, *b, *x, *precond, interrupter, state);
733 
734  // 6. Populate the output tree with values from the solution vector.
736  return createTreeFromVector<TreeValueT>(*x, *idxTree, /*background=*/zeroVal<TreeValueT>());
737 }
738 
739 } // namespace poisson
740 } // namespace tools
741 } // namespace OPENVDB_VERSION_NAME
742 } // namespace openvdb
743 
744 #endif // OPENVDB_TOOLS_POISSONSOLVER_HAS_BEEN_INCLUDED
745 
747 //
748 // Copyright (c) 2012-2015 DreamWorks Animation LLC
749 //
750 // All rights reserved. This software is distributed under the
751 // Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ )
752 //
753 // Redistributions of source code must retain the above copyright
754 // and license notice and the following restrictions and disclaimer.
755 //
756 // * Neither the name of DreamWorks Animation nor the names of
757 // its contributors may be used to endorse or promote products derived
758 // from this software without specific prior written permission.
759 //
760 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
761 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
762 // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
763 // A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
764 // OWNER OR CONTRIBUTORS BE LIABLE FOR ANY INDIRECT, INCIDENTAL,
765 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
766 // LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
767 // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
768 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
769 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
770 // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
771 // IN NO EVENT SHALL THE COPYRIGHT HOLDERS' AND CONTRIBUTORS' AGGREGATE
772 // LIABILITY FOR ALL CLAIMS REGARDLESS OF THEIR BASIS EXCEED US$250.00.
773 //
SourceTreeType::LeafNodeType LeafT
Definition: PoissonSolver.h:371
void operator()(LeafType &leaf, size_t leafIdx) const
Definition: PoissonSolver.h:297
SourceTreeType::ValueType TreeValueT
Definition: PoissonSolver.h:372
boost::shared_ptr< SparseStencilMatrix > Ptr
Definition: ConjGradient.h:271
Functor for use with LeafManager::foreach() to populate active leaf voxels with sequential indices...
Definition: PoissonSolver.h:293
TreeType::Ptr solveWithBoundaryConditions(const TreeType &, const BoundaryOp &, math::pcg::State &, Interrupter &)
Solve ∇2x = b for x with user-specified boundary conditions, where b is a vector comprising the valu...
Definition: PoissonSolver.h:674
VIndexTreeType::LeafNodeType VIdxLeafT
Definition: PoissonSolver.h:437
TreeType::template ValueConverter< VIndex >::Type::Ptr createIndexTree(const TreeType &)
Iterate over the active voxels of the input tree and for each one assign its index in the iteration s...
Definition: PoissonSolver.h:341
Read/write accessor to a row of this matrix.
Definition: ConjGradient.h:448
LeafCountOp(VIndex *count_)
Definition: PoissonSolver.h:283
Preconditioner using incomplete Cholesky factorization.
Definition: ConjGradient.h:71
VectorT * source
Definition: PoissonSolver.h:509
Functor for use with LeafManager::foreach() to populate a vector with the values of a tree's active v...
Definition: PoissonSolver.h:367
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
CopyFromVecOp(const VectorT &v, OutTreeT &t)
Definition: PoissonSolver.h:443
VectorT * vector
Definition: PoissonSolver.h:376
void operator()(const VIdxLeafT &idxLeaf, size_t) const
Definition: PoissonSolver.h:445
ValueType_ ValueType
Definition: ConjGradient.h:269
GridType::Ptr laplacian(const GridType &grid, bool threaded, InterruptT *interrupt)
Compute the Laplacian of the given scalar grid.
Definition: GridOperators.h:1016
Index32 SizeType
Definition: ConjGradient.h:61
Functor for use with LeafManager::foreach() to populate a sparse Laplacian matrix.
Definition: PoissonSolver.h:498
LaplacianMatrix::Ptr createISLaplacianWithBoundaryConditions(const typename BoolTreeType::template ValueConverter< VIndex >::Type &vectorIndexTree, const BoolTreeType &interiorMask, const BoundaryOp &boundaryOp, typename math::pcg::Vector< LaplacianMatrix::ValueType > &source)
Generate a sparse matrix of the index-space (Δx = 1) Laplacian operator with user-specified boundary ...
Definition: PoissonSolver.h:625
math::pcg::Vector< VectorValueType > VectorT
Definition: PoissonSolver.h:438
OPENVDB_STATIC_SPECIALIZATION void erodeVoxels(TreeType &tree, int iterations=1, NearestNeighbors nn=NN_FACE)
Topologically erode all leaf-level active voxels in the given tree.
Definition: Morphology.h:798
TreeType::Ptr solveWithBoundaryConditionsAndPreconditioner(const TreeType &, const DomainTreeType &, const BoundaryOp &, math::pcg::State &, Interrupter &)
Solve ∇2x = b for x with user-specified boundary conditions, where b is a vector comprising the valu...
Definition: PoissonSolver.h:695
math::pcg::Vector< ValueT > VectorT
Definition: PoissonSolver.h:503
const VIdxTreeT * idxTree
Definition: PoissonSolver.h:506
VIdxTreeT::LeafNodeType VIdxLeafT
Definition: PoissonSolver.h:501
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
LeafIndexOp(const VIndex *count_)
Definition: PoissonSolver.h:296
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
LaplacianMatrix * laplacian
Definition: PoissonSolver.h:505
CopyToVecOp(const SourceTreeType &t, VectorT &v)
Definition: PoissonSolver.h:378
bool isValueOn(const Coord &xyz) const
Return the active state of the voxel at the given coordinates.
Definition: ValueAccessor.h:263
ISLaplacianOp(LaplacianMatrix &m, const VIdxTreeT &idx, const BoolTreeType &mask, const BoundaryOp &op, VectorT &src)
Definition: PoissonSolver.h:511
Dummy NOOP interrupter class defining interface.
Definition: NullInterrupter.h:52
math::pcg::Vector< VectorValueType > VectorT
Definition: PoissonSolver.h:373
const BoundaryOp boundaryOp
Definition: PoissonSolver.h:508
VIndex * count
Definition: PoissonSolver.h:282
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
math::pcg::Vector< VectorValueType >::Ptr createVectorFromTree(const SourceTreeType &source, const typename SourceTreeType::template ValueConverter< VIndex >::Type &index)
Return a vector of the active voxel values of the scalar-valued source tree.
Definition: PoissonSolver.h:405
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: Types.h:419
Functor for use with LeafManager::foreach() to populate an array with per-leaf active voxel counts...
Definition: PoissonSolver.h:280
const BoolTreeType * interiorMask
Definition: PoissonSolver.h:507
Definition: Morphology.h:81
void populateIndexTree(VIndexTreeType &)
Overwrite each active voxel in the given scalar tree with a sequential index, starting from zero...
Definition: PoissonSolver.h:310
void operator()(const VIdxLeafT &idxLeaf, size_t) const
Definition: PoissonSolver.h:515
SourceTreeType::template ValueConverter< VIndex >::Type VIdxTreeT
Definition: PoissonSolver.h:369
LaplacianMatrix::ValueType ValueT
Definition: PoissonSolver.h:502
VIndexTreeType::template ValueConverter< TreeValueType >::Type OutTreeT
Definition: PoissonSolver.h:435
VIdxTreeT::LeafNodeType VIdxLeafT
Definition: PoissonSolver.h:370
math::pcg::SparseStencilMatrix< double, 7 > LaplacianMatrix
The type of a matrix used to represent a three-dimensional Laplacian operator.
Definition: PoissonSolver.h:113
OutTreeT * tree
Definition: PoissonSolver.h:441
void operator()(const Coord &, const Coord &, ValueType &, ValueType &diag) const
Definition: PoissonSolver.h:491
Definition: ValueAccessor.h:219
int32_t Int32
Definition: Types.h:60
const VectorT * vector
Definition: PoissonSolver.h:440
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
LaplacianMatrix::Ptr createISLaplacian(const typename BoolTreeType::template ValueConverter< VIndex >::Type &vectorIndexTree, const BoolTreeType &interiorMask)
Generate a sparse matrix of the index-space (Δx = 1) Laplacian operator using second-order finite dif...
Definition: PoissonSolver.h:612
VIndexTreeType::template ValueConverter< TreeValueType >::Type::Ptr createTreeFromVector(const math::pcg::Vector< VectorValueType > &values, const VIndexTreeType &index, const TreeValueType &background)
Return a tree with the same active voxel topology as the index tree but whose voxel values are taken ...
Definition: PoissonSolver.h:461
Int32 VIndex
Definition: PoissonSolver.h:110
Constant boundary condition functor.
Definition: PoissonSolver.h:490
OutTreeT::LeafNodeType OutLeafT
Definition: PoissonSolver.h:436
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:71
void operator()(const VIdxLeafT &idxLeaf, size_t) const
Definition: PoissonSolver.h:380
void operator()(const LeafType &leaf, size_t leafIdx) const
Definition: PoissonSolver.h:284
bool probeValue(const Coord &xyz, ValueType &value) const
Return the active state of the voxel as well as its value.
Definition: ValueAccessor.h:266
const SourceTreeType * tree
Definition: PoissonSolver.h:375
const VIndex * count
Definition: PoissonSolver.h:295
BoolTreeType::template ValueConverter< VIndex >::Type VIdxTreeT
Definition: PoissonSolver.h:500
TreeType::Ptr solve(const TreeType &, math::pcg::State &, Interrupter &)
Solve ∇2x = b for x, where b is a vector comprising the values of all of the active voxels in the in...
Definition: PoissonSolver.h:665
Functor for use with LeafManager::foreach() to populate a tree with values from a vector...
Definition: PoissonSolver.h:433