Liquid-rigid coupling

The next version of SpaceSim will have an improved handling of interactions between various objects in the simulation. How does it work?

Collision handling

One of the additions to this version is a solid boundary (previously only a periodic boundary was available). For this reason, the simulation has to handle collisions between:

• Two particles (handled by the IISPH solver)
• A rigid (non-deformable) object and a particle
• Two rigid objects
• A particle and a boundary
• A rigid object and a boundary

Importantly, there has to be a two-way coupling between the simulation objects that satisfies the principle of action and reaction. Rigid objects can exert a force on fluid particles and push them away, but they also experience a force and torque from the particles that push them in the opposite direction and change their rotation state.

There are different methods of dealing with object collisions. SpaceSim uses the Soft-sphere model, following this paper. Interactions are implemented as additional forces between objects, the magnitude of which depends on their overlap - the more they overlap, the larger the force that pushes them away from each other. Collisions happen over several time steps, as the added acceleration gradually changes the object's velocities.

In contrast, the Hard-sphere model does not allow any particle overlap and changes particle velocities instantaneously in a single time step. Both models have their advantages and disadvantages. I chose the soft-sphere model as it naturally conserves the total momentum and angular momentum of the system.

Arbitrary shapes

If all objects in the simulation were spheres, the collision handling would be straightforward. The tricky part is to make it work with objects of any shape, especially since the object can be defined by a custom triangle mesh that is not known in advance.

First, we need a unified way to describe an arbitrary shape in a way suitable for collision detection. For this, I used the signed-distance fields (SDF). An SDF of an object is a function that, given a point in space, returns the closest distance to the surface of the object, positive if the point lies outside of the object and negative if it lies inside. This is perfect for collisions - it provides a simple test to see if a particle lies outside or inside the object, but we also know the exact distance to the surface.

For most of the predefined shapes, the SDF is actually a simple analytical function. For example, a sphere of radius R can be described by the function SDF(p) = |p| - R, where |p| is the length of the vector. There exist similar (albeit more complex) functions for a lot of shapes - cube, cylinder, torus, etc. A large list of SDFs can be found on this website.

Naturally, there is no simple function describing an SDF for a generic triangular mesh. However, we can construct an approximate function by placing the mesh into a three-dimensional grid and calculating the distance to the mesh for each cell of the grid. If the cell size is small compared to the mesh, this approach will still give us an accurate result. Calculating the SDF can be quite slow, but it is done only once at the beginning of the simulation, so the actual collision handling is nearly as fast as using an analytical function.

The SDF is calculated in two steps:

Finding ray intersections

For each grid cell, the algorithm shoots several rays in different directions, starting in the center of the grid cell, and finds the first intersection of the ray with the mesh. The intersection is used to determine whether the cell lies 'inside' or 'outside' of the mesh - when the ray hits the back side of a triangle, it adds an 'inside' vote, if it hits the front side of the triangle or there is no intersection with the mesh, it adds an 'outside' vote. Finally, it marks the cell as 'inside' if there are more inside votes than 'outside' votes.

Furthermore, if the intersection is within a certain distance to the cell center, the algorithm stores the distance to the grid cell. This gives us a narrow band of cells with known distances to the mesh, with the rest of the grid having unknown distances. Although we could store the distances along the ray for all cells, this distance estimate is unreliable and sensitive to the choice of ray direction, and generally creates discontinuous SDFs, which is not ideal for collision handling. There is a better way to do this, explained below.

Solving the Eikonal equation

The second step is to extend the distance values from the input narrow band of cells close to the mesh where the distances are already known. This mathematical problem can be formulated as the Eikonal equation:

|∇u(x)| = f(x)

where the known narrow-band values are the boundary conditions of the equation. The equation, discretized on a grid using finite differences, can be efficiently solved using either the fast marching method or the fast sweeping method. Once the solution is calculated, we have the final SDF of the triangle mesh, which can be used for collision detection.

An example of the calculated SDF is shown below.

Input triangle mesh from which an SDF is calculated.

Animation showing layers of the signed distance field created from the boat mesh. Green and red colors mean positive and negative distances, respectively. Each animation frame corresponds to one slice through the volume at a fixed depth.

Collision of two rigid objects

Having the SDF of a rigid object, we can easily detect a collision with a spherical particle and calculate the overlap, used to derive the magnitude of the force. A collision of two rigid objects, however, is a bit more involved, the overlap could be calculated only for specific shapes, such as a sphere or a boundary plane. We need a general solution that works for all shapes.

To make this work in SpaceSim, I used a different representation for each of the colliding objects. One object is represented as the SDF, but the other object is simply a point cloud (a group of 3D points) randomly distributed on its surface. It is not feasible to use the mesh vertices directly, it would be too inaccurate and could even lead to objects passing through each other. Instead, the surface is sampled by 5000 points with uniform probability, as shown in the images below, ensuring it's accurately represented even for meshes containing triangles of various sizes.

Original mesh with many triangles of different sizes.

Point cloud of uniformly sampled surface points

The rest of the algorithm is the same as for particles. It goes over each point of the point cloud and calculates the overlap with the other rigid object, computing the necessary forces and torques. Even though 5000 points sounds like a lot of CPU work, the algorithm has a linear complexity (O(N), N being the number of points), so it can be run efficiently even for a large number of rigid objects.

A group of 9 rigid objects interacting gravitationally and colliding with each other.

The last part to solve is the collision with the boundary. Luckily, we already have everything we need - we can simply consider the boundary as a special rigid object and handle it the exact same way. There are only two things that make the boundary different from a rigid box. First, it is "inside out"; we have to prevent the particles from going outside of the box, not inside. Second, the boundary does not move, so it shouldn't change its position or rotation. Otherwise, it's the same.

That's all for now. Stay tuned for more updates!