archive/blender-archive
Archived
1
1
Fork 0
Code Pull Requests Packages Activity

Fluid: Updated Mantaflow source files

Browse Source 
This update introduces two improvements from the Mantaflow repository:

(1) Improved particle sampling:
- Liquid and secondary particles are sampled more predictably. With all parameters being equal, baked particles will be computed at the exact same position during every bake.
- Before, this was not guaranteed.

(2) Sparse grid caching:
- While saving grid data to disk, grids will from now on be saved in a sparse structure whenever possible (e.g. density, flame but not levelsets).
- With the sparse optimization grid cells with a value under the 'Empty Space' value (already present in domain settings) will not be cached.
- The main benefits of this optimization are: Smaller cache sizes and faster playback of simulation data in the viewport.
- This optimization works 'out-of-the-box'. There is no option in the UI to enable it.
- For now, only smoke simulation grids will take advantage of this optimization.
This commit is contained in:
Sebastián Barschkis
2020-11-25 23:17:47 +01:00
parent f7223d5f72
commit e09d0c0d07
13 changed files with 274 additions and 130 deletions
Download Patch File Download Diff File Expand all files Collapse all files

118
extern/mantaflow/preprocessed/fileio/iovdb.cpp vendored
View File

@@ -31,6 +31,8 @@
# include "openvdb/openvdb.h"
# include <openvdb/points/PointConversion.h>
# include <openvdb/points/PointCount.h>
# include <openvdb/tools/Clip.h>
# include <openvdb/tools/Dense.h>
#endif
#define POSITION_NAME "P"
@@ -45,15 +47,29 @@ namespace Manta {
template<class GridType, class T> void importVDB(typename GridType::Ptr from, Grid<T> *to)
{
using ValueT = typename GridType::ValueType;
typename GridType::Accessor accessor = from->getAccessor();
FOR_IJK(*to)
{
openvdb::Coord xyz(i, j, k);
ValueT vdbValue = accessor.getValue(xyz);
T toMantaValue;
convertFrom(vdbValue, &toMantaValue);
to->set(i, j, k, toMantaValue);
// Check if current grid is to be read as a sparse grid, active voxels (only) will be copied
if (to->saveSparse()) {
to->clear(); // Ensure that destination grid is empty before writing
for (typename GridType::ValueOnCIter iter = from->cbeginValueOn(); iter.test(); ++iter) {
ValueT vdbValue = *iter;
openvdb::Coord coord = iter.getCoord();
T toMantaValue;
convertFrom(vdbValue, &toMantaValue);
to->set(coord.x(), coord.y(), coord.z(), toMantaValue);
}
}
// When importing all grid cells, using a grid accessor is usually faster than a value iterator
else {
typename GridType::Accessor accessor = from->getAccessor();
FOR_IJK(*to)
{
openvdb::Coord xyz(i, j, k);
ValueT vdbValue = accessor.getValue(xyz);
T toMantaValue;
convertFrom(vdbValue, &toMantaValue);
to->set(i, j, k, toMantaValue);
}
}
}
@@ -173,19 +189,43 @@ static void setGridOptions(typename GridType::Ptr grid,
grid->setSaveFloatAsHalf(precision == PRECISION_MINI || precision == PRECISION_HALF);
}
template<class T, class GridType> typename GridType::Ptr exportVDB(Grid<T> *from)
template<class T, class GridType>
typename GridType::Ptr exportVDB(Grid<T> *from, float clip, openvdb::FloatGrid::Ptr clipGrid)
{
using ValueT = typename GridType::ValueType;
typename GridType::Ptr to = GridType::create();
typename GridType::Accessor accessor = to->getAccessor();
typename GridType::Ptr to = GridType::create(ValueT(0));
FOR_IJK(*from)
{
openvdb::Coord xyz(i, j, k);
T fromMantaValue = (*from)(i, j, k);
ValueT vdbValue;
convertTo(&vdbValue, fromMantaValue);
accessor.setValue(xyz, vdbValue);
// Copy data from grid by creating a vdb dense structure and then copying that into a vdb grid
// This is the fastest way to copy data for both dense and sparse grids -> if (true)
if (true) {
ValueT *data = (ValueT *)from->getData();
openvdb::math::CoordBBox bbox(
openvdb::Coord(0),
openvdb::Coord(from->getSizeX() - 1, from->getSizeY() - 1, from->getSizeZ() - 1));
openvdb::tools::Dense<ValueT, openvdb::tools::MemoryLayout::LayoutXYZ> dense(bbox, data);
// Trick: Set clip value to very small / negative value in order to copy all values of dense
// grids
float tmpClip = (from->saveSparse()) ? clip : -std::numeric_limits<Real>::max();
// Copy from dense to sparse grid structure considering clip value
openvdb::tools::copyFromDense(dense, *to, ValueT(tmpClip));
// If present, use clip grid to trim down current vdb grid even more
if (from->saveSparse() && clipGrid && !clipGrid->empty()) {
to = openvdb::tools::clip(*to, *clipGrid);
}
}
// Alternatively, reading all grid cells with an accessor (slightly slower) is possible like this
else {
typename GridType::Accessor accessor = to->getAccessor();
FOR_IJK(*from)
{
openvdb::Coord xyz(i, j, k);
T fromMantaValue = (*from)(i, j, k);
ValueT vdbValue;
convertTo(&vdbValue, fromMantaValue);
accessor.setValue(xyz, vdbValue);
}
}
return to;
}
@@ -346,7 +386,9 @@ int writeObjectsVDB(const string &filename,
float worldSize,
bool skipDeletedParts,
int compression,
int precision)
int precision,
float clip,
const Grid<Real> *clipGrid)
{
openvdb::initialize();
openvdb::io::File file(filename);
@@ -357,6 +399,19 @@ int writeObjectsVDB(const string &filename,
std::vector<ParticleDataBase *> pdbBuffer;
// Convert given clip grid to vdb clip grid
openvdb::FloatGrid::Ptr vdbClipGrid = nullptr;
if (clipGrid) {
vdbClipGrid = openvdb::FloatGrid::create();
Real *data = (Real *)clipGrid->getData();
openvdb::math::CoordBBox bbox(openvdb::Coord(0),
openvdb::Coord(clipGrid->getSizeX() - 1,
clipGrid->getSizeY() - 1,
clipGrid->getSizeZ() - 1));
openvdb::tools::Dense<float, openvdb::tools::MemoryLayout::LayoutXYZ> dense(bbox, data);
openvdb::tools::copyFromDense(dense, *vdbClipGrid, clip);
}
for (std::vector<PbClass *>::iterator iter = objects->begin(); iter != objects->end(); ++iter) {
openvdb::GridClass gClass = openvdb::GRID_UNKNOWN;
openvdb::GridBase::Ptr vdbGrid;
@@ -368,10 +423,14 @@ int writeObjectsVDB(const string &filename,
if (GridBase *mantaGrid = dynamic_cast<GridBase *>(*iter)) {
if (clipGrid) {
assertMsg(clipGrid->getSize() == mantaGrid->getSize(),
"writeObjectsVDB: Clip grid and exported grid must have the same size");
}
if (mantaGrid->getType() & GridBase::TypeInt) {
debMsg("Writing int grid '" << mantaGrid->getName() << "' to vdb file " << filename, 1);
Grid<int> *mantaIntGrid = (Grid<int> *)mantaGrid;
vdbGrid = exportVDB<int, openvdb::Int32Grid>(mantaIntGrid);
vdbGrid = exportVDB<int, openvdb::Int32Grid>(mantaIntGrid, clip, vdbClipGrid);
gridsVDB.push_back(vdbGrid);
}
else if (mantaGrid->getType() & GridBase::TypeReal) {
@@ -379,7 +438,9 @@ int writeObjectsVDB(const string &filename,
gClass = (mantaGrid->getType() & GridBase::TypeLevelset) ? openvdb::GRID_LEVEL_SET :
openvdb::GRID_FOG_VOLUME;
Grid<Real> *mantaRealGrid = (Grid<Real> *)mantaGrid;
vdbGrid = exportVDB<Real, openvdb::FloatGrid>(mantaRealGrid);
// Only supply clip grid if real grid is not equal to the clip grid
openvdb::FloatGrid::Ptr tmpClipGrid = (mantaRealGrid == clipGrid) ? nullptr : vdbClipGrid;
vdbGrid = exportVDB<Real, openvdb::FloatGrid>(mantaRealGrid, clip, tmpClipGrid);
gridsVDB.push_back(vdbGrid);
}
else if (mantaGrid->getType() & GridBase::TypeVec3) {
@@ -387,7 +448,7 @@ int writeObjectsVDB(const string &filename,
gClass = (mantaGrid->getType() & GridBase::TypeMAC) ? openvdb::GRID_STAGGERED :
openvdb::GRID_UNKNOWN;
Grid<Vec3> *mantaVec3Grid = (Grid<Vec3> *)mantaGrid;
vdbGrid = exportVDB<Vec3, openvdb::Vec3SGrid>(mantaVec3Grid);
vdbGrid = exportVDB<Vec3, openvdb::Vec3SGrid>(mantaVec3Grid, clip, vdbClipGrid);
gridsVDB.push_back(vdbGrid);
}
else {
@@ -620,9 +681,12 @@ template void importVDB<openvdb::Int32Grid, int>(openvdb::Int32Grid::Ptr from, G
template void importVDB<openvdb::FloatGrid, Real>(openvdb::FloatGrid::Ptr from, Grid<Real> *to);
template void importVDB<openvdb::Vec3SGrid, Vec3>(openvdb::Vec3SGrid::Ptr from, Grid<Vec3> *to);
template openvdb::Int32Grid::Ptr exportVDB<int, openvdb::Int32Grid>(Grid<int> *from);
template openvdb::FloatGrid::Ptr exportVDB<Real, openvdb::FloatGrid>(Grid<Real> *from);
template openvdb::Vec3SGrid::Ptr exportVDB<Vec3, openvdb::Vec3SGrid>(Grid<Vec3> *from);
template openvdb::Int32Grid::Ptr exportVDB<int, openvdb::Int32Grid>(
Grid<int> *from, float clip = 1e-4, openvdb::FloatGrid::Ptr clipGrid = nullptr);
template openvdb::FloatGrid::Ptr exportVDB<Real, openvdb::FloatGrid>(
Grid<Real> *from, float clip = 1e-4, openvdb::FloatGrid::Ptr clipGrid = nullptr);
template openvdb::Vec3SGrid::Ptr exportVDB<Vec3, openvdb::Vec3SGrid>(
Grid<Vec3> *from, float clip = 1e-4, openvdb::FloatGrid::Ptr clipGrid = nullptr);
openvdb::points::PointDataGrid::Ptr exportVDB(BasicParticleSystem *from,
std::vector<ParticleDataBase *> &fromPData,
@@ -652,7 +716,9 @@ int writeObjectsVDB(const string &filename,
float worldSize,
bool skipDeletedParts,
int compression,
int precision)
int precision,
float clip,
const Grid<Real> *clipGrid)
{
errMsg("Cannot save to .vdb file. Mantaflow has not been built with OpenVDB support.");
return 0;

20
extern/mantaflow/preprocessed/fileio/mantaio.cpp vendored
View File

@@ -83,7 +83,9 @@ int save(const string &name,
bool skipDeletedParts = false,
int compression = COMPRESSION_ZIP,
bool precisionHalf = true,
int precision = PRECISION_HALF)
int precision = PRECISION_HALF,
float clip = 1e-4,
const Grid<Real> *clipGrid = nullptr)
{
if (!precisionHalf) {
@@ -102,7 +104,8 @@ int save(const string &name,
else if (ext == ".vol")
return writeGridsVol(name, &objects);
if (ext == ".vdb")
return writeObjectsVDB(name, &objects, worldSize, skipDeletedParts, compression, precision);
return writeObjectsVDB(
name, &objects, worldSize, skipDeletedParts, compression, precision, clip, clipGrid);
else if (ext == ".npz")
return writeGridsNumpy(name, &objects);
else if (ext == ".txt")
@@ -129,8 +132,17 @@ static PyObject *_W_1(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
int compression = _args.getOpt<int>("compression", 4, COMPRESSION_ZIP, &_lock);
bool precisionHalf = _args.getOpt<bool>("precisionHalf", 5, true, &_lock);
int precision = _args.getOpt<int>("precision", 6, PRECISION_HALF, &_lock);
_retval = toPy(
save(name, objects, worldSize, skipDeletedParts, compression, precisionHalf, precision));
float clip = _args.getOpt<float>("clip", 7, 1e-4, &_lock);
const Grid<Real> *clipGrid = _args.getPtrOpt<Grid<Real>>("clipGrid", 8, nullptr, &_lock);
_retval = toPy(save(name,
objects,
worldSize,
skipDeletedParts,
compression,
precisionHalf,
precision,
clip,
clipGrid));
_args.check();
}
pbFinalizePlugin(parent, "save", !noTiming);

4
extern/mantaflow/preprocessed/fileio/mantaio.h vendored
View File

@@ -75,7 +75,9 @@ int writeObjectsVDB(const std::string &filename,
float scale = 1.0,
bool skipDeletedParts = false,
int compression = COMPRESSION_ZIP,
int precision = PRECISION_HALF);
int precision = PRECISION_HALF,
float clip = 1e-4,
const Grid<Real> *clipGrid = nullptr);
int readObjectsVDB(const std::string &filename,
std::vector<PbClass *> *objects,
float scale = 1.0);

2
extern/mantaflow/preprocessed/gitinfo.h vendored
View File

@@ -1,3 +1,3 @@
#define MANTA_GIT_VERSION "commit 73990a8a5b876e2b136a646258f714c08c5828da"
#define MANTA_GIT_VERSION "commit bb7cde47b6e04fa62815c70775dc70f02065599f"

29
extern/mantaflow/preprocessed/grid.cpp vendored
View File

@@ -60,7 +60,7 @@ template<> inline GridBase::GridType typeList<Vec3>()
}
template<class T>
Grid<T>::Grid(FluidSolver *parent, bool show) : GridBase(parent), externalData(false)
Grid<T>::Grid(FluidSolver *parent, bool show, bool sparse) : GridBase(parent), mExternalData(false)
{
mType = typeList<T>();
mSize = parent->getGridSize();
@@ -70,22 +70,23 @@ Grid<T>::Grid(FluidSolver *parent, bool show) : GridBase(parent), externalData(f
mDx = 1.0 / mSize.max();
clear();
setHidden(!show);
#if OPENVDB == 1
mSaveSparse = sparse;
#else
if (sparse)
debMsg("Cannot enable sparse save option without OpenVDB", 1);
mSaveSparse = false;
#endif
}
template<class T>
Grid<T>::Grid(FluidSolver *parent, T *data, bool show)
: GridBase(parent), mData(data), externalData(true)
Grid<T>::Grid(FluidSolver *parent, T *data, bool show) : Grid<T>::Grid(parent, show)
{
mType = typeList<T>();
mSize = parent->getGridSize();
mStrideZ = parent->is2D() ? 0 : (mSize.x * mSize.y);
mDx = 1.0 / mSize.max();
setHidden(!show);
mData = data;
}
template<class T> Grid<T>::Grid(const Grid<T> &a) : GridBase(a.getParent()), externalData(false)
template<class T> Grid<T>::Grid(const Grid<T> &a) : GridBase(a.getParent()), mExternalData(false)
{
mSize = a.mSize;
mType = a.mType;
@@ -98,7 +99,7 @@ template<class T> Grid<T>::Grid(const Grid<T> &a) : GridBase(a.getParent()), ext
template<class T> Grid<T>::~Grid()
{
if (!externalData) {
if (!mExternalData) {
mParent->freeGridPointer<T>(mData);
}
}
@@ -114,8 +115,8 @@ template<class T> void Grid<T>::swap(Grid<T> &other)
other.getSizeZ() != getSizeZ())
errMsg("Grid::swap(): Grid dimensions mismatch.");
if (externalData || other.externalData)
errMsg("Grid::swap(): Cannot swap if one grid stores externalData.");
if (mExternalData || other.mExternalData)
errMsg("Grid::swap(): Cannot swap if one grid stores mExternalData.");
T *dswap = other.mData;
other.mData = mData;

30
extern/mantaflow/preprocessed/grid.h vendored
View File

@@ -402,7 +402,7 @@ class GridBase : public PbClass {
template<class T> class Grid : public GridBase {
public:
//! init new grid, values are set to zero
Grid(FluidSolver *parent, bool show = true);
Grid(FluidSolver *parent, bool show = true, bool sparse = false);
static int _W_10(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
{
PbClass *obj = Pb::objFromPy(_self);
@@ -416,7 +416,8 @@ template<class T> class Grid : public GridBase {
ArgLocker _lock;
FluidSolver *parent = _args.getPtr<FluidSolver>("parent", 0, &_lock);
bool show = _args.getOpt<bool>("show", 1, true, &_lock);
obj = new Grid(parent, show);
bool sparse = _args.getOpt<bool>("sparse", 2, false, &_lock);
obj = new Grid(parent, show, sparse);
obj->registerObject(_self, &_args);
_args.check();
}
@@ -581,6 +582,16 @@ template<class T> class Grid : public GridBase {
DEBUG_ONLY(checkIndex(idx));
return mData[idx];
}
//! raw data access
inline T *getData() const
{
return mData;
}
//! query if this grid should be saved as a sparse grid
inline bool saveSparse()
{
return mSaveSparse;
}
//! set data
inline void set(int i, int j, int k, T &val)
@@ -1290,7 +1301,8 @@ template<class T> class Grid : public GridBase {
protected:
T *mData;
bool externalData; // True if mData is managed outside of the Fluidsolver
bool mExternalData; // True if mData is managed outside of the Fluidsolver
bool mSaveSparse; // True if this grid may be cached in a sparse structure
public:
PbArgs _args;
}
@@ -1302,7 +1314,8 @@ template<class T> class Grid : public GridBase {
//! Special function for staggered grids
class MACGrid : public Grid<Vec3> {
public:
MACGrid(FluidSolver *parent, bool show = true) : Grid<Vec3>(parent, show)
MACGrid(FluidSolver *parent, bool show = true, bool sparse = false)
: Grid<Vec3>(parent, show, sparse)
{
mType = (GridType)(TypeMAC | TypeVec3);
}
@@ -1319,7 +1332,8 @@ class MACGrid : public Grid<Vec3> {
ArgLocker _lock;
FluidSolver *parent = _args.getPtr<FluidSolver>("parent", 0, &_lock);
bool show = _args.getOpt<bool>("show", 1, true, &_lock);
obj = new MACGrid(parent, show);
bool sparse = _args.getOpt<bool>("sparse", 2, false, &_lock);
obj = new MACGrid(parent, show, sparse);
obj->registerObject(_self, &_args);
_args.check();
}
@@ -1425,7 +1439,8 @@ class MACGrid : public Grid<Vec3> {
//! Special functions for FlagGrid
class FlagGrid : public Grid<int> {
public:
FlagGrid(FluidSolver *parent, int dim = 3, bool show = true) : Grid<int>(parent, show)
FlagGrid(FluidSolver *parent, int dim = 3, bool show = true, bool sparse = false)
: Grid<int>(parent, show, sparse)
{
mType = (GridType)(TypeFlags | TypeInt);
}
@@ -1443,7 +1458,8 @@ class FlagGrid : public Grid<int> {
FluidSolver *parent = _args.getPtr<FluidSolver>("parent", 0, &_lock);
int dim = _args.getOpt<int>("dim", 1, 3, &_lock);
bool show = _args.getOpt<bool>("show", 2, true, &_lock);
obj = new FlagGrid(parent, dim, show);
bool sparse = _args.getOpt<bool>("sparse", 3, false, &_lock);
obj = new FlagGrid(parent, dim, show, sparse);
obj->registerObject(_self, &_args);
_args.check();
}

10
extern/mantaflow/preprocessed/particle.cpp vendored
View File

@@ -28,9 +28,15 @@
using namespace std;
namespace Manta {
ParticleBase::ParticleBase(FluidSolver *parent)
: PbClass(parent), mMaxParticles(0), mAllowCompress(true), mFreePdata(false)
int ParticleBase::globalSeed = 9832;
ParticleBase::ParticleBase(FluidSolver *parent, int fixedSeed)
: PbClass(parent), mMaxParticles(0), mAllowCompress(true), mFreePdata(false), mSeed(fixedSeed)
{
// use global random seed if none is given
if (fixedSeed == -1) {
mSeed = globalSeed;
}
}
ParticleBase::~ParticleBase()

40
extern/mantaflow/preprocessed/particle.h vendored
View File

@@ -47,7 +47,7 @@ class ParticleBase : public PbClass {
PINVALID = (1 << 30), // unused
};
ParticleBase(FluidSolver *parent);
ParticleBase(FluidSolver *parent, int fixedSeed = -1);
static int _W_0(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
{
PbClass *obj = Pb::objFromPy(_self);
@@ -60,7 +60,8 @@ class ParticleBase : public PbClass {
{
ArgLocker _lock;
FluidSolver *parent = _args.getPtr<FluidSolver>("parent", 0, &_lock);
obj = new ParticleBase(parent);
int fixedSeed = _args.getOpt<int>("fixedSeed", 1, -1, &_lock);
obj = new ParticleBase(parent, fixedSeed);
obj->registerObject(_self, &_args);
_args.check();
}
@@ -111,6 +112,11 @@ class ParticleBase : public PbClass {
return;
}
inline int getSeed()
{
return mSeed;
}
//! particle data functions
//! create a particle data object
@@ -192,9 +198,13 @@ class ParticleBase : public PbClass {
//! per particle)
std::vector<ParticleDataImpl<Real> *> mPdataReal;
std::vector<ParticleDataImpl<Vec3> *> mPdataVec3;
std::vector<ParticleDataImpl<int> *>
mPdataInt; //! indicate that pdata of this particle system is copied, and needs to be freed
std::vector<ParticleDataImpl<int> *> mPdataInt;
//! indicate that pdata of this particle system is copied, and needs to be freed
bool mFreePdata;
//! custom seed for particle systems, used by plugins
int mSeed; //! fix global random seed storage, used mainly by functions in this class
static int globalSeed;
public:
PbArgs _args;
}
@@ -2285,15 +2295,14 @@ void ParticleSystem<S>::advectInGrid(const FlagGrid &flags,
}
template<class S> struct KnProjectParticles : public KernelBase {
KnProjectParticles(ParticleSystem<S> &part, Grid<Vec3> &gradient)
: KernelBase(part.size()), part(part), gradient(gradient)
KnProjectParticles(ParticleSystem<S> &part, Grid<Vec3> &gradient, RandomStream &rand)
: KernelBase(part.size()), part(part), gradient(gradient), rand(rand)
{
runMessage();
run();
}
inline void op(IndexInt idx, ParticleSystem<S> &part, Grid<Vec3> &gradient)
inline void op(IndexInt idx, ParticleSystem<S> &part, Grid<Vec3> &gradient, RandomStream &rand)
{
static RandomStream rand(3123984);
const double jlen = 0.1;
if (part.isActive(idx)) {
@@ -2321,6 +2330,11 @@ template<class S> struct KnProjectParticles : public KernelBase {
return gradient;
}
typedef Grid<Vec3> type1;
inline RandomStream &getArg2()
{
return rand;
}
typedef RandomStream type2;
void runMessage()
{
debMsg("Executing kernel KnProjectParticles ", 3);
@@ -2332,15 +2346,17 @@ template<class S> struct KnProjectParticles : public KernelBase {
{
const IndexInt _sz = size;
for (IndexInt i = 0; i < _sz; i++)
op(i, part, gradient);
op(i, part, gradient, rand);
}
ParticleSystem<S> &part;
Grid<Vec3> &gradient;
RandomStream &rand;
};
template<class S> void ParticleSystem<S>::projectOutside(Grid<Vec3> &gradient)
{
KnProjectParticles<S>(*this, gradient);
RandomStream rand(globalSeed);
KnProjectParticles<S>(*this, gradient, rand);
}
template<class S> struct KnProjectOutOfBnd : public KernelBase {
@@ -2531,14 +2547,14 @@ template<class S> void ParticleSystem<S>::insertBufferedParticles()
int insertFlag;
Vec3 insertPos;
static RandomStream mRand(9832);
RandomStream rand(globalSeed);
for (IndexInt i = 0; i < numNewParts; ++i) {
// get random index in newBuffer vector
// we are inserting particles randomly so that they are sampled uniformly in the fluid region
// otherwise, regions of fluid can remain completely empty once mData.size() == maxParticles is
// reached.
int randIndex = floor(mRand.getReal() * mNewBufferPos.size());
int randIndex = floor(rand.getReal() * mNewBufferPos.size());
// get elements from new buffers with random index
std::swap(mNewBufferPos[randIndex], mNewBufferPos.back());

16
extern/mantaflow/preprocessed/plugin/flip.cpp vendored
View File

@@ -41,7 +41,7 @@ void sampleFlagsWithParticles(const FlagGrid &flags,
const bool is3D = flags.is3D();
const Real jlen = randomness / discretization;
const Vec3 disp(1.0 / discretization, 1.0 / discretization, 1.0 / discretization);
RandomStream mRand(9832);
RandomStream rand(parts.getSeed());
FOR_IJK_BND(flags, 0)
{
@@ -53,7 +53,7 @@ void sampleFlagsWithParticles(const FlagGrid &flags,
for (int dj = 0; dj < discretization; dj++)
for (int di = 0; di < discretization; di++) {
Vec3 subpos = pos + disp * Vec3(0.5 + di, 0.5 + dj, 0.5 + dk);
subpos += jlen * (Vec3(1, 1, 1) - 2.0 * mRand.getVec3());
subpos += jlen * (Vec3(1, 1, 1) - 2.0 * rand.getVec3());
if (!is3D)
subpos[2] = 0.5;
parts.addBuffered(subpos);
@@ -113,7 +113,7 @@ void sampleLevelsetWithParticles(const LevelsetGrid &phi,
const bool is3D = phi.is3D();
const Real jlen = randomness / discretization;
const Vec3 disp(1.0 / discretization, 1.0 / discretization, 1.0 / discretization);
RandomStream mRand(9832);
RandomStream rand(parts.getSeed());
if (reset) {
parts.clear();
@@ -132,7 +132,7 @@ void sampleLevelsetWithParticles(const LevelsetGrid &phi,
for (int dj = 0; dj < discretization; dj++)
for (int di = 0; di < discretization; di++) {
Vec3 subpos = pos + disp * Vec3(0.5 + di, 0.5 + dj, 0.5 + dk);
subpos += jlen * (Vec3(1, 1, 1) - 2.0 * mRand.getVec3());
subpos += jlen * (Vec3(1, 1, 1) - 2.0 * rand.getVec3());
if (!is3D)
subpos[2] = 0.5;
if (phi.getInterpolated(subpos) > 0.)
@@ -205,7 +205,7 @@ void sampleShapeWithParticles(const Shape &shape,
const bool is3D = flags.is3D();
const Real jlen = randomness / discretization;
const Vec3 disp(1.0 / discretization, 1.0 / discretization, 1.0 / discretization);
RandomStream mRand(9832);
RandomStream rand(parts.getSeed());
if (reset) {
parts.clear();
@@ -223,7 +223,7 @@ void sampleShapeWithParticles(const Shape &shape,
for (int dj = 0; dj < discretization; dj++)
for (int di = 0; di < discretization; di++) {
Vec3 subpos = pos + disp * Vec3(0.5 + di, 0.5 + dj, 0.5 + dk);
subpos += jlen * (Vec3(1, 1, 1) - 2.0 * mRand.getVec3());
subpos += jlen * (Vec3(1, 1, 1) - 2.0 * rand.getVec3());
if (!is3D)
subpos[2] = 0.5;
if (exclude && exclude->getInterpolated(subpos) <= 0.)
@@ -576,7 +576,7 @@ void adjustNumber(BasicParticleSystem &parts,
}
// seed new particles
RandomStream mRand(9832);
RandomStream rand(parts.getSeed());
FOR_IJK(tmp)
{
int cnt = tmp(i, j, k);
@@ -593,7 +593,7 @@ void adjustNumber(BasicParticleSystem &parts,
if (flags.isFluid(i, j, k) && cnt < minParticles) {
for (int m = cnt; m < minParticles; m++) {
Vec3 pos = Vec3(i, j, k) + mRand.getVec3();
Vec3 pos = Vec3(i, j, k) + rand.getVec3();
// Vec3 pos (i + 0.5, j + 0.5, k + 0.5); // cell center
parts.addBuffered(pos);
}

65
extern/mantaflow/preprocessed/plugin/secondaryparticles.cpp vendored
View File

@@ -479,7 +479,8 @@ struct knFlipSampleSecondaryParticlesMoreCylinders : public KernelBase {
const Real k_ta,
const Real k_wc,
const Real dt,
const int itype = FlagGrid::TypeFluid)
const int itype,
RandomStream &rand)
: KernelBase(&flags, 0),
flags(flags),
v(v),
@@ -497,7 +498,8 @@ struct knFlipSampleSecondaryParticlesMoreCylinders : public KernelBase {
k_ta(k_ta),
k_wc(k_wc),
dt(dt),
itype(itype)
itype(itype),
rand(rand)
{
runMessage();
run();
@@ -521,13 +523,13 @@ struct knFlipSampleSecondaryParticlesMoreCylinders : public KernelBase {
const Real k_ta,
const Real k_wc,
const Real dt,
const int itype = FlagGrid::TypeFluid)
const int itype,
RandomStream &rand)
{
if (!(flags(i, j, k) & itype))
return;
static RandomStream mRand(9832);
Real radius =
0.25; // diameter=0.5 => sampling with two cylinders in each dimension since cell size=1
for (Real x = i - radius; x <= i + radius; x += 2 * radius) {
@@ -549,9 +551,9 @@ struct knFlipSampleSecondaryParticlesMoreCylinders : public KernelBase {
cross(e1, dir)); // perpendicular to dir and e1, so e1 and e1 create reference plane
for (int di = 0; di < n; di++) {
const Real r = radius * sqrt(mRand.getReal()); // distance to cylinder axis
const Real theta = mRand.getReal() * Real(2) * M_PI; // azimuth
const Real h = mRand.getReal() * norm(dt * vi); // distance to reference plane
const Real r = radius * sqrt(rand.getReal()); // distance to cylinder axis
const Real theta = rand.getReal() * Real(2) * M_PI; // azimuth
const Real h = rand.getReal() * norm(dt * vi); // distance to reference plane
Vec3 xd = xi + r * cos(theta) * e1 + r * sin(theta) * e2 + h * getNormalized(vi);
if (!flags.is3D())
xd.z = 0;
@@ -561,7 +563,7 @@ struct knFlipSampleSecondaryParticlesMoreCylinders : public KernelBase {
vi; // init velocity of new particle
Real temp = (KE + TA + WC) / 3;
l_sec[l_sec.size() - 1] = ((lMax - lMin) * temp) + lMin +
mRand.getReal() * 0.1; // init lifetime of new particle
rand.getReal() * 0.1; // init lifetime of new particle
// init type of new particle
if (neighborRatio(i, j, k) < c_s) {
@@ -663,6 +665,11 @@ struct knFlipSampleSecondaryParticlesMoreCylinders : public KernelBase {
return itype;
}
typedef int type16;
inline RandomStream &getArg17()
{
return rand;
}
typedef RandomStream type17;
void runMessage()
{
debMsg("Executing kernel knFlipSampleSecondaryParticlesMoreCylinders ", 3);
@@ -696,7 +703,8 @@ struct knFlipSampleSecondaryParticlesMoreCylinders : public KernelBase {
k_ta,
k_wc,
dt,
itype);
itype,
rand);
}
const FlagGrid &flags;
const MACGrid &v;
@@ -715,6 +723,7 @@ struct knFlipSampleSecondaryParticlesMoreCylinders : public KernelBase {
const Real k_wc;
const Real dt;
const int itype;
RandomStream &rand;
};
// adds secondary particles to &pts_sec for every fluid cell in &flags according to the potential
@@ -738,7 +747,8 @@ struct knFlipSampleSecondaryParticles : public KernelBase {
const Real k_ta,
const Real k_wc,
const Real dt,
const int itype = FlagGrid::TypeFluid)
const int itype,
RandomStream &rand)
: KernelBase(&flags, 0),
flags(flags),
v(v),
@@ -756,7 +766,8 @@ struct knFlipSampleSecondaryParticles : public KernelBase {
k_ta(k_ta),
k_wc(k_wc),
dt(dt),
itype(itype)
itype(itype),
rand(rand)
{
runMessage();
run();
@@ -780,7 +791,8 @@ struct knFlipSampleSecondaryParticles : public KernelBase {
const Real k_ta,
const Real k_wc,
const Real dt,
const int itype = FlagGrid::TypeFluid)
const int itype,
RandomStream &rand)
{
if (!(flags(i, j, k) & itype))
@@ -793,9 +805,8 @@ struct knFlipSampleSecondaryParticles : public KernelBase {
const int n = KE * (k_ta * TA + k_wc * WC) * dt; // number of secondary particles
if (n == 0)
return;
static RandomStream mRand(9832);
Vec3 xi = Vec3(i, j, k) + mRand.getVec3(); // randomized offset uniform in cell
Vec3 xi = Vec3(i, j, k) + rand.getVec3(); // randomized offset uniform in cell
Vec3 vi = v.getInterpolated(xi);
Vec3 dir = dt * vi; // direction of movement of current particle
Vec3 e1 = getNormalized(Vec3(dir.z, 0, -dir.x)); // perpendicular to dir
@@ -803,9 +814,9 @@ struct knFlipSampleSecondaryParticles : public KernelBase {
cross(e1, dir)); // perpendicular to dir and e1, so e1 and e1 create reference plane
for (int di = 0; di < n; di++) {
const Real r = Real(0.5) * sqrt(mRand.getReal()); // distance to cylinder axis
const Real theta = mRand.getReal() * Real(2) * M_PI; // azimuth
const Real h = mRand.getReal() * norm(dt * vi); // distance to reference plane
const Real r = Real(0.5) * sqrt(rand.getReal()); // distance to cylinder axis
const Real theta = rand.getReal() * Real(2) * M_PI; // azimuth
const Real h = rand.getReal() * norm(dt * vi); // distance to reference plane
Vec3 xd = xi + r * cos(theta) * e1 + r * sin(theta) * e2 + h * getNormalized(vi);
if (!flags.is3D())
xd.z = 0;
@@ -815,7 +826,7 @@ struct knFlipSampleSecondaryParticles : public KernelBase {
vi; // init velocity of new particle
Real temp = (KE + TA + WC) / 3;
l_sec[l_sec.size() - 1] = ((lMax - lMin) * temp) + lMin +
mRand.getReal() * 0.1; // init lifetime of new particle
rand.getReal() * 0.1; // init lifetime of new particle
// init type of new particle
if (neighborRatio(i, j, k) < c_s) {
@@ -914,6 +925,11 @@ struct knFlipSampleSecondaryParticles : public KernelBase {
return itype;
}
typedef int type16;
inline RandomStream &getArg17()
{
return rand;
}
typedef RandomStream type17;
void runMessage()
{
debMsg("Executing kernel knFlipSampleSecondaryParticles ", 3);
@@ -947,7 +963,8 @@ struct knFlipSampleSecondaryParticles : public KernelBase {
k_ta,
k_wc,
dt,
itype);
itype,
rand);
}
const FlagGrid &flags;
const MACGrid &v;
@@ -966,6 +983,7 @@ struct knFlipSampleSecondaryParticles : public KernelBase {
const Real k_wc;
const Real dt;
const int itype;
RandomStream &rand;
};
void flipSampleSecondaryParticles(const std::string mode,
@@ -992,6 +1010,9 @@ void flipSampleSecondaryParticles(const std::string mode,
if (dt <= 0)
timestep = flags.getParent()->getDt();
/* Every particle needs to get a different random offset. */
RandomStream rand(pts_sec.getSeed());
if (mode == "single") {
knFlipSampleSecondaryParticles(flags,
v,
@@ -1009,7 +1030,8 @@ void flipSampleSecondaryParticles(const std::string mode,
k_ta,
k_wc,
timestep,
itype);
itype,
rand);
}
else if (mode == "multiple") {
knFlipSampleSecondaryParticlesMoreCylinders(flags,
@@ -1028,7 +1050,8 @@ void flipSampleSecondaryParticles(const std::string mode,
k_ta,
k_wc,
timestep,
itype);
itype,
rand);
}
else {
throw std::invalid_argument("Unknown mode: use \"single\" or \"multiple\" instead!");