some collision shapes fall through a flat terrain with contact filtering

[stephengold]

As discussed at the JME forum: https://hub.jmonkeyengine.org/t/rounded-object-falls-through-terrain/47584

Changes to the test app that prevent fallthrough:

1. disabling contact filtering on the terrain
2. changing the ball's collision shape to a SphereCollisionShape or MultiSphere
3. changing the ball's radius to 0.2 PSU or 10 PSU

Changes to the test app that don't prevent fallthrough:

1. changing the ball's collision shape to a GImpactCollisionShape
2. changing the ball's initial location to (0, 20, 0.001f) or (0.002f, 20, 0.001f) or (0.002f, 20.01f, 0.001f)

[stephengold]

Experiments using the DropTest app show that only a small fraction of randomized rigid bodies pass through filtered terrain. In particular, less than 1% of "frustum", "hemisphere", "hull", and "prism" drops pass through the "bedOfNails" and "dimples" platforms. (Note that "hemisphere" drops are custom convex shapes, not based on HullCollisionShape.) No passthru involving the "smooth" platform was observed.

Small shapes seem more likely to pass through than larger ones, which makes intuitive sense.

Initial focus will be on specific combinations of shapes known to be problematic, such as the one reported by ndebruyn.
Since the issue almost certainly originates in native code, I'll want a simple Libbulletjme app to reproduce it. I'll start by paring down ndebruyn's test to bare essentials.

[stephengold]

Changing the parameters of the Sphere constructor, from (20, 20, 1) to (20,19, 1) or (20, 21, 1) solves the issue. Changing them to (5, 7, 1) does not. That's good to know, since a hull shape with 23 vertices will be much easier to analyze than one with 384!

[stephengold]

Here's the pared-down test for Minie:

        flyCam.setEnabled(false);
        cam.setLocation(new Vector3f(-10.0f, 5.0f, 10.0f));
        cam.setRotation(new Quaternion(0.064f, 0.9106f, -0.156f, 0.377f));

        BulletAppState bulletAppState = new BulletAppState();
        bulletAppState.setDebugEnabled(true);
        stateManager.attach(bulletAppState);
        PhysicsSpace physicsSpace = bulletAppState.getPhysicsSpace();

        Sphere sphere = new Sphere(5, 7, 1f);
        HullCollisionShape hullShape = new HullCollisionShape(sphere);
        PhysicsRigidBody ballBody = new PhysicsRigidBody(hullShape);
        ballBody.setPhysicsLocation(new Vector3f(0f, 5f, 0f));
        physicsSpace.add(ballBody);

        CollisionShape heightShape = new HeightfieldCollisionShape(
                new float[9], new Vector3f(2f, 1f, 2f));
        PhysicsRigidBody terrainBody 
                = new PhysicsRigidBody(heightShape, PhysicsBody.massForStatic);
        //heightShape.setContactFilterEnabled(false);
        physicsSpace.add(terrainBody);

        Material solidGray = new Material(assetManager, Materials.UNSHADED);
        solidGray.setColor("Color", ColorRGBA.DarkGray);
        terrainBody.setDebugMaterial(solidGray);

[stephengold]

Activating CCD solves the issue, which I find surprising since the ball doesn't move very fast. During step 52, just before making contact, the ball moves about 0.145 PSU, about 7% of its diameter. At such low speeds, CCD shouldn't be necessary.

Without CCD:
...
tick 51 y=1.3866501 vy=-8.338497
tick 52 y=1.2449502 vy=-8.501997
tick 53 y=1.1005253 vy=-8.665497
tick 54 y=0.95337534 vy=-8.828997
tick 55 y=0.8035004 vy=-8.9924965
tick 56 y=0.6509005 vy=-9.155996

With CCD:

        ballBody.setCcdMotionThreshold(0.1f);
        ballBody.setCcdSweptSphereRadius(1f);

...
tick 51 y=1.3866501 vy=-8.338497
tick 52 y=1.2449502 vy=-8.501997
tick 53 y=1.1005253 vy=-8.665497
tick 54 y=1.04 vy=-3.6315176
tick 55 y=1.0412261 vy=0.07357079
tick 56 y=1.0417278 vy=0.03009659

[stephengold]

Debugging in Java as long as possible (to postpone the use of GDB) ...

With both contact filtering and CCD disabled, the maximum number of contact manifolds is 1.
There are ContactListener callbacks for the manifold, albeit a couple ticks later than I would expect:

tick 51  y=1.3866501  vy=-8.338497 numManifolds=0
tick 52  y=1.2449502  vy=-8.501997 numManifolds=0
tick 53  y=1.1005253  vy=-8.665497 numManifolds=1
tick 54  y=0.95337534  vy=-8.828997 numManifolds=1
  created manifold 140193185722400
  processed point 140193185722616
  processed point 140193185722408
tick 55  y=0.9631649  vy=-0.47508544 numManifolds=1
  processed point 140193185722616
  processed point 140193185722408
tick 56  y=0.97778815  vy=0.03134674 numManifolds=1
  processed point 140193185722616
  processed point 140193185722408

---

With contact filtering enabled (and CCD still disabled), the maximum number of contact manifolds is still 1:

tick 51  y=1.3866501  vy=-8.338497 numManifolds=0
tick 52  y=1.2449502  vy=-8.501997 numManifolds=0
tick 53  y=1.1005253  vy=-8.665497 numManifolds=1
tick 54  y=0.95337534  vy=-8.828997 numManifolds=1
tick 55  y=0.8035004  vy=-8.9924965 numManifolds=1
tick 56  y=0.6509005  vy=-9.155996 numManifolds=1
tick 57  y=0.49557555  vy=-9.319496 numManifolds=1
tick 58  y=0.3375256  vy=-9.482996 numManifolds=1
tick 59  y=0.17675067  vy=-9.646496 numManifolds=1
tick 60  y=0.013250738  vy=-9.809996 numManifolds=1
tick 61  y=-0.15297419  vy=-9.9734955 numManifolds=1
tick 62  y=-0.32192412  vy=-10.136995 numManifolds=1
tick 63  y=-0.49359906  vy=-10.300495 numManifolds=1
tick 64  y=-0.66799897  vy=-10.463995 numManifolds=1
tick 65  y=-0.8451239  vy=-10.627495 numManifolds=1
tick 66  y=-1.0249739  vy=-10.790995 numManifolds=1
  created manifold 140535113322528
  processed point 140535113322536
tick 67  y=-1.2075487  vy=-10.954494 numManifolds=1
  removed manifold 140535113322528
tick 68  y=-1.3928486  vy=-11.117994 numManifolds=1
tick 69  y=-1.5808735  vy=-11.281494 numManifolds=0

The original manifold (added during tick 52) has little or no effect on the ball's motion. Also, there seem to be no ContactListener callbacks for it. It's a good bet that these differences result from contact filtering.

A second manifold is created as the ball emerges out the bottom of the terrain, but by then it's too late to stop the ball's descent. (Minie's contact bookkeeping is confusing!)

---

With both contact filtering and CCD enabled, the maximum number of contact manifolds increases to 2:

tick 51  y=1.3866501  vy=-8.338497 numManifolds=0
tick 52  y=1.2449502  vy=-8.501997 numManifolds=0
tick 53  y=1.1005253  vy=-8.665497 numManifolds=1
tick 54  y=1.04  vy=-3.6315176 numManifolds=2
  removed manifold 140034432357264
  created manifold 140034432356384
  processed point 140034432356392
tick 55  y=1.0412261  vy=0.07357079 numManifolds=1
  processed point 140034432356392
tick 56  y=1.0417278  vy=0.03009659 numManifolds=1
  processed point 140034432356392

It's a good bet that the new manifold (added during tick 53) originates from CCD.

For debugging, it would nice to visualize HeightfieldCollisionShape in a way that shows the thickness of each triangle.

It also would be nice to visualize all contact points [as Bullet does ... see btCollisionWorld::debugDrawWorld()]. The first step would be to expose btDispatcher::getManifoldByIndexInternal() via JNI. That functionality might also help clarify the bookkeeping discrepancies.

[stephengold]

Actually, getManifoldByIndexInternal() is already exposed by PhysicsSpace.listManifoldIds().

[stephengold]

Enabling CCD doesn't solve the issue with ndebruyn's test. The ball still falls through, albeit more slowly.
So disabling contact filtering is a much better workaround.

[stephengold]

With contact filtering enabled and CCD disabled, the 1st persistent manifold is created without any contact points, which may explain why it's ineffective.

[stephengold]

Contact filtering is mainly implemented in btManifoldResult::addContactPoint(), which is first invoked during tick 54, 2 ticks after the btPersistentManifold was created. Comments in the Bullet source suggest why the btPersistentManifold might be created before there are any points to add.
Call stack:

#0  btManifoldResult::addContactPoint
#1  btGjkPairDetector::getClosestPointsNonVirtual
#2  btGjkPairDetector::getClosestPoints
#3  btConvexConvexAlgorithm::processCollision 
#4  btConvexTriangleCallback::processTriangle
#5  btHeightfieldTerrainShape::processAllTriangles
#6  btConvexConcaveCollisionAlgorithm::processCollision
#7  btCollisionDispatcher::defaultNearCallback
#8  btCollisionPairCallback::processOverlap
#9  btHashedOverlappingPairCache::processAllOverlappingPairs
#10 btHashedOverlappingPairCache::processAllOverlappingPairs
#11 btCollisionDispatcher::dispatchAllCollisionPairs
#12 btCollisionWorld::performDiscreteCollisionDetection
#13 btDiscreteDynamicsWorld::internalSingleStepSimulation
#14 btDiscreteDynamicsWorld::stepSimulation

isSwapped is false
isNewCollision is true
pcoA is the btConvexHullShape (type=4)
isValidContact(localA, -1, -1) returns true
pcoB is the HeightfieldShape (type=24)
localB is (-0.222395927, 0.0381531119, 0.00952136237)
m_partId1 and m_index1 are both 0
isValidContact(localB, 0, 0) returns false

Stepping through isValidContact(localB, 0, 0):
margin set to 0.0399600007
aabbMin set to (-0.262355924, -0.00180688873, -0.0304386392)
aabbMax set to (-0.18243593, 0.0781131089, 0.0494813621)

Stepping through btHeightfieldTerrainShape::processAllTriangles():
m_localScaling is (2, 1, 2) (could scaling of margin be handled better?)
m_localOrigin is (1, 0, 1)
startJ and startX set to 0
endJ and endZ set to 2

All 9 heightfield squares will be processed. I set a breakpoint on btTriangleCallback.h line 51...
1st time: partId=0, triangleIndex=0 -> early return from line 51
2nd time: partId=1, triangleIndex=0 -> isInside() returns false
partId=2, triangleIndex=0 -> false
partId=3, triangleIndex=0 -> false
partId=0, triangleIndex=1 -> false
partId=1, triangleIndex=1 -> true (could we implement early return from processAllTriangles()?)
partId=2, triangleIndex=1 -> false
partId=3, triangleIndex=1 -> false

After that, the addContactPoint() breakpoint was hit, so only 8 heightfield triangles were processed. I guess the other 10 triangles failed the bounding-box test in btHeightfieldTerrainShape::processAllTriangles().

This time:
isSwapped is false
isNewCollision is true
pcoA is the btConvexHullShape (type=4)
isValidContact(localA, -1, -1) returns true
pcoB is the HeightfieldShape (type=24)
localB is (-0.014738081, 0.0380637944, 0.00264857057)
m_partId1 is 1 and m_index1 is 0

The next addContactPoint() has m_partId1=2 and mIndex1=0. Again, it is isValidContact(localB, ...) that returns false

Next step: instrument the code to get a more complete picture of which triangles invalidate which contact points.

[stephengold]

Actually, there are only 4 heightfield squares, thus 8 triangles, none of which failed the bounding-box test.

Note btHeightfieldTerrainShape::getVertex() factors in m_localScaling, so the btTriangleShape is defined in the (non-uniformly) scaled local coordinates of the heightfield.

[stephengold]

get a more complete picture of which triangles invalidate which contact points

tick 52  y=1.2449502  vy=-8.501997  manifolds={}
tick 53  y=1.1005253  vy=-8.665497  manifolds={7f1359960020<0>}
tick 54  y=0.95337534  vy=-8.828997  manifolds={7f1359960020<0>}
Contact at (-0.222396,0.0381531,0.00952136) on tri(0,0) inside tri(1,1)
 coords (-2,0,0)(0,0,2)(0,0,0)
 marg=0.039960 oldCnt=0
Contact at (-0.0147381,0.0380638,0.00264857) on tri(1,0) inside tri(0,0)
 coords (-2,0,-2)(-2,0,0)(0,0,0)
 marg=0.039960 oldCnt=0
Contact at (-0.0147381,0.0380638,0.00264857) on tri(1,0) inside tri(1,1)
 coords (-2,0,0)(0,0,2)(0,0,0)
 marg=0.039960 oldCnt=1
Contact at (-0.00890616,0.0389959,0) on tri(2,0) inside tri(0,0)
 coords (-2,0,-2)(-2,0,0)(0,0,0)
 marg=0.039960 oldCnt=0
Contact at (-0.00890616,0.0389959,0) on tri(2,0) inside tri(1,0)
 coords (-2,0,-2)(0,0,0)(0,0,-2)
 marg=0.039960 oldCnt=1
Contact at (-0.00890616,0.0389959,0) on tri(2,0) inside tri(1,1)
 coords (-2,0,0)(0,0,2)(0,0,0)
 marg=0.039960 oldCnt=2
Contact at (-0.2224,0.0381437,-0.00956095) on tri(1,1) inside tri(0,0)
 coords (-2,0,-2)(-2,0,0)(0,0,0)
 marg=0.039960 oldCnt=0
Contact at (-0.00890623,0.0389959,0) on tri(2,1) inside tri(0,0)
 coords (-2,0,-2)(-2,0,0)(0,0,0)
 marg=0.039960 oldCnt=0
Contact at (-0.00890623,0.0389959,0) on tri(2,1) inside tri(1,0)
 coords (-2,0,-2)(0,0,0)(0,0,-2)
 marg=0.039960 oldCnt=1
Contact at (-0.00890623,0.0389959,0) on tri(2,1) inside tri(1,1)
 coords (-2,0,0)(0,0,2)(0,0,0)
 marg=0.039960 oldCnt=2
Contact at (-0.00890623,0.0389959,0) on tri(3,1) inside tri(0,0)
 coords (-2,0,-2)(-2,0,0)(0,0,0)
 marg=0.039960 oldCnt=0
Contact at (-0.00890623,0.0389959,0) on tri(3,1) inside tri(1,0)
 coords (-2,0,-2)(0,0,0)(0,0,-2)
 marg=0.039960 oldCnt=1
Contact at (-0.00890623,0.0389959,0) on tri(3,1) inside tri(1,1)
 coords (-2,0,0)(0,0,2)(0,0,0)
 marg=0.039960 oldCnt=2
tick 55  y=0.8035004  vy=-8.9924965  manifolds={7f1359960020<0>}

[stephengold]

Surprise: 6 contacts, all with Y < 0.03996
To sort this out, I'll need a diagram.

[stephengold]

All 6 contacts lie close together near the center of the grid, where 6 triangles meet...

The 3 contacts with Z=0 lie on the -X side of the grid, on the edge shared by tri(0,0) and tri(1,1):

• at (-0.00890616,0.0389959,0) on tri(2,0)
• at (-0.00890623,0.0389959,0) on tri(2,1)
• at (-0.00890623,0.0389959,0) on tri(3,1)

The other 3 contacts (with non-zero Z) also lie very near that same edge:

• at (-0.222396,0.0381531,0.00952136) on tri(0,0)
• at (-0.0147381,0.0380638,0.00264857) on tri(1,0)
• at (-0.2224,0.0381437,-0.00956095) on tri(1,1)

I still need to understand why none of the Y values is 0.040, which probably means understanding the Gilbert-Johnson-Keerthi distance algorithm.

[stephengold]

In getClosestPointsNonVirtual(), pointInWorld is calculated as pointOnB + positionOffset.

For the first invocation of addContactPoint() that is (-0.222395927, -0.438534558, 0.00952136237) + (0, 0.47668767, 0). The Y component is the sum of positive and negative scalars that have almost the same magnitude, making it vulnerable to rounding error. I re-ran the test with "DebugDp" natives, but the ball still falls through.

If PointOnB is inaccurate, perhaps it's because the G-J-K loop was terminated too soon?

[stephengold]

Still using "DebugDp" natives ...

The first loop in getClosestPointsNonVirtual exits from line 794 with iterations == 4.

The second loop exits from line 911 with m_cachedSeparatingAxis == (0,0,0).

m_lastUsedMethod gets set to 2.

btGjkEpaPenetrationDepthSolver::calcPenDepth is invoked at line 1022 (see NarrowPhaseCollision/btGjkEpaPenetrationDepthSolver.cpp). It returns true (valid penetration) from line 62.

isValid changes from false to true at line 1053, with the following results:

• distance = -0.099983147638179837
• pointOnA = (-0.22165398902262695, -0.53229978103984954, -0.022477081784502237)
• pointOnB = (-0.22165398902262751, -0.43960173343572007, 0.014988926570013071)
• normalInB = (-5.5520507748107931e-15, 0.92713672047599693, 0.37472323326024631)

Note the normal isn't vertical, and it's the same as the normal passed to addContactPoint().
Note positionOffset is the midpoint between the centers of the 2 shapes.

[stephengold]

Reading this article has me wondering whether calcPenDepth actually implements the Expanding Polytope Algorithm. The code in btGjkEpaPenetrationDepthSolver.cpp sure looks like it is guessing.

[stephengold]

At the start of tick 54, the lowest vertex of the hull shape is vertex[0] with world location (-0.22234385, -0.021553636, 0.0)

[stephengold]

Stepping through btGjkPairDetector::getClosestPointsNonVirtual() again, in double precision.
Note: A is the hull and B is the heightfield triangle.
No notable progress.

Changing btDefaultCollisionConstructionInfo so that the default is m_useEpaPenetrationAlgorithm(false) causes the btMinkowskiPenetrationDepthSolver to be used in place of btGjkEpaPenetrationDepthSolver. That solves this issue for the 23-vertex test case, both with double and single precision. It also solves the issue for ndebruyn's test app.

Next step will be to add mechanisms to configure the depth solver at runtime and compare performance.

[stephengold]

After many digressions, I got the Minkowski PDS working in Minie. Using DropTest, I dropped "prism" drops onto the "dimples" platform (with contact filtering enabled) and observed >1% fallthrough. So Minkowski PDS doesn't solve the general issue.

[stephengold]

I've been busy with other issues.

Unless this issue is solved before the next Minie release, I'll probably disable contact filtering by default. (The option has been enabled by default since it was introduced in 2021.)

[stephengold]

I disabled contact filtering by default in Libbulletjme v21.2.1, but now I'm having second thoughts about that decision.

In DropTest of "prism" drops on the "dimples" platform, I'm seeing occasional fallthrough, even with filtering disabled and the CCD motion thresholds reduced to 1. The "candyDish" platform (a MeshCollisionShape) and the "smooth" platform (a HeightfieldCollisionShape) don't seem to have this issue. The "bedOfNails" has it, but only very rarely.

The next step would be to create more simple testcases for detailed study, including some that fail even with filtering disabled.

[stephengold]

Simple testcases aren't easy to find.
The prism-on-dimples fallthroughs without filtering are easily produced with 80 dynamic bodies, but difficult to produce using less than 12.

[stephengold]

Focusing on drops that fall through: after studying more than a dozen instances, all were among the first 16 drops created after restarting the scenario.
I have a hunch that those drops are getting hammered through the platform when later drops land on top of them.
If that's true, I might be able to increase the rate of fallthrough by concentrating drops in a smaller area.

[stephengold]

Concentrating the drops helps a little, but probably not enough to "crack" this issue.

It occurs to me that HeightfieldCollisionShape, being a thin, two-sided shape (typically 0.08 psu thick) isn't a good fit for typical (infinitely thick) use cases. Increasing its thickness might make it more bullet-proof. (Pardon the choice of words!) Currently, thickness can be enhanced by increasing the collision margin, but this can only be taken so far before margin starts destroying surface details.

[stephengold]

Bullet collision shapes for rigid bodies fall into 3 categories: convex (defined by supporting vertices), concave (composed of thin triangles), and compound (composed of sub-shapes). Perhaps there's a need for a new terrain shape, one with greater thickness, implemented more like a compound shape than a concave.