Fluid: Added APIC simulation method

Basic support for velocity updates with the APIC method.

This commit adds APIC to the already existing dropdown menu for the simulation method. The APIC plugin within Mantaflow has been updated to the latest version.
This commit is contained in:
Sebastián Barschkis 2020-10-30 00:09:03 +01:00
parent e3e5d595f6
commit 8bdf191461
6 changed files with 279 additions and 154 deletions

View File

@ -6,7 +6,7 @@
// ----------------------------------------------------------------------------
//
// MantaFlow fluid solver framework
// Copyright 2016-2017 Kiwon Um, Nils Thuerey
// Copyright 2016-2020 Kiwon Um, Nils Thuerey
//
// This program is free software, distributed under the terms of the
// Apache License, Version 2.0
@ -21,6 +21,54 @@
namespace Manta {
#define FOR_INT_IJK(num) \
for (int i = 0; i < num; ++i) \
for (int j = 0; j < num; ++j) \
for (int k = 0; k < num; ++k)
static inline IndexInt indexUFace(const Vec3 &pos, const MACGrid &ref)
{
const Vec3i f = toVec3i(pos), c = toVec3i(pos - 0.5);
const IndexInt index = f.x * ref.getStrideX() + c.y * ref.getStrideY() + c.z * ref.getStrideZ();
assertDeb(ref.isInBounds(index),
"Grid index out of bounds for particle position [" << pos.x << ", " << pos.y << ", "
<< pos.z << "]");
return (ref.isInBounds(index)) ? index : -1;
}
static inline IndexInt indexVFace(const Vec3 &pos, const MACGrid &ref)
{
const Vec3i f = toVec3i(pos), c = toVec3i(pos - 0.5);
const IndexInt index = c.x * ref.getStrideX() + f.y * ref.getStrideY() + c.z * ref.getStrideZ();
assertDeb(ref.isInBounds(index),
"Grid index out of bounds for particle position [" << pos.x << ", " << pos.y << ", "
<< pos.z << "]");
return (ref.isInBounds(index)) ? index : -1;
}
static inline IndexInt indexWFace(const Vec3 &pos, const MACGrid &ref)
{
const Vec3i f = toVec3i(pos), c = toVec3i(pos - 0.5);
const IndexInt index = c.x * ref.getStrideX() + c.y * ref.getStrideY() + f.z * ref.getStrideZ();
assertDeb(ref.isInBounds(index),
"Grid index out of bounds for particle position [" << pos.x << ", " << pos.y << ", "
<< pos.z << "]");
return (ref.isInBounds(index)) ? index : -1;
}
static inline IndexInt indexOffset(
const IndexInt gidx, const int i, const int j, const int k, const MACGrid &ref)
{
const IndexInt dX[2] = {0, ref.getStrideX()};
const IndexInt dY[2] = {0, ref.getStrideY()};
const IndexInt dZ[2] = {0, ref.getStrideZ()};
const IndexInt index = gidx + dX[i] + dY[j] + dZ[k];
assertDeb(ref.isInBounds(index),
"Grid index out of bounds for particle position [" << pos.x << ", " << pos.y << ", "
<< pos.z << "]");
return (ref.isInBounds(index)) ? index : -1;
}
struct knApicMapLinearVec3ToMACGrid : public KernelBase {
knApicMapLinearVec3ToMACGrid(const BasicParticleSystem &p,
MACGrid &mg,
@ -30,7 +78,8 @@ struct knApicMapLinearVec3ToMACGrid : public KernelBase {
const ParticleDataImpl<Vec3> &cpy,
const ParticleDataImpl<Vec3> &cpz,
const ParticleDataImpl<int> *ptype,
const int exclude)
const int exclude,
const int boundaryWidth)
: KernelBase(p.size()),
p(p),
mg(mg),
@ -40,7 +89,8 @@ struct knApicMapLinearVec3ToMACGrid : public KernelBase {
cpy(cpy),
cpz(cpz),
ptype(ptype),
exclude(exclude)
exclude(exclude),
boundaryWidth(boundaryWidth)
{
runMessage();
run();
@ -54,73 +104,91 @@ struct knApicMapLinearVec3ToMACGrid : public KernelBase {
const ParticleDataImpl<Vec3> &cpy,
const ParticleDataImpl<Vec3> &cpz,
const ParticleDataImpl<int> *ptype,
const int exclude)
const int exclude,
const int boundaryWidth)
{
if (!p.isActive(idx) || (ptype && ((*ptype)[idx] & exclude)))
return;
const IndexInt dX[2] = {0, vg.getStrideX()};
const IndexInt dY[2] = {0, vg.getStrideY()};
const IndexInt dZ[2] = {0, vg.getStrideZ()};
if (!vg.isInBounds(p.getPos(idx), boundaryWidth)) {
debMsg("Skipping particle at index " << idx
<< ". Is out of bounds and cannot be applied to grid.",
1);
return;
}
const Vec3 &pos = p.getPos(idx), &vel = vp[idx];
const Vec3i f = toVec3i(pos);
const Vec3i c = toVec3i(pos - 0.5);
const Vec3 wf = clamp(pos - toVec3(f), Vec3(0.), Vec3(1.));
const Vec3 wc = clamp(pos - toVec3(c) - 0.5, Vec3(0.), Vec3(1.));
const Vec3 &pos = p[idx].pos, &vel = vp[idx];
const IndexInt fi = static_cast<IndexInt>(pos.x), fj = static_cast<IndexInt>(pos.y),
fk = static_cast<IndexInt>(pos.z);
const IndexInt ci = static_cast<IndexInt>(pos.x - 0.5),
cj = static_cast<IndexInt>(pos.y - 0.5),
ck = static_cast<IndexInt>(pos.z - 0.5);
const Real wfi = clamp(pos.x - fi, Real(0), Real(1)),
wfj = clamp(pos.y - fj, Real(0), Real(1)),
wfk = clamp(pos.z - fk, Real(0), Real(1));
const Real wci = clamp(Real(pos.x - ci - 0.5), Real(0), Real(1)),
wcj = clamp(Real(pos.y - cj - 0.5), Real(0), Real(1)),
wck = clamp(Real(pos.z - ck - 0.5), Real(0), Real(1));
// TODO: check index for safety
{ // u-face
const IndexInt gidx = fi * dX[1] + cj * dY[1] + ck * dZ[1];
const Vec3 gpos(fi, cj + 0.5, ck + 0.5);
const Real wi[2] = {Real(1) - wfi, wfi};
const Real wj[2] = {Real(1) - wcj, wcj};
const Real wk[2] = {Real(1) - wck, wck};
for (int i = 0; i < 2; ++i)
for (int j = 0; j < 2; ++j)
for (int k = 0; k < 2; ++k) {
const Real w = wi[i] * wj[j] * wk[k];
mg[gidx + dX[i] + dY[j] + dZ[k]].x += w;
vg[gidx + dX[i] + dY[j] + dZ[k]].x += w * vel.x;
vg[gidx + dX[i] + dY[j] + dZ[k]].x += w * dot(cpx[idx], gpos + Vec3(i, j, k) - pos);
}
const IndexInt gidx = indexUFace(pos, vg);
if (gidx < 0)
return; // debug will fail
const Vec3 gpos(f.x, c.y + 0.5, c.z + 0.5);
const Real wi[2] = {Real(1) - wf.x, wf.x};
const Real wj[2] = {Real(1) - wc.y, wc.y};
const Real wk[2] = {Real(1) - wc.z, wc.z};
FOR_INT_IJK(2)
{
const Real w = wi[i] * wj[j] * wk[k];
const IndexInt vidx = indexOffset(gidx, i, j, k, vg);
if (vidx < 0)
continue; // debug will fail
mg[vidx].x += w;
vg[vidx].x += w * vel.x;
vg[vidx].x += w * dot(cpx[idx], gpos + Vec3(i, j, k) - pos);
}
}
{ // v-face
const IndexInt gidx = ci * dX[1] + fj * dY[1] + ck * dZ[1];
const Vec3 gpos(ci + 0.5, fj, ck + 0.5);
const Real wi[2] = {Real(1) - wci, wci};
const Real wj[2] = {Real(1) - wfj, wfj};
const Real wk[2] = {Real(1) - wck, wck};
for (int i = 0; i < 2; ++i)
for (int j = 0; j < 2; ++j)
for (int k = 0; k < 2; ++k) {
const Real w = wi[i] * wj[j] * wk[k];
mg[gidx + dX[i] + dY[j] + dZ[k]].y += w;
vg[gidx + dX[i] + dY[j] + dZ[k]].y += w * vel.y;
vg[gidx + dX[i] + dY[j] + dZ[k]].y += w * dot(cpy[idx], gpos + Vec3(i, j, k) - pos);
}
const IndexInt gidx = indexVFace(pos, vg);
if (gidx < 0)
return;
const Vec3 gpos(c.x + 0.5, f.y, c.z + 0.5);
const Real wi[2] = {Real(1) - wc.x, wc.x};
const Real wj[2] = {Real(1) - wf.y, wf.y};
const Real wk[2] = {Real(1) - wc.z, wc.z};
FOR_INT_IJK(2)
{
const Real w = wi[i] * wj[j] * wk[k];
const IndexInt vidx = indexOffset(gidx, i, j, k, vg);
if (vidx < 0)
continue;
mg[vidx].y += w;
vg[vidx].y += w * vel.y;
vg[vidx].y += w * dot(cpy[idx], gpos + Vec3(i, j, k) - pos);
}
}
if (!vg.is3D())
return;
{ // w-face
const IndexInt gidx = ci * dX[1] + cj * dY[1] + fk * dZ[1];
const Vec3 gpos(ci + 0.5, cj + 0.5, fk);
const Real wi[2] = {Real(1) - wci, wci};
const Real wj[2] = {Real(1) - wcj, wcj};
const Real wk[2] = {Real(1) - wfk, wfk};
for (int i = 0; i < 2; ++i)
for (int j = 0; j < 2; ++j)
for (int k = 0; k < 2; ++k) {
const Real w = wi[i] * wj[j] * wk[k];
mg[gidx + dX[i] + dY[j] + dZ[k]].z += w;
vg[gidx + dX[i] + dY[j] + dZ[k]].z += w * vel.z;
vg[gidx + dX[i] + dY[j] + dZ[k]].z += w * dot(cpz[idx], gpos + Vec3(i, j, k) - pos);
}
const IndexInt gidx = indexWFace(pos, vg);
if (gidx < 0)
return;
const Vec3 gpos(c.x + 0.5, c.y + 0.5, f.z);
const Real wi[2] = {Real(1) - wc.x, wc.x};
const Real wj[2] = {Real(1) - wc.y, wc.y};
const Real wk[2] = {Real(1) - wf.z, wf.z};
FOR_INT_IJK(2)
{
const Real w = wi[i] * wj[j] * wk[k];
const IndexInt vidx = indexOffset(gidx, i, j, k, vg);
if (vidx < 0)
continue;
mg[vidx].z += w;
vg[vidx].z += w * vel.z;
vg[vidx].z += w * dot(cpz[idx], gpos + Vec3(i, j, k) - pos);
}
}
}
inline const BasicParticleSystem &getArg0()
@ -168,6 +236,11 @@ struct knApicMapLinearVec3ToMACGrid : public KernelBase {
return exclude;
}
typedef int type8;
inline const int &getArg9()
{
return boundaryWidth;
}
typedef int type9;
void runMessage()
{
debMsg("Executing kernel knApicMapLinearVec3ToMACGrid ", 3);
@ -179,7 +252,7 @@ struct knApicMapLinearVec3ToMACGrid : public KernelBase {
{
const IndexInt _sz = size;
for (IndexInt i = 0; i < _sz; i++)
op(i, p, mg, vg, vp, cpx, cpy, cpz, ptype, exclude);
op(i, p, mg, vg, vp, cpx, cpy, cpz, ptype, exclude, boundaryWidth);
}
const BasicParticleSystem &p;
MACGrid &mg;
@ -190,6 +263,7 @@ struct knApicMapLinearVec3ToMACGrid : public KernelBase {
const ParticleDataImpl<Vec3> &cpz;
const ParticleDataImpl<int> *ptype;
const int exclude;
const int boundaryWidth;
};
void apicMapPartsToMAC(const FlagGrid &flags,
@ -201,23 +275,22 @@ void apicMapPartsToMAC(const FlagGrid &flags,
const ParticleDataImpl<Vec3> &cpz,
MACGrid *mass = NULL,
const ParticleDataImpl<int> *ptype = NULL,
const int exclude = 0)
const int exclude = 0,
const int boundaryWidth = 0)
{
// affine map
// let's assume that the particle mass is constant, 1.0
const bool freeMass = !mass;
if (!mass)
mass = new MACGrid(flags.getParent());
else
mass->clear();
// affine map: let's assume that the particle mass is constant, 1.0
if (!mass) {
MACGrid tmpmass(vel.getParent());
mass = &tmpmass;
}
mass->clear();
vel.clear();
knApicMapLinearVec3ToMACGrid(parts, *mass, vel, partVel, cpx, cpy, cpz, ptype, exclude);
knApicMapLinearVec3ToMACGrid(
parts, *mass, vel, partVel, cpx, cpy, cpz, ptype, exclude, boundaryWidth);
mass->stomp(VECTOR_EPSILON);
vel.safeDivide(*mass);
if (freeMass)
delete mass;
}
static PyObject *_W_0(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
{
@ -241,8 +314,10 @@ static PyObject *_W_0(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
const ParticleDataImpl<int> *ptype = _args.getPtrOpt<ParticleDataImpl<int>>(
"ptype", 8, NULL, &_lock);
const int exclude = _args.getOpt<int>("exclude", 9, 0, &_lock);
const int boundaryWidth = _args.getOpt<int>("boundaryWidth", 10, 0, &_lock);
_retval = getPyNone();
apicMapPartsToMAC(flags, vel, parts, partVel, cpx, cpy, cpz, mass, ptype, exclude);
apicMapPartsToMAC(
flags, vel, parts, partVel, cpx, cpy, cpz, mass, ptype, exclude, boundaryWidth);
_args.check();
}
pbFinalizePlugin(parent, "apicMapPartsToMAC", !noTiming);
@ -270,7 +345,8 @@ struct knApicMapLinearMACGridToVec3 : public KernelBase {
const MACGrid &vg,
const FlagGrid &flags,
const ParticleDataImpl<int> *ptype,
const int exclude)
const int exclude,
const int boundaryWidth)
: KernelBase(vp.size()),
vp(vp),
cpx(cpx),
@ -280,7 +356,8 @@ struct knApicMapLinearMACGridToVec3 : public KernelBase {
vg(vg),
flags(flags),
ptype(ptype),
exclude(exclude)
exclude(exclude),
boundaryWidth(boundaryWidth)
{
runMessage();
run();
@ -294,78 +371,94 @@ struct knApicMapLinearMACGridToVec3 : public KernelBase {
const MACGrid &vg,
const FlagGrid &flags,
const ParticleDataImpl<int> *ptype,
const int exclude) const
const int exclude,
const int boundaryWidth) const
{
if (!p.isActive(idx) || (ptype && ((*ptype)[idx] & exclude)))
return;
if (!vg.isInBounds(p.getPos(idx), boundaryWidth)) {
debMsg("Skipping particle at index " << idx
<< ". Is out of bounds and cannot get value from grid.",
1);
return;
}
vp[idx] = cpx[idx] = cpy[idx] = cpz[idx] = Vec3(Real(0));
const IndexInt dX[2] = {0, vg.getStrideX()}, dY[2] = {0, vg.getStrideY()},
dZ[2] = {0, vg.getStrideZ()};
const Real gw[2] = {-Real(1), Real(1)};
const Vec3 &pos = p[idx].pos;
const IndexInt fi = static_cast<IndexInt>(pos.x), fj = static_cast<IndexInt>(pos.y),
fk = static_cast<IndexInt>(pos.z);
const IndexInt ci = static_cast<IndexInt>(pos.x - 0.5),
cj = static_cast<IndexInt>(pos.y - 0.5),
ck = static_cast<IndexInt>(pos.z - 0.5);
const Real wfi = clamp(pos.x - fi, Real(0), Real(1)),
wfj = clamp(pos.y - fj, Real(0), Real(1)),
wfk = clamp(pos.z - fk, Real(0), Real(1));
const Real wci = clamp(Real(pos.x - ci - 0.5), Real(0), Real(1)),
wcj = clamp(Real(pos.y - cj - 0.5), Real(0), Real(1)),
wck = clamp(Real(pos.z - ck - 0.5), Real(0), Real(1));
// TODO: check index for safety
{ // u
const IndexInt gidx = fi * dX[1] + cj * dY[1] + ck * dZ[1];
const Real wx[2] = {Real(1) - wfi, wfi};
const Real wy[2] = {Real(1) - wcj, wcj};
const Real wz[2] = {Real(1) - wck, wck};
for (int i = 0; i < 2; ++i)
for (int j = 0; j < 2; ++j)
for (int k = 0; k < 2; ++k) {
const IndexInt vidx = gidx + dX[i] + dY[j] + dZ[k];
Real vgx = vg[vidx].x;
vp[idx].x += wx[i] * wy[j] * wz[k] * vgx;
cpx[idx].x += gw[i] * wy[j] * wz[k] * vgx;
cpx[idx].y += wx[i] * gw[j] * wz[k] * vgx;
cpx[idx].z += wx[i] * wy[j] * gw[k] * vgx;
}
const Vec3 &pos = p.getPos(idx);
const Vec3i f = toVec3i(pos);
const Vec3i c = toVec3i(pos - 0.5);
const Vec3 wf = clamp(pos - toVec3(f), Vec3(0.), Vec3(1.));
const Vec3 wc = clamp(pos - toVec3(c) - 0.5, Vec3(0.), Vec3(1.));
{ // u-face
const IndexInt gidx = indexUFace(pos, vg);
if (gidx < 0)
return; // debug will fail
const Real wx[2] = {Real(1) - wf.x, wf.x};
const Real wy[2] = {Real(1) - wc.y, wc.y};
const Real wz[2] = {Real(1) - wc.z, wc.z};
FOR_INT_IJK(2)
{
const IndexInt vidx = indexOffset(gidx, i, j, k, vg);
if (vidx < 0)
continue; // debug will fail
const Real vgx = vg[vidx].x;
vp[idx].x += wx[i] * wy[j] * wz[k] * vgx;
cpx[idx].x += gw[i] * wy[j] * wz[k] * vgx;
cpx[idx].y += wx[i] * gw[j] * wz[k] * vgx;
cpx[idx].z += wx[i] * wy[j] * gw[k] * vgx;
}
}
{ // v
const IndexInt gidx = ci * dX[1] + fj * dY[1] + ck * dZ[1];
const Real wx[2] = {Real(1) - wci, wci};
const Real wy[2] = {Real(1) - wfj, wfj};
const Real wz[2] = {Real(1) - wck, wck};
for (int i = 0; i < 2; ++i)
for (int j = 0; j < 2; ++j)
for (int k = 0; k < 2; ++k) {
const IndexInt vidx = gidx + dX[i] + dY[j] + dZ[k];
Real vgy = vg[vidx].y;
vp[idx].y += wx[i] * wy[j] * wz[k] * vgy;
cpy[idx].x += gw[i] * wy[j] * wz[k] * vgy;
cpy[idx].y += wx[i] * gw[j] * wz[k] * vgy;
cpy[idx].z += wx[i] * wy[j] * gw[k] * vgy;
}
{ // v-face
const IndexInt gidx = indexVFace(pos, vg);
if (gidx < 0)
return;
const Real wx[2] = {Real(1) - wc.x, wc.x};
const Real wy[2] = {Real(1) - wf.y, wf.y};
const Real wz[2] = {Real(1) - wc.z, wc.z};
FOR_INT_IJK(2)
{
const IndexInt vidx = indexOffset(gidx, i, j, k, vg);
if (vidx < 0)
continue;
const Real vgy = vg[vidx].y;
vp[idx].y += wx[i] * wy[j] * wz[k] * vgy;
cpy[idx].x += gw[i] * wy[j] * wz[k] * vgy;
cpy[idx].y += wx[i] * gw[j] * wz[k] * vgy;
cpy[idx].z += wx[i] * wy[j] * gw[k] * vgy;
}
}
if (!vg.is3D())
return;
{ // w
const IndexInt gidx = ci * dX[1] + cj * dY[1] + fk * dZ[1];
const Real wx[2] = {Real(1) - wci, wci};
const Real wy[2] = {Real(1) - wcj, wcj};
const Real wz[2] = {Real(1) - wfk, wfk};
for (int i = 0; i < 2; ++i)
for (int j = 0; j < 2; ++j)
for (int k = 0; k < 2; ++k) {
const IndexInt vidx = gidx + dX[i] + dY[j] + dZ[k];
Real vgz = vg[vidx].z;
vp[idx].z += wx[i] * wy[j] * wz[k] * vgz;
cpz[idx].x += gw[i] * wy[j] * wz[k] * vgz;
cpz[idx].y += wx[i] * gw[j] * wz[k] * vgz;
cpz[idx].z += wx[i] * wy[j] * gw[k] * vgz;
}
{ // w-face
const IndexInt gidx = indexWFace(pos, vg);
if (gidx < 0)
return;
const Real wx[2] = {Real(1) - wc.x, wc.x};
const Real wy[2] = {Real(1) - wc.y, wc.y};
const Real wz[2] = {Real(1) - wf.z, wf.z};
FOR_INT_IJK(2)
{
const IndexInt vidx = indexOffset(gidx, i, j, k, vg);
if (vidx < 0)
continue;
const Real vgz = vg[vidx].z;
vp[idx].z += wx[i] * wy[j] * wz[k] * vgz;
cpz[idx].x += gw[i] * wy[j] * wz[k] * vgz;
cpz[idx].y += wx[i] * gw[j] * wz[k] * vgz;
cpz[idx].z += wx[i] * wy[j] * gw[k] * vgz;
}
}
}
inline ParticleDataImpl<Vec3> &getArg0()
@ -413,6 +506,11 @@ struct knApicMapLinearMACGridToVec3 : public KernelBase {
return exclude;
}
typedef int type8;
inline const int &getArg9()
{
return boundaryWidth;
}
typedef int type9;
void runMessage()
{
debMsg("Executing kernel knApicMapLinearMACGridToVec3 ", 3);
@ -423,7 +521,7 @@ struct knApicMapLinearMACGridToVec3 : public KernelBase {
void operator()(const tbb::blocked_range<IndexInt> &__r) const
{
for (IndexInt idx = __r.begin(); idx != (IndexInt)__r.end(); idx++)
op(idx, vp, cpx, cpy, cpz, p, vg, flags, ptype, exclude);
op(idx, vp, cpx, cpy, cpz, p, vg, flags, ptype, exclude, boundaryWidth);
}
void run()
{
@ -438,6 +536,7 @@ struct knApicMapLinearMACGridToVec3 : public KernelBase {
const FlagGrid &flags;
const ParticleDataImpl<int> *ptype;
const int exclude;
const int boundaryWidth;
};
void apicMapMACGridToParts(ParticleDataImpl<Vec3> &partVel,
@ -448,9 +547,11 @@ void apicMapMACGridToParts(ParticleDataImpl<Vec3> &partVel,
const MACGrid &vel,
const FlagGrid &flags,
const ParticleDataImpl<int> *ptype = NULL,
const int exclude = 0)
const int exclude = 0,
const int boundaryWidth = 0)
{
knApicMapLinearMACGridToVec3(partVel, cpx, cpy, cpz, parts, vel, flags, ptype, exclude);
knApicMapLinearMACGridToVec3(
partVel, cpx, cpy, cpz, parts, vel, flags, ptype, exclude, boundaryWidth);
}
static PyObject *_W_1(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
{
@ -473,8 +574,10 @@ static PyObject *_W_1(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
const ParticleDataImpl<int> *ptype = _args.getPtrOpt<ParticleDataImpl<int>>(
"ptype", 7, NULL, &_lock);
const int exclude = _args.getOpt<int>("exclude", 8, 0, &_lock);
const int boundaryWidth = _args.getOpt<int>("boundaryWidth", 9, 0, &_lock);
_retval = getPyNone();
apicMapMACGridToParts(partVel, cpx, cpy, cpz, parts, vel, flags, ptype, exclude);
apicMapMACGridToParts(
partVel, cpx, cpy, cpz, parts, vel, flags, ptype, exclude, boundaryWidth);
_args.check();
}
pbFinalizePlugin(parent, "apicMapMACGridToParts", !noTiming);

View File

@ -697,12 +697,6 @@ void MANTA::initializeRNAMap(FluidModifierData *fmd)
if ((fds->border_collisions & FLUID_DOMAIN_BORDER_TOP) == 0)
borderCollisions += "Z";
string simulationMethod = "";
if (fds->simulation_method & FLUID_DOMAIN_METHOD_FLIP)
simulationMethod += "'FLIP'";
else if (fds->simulation_method & FLUID_DOMAIN_METHOD_APIC)
simulationMethod += "'APIC'";
string particleTypesStr = "";
if (fds->particle_type & FLUID_DOMAIN_PARTICLE_SPRAY)
particleTypesStr += "PtypeSpray";
@ -837,7 +831,7 @@ void MANTA::initializeRNAMap(FluidModifierData *fmd)
mRNAMap["CACHE_MESH_FORMAT"] = getCacheFileEnding(fds->cache_mesh_format);
mRNAMap["CACHE_NOISE_FORMAT"] = getCacheFileEnding(fds->cache_noise_format);
mRNAMap["CACHE_PARTICLE_FORMAT"] = getCacheFileEnding(fds->cache_particle_format);
mRNAMap["SIMULATION_METHOD"] = simulationMethod;
mRNAMap["USING_APIC"] = getBooleanString(fds->simulation_method == FLUID_DOMAIN_METHOD_APIC);
mRNAMap["FLIP_RATIO"] = to_string(fds->flip_ratio);
mRNAMap["PARTICLE_RANDOMNESS"] = to_string(fds->particle_randomness);
mRNAMap["PARTICLE_NUMBER"] = to_string(fds->particle_number);

View File

@ -50,7 +50,8 @@ smoothenPos_s$ID$ = $MESH_SMOOTHEN_POS$\n\
smoothenNeg_s$ID$ = $MESH_SMOOTHEN_NEG$\n\
randomness_s$ID$ = $PARTICLE_RANDOMNESS$\n\
surfaceTension_s$ID$ = $LIQUID_SURFACE_TENSION$\n\
maxSysParticles_s$ID$ = $PP_PARTICLE_MAXIMUM$\n";
maxSysParticles_s$ID$ = $PP_PARTICLE_MAXIMUM$\n\
using_apic_s$ID$ = $USING_APIC$\n";
const std::string liquid_variables_particles =
"\n\
@ -91,6 +92,14 @@ curvature_s$ID$ = None\n\
pp_s$ID$ = s$ID$.create(BasicParticleSystem, name='$NAME_PARTS$')\n\
pVel_pp$ID$ = pp_s$ID$.create(PdataVec3, name='$NAME_PARTSVELOCITY$')\n\
\n\
pCx_pp$ID$ = None\n\
pCy_pp$ID$ = None\n\
pCz_pp$ID$ = None\n\
if using_apic_s$ID$:\n\
pCx_pp$ID$ = pp_s$ID$.create(PdataVec3)\n\
pCy_pp$ID$ = pp_s$ID$.create(PdataVec3)\n\
pCz_pp$ID$ = pp_s$ID$.create(PdataVec3)\n\
\n\
# Acceleration data for particle nbs\n\
pindex_s$ID$ = s$ID$.create(ParticleIndexSystem, name='$NAME_PINDEX$')\n\
gpi_s$ID$ = s$ID$.create(IntGrid, name='$NAME_GPI$')\n\
@ -275,8 +284,12 @@ def liquid_step_$ID$():\n\
resetOutflow(flags=flags_s$ID$, phi=phi_s$ID$, parts=pp_s$ID$, index=gpi_s$ID$, indexSys=pindex_s$ID$)\n\
flags_s$ID$.updateFromLevelset(phi_s$ID$)\n\
\n\
# combine particles velocities with advected grid velocities\n\
mapPartsToMAC(vel=velParts_s$ID$, flags=flags_s$ID$, velOld=velOld_s$ID$, parts=pp_s$ID$, partVel=pVel_pp$ID$, weight=mapWeights_s$ID$)\n\
# combine particle velocities with advected grid velocities\n\
if using_apic_s$ID$:\n\
apicMapPartsToMAC(flags=flags_s$ID$, vel=vel_s$ID$, parts=pp_s$ID$, partVel=pVel_pp$ID$, cpx=pCx_pp$ID$, cpy=pCy_pp$ID$, cpz=pCz_pp$ID$)\n\
else:\n\
mapPartsToMAC(vel=velParts_s$ID$, flags=flags_s$ID$, velOld=velOld_s$ID$, parts=pp_s$ID$, partVel=pVel_pp$ID$, weight=mapWeights_s$ID$)\n\
\n\
extrapolateMACFromWeight(vel=velParts_s$ID$, distance=2, weight=mapWeights_s$ID$)\n\
combineGridVel(vel=velParts_s$ID$, weight=mapWeights_s$ID$, combineVel=vel_s$ID$, phi=phi_s$ID$, narrowBand=combineBandWidth_s$ID$, thresh=0)\n\
velOld_s$ID$.copyFrom(vel_s$ID$)\n\
@ -319,7 +332,11 @@ def liquid_step_$ID$():\n\
# set source grids for resampling, used in adjustNumber!\n\
pVel_pp$ID$.setSource(grid=vel_s$ID$, isMAC=True)\n\
adjustNumber(parts=pp_s$ID$, vel=vel_s$ID$, flags=flags_s$ID$, minParticles=minParticles_s$ID$, maxParticles=maxParticles_s$ID$, phi=phi_s$ID$, exclude=phiObs_s$ID$, radiusFactor=radiusFactor_s$ID$, narrowBand=adjustedNarrowBandWidth_s$ID$)\n\
flipVelocityUpdate(vel=vel_s$ID$, velOld=velOld_s$ID$, flags=flags_s$ID$, parts=pp_s$ID$, partVel=pVel_pp$ID$, flipRatio=flipRatio_s$ID$)\n";
\n\
if using_apic_s$ID$:\n\
apicMapMACGridToParts(partVel=pVel_pp$ID$, cpx=pCx_pp$ID$, cpy=pCy_pp$ID$, cpz=pCz_pp$ID$, parts=pp_s$ID$, vel=vel_s$ID$, flags=flags_s$ID$)\n\
else:\n\
flipVelocityUpdate(vel=vel_s$ID$, velOld=velOld_s$ID$, flags=flags_s$ID$, parts=pp_s$ID$, partVel=pVel_pp$ID$, flipRatio=flipRatio_s$ID$)\n";
const std::string liquid_step_mesh =
"\n\

View File

@ -486,7 +486,8 @@ class PHYSICS_PT_liquid(PhysicButtonsPanel, Panel):
col = flow.column()
col.prop(domain, "simulation_method", expand=False)
col.prop(domain, "flip_ratio", text="FLIP Ratio")
if domain.simulation_method == 'FLIP':
col.prop(domain, "flip_ratio", text="FLIP Ratio")
col.prop(domain, "sys_particle_maximum", text="System Maximum")
col = col.column(align=True)
col.prop(domain, "particle_radius", text="Particle Radius")

View File

@ -4993,6 +4993,7 @@ void BKE_fluid_modifier_copy(const struct FluidModifierData *fmd,
tfds->fractions_threshold = fds->fractions_threshold;
tfds->fractions_distance = fds->fractions_distance;
tfds->sys_particle_maximum = fds->sys_particle_maximum;
tfds->simulation_method = fds->simulation_method;
/* diffusion options*/
tfds->surface_tension = fds->surface_tension;

View File

@ -1456,8 +1456,16 @@ static void rna_def_fluid_domain_settings(BlenderRNA *brna)
{0, NULL, 0, NULL, NULL}};
static EnumPropertyItem simulation_methods[] = {
{FLUID_DOMAIN_METHOD_FLIP, "FLIP", 0, "FLIP", "Use FLIP as the simulation method"},
/*{FLUID_DOMAIN_METHOD_APIC, "APIC", 0, "APIC", "Use APIC as the simulation method"},*/
{FLUID_DOMAIN_METHOD_FLIP,
"FLIP",
0,
"FLIP",
"Use FLIP as the simulation method (more splashy behavior)"},
{FLUID_DOMAIN_METHOD_APIC,
"APIC",
0,
"APIC",
"Use APIC as the simulation method (more energetic and stable behavior)"},
{0, NULL, 0, NULL, NULL},
};
@ -1820,6 +1828,7 @@ static void rna_def_fluid_domain_settings(BlenderRNA *brna)
RNA_def_property_enum_sdna(prop, NULL, "simulation_method");
RNA_def_property_enum_items(prop, simulation_methods);
RNA_def_property_ui_text(prop, "Simulation Method", "Change the underlying simulation method");
RNA_def_property_clear_flag(prop, PROP_ANIMATABLE);
RNA_def_property_update(prop, NC_OBJECT | ND_MODIFIER, "rna_Fluid_domain_data_reset");
prop = RNA_def_property(srna, "flip_ratio", PROP_FLOAT, PROP_NONE);