Hi all,

I’m currently working on a physics implementation for an icosahedron shape in a project I’m developing. The physics are applied when one of the vertices collides with the ground, and it works well if a single vertex is hitting the ground.

However, I’m encountering an issue when multiple vertices hit the ground at the same time. In this scenario, the physics stops applying any forces to the object. The first initial contact is fine, but when multiple vertices are on the floor, the system doesn’t seem to know what to do and the object sinks into the floor.

I’ve been trying to troubleshoot this issue but haven’t had much luck so far. I would greatly appreciate any insights or suggestions on how to handle this situation. Has anyone encountered a similar issue before, or does anyone have any ideas on how to prevent the object from sinking when multiple vertices are in contact with the ground?

vec3 Engine::DetectCollision(Surface& object, float height)
{
    for(int i = 0; i < object.vertices.size(); i++)
    {
        vec3 vert = object.vertices[i];
        vec4 vertex = vec4(vert.x, vert.y, vert.z, 1.0f);

        Matrix4 rotateX = Matrix4::RotateX(angularVelocity.x);
        Matrix4 rotateY = Matrix4::RotateY(angularVelocity.y);
        Matrix4 rotateZ = Matrix4::RotateZ(angularVelocity.z);
        Matrix4 rotation = Matrix4::RotateX(-90.0f) * rotateX * rotateY * rotateZ;
        Matrix4 translate = Matrix4::Translate(object_translation);
        vertex = translate * rotation * vertex;
        if(vertex.z <= height)
        {
            isHittingFloor = true;
            return vertex.Point();            
        }
    }
}

bool Engine::isCollidingFloor()
{
  return isHittingFloor;
}

void Engine::Physics()
{
    velocity = velocity + gravity * frameTime;
    object_translation = object_translation + velocity * frameTime;
    float height = rollingLandModel.getHeight(object_translation.x, object_translation.y);

    if(isStatic)
    {
        object_translation.z = height + collisionRadius;

    } else 
    {
        // Calculate the position of the center of mass after transformations
        center_of_mass_update = vec4(center_of_mass.x, center_of_mass.y, center_of_mass.z, 1.0f);
        Matrix4 translate = Matrix4::Translate(vec3(object_translation.x, object_translation.y, object_translation.z));
        center_of_mass_update = translate * center_of_mass_update;
        
        normal = rollingLandModel.getNormal(object_translation.x, object_translation.y);
        normal = normal.unit();
        vec3 point = DetectCollision(dodecahedronMesh, height); // find colliding vertex with ground
        if(isCollidingFloor())
        {
            ApplyLinearImpulse(normal, point, height, 0);
            isHittingFloor = false;
        }
    }
}

// Applying Linear impulse and angular velocity 
void Engine:ApplyLinearImpulse(const vec3& normal, vec3&vert, float height, int c)
{
	float relativeVelo = velocity.dot(normal);

	float e = coefficient_of_restitution;

	float invMass1 = 1 / mass;
	float invMass2 = 1 / groundMass;
	float invMassSum = invMass1 + invMass2;

	// j is the magnitude of the impulse
	// j is calculated per contact, we should divide j per contact
	auto j = -(1 + e) * relativeVelo / (invMassSum);

	// multiply the collision normal by the magnitude of the impulse and apply
	// resulting vector to the velocity of the object. this is how we apply linear impulse
	vec3 impulse = normal * j;

	vec3 updateVelocity = velocity + impulse * invMass1;
	if(updateVelocity.length() <= 0.1)
	{
		velocity = vec3(0,0,0);
		isStatic = true;
		return;
	}

	velocity = velocity + impulse * invMass1;
	
	// After applying linear impulse, we must apply friction.
	// to apply friction we find a vector tangent to the collision normal
	float frictionCoefficient = 0.1; 
	vec3 tangent = velocity - normal * velocity.dot(normal);
	tangent = tangent.unit();

	vec3 frictionImpulse = -tangent * frictionCoefficient * j;
	velocity = velocity + frictionImpulse / mass;

	vec3 radius = vert - center_of_mass_update.Point();
	vec3 torque = radius.cross(impulse);

	Matrix3 inertialTensor = isocahedronMesh.InertialTensor();

	vec3 angularImpulse = inertialTensor.Inverse() * torque;

	angularVelocity = angularVelocity + angularImpulse;	
}

snoken said:
Matrix4 rotateX = Matrix4::RotateX(angularVelocity.x);

Matrix4 rotateY = Matrix4::RotateY(angularVelocity.y);

Matrix4 rotateZ = Matrix4::RotateZ(angularVelocity.z);

Matrix4 rotation = Matrix4::RotateX(-90.0f) * rotateX * rotateY * rotateZ;

Omg, You use Euler angles to convert angular velocity to a 3D rotation?
That's totally wrong. Real world physics does not compose rotations in some given order. It is always just one single rotation.
Fun fact: I have seen this terrible mistake even in physics book for game dev.

It's probably not the reason for your problem, since your object does not start to rotate (although it should), but you should fix this. Otherwise you only get some fake physics in the best case, but nothing realistic, because you ignore the laws of physics.
You should read a book, look at code from physics engines, etc.

When i've personally worked on rigid body physics 30 years ago, i was following those articles here:
https://chrishecker.com/Rigid_Body_Dynamics

(new post because forum hangs when posting links)

The article explains the mass properties if a rigid body, including moment of inertia. We use this to model ‘rotation resistance’ of some object. A stick is easy to rotate if we rotate it around it's longest axis, but hard to rotate if we do so around a direction perpendicular to that longest axis. MOI gives us only an approximation of this behavior, but everybody is fine with this level of approximation. Many don't even know it's actually an approximation.

However, since modeling this adds some complexity, you can decide to ignore it. Then every object behaves like a ball, having the same rotation resistance on any direction. That's not realistic, but technically it works.

The articles also go into detail about how to calculate collision response.

And it explains how to integrate angular velocity.
In short: Angular velocity is a vector quantity, describing how much we rotate around x,y and z. But it's not 3 rotation in some order - they happen simultaneously. Euler rotations are useless here, and simply the wrong concept.

The proper way to think about 3D rotations is ‘axis and angle rotation’.
Imagine some axis, a unit vector describing a direction. Our rotation happens around this direction, so the direction does not change while the abject rotates around it.
This, and only this is how you should think about 3D rotations in general. No more need for composing multiple rotations in some order.

Knowing our direction, we also want to know how much to rotate. So we add a rotation angle to our definition of rotation axis, given in radians.
With this in mind, we can understand angular velocity, which simply is axis * angle. Instead using a unit vector for axis, we now have the angle in its magnitude.
We now can describe rotations with only 3 numbers, unlike quaternion or 3x3 matrices.
And those 3 numbers even make sense for our 3 primary axis of defining the 3D space, unlike Euler angles.
We get an advantage: We can now simply add or subtract two angular velocity vectors, and it works just a swell as with adding / subtracting position coordinates. We can do basic linear algebra with 3D rotations! We could not do this with any other representation of 3D rotations. This is essential if we want to simulate physics.

Knowing this, we can now derive our method to integrate angular velocity:

vec3 angVel (0.1f, 0, 0); // means: angular velocity which rotates 0.1 radians per second around x axis in world space
vec3 axis = angVel.Unit(); // unit vector representing axis of rotation (1,0,0)
float angle = angVel.Length(); // 0.1 radians of rotation per second

quaternion rotation; rotation.FromAxisAndAngle(axis, angle * timestep); // if our timestep is 1/60, we will rotate 0.1 / 60 radians per frame.

rigidBody.orientation *= rotation; // apply the rotation to our object, aSSUMING YOU STORE ORIENTATION AS QUATERNION. Otherwise convert the rotation quaternion to matrix to multiply

// or: rigidBody.orientation = rotation * rigidBody.orientation; depending on convention

That's correct math and it should be good enough to get started.
But there is a problem: We need to normalize to get our axis unit vector from angular velocity.
And if our angular velocity is very small, we can't do the division with a close to zero magnitude.
To solve this, we can use ‘small angle approximations’, which isn't exactly correct but works well in practice, and avoids issues from clipping small velocities to zero.
But just to mention - remember there is a solution if you one day face related problems.

Reading the rest of your code i see you actually did look up proper resources.
But the way you use angular velocity in the function on top surely is wrong and part of the problem.
I also miss any use of delta time, and i don't see any integration.
So i guess you want to learn about integration first, which usually is the easiest part.
Make sure bodies move at realistic parabola trajectories under gravity in free space, and they rotate properly if you give them some initial angular velocity.
Only after this is understood and works, you're ready to continue with collision response.