Go to the documentation of this file.
10 #ifndef OPENVDB_TOOLS_POTENTIAL_FLOW_HAS_BEEN_INCLUDED
11 #define OPENVDB_TOOLS_POTENTIAL_FLOW_HAS_BEEN_INCLUDED
28 template<
typename VecGr
idT>
31 typename VecGridT::template ValueConverter<typename VecGridT::ValueType::value_type>::Type;
32 using Ptr =
typename Type::Ptr;
42 template<typename GridT, typename MaskT = typename GridT::template ValueConverter<ValueMask>::Type>
43 inline typename MaskT::Ptr
56 template<
typename Vec3T,
typename Gr
idT,
typename MaskT>
57 inline typename GridT::template ValueConverter<Vec3T>::Type::Ptr
59 const typename GridT::template ValueConverter<Vec3T>::Type::ConstPtr boundaryVelocity,
60 const Vec3T& backgroundVelocity);
74 template<
typename Vec3Gr
idT,
typename MaskT,
typename InterrupterT = util::NullInterrupter>
77 InterrupterT* interrupter =
nullptr);
86 template<
typename Vec3Gr
idT>
87 inline typename Vec3GridT::Ptr
89 const Vec3GridT& neumann,
90 const typename Vec3GridT::ValueType backgroundVelocity =
91 zeroVal<typename Vec3GridT::TreeType::ValueType>());
97 namespace potential_flow_internal {
102 template<
typename Gr
idT>
103 inline typename GridT::TreeType::template ValueConverter<ValueMask>::Type::Ptr
104 extractOuterVoxelMask(GridT& inGrid)
106 using MaskTreeT =
typename GridT::TreeType::template ValueConverter<ValueMask>::Type;
108 typename MaskTreeT::Ptr boundaryMask(
new MaskTreeT(inGrid.tree(),
false,
TopologyCopy()));
117 template<
typename Vec3Gr
idT,
typename GradientT>
120 using ValueT =
typename Vec3GridT::ValueType;
123 typename Vec3GridT::ConstAccessor,
BoxSampler>;
127 const Vec3GridT& velocity,
128 const ValueT& backgroundVelocity)
130 , mVelocity(&velocity)
131 , mBackgroundVelocity(backgroundVelocity) { }
134 const ValueT& backgroundVelocity)
136 , mBackgroundVelocity(backgroundVelocity) { }
138 void operator()(
typename Vec3GridT::TreeType::LeafNodeType& leaf,
size_t)
const {
139 auto gradientAccessor = mGradient.getConstAccessor();
141 std::unique_ptr<VelocityAccessor> velocityAccessor;
142 std::unique_ptr<VelocitySamplerT> velocitySampler;
144 velocityAccessor.reset(
new VelocityAccessor(mVelocity->getConstAccessor()));
145 velocitySampler.reset(
new VelocitySamplerT(*velocityAccessor, mVelocity->transform()));
148 for (
auto it = leaf.beginValueOn(); it; ++it) {
149 Coord ijk = it.getCoord();
150 auto gradient = gradientAccessor.getValue(ijk);
152 const Vec3d xyz = mGradient.transform().indexToWorld(ijk);
153 const ValueT sampledVelocity = velocitySampler ?
154 velocitySampler->wsSample(xyz) : zeroVal<ValueT>();
155 auto velocity = sampledVelocity + mBackgroundVelocity;
166 const GradientT& mGradient;
167 const Vec3GridT* mVelocity =
nullptr;
168 const ValueT& mBackgroundVelocity;
173 template<
typename Vec3Gr
idT,
typename MaskT>
177 : mVoxelSize(domainGrid.voxelSize()[0])
179 , mDomainGrid(domainGrid)
183 double& source,
double& diagonal)
const {
185 typename Vec3GridT::ConstAccessor velGridAccessor = mVelGrid.getAccessor();
186 const Coord diff = (ijk - neighbor);
188 if (velGridAccessor.isValueOn(ijk)) {
189 const typename Vec3GridT::ValueType& sampleVel = velGridAccessor.getValue(ijk);
190 source += mVoxelSize*diff[0]*sampleVel[0];
191 source += mVoxelSize*diff[1]*sampleVel[1];
192 source += mVoxelSize*diff[2]*sampleVel[2];
210 template<
typename Gr
idT,
typename MaskT>
211 inline typename MaskT::Ptr
214 using MaskTreeT =
typename MaskT::TreeType;
216 if (!grid.hasUniformVoxels()) {
224 typename MaskTreeT::Ptr maskTree(
new MaskTreeT(interior->tree(),
false,
TopologyCopy()));
225 typename MaskT::Ptr mask = MaskT::create(maskTree);
226 mask->setTransform(grid.transform().copy());
231 mask->tree().topologyDifference(interior->tree());
237 template<
typename Vec3T,
typename Gr
idT,
typename MaskT>
239 const GridT& collider,
241 const typename GridT::template ValueConverter<Vec3T>::Type::ConstPtr boundaryVelocity,
242 const Vec3T& backgroundVelocity)
244 using Vec3GridT =
typename GridT::template ValueConverter<Vec3T>::Type;
245 using TreeT =
typename Vec3GridT::TreeType;
246 using ValueT =
typename TreeT::ValueType;
254 !std::is_floating_point<typename GridT::TreeType::ValueType>::value) {
259 if (backgroundVelocity == zeroVal<Vec3T>() &&
260 (!boundaryVelocity || boundaryVelocity->empty())) {
261 auto neumann = Vec3GridT::create();
262 neumann->setTransform(collider.transform().copy());
267 using MaskTreeT =
typename GridT::TreeType::template ValueConverter<ValueMask>::Type;
268 typename MaskTreeT::Ptr boundary(
new MaskTreeT(domain.tree(),
false,
TopologyCopy()));
269 boundary->topologyIntersection(collider.tree());
271 typename TreeT::Ptr neumannTree(
new TreeT(*boundary, zeroVal<ValueT>(),
TopologyCopy()));
272 neumannTree->voxelizeActiveTiles();
279 if (boundaryVelocity && !boundaryVelocity->empty()) {
280 ComputeNeumannVelocityOp<Vec3GridT, GradientT>
281 neumannOp(*
gradient, *boundaryVelocity, backgroundVelocity);
282 leafManager.
foreach(neumannOp,
false);
285 ComputeNeumannVelocityOp<Vec3GridT, GradientT>
286 neumannOp(*
gradient, backgroundVelocity);
287 leafManager.
foreach(neumannOp,
false);
293 typename Vec3GridT::Ptr neumann(Vec3GridT::create(neumannTree));
294 neumann->setTransform(collider.transform().copy());
300 template<
typename Vec3Gr
idT,
typename MaskT,
typename InterrupterT>
301 inline typename VectorToScalarGrid<Vec3GridT>::Ptr
305 using ScalarT =
typename Vec3GridT::ValueType::value_type;
306 using ScalarTreeT =
typename Vec3GridT::TreeType::template ValueConverter<ScalarT>::Type;
307 using ScalarGridT =
typename Vec3GridT::template ValueConverter<ScalarT>::Type;
312 ScalarTreeT solveTree(domain.tree(), zeroVal<ScalarT>(),
TopologyCopy());
313 solveTree.voxelizeActiveTiles();
316 if (!interrupter) interrupter = &nullInterrupt;
319 SolveBoundaryOp<Vec3GridT, MaskT>
solve(neumann, domain);
320 typename ScalarTreeT::Ptr potentialTree =
323 auto potential = ScalarGridT::create(potentialTree);
324 potential->setTransform(domain.transform().copy());
330 template<
typename Vec3Gr
idT>
331 inline typename Vec3GridT::Ptr
333 const Vec3GridT& neumann,
334 const typename Vec3GridT::ValueType backgroundVelocity)
336 using Vec3T =
const typename Vec3GridT::ValueType;
337 using potential_flow_internal::extractOuterVoxelMask;
352 auto applyNeumann = [&
gradient, &neumann] (
353 const MaskGrid::TreeType::LeafNodeType& leaf, size_t)
355 typename Vec3GridT::Accessor gradientAccessor =
gradient->getAccessor();
356 typename Vec3GridT::ConstAccessor neumannAccessor = neumann.getAccessor();
357 for (
auto it = leaf.beginValueOn(); it; ++it) {
358 const Coord ijk = it.getCoord();
359 typename Vec3GridT::ValueType value;
360 if (neumannAccessor.probeValue(ijk, value)) {
361 gradientAccessor.setValue(ijk, value);
368 leafManager.
foreach(applyNeumann);
372 if (backgroundVelocity != zeroVal<Vec3T>()) {
373 auto applyBackgroundVelocity = [&backgroundVelocity] (
374 typename Vec3GridT::TreeType::LeafNodeType& leaf, size_t)
376 for (
auto it = leaf.beginValueOn(); it; ++it) {
377 it.setValue(it.getValue() - backgroundVelocity);
382 leafManager2.
foreach(applyBackgroundVelocity);
395 #endif // OPENVDB_TOOLS_POTENTIAL_FLOW_HAS_BEEN_INCLUDED
Information about the state of a conjugate gradient solution.
Definition: ConjGradient.h:46
void foreach(const LeafOp &op, bool threaded=true, size_t grainSize=1)
Threaded method that applies a user-supplied functor to each leaf node in the LeafManager.
Definition: LeafManager.h:496
Signed (x, y, z) 32-bit integer coordinates.
Definition: Coord.h:26
Solve Poisson's equation ∇2x = b for x, where b is a vector comprising the values of all of the activ...
SharedPtr< Grid > Ptr
Definition: Grid.h:574
Vec3< double > Vec3d
Definition: Vec3.h:662
Dummy NOOP interrupter class defining interface.
Definition: NullInterrupter.h:26
State solve(const PositiveDefMatrix &A, const Vector< typename PositiveDefMatrix::ValueType > &b, Vector< typename PositiveDefMatrix::ValueType > &x, Preconditioner< typename PositiveDefMatrix::ValueType > &preconditioner, const State &termination=terminationDefaults< typename PositiveDefMatrix::ValueType >())
Solve Ax = b via the preconditioned conjugate gradient method.
Definition: ConjGradient.h:1612
Apply an operator to an input grid to produce an output grid with the same active voxel topology but ...
Definition: Exceptions.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:83
Tag dispatch class that distinguishes topology copy constructors from deep copy constructors.
Definition: Types.h:681
Construct boolean mask grids from grids of arbitrary type.
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:154
Implementation of morphological dilation and erosion.
#define OPENVDB_VERSION_NAME
The version namespace name for this library version.
Definition: version.h:102
Definition: Exceptions.h:65
@ GRID_LEVEL_SET
Definition: Types.h:454
Definition: Exceptions.h:13
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:82