2b:null 2c:["$","title",null,{"children":"2D GJK and EPA algorithms | Claudio Esperança | Observable"}] 31:Ta42,md`# 2D GJK and EPA algorithms The GJK (Gilbert-Johnson-Keerthi) algorithm is used to test whether two convex shapes intersect. It is frequently used in physical simulations to test for collisions between objects. In this implementation, we are focusing on the bidimensional case and polygons. A pretty good explanation of the algorithm and its implementation can be found in the [Wikipedia page](https://en.wikipedia.org/wiki/Gilbert%E2%80%93Johnson%E2%80%93Keerthi_distance_algorithm), which referred me to this really cool [video](https://caseymuratori.com/blog_0003) by Casey Muratori. Since we are working in 2D, the implementation is considerably simpler. The GJK algorithm is based on the concept of [Minkowski differences](http://twistedoakstudios.com/blog/Post554_minkowski-sums-and-differences). Let ${tex`A`} and ${tex`B`} be convex sets of vectors. Then their Minkowski difference is ${tex.block`A \ominus B=\{ a-b ~|~ a \in A, b \in B\}.`} We can show that ${tex`A`} and ${tex`B`} intersect if and only if ${tex`A \ominus B`} contains the origin. The [Expanding Polytope Algorithm (EPA)](https://graphics.stanford.edu/courses/cs468-01-fall/Papers/van-den-bergen.pdf) can be thought as a companion for the GJK, and is used to compute the penetration depth and direction of the collision. This information can then be used by the physics engine to compute a suitable response to the collision, for instance, to make the two shapes bounce from each other. Our implementation was strongly influenced by [this post](http://www.dyn4j.org/2010/05/epa-expanding-polytope-algorithm/). The demo below (left) shows two shapes which you can drag with the mouse. Watch the GJK algorithm detect the intersection and the EPA compute the vertices involved in the collision. The panel on the right shows what's being computed: If ${tex` \partial A`} is the set of vertices of polygon ${tex`A`} and similarly for ${tex`B`}, then the green points are ${tex`\{ a-b ~|~ a \in \partial A, b \in \partial B\}`}. The green polygon is the convex hull of this set. If the hull contains the origin (blue point), then the shapes intersect. The objective of the GJK algorithm is to efficiently build a triangle (shown in black) that contains the origin. If successful, this triangle is fed to the EPA algorithm that will then locate the hull edge closest to the origin (shown in yellow). This edge is important because it contains information about the vertices of ${tex`A`} and ${tex`B`} more closely involved in the penetration, which will help in the collision response. These are shown as black dots in the canvas on the left. `32:T583,demo = { let dsp = new Display(width / 2 - 2, width / 2 - 2); let mouse = null; let selected = null; let iSelected = -1; function refresh() { dsp.clear(); let hit = gjk(...polygons); let color = hit ? 'red' : 'black'; for (let poly of polygons) { dsp.drawPolygon(poly, { stroke: color }); } if (hit) { let { p, q, n } = epa(...polygons, ...hit); let pts = supportEdge(polygons[0], n).concat( supportEdge(polygons[1], [-n[0], -n[1]]) ); let [cx, cy] = centerOfMass(pts); let sep = [ [cx - n[1] * 1000, cy + n[0] * 1000], [cx + n[1] * 1000, cy - n[0] * 1000] ]; dsp.drawPoints(pts); dsp.drawPolygon(sep, { stroke: 'blue', close: false, dash: [5, 5] }); } } dsp.canvas.onmousedown = e => { mouse = [e.offsetX, e.offsetY]; [selected, iSelected] = [null, -1]; polygons.forEach((poly, i) => { if (pointInPolygon(mouse, poly)) [selected, iSelected] = [poly, i]; }); refresh(); }; dsp.canvas.onmousemove = e => { if (selected) { let newMouse = [e.offsetX, e.offsetY]; let v = vec2.sub([], newMouse, mouse); selected.forEach(p => vec2.add(p, p, v)); mutable changedPolygon++; refresh(); mouse = newMouse; } }; dsp.canvas.onmouseup = e => { [selected, iSelected] = [null, -1]; }; refresh(); return dsp.canvas; }33:T6db,minkDisplay = { changedPolygon; let dsp = new Display(width / 2 - 2, width / 2 - 2); // Bounding box of the Minkowski difference let minp = [Number.POSITIVE_INFINITY, Number.POSITIVE_INFINITY]; let maxp = [Number.NEGATIVE_INFINITY, Number.NEGATIVE_INFINITY]; for (let p of minkDiff(...polygons)) { minp = minp.map((x, i) => Math.min(p[i], x)); maxp = maxp.map((x, i) => Math.max(p[i], x)); } // Compute an affine transformation to show the Minkowski difference and its convex hull // centered inside the display let center = vec2.scale([], vec2.add([], minp, maxp), 0.5); let size = vec2.sub([], maxp, minp); let ccenter = [dsp.canvas.width / 2, dsp.canvas.height / 2]; let scale = Math.min(...vec2.div([], ccenter, vec2.scale([], size, 0.5))) * 0.8; let transf = p => p.map((x, i) => (x - center[i]) * scale + ccenter[i]); // Draw the origin in blue dsp.drawPoints([[0, 0]].map(transf), { fill: 'blue' }); // Draw the Minkowski difference and its hull let diff = minkDiff(...polygons).map(transf); dsp.drawPoints(diff, { radius: 4, fill: null, stroke: 'green' }); dsp.drawPolygon(minkPolygon(...polygons).map(transf), { stroke: 'green' }); // The result of the gjk algorithm let hit = gjk(...polygons); if (hit) { let [a, b, c] = hit; dsp.drawPolygon(hit.map(transf), { lineWidth: 3, stroke: 'black' }); // The EPA edge let { p, q } = epa(...polygons, a, b, c); dsp.drawPolygon([p, q].map(transf), { lineWidth: 4, stroke: 'yellow' }); // The hull indices // let l = polygons[1].length; // diff.forEach(([x, y], i) => // dsp.ctx.fillText(~~(i / l) + "," + (i % l), x, y) // ); } return dsp.canvas; }34:T4cb,// The Gilbert-Johnson-Keerthi algorithm. Given two polygons (pts1, pts2), visits the // convex hull of the Minkowski difference between them and locates a simplex (a triangle) within the hull // that contains the origin. Returns an array [a,b,c] with the vertices of that simplex or false if // the simplex can't be found, i.e., if the hull does not contain the origin. function gjk(pts1, pts2) { // First point of pts1-pts2 let a = support2(pts1, pts2, [1, 1]); // First direction let v = [-a[0], -a[1]]; // Second point let b = support2(pts1, pts2, v); if (vec2.dot(b, v) <= 0) return false; // Second point fails // Second direction let ab = vec2.sub([], b, a); v = tripleProduct(ab, [-a[0], -a[1]], ab); for (;;) { // Third point let c = support2(pts1, pts2, v); if (vec2.dot(c, v) <= 0) return false; // Third point fails let c0 = [-c[0], -c[1]]; let cb = vec2.sub([], b, c); let ca = vec2.sub([], a, c); let cbPerp = tripleProduct(ca, cb, cb); let caPerp = tripleProduct(cb, ca, ca); if (vec2.dot(caPerp, c0) > 0) { b = c; v = caPerp; } else if (vec2.dot(cbPerp, c0) > 0) { a = c; v = cbPerp; } else return [a, b, c]; } }35:T5f5,// A canvas display helper class class Display { constructor(w = 640, h = 480) { [this.w, this.h] = [w, h]; this.canvas = html`