31:T47a,buttons = { let div = d3 .create('div') .attr('class', 'pure-button-group') .style('margin', '10px') .style('padding', '10px') .style('width', width) .style('text-align', 'center'); for (let i = 0; i < aps.length; i++) { let s = 50; let canvas = d3 .create('canvas') .attr('width', s) .attr('height', s); draw(...aps[i], { canvas: canvas.node(), show_curvatures: false, show_more: false, show_regions: false, interactive: false }); let button = div .append('button') .attr('class', function() { if (i == 0) { return 'pure-button pure-button-active'; } else { return 'pure-button'; } }) .style('padding', '10px') .on('click', function() { div.selectAll('button').attr('class', 'pure-button'); d3.select(this).attr('class', 'pure-button pure-button-active'); d3.select(bowl) .selectAll('canvas') .remove(); d3.select(bowl).append(() => draw(...aps[i])); }) .append(() => canvas.node()); } return div.node(); }32:T4a2,comments = md`## Mathematical notes An [Apollonian packing](https://en.wikipedia.org/wiki/Apollonian_gasket) is generated from a quadruple of mutually tangent circles by recursively filling in more circles - each tangent to three of the circles in the configuration. In the special cases shown here the reciprocal of the radius of each circle is an integer - a configuration sometimes called a *bowl of integers*. This notebook displays several such packings and allows you to zoom in for a closer look. Of course, zooming in will eventually produce numerical error. ### Curvature and Descarte's circle theorem The reciprocal of the radius of a circle is known more succinctly as it's *curvature*, which are shown by the labels in the diagram. If the curvature of four mutually tangent circles are ${tex`\kappa_1`}, ${tex`\kappa_2`}, ${tex`\kappa_3`}, and ${tex`\kappa_4`}, then they must satisfy [Descartes' circle theorem](https://en.wikipedia.org/wiki/Descartes%27_theorem), namely: ${tex.block`(\kappa_1+\kappa_2+\kappa_3+\kappa_4)^2 = 2(\kappa_1^2+\kappa_2^2+\kappa_3^2+\kappa_4^2).`} This is true even for the exterior circle, assuming we take it's curvature to be negative.`33:T1b82,function draw(c0, c1, c2, c3, opts = {}) { let canvas_width, canvas_height, d3Canvas, context; // Create the canvas and and set the bounds automatically // when called the first time let { canvas = "new", bounds = "auto", interactive = true, show_curvatures = true, show_more = true, show_regions = false, n = 0 } = opts; if (canvas == "new") { if (width <= 720) { canvas_width = width; canvas_height = width; } else { canvas_width = 0.7 * width; canvas_height = 0.7 * width; } d3Canvas = d3 .create("canvas") .attr("width", canvas_width) .attr("height", canvas_height); context = d3Canvas.node().getContext("2d"); context.textAlign = "center"; context.textBaseline = "middle"; } // Redraw to the canvas when called later. else { d3Canvas = d3.select(opts.canvas); canvas_width = opts.canvas.width; canvas_height = opts.canvas.height; context = opts.canvas.getContext("2d"); } // Set the size context.clearRect(0, 0, canvas_width, canvas_height); let r = -1 / c0[2]; let R = 1.1 * r; // Set the view rectangle - either automatically // or based on specified bounds let xmin, xmax, ymin, ymax, rect; if (bounds == "auto") { xmin = -R; xmax = R; ymin = -R; ymax = R; rect = [ [xmin, ymin], [xmax, ymin], [xmax, ymax], [xmin, ymax] ]; } else { xmin = opts.xmin; ymin = opts.ymin; xmax = xmin + opts.dd; ymax = ymin + opts.dd; rect = [ [xmin, ymin], [xmax, ymin], [xmax, ymax], [xmin, ymax] ]; } // Set scale functions and path generator let x_scale = d3.scaleLinear().domain([xmin, xmax]).range([0, canvas_width]); let y_scale = d3.scaleLinear().domain([ymax, ymin]).range([0, canvas_height]); let r_scale = d3 .scaleLinear() .domain([0, (xmax - xmin) / 2]) .range([0, canvas_width / 2]); let path = d3 .line() .x((d) => x_scale(d[0])) .y((d) => y_scale(d[1])); // This option is used to generate an explantory image on zooming // later in this notebook. if (show_regions) { context.beginPath(); let r = -1 / c0[2]; context.moveTo(x_scale(r), y_scale(0)); context.arc(x_scale(0), y_scale(0), r_scale(r), 0, 2 * Math.PI); context.fillStyle = "#dddd00"; context.fill(); let [v1, v2, v3] = inner_triangle(c1, c2, c3); let x1 = x_scale(v1[0]); let y1 = y_scale(v1[1]); let x2 = x_scale(v2[0]); let y2 = y_scale(v2[1]); let x3 = x_scale(v3[0]); let y3 = y_scale(v3[1]); context.beginPath(); context.moveTo(x1, y1); context.lineTo(x2, y2); context.lineTo(x3, y3); context.closePath(); context.fillStyle = "#0000dd"; context.fill(); } // Generate the gasket and store let all_circles; if (show_more) { all_circles = gasket(c0, c1, c2, c3, rect); } else { all_circles = [c0, c1, c2, c3]; } // Draw most of the circles for (let [x, y, c] of all_circles.slice(4)) { context.beginPath(); let r = 1 / c; context.moveTo(x_scale(x + r), y_scale(y)); context.arc(x_scale(x), y_scale(y), r_scale(r), 0, 2 * Math.PI); context.stroke(); if (show_curvatures) { let s = r_scale(1 / (2 * c)); context.font = `${s}px sans-serif`; context.fillText(c, x_scale(x), y_scale(y)); context.closePath(); } } // Draw and shade the first three interior circles for (let [x, y, c] of [c1, c2, c3]) { context.beginPath(); let r = 1 / c; context.moveTo(x_scale(x + r), y_scale(y)); context.arc(x_scale(x), y_scale(y), r_scale(r), 0, 2 * Math.PI); context.fillStyle = "#dddddd"; context.fill(); context.stroke(); if (show_curvatures) { let s = r_scale(1 / (2 * c)); context.font = `${s}px sans-serif`; context.fillStyle = "#000000"; context.fillText(c, x_scale(x), y_scale(y)); context.closePath(); } } // Draw the exterior circle let [x, y, c] = c0; context.beginPath(); r = Math.abs(1 / c); context.lineWidth = 2; context.moveTo(x_scale(x + r), y_scale(y)); context.arc(x_scale(x), y_scale(y), r_scale(r), 0, 2 * Math.PI); context.stroke(); context.lineWidth = 1; context.closePath(); if (show_curvatures) { let s = 50; context.font = `${s}px sans-serif`; context.fillStyle = "#000000"; context.fillText(c, x_scale(-0.75 * R), y_scale(0.75 * R)); context.stroke(); context.closePath(); } // Set up the zooming behavior if (interactive) { d3Canvas.call( d3 .zoom() .scaleExtent([0, 8]) // On zoom, we'll simply draw the circles that have already been generated .on("zoom", function () { context.save(); context.clearRect(0, 0, canvas_width, canvas_height); context.beginPath(); for (const [x, y, c] of all_circles) { context.beginPath(); const r = Math.abs(1 / c); let [xxr, yyr] = d3.event.transform.apply([ x_scale(x + r), y_scale(y) ]); context.moveTo(xxr, yyr); let [xx, yy] = d3.event.transform.apply([x_scale(x), y_scale(y)]); context.arc(xx, yy, xxr - xx, 0, 2 * Math.PI); context.stroke(); context.closePath(); } for (let [x, y, c] of [c1, c2, c3]) { context.beginPath(); const r = Math.abs(1 / c); let [xxr, yyr] = d3.event.transform.apply([ x_scale(x + r), y_scale(y) ]); context.moveTo(xxr, yyr); let [xx, yy] = d3.event.transform.apply([x_scale(x), y_scale(y)]); context.arc(xx, yy, xxr - xx, 0, 2 * Math.PI); context.fillStyle = "#dddddd"; context.fill(); context.stroke(); context.closePath(); } context.restore(); }) // On zoom-end, we'll refine, regenerate, and redraw .on("end", function () { // Reset the window based on the transform. let [new_xmin, new_ymin] = d3.event.transform.invert([ 0, canvas_height ]); let [new_xmax, new_ymax] = d3.event.transform.invert([ canvas_width, 0 ]); new_xmin = x_scale.invert(new_xmin); new_ymin = y_scale.invert(new_ymin); new_xmax = x_scale.invert(new_xmax); let dd = new_xmax - new_xmin; // Draw to the existing canvas using the new view window draw(c0, c1, c2, c3, { canvas: d3Canvas.node(), xmin: new_xmin, ymin: new_ymin, dd: dd, bounds: "set" }); // Reset the canvas transform since drawing is based on // the window we just computed. d3.zoom().transform(d3Canvas, d3.zoomIdentity); }) ); } return d3Canvas.node(); }34:T472,function* circles(c0, c1, c2, c3, B, rect) { // This function has substantial additions const c = circle(c0, c1, c2, c3); // If the curvature of c is too large, // then we can stop. if (!(Math.abs(c[2]) < B)) { return; } // We now determine whether our viewing rectangle // intersects the region in which c lies. If not // we can stop. let yield_more = false; if (c0[2] > 0 && c1[2] > 0 && c2[2] > 0) { let T = inner_triangle(c0, c1, c2); if (intersect(T, rect)) { yield_more = true; } } else { if (c0[2] < 0) { let A = c1; let B = c2; if (some_sameSide(A, B, c, rect)) { yield_more = true; } } else if (c1[2] < 0) { let A = c0; let B = c2; if (some_sameSide(A, B, c, rect)) { yield_more = true; } } else if (c2[2] < 0) { let A = c0; let B = c1; if (some_sameSide(A, B, c, rect)) { yield_more = true; } } } if (yield_more) { yield c; yield* circles(c0, c1, c, c2, B, rect); yield* circles(c1, c, c2, c0, B, rect); yield* circles(c2, c, c0, c1, B, rect); } }35:[] 2f:{"title":"Integer Apollonian Packings","cells":[{"id":0,"value":"title = md`# Integer Apollonian Packings`","pinned":false,"hidden":false,"mode":"ojs"},{"id":1038,"value":"$31","pinned":false,"hidden":false,"mode":"ojs"},{"id":190,"value":"bowl = {\n let w, h;\n if (width <= 720) {\n w = width;\n h = width;\n } else {\n w = width;\n h = 0.7 * width;\n }\n let div = d3\n .create('div')\n .attr('class', 'bowl')\n .style('width', `${w}px`)\n .style('height', `${h}px`)\n .style('text-align', 'center');\n\n div.append(() => draw(...aps[0]));\n return div.node();\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":1089,"value":"$32","pinned":false,"hidden":false,"mode":"ojs"},{"id":791,"value":"md`## Code\n### The main draw function\n\nThe main draw function is kinda long but pretty straightforward. It's got some options so that it can display several types of configurations. For example, the button labels are drawn using this function.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":35,"value":"$33","pinned":true,"hidden":false,"mode":"ojs"},{"id":796,"value":"md`### Compute the gasket\n\nThese functions are based largely on Mike's [Appollonian Gasket notebook](https://observablehq.com/@mbostock/apollonian-gasket), though there are significant additions to generate circles within the zoom rectangle.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":17,"value":"function gasket(c0, c1, c2, c3, rect) {\n // A lot like Mike's function\n let B = (800 * (-c0[2]) ** 0.3) / (rect[1][0] - rect[0][0]);\n return [\n c0,\n c1,\n c2,\n c3,\n ...circles(c1, c2, c3, c0, B, rect),\n ...circles(c2, c3, c0, c1, B, rect),\n ...circles(c3, c0, c1, c2, B, rect),\n ...circles(c0, c1, c2, c3, B, rect)\n ];\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":15,"value":"$34","pinned":true,"hidden":false,"mode":"ojs"},{"id":12,"value":"function circle([x0, y0, b0], [x1, y1, b1], [x2, y2, b2], [x3, y3, b3]) {\n // Exactly Mike's function. Looks like it's based on Descartes' theorem\n const b = 2 * b0 + 2 * b1 + 2 * b2 - b3;\n return [\n (2 * b0 * x0 + 2 * b1 * x1 + 2 * b2 * x2 - b3 * x3) / b,\n (2 * b0 * y0 + 2 * b1 * y1 + 2 * b2 * y2 - b3 * y3) / b,\n b\n ];\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":817,"value":"md`### Geometric computations`","pinned":false,"hidden":false,"mode":"ojs"},{"id":1074,"value":"md`For effective zooming, we don't want to draw too many circles outside of our viewing rectangle. For example, when extending the initial configuration shown below into the superimposed viewing rectangle, we'd certainly want to extend into the central blue region and the yellow region on the lower right but not into the other two yellow regions. Checking for intersections between the viewing rectangle and these regions is the realm of compoutational geometry. In the code, the regions are replaced with polygons so that we end up including more circles than necessary. Nonetheless, we exclude enough so that the code runs smoothly.`","pinned":false,"hidden":false,"mode":"ojs"},{"id":1046,"value":"{\n let s = 300;\n let canvas = d3\n .create('canvas')\n .attr('width', s)\n .attr('height', s);\n draw(...aps[0], {\n canvas: canvas.node(),\n show_curvatures: false,\n show_more: false,\n show_regions: true,\n interactive: false\n });\n let context = canvas.node().getContext('2d');\n context.moveTo(0.525 * s, 0.525 * s);\n context.lineTo(0.8 * s, 0.525 * s);\n context.lineTo(0.8 * s, 0.8 * s);\n context.lineTo(0.525 * s, 0.8 * s);\n context.closePath();\n context.stroke();\n return canvas.node();\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":363,"value":"function tangency(c1, c2) {\n let [x1, y1, k1] = c1;\n let [x2, y2, k2] = c2;\n let r1 = 1 / k1;\n let r2 = 1 / k2;\n let d = r1 + r2;\n let x = (x1 * r2 + x2 * r1) / d;\n let y = (y1 * r2 + y2 * r1) / d;\n return [x, y];\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":376,"value":"function inner_triangle(c1, c2, c3) {\n let p = tangency(c1, c2);\n let q = tangency(c2, c3);\n let r = tangency(c3, c1);\n return [p, q, r];\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":374,"value":"function intersect(Tverts, Rverts) {\n let T = new Flatten.Polygon(Tverts);\n let R = new Flatten.Polygon(Rverts);\n\n // Edges might intersect\n if (T.intersect(R).length > 0) {\n return true;\n }\n\n // R might fully contain T\n let Tvs = T.vertices;\n for (let i = 0; i < Tvs.length; i++) {\n if (R.contains(Tvs[i])) {\n return true;\n }\n }\n\n // T might fully contain R\n let Rvs = R.vertices;\n for (let i = 0; i < Rvs.length; i++) {\n if (T.contains(Rvs[i])) {\n return true;\n }\n }\n\n // Or they might not intersect at all\n return false;\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":448,"value":"function some_sameSide([x0, y0], [x1, y1], [x, y], rect) {\n for (let i = 0; i < rect.length; i++) {\n if (sameSide([x0, y0], [x1, y1], [x, y], rect[i])) {\n return true;\n }\n }\n return false;\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":431,"value":"function sameSide([x0, y0], [x1, y1], [a, b], [c, d]) {\n let ab = leftOf([x0, y0], [x1, y1], [a, b]);\n let cd = leftOf([x0, y0], [x1, y1], [c, d]);\n if ((ab && cd) || (!ab && !cd)) {\n return true;\n } else {\n return false;\n }\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":429,"value":"function leftOf([x0, y0], [x1, y1], [x, y]) {\n return x1 * (y - y0) + x * (y0 - y1) + x0 * (y1 - y) > 0;\n}","pinned":true,"hidden":false,"mode":"ojs"},{"id":1077,"value":"md`### Initial Appollonian packings\n\nThere are infinitely many integer Apollonian packings but they are not necessarily easy to find. Here's a list of the six that are displayed above.\n`","pinned":false,"hidden":false,"mode":"ojs"},{"id":826,"value":"aps = [\n [[0, 0, -1], [-1 / 2, 0, 2], [1 / 2, 0, 2], [0, -2 / 3, 3]],\n [[0, 0, -2], [-1 / 6, 0, 3], [1 / 3, 0, 6], [3 / 14, -2 / 7, 7]],\n [[0, 0, -3], [-2 / 15, 0, 5], [1 / 6, -1 / 8, 8], [1 / 6, 1 / 8, 8]],\n [[0, 0, -6], [-5 / 66, 0, 11], [8 / 105, -2 / 35, 14], [3 / 50, 4 / 50, 15]],\n [\n [0, 0, -13],\n [-10 / 299, 0, 23],\n [14 / 325, -2 / 300, 30],\n [22 / 1235, 9 / 190, 38]\n ],\n [\n [0, 0, -76],\n [-1 / 5852, 0, 77],\n [1 / 77, 0, 5852],\n [2000 / 154053, -2 / 5853, 5853]\n ]\n]","pinned":true,"hidden":false,"mode":"ojs"},{"id":1080,"value":"md`### Imports`","pinned":false,"hidden":false,"mode":"ojs"},{"id":372,"value":"Flatten = require('@flatten-js/core')","pinned":true,"hidden":false,"mode":"ojs"},{"id":823,"value":"import { slider } from \"@jashkenas/inputs\"","pinned":true,"hidden":false,"mode":"ojs"},{"id":4,"value":"d3 = require('d3@5')","pinned":true,"hidden":false,"mode":"ojs"},{"id":1144,"value":"// html`\n// `","pinned":true,"hidden":false,"mode":"ojs"}],"files":"$W35","stdlib":"1","titling":"auto"}