## Billiard simulation–part 8: double collision

June 30, 2013 at 4:47 AM

As explained previously, the current algorithm does not handle collisions of multiple balls at the same time. I found a solution for a double collision, and I am pretty sure that there is no single solution for more complex cases using only conservation laws and the coefficient of restitution formula – but I still have to prove it. Of course there would be a solution modeling the actual physics of the collisions – a problem for another day.

### Double ball collision

..or more formally: elastic collision of three spheres in two dimensions

Ball 0 – with velocity $$(v_0,w_0)$$, hits at the same time ball 1 and ball 2 –  that have velocities $$(v_1, w_1)$$ and $$(v_2, w_2)$$ respectively.

Projecting the velocity $$(v_0,w_0)$$ onto the lines connecting the center of ball 0 with the centers of balls 1 and 2 gives the velocity components $$q_{01}$$ and $$q_{02}$$. Projecting the velocity $$(v_1,w_1)$$ onto the line connecting the center of ball 0 with the center of ball 1 and onto its perpendicular through the center of ball 1 gives the velocity components $$q_1$$ and $$p_1$$. Doing the same thing for ball 2 gives the velocity components $$q_2$$ and $$p_2$$.

Here is the drawing:

As usual, balls are not spinning and there is no ball-ball friction, so the forces between balls are only along the lines connecting the center of ball 0 with the centers of balls 1 and 2. This means that the velocity components perpendicular to those lines are the same before and after the collision:

$$p’_1 = p_1$$

$$p’_2 = o_2$$

(indicating with a prime the velocities after the collision), whereas the velocity components along those lines are related by the coefficient of restitution formula:

$$q’_{01} – q’_1 = c(q_1 – q_{01})$$

$$q’_{02} – q’_2 = c(q_2 – q_{02})$$

Let’s not use trigonometry this time. The projection of a velocity onto the line connecting the center of two ball is the scalar product of the velocity with the unit vector pointing from one center to the other. For balls 0 and 1 this unit vector is:

$$\left( \frac{x_{01}}{\sqrt{x_{01}^2 + y_{01}^2}}, \frac{y_{01}}{\sqrt{x_{01}^2 + y_{01}^2}} \right)$$

where $$x_{01} = x_1 – x_0$$ and $$y_{01}=y_1-y_0$$, hence:

$$q_1 = \frac{v_1x_{01} + w_1y_{01}}{\sqrt{x_{01}^2 + y_{01}^2}}$$

$$q_01 = \frac{v_0x_{01} + w_0y_{01}}{\sqrt{x_{01}^2 + y_{01}^2}}$$

Similarly, the projection on the perpendicular of that line is the scalar product of the velocity with the perpendicular unit vector:

$$\left( -\frac{y_{01}}{\sqrt{x_{01}^2 + y_{01}^2}} , \frac{x_{01}}{\sqrt{x_{01}^2 + y_{01}^2}} \right)$$

hence:

$$p_1 = \frac{-v_1y_{01} + w_1x_{01}}{\sqrt{x_{01}^2 + y_{01}^2}}$$

Replacing these values in the four initial equations gives:

$$–v’_1y_{01} + w’_1x_{01} = -v_1y_{01} + w_1x_{01}$$

$$–v’_2y_{02} + w’_2x_{02} = -v_2y_{02} + w_2x_{02}$$

$$v’_0x_{01} + w’_0y_{01} – v’_1x_{01} – w’_1y_{01} = c(v_1x_{01} + w_1y_{01} - v_0x_{01} - w_0y_{01})$$

$$v’_0x_{02} + w’_0y_{02} – v’_1x_{02} – w’_2y_{02} = c(v_2x_{02} + w_2y_{02} - v_0x_{02} - w_0y_{02})$$

The conservation of momentum gives two additional equations:

$$m_0v’_0 + m_1v’_1 + m_2v’_2 = m_0v_0 + m_1v_1 + m_2v_2$$

$$m_0w’_0 + m_1w’_1 + m_2w’_2 = m_0w_0 + m_1w_1 + m_2w_2$$

- where $$m_0$$, $$m_1$$ and $$m_2$$ are the masses of the balls, resulting in a system of 6 linear equations for the 6 unknowns $$v’_0, w’_0, v’_1, w’_1, v’_2, w’_2$$. Introducing new unknowns $$u_j$$ defined as:

$$u_0 = v’_0 – v_0,\quad u_1 = w’_0 – w_0,\quad u_2 = v’_1 – v_1,\quad u_3 = w’_1 – w_1,\quad u_4 = v’_2 – v_2,\quad u_5 = w’_2 – w_2$$

the system becomes a bit simpler:

$$x_{01} u_3 – y_{01} u_2 = 0$$

$$x_{02} u_5 – y_{02} u_4 = 0$$

$$x_{01} (u_0-u_2) + y_{01}(u_1-u_3) = B_{01}$$

$$x_{02} (u_0-u_4) + y_{02}(u_1-u_5) = B_{02}$$

$$m_0u_0 + m_1u_2 + m_2u_4 = 0$$

$$m_0u_1 + m_1u_3 + m_2u_5 = 0$$

with

$$B_{01} = (1+c)[x_{01}(v_1 – v_0) + y_{01}(w_1-w_0)]$$

$$B_{02} = (1+c)[x_{02}(v_2 – v_0) + y_{02}(w_2-w_0)]$$

After quite a lot of substitutions and simplification the result is

$$u_0 = v’_0 – v_0 = –\frac{1}{\Delta}[B_{01}m_1(x_{02}m_2P – x_{01}(m_0+m_2)S_{02}) + B_{02}m_2(x_{01}m_1P-x_{02}(m_0+m_1)S_{01})]$$

$$u_1= w’_0 – w_0 = –\frac{1}{\Delta}[B_{01}m_1(y_{02}m_2P – y_{01}(m_0+m_2)S_{02}) + B_{02}m_2(y_{01}m_1P-y_{02}(m_0+m_1)S_{01})]$$

$$u_2 = v’_1 – v_1 = \frac{x_{01}m_0}{\Delta}[B_{02}m_2P-B_{01}(m_0+m_2)S_{02})]$$

$$u_3= w’_1 – w_1 = \frac{y_{01}m_0}{\Delta}[B_{02}m_2P-B_{01}(m_0+m_2)S_{02})]$$

$$u_4 = v’_2 – v_2 = \frac{x_{02}m_0}{\Delta}[B_{01}m_1P-B_{02}(m_0+m_1)S_{01})]$$

$$u_5 = w’_2 – w_2 = \frac{y_{02}m_0}{\Delta}[B_{01}m_1P-B_{02}(m_0+m_1)S_{01})]$$

with

$$P = x_{01}x_{02}+y_{01}y_{02}$$

$$S_{01} = x_{01}^2 + y_{01}^2$$

$$S_{02} = x_{02}^2 + y_{02}^2$$

$$\Delta = (m_0 + m_1)(m_0 + m_2) S_{01}S_{02}-m_1m_2P^2$$

The corresponding code is:

/**
* Updates the velocities of this ball and two other balls after a double collisions
* with both balls at the exact same time
* The coordinate of the balls must be at the collision point.
* @param otherBall1 second colliding ball
* @param otherBall2 third colliding ball
* @param restitution coefficient of restitution for a ball-ball collision
*/
collide2(otherBall1: Ball, otherBall2: Ball, restitution: number) {
var x01 = otherBall1.x - this.x;
var y01 = otherBall1.y - this.y;
var x02 = otherBall2.x - this.x;
var y02 = otherBall2.y - this.y;
var p = x01 * x02 + y01 * y02;
var s01 = x01 * x01 + y01 * y01;
var s02 = x02 * x02 + y02 * y02;
var delta = (m0 + m1) * (m0 + m2) * s01 * s02 - m1 * m2 * p * p;
var v01 = otherBall1.v - this.v;
var w01 = otherBall1.w - this.w;
var b01 = (1 + restitution) * (x01 * v01 + y01 * w01);
var v02 = otherBall2.v - this.v;
var w02 = otherBall2.w - this.w;
var b02 = (1 + restitution) * (x02 * v02 + y02 * w02);
this.v = -(b01 * m1 * (x02 * m2 * p - x01 * (m0 + m2) * s02) + b02 * m2 * (x01 * m1 * p - x02 * (m0 + m1) * s01)) / delta + this.v;
this.w = -(b01 * m1 * (y02 * m2 * p - y01 * (m0 + m2) * s02) + b02 * m2 * (y01 * m1 * p - y02 * (m0 + m1) * s01)) / delta + this.w;
var r1 = m0 * (b02 * m2 * p - b01 * (m0 + m2) * s02) / delta;
otherBall1.v = x01 * r1 + otherBall1.v;
otherBall1.w = y01 * r1 + otherBall1.w;
var r2 = m0 * (b01 * m1 * p - b02 * (m0 + m1) * s01) / delta;
otherBall2.v = x02 * r2 + otherBall2.v;
otherBall2.w = y02 * r2 + otherBall2.w;
} // collide2

As a verification, when $$m_2 = 0$$ the first four equations become:

$$u_0 = v’_0 – v_0 = x_{01} \frac{m_1}{m_0+m_1} \frac{B_{01}}{S_{01}}$$

$$u_1 = w’_0 – w_0 = y_{01} \frac{m_1}{m_0+m_1} \frac{B_{01}}{S_{01}}$$

$$u_2 = v’_1 – v_1 = -x_{01} \frac{m_0}{m_0+m_1} \frac{B_{01}}{S_{01}}$$

$$u_3 = w’_1 – w_1 = -y_{01} \frac{m_0}{m_0+m_1} \frac{B_{01}}{S_{01}}$$

if $$\alpha$$ is the angle between the vector from ball 0 to ball 1 and the horizontal, then:

$$\cos\alpha = \frac{x_{01}}{\sqrt{x_{01}^2 + y_{01}^2}}$$

$$\sin\alpha = \frac{y_{01}}{\sqrt{x_{01}^2 + y_{01}^2}}$$

hence

$$\frac{B_{01}}{S_{01}} = -\frac{1+c}{\sqrt{x_{01}^2 + y_{01}^2}}[\cos\alpha(v_0 – v_1) + \sin\alpha(w_0-w_1)]$$

and the result is the same as in the two balls case.

Posted in: Physics | Programming