 changed status to resolved
Hello, I'm working on deep learning applied to fluid mechanics following the approach of Tompson et al.  Accelerating Eulerian Fluid Simulation With Convolutional Networks. They generated their dataset using Mantaflow using a timestep of 0.1 for their solver and I have one question about the implementation of your code. In your plume_2d example you choose your timestep = 1, but what if I modify this value (e.g dt = 0.1). I've seen that this time step is naturally taken into account in the SemiLagrangian advection but NOT in the resolution of the Poisson equation (e.g. coefficient in the diagonal are 1 whereas they should be dt/ (rhodxdx) and knowing that dx=1 and rho=1 in all your code, the dt should still remain there.), NEITHER in the velocity correction (source/plugin/pressure.cpp, CorrectVelocity): there, you got
//! Kernel: Apply velocity update from poisson equation KERNEL(bnd=1) void CorrectVelocity(FlagGrid& flags, MACGrid& vel, Grid<Real>& pressure) { int idx = flags.index(i,j,k); if (flags.isEmpty(idx)) { // This was originally a debug print, but now I'll leave it in here as a // way to make sure we're never marking flags as empty (otherwise our own // velocity update code will break). (i.e. We DO NOT handle the empty // cell case). printf("isEmpty (i, j, k) = (%d, %d, %d)  NOT EXPECTED!\n", i, j, k); } if (flags.isFluid(idx)) { if (flags.isFluid(i1,j,k)) vel[idx].x = (pressure[idx]  pressure(i1,j,k)); if (flags.isFluid(i,j1,k)) vel[idx].y = (pressure[idx]  pressure(i,j1,k)); if (flags.is3D() && flags.isFluid(i,j,k1)) vel[idx].z = (pressure[idx]  pressure(i,j,k1)); if (flags.isEmpty(i1,j,k)) vel[idx].x = pressure[idx]; if (flags.isEmpty(i,j1,k)) vel[idx].y = pressure[idx]; if (flags.is3D() && flags.isEmpty(i,j,k1)) vel[idx].z = pressure[idx]; } else if (flags.isEmpty(idx)&&!flags.isOutflow(idx)) // don't change velocities in outflow cells { if (flags.isFluid(i1,j,k)) vel[idx].x += pressure(i1,j,k); else vel[idx].x = 0.f; if (flags.isFluid(i,j1,k)) vel[idx].y += pressure(i,j1,k); else vel[idx].y = 0.f; if (flags.is3D() ) { if (flags.isFluid(i,j,k1)) vel[idx].z += pressure(i,j,k1); else vel[idx].z = 0.f; } } }
and I would have expected a dt multiplying the discrete pressure gradient (p(i,j,k)  p(i1,j,k)). Could you confirm that you have make the assumption of dt=1? This timestep is used in the advection and buoyancy/gravity algorithms but not in pressure projection/velocity correction.
Thank you very much for your time! Antonio
Comments (1)

reporter  Log in to comment
I'm going to leave the response that I kindly got from Nils Thuerey: no time step multiplication is needed in the projection step as if you divide the rhs of the system by dt, you also multiply the gradient of pressure by the same amount, effectively cancelling both contributions. Therefore, it is better to leave it as in Mantaflow in order to obtain a well conditionned Laplacian matrix, Thank you very much for your work!!