A physics engine is a loop. Every rigid-body engine, from Bullet to MuJoCo to PhysX, advances the world by repeating three steps at a fixed timestep $h$ :

Collide. Find which bodies touch, and where. Output: contact points, normals, penetration depths.
Solve. Compute the forces or impulses that stop interpenetration and respect joints and friction. Output: corrected velocities.
Integrate. Advance positions and velocities forward by $h$ .

Step through one tick, one phase at a time:

One simulation tickinteractive

statey = 2.20v = 0.00

1 collide() ← next

2 solve()

3 integrate()

press Step

Step through one tick to see the numbers.

Collide detects the floor contact, solve reflects the velocity (restitution 0.72) and depenetrates, integrate advances under gravity. Then it repeats.

Where the difficulty lives

The integrate step is the easy one: given a velocity, move. Almost the entire research literature is about the other two, and about one tension between them.

Collision detection has to be fast (it runs on every pair, every step) and robust (a missed contact means bodies pass through each other). That is a computational-geometry problem: broad-phase culling, then narrow-phase tests like GJK. Part 4.

Solving is where the physics and the numerics fight. Contacts are inequality constraints (a body can push but not pull), friction is a cone constraint, and both switch on and off as bodies touch and separate. Writing that down honestly gives a complementarity problem; making it solvable in a few milliseconds is why engines differ so much. Parts 5 through 7.

The map onto this series

Everything ahead is a piece of this loop:

Parts 1 to 2 set up the state the loop advances: how you write a configuration (minimal vs maximal coordinates) and the equations of motion $M(q)\ddot q + C\dot q + g = \tau$ that integrate obeys.
Part 3 is the integrator itself, and why the obvious choice leaks energy.
Part 4 is collide.
Parts 5 to 7 are solve: complementarity and the friction cone, the iterative solvers (PGS, TGS) that crank through it, and the position-based shortcut (XPBD) that reframes the whole thing.
Parts 8 to 9 are what happens when you run thousands of these loops in parallel on a GPU, and when you ask for the derivative of the whole loop.

If you have not yet, the engine landscape is the companion map: same loop, nine different sets of choices for steps 2 and 3.