Any Plans for Liquids?

[Leadwerks]

What is so special about your algorithm, other than the fact you put everything into a 3D grid to reduce the number of particle combinations you have iterate through? Is there any other major optimization that can be done?

[JoeJ]

What is so special about your algorithm, other than the fact you put everything into a 3D grid to reduce the number of particle combinations you have iterate through? Is there any other major optimization that can be done?

As said, the only difference to state of the art methods from papers is that i trade a need for atomics or spin locks against a coarse division of work so it can be parallelized lock free *). But that's nothing special and does not explain why i seemingly beat mGPU paper results.

My 3D grid is certainly not special - instead i only adopted this idea from said papers. But it's not about particle combinations at all, as there are no direct particle interactions in MPM (unlike SPH).
Traditionally, MPM is simply like so:
'G2P': Particles scatter their velocity and mass to a 3^2 neighborhood of grid cells.
'P2G': Particles gather updated velocity from the grid in the same way.
(Latest paper proposes to fuse both loops into one, halfing particle BW and also halfing the need for expensive SVD matrix decompositions. Very promising, but i have no tried this yet.)

So you see the grid acts to smooth out velocities which we want for fluid.
To make this realistic, they use 3x3 matrices per particle to model things like deformation of space around them, and those matrices 'skew' the particle velocity as injected to adjacent grid cells, allowing to model elastic collisions, or any kind of shrinking / expanding of the surroundings, e.g. to preserve volume. This way a particle can also track state over time, e.g. 'I'm compressed to 140 % of density, and i want to expand to get back to 100 %.

Beside that, a performance relation of the grid (or sparse grid, etc.) has motivations like this:
Make groups for parallel processing and coherent memory access.
Limit memory hazards (particles may write to the same cells concurrently)
Allow larger domains.



*) I just came here to point out that this is no good idea for games either.
I have disabled MT for debugging and noticed on a single core the 10k particle fluid experiment goes from 6ms to only 10ms. So, for a low particle count my method performs poorly, even if i would use a simple regular grid. The workloads become too small, and threads spend too much time on waiting on the last one.

This raises the question if particle binning and reordering in memory to grid is worth it at all. Without that, particles memory access to the momentum grid becomes random over time, but i can imagine that's the lesser evil.
If we adopt the fused 'P2G2P' idea from the paper, we will access each particle only once during simulation. No need to reorder stuff just for that, i guess.

[JoshKlint]

Okay, so if I understand correctly, it works like this:

1. Iterate through all particles and add their properties into a grid. Each grid cell gets a velocity property that is the average of all the particles that contribute to that grid cell.

2. Iterate through all particles and apply the properties of the grid cell(s) they occupy. If the grid cell has a velocity property, that gets mixed with the particle's current velocity, with the idea that other particles in that cell are pushing it in that direction.

This makes for a course estimate of particle interactions, so your loop is n * 2 instead of n * n. If you just use a thread ID and an array for each property you can easily eliminate the need for any mutex locks.

[FSA]

Is this feature currently in development, or are other thing of higher priority in the making?

[Julio Jerez]

I need to get the opencl integration, so that at least we have one way to get a decent usable feature.

Using cpus, can probably yield 8 to 10 thousand particle simulation, but it will tax the system with a mindless useful demo with not practical application, and I am tire of that.

[JoshKlint]

I'm sure what you come up with will be great.

[JoeJ]

2. Iterate through all particles and apply the properties of the grid cell(s) they occupy. If the grid cell has a velocity property, that gets mixed with the particle's current velocity

Slight mistake here: There is no mixing - particles get their entire velocity from sampling the grid. Standard is a 3^3 cubic filter (did not try how bilinear would do). Weights come from particle mass.

This gives some limitations: While particles can have different materials, e.g. water and sand, behavior blends a lot due to shared velocity grid, so this needs very high resolution to support true multiple materials.
Also particles are not allowed to move faster than one cell per step.

[JoeJ]

BTW, why not some 2D shallow water sim? Much more practical for games?
Random video:

[JoshKlint]

I would certainly be happy with a limited capability like that. You could enable and disable areas of water, or run the simulation in not-quite-realtime and record the water flow, then play it back as a looping animation.

[JoeJ]

I have implemented shallow water sim along working on erosion of 2d terrains. Like it a lot. For example i see basins of water forming on different heights of the terrain, with waterfalls connecting them. Swimming around there really would be fun.
I guess games use it a lot already, e.g. when there are water ripples in puddles. But could be used for global water at lower res. as well. Global water levels usually remain constant, so simulation could be reduced to area around player.
Here's a nice WebGL demo running it on GPU. Press 'c' to drop water: https://lanlou123.github.io/Webgl-Erosion/

[JoeJ]

I have now added the option to do SPH to my simulator, and i see i was wrong with assumption MPM would be much faster.

To find close particles, i did those changes:
For MPM, particles are binned to sparse grid blocks, each block has 4^3 volume cells.
I have now refined this so particles are binned to each individual cell, so i can iterate particles per cell, to reduce the number of potential neighbors.

Actually i use it only to resolve collisions. It's not a full SPH implementation modeling things like viscosity. But we can see the cost of traversing neighbors.

I have 216k particles, and initially there are about 8 particles in each cell on average. MPM is still running as well, and i get this timings:

P2G and G2P (MPM): 37ms
SPH collisions with a search radius of 0.5 times the cell size: 20ms
Sparse Grid update: 37ms (forgot to disable debug checks)
Extra binning per cell for SPH: 4ms

So that's pretty good. But each particle interacts only with 8 neighbors on average this way. Likely not enough for fluid. If i increase the collision radius to 1 cell size, average interaction becomes something like 64 neighbors, and timing increases from 20 to 40ms.

I still think MPM is faster, but not so much.
And if SPH allows larger timesteps (can't say), it would be the better option from perf. perspective.

[Julio Jerez]

I did not think that this MPM method was the breakthrough you'd mentioned.
I read some of the papers and after seen all the limitations, I soon lost the interest,
Stuff like not really using a kernel like sph, where you can make it as large as you want trading accuracy and correctness for speed. Or the physic limits, like impossibility of stuff like turbulence, surface tension, and so on.
Those things convenced me that mpm is more like a game thing and not a real physic based approach.
SPH in the limit converges to a full Navier Stock model if you make your kernel wide enought and you can get a very wide range of physic vehavior.

I have stop the particle at the moment in order to advance the core features of 4.00 and bring it to be in part with 3.14

But after I get that, I will resume it.

Now I am completing the vehicle model, which is comming out quite good. I am very happy with the result so far.

After that I will give it another shot at the walking character.
Newto 4.0 fixes some of the core intrinsic problems I had in 3.14 that I had to spend most my time finding work around.
4.0 is more geared to physic modeling with most emphasis in making it user friendly.

The after that I will resume the fluid with SPH and add open cl support.

[JoeJ]

I did not think that this MPM method was the breakthrough you'd mentioned.
I read some of the papers and after seen all the limitations, I soon lost the interest,
Stuff like not really using a kernel like sph, where you can make it as large as you want trading accuracy and correctness for speed. Or the physic limits, like impossibility of stuff like turbulence, surface tension, and so on.

Yeah, you're right. Those are the reasons i now added SPH support.
I was able to do surface tension with MPM for example, but i need high resolutions then, because 'kernel size' is fixed to the grid size with MPM. It's also hard to process at variable frequencies because of that.

Likely another SPH advantage also is it needs less particle memory. Traditional MPM needs one 3x3 matrix, and MLS-MPM needs two of them. So my particles are fat and that's why binning/reordering takes 4ms due to bandwidth.

[JoeJ]

Ok, just started to work on SPH for real.
First impression is MPM seems overall better, both regarding realism and performance, but just worked on it one day. I want it to model sediment of landslides, and for that it seems nice because particles act more individually and less like a blob.

However, i came up with an interesting optimization, which would work for you too.
Actually i have two passes: 1. calculate density, 2. viscosity force and integration.
So we need to iterate adjacent articles two times per step.

I just tried to use density from the previous step, to iterate just once.
Simulating a cube of particles which flows over the floor for 2000 steps. And i can't see any difference in the output, while i get a nice speedup of two.

I saw you cache the neighbour list. But that's a lot of memory, bandwidth, and errors if count exceeds 32.

What do you think about such trick?

[Julio Jerez]

difficult to say, since I have not use MPM.
and yes I precompute the hash in an array yes. I set it to 32 max.
for a grid size that is twice as large as the radius of the particle it will never generate more than 32 pairs. in fact some of the demos you see out there cut the max to 16.
but you are correct, if is ever produce more than 32 pair it will fail.
I see that some people us hash map, but I am going for the simple approach, and I am counting on the large memory bandwidth of GPUs. Hash maps algorithm are no very good for gpu, but that is Indeed an option.

but I am not going to leave it be like that.
In fact, when I come back to it, the number or pair will have not upper bound.
this is why I am spending so much time on a very efficient prefix scan and sorting algorithm.
we an efficient scan, you get an array of start and number of elements, so not need to max limit.
this will be important for introducing rigid body solids or different density fluids.

the way I am to add rigid solid is, as follows.
-voxelize the solid with a flood field algorithm.
-each voxel will be a particle, with the rigid body ID
-add the clump of voxel to the mix.
-let the algorisms run, when the acceleration on each particle is calculated.
-is a particle is solid it invmass is zero so it just bounce of the solid.
-is tow particle are solid, they accel is zero.
-after the fluids are all calculated, them for each solid it iterated summing all the force and torque contribution of each particle. and add that to the body.

this will apply to all boundaries as well.
and as you can see is may be that the voxelization requires a difference particle radius.