This applet simulates rectangular shaped objects colliding in 2 dimensions. You can set gravity, elasticity (bounciness), and damping (friction) using the slider controls. You can choose from one to six objects. The mass of the objects is initially set to 1.0, but the mass of the green object is adjustable.

Click near an object to exert a rubber band force with your mouse. With the keyboard you can control four "thrusters". The keys S,D,F,E control the blue object. The keys J,K,L,I (and also the arrow keys) control the green object. (If the keys don't work, try clicking near an object first - this ensures that key strokes are passed to the simulation instead of the buttons or scrollbars).

If elasticity is less than 1, and gravity is greater than zero, then the objects eventually settle onto the floor. Because this simulation does not implement resting contact the simulation will get 'stuck'.

To check the correctness of the simulation, you can click the "show energy" checkbox and look at the display of the pre- and post-collision momentum and energy. If damping = 0 and elasticity = 1, then the energy should not change. The momentum should not change when objects collide, however it will change when an object collides with a wall (walls are assumed to be infinitely massive, so we don't model their movement).

The bar graph shows the gravitational, rotational and translational energy. Note that if gravity is zero, then there is no gravitational potential energy.

See the section on Energy and Momentum for how these quantities are calculated.

Ignoring collisions for the moment, we develop the equations of motion for the objects. The forces on the objects are

- Gravity
- Friction
- Thrust (keyboard controls)
- Rubber band (mouse control)

For each body, there are three variables to specify position:

x = horizontal position of center of massy = vertical position of center of massθ = angle of rotation

In addition, each body has a velocity corresponding to each of these positions:

x' = v horizontal velocity of center of mass_{x}=y' = v vertical velocity of center of mass_{y}=θ' = ω = angular velocity

The equations of motion for the body involve the total force
**F** = (F_{x}, F_{y})

_{x} = m x''

F_{y} = m y''

τ = I θ''

where **F**

We now build up the equations of motion by adding in one force at a time, beginning with the thrust force. Let **T** = (T_{x}, T_{y})

_{x}

m y'' = T_{y}

In these equations it doesn't matter where on the body the thrust force is applied. The point

**Moment of inertia** is the equivalent of mass for rotational physics. It measures how difficult it is to rotate a body about a given point. Since our rectangle bodies rotate freely about their center of mass, we use center of mass as the point for calculating moment of inertia. From a physics textbook you can find the equation for the moment of inertia about the center of a thin rectangular plate. It is given by

^{2} + height^{2})/12

Let **R** = (R_{x}, R_{y}) =**torque** at the point **R** × **T** = R_{x} T_{y} - R_{y} T_{x}

**R** × **T**

Actually torque is a vector, but since we are working in 2 dimensions we know that torque is always perpendicular to the plane and so we aren't using vector notation for it. The true vector cross product results in a vector. So the above cross product corresponds to the following 3 dimensional calculation:

_{x}, R_{y}, 0) × (T_{x}, T_{y}, 0) = (0, 0, R_{x} T_{y} - R_{y} T_{x})

You can see that the result is always perpendicular to the plane, with zero x and y components. The general vector cross product of two vectors is defined as:

The rubber band force in the simulation operates along the vector from point **L****B** = (B_{x}, B_{y}) = s **L**

_{x} + B_{x}

m y'' = T_{y} + B_{y}

I θ'' = **R** × **T** + **R** × **B**

Let

_{y} + B_{y} - m g

Because gravity works on all points of the body, no torque is generated by gravity.

Damping (friction) causes a force opposite to the direction of motion. The faster you go, the more friction resists your motion. So the magnitude of the damping force is proportional to the velocity. Let

_{x} + B_{x} - k x' (eqn 1)

m y'' = T_{y} + B_{y} - m g - k y' (eqn 2)

I θ'' = **R** × **T** + **R** × **B** - k θ' (eqn 3)

In a more realistic simulation, there may be different damping constants for rotational versus translational motion. But here, we use the same constant for both.

Equations (1-3) are the equations of motion for one of our rectangular bodies. To use the Runge-Kutta method for solving a set of differential equations we need to convert these 3 second order equations into 6 first order equations. Define the velocity variables

v horizontal velocity_{x}= x' =v vertical velocity_{y}= y' =ω = θ' = angular velocity

Then we have the 6 first order differential equations

_{x}

y' = v_{y}

θ' = ω

v_{x}' = (T_{x} + B_{x} - k v_{x})/m

v_{y}' = (T_{y} + B_{y} - m g - k v_{y})/m

ω' = (**R** × **T** + **R** × **B** - k ω)/I

These equations are now in the form needed for solving numerically with the Runge-Kutta method. Each rectangle body will have its own set of 6 variables (3 positions _{x},v_{y},ω

This section explains how the energy and momentum of the objects are calculated. See the description of the energy bar graph for how to observe these quantities in the simulation.

If there is no loss of energy to friction (damping

A collision between objects should not change the angular momentum. However, a collision with a wall will not preserve angular momentum because the super massive walls are not included in the calculations of the momentum of the objects.

Gravitational energy is given by

Translational energy is **v**^{2}**v**

Rotational energy is ^{2}

Linear momentum in the horizontal direction is given by _{x}_{y}

Angular momentum is measured with regard to a particular point in space, for example the origin. It is given by:

**k** + m **r** × **v**

where **r****k****k****r** × **v**

At each step in the simulation, we check to see if there is a collision. The bodies can collide with each other or with a boundary wall. For the rectangular shapes we are using it is simple geometry to determine if a collision has occurred by checking if any vertex is within a wall or foreign body. If so, we use a binary search to back up to an earlier time just shortly before the collision occurred. We then make the approximation that the collision takes place at this exact time, and calculate the resulting velocities as described below. The Colliding Blocks simulation further describes these aspects of collision handling.

Handling collisions is the most challenging part of this simulation. The explanation here is fairly condensed, so you may want to read some other descriptions as well.

- Chris Hecker's Rigid Body Dynamics Information. See especially his article for Game Developer Magazine on Collision Response. Easy to read, with many other references listed on the web site. Note that some of Mr. Hecker's PDF documents do not print their math correctly.
- Physically Based Modeling is a set of course notes from Siggraph '97 by Andrew Witkin and David Baraff. See especially the sections on Rigid Body Dynamics. Fairly complete and rigorous explanations, but still quite readable. Covers the 3 dimensional case and also resting contact forces.
- Physics For Game Developers by David M. Bourg (O'Reilly Press) has a section on linear and angular impulse on page 95, and resting contact forces on page 258. This book has good explanations, but often relies on hard-to-follow programming code to illustrate the physics.
- The textbook Engineering Mechanics Dynamics by Robert W. Soutas-Little and Daniel J. Inman has relevant sections on Eccentric Impact.

Suppose a vertex on body A is colliding into an edge of body B at the point P. Define the following variables

m mass of bodies A, B_{a}, m_{b}= distance vector from center of mass of body A to point P**r**_{ap}= distance vector from center of mass of body B to point P**r**_{bp}=ω initial pre-collision angular velocity of bodies A, B_{a1}, ω_{b1}=ω final post-collision angular velocity of bodies A, B_{a2}, ω_{b2}= initial pre-collision velocities of center of mass bodies A, B**v**_{a1},**v**_{b1}= final post-collision velocities of center of mass bodies A, B**v**_{a2},**v**_{b2}= initial pre-collision velocity of impact point on body A**v**_{ap1}= initial pre-collision velocity of impact point on body B**v**_{bp1}= normal (perpendicular) vector to edge of body B**n**=

We now use a standard formula for the velocity of an arbitrary point on a rotating and translating rigid body to get the pre-collision velocities of the points of collision (which is the point

**v**_{ap1} = **v**_{a1} + ω_{a1} × **r**_{ap}

**v**_{bp1} = **v**_{b1} + ω_{b1} × **r**_{bp}

Similarly we have the final post-collision velocities **v**_{ap2}**v**_{bp2}

**v**_{ap2} = **v**_{a2} + ω_{a2} × **r**_{ap}

**v**_{bp2} = **v**_{b2} + ω_{b2} × **r**_{bp}

Here we are regarding the angular velocity as a 3 dimensional vector perpendicular to the plane, so that the cross product is calculated as

**r** = (0, 0, ω) × (r_{x}, r_{y}, 0) = (-ω r_{x}, ω r_{y}, 0)

Now we can find an expression for the velocity with which the colliding points are approching each other. We call this the relative velocity. Let **v**_{ab1}**v**_{ab2}

**v**_{ab1} = **v**_{ap1} - **v**_{bp1}

**v**_{ab2} = **v**_{ap2} - **v**_{bp2}

Using the formulas given above for velocity of a point on a rigid body we can expand these to

**v**_{ab1} = **v**_{a1} + ω_{a1} × **r**_{ap} - **v**_{b1} - ω_{b1} × **r**_{bp} (eqn 4)

**v**_{ab2} = **v**_{a2} + ω_{a2} × **r**_{ap} - **v**_{b2} - ω_{b2} × **r**_{bp} (eqn 5)

Let the vector **n****n****n**

**v**_{ab1} ∙ **n** = **v**_{ab1x} **n**_{x} + **v**_{ab1y} **n**_{y}

Note that for a collision to occur this relative normal velocity must be negative (that is, the objects must be approaching each other).
Let

**v**_{ab2} ∙ **n** = -e **v**_{ab1} ∙ **n** (eqn 6)

This says that the velocity at which the objects fly apart is proportional to the velocity with which they were coming together. The proportionality factor is the elasticity

To resolve the collision, we will use the concept of an *impulse*. An impulse is the change in momentum of an object when a large force is applied over a very brief period of time. We imagine that during the collision there is a very large force acting for a very brief period of time. If you integrate (sum) that force over that brief time, you get the impulse.

Why do we need this strange concept of an impulse? Why not just use the familiar concept of force as in

This is far beyond what our simple simulation can do. Luckily, we can assume that the collision happens so quickly that the position and orientation of the bodies do not change during the collision. Instead, all that changes is the velocities of the bodies. Since a change in velocity is a change in momentum (remember momentum = velocity times mass) we have the concept of an impulse.

We are assuming no friction for our collision, so the only force during the collision is in the direction perpendicular to the edge, which is given by the vector **n****n****n****n**

**v**_{a2} = **v**_{a1} + j **n** / m_{a} (eqn 7)

**v**_{b2} = **v**_{b1} - j **n** / m_{b} (eqn 8)

The change in angular momentum of body A from the impulse **n****r**_{ap} × j **n**

_{a2} = ω_{a1} + (**r**_{ap} × j **n**) / I_{a} (eqn 9)

ω_{b2} = ω_{b1} - (**r**_{bp} × j **n**) / I_{b} (eqn 10)

Now we can put all these various equations together and solve for the impulse parameter

We start with equation (6), and then expand using our definition of relative velocity in equations (4-5).

**v**_{ab2} ∙ **n** = -e **v**_{ab1} ∙ **n**

(**v**_{ap2} - **v**_{bp2}) ∙ **n** = -e **v**_{ab1} ∙ **n**

[**v**_{a2} + ω_{a2} × **r**_{ap} - **v**_{b2} - ω_{b2} × **r**_{bp}] ∙ **n** = -e **v**_{ab1} ∙ **n**

Let's expand the left hand side, using the impulse relationships in equations (7-10).

**v**_{a1} + j **n** / m_{a}) + (ω_{a1} + (**r**_{ap} × j **n**) / I_{a}) × **r**_{ap} **v**_{b1} - j **n** / m_{b}) _{b1} - (**r**_{bp} × j **n**) / I_{b}) × **r**_{bp}] ∙ **n** **v**_{ab1} ∙ **n**

We recognize that the left hand side contains the quantity **v**_{ab1} ∙ **n**

**n** / m_{a} + (**r**_{ap} × j **n**) × **r**_{ap} / I_{a} **n** / m_{b} **r**_{bp} × j **n**) ×**r**_{bp} / I_{b}] ∙ **n** **v**_{ab1} ∙ **n**

Note that **n****n** ∙ **n** = 1

**A** × **B**) ∙ **C** = (**B** × **C**) ∙ **A**

to derive the following identity

**A** × **B**) × **A** ∙ **B** = (**A** × **B**) ∙ (**A** × **B**) = (**A** × **B**)^{2}

(Here squaring a vector means taking its dot product with itself.) We can then simplify to

_{a} + 1 / m_{b} + (**r**_{ap} × **n**)^{2} / I_{a} **r**_{bp} × **n**)^{2} / I_{b})**v**_{ab1} ∙ **n**

Dividing leads to our final expression for

-(1 + e) |

We can now calculate

We can use this same expression for calculating collisions with a wall by assuming that the mass of the wall is infinite. So let _{b} → ∞_{b} → ∞

-(1 + e) |

The above analysis only handles the case of a single corner and edge impact. There are several other simultaneous multiple impact cases, such as:

- Two corners of object A impact a single wall
- Two corners of object A impact object B
- Two corners of object A impact different walls
- Two corners of object A impact different objects

If there are simultaneous impacts among unrelated objects, these can be easily handled separately. For example, suppose objects A and B collide at the same time that objects C and D collide, then we can deal with each collision separately because they don't affect each other.

For the case of two adjacent corners of object A impacting either a single wall or another object's edge, we change the impact point to be the midpoint between the two corners and then handle it as other collisions.

The more complicated cases with corners of object A impacting different walls or multiple objects is not handled correctly by this simulation. What the code does is handle each impact separately, which will probably be wrong because these collisions are not independent.

You will notice that certain settings of gravity and elasticity cause the simulation to eventually get stuck. Specifically, this happens when gravity

The right thing to do at this point (which this simulation does not do) is recognize that the object is in *resting contact* with the floor. This would change the equations of motion for the object to include a contact force from the floor, which prevents the penetration of the object. Some of the references cited earlier have sections on calculating resting contact forces.