I want to make more demos and games with rigid body physics, but I don't like how big and clunky most engines are. Unfortunately I can rarely find complete and standalone versions of the algorithms that these physics engines are composed of. This article is one such stand-alone implementation.
In this article, I demonstrate a concise convex hull to convex hull collision implementation using the separating axis theorem. When a collision is detected, the implementation also generates a contact point and a minimum translation vector to push the two hulls apart.
The code in this article is designed to avoid allocations in typical hot code paths for a physics
engine by using static memory and output memory locations that were allocated outside of the
routines and passed in as arguments. The code also makes use of the
vec3 routines from
glmatrix (also available as a standalone package).
Read on or jump to the final result: hull.js
Or jump ahead to the interactive demo.
The method of separating axes compares the interval of the minimum and maximum values for each vertex in each hull projected onto each separating axis. If a separation exists, then the hulls are not colliding. If no separation exists, then the interval with the minimum overlap is returned as the minimum translation vector.
The contact point for a face separation is the support point on the opposite hull and the contact point for edges is a point on one of the hulls that minimizes the distance between the two edges.
Each face normal must be tested as a separating axis. Each edge must be tested against all the non-parallel edges in the other hull with the vector orthogonal to edge pair as the separating axis to test.
To test by separating axes, project each vertex from each hull onto the separating axis with the dot product.
The dot product on its own will project a point onto a separating axis as a single scalar value.
It's also sometimes helpful to subtract the projection of a point Q that lies on the
face that the normal comes from, so that the sign of the result determines which side of the plane a
point lies on.
That is, to project a point P onto a separating axis N from a plane
with normal N (the separating axis) and a point on the plane Q:
vec3.dot(N, P) - vec3.dot(N, Q)A convex hull can be described from its vertices alone, but our implementation requires a few more supplementary pieces of information. You could choose to arrange these differently in your own implementation since they can all be derived from the set of verticies.
Each hull in this implementation stores:
[x,y,z] array.
Now we can define a Hull class with some static fields we will use in multiple places.
class Hull {
static #tv0 = [0,0,0]; static #tv1 = [0,0,0]; static #tv2 = [0,0,0]
static #supports = [0,0,0,0] // used by queryFaceDirections and queryEdgeDirections
constructor({ positions, edges, faces, normals }) {
this.faces = faces
this.edges = edges
this.normals = normals
this.positions = positions
}
// ...
}
First we need a function that takes a plane normal (and separating axis) N and point on
the plane Q and projects all the points in the hull onto the plane normal, where points
on the "inside" (with respect to the hull center) will have a negative projection and "outside" will
be negative.
This function writes to out with the index and projected value for the vertices with
the minimum and maximum projected values.
getSupports(out, Q, N) { // out: [minI,min,maxI,max]
let dNQ = vec3.dot(N, Q)
out[0] = -1
out[1] = Infinity
out[2] = -1
out[3] = -Infinity
for (let i = 0; i < this.positions.length; i++) {
let p = this.positions[i]
let d = vec3.dot(N, p) - dNQ
if (d < out[1]) {
out[0] = i
out[1] = d
}
if (d > out[3]) {
out[2] = i
out[3] = d
}
}
return out
}
Now we can implement the first type of separating axis test by looping over each face and its corresponding normal and finding support points for that separating axis on the other hull.
If both the minimum and maximum support projections are positive, the hull is completely on the other side of the current face and we've found a separation and so the hulls must not be collding.
Otherwise if the minimum is negative and the maximum is positive, the minimum support could be the minimum separation between the hulls if it has the greatest value (because overlapping separations here are negative).
static queryFaceDirections(out, A, B) {
out.face = -1
out.support = -1
out.separation = -Infinity
let min = Infinity, max = -Infinity
for (let i = 0; i < A.faces.length; i++) {
let aN = A.normals[i]
let aQ = A.positions[A.faces[i][0]]
let [minI,min,maxI,max] = B.getSupports(this.#supports, aQ, aN)
if (min > 0.0 && max > 0.0) return false
if (min < 0.0 && max < 0.0) continue
if (min > out.separation) {
out.separation = min
out.support = minI
out.face = i
}
}
return true
}
Similar to queryFaceDirections, we can search for minimum separation candidates or rule
out any separation by testing separation axes orthogonal to each pair of edges from each hull.
This version calculates supports for both hulls and compares the overlap intervals, ruling out axes that show a separation and otherwise saving the maximum negative separation.
static queryEdgeDirections(out, A, B) {
out.index[0] = -1
out.index[1] = -1
out.separation = -Infinity
for (let i = 0; i < A.edges.length; i++) {
let ea = A.edges[i]
let a0 = A.positions[ea[0]]
let a1 = A.positions[ea[1]]
vec3.subtract(Hull.#tv0, a0, a1)
vec3.normalize(Hull.#tv0, Hull.#tv0)
for (let j = 0; j < B.edges.length; j++) {
let eb = B.edges[j]
let b0 = B.positions[eb[0]]
let b1 = B.positions[eb[1]]
vec3.subtract(Hull.#tv1, b0, b1)
vec3.normalize(Hull.#tv1, Hull.#tv1)
let N = vec3.cross(Hull.#tv2, Hull.#tv0, Hull.#tv1)
if (vec3.length(N) < 0.5) continue
vec3.normalize(N, N)
// flip edges that are pointing the wrong way
if (vec3.dot(N, b0) > vec3.dot(N, a0)) {
vec3.negate(N, N)
}
let [minIA,minA,maxIA,maxA] = A.getSupports(this.#supports, a0, N)
let [minIB,minB,maxIB,maxB] = B.getSupports(this.#supports, a0, N)
if (minA >= maxB) return false
let s = minB - maxA
if (s > out.separation) {
out.index[0] = i
out.index[1] = j
out.separation = s
vec3.copy(out.normal, N)
}
}
}
return true
}
The final piece required before collision can be implemented is a routine to calculate the closest point on each line between two lines. This will be needed to calculate the contact point for edge/edge intersections.
// based on https://paulbourke.net/geometry/pointlineplane/
static #p13 = [0,0,0]; static #p43 = [0,0,0]; static #p21 = [0,0,0]
static #lineLine(out0, out1, p1, p2, p3, p4, E=1e-8) {
let p13 = vec3.subtract(this.#p13, p1, p3)
let p43 = vec3.subtract(this.#p43, p4, p3)
if (vec3.length(p43) < E) return false
let p21 = vec3.subtract(this.#p21, p2, p1)
if (vec3.length(p21) < E) return false
let d1343 = vec3.dot(p13, p43)
let d4321 = vec3.dot(p43, p21)
let d1321 = vec3.dot(p13, p21)
let d4343 = vec3.dot(p43, p43)
let d2121 = vec3.dot(p21, p21)
let denom = d2121 * d4343 - d4321 * d4321
if (Math.abs(denom) < E) return false
let mua = (d1343 * d4321 - d1321 * d4343) / denom
let mub = (d1343 + d4321 * mua) / d4343
vec3.scaleAndAdd(out0, p1, p21, mua)
vec3.scaleAndAdd(out1, p3, p43, mub)
return true
}
Now that all the required functions have been specified, the collide function can use them to test for collisions and calculate a contact point and minimum translation vector. If any of the functions finds a positive separation along a separating axis, the hulls are not colliding so the rest of the collision checks can be skipped.
static TYPE_A = 0; static TYPE_B = 1; static TYPE_E = 2
static #qA = { type: this.TYPE_A, face: -1, support: -1, separation: 0 }
static #qB = { type: this.TYPE_B, face: -1, support: -1, separation: 0 }
static #qE = { type: this.TYPE_E, index: [-1,-1], separation: -1, normal: [0,0,0] }
static collide(out, A, B) {
if (!this.queryFaceDirections(this.#qA, A, B)) return null
if (!this.queryFaceDirections(this.#qB, B, A)) return null
if (!this.queryEdgeDirections(this.#qE, A, B)) return null
out.info = this.#qA
if (this.#qB.separation > out.info.separation) {
out.info = this.#qB
}
if (this.#qE.separation > out.info.separation) {
out.info = this.#qE
}
if (out.info.type === this.TYPE_A) {
vec3.copy(out.point, B.positions[this.#qA.support])
let N = A.normals[this.#qA.face]
vec3.scale(out.mtv, N, -this.#qA.separation)
} else if (out.info.type === this.TYPE_B) {
vec3.copy(out.point, A.positions[this.#qB.support])
let N = B.normals[this.#qB.face]
vec3.scale(out.mtv, N, this.#qB.separation)
} else if (out.info.type === this.TYPE_E) {
let ea = A.edges[this.#qE.index[0]]
let eb = B.edges[this.#qE.index[1]]
let a0 = A.positions[ea[0]]
let a1 = A.positions[ea[1]]
let b0 = B.positions[eb[0]]
let b1 = B.positions[eb[1]]
vec3.scale(out.mtv, this.#qE.normal, -this.#qE.separation)
this.#lineLine(Hull.#tv0, Hull.#tv1, a0, a1, b0, b1)
vec3.copy(out.point, Hull.#tv0)
}
return out
}
The out argument to give to collide() to write its output should look
like:
{ mtv: [0,0,0], point: [0,0,0], info: null }Here is some example data with an example usage of a tetrahedron (A) colliding with a cube (B).
let A = new Hull({ // tetrahedron
faces: [[0,1,2],[0,2,3],[0,3,1],[1,2,3]],
edges: [[0,1],[0,2],[0,3],[1,2],[2,3],[3,1]],
normals: [
[0.38748992453971465,0.004074648544093487,0.9218649335013498],
[-0.6349199599421436,0.7711759948353863,0.046521279611796626],
[0.7500043993541771,0.3903606074731003,-0.5339587971769095],
[-0.3746117365449317,-0.8688005053127466,-0.3238081666837741]
],
positions:[
[-0.015474963193511798,0.33537417751658094,-0.1404996212167574],
[0.5141696198411252,-0.9568235666566539,-0.35124017154634457],
[-1.0335072548392017,-0.5373752088540689,0.30094133605231416],
[-0.6342374754213884,-0.1119237570390201,-1.3282023008556973]
],
})
let B = new Hull({ // cube
faces: [[0,1,2,3],[4,5,6,7],[0,1,4,5],[1,2,5,6],[2,3,6,7],[3,0,7,4]],
edges: [[0,1],[1,2],[2,3],[3,0],[4,5],[5,6],[6,7],[7,4],[0,4],[1,5],[2,6],[3,7]],
normals: [
[-0.2016831109350874,-0.9760273663061185,-0.08181994246567997],
[0.2016831109350874,0.9760273663061185,0.08181994246567997],
[-0.688355761979541,0.08181994246567997,0.720743950348859],
[0.6967713166549776,-0.2016831109350874,0.6883557619795408],
[0.688355761979541,-0.08181994246567997,-0.720743950348859],
[-0.6967713166549776,0.2016831109350874,-0.6883557619795408]
],
positions: [
[-1.2969646066524665,-0.14439705019398202,0.08105014237762323],
[-0.5917402473790571,-0.34852692650445355,0.7656251200986752],
[0.09283473034199496,-0.43133948655379123,0.0361372975730391],
[-0.6123896289314145,-0.2272096102433197,-0.648437680148013],
[-1.092834730341995,0.8313394865537913,0.1638627024269609],
[-0.38761037106858554,0.6272096102433198,0.8484376801480129],
[0.2969646066524665,0.5443970501939821,0.11894985762237678],
[-0.40825975262094294,0.7485269265044536,-0.5656251200986753]
],
})
let manifold = { mtv: [0,0,0], point: [0,0,0], info: null }
if (Hull.collide(manifold, A, B)) {
console.log(manifold)
} else {
console.log('no collision')
}
/* output:
{
mtv: [ -0.17742709499097423, 0.3525922983309437, 0.2879166461579119 ],
point: [ -0.39404904235014193, 0.010826225559499636, 0.023658338478221408 ],
info: {
type: 2,
index: [ 1, 2 ],
separation: -0.48856698416292277,
normal: [ -0.363158176344981, 0.7216867077808199, 0.5893084377185425 ]
}
}
*/
The result is an edge-edge collision (type 2) with a penetration depth of 0.4886 between edge index 1 on hull A and edge index 2 on hull B.
Check out the interactive demo which features an interactive editor and viewer along with preset geometries.
The demo uses these additional values for positioning and rotation: