### Aperiodic Hat Monotiles

Code adapted for Turtletoy from github.com/isohedral/hatviz

An aperiodic monotile
David Smith, Joseph Samuel Myers, Craig S. Kaplan, and Chaim Goodman-Strauss, 2023
cs.uwaterloo.ca/~csk/hat/

```// Aperiodic Hat Monotiles. Created by Reinder Nijhoff 2023
// @reindernijhoff
//
// https://turtletoy.net/turtle/61a52c764a
//
// Copyright (c) 2023, Craig S. Kaplan
//

const turtle = new Turtle();

const levels = 3; // min=1, max=6, step=1
const scale = 6; // min=0.5, max=10, step=0.001
const draw  = 1; // min=1, max=3, step=1 (Draw Hats, Draw Supertiles, Draw All)
const hatching = 1; // min=0, max=1, step=1 (No, Yes)

const polygons = new Polygons();
const iterator = createEinsteinTiling(levels, draw, scale);

function walk(i) {
const d = iterator.next();
if (!d.done) {
const v = d.value;

const poly = polygons.create();
if (hatching && v.geom.children.length === 0) {
switch(v.geom.fill) {
break;
break;
break;
}
}
polygons.draw(turtle, poly, false);

return true;
} else {
return false;
}
}

//
// An aperiodic monotile
// David Smith, Joseph Samuel Myers, Craig S. Kaplan, and Chaim Goodman-Strauss, 2023
// https://cs.uwaterloo.ca/~csk/hat/
//
// Copyright (c) 2023, Craig S. Kaplan
//

function* createEinsteinTiling(levels, draw = 1, scale = 1) {
const r3 = Math.sqrt(3);
const hr3 = r3 / 2;
const ident = [1, 0, 0, 0, 1, 0];

const scl2   = (a,b)   => [a[0]*b, a[1]*b];
const add2   = (a,b)   => [a[0]+b[0], a[1]+b[1]];
const sub2   = (a,b)   => [a[0]-b[0], a[1]-b[1]];
const rot2   = (a)     => [Math.cos(a), -Math.sin(a), Math.sin(a), Math.cos(a)];
const hexPt  = (x, y)  => [x + 0.5 * y, hr3 * y];
const intersect = (p1, q1, p2, q2) => {
const d = (q2[1] - p2[1]) * (q1[0] - p1[0]) - (q2[0] - p2[0]) * (q1[1] - p1[1]);
const uA = ((q2[0] - p2[0]) * (p1[1] - p2[1]) - (q2[1] - p2[1]) * (p1[0] - p2[0])) / d;
return [p1[0] + uA * (q1[0] - p1[0]), p1[1] + uA * (q1[1] - p1[1])];
}

// Affine matrix inverse
const inv = (T) => {
const det = T[0] * T[4] - T[1] * T[3];
return [ T[4] / det, -T[1] / det, (T[1] * T[5] - T[2] * T[4]) / det,
-T[3] / det, T[0] / det, (T[2] * T[3] - T[0] * T[5]) / det ];
}
// Affine matrix multiply
const mul = (A, B) => [
A[0] * B[0] + A[1] * B[3],
A[0] * B[1] + A[1] * B[4],
A[0] * B[2] + A[1] * B[5] + A[2],

A[3] * B[0] + A[4] * B[3],
A[3] * B[1] + A[4] * B[4],
A[3] * B[2] + A[4] * B[5] + A[5]
];
// Rotation matrix
const trot = (a) => [Math.cos(a), -Math.sin(a), 0, Math.sin(a), Math.cos(a), 0];
// Translation matrix
const ttrans = (p) => [1, 0, p[0], 0, 1, p[1]];
const rotAbout = (p, ang) => mul(ttrans(p), mul(trot(ang), ttrans([-p[0], -p[1]])));
// Matrix * point
const transPt = (M, P) => [M[0] * P[0] + M[1] * P[1] + M[2], M[3] * P[0] + M[4] * P[1] + M[5]];
// Match unit interval to line segment p->q
const matchSeg = (p, q) => [q[0] - p[0], p[1] - q[1], p[0], q[1] - p[1], q[0] - p[0], p[1]];
// Match line segment p1->q1 to line segment p2->q2
const matchTwo = (p1, q1, p2, q2) => mul(matchSeg(p2, q2), inv(matchSeg(p1, q1)));

class Geom {
constructor(pgon, fill = 0) {
this.shape = pgon;
this.fill = fill
this.width = 1.0;
this.children = [];
}

this.children.push({
T: T,
geom: geom
});
}

evalChild(n, i) {
return transPt(this.children[n].T, this.children[n].geom.shape[i]);
}

recentre() {
let tr = [0, 0];
this.shape.forEach( p => tr = add2(tr, p));
tr = scl2( tr, -1 / this.shape.length);

this.shape = this.shape.map (p => add2(p, tr));

const M = ttrans(tr);
for (let ch of this.children) {
ch.T = mul(M, ch.T);
}
}
}

const hat_outline = [
hexPt(0, 0), hexPt(-1, -1), hexPt(0, -2), hexPt(2, -2),
hexPt(2, -1), hexPt(4, -2), hexPt(5, -1), hexPt(4, 0),
hexPt(3, 0), hexPt(2, 2), hexPt(0, 3), hexPt(0, 2),
hexPt(-1, 2)
];

const H1_hat = new Geom(hat_outline, 0);
const H_hat = new Geom(hat_outline, 1);
const T_hat = new Geom(hat_outline, 2);
const P_hat = new Geom(hat_outline, 3);
const F_hat = new Geom(hat_outline, 4);

function constructPatch(H, T, P, F) {
const rules = [
['H'],
[0, 0, 'P', 2],
[1, 0, 'H', 2],
[2, 0, 'P', 2],
[3, 0, 'H', 2],
[4, 4, 'P', 2],
[0, 4, 'F', 3],
[2, 4, 'F', 3],
[4, 1, 3, 2, 'F', 0],
[8, 3, 'H', 0],
[9, 2, 'P', 0],
[10, 2, 'H', 0],
[11, 4, 'P', 2],
[12, 0, 'H', 2],
[13, 0, 'F', 3],
[14, 2, 'F', 1],
[15, 3, 'H', 4],
[8, 2, 'F', 1],
[17, 3, 'H', 0],
[18, 2, 'P', 0],
[19, 2, 'H', 2],
[20, 4, 'F', 3],
[20, 0, 'P', 2],
[22, 0, 'H', 2],
[23, 4, 'F', 3],
[23, 0, 'F', 3],
[16, 0, 'P', 2],
[9, 4, 0, 2, 'T', 2],
[4, 0, 'F', 3]
];

ret = new Geom([], null, null);
ret.width = H.width;
shapes = {
'H': H,
'T': T,
'P': P,
'F': F
};

for (let r of rules) {
if (r.length == 1) {
} else if (r.length == 4) {
const poly = ret.children[r[0]].geom.shape;
const T = ret.children[r[0]].T;
const P = transPt(T, poly[(r[1] + 1) % poly.length]);
const Q = transPt(T, poly[r[1]]);
const nshp = shapes[r[2]];
const npoly = nshp.shape;

ret.addChild(matchTwo(npoly[r[3]], npoly[(r[3] + 1) % npoly.length], P, Q), nshp);
} else {
const chP = ret.children[r[0]];
const chQ = ret.children[r[2]];

const P = transPt(chQ.T, chQ.geom.shape[r[3]]);
const Q = transPt(chP.T, chP.geom.shape[r[1]]);
const nshp = shapes[r[4]];
const npoly = nshp.shape;

ret.addChild(matchTwo(npoly[r[5]], npoly[(r[5] + 1) % npoly.length], P, Q), nshp);
}
}

return ret;
}

function constructMetatiles(patch) {
const bps1 = patch.evalChild(8, 2);
const bps2 = patch.evalChild(21, 2);
const rbps = transPt(rotAbout(bps1, -2.0 * Math.PI / 3.0), bps2);

const p72 = patch.evalChild(7, 2);
const p252 = patch.evalChild(25, 2);

const llc = intersect(bps1, rbps, patch.evalChild(6, 2), p72);
let w = sub2(patch.evalChild(6, 2), llc);

const new_H_outline = [llc, bps1];
w = transPt(trot(-Math.PI / 3), w);
new_H_outline.push(patch.evalChild(14, 2));
w = transPt(trot(-Math.PI / 3), w);
new_H_outline.push(sub2(new_H_outline[3], w));
new_H_outline.push(patch.evalChild(6, 2));

const new_H = new Geom(new_H_outline);
new_H.width = patch.width * 2;
for (let ch of [0, 9, 16, 27, 26, 6, 1, 8, 10, 15]) {
}

const new_P_outline = [p72, add2(p72, sub2(bps1, llc)), bps1, llc];
const new_P = new Geom(new_P_outline);
new_P.width = patch.width * 2;
for (let ch of [7, 2, 3, 4, 28]) {
}

const new_F_outline = [bps2, patch.evalChild(24, 2), patch.evalChild(25, 0), p252, add2(p252, sub2(llc, bps1))];
const new_F = new Geom(new_F_outline);
new_F.width = patch.width * 2;
for (let ch of [21, 20, 22, 23, 24, 25]) {
}

const AAA = new_H_outline[2];
const BBB = add2(new_H_outline[1], sub2(new_H_outline[4], new_H_outline[5]));
const CCC = transPt(rotAbout(BBB, -Math.PI / 3), AAA);
const new_T_outline = [BBB, CCC, AAA];
const new_T = new Geom(new_T_outline);
new_T.width = patch.width * 2;

new_H.recentre();
new_P.recentre();
new_F.recentre();
new_T.recentre();

return [new_H, new_T, new_P, new_F]
}

// init start tiles
const H_outline = [[0, 0], [4, 0], [4.5, hr3], [2.5, 5 * hr3], [1.5, 5 * hr3], [-0.5, hr3]];
const H_init = new Geom(H_outline);
H_init.width = 2;
H_init.addChild(mul(ttrans([2.5, hr3]), mul([-0.5, -hr3, 0, hr3, -0.5, 0], [0.5, 0, 0, 0, -0.5, 0])), H1_hat);

const T_outline = [[0, 0], [3, 0], [1.5, 3 * hr3]];
const T_init = new Geom(T_outline);
T_init.width = 2;
T_init.addChild([0.5, 0, 0.5, 0, 0.5, hr3], T_hat);

const P_outline = [[0, 0], [4, 0], [3, 2 * hr3], [-1, 2 * hr3]];
const P_init = new Geom(P_outline);
P_init.width = 2;
P_init.addChild([0.5, 0, 1.5, 0, 0.5, hr3], P_hat);
P_init.addChild(mul(ttrans([0, 2 * hr3]), mul([0.5, hr3, 0, -hr3, 0.5, 0], [0.5, 0.0, 0.0, 0.0, 0.5, 0.0])), P_hat);

const F_outline = [[0, 0], [3, 0], [3.5, hr3], [3, 2 * hr3], [-1, 2 * hr3]];
const F_init = new Geom(F_outline);
F_init.width = 2;
F_init.addChild([0.5, 0, 1.5, 0, 0.5, hr3], F_hat);
F_init.addChild(mul(ttrans([0, 2 * hr3]), mul([0.5, hr3, 0, -hr3, 0.5, 0], [0.5, 0.0, 0.0, 0.0, 0.5, 0.0])), F_hat);

let tiles = [H_init, T_init, P_init, F_init];

// create all tiles for level
for (let i = 0; i < levels; i++) {
const patch = constructPatch(...tiles);
tiles = constructMetatiles(patch);
}

// generator code
const queue = [{
T: [scale, 0, 0, 0, scale, 0],
geom: tiles[0],
level: levels
}];
do {
const t = queue.pop();
if (t.level >= 0) {
for (let g of t.geom.children) {
queue.push({
T: mul(t.T, g.T),
geom: g.geom,
level: t.level - 1
});
}
}
if ( (t.level < 0 && draw == 1) || (t.level >= 0 && draw == 2) || draw == 3) {
yield {
geom: t.geom,
polygon: t.geom.shape.map(p => transPt(t.T, p))
};
}
} while (queue.length > 0);
}

////////////////////////////////////////////////////////////////
// Polygon Clipping utility code - Created by Reinder Nijhoff 2019
// (Polygon binning by Lionel Lemarie 2021)
// https://turtletoy.net/turtle/a5befa1f8d
////////////////////////////////////////////////////////////////
function Polygons(){const t=[],s=25,e=Array.from({length:s**2},t=>[]),n=class{constructor(){this.cp=[],this.dp=[],this.aabb=[]}addPoints(...t){let s=1e5,e=-1e5,n=1e5,h=-1e5;(this.cp=[...this.cp,...t]).forEach(t=>{s=Math.min(s,t[0]),e=Math.max(e,t[0]),n=Math.min(n,t[1]),h=Math.max(h,t[1])}),this.aabb=[s,n,e,h]}addSegments(...t){t.forEach(t=>this.dp.push(t))}addOutline(){for(let t=0,s=this.cp.length;t<s;t++)this.dp.push(this.cp[t],this.cp[(t+1)%s])}draw(t){for(let s=0,e=this.dp.length;s<e;s+=2)t.jump(this.dp[s]),t.goto(this.dp[s+1])}addHatching(t,s){const e=new n;e.cp.push([-1e5,-1e5],[1e5,-1e5],[1e5,1e5],[-1e5,1e5]);const h=Math.sin(t)*s,o=Math.cos(t)*s,a=200*Math.sin(t),i=200*Math.cos(t);for(let t=.5;t<150/s;t++)e.dp.push([h*t+i,o*t-a],[h*t-i,o*t+a]),e.dp.push([-h*t+i,-o*t-a],[-h*t-i,-o*t+a]);e.boolean(this,!1),this.dp=[...this.dp,...e.dp]}inside(t){let s=0;for(let e=0,n=this.cp.length;e<n;e++)this.segment_intersect(t,[.1,-1e3],this.cp[e],this.cp[(e+1)%n])&&s++;return 1&s}boolean(t,s=!0){const e=[];for(let n=0,h=this.dp.length;n<h;n+=2){const h=this.dp[n],o=this.dp[n+1],a=[];for(let s=0,e=t.cp.length;s<e;s++){const n=this.segment_intersect(h,o,t.cp[s],t.cp[(s+1)%e]);!1!==n&&a.push(n)}if(0===a.length)s===!t.inside(h)&&e.push(h,o);else{a.push(h,o);const n=o[0]-h[0],i=o[1]-h[1];a.sort((t,s)=>(t[0]-h[0])*n+(t[1]-h[1])*i-(s[0]-h[0])*n-(s[1]-h[1])*i);for(let n=0;n<a.length-1;n++)(a[n][0]-a[n+1][0])**2+(a[n][1]-a[n+1][1])**2>=.001&&s===!t.inside([(a[n][0]+a[n+1][0])/2,(a[n][1]+a[n+1][1])/2])&&e.push(a[n],a[n+1])}}return(this.dp=e).length>0}segment_intersect(t,s,e,n){const h=(n[1]-e[1])*(s[0]-t[0])-(n[0]-e[0])*(s[1]-t[1]);if(0===h)return!1;const o=((n[0]-e[0])*(t[1]-e[1])-(n[1]-e[1])*(t[0]-e[0]))/h,a=((s[0]-t[0])*(t[1]-e[1])-(s[1]-t[1])*(t[0]-e[0]))/h;return o>=0&&o<=1&&a>=0&&a<=1&&[t[0]+o*(s[0]-t[0]),t[1]+o*(s[1]-t[1])]}};return{list:()=>t,create:()=>new n,draw:(n,h,o=!0)=>{reducedPolygonList=function(n){const h={},o=200/s;for(var a=0;a<s;a++){const c=a*o-100,r=[0,c,200,c+o];if(!(n[3]<r[1]||n[1]>r[3]))for(var i=0;i<s;i++){const c=i*o-100;r[0]=c,r[2]=c+o,n[0]>r[2]||n[2]<r[0]||e[i+a*s].forEach(s=>{const e=t[s];n[3]<e.aabb[1]||n[1]>e.aabb[3]||n[0]>e.aabb[2]||n[2]<e.aabb[0]||(h[s]=1)})}}return Array.from(Object.keys(h),s=>t[s])}(h.aabb);for(let t=0;t<reducedPolygonList.length&&h.boolean(reducedPolygonList[t]);t++);h.draw(n),o&&function(n){t.push(n);const h=t.length-1,o=200/s;e.forEach((t,e)=>{const a=e%s*o-100,i=(e/s|0)*o-100,c=[a,i,a+o,i+o];c[3]<n.aabb[1]||c[1]>n.aabb[3]||c[0]>n.aabb[2]||c[2]<n.aabb[0]||t.push(h)})}(h)}}}

```