I've been developing a simple physics engine as a side project for a while, but I've been stuck on an issue when incorporating large amounts of bodies. When I start using around 50+ rigid bodies in my scene, some bodies seem to gain velocity from nowhere, and get shot off to the side, and eventually all of the bodies disappear.

For collision, I'm using a simple AABB:

struct AABB
{
Math::Vector2 center;
Math::Vector2 hSize;

AABB() = default;
AABB(const Math::Vector2& _center,
const Math::Vector2& _hSize)
:
center{_center},
hSize{_hSize}
{}
};

My broadphase is a quadtree, and I'm using impulses for collision resolution:

Math::Vector2 relativeVelocity {
contact->bodies[1]->linearVelocity - contact->bodies[0]->linearVelocity
};

float normalVelocity { relativeVelocity * contact->normal };
if (normalVelocity <= 0) return;

float e { std::min(contact->bodies[0]->restitution, contact->bodies[1]->restitution) };
float totalInvMass { contact->bodies[0]->invMass + contact->bodies[1]->invMass };

float j { -(1.0f + e) * normalVelocity / totalInvMass };

Math::Vector2 impulse { contact->normal * j };
contact->bodies[0]->linearVelocity -= impulse * contact->bodies[0]->invMass;
contact->bodies[1]->linearVelocity += impulse * contact->bodies[1]->invMass;

And a position solver:

constexpr float slop { 0.01f };
constexpr float percent { 0.8f };

Math::Vector2 correction {
contact->normal * percent * std::max(contact->intersection - slop, 0.0f) /
(contact->bodies[0]->invMass + contact->bodies[1]->invMass)
};

contact->bodies[0]->position -= correction * contact->bodies[0]->invMass;
contact->bodies[1]->position += correction * contact->bodies[1]->invMass;

My physics engine step is applying gravity to bodies that are affected, integrate the bodies, collision/detection+resolution, then clear forces.

Is there any reason I would be having an issue?

@bryant5 I have experienced similar issues when acceleration is high or the step time is high, or both. In both situations, there is a big penetration between objects and the correction force is too high.

My fisrt recommendation is reducing the step time between two frames to a small value compared with the one you have now (Or alternatively reduce the gravity).

Where did you get your physics relations from ?

@Gotanod I'm currently just using the time between frames, I don't have a fixed step atm. I've tried reducing gravity before too

@bryant5 Instead of calling the physics update with the time between frames (dt). Create a loop to call the physics n times with a fraction of the time (dt/n) (Sub-stepping). This reduces issues like bullet-through-paper issue. Because it moves the physics simulation slower.
But you multiply the CPU cycles that it needs by the sub-stepping number that you have selected.
Even if you don't use it in production, it helps to identify that you have a problem with large velocities/acceleration/delta times/floating point margins…

physics.update(dt);

Replace with sub-stepping (n steps)

n=10; // for example
subdt = dt/n;
for (int i=0;i<n;++i) {
physics.update(subdt);
}

or if you want a fixed dt for physics

subdt = 0.01; // fixed
int n = (dt + remaining) / subdt; // integer division or floor
remaining = dt - n * subdt; // keep that time for next update
for (int i=0;i<n;++i) {
physics.update(subdt);
}

Reference to Unreal Engine Sub-stepping
https://docs.unrealengine.com/5.1/en-US/physics-sub-stepping-in-unreal-engine/