diff options
Diffstat (limited to 'extern/mantaflow')
23 files changed, 615 insertions, 169 deletions
diff --git a/extern/mantaflow/CMakeLists.txt b/extern/mantaflow/CMakeLists.txt index d2f78f7ab6c..bdee06349d2 100644 --- a/extern/mantaflow/CMakeLists.txt +++ b/extern/mantaflow/CMakeLists.txt @@ -54,6 +54,10 @@ if(WITH_OPENVDB) add_definitions(-DOPENVDB_STATICLIB) endif() +if(WIN32) + add_definitions(-D_USE_MATH_DEFINES) +endif() + set(INC ${MANTA_PP} ${MANTA_PP}/fileio @@ -80,22 +84,11 @@ endif() if(WITH_OPENVDB) list(APPEND INC_SYS - ${BOOST_INCLUDE_DIR} - ${OPENEXR_INCLUDE_DIRS} ${OPENVDB_INCLUDE_DIRS} ) list(APPEND LIB ${OPENVDB_LIBRARIES} - ${OPENEXR_LIBRARIES} - ${ZLIB_LIBRARIES} - ${BOOST_LIBRARIES} ) - if(WITH_OPENVDB_BLOSC) - list(APPEND LIB - ${BLOSC_LIBRARIES} - ${ZLIB_LIBRARIES} - ) - endif() endif() set(SRC @@ -113,6 +106,7 @@ set(SRC ${MANTA_PP}/fastmarch.cpp ${MANTA_PP}/fastmarch.h ${MANTA_PP}/fastmarch.h.reg.cpp + ${MANTA_PP}/fileio/ioutil.cpp ${MANTA_PP}/fileio/iogrids.cpp ${MANTA_PP}/fileio/iomeshes.cpp ${MANTA_PP}/fileio/ioparticles.cpp diff --git a/extern/mantaflow/UPDATE.sh b/extern/mantaflow/UPDATE.sh index 28d96cdf6d8..0db18d03094 100644 --- a/extern/mantaflow/UPDATE.sh +++ b/extern/mantaflow/UPDATE.sh @@ -47,17 +47,17 @@ if [[ "$CLEAN_REPOSITORY" -eq "1" ]]; then git checkout develop fi -MANTA_BUILD_PATH=$MANTA_INSTALLATION/mantaflowgit/build_blender/ +MANTA_BUILD_PATH=$MANTA_INSTALLATION/build_blender/ mkdir -p $MANTA_BUILD_PATH cd $MANTA_BUILD_PATH -cmake .. -DGUI=OFF -DOPENMP=$USE_OMP -DTBB=$USE_TBB -DBLENDER=ON -DPREPDEBUG=ON && make -j8 +cmake ../mantaflowgit -DGUI=OFF -DOPENMP=$USE_OMP -DTBB=$USE_TBB -DBLENDER=ON -DPREPDEBUG=ON && make -j8 # ==================== 3) COPY MANTAFLOW FILES TO BLENDER ROOT =========================== mkdir -p $BLENDER_INSTALLATION/blender/tmp/dependencies/ && cp -Rf $MANTA_INSTALLATION/mantaflowgit/dependencies/cnpy "$_" mkdir -p $BLENDER_INSTALLATION/blender/tmp/helper/ && cp -Rf $MANTA_INSTALLATION/mantaflowgit/source/util "$_" mkdir -p $BLENDER_INSTALLATION/blender/tmp/helper/ && cp -Rf $MANTA_INSTALLATION/mantaflowgit/source/pwrapper "$_" -mkdir -p $BLENDER_INSTALLATION/blender/tmp/preprocessed/ && cp -Rf $MANTA_INSTALLATION/mantaflowgit/build_blender/pp/source/. "$_" +mkdir -p $BLENDER_INSTALLATION/blender/tmp/preprocessed/ && cp -Rf $MANTA_INSTALLATION/build_blender/pp/source/. "$_" # Remove some files that are not need in Blender rm $BLENDER_INSTALLATION/blender/tmp/dependencies/cnpy/example1.cpp diff --git a/extern/mantaflow/helper/pwrapper/registry.cpp b/extern/mantaflow/helper/pwrapper/registry.cpp index 3cdc2248b98..5b313e52fb4 100644 --- a/extern/mantaflow/helper/pwrapper/registry.cpp +++ b/extern/mantaflow/helper/pwrapper/registry.cpp @@ -192,7 +192,7 @@ int cbDisableConstructor(PyObject *self, PyObject *args, PyObject *kwds) return -1; } -PyMODINIT_FUNC PyInit_Main(void) +PyMODINIT_FUNC PyInit_manta_main(void) { MantaEnsureRegistration(); #if PY_MAJOR_VERSION >= 3 @@ -567,7 +567,7 @@ void WrapperRegistry::construct(const string &scriptname, const vector<string> & registerDummyTypes(); // work around for certain gcc versions, cast to char* - PyImport_AppendInittab((char *)gDefaultModuleName.c_str(), PyInit_Main); + PyImport_AppendInittab((char *)gDefaultModuleName.c_str(), PyInit_manta_main); } inline PyObject *castPy(PyTypeObject *p) diff --git a/extern/mantaflow/helper/pwrapper/registry.h b/extern/mantaflow/helper/pwrapper/registry.h index 139863df85d..d9d2bbb624b 100644 --- a/extern/mantaflow/helper/pwrapper/registry.h +++ b/extern/mantaflow/helper/pwrapper/registry.h @@ -62,7 +62,7 @@ void MantaEnsureRegistration(); #ifdef BLENDER # ifdef PyMODINIT_FUNC -PyMODINIT_FUNC PyInit_Main(void); +PyMODINIT_FUNC PyInit_manta_main(void); # endif #endif diff --git a/extern/mantaflow/helper/util/rcmatrix.h b/extern/mantaflow/helper/util/rcmatrix.h index 39951cece2e..3a466952ac2 100644 --- a/extern/mantaflow/helper/util/rcmatrix.h +++ b/extern/mantaflow/helper/util/rcmatrix.h @@ -17,16 +17,13 @@ // link to omp & tbb for now #if OPENMP == 1 || TBB == 1 -# define MANTA_ENABLE_PARALLEL 0 +# define MANTA_ENABLE_PARALLEL 1 // allow the preconditioner to be computed in parallel? (can lead to slightly non-deterministic // results) # define MANTA_ENABLE_PARALLEL_PC 0 -// use c++11 code? -# define MANTA_USE_CPP11 1 #else # define MANTA_ENABLE_PARALLEL 0 # define MANTA_ENABLE_PARALLEL_PC 0 -# define MANTA_USE_CPP11 0 #endif #if MANTA_ENABLE_PARALLEL == 1 @@ -503,11 +500,7 @@ template<class N, class T> struct RCMatrix { for (Iterator it = row_begin(i); it; ++it) result.fix_element(it.index(), i, it.value()); } -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else return result; -#endif } RCMatrix getKtK() const @@ -532,12 +525,7 @@ template<class N, class T> struct RCMatrix { } } } - parallel_end -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else - return result; -#endif + parallel_end return result; } RCMatrix operator*(const RCMatrix &m) const @@ -561,12 +549,7 @@ template<class N, class T> struct RCMatrix { } } } - parallel_end -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else - return result; -#endif + parallel_end return result; } RCMatrix sqrt() const @@ -581,12 +564,7 @@ template<class N, class T> struct RCMatrix { result.set_element(i, j, std::sqrt(it_A.value())); } } - parallel_end -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else - return result; -#endif + parallel_end return result; } RCMatrix operator*(const double k) const @@ -601,12 +579,7 @@ template<class N, class T> struct RCMatrix { result.add_to_element(i, j, it_A.value() * k); } } - parallel_end -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else - return result; -#endif + parallel_end return result; } RCMatrix applyKernel(const RCMatrix &kernel, const int nx, const int ny) const @@ -640,12 +613,7 @@ template<class N, class T> struct RCMatrix { } } } - parallel_end -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else - return result; -#endif + parallel_end return result; } RCMatrix applyHorizontalKernel(const RCMatrix &kernel, const int nx, const int ny) const @@ -679,12 +647,7 @@ template<class N, class T> struct RCMatrix { } } } - parallel_end -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else - return result; -#endif + parallel_end return result; } RCMatrix applyVerticalKernel(const RCMatrix &kernel, const int nx, const int ny) const @@ -718,12 +681,7 @@ template<class N, class T> struct RCMatrix { } } } - parallel_end -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else - return result; -#endif + parallel_end return result; } RCMatrix applySeparableKernel(const RCMatrix &kernelH, @@ -747,11 +705,7 @@ template<class N, class T> struct RCMatrix { { std::vector<T> result(n, 0.0); multiply(rhs, result); -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else return result; -#endif } void multiply(const std::vector<T> &rhs, std::vector<T> &result) const { @@ -832,11 +786,7 @@ template<class N, class T> struct RCMatrix { for (N i = 0; i < result.size(); i++) { result[i] = std::abs(result[i] - rhs[i]); } -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else return result; -#endif } T norm() const { @@ -1024,11 +974,7 @@ template<class N, class T> struct RCFixedMatrix { { std::vector<T> result(n, 0.0); multiply(rhs, result); -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else return result; -#endif } void multiply(const std::vector<T> &rhs, std::vector<T> &result) const { @@ -1064,12 +1010,7 @@ template<class N, class T> struct RCFixedMatrix { } } } - parallel_end -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else - return result; -#endif + parallel_end return result; } RCMatrix<N, T> toRCMatrix() const @@ -1087,12 +1028,7 @@ template<class N, class T> struct RCFixedMatrix { result.matrix[i]->value[j] = value[rowstart[i] + j]; } } - parallel_end -#if MANTA_USE_CPP11 == 1 - return std::move(result); -#else - return result; -#endif + parallel_end return result; } }; diff --git a/extern/mantaflow/helper/util/vectorbase.h b/extern/mantaflow/helper/util/vectorbase.h index a3135431eb3..41584663a0f 100644 --- a/extern/mantaflow/helper/util/vectorbase.h +++ b/extern/mantaflow/helper/util/vectorbase.h @@ -15,7 +15,7 @@ #define _VECTORBASE_H // get rid of windos min/max defines -#if defined(WIN32) || defined(_WIN32) +#if (defined(WIN32) || defined(_WIN32)) && !defined(NOMINMAX) # define NOMINMAX #endif diff --git a/extern/mantaflow/preprocessed/conjugategrad.cpp b/extern/mantaflow/preprocessed/conjugategrad.cpp index ac317402a49..de522d0b8f4 100644 --- a/extern/mantaflow/preprocessed/conjugategrad.cpp +++ b/extern/mantaflow/preprocessed/conjugategrad.cpp @@ -650,7 +650,7 @@ void cgSolveDiffusion(const FlagGrid &flags, << ", res:" << gcg->getSigma(), CG_DEBUGLEVEL); } - else if ((grid.getType() & GridBase::TypeVec3) || (grid.getType() & GridBase::TypeVec3)) { + else if ((grid.getType() & GridBase::TypeVec3) || (grid.getType() & GridBase::TypeMAC)) { Grid<Vec3> &vec = ((Grid<Vec3> &)grid); Grid<Real> u(parent); diff --git a/extern/mantaflow/preprocessed/fileio/iogrids.cpp b/extern/mantaflow/preprocessed/fileio/iogrids.cpp index 2f6cdaa6209..acd1bda5174 100644 --- a/extern/mantaflow/preprocessed/fileio/iogrids.cpp +++ b/extern/mantaflow/preprocessed/fileio/iogrids.cpp @@ -298,7 +298,7 @@ template<class T> void writeGridRaw(const string &name, Grid<T> *grid) debMsg("writing grid " << grid->getName() << " to raw file " << name, 1); #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "wb1"); // do some compression + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); // do some compression if (!gzf) errMsg("writeGridRaw: can't open file " << name); gzwrite(gzf, &((*grid)[0]), sizeof(T) * grid->getSizeX() * grid->getSizeY() * grid->getSizeZ()); @@ -313,7 +313,7 @@ template<class T> void readGridRaw(const string &name, Grid<T> *grid) debMsg("reading grid " << grid->getName() << " from raw file " << name, 1); #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "rb"); + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "rb"); if (!gzf) errMsg("readGridRaw: can't open file " << name); @@ -350,7 +350,7 @@ void getUniFileSize(const string &name, int &x, int &y, int &z, int *t, std::str { x = y = z = 0; #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "rb"); + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "rb"); if (gzf) { char ID[5] = {0, 0, 0, 0, 0}; gzread(gzf, ID, 4); @@ -499,7 +499,7 @@ template<class T> void writeGridUni(const string &name, Grid<T> *grid) else errMsg("writeGridUni: unknown element type"); - gzFile gzf = gzopen(name.c_str(), "wb1"); // do some compression + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); // do some compression if (!gzf) errMsg("writeGridUni: can't open file " << name); @@ -527,7 +527,7 @@ template<class T> void readGridUni(const string &name, Grid<T> *grid) debMsg("Reading grid " << grid->getName() << " from uni file " << name, 1); #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "rb"); + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "rb"); if (!gzf) errMsg("readGridUni: can't open file " << name); @@ -736,7 +736,7 @@ template<class T> void writeGrid4dUni(const string &name, Grid4d<T> *grid) else errMsg("writeGrid4dUni: unknown element type"); - gzFile gzf = gzopen(name.c_str(), "wb1"); // do some compression + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); // do some compression if (!gzf) errMsg("writeGrid4dUni: can't open file " << name); @@ -778,7 +778,7 @@ void readGrid4dUni( // optionally - reuse file handle, if valid one is passed in fileHandle pointer... if ((!fileHandle) || (fileHandle && (*fileHandle == NULL))) { - gzf = gzopen(name.c_str(), "rb"); + gzf = (gzFile)safeGzopen(name.c_str(), "rb"); if (!gzf) errMsg("readGrid4dUni: can't open file " << name); @@ -905,7 +905,7 @@ template<class T> void writeGrid4dRaw(const string &name, Grid4d<T> *grid) debMsg("writing grid4d " << grid->getName() << " to raw file " << name, 1); #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "wb1"); // do some compression + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); // do some compression if (!gzf) errMsg("writeGrid4dRaw: can't open file " << name); gzwrite(gzf, @@ -922,7 +922,7 @@ template<class T> void readGrid4dRaw(const string &name, Grid4d<T> *grid) debMsg("reading grid4d " << grid->getName() << " from raw file " << name, 1); #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "rb"); + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "rb"); if (!gzf) errMsg("readGrid4dRaw: can't open file " << name); @@ -954,6 +954,79 @@ template<class T> void readGridVDB(const string &name, Grid<T> *grid) 1); } +template<> void writeGridVDB(const string &name, Grid<int> *grid) +{ + debMsg("Writing int grid " << grid->getName() << " to vdb file " << name, 1); + + // Create an empty int32-point grid with background value 0. + openvdb::initialize(); + openvdb::Int32Grid::Ptr gridVDB = openvdb::Int32Grid::create(); + gridVDB->setTransform( + openvdb::math::Transform::createLinearTransform(1. / grid->getSizeX())); // voxel size + + // Get an accessor for coordinate-based access to voxels. + openvdb::Int32Grid::Accessor accessor = gridVDB->getAccessor(); + + gridVDB->setGridClass(openvdb::GRID_UNKNOWN); + + // Name the grid "density". + gridVDB->setName(grid->getName()); + + openvdb::io::File file(name); + + FOR_IJK(*grid) + { + openvdb::Coord xyz(i, j, k); + accessor.setValue(xyz, (*grid)(i, j, k)); + } + + // Add the grid pointer to a container. + openvdb::GridPtrVec gridsVDB; + gridsVDB.push_back(gridVDB); + + // Write out the contents of the container. + file.write(gridsVDB); + file.close(); +} + +template<> void readGridVDB(const string &name, Grid<int> *grid) +{ + debMsg("Reading int grid " << grid->getName() << " from vdb file " << name, 1); + + openvdb::initialize(); + openvdb::io::File file(name); + file.open(); + + openvdb::GridBase::Ptr baseGrid; + for (openvdb::io::File::NameIterator nameIter = file.beginName(); nameIter != file.endName(); + ++nameIter) { +# ifndef BLENDER + // Read in only the grid we are interested in. + if (nameIter.gridName() == grid->getName()) { + baseGrid = file.readGrid(nameIter.gridName()); + } + else { + debMsg("skipping grid " << nameIter.gridName(), 1); + } +# else + // For Blender, skip name check and pick first grid from loop + baseGrid = file.readGrid(nameIter.gridName()); + break; +# endif + } + file.close(); + openvdb::Int32Grid::Ptr gridVDB = openvdb::gridPtrCast<openvdb::Int32Grid>(baseGrid); + + openvdb::Int32Grid::Accessor accessor = gridVDB->getAccessor(); + + FOR_IJK(*grid) + { + openvdb::Coord xyz(i, j, k); + int v = accessor.getValue(xyz); + (*grid)(i, j, k) = v; + } +} + template<> void writeGridVDB(const string &name, Grid<Real> *grid) { debMsg("Writing real grid " << grid->getName() << " to vdb file " << name, 1); diff --git a/extern/mantaflow/preprocessed/fileio/iomeshes.cpp b/extern/mantaflow/preprocessed/fileio/iomeshes.cpp index 79a9e76da3f..906b849fffb 100644 --- a/extern/mantaflow/preprocessed/fileio/iomeshes.cpp +++ b/extern/mantaflow/preprocessed/fileio/iomeshes.cpp @@ -158,7 +158,7 @@ void readBobjFile(const string &name, Mesh *mesh, bool append) const Real dx = mesh->getParent()->getDx(); const Vec3 gs = toVec3(mesh->getParent()->getGridSize()); - gzFile gzf = gzopen(name.c_str(), "rb1"); // do some compression + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "rb1"); // do some compression if (!gzf) errMsg("readBobj: unable to open file"); @@ -213,7 +213,7 @@ void writeBobjFile(const string &name, Mesh *mesh) const Real dx = mesh->getParent()->getDx(); const Vec3i gs = mesh->getParent()->getGridSize(); - gzFile gzf = gzopen(name.c_str(), "wb1"); // do some compression + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); // do some compression if (!gzf) errMsg("writeBobj: unable to open file"); @@ -412,7 +412,7 @@ template<class T> void readMdataUni(const std::string &name, MeshDataImpl<T> *md debMsg("reading mesh data " << mdata->getName() << " from uni file " << name, 1); #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "rb"); + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "rb"); if (!gzf) errMsg("can't open file " << name); @@ -460,7 +460,7 @@ template<class T> void writeMdataUni(const std::string &name, MeshDataImpl<T> *m MuTime stamp; head.timestamp = stamp.time; - gzFile gzf = gzopen(name.c_str(), "wb1"); // do some compression + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); // do some compression if (!gzf) errMsg("can't open file " << name); gzwrite(gzf, ID, 4); diff --git a/extern/mantaflow/preprocessed/fileio/ioparticles.cpp b/extern/mantaflow/preprocessed/fileio/ioparticles.cpp index 432cbc9f100..2eab485beb3 100644 --- a/extern/mantaflow/preprocessed/fileio/ioparticles.cpp +++ b/extern/mantaflow/preprocessed/fileio/ioparticles.cpp @@ -176,7 +176,7 @@ void writeParticlesUni(const std::string &name, const BasicParticleSystem *parts MuTime stamp; head.timestamp = stamp.time; - gzFile gzf = gzopen(name.c_str(), "wb1"); // do some compression + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); // do some compression if (!gzf) errMsg("can't open file " << name); @@ -206,7 +206,7 @@ void readParticlesUni(const std::string &name, BasicParticleSystem *parts) debMsg("reading particles " << parts->getName() << " from uni file " << name, 1); #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "rb"); + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "rb"); if (!gzf) errMsg("can't open file " << name); @@ -273,7 +273,7 @@ template<class T> void writePdataUni(const std::string &name, ParticleDataImpl<T MuTime stamp; head.timestamp = stamp.time; - gzFile gzf = gzopen(name.c_str(), "wb1"); // do some compression + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); // do some compression if (!gzf) errMsg("can't open file " << name); gzwrite(gzf, ID, 4); @@ -299,7 +299,7 @@ template<class T> void readPdataUni(const std::string &name, ParticleDataImpl<T> debMsg("reading particle data " << pdata->getName() << " from uni file " << name, 1); #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "rb"); + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "rb"); if (!gzf) errMsg("can't open file " << name); @@ -310,6 +310,8 @@ template<class T> void readPdataUni(const std::string &name, ParticleDataImpl<T> UniPartHeader head; assertMsg(gzread(gzf, &head, sizeof(UniPartHeader)) == sizeof(UniPartHeader), "can't read file, no header present"); + pdata->resize(head.dim); + assertMsg(head.dim == pdata->size(), "pdata size doesn't match"); # if FLOATINGPOINT_PRECISION != 1 ParticleDataImpl<T> temp(pdata->getParent()); diff --git a/extern/mantaflow/preprocessed/fileio/ioutil.cpp b/extern/mantaflow/preprocessed/fileio/ioutil.cpp new file mode 100644 index 00000000000..e04633c5634 --- /dev/null +++ b/extern/mantaflow/preprocessed/fileio/ioutil.cpp @@ -0,0 +1,60 @@ + + +// DO NOT EDIT ! +// This file is generated using the MantaFlow preprocessor (prep generate). + +/****************************************************************************** + * + * MantaFlow fluid solver framework + * Copyright 2011-2020 Tobias Pfaff, Nils Thuerey + * + * This program is free software, distributed under the terms of the + * Apache License, Version 2.0 + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Helper functions to handle file IO + * + ******************************************************************************/ + +#include "mantaio.h" + +#if NO_ZLIB != 1 +extern "C" { +# include <zlib.h> +} + +# if defined(WIN32) || defined(_WIN32) +# include <windows.h> +# include <string> +# endif + +using namespace std; + +namespace Manta { + +# if defined(WIN32) || defined(_WIN32) +static wstring stringToWstring(const char *str) +{ + const int length_wc = MultiByteToWideChar(CP_UTF8, 0, str, strlen(str), NULL, 0); + wstring strWide(length_wc, 0); + MultiByteToWideChar(CP_UTF8, 0, str, strlen(str), &strWide[0], length_wc); + return strWide; +} +# endif + +void *safeGzopen(const char *filename, const char *mode) +{ + gzFile gzfile; + +# if defined(WIN32) || defined(_WIN32) + wstring filenameWide = stringToWstring(filename); + gzfile = gzopen_w(filenameWide.c_str(), mode); +# else + gzfile = gzopen(filename, mode); +# endif + + return gzfile; +} +#endif + +} // namespace diff --git a/extern/mantaflow/preprocessed/fileio/mantaio.h b/extern/mantaflow/preprocessed/fileio/mantaio.h index 8bb0a5af6a4..fbfe4bdd5d4 100644 --- a/extern/mantaflow/preprocessed/fileio/mantaio.h +++ b/extern/mantaflow/preprocessed/fileio/mantaio.h @@ -76,6 +76,8 @@ template<class T> void readMdataUni(const std::string &name, MeshDataImpl<T> *md void getUniFileSize( const std::string &name, int &x, int &y, int &z, int *t = NULL, std::string *info = NULL); +void *safeGzopen(const char *filename, const char *mode); + } // namespace Manta #endif diff --git a/extern/mantaflow/preprocessed/fluidsolver.cpp b/extern/mantaflow/preprocessed/fluidsolver.cpp index 814d5444b15..8c48e60760b 100644 --- a/extern/mantaflow/preprocessed/fluidsolver.cpp +++ b/extern/mantaflow/preprocessed/fluidsolver.cpp @@ -136,9 +136,9 @@ FluidSolver::FluidSolver(Vec3i gridsize, int dim, int fourthDim) mDtMin(1.), mDtMax(1.), mFrameLength(1.), + mTimePerFrame(0.), mGridSize(gridsize), mDim(dim), - mTimePerFrame(0.), mLockDt(false), mFourthDim(fourthDim) { diff --git a/extern/mantaflow/preprocessed/gitinfo.h b/extern/mantaflow/preprocessed/gitinfo.h index 216aaf4e2bc..208d8008a7e 100644 --- a/extern/mantaflow/preprocessed/gitinfo.h +++ b/extern/mantaflow/preprocessed/gitinfo.h @@ -1,3 +1,3 @@ -#define MANTA_GIT_VERSION "commit 8d19f1096a4d8e115d9d1ec6024c65d53a94f47b" +#define MANTA_GIT_VERSION "commit 5fbd3d04381b21afce4a593d1fe2d9bc7bef5424" diff --git a/extern/mantaflow/preprocessed/grid.cpp b/extern/mantaflow/preprocessed/grid.cpp index c21d56d8879..f10052349d5 100644 --- a/extern/mantaflow/preprocessed/grid.cpp +++ b/extern/mantaflow/preprocessed/grid.cpp @@ -1244,15 +1244,67 @@ void PbRegister_gridMaxDiffVec3() } } +struct knCopyMacToVec3 : public KernelBase { + knCopyMacToVec3(MACGrid &source, Grid<Vec3> &target) + : KernelBase(&source, 0), source(source), target(target) + { + runMessage(); + run(); + } + inline void op(int i, int j, int k, MACGrid &source, Grid<Vec3> &target) const + { + target(i, j, k) = source(i, j, k); + } + inline MACGrid &getArg0() + { + return source; + } + typedef MACGrid type0; + inline Grid<Vec3> &getArg1() + { + return target; + } + typedef Grid<Vec3> type1; + void runMessage() + { + debMsg("Executing kernel knCopyMacToVec3 ", 3); + debMsg("Kernel range" + << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ", + 4); + }; + void operator()(const tbb::blocked_range<IndexInt> &__r) const + { + const int _maxX = maxX; + const int _maxY = maxY; + if (maxZ > 1) { + for (int k = __r.begin(); k != (int)__r.end(); k++) + for (int j = 0; j < _maxY; j++) + for (int i = 0; i < _maxX; i++) + op(i, j, k, source, target); + } + else { + const int k = 0; + for (int j = __r.begin(); j != (int)__r.end(); j++) + for (int i = 0; i < _maxX; i++) + op(i, j, k, source, target); + } + } + void run() + { + if (maxZ > 1) + tbb::parallel_for(tbb::blocked_range<IndexInt>(minZ, maxZ), *this); + else + tbb::parallel_for(tbb::blocked_range<IndexInt>(0, maxY), *this); + } + MACGrid &source; + Grid<Vec3> ⌖ +}; // simple helper functions to copy (convert) mac to vec3 , and levelset to real grids // (are assumed to be the same for running the test cases - in general they're not!) void copyMacToVec3(MACGrid &source, Grid<Vec3> &target) { - FOR_IJK(target) - { - target(i, j, k) = source(i, j, k); - } + knCopyMacToVec3(source, target); } static PyObject *_W_3(PyObject *_self, PyObject *_linargs, PyObject *_kwds) { @@ -1323,10 +1375,14 @@ void PbRegister_convertMacToVec3() } } -//! vec3->mac grid conversion , but with full resampling -void resampleVec3ToMac(Grid<Vec3> &source, MACGrid &target) -{ - FOR_IJK_BND(target, 1) +struct knResampleVec3ToMac : public KernelBase { + knResampleVec3ToMac(Grid<Vec3> &source, MACGrid &target) + : KernelBase(&source, 1), source(source), target(target) + { + runMessage(); + run(); + } + inline void op(int i, int j, int k, Grid<Vec3> &source, MACGrid &target) const { target(i, j, k)[0] = 0.5 * (source(i - 1, j, k)[0] + source(i, j, k))[0]; target(i, j, k)[1] = 0.5 * (source(i, j - 1, k)[1] + source(i, j, k))[1]; @@ -1334,6 +1390,55 @@ void resampleVec3ToMac(Grid<Vec3> &source, MACGrid &target) target(i, j, k)[2] = 0.5 * (source(i, j, k - 1)[2] + source(i, j, k))[2]; } } + inline Grid<Vec3> &getArg0() + { + return source; + } + typedef Grid<Vec3> type0; + inline MACGrid &getArg1() + { + return target; + } + typedef MACGrid type1; + void runMessage() + { + debMsg("Executing kernel knResampleVec3ToMac ", 3); + debMsg("Kernel range" + << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ", + 4); + }; + void operator()(const tbb::blocked_range<IndexInt> &__r) const + { + const int _maxX = maxX; + const int _maxY = maxY; + if (maxZ > 1) { + for (int k = __r.begin(); k != (int)__r.end(); k++) + for (int j = 1; j < _maxY; j++) + for (int i = 1; i < _maxX; i++) + op(i, j, k, source, target); + } + else { + const int k = 0; + for (int j = __r.begin(); j != (int)__r.end(); j++) + for (int i = 1; i < _maxX; i++) + op(i, j, k, source, target); + } + } + void run() + { + if (maxZ > 1) + tbb::parallel_for(tbb::blocked_range<IndexInt>(minZ, maxZ), *this); + else + tbb::parallel_for(tbb::blocked_range<IndexInt>(1, maxY), *this); + } + Grid<Vec3> &source; + MACGrid ⌖ +}; +//! vec3->mac grid conversion , but with full resampling + +void resampleVec3ToMac(Grid<Vec3> &source, MACGrid &target) +{ + knResampleVec3ToMac(source, target); } static PyObject *_W_5(PyObject *_self, PyObject *_linargs, PyObject *_kwds) { @@ -1367,13 +1472,66 @@ void PbRegister_resampleVec3ToMac() } } -//! mac->vec3 grid conversion , with full resampling -void resampleMacToVec3(MACGrid &source, Grid<Vec3> &target) -{ - FOR_IJK_BND(target, 1) +struct knResampleMacToVec3 : public KernelBase { + knResampleMacToVec3(MACGrid &source, Grid<Vec3> &target) + : KernelBase(&source, 1), source(source), target(target) + { + runMessage(); + run(); + } + inline void op(int i, int j, int k, MACGrid &source, Grid<Vec3> &target) const { target(i, j, k) = source.getCentered(i, j, k); } + inline MACGrid &getArg0() + { + return source; + } + typedef MACGrid type0; + inline Grid<Vec3> &getArg1() + { + return target; + } + typedef Grid<Vec3> type1; + void runMessage() + { + debMsg("Executing kernel knResampleMacToVec3 ", 3); + debMsg("Kernel range" + << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ", + 4); + }; + void operator()(const tbb::blocked_range<IndexInt> &__r) const + { + const int _maxX = maxX; + const int _maxY = maxY; + if (maxZ > 1) { + for (int k = __r.begin(); k != (int)__r.end(); k++) + for (int j = 1; j < _maxY; j++) + for (int i = 1; i < _maxX; i++) + op(i, j, k, source, target); + } + else { + const int k = 0; + for (int j = __r.begin(); j != (int)__r.end(); j++) + for (int i = 1; i < _maxX; i++) + op(i, j, k, source, target); + } + } + void run() + { + if (maxZ > 1) + tbb::parallel_for(tbb::blocked_range<IndexInt>(minZ, maxZ), *this); + else + tbb::parallel_for(tbb::blocked_range<IndexInt>(1, maxY), *this); + } + MACGrid &source; + Grid<Vec3> ⌖ +}; +//! mac->vec3 grid conversion , with full resampling + +void resampleMacToVec3(MACGrid &source, Grid<Vec3> &target) +{ + knResampleMacToVec3(source, target); } static PyObject *_W_6(PyObject *_self, PyObject *_linargs, PyObject *_kwds) { @@ -1407,12 +1565,65 @@ void PbRegister_resampleMacToVec3() } } -void copyLevelsetToReal(LevelsetGrid &source, Grid<Real> &target) -{ - FOR_IJK(target) +struct knCopyLevelsetToReal : public KernelBase { + knCopyLevelsetToReal(LevelsetGrid &source, Grid<Real> &target) + : KernelBase(&source, 0), source(source), target(target) + { + runMessage(); + run(); + } + inline void op(int i, int j, int k, LevelsetGrid &source, Grid<Real> &target) const { target(i, j, k) = source(i, j, k); } + inline LevelsetGrid &getArg0() + { + return source; + } + typedef LevelsetGrid type0; + inline Grid<Real> &getArg1() + { + return target; + } + typedef Grid<Real> type1; + void runMessage() + { + debMsg("Executing kernel knCopyLevelsetToReal ", 3); + debMsg("Kernel range" + << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ", + 4); + }; + void operator()(const tbb::blocked_range<IndexInt> &__r) const + { + const int _maxX = maxX; + const int _maxY = maxY; + if (maxZ > 1) { + for (int k = __r.begin(); k != (int)__r.end(); k++) + for (int j = 0; j < _maxY; j++) + for (int i = 0; i < _maxX; i++) + op(i, j, k, source, target); + } + else { + const int k = 0; + for (int j = __r.begin(); j != (int)__r.end(); j++) + for (int i = 0; i < _maxX; i++) + op(i, j, k, source, target); + } + } + void run() + { + if (maxZ > 1) + tbb::parallel_for(tbb::blocked_range<IndexInt>(minZ, maxZ), *this); + else + tbb::parallel_for(tbb::blocked_range<IndexInt>(0, maxY), *this); + } + LevelsetGrid &source; + Grid<Real> ⌖ +}; + +void copyLevelsetToReal(LevelsetGrid &source, Grid<Real> &target) +{ + knCopyLevelsetToReal(source, target); } static PyObject *_W_7(PyObject *_self, PyObject *_linargs, PyObject *_kwds) { @@ -1446,17 +1657,95 @@ void PbRegister_copyLevelsetToReal() } } -void copyVec3ToReal(Grid<Vec3> &source, - Grid<Real> &targetX, - Grid<Real> &targetY, - Grid<Real> &targetZ) -{ - FOR_IJK(source) +struct knCopyVec3ToReal : public KernelBase { + knCopyVec3ToReal(Grid<Vec3> &source, + Grid<Real> &targetX, + Grid<Real> &targetY, + Grid<Real> &targetZ) + : KernelBase(&source, 0), + source(source), + targetX(targetX), + targetY(targetY), + targetZ(targetZ) + { + runMessage(); + run(); + } + inline void op(int i, + int j, + int k, + Grid<Vec3> &source, + Grid<Real> &targetX, + Grid<Real> &targetY, + Grid<Real> &targetZ) const { targetX(i, j, k) = source(i, j, k).x; targetY(i, j, k) = source(i, j, k).y; targetZ(i, j, k) = source(i, j, k).z; } + inline Grid<Vec3> &getArg0() + { + return source; + } + typedef Grid<Vec3> type0; + inline Grid<Real> &getArg1() + { + return targetX; + } + typedef Grid<Real> type1; + inline Grid<Real> &getArg2() + { + return targetY; + } + typedef Grid<Real> type2; + inline Grid<Real> &getArg3() + { + return targetZ; + } + typedef Grid<Real> type3; + void runMessage() + { + debMsg("Executing kernel knCopyVec3ToReal ", 3); + debMsg("Kernel range" + << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ", + 4); + }; + void operator()(const tbb::blocked_range<IndexInt> &__r) const + { + const int _maxX = maxX; + const int _maxY = maxY; + if (maxZ > 1) { + for (int k = __r.begin(); k != (int)__r.end(); k++) + for (int j = 0; j < _maxY; j++) + for (int i = 0; i < _maxX; i++) + op(i, j, k, source, targetX, targetY, targetZ); + } + else { + const int k = 0; + for (int j = __r.begin(); j != (int)__r.end(); j++) + for (int i = 0; i < _maxX; i++) + op(i, j, k, source, targetX, targetY, targetZ); + } + } + void run() + { + if (maxZ > 1) + tbb::parallel_for(tbb::blocked_range<IndexInt>(minZ, maxZ), *this); + else + tbb::parallel_for(tbb::blocked_range<IndexInt>(0, maxY), *this); + } + Grid<Vec3> &source; + Grid<Real> &targetX; + Grid<Real> &targetY; + Grid<Real> &targetZ; +}; + +void copyVec3ToReal(Grid<Vec3> &source, + Grid<Real> &targetX, + Grid<Real> &targetY, + Grid<Real> &targetZ) +{ + knCopyVec3ToReal(source, targetX, targetY, targetZ); } static PyObject *_W_8(PyObject *_self, PyObject *_linargs, PyObject *_kwds) { @@ -1492,17 +1781,95 @@ void PbRegister_copyVec3ToReal() } } -void copyRealToVec3(Grid<Real> &sourceX, - Grid<Real> &sourceY, - Grid<Real> &sourceZ, - Grid<Vec3> &target) -{ - FOR_IJK(target) +struct knCopyRealToVec3 : public KernelBase { + knCopyRealToVec3(Grid<Real> &sourceX, + Grid<Real> &sourceY, + Grid<Real> &sourceZ, + Grid<Vec3> &target) + : KernelBase(&sourceX, 0), + sourceX(sourceX), + sourceY(sourceY), + sourceZ(sourceZ), + target(target) + { + runMessage(); + run(); + } + inline void op(int i, + int j, + int k, + Grid<Real> &sourceX, + Grid<Real> &sourceY, + Grid<Real> &sourceZ, + Grid<Vec3> &target) const { target(i, j, k).x = sourceX(i, j, k); target(i, j, k).y = sourceY(i, j, k); target(i, j, k).z = sourceZ(i, j, k); } + inline Grid<Real> &getArg0() + { + return sourceX; + } + typedef Grid<Real> type0; + inline Grid<Real> &getArg1() + { + return sourceY; + } + typedef Grid<Real> type1; + inline Grid<Real> &getArg2() + { + return sourceZ; + } + typedef Grid<Real> type2; + inline Grid<Vec3> &getArg3() + { + return target; + } + typedef Grid<Vec3> type3; + void runMessage() + { + debMsg("Executing kernel knCopyRealToVec3 ", 3); + debMsg("Kernel range" + << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ", + 4); + }; + void operator()(const tbb::blocked_range<IndexInt> &__r) const + { + const int _maxX = maxX; + const int _maxY = maxY; + if (maxZ > 1) { + for (int k = __r.begin(); k != (int)__r.end(); k++) + for (int j = 0; j < _maxY; j++) + for (int i = 0; i < _maxX; i++) + op(i, j, k, sourceX, sourceY, sourceZ, target); + } + else { + const int k = 0; + for (int j = __r.begin(); j != (int)__r.end(); j++) + for (int i = 0; i < _maxX; i++) + op(i, j, k, sourceX, sourceY, sourceZ, target); + } + } + void run() + { + if (maxZ > 1) + tbb::parallel_for(tbb::blocked_range<IndexInt>(minZ, maxZ), *this); + else + tbb::parallel_for(tbb::blocked_range<IndexInt>(0, maxY), *this); + } + Grid<Real> &sourceX; + Grid<Real> &sourceY; + Grid<Real> &sourceZ; + Grid<Vec3> ⌖ +}; + +void copyRealToVec3(Grid<Real> &sourceX, + Grid<Real> &sourceY, + Grid<Real> &sourceZ, + Grid<Vec3> &target) +{ + knCopyRealToVec3(sourceX, sourceY, sourceZ, target); } static PyObject *_W_9(PyObject *_self, PyObject *_linargs, PyObject *_kwds) { diff --git a/extern/mantaflow/preprocessed/levelset.cpp b/extern/mantaflow/preprocessed/levelset.cpp index dcc10718d71..a6f9fb40f3f 100644 --- a/extern/mantaflow/preprocessed/levelset.cpp +++ b/extern/mantaflow/preprocessed/levelset.cpp @@ -773,9 +773,9 @@ void LevelsetGrid::createMesh(Mesh &mesh) Grid<int> edgeVY(mParent); Grid<int> edgeVZ(mParent); - for (int i = 0; i < mSize.x - 1; i++) + for (int k = 0; k < mSize.z - 1; k++) for (int j = 0; j < mSize.y - 1; j++) - for (int k = 0; k < mSize.z - 1; k++) { + for (int i = 0; i < mSize.x - 1; i++) { Real value[8] = {get(i, j, k), get(i + 1, j, k), get(i + 1, j + 1, k), diff --git a/extern/mantaflow/preprocessed/particle.cpp b/extern/mantaflow/preprocessed/particle.cpp index 478f1417109..41af5d3ba81 100644 --- a/extern/mantaflow/preprocessed/particle.cpp +++ b/extern/mantaflow/preprocessed/particle.cpp @@ -182,7 +182,7 @@ void BasicParticleSystem::writeParticlesText(const string name) const void BasicParticleSystem::writeParticlesRawPositionsGz(const string name) const { #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "wb1"); + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); if (!gzf) errMsg("can't open file " << name); for (IndexInt i = 0; i < this->size(); ++i) { @@ -198,7 +198,7 @@ void BasicParticleSystem::writeParticlesRawPositionsGz(const string name) const void BasicParticleSystem::writeParticlesRawVelocityGz(const string name) const { #if NO_ZLIB != 1 - gzFile gzf = gzopen(name.c_str(), "wb1"); + gzFile gzf = (gzFile)safeGzopen(name.c_str(), "wb1"); if (!gzf) errMsg("can't open file " << name); if (mPdataVec3.size() < 1) diff --git a/extern/mantaflow/preprocessed/plugin/extforces.cpp b/extern/mantaflow/preprocessed/plugin/extforces.cpp index df0ddb15b33..36221fbbc10 100644 --- a/extern/mantaflow/preprocessed/plugin/extforces.cpp +++ b/extern/mantaflow/preprocessed/plugin/extforces.cpp @@ -1335,7 +1335,7 @@ struct KnConfForce : public KernelBase { void vorticityConfinement(MACGrid &vel, const FlagGrid &flags, - Real strengthGlobal = 0, + Real strength = 0, const Grid<Real> *strengthCell = NULL) { Grid<Vec3> velCenter(flags.getParent()), curl(flags.getParent()), force(flags.getParent()); @@ -1344,7 +1344,7 @@ void vorticityConfinement(MACGrid &vel, GetCentered(velCenter, vel); CurlOp(velCenter, curl); GridNorm(norm, curl); - KnConfForce(force, norm, curl, strengthGlobal, strengthCell); + KnConfForce(force, norm, curl, strength, strengthCell); KnApplyForceField(flags, vel, force, NULL, true, false); } static PyObject *_W_8(PyObject *_self, PyObject *_linargs, PyObject *_kwds) @@ -1359,11 +1359,11 @@ static PyObject *_W_8(PyObject *_self, PyObject *_linargs, PyObject *_kwds) ArgLocker _lock; MACGrid &vel = *_args.getPtr<MACGrid>("vel", 0, &_lock); const FlagGrid &flags = *_args.getPtr<FlagGrid>("flags", 1, &_lock); - Real strengthGlobal = _args.getOpt<Real>("strengthGlobal", 2, 0, &_lock); + Real strength = _args.getOpt<Real>("strength", 2, 0, &_lock); const Grid<Real> *strengthCell = _args.getPtrOpt<Grid<Real>>( "strengthCell", 3, NULL, &_lock); _retval = getPyNone(); - vorticityConfinement(vel, flags, strengthGlobal, strengthCell); + vorticityConfinement(vel, flags, strength, strengthCell); _args.check(); } pbFinalizePlugin(parent, "vorticityConfinement", !noTiming); diff --git a/extern/mantaflow/preprocessed/plugin/flip.cpp b/extern/mantaflow/preprocessed/plugin/flip.cpp index f6d082900b5..8ac167c1166 100644 --- a/extern/mantaflow/preprocessed/plugin/flip.cpp +++ b/extern/mantaflow/preprocessed/plugin/flip.cpp @@ -18,6 +18,7 @@ ******************************************************************************/ #include "particle.h" +#include "general.h" #include "grid.h" #include "commonkernels.h" #include "randomstream.h" @@ -1407,7 +1408,6 @@ struct correctLevelset : public KernelBase { { if (rAcc(i, j, k) <= VECTOR_EPSILON) return; // outside nothing happens - Real x = pAcc(i, j, k).x; // create jacobian of pAcc via central differences Matrix3x3f jacobian = Matrix3x3f(0.5 * (pAcc(i + 1, j, k).x - pAcc(i - 1, j, k).x), @@ -1430,9 +1430,9 @@ struct correctLevelset : public KernelBase { Real t = (t_high - maxEV) / (t_high - t_low); correction = t * t * t - 3 * t * t + 3 * t; } - correction = (correction < 0) ? - 0 : - correction; // enforce correction factor to [0,1] (not explicitly in paper) + correction = clamp(correction, + Real(0), + Real(1)); // enforce correction factor to [0,1] (not explicitly in paper) const Vec3 gridPos = Vec3(i, j, k) + Vec3(0.5); // shifted by half cell const Real correctedPhi = fabs(norm(gridPos - pAcc(i, j, k))) - rAcc(i, j, k) * correction; diff --git a/extern/mantaflow/preprocessed/plugin/fluidguiding.cpp b/extern/mantaflow/preprocessed/plugin/fluidguiding.cpp index 13383581123..18a5a37771f 100644 --- a/extern/mantaflow/preprocessed/plugin/fluidguiding.cpp +++ b/extern/mantaflow/preprocessed/plugin/fluidguiding.cpp @@ -381,7 +381,6 @@ void getSpiralVelocity(const FlagGrid &flags, nz = flags.getSizeZ(); Real midX = 0.5 * (Real)(nx - 1); Real midY = 0.5 * (Real)(ny - 1); - Real midZ = 0.5 * (Real)(nz - 1); for (int i = 0; i < nx; i++) { for (int j = 0; j < ny; j++) { for (int k = 0; k < nz; k++) { diff --git a/extern/mantaflow/preprocessed/plugin/initplugins.cpp b/extern/mantaflow/preprocessed/plugin/initplugins.cpp index 519aa41a08d..7a765813f9f 100644 --- a/extern/mantaflow/preprocessed/plugin/initplugins.cpp +++ b/extern/mantaflow/preprocessed/plugin/initplugins.cpp @@ -1829,13 +1829,15 @@ struct KnUpdateFlagsObs : public KernelBase { const MACGrid *fractions, const Grid<Real> &phiObs, const Grid<Real> *phiOut, - const Grid<Real> *phiIn) - : KernelBase(&flags, 1), + const Grid<Real> *phiIn, + int boundaryWidth) + : KernelBase(&flags, boundaryWidth), flags(flags), fractions(fractions), phiObs(phiObs), phiOut(phiOut), - phiIn(phiIn) + phiIn(phiIn), + boundaryWidth(boundaryWidth) { runMessage(); run(); @@ -1847,7 +1849,8 @@ struct KnUpdateFlagsObs : public KernelBase { const MACGrid *fractions, const Grid<Real> &phiObs, const Grid<Real> *phiOut, - const Grid<Real> *phiIn) const + const Grid<Real> *phiIn, + int boundaryWidth) const { bool isObs = false; @@ -1910,6 +1913,11 @@ struct KnUpdateFlagsObs : public KernelBase { return phiIn; } typedef Grid<Real> type4; + inline int &getArg5() + { + return boundaryWidth; + } + typedef int type5; void runMessage() { debMsg("Executing kernel KnUpdateFlagsObs ", 3); @@ -1923,15 +1931,15 @@ struct KnUpdateFlagsObs : public KernelBase { const int _maxY = maxY; if (maxZ > 1) { for (int k = __r.begin(); k != (int)__r.end(); k++) - for (int j = 1; j < _maxY; j++) - for (int i = 1; i < _maxX; i++) - op(i, j, k, flags, fractions, phiObs, phiOut, phiIn); + for (int j = boundaryWidth; j < _maxY; j++) + for (int i = boundaryWidth; i < _maxX; i++) + op(i, j, k, flags, fractions, phiObs, phiOut, phiIn, boundaryWidth); } else { const int k = 0; for (int j = __r.begin(); j != (int)__r.end(); j++) - for (int i = 1; i < _maxX; i++) - op(i, j, k, flags, fractions, phiObs, phiOut, phiIn); + for (int i = boundaryWidth; i < _maxX; i++) + op(i, j, k, flags, fractions, phiObs, phiOut, phiIn, boundaryWidth); } } void run() @@ -1939,13 +1947,14 @@ struct KnUpdateFlagsObs : public KernelBase { if (maxZ > 1) tbb::parallel_for(tbb::blocked_range<IndexInt>(minZ, maxZ), *this); else - tbb::parallel_for(tbb::blocked_range<IndexInt>(1, maxY), *this); + tbb::parallel_for(tbb::blocked_range<IndexInt>(boundaryWidth, maxY), *this); } FlagGrid &flags; const MACGrid *fractions; const Grid<Real> &phiObs; const Grid<Real> *phiOut; const Grid<Real> *phiIn; + int boundaryWidth; }; //! update obstacle and outflow flags from levelsets @@ -1954,9 +1963,10 @@ void setObstacleFlags(FlagGrid &flags, const Grid<Real> &phiObs, const MACGrid *fractions = NULL, const Grid<Real> *phiOut = NULL, - const Grid<Real> *phiIn = NULL) + const Grid<Real> *phiIn = NULL, + int boundaryWidth = 1) { - KnUpdateFlagsObs(flags, fractions, phiObs, phiOut, phiIn); + KnUpdateFlagsObs(flags, fractions, phiObs, phiOut, phiIn, boundaryWidth); } static PyObject *_W_18(PyObject *_self, PyObject *_linargs, PyObject *_kwds) { @@ -1973,8 +1983,9 @@ static PyObject *_W_18(PyObject *_self, PyObject *_linargs, PyObject *_kwds) const MACGrid *fractions = _args.getPtrOpt<MACGrid>("fractions", 2, NULL, &_lock); const Grid<Real> *phiOut = _args.getPtrOpt<Grid<Real>>("phiOut", 3, NULL, &_lock); const Grid<Real> *phiIn = _args.getPtrOpt<Grid<Real>>("phiIn", 4, NULL, &_lock); + int boundaryWidth = _args.getOpt<int>("boundaryWidth", 5, 1, &_lock); _retval = getPyNone(); - setObstacleFlags(flags, phiObs, fractions, phiOut, phiIn); + setObstacleFlags(flags, phiObs, fractions, phiOut, phiIn, boundaryWidth); _args.check(); } pbFinalizePlugin(parent, "setObstacleFlags", !noTiming); diff --git a/extern/mantaflow/preprocessed/plugin/pressure.cpp b/extern/mantaflow/preprocessed/plugin/pressure.cpp index 7def2669e36..780ba44a2b5 100644 --- a/extern/mantaflow/preprocessed/plugin/pressure.cpp +++ b/extern/mantaflow/preprocessed/plugin/pressure.cpp @@ -1171,6 +1171,11 @@ void solvePressureSystem(Grid<Real> &rhs, maxIter = 100; pmg = gMapMG[parent]; + // Release MG from previous step if present (e.g. if previous solve was with MGStatic) + if (pmg && preconditioner == PcMGDynamic) { + releaseMG(parent); + pmg = nullptr; + } if (!pmg) { pmg = new GridMg(pressure.getSize()); gMapMG[parent] = pmg; diff --git a/extern/mantaflow/preprocessed/plugin/secondaryparticles.cpp b/extern/mantaflow/preprocessed/plugin/secondaryparticles.cpp index 281e12ef04b..18582d57e64 100644 --- a/extern/mantaflow/preprocessed/plugin/secondaryparticles.cpp +++ b/extern/mantaflow/preprocessed/plugin/secondaryparticles.cpp @@ -1099,7 +1099,7 @@ void PbRegister_flipSampleSecondaryParticles() // evaluates cubic spline with radius h and distance l in dim dimensions Real cubicSpline(const Real h, const Real l, const int dim) { - const Real h2 = square(h), h3 = h2 * h, h4 = h3 * h, h5 = h4 * h; + const Real h2 = square(h), h3 = h2 * h; const Real c[] = { Real(2e0 / (3e0 * h)), Real(10e0 / (7e0 * M_PI * h2)), Real(1e0 / (M_PI * h3))}; const Real q = l / h; @@ -1175,9 +1175,6 @@ struct knFlipUpdateSecondaryParticlesLinear : public KernelBase { } Vec3i gridpos = toVec3i(pts_sec[idx].pos); - int i = gridpos.x; - int j = gridpos.y; - int k = gridpos.z; // spray particle if (neighborRatio(gridpos) < c_s) { |