Unstable hydroelastic cube contact

[amcastro-tri]

What happened?

This was first reported internally by @AdityaBhat-TRI when manipulating a cube in simulation.
Here I provide a minimum reproduction with a cube on top of the ground. The cube has 10 cm sides and it has a 1 kg mass.

We observe rattling with either zero margin or the recommended 1e-4 m margin.
Increasing margin to a large value of 1 cm should remove any rattling. Rattling is removed, however we observe the box leaning towards one of its corners. This is most unexpected.

Finally, @joemasterjohn observed that if we perturbe one of the lenghts of the cubes by something as small as 1e-4m, the problem goes away. Also most unexpected. The problem seems to be for perfect or close to perfect cubes?

Version

1.35.0

[jbinagia]

Just to add some additional data points, I also am observing this (I believe same) issue, reproduced in deepnote here. Some notes:

• In my case, I also observe the rattling for non-cubes (general cuboids)
• I do not observe the rattling anymore if I switch discrete_contact_approximation to sap (as opposed to lagged or similar)

[amcastro-tri]

Below is an animation of the minimum reproduction case I provided above:

And an animation of the contact surface:

What happens here is that our custom hydro box meshes are as coarse as possible while still faitfully reproducting the pressure field. However, this causes trouble in discrete mode. As you can see, we get an oscillatory motion along an axis pependicular to the line along which we (arbitrarly) decide to split the box's face (while I don't show this directly, it's easy to infer from the contact surface animation).
The oscillations happen because we apply a single contact force per polygon (in this case mostly the two big triangles) at their geometric centroid. In equilibrium, the centroid would be right in the middle. Small deviations from the equilibrium cause the center of pressure to move back and forth. Since this is treated explicitly (i.e. geometry is "frozen" during the SAP resolution phase), the two only contact points oscillate back and forth.

TLDR; this effectively behaves as if we placed the box balancing on top of a blade along one of its diagonals. Obviously pretty unstable.

[amcastro-tri]

The next thing I tried, is to implent the Center of Pressure (COP) correction suggested separately in #22440.
I developed the analytical equations in dyadic form and implemented them in this branch.

Below is the meshcat animation:

And the contact surface animation:

The problem persists, basically for the same reason stated above: once the contact point is computed (whether centroid or COP), it is "frozen" (explicit) during the SAP (implicit) contact resolution phase.

Now, we do see this goes in the right direction, with the magnitude of the oscillations actually smaller (we see that the patches stay triangular most of the time, i.e. they do not escape the ground box as before) and I also observed the frequency of the oscillations to go higher (since the COP reacts more quickly to angular changes than the centroid).

Still, the oscillations persists given the explicit treatment of the position of the contact point.

TLDR; This does not solver our problem. I strongly discourage an implicit treatment of the COP or the torques produced by each patch. While this should solve the problem, it'd lead to a significatn the complexity of the code, increase runtime cost (more expensive constraints), and at little if no gain at all for the general hydro case.

[amcastro-tri]

The last thing I tried is a custom box mesh for discrete contact, @xuchenhan-tri suggested this idea.

In the same way we do have a custom box mesh today, it is natural to have a custom box mesh that solves this particular problem. The simplest thing we can try, is a mesh that splits each face of the box into four triangles. I implemnted this idea in this branch (along with CoP calculation, though obviously located in completely separate files in the codebase).

Here is a meshcat animation:

And the contact surface animation:

NOTE: Look carefully, those animations are very stable and only show a quick transient at the beginning, so do not think they are still images.

From the contact surface animation, you can see that now the patch is symmetrically, with four big triangles covering most of the surface. This is enough to impliictly resolve the restoring torques produced by hydro, i.e. the instability is now removed.

TLDR;. This is the best solution that trades off complexity, cost, and still solves very effectively. My plan is to open a PR with the cleaned up solution from my prototype.

[sherm1]

For the final images above did you do both the remeshing and the COP change? Or is that just remeshing?

[amcastro-tri]

did you do both the remeshing and the COP change?

One change at a time:
For CoP, just CoP.
For remeshing, just remeshing.