stephengold / Minie

Integrate Bullet Physics and V-HACD into jMonkeyEngine projects (code has New BSD license)
https://stephengold.github.io/Minie/minie/overview.html
Other
122 stars 22 forks source link

BetterCharacterControl hops across seams #18

Closed stephengold closed 2 years ago

stephengold commented 2 years ago

When a BetterCharacterController moves back and forth across the seam between 2 triangles in a flat MeshCollisionShape, it sometimes makes an unexpected vertical hop.

This issue was reported (in JMonkeyEngine) as long ago as December 2013: https://hub.jmonkeyengine.org/t/bettercharactercontrol-acts-like-trampoline/28186 and also as recently as October 2021: https://hub.jmonkeyengine.org/t/bettercharactercontroller-bouncing-off-of-the-ground-while-moving/44991.

Thanks to @jcfandino, we now have a simple test that reproduces the issue. Note that the pauses (desiredVelocity.zero()) are necessary to reproduce the issue.

package jme3utilities.minie.test.issue;

import com.jme3.app.SimpleApplication;
import com.jme3.bullet.BulletAppState;
import com.jme3.bullet.PhysicsSpace;
import com.jme3.bullet.PhysicsTickListener;
import com.jme3.bullet.collision.shapes.CollisionShape;
import com.jme3.bullet.control.BetterCharacterControl;
import com.jme3.bullet.control.RigidBodyControl;
import com.jme3.bullet.objects.PhysicsRigidBody;
import com.jme3.bullet.util.CollisionShapeFactory;
import com.jme3.math.FastMath;
import com.jme3.math.Vector3f;
import com.jme3.scene.Geometry;
import com.jme3.scene.Mesh;
import com.jme3.scene.Node;
import com.jme3.scene.Spatial;
import com.jme3.scene.shape.Quad;

/**
 * Reproduce the hopping issue in BetterCharacterControl. If the issue is
 * present, numeric data will be printed to System.out .
 *
 * @author Stephen Gold sgold@sonic.net
 */
public class TestIssueXX
        extends SimpleApplication
        implements PhysicsTickListener {
    // *************************************************************************
    // fields

    /**
     * control under test
     */
    private BetterCharacterControl bcc;
    /**
     * true if character will move toward +X, false if it will move toward -X
     */
    private boolean increasingX;
    /**
     * largest Y value seen so far: anything larger than 0.05 is an issue
     */
    private float maxHeight = 0.05f;
    /**
     * count of physics timesteps simulated
     */
    private int tickCount = 0;
    // *************************************************************************
    // new methods exposed

    public static void main(String[] ignored) {
        TestIssueXX application = new TestIssueXX();
        application.start();
    }
    // *************************************************************************
    // SimpleApplication methods

    @Override
    public void simpleInitApp() {
        PhysicsSpace physicsSpace = configurePhysics();
        addGround(physicsSpace);

        Node controlledNode = new Node("controlled node");
        rootNode.attachChild(controlledNode);

        float characterRadius = 1f;
        float characterHeight = 4f;
        float characterMass = 1f;
        bcc = new BetterCharacterControl(characterRadius, characterHeight,
                characterMass);
        controlledNode.addControl(bcc);
        physicsSpace.add(bcc);
    }

    @Override
    public void simpleUpdate(float tpf) {
        /*
         * Terminate the test after 200 time steps.
         */
        if (tickCount > 200) {
            stop();
        }
    }
    // *************************************************************************
    // PhysicsTickListener methods

    @Override
    public void physicsTick(PhysicsSpace space, float timeStep) {
        /*
         * Determine the character's height and print it if it's a new high.
         */
        PhysicsRigidBody body = bcc.getRigidBody();
        Vector3f location = body.getPhysicsLocation();
        if (location.y > maxHeight) {
            maxHeight = location.y;
            System.out.println(tickCount + ": " + location);
        }
    }

    @Override
    public void prePhysicsTick(PhysicsSpace space, float timeStep) {
        ++tickCount;
        /*
         * Walk rapidly back and forth across the seam between the 2 triangles.
         */
        Vector3f desiredVelocity = new Vector3f();
        float walkSpeed = 99f;
        if (increasingX) {
            desiredVelocity.x = walkSpeed;
        } else {
            desiredVelocity.x = -walkSpeed;
        }

        Vector3f location = bcc.getRigidBody().getPhysicsLocation();
        if (increasingX && location.x > 7f) {
            // stop and reverse direction
            desiredVelocity.zero();
            increasingX = false;
        } else if (!increasingX && location.x < -7f) {
            // stop and reverse direction
            desiredVelocity.zero();
            increasingX = true;
        }

        bcc.setWalkDirection(desiredVelocity);
    }
    // *************************************************************************
    // private methods

    /**
     * Add a ground body to the specified PhysicsSpace.
     *
     * @param physicsSpace (not null)
     */
    private void addGround(PhysicsSpace physicsSpace) {
        Mesh quad = new Quad(1000f, 1000f);
        Spatial ground = new Geometry("ground", quad);
        ground.move(-500f, 0f, 500f);
        ground.rotate(-FastMath.HALF_PI, 0f, 0f);

        CollisionShape shape = CollisionShapeFactory.createMeshShape(ground);
        RigidBodyControl rbc = new RigidBodyControl(shape, 0f);
        rbc.setPhysicsSpace(physicsSpace);
        ground.addControl(rbc);
    }

    /**
     * Configure physics during startup.
     */
    private PhysicsSpace configurePhysics() {
        BulletAppState bulletAppState = new BulletAppState();
        stateManager.attach(bulletAppState);

        PhysicsSpace result = bulletAppState.getPhysicsSpace();
        result.addTickListener(this);

        return result;
    }
}
stephengold commented 2 years ago

Note that the test app relies on the app's tick listener being added before the BCC's tick listener!

The first anomaly detected is in TestIssue18.physicsTick() with tickCount=12:

Libbulletjme version 12.2.2 initializing
12: (1.4449339, 0.15935323, -0.009808682)

This is highly reproducible. The character's post-tick location with tickCount=11 was (-0.16318321, 0.036596943, -2.5927548E-7), so the anomaly happens when character crosses the seam for the first time.

During tickCount=12, onGround gets set to true, so the issue isn't with the raycast. May need to dive into native code.

stephengold commented 2 years ago

On a hunch, I replaced the CapsuleCollisionShape in BCC with a MultiSphere, and that solved the issue. So probably the root cause is in Bullet's btCapsuleShape class.

Before implementing this workaround, try to pin down the root cause, since it presumably affects non-BCC rigid bodies as well.

stephengold commented 2 years ago

Today the MultiSphere workaround isn't working for me. Not sure what happened. Perhaps I was drowsy last night and got confused.

stephengold commented 2 years ago

Shortening the physics timestep slightly (from 0.0166667 to 0.013) resolves the issue, at least for TestIssue18. However, doubling the walk speed (from 99 to 198) causes the issue to return. From there, further shortening of the timestep doesn't seem to help

stephengold commented 2 years ago

Another symptom of this issue is that when tickCount=12, the applied impulse of the new contact (obtained using addCollisionListener() and getAppliedImpulse()) is about 50x larger than expected:

Debug_Libbulletjme version 12.2.2 initializing
1: (-1.6473188, 0.007999993, 0.0)
2: (-3.2970002, 0.014549839, -7.7255095E-8)
3: (-4.946684, 0.019639885, -1.0401872E-7)
4: (-6.5962634, 0.023687199, -1.066944E-7)
index = 1, impulse = 0.100936
index = 1, impulse = 0.100936
index = 1, impulse = 0.100936
index = 1, impulse = 0.100936
5: (-8.245772, 0.026974354, -1.0696189E-7)
6: (-8.410214, 0.029572152, -7.671717E-8)
7: (-6.7607675, 0.031657696, -1.0648709E-7)
8: (-5.111332, 0.033277355, -1.4287345E-7)
index = 1, impulse = 0.135546
index = 1, impulse = 0.135546
index = 1, impulse = 0.135546
index = 1, impulse = 0.135546
9: (-3.4619443, 0.034621917, -1.8268327E-7)
index = 1, impulse = 0.146994
10: (-1.8125705, 0.03574617, -1.8666273E-7)
index = 1, impulse = 0.150283
11: (-0.16318321, 0.036596943, -2.5927548E-7)
index = 1, impulse = 0.147091
12: (1.4449339, 0.15935323, -0.009808682)
index = 0, impulse = 7.674518

Need to look into how the impulse is calculated.

stephengold commented 2 years ago

With addOngoingCollisionListener() the tickCount=12 collision events look even stranger:

Debug_Libbulletjme version 12.2.2 initializing
1: (-1.6473188, 0.007999993, 0.0)
2: (-3.2970002, 0.014549839, -7.7255095E-8)
3: (-4.946684, 0.019639885, -1.0401872E-7)
4: (-6.5962634, 0.023687199, -1.066944E-7)
index = 1, impulse = 0.100936
index = 1, impulse = 0.100936
index = 1, impulse = 0.100936
index = 1, impulse = 0.100936
5: (-8.245772, 0.026974354, -1.0696189E-7)
6: (-8.410214, 0.029572152, -7.671717E-8)
7: (-6.7607675, 0.031657696, -1.0648709E-7)
8: (-5.111332, 0.033277355, -1.4287345E-7)
index = 1, impulse = 0.135546
index = 1, impulse = 0.135546
index = 1, impulse = 0.135546
index = 1, impulse = 0.135546
9: (-3.4619443, 0.034621917, -1.8268327E-7)
index = 1, impulse = 0.146994
10: (-1.8125705, 0.03574617, -1.8666273E-7)
11: (-0.16318321, 0.036596943, -2.5927548E-7)
index = 1, impulse = 0.147091
index = 1, impulse = 0.147091
12: (1.4449339, 0.15935323, -0.009808682)
index = 1, impulse = 0.000000
index = 0, impulse = 7.674518
13: (3.094934, 0.27938452, -0.010789524)
14: (4.744934, 0.39669082, -0.010887608)
stephengold commented 2 years ago

The m_appliedImpulse of a btManifoldPoint is apparently set at btSequentialImpulseConstraintSolver.cpp:1766, copied from a btSolverConstraint in the contact constraint pool. In practice, that field seems to be set 5 places in the same file, in lines 61-95.

stephengold commented 2 years ago

By adding printfs to the btSISC.cpp, I convinced myself that the return value of getAppliedImpulse() is the value from the previous timestep, kept around for warm starting.

The mysterious 7.674518 impulse is set at line 95:

LowerLimit95** deltaImpulse=7.674517 m_appliedImpulse=0.000000  m_rhs=7.674517 m_cfm=0.000000 deltaVel1Dotn=0.000000 deltaVel2Dotn=0.000000 m_jacDiagABInv=1.000000

7.674517 is an abnormally large value for m_rhs. There are only a few places where m_rhs could be set. I'll pursue that next.

stephengold commented 2 years ago

It's set in btSequentialImpulseConstraintSolver::setupContactConstraint(), specifically the if clause:

        if (!infoGlobal.m_splitImpulse || (penetration > infoGlobal.m_splitImpulsePenetrationThreshold))
        {
            //combine position and velocity into rhs
            solverConstraint.m_rhs = penetrationImpulse + velocityImpulse;  //-solverConstraint.m_contactNormal1.dot(bodyA->m_externalForce*bodyA->m_invMass-bodyB->m_externalForce/bodyB->m_invMass)*solverConstraint.m_jacDiagABInv;
stephengold commented 2 years ago

velocityImpulse=7.674517, which results from velocityError=7.674517, which results from rel_vel=-7.854604 (on earlier timesteps, rel_vel was more like -0.1) which results from vel1Dotn=-7.854604 (previously -0.096045), which results from a slight -X displacement of m_contactNormal1.

Step 10 had:

948     m_contactNormal1=(0.000000 1.000000 -0.000029) m_linearVelocity=(99.000000 0.067455 -0.000000) externalForceImpulseA=(0.000000 -0.163500 0.000000)
        m_relpos1CrossNormal=(-0.000046 -0.000000 -0.000021) m_angularVelocity=(0.000000 0.000000 0.000000) externalTorqueImpulseA=(0.000000 0.000000 0.000000)

whereas step 11 has:

948     m_contactNormal1=(-0.078211 0.993862 -0.078240) m_linearVelocity=(99.000000 0.051046 -0.000000) externalForceImpulseA=(0.000000 -0.163500 0.000000)
        m_relpos1CrossNormal=(-0.078242 0.000000 0.078217) m_angularVelocity=(0.000000 0.000000 0.000000) externalTorqueImpulseA=(0.000000 0.000000 0.000000)

That partly explains why the issue is more pronounced at high speeds.

stephengold commented 2 years ago

On both step 10 and step 11, m_contactNormal1 is set from cp.m_normalWorldOnB, where cp is a function argument. Here's the call stack:

(gdb) bt
#0  btSequentialImpulseConstraintSolver::setupContactConstraint(btSolverConstraint&, int, int, btManifoldPoint&, btContactSolverInfo const&, float&, btVector3 const&, btVector3 const&)
    (this=0x7f512500fa30, solverConstraint=..., solverBodyIdA=0, solverBodyIdB=1, cp=..., infoGlobal=..., relaxation=@0x7f512993dffc: 1, rel_pos1=..., rel_pos2=...)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletDynamics/ConstraintSolver/btSequentialImpulseConstraintSolver.cpp:887
#1  0x00007f508ecbc33e in btSequentialImpulseConstraintSolver::convertContact(btPersistentManifold*, btContactSolverInfo const&)
    (this=0x7f512500fa30, manifold=0x7f50782dc020, infoGlobal=...)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletDynamics/ConstraintSolver/btSequentialImpulseConstraintSolver.cpp:1094
#2  0x00007f508ecbcd08 in btSequentialImpulseConstraintSolver::convertContacts(btPersistentManifold**, int, btContactSolverInfo const&)
    (this=0x7f512500fa30, manifoldPtr=0x7f5124631aa0, numManifolds=1, infoGlobal=...)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletDynamics/ConstraintSolver/btSequentialImpulseConstraintSolver.cpp:1200
#3  0x00007f508ecbe512 in btSequentialImpulseConstraintSolver::solveGroupCacheFriendlySetup(btCollisionObject**, int, btPersistentManifold**, int, btTypedConstraint**, int, btContactSolverInfo const&, btIDebugDraw*)
    (this=0x7f512500fa30, bodies=0x7f5124631a70, numBodies=1, manifoldPtr=0x7f5124631aa0, numManifolds=1, constraints=0x0, numConstraints=0, infoGlobal=..., debugDrawer=0x0)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletDynamics/ConstraintSolver/btSequentialImpulseConstraintSolver.cpp:1533
#4  0x00007f508ecc02c9 in btSequentialImpulseConstraintSolver::solveGroup(btCollisionObject**, int, btPersistentManifold**, int, btTypedConstraint**, int, btContactSolverInfo const&, btIDebugDraw*, btDispatcher*) (this=0x7f512500fa30, bodies=0x7f5124631a70, numBodies=1, manifoldPtr=0x7f5124631aa0, numManifolds=1, constraints=0x0, numConstraints=0, infoGlobal=..., debugDrawer=0x0)
     at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletDynamics/ConstraintSolver/btSequentialImpulseConstraintSolver.cpp:1905
#5  0x00007f508edf3230 in InplaceSolverIslandCallback::processConstraints() (this=0x7f512500fc80) at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletDynamics/Dynamics/btDiscreteDynamicsWorld.cpp:184
#6  0x00007f508edee1b8 in btDiscreteDynamicsWorld::solveConstraints(btContactSolverInfo&) (this=0x7f5124927900, solverInfo=...) at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletDynamics/Dynamics/btDiscreteDynamicsWorld.cpp:704
#7  0x00007f508eded50f in btDiscreteDynamicsWorld::internalSingleStepSimulation(float) (this=0x7f5124927900, timeStep=0.0166666675) at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletDynamics/Dynamics/btDiscreteDynamicsWorld.cpp:480
stephengold commented 2 years ago

cp was set by

btManifoldPoint& cp = manifold->getContactPoint(j);

in convertContact() with j=0. Step 11 is the first timestep in which the manifold has 2 contacts.

Now to figure out where that m_normalWorldOnB X-value got set. I suspect it's in the btManifoldPoint constructor, but that needs to be verified.

stephengold commented 2 years ago

It's the 4-arg constructor invoked at btManifoldResult.cpp:131. Here's the call stack:

#0  btManifoldPoint::btManifoldPoint(btVector3 const&, btVector3 const&, btVector3 const&, float)
    (this=0x7fbf3400ac70, pointA=..., pointB=..., normal=..., distance=0.00299906731)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/NarrowPhaseCollision/btManifoldPoint.h:104
#1  0x00007fbeb8378c61 in btManifoldResult::addContactPoint(btVector3 const&, btVector3 const&, float)
    (this=0x7fbf3400c120, normalOnBInWorld=..., pointInWorld=..., depth=0.00299906731)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btManifoldResult.cpp:131
#2  0x00007fbeb83f1c8d in btGjkPairDetector::getClosestPointsNonVirtual(btDiscreteCollisionDetectorInterface::ClosestPointInput const&, btDiscreteCollisionDetectorInterface::Result&, btIDebugDraw*)
    (this=0x7fbf3400b380, input=..., output=..., debugDraw=0x0)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/NarrowPhaseCollision/btGjkPairDetector.cpp:1174
#3  0x00007fbeb83eea5e in btGjkPairDetector::getClosestPoints(btDiscreteCollisionDetectorInterface::ClosestPointInput const&, btDiscreteCollisionDetectorInterface::Result&, btIDebugDraw*, bool)
    (this=0x7fbf3400b380, input=..., output=..., debugDraw=0x0, swapResults=false)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/NarrowPhaseCollision/btGjkPairDetector.cpp:77
#4  0x00007fbeb837ccfe in btConvexConvexAlgorithm::processCollision(btCollisionObjectWrapper const*, btCollisionObjectWrapper const*, btDispatcherInfo const&, btManifoldResult*)
    (this=0x7fbf2d060d60, body0Wrap=0x7fbf3400bc20, body1Wrap=0x7fbf3400b790, dispatchInfo=..., resultOut=0x7fbf3400c120)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btConvexConvexAlgorithm.cpp:699
#5  0x00007fbeb82f2e8a in btConvexTriangleCallback::processTriangle(btVector3*, int, int)
   (this=0x7fbf2d060cb0, triangle=0x7fbf3400ba18, partId=0, triangleIndex=0)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btConvexConcaveCollisionAlgorithm.cpp:135
#6  0x00007fbeb84a8138 in btBvhTriangleMeshShape::MyNodeOverlapCallback::processNode(int, int)
    (this=0x7fbf3400ba00, nodeSubPart=0, nodeTriangleIndex=0)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionShapes/btBvhTriangleMeshShape.cpp:319
#7  0x00007fbeb826ecd0 in btQuantizedBvh::walkStacklessQuantizedTree(btNodeOverlapCallback*, unsigned short*, unsigned short*, int, int) const
    (this=0x7fbf2c62a460, nodeCallback=0x7fbf3400ba00, quantizedQueryAabbMin=0x7fbf3400b9bc, quantizedQueryAabbMax=0x7fbf3400b9c2, startNodeIndex=0, endNodeIndex=3)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/BroadphaseCollision/btQuantizedBvh.cpp:697
#8  0x00007fbeb826dbeb in btQuantizedBvh::reportAabbOverlappingNodex(btNodeOverlapCallback*, btVector3 const&, btVector3 const&) const
    (this=0x7fbf2c62a460, nodeCallback=0x7fbf3400ba00, aabbMin=..., aabbMax=...) at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/BroadphaseCollision/btQuantizedBvh.cpp:327
#9  0x00007fbeb84a81c5 in btBvhTriangleMeshShape::processAllTriangles(btTriangleCallback*, btVector3 const&, btVector3 const&) const
    (this=0x7fbf2c62a3d0, callback=0x7fbf2d060cb0, aabbMin=..., aabbMax=...)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionShapes/btBvhTriangleMeshShape.cpp:326
#10 0x00007fbeb82f38b0 in btConvexConcaveCollisionAlgorithm::processCollision(btCollisionObjectWrapper const*, btCollisionObjectWrapper const*, btDispatcherInfo const&, btManifoldResult*)
    (this=0x7fbf2d060ca0, body0Wrap=0x7fbf3400bc20, body1Wrap=0x7fbf3400c0c0, dispatchInfo=..., resultOut=0x7fbf3400c120)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btConvexConcaveCollisionAlgorithm.cpp:262
#11 0x00007fbeb831a0ea in btCompoundLeafCallback::ProcessChildShape(btCollisionShape const*, int)
    (this=0x7fbf3400bf10, childShape=0x7fbf2cff54f0, index=0)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btCompoundCollisionAlgorithm.cpp:183
#12 0x00007fbeb831a227 in btCompoundLeafCallback::Process(btDbvtNode const*)
    (this=0x7fbf3400bf10, leaf=0x7fbf2cff64d0)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btCompoundCollisionAlgorithm.cpp:226
#13 0x00007fbeb8319919 in btDbvt::collideTVNoStackAlloc(btDbvtNode const*, btDbvtAabbMm const&, btAlignedObjectArray<btDbvtNode const*>&, btDbvt::ICollide&) const
    (this=0x7fbf2cff7f30, root=0x7fbf2cff64d0, vol=..., stack=..., policy=...)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/BroadphaseCollision/btDbvt.h:1215
#14 0x00007fbeb831939c in btCompoundCollisionAlgorithm::processCollision(btCollisionObjectWrapper const*, btCollisionObjectWrapper const*, btDispatcherInfo const&, btManifoldResult*)
    (this=0x7fbf2d060be0, body0Wrap=0x7fbf3400c0c0, body1Wrap=0x7fbf3400c0f0, dispatchInfo=..., resultOut=0x7fbf3400c120)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btCompoundCollisionAlgorithm.cpp:293
#15 0x00007fbeb8369eba in btCollisionDispatcher::defaultNearCallback(btBroadphasePair&, btCollisionDispatcher&, btDispatcherInfo const&)
    (collisionPair=..., dispatcher=..., dispatchInfo=...)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btCollisionDispatcher.cpp:256
#16 0x00007fbeb836a607 in btCollisionPairCallback::processOverlap(btBroadphasePair&)
    (this=0x7fbf3400c2e0, pair=...)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btCollisionDispatcher.cpp:212
#17 0x00007fbeb8269df5 in btHashedOverlappingPairCache::processAllOverlappingPairs(btOverlapCallback*, btDispatcher*)
    (this=0x7fbf2d014960, callback=0x7fbf3400c2e0, dispatcher=0x7fbf2d120bf0)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/BroadphaseCollision/btOverlappingPairCache.cpp:341
#18 0x00007fbeb826a128 in btHashedOverlappingPairCache::processAllOverlappingPairs(btOverlapCallback*, btDispatcher*, btDispatcherInfo const&)
    (this=0x7fbf2d014960, callback=0x7fbf3400c2e0, dispatcher=0x7fbf2d120bf0, dispatchInfo=...)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/BroadphaseCollision/btOverlappingPairCache.cpp:412
#19 0x00007fbeb8369c44 in btCollisionDispatcher::dispatchAllCollisionPairs(btOverlappingPairCache*, btDispatcherInfo const&, btDispatcher*)
    (this=0x7fbf2d120bf0, pairCache=0x7fbf2d014960, dispatchInfo=..., dispatcher=0x7fbf2d120bf0)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btCollisionDispatcher.cpp:225
#20 0x00007fbeb825ae92 in btCollisionWorld::performDiscreteCollisionDetection()
    (this=0x7fbf2c8ffbb0)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletCollision/CollisionDispatch/btCollisionWorld.cpp:235
#21 0x00007fbeb8440479 in btDiscreteDynamicsWorld::internalSingleStepSimulation(float)
    (this=0x7fbf2c8ffbb0, timeStep=0.0166666675)
    at /home/sgold/Git/Libbulletjme/src/main/native/bullet3/BulletDynamics/Dynamics/btDiscreteDynamicsWorld.cpp:473
stephengold commented 2 years ago

Looks like it comes from btGjkPairDetector::getClosestPointsNonVirtual(). Specifically, in the if (checkSimplex) clause:

                normalInB = m_cachedSeparatingAxis;

and

                    btScalar rlen = btScalar(1.) / btSqrt(lenSqr);
                    normalInB *= rlen;  //normalize

There oughta be a law against functions that are 500 lines long!

stephengold commented 2 years ago

Here are the normalInB values on step 11:

974    (-0.081573 1.036597 -0.081604)
986    (-0.078211 0.993862 -0.078240)
974    (0.000092 1.036597 -0.000046)
986    (0.000088 1.000000 -0.000044)

Compare with step 10:

974    (-0.906296 1.035746 -0.906281)
986    (-0.549979 0.628535 -0.549970)
974    (0.000000 1.035746 -0.000031)
986    (0.000000 1.000000 -0.000029)

So now looking for where m_cachedSeparatingAxis got set to (-0.081573 1.036597 -0.081604).

stephengold commented 2 years ago

There's an ominous comment at line 735:

    //note, for large differences in shapes, use double precision build!

I tried using double-precision, and it didn't seem to help.

The history of m_cachedSeparatingAxis is interesting. During the calculation that leads to (-0.081573 1.036597 -0.081604), it gets set twice:

11: (-0.16318321, 0.036596943, -2.5927548E-7)
index = 1, impulse = 0.147091
939    (499.836823 1.036597 -500.000000)
939    (-0.081573 1.036597 -0.081604)
974    (-0.081573 1.036597 -0.081604)

The code iterates, in a loop extending from line 839 to line 969, with 8 possible exits from the loop. (There oughta be a law!) I wonder which exit is taken here...

The main exits seem to be the first 2, which I've designated A and B. The one leading to (-0.081573 1.036597 -0.081604) is exit B:

//exit 0: the new point is already in the simplex, or we didn't come any closer
                if (m_simplexSolver->inSimplex(w))
                {
                    m_degenerateSimplex = 1;
                    checkSimplex = true;
                    break;
                }

Exit B was also taken on the previous timestep, only there the loop went around 3 times:

939    (498.187439 1.035746 -500.000000)
939    (-0.906296 1.035746 -0.906281)
939    (0.000000 1.035746 -0.000031)
exit B

Perhaps the loop terminated too soon.

stephengold commented 2 years ago

Despite comments about the Johnson algorithm and alternative approaches, the only implementation of btSimplexSolverInterface is btVoronoiSimplexSolver.

In btVoronoiSimplexSolver::inSimplex(), BT_USE_EQUAL_VERTEX_THRESHOLD is defined, and m_equalVertexThreshold = 9.99999975e-05.

Adding "" to "btVoronoiSimplexSolver.cpp" broke raycasting, so I stepped through in GDB.

On the 2nd iteration of step 11, GDB sees:

(gdb) p numverts
$6 = 2
(gdb) p m_simplexVectorW[0]
$7 = {m_floats = {499.836823, 1.03659701, -500, 0}}
(gdb) p m_simplexVectorW[1]
$8 = {m_floats = {-500.163208, 1.03659701, 500, 0}}
(gdb) p w
$9 = (const btVector3 &) @0x7f944de45f10: {m_floats = {-500.163208, 
    1.03659701, 500, 0}}

The function returns true, and exit B is taken, with m_cachedSeparatingAxis = (-0.0815734863, 1.03659701, -0.0816040039). Given the exact match for w, a lower threshold probably won't prolong the loop.

GJK seems to believe that the contact point lies exactly on the edge of the mesh triangle, so the normal is undefined. I need to study the code some more.

stephengold commented 2 years ago

I think what goes on in step 11 is, there are 2 contact points calculated, one for each triangle, and the SI solver sets up a contact constraint for each point. But only one of the constraints is really valid. For the triangle with index=0, m_distance1=0.003001 and the normal is inaccurate. For the triangle with index=1, m_distance1=-0.003403 and the normal is good. Maybe we can sort this out based on the separation distances, giving priority to the point with the smallest distance...

stephengold commented 2 years ago

I now understand this issue well enough to contemplate solving it. Some approaches I'm considering:

  1. modify the normal returned by GJK for triangle edges, on a per-simplex basis
  2. modify the convex-concave collision algorithm to skip contacts that lie within the margin of another triangle of the concave
  3. add a mechanism to override normals in the sequential-impulse solver, on a per-body basis

1 and 3 seem straightforward to implement, but introduce approximations and might break assumptions made elsewhere. 2 seems the most logical approach, but might prove inefficient and/or tricky to implement in practice.

stephengold commented 2 years ago

Found btInternalEdgeUtility.h which pointed me to http://code.google.com/p/bullet/issues/detail?id=27 which pointed me to https://pybullet.org/Bullet/phpBB3/viewtopic.php?f=9&t=1814

stephengold commented 2 years ago

The intent of btInternalEdgeUtility seems similar to approach 3. Unfortunately, "ConcaveDemo" no longer exists in the Bullet repo, and the links to the GoogleCode attachments are dead, so it may be difficult to find a complete example app to study. Still, the instructions for use look straightforward.

stephengold commented 2 years ago

Pursuing approach 2, we can use btTriangleShape::isInside() to detect internal contacts. It would be nice to expose isInside() in the Java API. But although it's declared as part of btPolyhedralConvexShape, isInside() isn't yet implemented for btConvexHullShape or btBU_Simplex1to4. It ought to be implemented for all collision shapes, but that's a hobby horse for another day.

stephengold commented 2 years ago

The plan is to modify btManifoldResult::addContactPoint() to perform shape-dependent contact filtering. That way we can solve similar issues with gimpact, heightfield, and compound shapes in the future.

btTriangleShape::isInside() has a couple bugs, previously not noticed because it was unused.

stephengold commented 2 years ago

Solved in Libbulletjme for MeshCollisionShape at commit a67f3bc8. I hope to confirm that there's a similar bug in HeightfieldCollisionShape and fix that too.

stephengold commented 2 years ago

There was a bug when the capsule rests exactly on the seam, in which case both contacts get filtered out. That was solved by reducing the tolerance by 0.1%, for "less aggressive" filtering.

Now there's another bug, which seems to affect sphere shapes only---it shows up in PoolDemo but not when the ball shape is changed to MultiSphere. Apparently SphereTriangleDetector::getClosestPoints() reports the contact locations differently than btGjkPairDetector::getClosestPointsNonVirtual(), causing the cue ball to fall through the table.

stephengold commented 2 years ago

An easy solution to the latest bug would be to avoid SphereTriangleDetector by commenting out 4 return statements in "btDefaultCollisionConfiguration.cpp". However, that might be inefficient.

stephengold commented 2 years ago

Apparently SphereTriangleDetector doesn't take the triangle's margin into account.

stephengold commented 2 years ago

The SphereTriangleDetector bug was corrected in Libbulletjme v12.4.1.

I have concerns about the performance impact of contact filtering and also the likelihood that these changes to Bullet might expose more bugs, or introduce new ones. Due to these concerns, the fix will appear in a test release (Minie v4.5.0-test1) before being put into production.

stephengold commented 2 years ago

Consider this issue solved with the release of Minie v4.5.0-test1. Future work may include adding contact filtering to compound and GImpact shapes.

stephengold commented 2 years ago

I added contact filtering for Gimpact shapes to Libbulletjme at af2752b033a34a76b38336b173a3ae12858af04b. That fix should be included in the stable release of Minie v4.5 . Contact filtering for compound shapes will require a lot of work, and I don't foresee it happening in 2021.

stephengold commented 5 months ago

Another issue related to contact filtering: #40