Colliding a box with a plane - strange rotation

Question

I am following the Game Physics Engine Development book. I understand splitting up a collision: generating a contact first and then resolving it based on its direction, position, penetation, etc. I already have sphere - plane and sphere - sphere collisions working correctly.

I am having major problems with box - plane however. Here is a video of how it currently looks: https://imgur.com/a/YqEjCoK

Notice that the cube has some very questionable behaviour. It never settles on one side and is forever stuck laying on a single vertex or a single side.

I want to understand the first time an impact occurs between the ground and the box. Notice that although the vertices of the bottom face lie on the same plane, travel with the same velocity and hit the ground uniformly at the same time, different impulses end up being applied to each one and rotation is wrongly produced.

The book I follow finds the contact with greatest penetration, applies rotational / linear changes to it and then adjusts the penetration of all the other contacts based on the angular / linear component of the contact with greatest penetration:

  var velocityChange = [float3.zero, float3.zero]
  var rotationChange = [float3.zero, float3.zero]

  velocityIterationsUsed = 0

  while (velocityIterationsUsed < velocityIterations) {
    var max = velocityEpsilon
    var index = contacts.count

    // find contact with greatest "penetration"
    for i in 0 ..< contacts.count {
      if (contacts[i].desiredDeltaVelocity > max) {
        max = contacts[i].desiredDeltaVelocity
        index = i
      }
    }

    if (index == contacts.count) {
      break
    }

    // Match the awake state at the contact
    contacts[index].matchAwakeState()

    // Apply velocity change to the greatest penetration contact
    contacts[index].applyVelocityChange(
      velocityChange: &velocityChange,
      rotationChange: &rotationChange
    )

    // Loop all generated contacts
    for i in 0 ..< contacts.count {
      // Check each body in the contact
      for b in 0 ..< 2 {
        if (contacts[i].bodies[b] != nil) {
          // Check for a match with each body in the newly resolved contact
          for d in 0 ..< 2 {
            if (contacts[i].bodies[b] == contacts[index].bodies[d]) {
              let deltaVel = velocityChange[d] + cross(rotationChange[d], contacts[i].relativeContactPosition[b])

              print("""
                    =======================
                    contact idx \(i):
                    velocityChange: \(velocityChange[d])
                    rotationChange: \(rotationChange[d])
                    relativeContactPosition: (\(contacts[i].relativeContactPosition[b].x) \(contacts[i].relativeContactPosition[b].y) \(contacts[i].relativeContactPosition[b].z))
                    deltaVel: (\(deltaVel.x) \(deltaVel.y) \(deltaVel.z))
              """)

              // deltaVel does not make sense on first impact

            }
          }
        }
      }
    }

  }

  velocityIterationsUsed += 1

}

Notice the deltaVel variable that reconstructs each contact velocity offset based on the linear and angular components. On first box impact, as the 4 vertices of the bottom face hit the ground uniformly, at the same time and with the same velocity, deltaVel reports really strange values:

=======================
contact idx 0:
velocityChange: SIMD3<Float>(0.0, 3.1689641, 0.0)
rotationChange: SIMD3<Float>(0.26407933, 0.0, 0.2640793)
relativeContactPosition: (0.5 -0.5010185 -0.50000006)
deltaVel: (0.13230862 3.4330435 -0.13230863)
=======================
contact idx 1:
velocityChange: SIMD3<Float>(0.0, 3.1689641, 0.0)
rotationChange: SIMD3<Float>(0.26407933, 0.0, 0.2640793)
relativeContactPosition: (-0.5 -0.5010185 -0.50000006)
deltaVel: (0.13230862 3.1689641 -0.13230863)
=======================
contact idx 2:
velocityChange: SIMD3<Float>(0.0, 3.1689641, 0.0)
rotationChange: SIMD3<Float>(0.26407933, 0.0, 0.2640793)
relativeContactPosition: (0.5 -0.5010185 0.4999999)
deltaVel: (0.13230862 3.1689641 -0.13230863)
=======================
contact idx 3:
velocityChange: SIMD3<Float>(0.0, 3.1689641, 0.0)
rotationChange: SIMD3<Float>(0.26407933, 0.0, 0.2640793)
relativeContactPosition: (-0.5 -0.5010185 0.4999999)
deltaVel: (0.13230862 2.9048848 -0.13230863)

Notice the line

let deltaVel = velocityChange[d] + cross(rotationChange[d], contacts[i].relativeContactPosition[b])

And how it produces almost the same deltaVel for each vertex, with the exception of the Y component which gradually decreases with each next vertex? But why does it decrease? I did this equation for each vertex in matlab and it is indeed correct. It is my understanding that it should be equal along all vertices (as they lay on the same plane)? This way they will all bounce back with the same Y velocity and no rotation will happen.

As you can see, the only different values between the iterations are the contact[i].relativeContactPosition. As the box is with width, height and depth of 1, these values are around +/- 0.5 in each axis.

// vertex #1
let linearChange0 = float3(0, 3.1854, 0.0)
let angularChange0 = float3(0.265, 0, 0.265)
let relativeContactPos0 = float3(-0.5, -0.5265012, -0.50000006)
linearChange0 + cross(angularChange0, relativeContactPos0) = (0.139522806, 3.18540001, -0.139522806)

// vertex #2
let linearChange1 = float3(0, 3.1854, 0.0)
let angularChange1 = float3(0.265, 0, 0.265)
let relativeContactPos1 = float3(-0.5, -0.5265012, 0.4999999)
linearChange1 + cross(angularChange1, relativeContactPos1) = (0.139522806, 2.92040014, -0.139522806)

Given that the only difference between the two equation is the Z component of the relativeContactPos value, why does the resulting Y value differ? How to make them uniform?

TLDR: Each vertex in the attached above IMGUR video gets different Y velocity offset, even though they are all positioned on the same plane, traveling at the same velocities. What to do?

Answer 1

Given that the only difference between the two equation is the Z component of the relativeContactPos value, why does the resulting Y value differ? How to make them uniform?

Because you are taking a cross-product. Only the off-axis terms contribute to the resulting axis term:

That you are taking a cross product is correct -- torque is a cross product, after all -- but clearly the wrong cross product is being taken.

The values calculated for change in velocity (velocityChange) appear correct, as do the position vectors (relativeContactPos).

What we need is the impulse due to torque on each vertex which we know should take the form

Since torque is given by

And since impulse is equal to force multiplied by time, we can construct this equation for the impulse due to torque:

but since acceleration is equal to change in velocity over time, this becomes

Hence, the cross product should be

cross(contacts[i].relativeContactPosition[b], velocityChange[d])

Which will result in the completely offsetting vectors shown below (top-down view)

Note, however, that this does not handle the impulse aspect, only the torque aspect. The magnitude of the torques will be scaled by the mass and the time-delta. In this specific case of face-on collision it doesn't matter because the resulting torques cancel out, but other cases will see their impulse calculations incorrectly scaled.

Since you haven't included the necessary information needed to complete the impulse calculation, it can't be said for certain what you need to do, but in broad strokes the complete contribution to change in velocity from the torque will be calculated as:

cross(contacts[i].relativeContactPosition[b], velocityChange[d]) * deltaT * mass