April 2017 – Dylan's AMCS Blog

Notes:

As a final project for one of my Computer Science courses I made a Google App Engine game of billiards where you can interact with the mathematical derivations from this post. It requires two players, so you’ll need to navigate to the web page in two separate browser tabs if you don’t have an opponent. Visit http://our-mediator-166519.appspot.com/ to play.
Additionally, all of the images from this post were rendered using some of the components I developed for this project. You can view the full source on my Github repository.

I’ve been allowed some creative freedom with a final project for one of my computer science courses. Naturally, I’m going to make a video game. This has led me to think lately about some of the ensuing technical hurdles and mathematical tools I’ll need to develop in order to get the project underway.

Matrix determinant and parallelogram area

P5Image1

I’ll begin by defining a vector operation, $\circ$ , where

$a \circ b = a_x b_y - a_y b_x$ .

The operation is just a 2×2 matrix determinant, but without getting too technical, I’ll just think of it as an operation on two vectors that gives a number, the absolute value of which is the area of the parallelogram that they subtend.

P5Image2

Note that the operation is not commutative, it is anti-commutative.

$\begin{aligned} a \circ b &= a_x b_y - a_y b_x \\ &= -(a_y b_x - a_x b_y) \\ &= -(b_x a_y - b_y a_x) \\ &= - (b \circ a) \end{aligned}$

This is a useful property of the operation. if $a \circ b$ is positive, then I know the area of the parallelogram that they subtend, but also that vector b is oriented counter-clockwise with respect to vector a. On the other hand, if $a \circ b$ is negative, then I know that vector b is oriented clockwise with respect to vector a.

Also useful to know, the operation is left and right distributive over addition/subtraction:

$\begin{aligned} a \circ (b + c) &= a_x (b_y + c_y) - a_y (b_x + c_x) \\ &= a_x b_y + a_x c_y- a_y b_x - a_y c_x \\ &= a_x b_y - a_y b_x + a_x c_y - a_y c_x\\ &= a \circ b + a \circ c \\ \end{aligned}$

also,

$(a+b) \circ c = a \circ c + b \circ c$

I will be using a $\circ$ b notation throughout this post to save typing. In addition to this, I will use more familiar notation for other common operations, such as:

Dot Product

$a \cdot b = a_x b_x + a_y b_y$

and Euclidean norm

$|a| = \sqrt{a_x ^ 2 + a_y ^ 2}$

What is the height of the blue vector over the red vector?

P5ImageX1

To find the length of the yellow line, I’ll begin by finding the area of the parallelogram subtended by the blue and red vectors.

$\begin{aligned} A_p &= a_x b_y - a_y b_x \\ &= a \circ b \end{aligned}$

P5ImageX2

Half of this area gives the area of the triangle defined by the blue and red vectors.

$\begin{aligned} A_t &= \frac{a_x b_y - a_y b_x}{2} \\ &= \frac{a \circ b}{2} \end{aligned}$

P5ImageX3

The length of the yellow line is the height of this triangle, with base equal to the length of the red vector. I can find the height of the yellow line with the familiar formula,

$\begin{aligned} A_t &= \frac{(base)(height)}{2}\\ height &= \frac{2 A_t}{base}. \end{aligned}$

Substituting for $A_t$ and the base of the triangle gives

$\begin{aligned} h &= \frac{2[\frac{a_x b_y - a_y b_x}{2}]}{[\sqrt{b_x^2 + b_y^2}]}\\ &= \frac{a_x b_y - a_y b_x}{\sqrt{b_x^2 + b_y^2}} \\ &= \frac{a \circ b}{|b|} \end{aligned}$

Note that this operation is related to the formula for vector projection,

$\begin{aligned} aProjB &= \frac{a_x b_x + a_y b_y}{\sqrt{b_x^2 + b_y^2}}\\ &= \frac{a \cdot b}{|b|} \end{aligned}$

When paired together, the two operations satisfy the Pythagorean theorem:

$\begin{aligned} |a|^2 &= (\frac{a \cdot b}{|b|})^2 + (\frac{a \circ b}{|b|})^2\\ &= (\frac{a_x b_x + a_y b_y}{\sqrt{b_x^2 + b_y^2}})^2 + (\frac{a_x b_y - a_y b_x}{\sqrt{b_x^2 + b_y^2}})^2\\ &= \frac{a_x^2 b_x^2 + 2a_x b_x a_y b_y + a_y^2 b_y^2 + a_x^2 b_y^2 - 2a_x b_y a_y b_x+ a_y^2 b_x^2}{b_x^2 + b_y^2}\\ &= \frac{a_x^2 b_x^2 + a_y^2 b_y^2 + a_x^2 b_y^2 + a_y^2 b_x^2}{b_x^2 + b_y^2}\\ &= \frac{a_x^2(b_x^2 + b_y^2) + a_y^2(b_x^2 + b_y^2)}{b_x^2 + b_y^2}\\ &= \frac{(a_x^2+a_y^2)(b_x^2 + b_y^2)}{b_x^2 + b_y^2}\\ &= a_x^2+a_y^2\\ &= |a|^2 \end{aligned}$

Together they give the base and height, respectively, of the generalized right triangle formed by two arbitrary vectors.

$\begin{aligned} base &= \frac{a \cdot b}{|b|} & height &= \frac{a \circ b}{|b|} \end{aligned}$

P5ImageX4

When will the circle run into the edge?

P5Image4

The circle will run into the edge when a line (yellow) drawn perpendicular from the edge to the center of the circle has a length that is equal to the circle’s radius. P5Image5

I’ll use the operation defined earlier to get the height of the circle, c, above the edge, e

$\begin{aligned} h &= \frac{c \circ e}{|e|}\\ &= \frac{c_x e_y - c_y e_x}{\sqrt{e_x^2 + e_y^2}} \end{aligned}$

P5Image7

As circle, c, travels along its velocity vector, v, through time,

$\begin{aligned} c_x &= v_x t + c_{x0}\\ c_y &= v_y t + c_{y0} \end{aligned}$

I’ll make this substitution for $c_x$ and $c_y$ and solve for time t, when the height of the yellow line is equal to the radius of the circle, $c_r$ .

$h = \pm c_r = \frac{[v_x t + c_{x0}]e_y-[v_y t + c_{y0}]e_x}{\sqrt{e_x^2 + e_y^2}}$

With a little algebra, the time at which the circle runs into the edge,

$t = \frac{c_{x0}e_y-c_{y0}e_x \pm c_r\sqrt{e_x^2 + e_y^2}}{e_x v_y - e_y v_x}$

This can be written more compactly:

$t = \frac{c \circ e \pm c_r |e|}{e \circ v}$

The plus or minus symbol in the formula arises because the circle can be in contact with the line segment from both above and below, or more accurately, when the circle is oriented counter-clockwise and when it is oriented clockwise with respect to the line segment.

In what direction should the circle be deflected off of the edge?
P5Image15
The circle should be deflected in the direction of a new velocity vector whose component that is perpendicular to the edge has been negated, but whose component that is parallel to the edge has remained unchanged.
P5Image19
With respect to the edge, e, the perpendicular and parallel components of the velocity vector before deflection, $v_i$ , are respectively

$\begin{aligned} &\frac{v_i \circ e}{|e|}, &\frac{v_i \cdot e}{|e|} \end{aligned}$

P5Image16
And so the logic above gives the following system of equations for the velocity after deflection, $v_f$ :

$\begin{cases} \frac{v_f \circ e}{|e|} &= -\frac{v_i \circ e}{|e|}\\ \frac{v_f \cdot e}{|e|} &= \frac{v_i \cdot e}{|e|}\\ \end{cases}$

i.e.

$\begin{cases} v_{fx}e_y - v_{fy}e_x &= -(v_{ix}e_y - v_{iy}e_x)\\ v_{fx}e_x + v_{fy}e_y &= v_{ix}e_x + v_{iy}e_y \end{cases}$

Because it will appear in other calculations, for the sake of reuse I’ll state that the coefficient matrix of this system has inverse,

$\frac{1}{|e|^2} \begin{bmatrix} e_y & e_x\\ -e_x & e_y \end{bmatrix}$

Multiplying both sides of the system gives the solution,

$\begin{cases} v_{fx} &= \frac{-e_y(v_{ix}e_y - v_{iy}e_x) + e_x(v_{ix}e_x + v_{iy}e_y)}{e_x^2 + e_y^2}\\ v_{fy} &= \frac{e_x(v_{ix}e_y - v_{iy}e_x) + e_y(v_{ix}e_x + v_{iy}e_y)}{e_x^2 + e_y^2} \end{cases}$

i.e.

$\begin{cases} v_{fx} &= \frac{-e_y(v_i \circ e) + e_x(v_i \cdot e)}{|e|^2}\\ v_{fy} &= \frac{e_x(v_i \circ e) + e_y(v_i \cdot e)}{|e|^2} \end{cases}$

P5Image17

When will the two circles collide with one another?

P5Image11
The two circles will collide with one another when the distance between the centers of the two circles is equal to the sum of their radii.

$\sqrt{(A_x - B_x)^2 + (A_y - B_y)^2} = A_r + B_r$

As circles A and B move in accordance of their velocity vectors,

$\begin{cases} A_x &= A_{x0} + A_{vx}t \\ A_y &= A_{y0} + A_{vy}t \\ \end{cases} \\ \begin{cases} B_x &= B_{x0} + B_{vx}t \\ B_y &= B_{y0} + B_{vy}t \\ \end{cases}$

Then

$\begin{aligned} ([A_{x0} + A_{vx}t] - [B_{x0} + B_{vx}t])^2 + ([A_{y0} + A_{vy}t] - [B_{y0} + B_{vy}t])^2 &= (A_r + B_r)^2 \\ ([A_{x0} - B_{x0}] + [A_{vx}-B_{vx}]t])^2 + ([A_{y0} - B_{y0}] + [A_{vy}-B_{vy}]t)^2 &= (A_r + B_r)^2 \\ \end{aligned}$

gives

$\begin{aligned} &(A_{x0} - B_{x0})^2 + 2*(A_{x0} - B_{x0})(A_{vx}-B_{vx})t + (A_{vx}-B_{vx})^2 t^2 \\ + &(A_{y0} - B_{y0})^2 + 2*(A_{y0} - B_{y0})(A_{vy}-B_{vy})t + (A_{vy}-B_{vy})^2 t^2\\ - &(A_r + B_r)^2 \\= &0 \end{aligned}$

let

$\begin{cases} \alpha &= (A_{vx}-B_{vx})^2 + (A_{vy}-B_{vy})^2 \\ \beta &= 2[(A_{x0} - B_{x0})(A_{vx}-B_{vx}) + (A_{y0} - B_{y0})(A_{vy}-B_{vy})]\\ \gamma &= (A_{x0} - B_{x0})^2 + (A_{y0} - B_{y0})^2 - (A_r + B_r)^2 \end{cases}$

Solving for t with the quadratic formula gives two solutions.

$\begin{aligned} t = \frac{-\beta \pm \sqrt{\beta^2 - 4\alpha\gamma}}{2\alpha} \end{aligned}$

P5Image14
The solution that is the smallest positive value gives the time at which the soonest collision occurs.
P5Image12
The solution that is the largest positive value gives the time at which another collision would occur should the circles pass through one another, rather than bounce away from each other.
P5Image13
If the two circles are not going to collide, t will either be a complex number or infinity.

How should the two circles be deflected away from one another?

P5Image20

The component of the circle’s velocity vector that is perpendicular to an edge connecting the center of the two circles should remain constant, while the component of the circle’s velocity vector that is parallel to the edge connecting the centers should be the result of a one-dimensional elastic collision.

P5Image26

Let a vector, e, extend from the center of circle A, to the center of circle B.

P5Image21

Obtain the projection, $A_{pi}$ , of A’s velocity vector on vector e.

P5Image22

Also, obtain the projection, $B_{pi}$ , of B’s velocity vector onto vector e.

P5Image24

In an elastic collision, where circles have mass, $A_m$ and $B_m$ the final velocities will have projections onto edge, e,

$\begin{cases} A_{pf} &= \frac{ 2B_m B_{pi} + A_{pi}(A_m - B_m)}{A_m + B_m} \\ B_{pf} &= \frac{ 2A_m A_{pi} + A_{pi}(B_m - A_m)}{A_m + B_m} \end{cases}$

Now Circle A’s velocity after deflection can be described with the system,

$\begin{cases} \frac{A_{vf} \circ e}{|e|} &= \frac{A_{vi} \circ e}{|e|} \\ \frac{A_{vf} \cdot e}{|e|} &= A_{pf} \end{cases}$

i.e.

$\begin{cases} A_{vfx}e_y - A_{vfy}e_x &= A_{vix}e_y - A_{viy}e_x\\ A_{vfx}e_x + A_{vfy}e_y &= A_{pf}|e| \end{cases}$

Again, the coefficient matrix of this system has inverse,

$\frac{1}{|e|^2} \begin{bmatrix} e_y & e_x\\ -e_x & e_y \end{bmatrix}$

multiplying both sides gives

$\begin{cases} A_{vfx} &= \frac{e_y(A_{vix}e_y - A_{viy}e_x) + e_x(A_{pf}|e|)}{e_x^2 + e_y^2}\\ A_{vfy} &= \frac{-e_x(A_{vix}e_y - A_{viy}e_x) + e_y(A_{pf}|e|)}{e_x^2 + e_y^2} \end{cases}$

i.e.

$\begin{cases} A_{vfx} &= \frac{e_y(A_{vi} \circ e) + e_x(A_{pf}|e|)}{|e|^2}\\ A_{vfy} &= \frac{-e_x(A_{vi} \circ e) + e_y(A_{pf}|e|)}{|e|^2} \end{cases}$

Use the same process to determine the velocity of the other circle after deflection.

Parametric Lines
A useful parameterization gives points along a line defined by endpoints $(x_1, y_1)$ and $(x_2, y_2)$

$\begin{aligned} x(t) &= x_1 + t(x_2 - x_1)\\ y(t) &= y_1 + t(y_2 - y_1) \end{aligned}$

This parameterization has the property that for values $0 \leq t \leq 1$ , $(x(t), y(t))$ is on the line segment, between the two endpoints. On the other hand, if t is less than zero, or t is greater than one, $(x(t), y(t))$ is outside of the two endpoints.

Below are points (from left to right) for t = -0.5 in red, t = 0, 0.25, 0.5, 0.75, and 1 in green, and for t = 1.5 in red again:
P5Image8

At what point do two line segments intersect?

P5Image9

Lines A and B intersect for parametric values s and t, such that

$\begin{aligned} A_x(s) &= B_x(t)\\ A_y(s) &= B_y(t) \end{aligned}$

I’ll use the parametric equation of lines A and B and solve for the values s and t for which the equality is satisfied.

$\begin{aligned} A_{x1} + s(A_{x2} - A_{x1}) &= B_{x1} + t(B_{x2} - B_{x1})\\ A_{y1} + s(A_{y2} - A_{y1}) &= B_{y1} + t(B_{y2} - B_{y1}) \end{aligned}$

then

$\begin{aligned} s(A_{x2} - A_{x1}) - t(B_{x2} - B_{x1}) &= B_{x1} - A_{x1}\\ s(A_{y2} - A_{y1}) - t(B_{y2} - B_{y1}) &= B_{y1} - A_{y1} \end{aligned}$

In matrix form

$\begin{bmatrix} (A_{x2} - A_{x1}) & -(B_{x2} - B_{x1}) \\ (A_{y2} - A_{y1}) & -(B_{y2} - B_{y1}) \end{bmatrix} \begin{bmatrix} s\\ t \end{bmatrix} = \begin{bmatrix} B_{x1} - A_{x1}\\ B_{y1} - A_{y1} \end{bmatrix}$

If I define vectors

$\begin{aligned} a &= ( A_{x2} - A_{x1}, A_{y2} - A_{y1} )\\ b &= ( B_{x2} - B_{x1}, B_{y2} - B_{y1} )\\ c &= ( B_{x1} - A_{x1}, B_{y1} - A_{y1} ) \end{aligned}$

The solutions can be expressed as simple ratios of determinants.

$\begin{array}{cc} s = \frac{b \circ c}{b \circ a} & t = \frac{a \circ c}{b \circ a} \end{array}$

If the line segments do intersect, they do at point

$A(\frac{b \circ c}{b \circ a}) = B(\frac{a \circ c}{b \circ a})$

If either s or t are not in the range of 0 to 1, then the line segments do not cross. If the line segments do not cross because they run parallel to one another, the determinant $a \circ b$ must equal zero.
P5Image10

How can I tell if an object is enclosed within a region?

P5Image3

If, by protocol, regions are polygons defined with counter-clockwise directed edges, an object inside the region is then oriented counter-clockwise with respect to all the edges that define the region, whereas an object that is outside of the region must be oriented clockwise with at least one of the region’s edges.

With the operation defined earlier, it is possible to determine this property. While a directed edges aren’t vectors, the coordinates of the involved objects can be transformed to the origin and thought of as such.

bool aIsCounterClockwiseToB(point a, edge b) {

    int u_x <- b.x2 - b.x1
    int u_y <- b.y2 - b.y1
    int v_x <- a.x - b.x1
    int v_y <- a.y - b.y1
    int A <- u_x*v_y - u_y*v_x     return A > 0
}

bool aIsWithinRegionB(point a, edge[] b) {

    bool isInRegion <- true
    for edge in b {
        if !aIsCounterClockwiseToB(a, edge) {
            isInRegion <- false
        }
    }

    return isInRegion
}

Of course, regions could be defined with clockwise oriented edges and the reverse of the above algorithms could determine if an object is enclosed or not, just as well. But note that this algorithm will only work for convex polygonal regions. In order to determine if a point is inside of a concave polygonal region, it would have to first be broken down into several convex polygonal regions.

More to follow…

Month: April 2017

Assorted Derivations in 2D