This blog is going to talk about a section of my Master’s Thesis. I applied the Parker-Sochacki Method to a section of the FLIP solver in Houdini and will discuss the results.
I created a simulation of a Rigid Body sphere falling into a fluid tank. The fluid is a FLIP simulation. I’ll start with showing what I mean:
Standard FLIP Method:
Notice that the splash does not go straight up and down. It moves to the side. This is because, under the standard settings, it is not very accurate.
Part of my Master’s Thesis was applying the Parker-Sochacki method to part of the FLIP method, a hybrid of Eulerian and Lagrangian based fluids, to see if I could speed up a part of the simulation.
I will be referencing Houdini 13 so some of this might be outdated since Houdini 15 is now available. The FLIP solver inside Houdini is split into 5 sections:
I modified section 5
First, I’m going to show a slight modification I did that made the simulation more accurate. This is what the Particle Motion section looks like in Houdini 13:
The Advect Particles node, second from the left, is a microsolver that advects fields and geometry by a velocity field. I removed it and replaced it with the following:
This was replaced with the Euler step. The Euler step was simply taking the velocity, multiplying it by a time step, and adding it to the current position. That is the first video you saw in this blog. The splash moved sideways showing that it wasn’t accurate. In this version, we multiplied velocity by the time step divided by 2 and added that to the current position. Then that was multiplied by the time step and added to the current position for the position update. Look at the simulation and see that this splash is more straight up and down than the previous video. I got the idea from reading Professor Egbert’s Fluid Flow for the Rest of Us paper.
Slightly Modified FLIP Method:
I’m showing that accuracy matters.
Here is the velocity equation that comes from FLIP:
The velocity is calculated by taking the current velocity and subtracting by the time step divided by density times grid size h, multiplied by the pressure gradient. Density and grid size are both constants that can be adjusted in Houdini. Professor Seth Holladay previously wrote a microsolver to compute the Pressure Gradient in Houdini and allowed me to use his microsolver.
To do a Parker-Sochacki type implementation, I let the right hand side(without the time step) be acceleration:
Acceleration is the same as v'(t). Density and grid size are constants so I let the pressure gradient be an auxiliary variable. We converted this into a tenth order method. The pressure gradient at order m was created, then that was used to find the velocity at order m+1 and then that was used to find the position at order m+1:
The final result is a Maclaurin Series that updates Position and Velocity:
To help out, here is an example of how I find the m and m+1 terms at order 7 and 8:
I would like to note that I normally use “i” and “i+1” rather than “m” and “m-1” for Parker-Sochacki. You can use whatever letters you like.
This simulation is different from the videos I have shown(the ball was higher) but the top image was a default setting and the bottom was the Parker-Sochacki implementation:
Here is a video where the ball was dropped from a lower height. In the videos, the ball was dropped from a height of 7. In the following simulation, the ball was dropped from a height of 5. In the images it was dropped from a heighr of 10.
Parker-Sochacki FLIP Method:
I need to note that while I did find a new way to update the velocity and position for a FLIP simulation, it is not an effective one compared to what is available. It was significantly slower.
First, this image might be hard to see but look how much larger the node network is in section 5 compared to the what the original setting is:
Most of the nodes on the left are using the gasmatchfield node in Houdini. I also have a special subnetwork for the Parker-Sochacki method:
Here is the problem with this implementation. The Performance monitor in Houdini measures how much time it takes a specific node. Having more nodes in a simulation does not mean that the simulation will take longer. However, the slowest part of a simulation is red. Take a look at what is red:
Under the default settings, the popsolver is red, the gasparticletofield_vel is red, and the advenct_vel_by_oldvel is red. So this a problem. Here is why Parker-Sochacki is slow for this implementation:
This was hard coded using Geometry VOPs. I have to go back and forth between the Eulerian and Lagrangian parts of FLIP just to compute the next iteration m. This takes too much time. I solve for the pressure gradient and then a “special auxiliary variable pressure gradient” over and over again recomputing the same terms. The entire simulation became significantly slower.