Loopy / grid.c: support the new Spectre monotiling.

This uses a tile shape very similar to the hat, but the tiling
_structure_ is totally changed so that there aren't any reflected
copies of the tile.

I'm not sure how much difference this makes to gameplay: the two
tilings are very similar for Loopy purposes. But the code was fun to
write, and I think the Spectre shape is noticeably prettier, so I'm
adding this to the collection anyway.

The test programs also generate a pile of SVG images used in the
companion article on my website.

CMakeLists.txt

@@ -10,7 +10,7 @@ add_library(core_obj OBJECT
   laydomino.c loopgen.c malloc.c matching.c midend.c misc.c penrose.c
   ps.c random.c sort.c tdq.c tree234.c version.c ${platform_common_sources})
 add_library(core $<TARGET_OBJECTS:core_obj>)
-add_library(common $<TARGET_OBJECTS:core_obj> hat.c)
+add_library(common $<TARGET_OBJECTS:core_obj> hat.c spectre.c)
 
 include_directories(${CMAKE_CURRENT_SOURCE_DIR})
 

auxiliary/CMakeLists.txt

@@ -7,4 +7,6 @@ cliprogram(matching matching.c)
 cliprogram(obfusc obfusc.c)
 cliprogram(penrose-test penrose-test.c)
 cliprogram(sort-test sort-test.c)
+cliprogram(spectre-gen spectre-gen.c spectre-help.c CORE_LIB)
+cliprogram(spectre-test spectre-test.c spectre-help.c)
 cliprogram(tree234-test tree234-test.c)

auxiliary/spectre-gen.c

@@ -0,0 +1,709 @@
+/*
+ * Generate the lookup tables used by the Spectre tiling.
+ */
+
+#include <assert.h>
+#include <errno.h>
+#include <math.h>
+#include <stdio.h>
+#include <string.h>
+
+#include "puzzles.h"
+#include "tree234.h"
+#include "spectre-internal.h"
+#include "spectre-tables-manual.h"
+#include "spectre-tables-extra.h"
+#include "spectre-help.h"
+
+struct HexData {
+    const Hex *subhexes;
+    const unsigned *orientations;
+    const int *edges;
+    Point hex_outline_start, hex_outline_direction;
+    unsigned spectre_outline_start_spec, spectre_outline_start_vertex;
+};
+
+static const struct HexData hexdata[] = {
+    #define HEXDATA_ENTRY(x) { subhexes_##x, orientations_##x, edges_##x, \
+            HEX_OUTLINE_START_##x, SPEC_OUTLINE_START_##x },
+    HEX_LETTERS(HEXDATA_ENTRY)
+    #undef HEXDATA_ENTRY
+};
+
+/*
+ * Store information about an edge of the hexagonal tiling.
+ */
+typedef struct EdgeData {
+    /* Edges are regarded as directed, so that we can store
+     * information separately about what's on each side of one. The
+     * names 'start' and 'finish' indicate a direction of travel,
+     * which is taken to be anticlockwise around a hexagon, i.e. if
+     * you walk from 'start' to 'finish' then the hexagon in question
+     * is the one on your left. */
+    Point start, finish;
+
+    /* Whether this edge is internal (i.e. owned by a hexagon). */
+    bool internal;
+
+    /*
+     * High- and low-order parts of the edge identity.
+     *
+     * If the edge is internal, then 'hi' indexes the hexagon it's an
+     * edge of, and 'lo' identifies one of its edges.
+     *
+     * If it's external, then 'hi' is the index of the edge segment
+     * corresponding to a particular edge of the superhex, and 'lo'
+     * the sub-index within that segment.
+     */
+    unsigned hi, lo;
+} EdgeData;
+
+static int edge_cmp(void *av, void *bv)
+{
+    const EdgeData *a = (const EdgeData *)av;
+    const EdgeData *b = (const EdgeData *)bv;
+    size_t i;
+
+    for (i = 0; i < 4; i++) {
+        if (a->start.coeffs[i] < b->start.coeffs[i])
+            return -1;
+        if (a->start.coeffs[i] > b->start.coeffs[i])
+            return +1;
+    }
+    for (i = 0; i < 4; i++) {
+        if (a->finish.coeffs[i] < b->finish.coeffs[i])
+            return -1;
+        if (a->finish.coeffs[i] > b->finish.coeffs[i])
+            return +1;
+    }
+    return 0;
+}
+
+static void lay_out_hexagons(Hex h, Graphics *gr, FILE *hdr)
+{
+    size_t i, j;
+    tree234 *edge_map = newtree234(edge_cmp);
+    EdgeData *edge;
+    EdgeData *intmap[48], *extmap[22];
+    unsigned edgestarts[7];
+    const struct HexData *hd = h == NO_HEX ? NULL : &hexdata[h];
+
+    /*
+     * Iterate over all hexagons and enter their edges into the edge
+     * map.
+     */
+    for (i = 0; i < (h == NO_HEX ? 8 : num_subhexes(h)); i++) {
+        Point centre = hex_centres[i];
+        Point vrel = {{ -2, 0, 4, 0 }};
+        Point vertices[6];
+
+        if (hd)
+            vrel = point_mul(vrel, point_rot(2*hd->orientations[i]));
+        for (j = 0; j < 6; j++) {
+            Point vrelnext = point_mul(vrel, point_rot(2));
+
+            edge = snew(EdgeData);
+            edge->start = point_add(centre, vrel);
+            edge->finish = point_add(centre, vrelnext);
+            edge->internal = true;
+            edge->hi = i;
+            edge->lo = j;
+            add234(edge_map, edge);
+            intmap[6*i + j] = edge;
+
+            vertices[j] = edge->start;
+
+            vrel = vrelnext;
+        }
+
+        gr_draw_hex(gr, gr->jigsaw_mode ? -1 : i,
+                    hd ? hd->subhexes[i] : NO_HEX, vertices);
+    }
+
+    /*
+     * Trace round the exterior outline of the hex expansion,
+     * following the list of edge types.
+     */
+    if (hd) {
+        Point pos, dir;
+        size_t mappos = 0;
+
+        pos = hd->hex_outline_start;
+        dir = hd->hex_outline_direction;
+
+        for (i = 0; i < 6; i++) {
+            int edge_type = hd->edges[i];
+            int sign = edge_type < 0 ? -1 : +1;
+            const int *edge_shape = hex_edge_shapes[abs(edge_type)];
+            size_t len = hex_edge_lengths[abs(edge_type)];
+            size_t index = sign < 0 ? len-2 : 0;
+
+            if (gr->vertex_blobs)
+                gr_draw_blob(gr, (i == 0 ? "startpoint" : "edgesep"),
+                             gr_logcoords(pos), (i == 0 ? 0.6 : 0.3));
+
+            edgestarts[i] = mappos;
+
+            for (j = 0; j < len; j++) {
+                Point posnext = point_add(pos, dir);
+                if (j < len-1) {
+                    dir = point_mul(dir, point_rot(sign * edge_shape[index]));
+                    index += sign;
+                }
+
+                edge = snew(EdgeData);
+                edge->start = pos;
+                edge->finish = posnext;
+                edge->internal = false;
+                edge->hi = i;
+                edge->lo = j;
+                add234(edge_map, edge);
+
+                assert(mappos < lenof(extmap));
+                extmap[mappos++] = edge;
+
+                pos = posnext;
+            }
+
+            /*
+             * In the hex expansion, every pair of edges meet at a
+             * 60-degree left turn.
+             */
+            dir = point_mul(dir, point_rot(-2));
+        }
+
+        edgestarts[i] = mappos;        /* record end position */
+
+        for (i = 0; i < 4; i++)
+            assert(pos.coeffs[i] == hd->hex_outline_start.coeffs[i]);
+    }
+
+    /*
+     * Draw the labels on the edges.
+     */
+    if (gr->number_edges) {
+        for (i = 0; (edge = index234(edge_map, i)) != NULL; i++) {
+            char buf[64];
+            double textheight = 0.8, offset = textheight * 0.2;
+            GrCoords start = gr_logcoords(edge->start);
+            GrCoords finish = gr_logcoords(edge->finish);
+            GrCoords len = { finish.x - start.x, finish.y - start.y };
+            GrCoords perp = { -len.y, +len.x };
+            GrCoords mid = { (start.x+finish.x)/2, (start.y+finish.y)/2 };
+
+            if (edge->internal) {
+                sprintf(buf, "%u", edge->lo);
+            } else {
+                sprintf(buf, "%u.%u", edge->lo, edge->hi);
+                offset = textheight * 0.3;
+            }
+
+            {
+                GrCoords pos = {
+                    mid.x + offset * perp.x,
+                    mid.y + offset * perp.y,
+                };
+                gr_draw_text(gr, pos, textheight, buf);
+            }
+        }
+    }
+
+    /*
+     * Write out C array declarations for the machine-readable version
+     * of the maps we just generated.
+     */
+    if (hdr) {
+        fprintf(hdr, "static const struct MapEntry hexmap_%s[] = {\n",
+                hex_names[h]);
+        for (i = 0; i < 6 * num_subhexes(h); i++) {
+            EdgeData *our_edge = intmap[i];
+            EdgeData key, *rev_edge;
+            key.finish = our_edge->start;
+            key.start = our_edge->finish;
+            rev_edge = find234(edge_map, &key, NULL);
+            assert(rev_edge);
+            fprintf(hdr, "    { %-6s %u, %u }, /* edge %u of hex %u (%s) */\n",
+                    rev_edge->internal ? "true," : "false,",
+                    rev_edge->hi, rev_edge->lo,
+                    our_edge->lo, our_edge->hi,
+                    hex_names[hd->subhexes[our_edge->hi]]);
+        }
+        fprintf(hdr, "};\n");
+
+        fprintf(hdr, "static const struct MapEdge hexedges_%s[] = {\n",
+                hex_names[h]);
+        for (i = 0; i < 6; i++)
+            fprintf(hdr, "    { %2u, %u },\n", edgestarts[i],
+                    edgestarts[i+1] - edgestarts[i]);
+        fprintf(hdr, "};\n");
+
+        fprintf(hdr, "static const struct MapEntry hexin_%s[] = {\n",
+                hex_names[h]);
+        for (i = 0; i < edgestarts[6]; i++) {
+            EdgeData *our_edge = extmap[i];
+            EdgeData key, *rev_edge;
+            key.finish = our_edge->start;
+            key.start = our_edge->finish;
+            rev_edge = find234(edge_map, &key, NULL);
+            assert(rev_edge);
+            fprintf(hdr, "    { %-6s %u, %u }, /* subedge %u of edge %u */\n",
+                    rev_edge->internal ? "true," : "false,",
+                    rev_edge->hi, rev_edge->lo,
+                    our_edge->lo, our_edge->hi);
+        }
+        fprintf(hdr, "};\n");
+    }
+
+    while ((edge = delpos234(edge_map, 0)) != NULL)
+        sfree(edge);
+    freetree234(edge_map);
+}
+
+static void lay_out_spectres(Hex h, Graphics *gr, FILE *hdr)
+{
+    size_t i, j;
+    tree234 *edge_map = newtree234(edge_cmp);
+    EdgeData *edge;
+    EdgeData *intmap[28], *extmap[24];
+    Point vertices[28];
+    unsigned edgestarts[7];
+    const struct HexData *hd = (h == NO_HEX ? NULL : &hexdata[h]);
+
+    /*
+     * Iterate over the Spectres in a hex (usually only one), and enter
+     * their edges into the edge map.
+     */
+    for (i = 0; i < (h == NO_HEX ? 2 : num_spectres(h)); i++) {
+        Point start = {{ 0, 0, 0, 0 }};
+        Point pos = start;
+        Point diag = {{ 2, 0, 0, 2 }};
+        Point dir = point_mul(diag, point_rot(5));
+
+        /*
+         * Usually the single Spectre in each map is oriented in the
+         * same place. For spectre #1 in the G map, however, we orient
+         * it manually in a different location. (There's no point
+         * making an organised lookup table for just this one
+         * exceptional case.)
+         */
+        if (i == 1) {
+            Point unusual_start = {{ 2, 6, 2, 0 }};
+            pos = unusual_start;
+            dir = point_mul(dir, point_rot(+1));
+        }
+
+        for (j = 0; j < 14; j++) {
+            edge = snew(EdgeData);
+            edge->start = pos;
+            edge->finish = point_add(pos, dir);
+            edge->internal = true;
+            edge->hi = i;
+            edge->lo = j;
+            add234(edge_map, edge);
+            intmap[14*i + j] = edge;
+
+            vertices[14*i + j] = edge->start;
+
+            pos = edge->finish;
+            dir = point_mul(dir, point_rot(spectre_angles[(j+1) % 14]));
+        }
+
+        gr_draw_spectre(gr, h, i, vertices + 14*i);
+    }
+
+    /*
+     * Trace round the exterior outline of the hex expansion,
+     * following the list of edge types. Due to the confusing
+     * reflection of all the expansions, we end up doing this in the
+     * reverse order to the hexes code above.
+     */
+    if (hd) {
+        Point start, pos, dir;
+        size_t mappos = lenof(extmap);
+
+        start = pos = vertices[14 * hd->spectre_outline_start_spec +
+                               hd->spectre_outline_start_vertex];
+
+        edgestarts[6] = mappos;
+
+        for (i = 0; i < 6; i++) {
+            int edge_type = hd->edges[5-i];
+            int sign = edge_type < 0 ? -1 : +1;
+            const int *edge_shape = spec_edge_shapes[abs(edge_type)];
+            size_t len = spec_edge_lengths[abs(edge_type)];
+            size_t index = sign < 0 ? len-2 : 0;
+
+            if (gr->vertex_blobs)
+                gr_draw_blob(gr, (i == 0 ? "startpoint" : "edgesep"),
+                             gr_logcoords(pos), (i == 0 ? 0.6 : 0.3));
+
+            if (h == HEX_S && i >= 4) {
+                /*
+                 * Two special cases
+                 */
+                if (i == 4)
+                    /* leave dir from last time */;
+                else
+                    dir = point_mul(dir, point_rot(6)); /* reverse */
+            } else {
+                /*
+                 * Determine the direction of the first sub-edge of
+                 * this edge expansion, by iterating over all the
+                 * edges in edge_map starting at this point and
+                 * finding one whose reverse isn't in the map (hence,
+                 * it's an exterior edge).
+                 */
+                EdgeData dummy, *iter, *found = NULL;
+                dummy.start = pos;
+                for (j = 0; j < 4; j++)
+                    dummy.finish.coeffs[j] = INT_MIN;
+                for (iter = findrel234(edge_map, &dummy, NULL, REL234_GE);
+                     iter != NULL && point_equal(iter->start, pos);
+                     iter = findrel234(edge_map, iter, NULL, REL234_GT)) {
+                    EdgeData *rev;
+
+                    dummy.finish = iter->start;
+                    dummy.start = iter->finish;
+                    rev = find234(edge_map, &dummy, NULL);
+                    if (!rev) {
+                        found = iter;
+                        break;
+                    }
+                }
+
+                assert(found);
+                dir = point_sub(found->finish, found->start);
+            }
+
+            for (j = 0; j < len; j++) {
+                Point posnext = point_add(pos, dir);
+                if (j < len-1) {
+                    dir = point_mul(dir, point_rot(sign * edge_shape[index]));
+                    index += sign;
+                }
+
+                edge = snew(EdgeData);
+                edge->start = posnext;
+                edge->finish = pos;
+                edge->internal = false;
+                edge->hi = 5-i;
+                edge->lo = len-1-j;
+                add234(edge_map, edge);
+
+                assert(mappos > 0);
+                extmap[--mappos] = edge;
+
+                pos = posnext;
+            }
+
+            edgestarts[5-i] = mappos;
+        }
+
+        assert(point_equal(pos, start));
+    }
+
+    /*
+     * Draw the labels on the edges.
+     */
+    if (gr->number_edges) {
+        for (i = 0; (edge = index234(edge_map, i)) != NULL; i++) {
+            char buf[64];
+            double textheight = 0.8, offset = textheight * 0.2;
+            GrCoords start = gr_logcoords(edge->start);
+            GrCoords finish = gr_logcoords(edge->finish);
+            GrCoords len = { finish.x - start.x, finish.y - start.y };
+            GrCoords perp = { +len.y, -len.x };
+            GrCoords mid = { (start.x+finish.x)/2, (start.y+finish.y)/2 };
+
+            if (edge->internal) {
+                sprintf(buf, "%u", edge->lo);
+            } else {
+                sprintf(buf, "%u.%u", edge->lo, edge->hi);
+                textheight = 0.6;
+            }
+            if (strlen(buf) > 1)
+                offset = textheight * 0.35;
+
+            {
+                GrCoords pos = {
+                    mid.x + offset * perp.x,
+                    mid.y + offset * perp.y,
+                };
+                gr_draw_text(gr, pos, textheight, buf);
+            }
+        }
+    }
+
+    /*
+     * Write out C array declarations for the machine-readable version
+     * of the maps we just generated.
+     *
+     * Also, because it's easier than having a whole extra iteration,
+     * draw lines for the extraordinary edges outside the S diagram.
+     */
+    if (hdr) {
+        fprintf(hdr, "static const struct MapEntry specmap_%s[] = {\n",
+                hex_names[h]);
+        for (i = 0; i < 14 * num_spectres(h); i++) {
+            EdgeData *our_edge = intmap[i];
+            EdgeData key, *rev_edge;
+            key.finish = our_edge->start;
+            key.start = our_edge->finish;
+            rev_edge = find234(edge_map, &key, NULL);
+            assert(rev_edge);
+            fprintf(hdr, "    { %-6s %u, %2u }, /* edge %2u of Spectre %u */\n",
+                    rev_edge->internal ? "true," : "false,",
+                    rev_edge->hi, rev_edge->lo,
+                    our_edge->lo, our_edge->hi);
+        }
+        fprintf(hdr, "};\n");
+
+        fprintf(hdr, "static const struct MapEdge specedges_%s[] = {\n",
+                hex_names[h]);
+        for (i = 0; i < 6; i++)
+            fprintf(hdr, "    { %2u, %u },\n", edgestarts[i] - edgestarts[0],
+                    edgestarts[i+1] - edgestarts[i]);
+        fprintf(hdr, "};\n");
+
+        fprintf(hdr, "static const struct MapEntry specin_%s[] = {\n",
+                hex_names[h]);
+        for (i = edgestarts[0]; i < edgestarts[6]; i++) {
+            EdgeData *our_edge = extmap[i];
+            EdgeData key, *rev_edge;
+            key.finish = our_edge->start;
+            key.start = our_edge->finish;
+            rev_edge = find234(edge_map, &key, NULL);
+            assert(rev_edge);
+            fprintf(hdr, "    { %-6s %u, %2u }, /* subedge %u of edge %u */\n",
+                    rev_edge->internal ? "true," : "false,",
+                    rev_edge->hi, rev_edge->lo,
+                    our_edge->lo, our_edge->hi);
+
+            if (!our_edge->internal && !rev_edge->internal)
+                gr_draw_extra_edge(gr, key.start, key.finish);
+        }
+        fprintf(hdr, "};\n");
+    }
+
+    while ((edge = delpos234(edge_map, 0)) != NULL)
+        sfree(edge);
+    freetree234(edge_map);
+}
+
+static void draw_base_hex(Hex h, Graphics *gr)
+{
+    size_t i;
+    Point vertices[6];
+
+    /*
+     * Plot the points of the hex.
+     */
+    for (i = 0; i < 6; i++) {
+        Point startvertex = {{ -2, 0, 4, 0 }};
+        vertices[i] = point_mul(startvertex, point_rot(2*i));
+    }
+
+    /*
+     * Draw the hex itself.
+     */
+    gr_draw_hex(gr, -1, h, vertices);
+
+    if (gr->vertex_blobs) {
+        /*
+         * Draw edge-division blobs on all vertices, to match the ones on
+         * the expansion diagrams.
+         */
+        for (i = 0; i < 6; i++) {
+            gr_draw_blob(gr, (i == 0 ? "startpoint" : "edgesep"),
+                         gr_logcoords(vertices[i]), (i == 0 ? 0.6 : 0.3));
+        }
+    }
+
+    if (gr->number_edges) {
+        /*
+         * Draw the labels on its edges.
+         */
+        for (i = 0; i < 6; i++) {
+            char buf[64];
+            double textheight = 0.8, offset = textheight * 0.2;
+            GrCoords start = gr_logcoords(vertices[i]);
+            GrCoords finish = gr_logcoords(vertices[(i+1) % 6]);
+            GrCoords len = { finish.x - start.x, finish.y - start.y };
+            GrCoords perp = { -len.y, +len.x };
+            GrCoords mid = { (start.x+finish.x)/2, (start.y+finish.y)/2 };
+
+            sprintf(buf, "%zu", i);
+
+            {
+                GrCoords pos = {
+                    mid.x + offset * perp.x,
+                    mid.y + offset * perp.y,
+                };
+                gr_draw_text(gr, pos, textheight, buf);
+            }
+        }
+    }
+}
+
+static void draw_one_spectre(Graphics *gr)
+{
+    size_t i, j;
+    Point vertices[14];
+
+    {
+        Point start = {{ 0, 0, 0, 0 }};
+        Point pos = start;
+        Point diag = {{ 2, 0, 0, 2 }};
+        Point dir = point_mul(diag, point_rot(9));
+
+        for (j = 0; j < 14; j++) {
+            vertices[j] = pos;
+            pos = point_add(pos, dir);
+            dir = point_mul(dir, point_rot(spectre_angles[(j+1) % 14]));
+        }
+
+        gr_draw_spectre(gr, NO_HEX, -1, vertices);
+    }
+
+    /*
+     * Draw the labels on the edges.
+     */
+    if (gr->number_edges) {
+        for (i = 0; i < 14; i++) {
+            char buf[64];
+            double textheight = 0.8, offset = textheight * 0.2;
+            GrCoords start = gr_logcoords(vertices[i]);
+            GrCoords finish = gr_logcoords(vertices[(i+1) % 14]);
+            GrCoords len = { finish.x - start.x, finish.y - start.y };
+            GrCoords perp = { +len.y, -len.x };
+            GrCoords mid = { (start.x+finish.x)/2, (start.y+finish.y)/2 };
+
+            sprintf(buf, "%zu", i);
+            if (strlen(buf) > 1)
+                offset = textheight * 0.35;
+
+            {
+                GrCoords pos = {
+                    mid.x + offset * perp.x,
+                    mid.y + offset * perp.y,
+                };
+                gr_draw_text(gr, pos, textheight, buf);
+            }
+        }
+    }
+
+}
+
+static void make_parent_tables(FILE *fp)
+{
+    size_t i, j, k;
+
+    for (i = 0; i < 9; i++) {
+        fprintf(fp, "static const struct Possibility poss_%s[] = {\n",
+                hex_names[i]);
+        for (j = 0; j < 9; j++) {
+            for (k = 0; k < num_subhexes(j); k++) {
+                if (hexdata[j].subhexes[k] == i) {
+                    fprintf(fp, "    { HEX_%s, %zu, PROB_%s },\n",
+                            hex_names[j], k, hex_names[j]);
+                }
+            }
+        }
+        fprintf(fp, "};\n");
+    }
+
+    fprintf(fp, "static const struct Possibility poss_spectre[] = {\n");
+    for (j = 0; j < 9; j++) {
+        for (k = 0; k < num_spectres(j); k++) {
+            fprintf(fp, "    { HEX_%s, %zu, PROB_%s },\n",
+                    hex_names[j], k, hex_names[j]);
+        }
+    }
+    fprintf(fp, "};\n");
+}
+
+int main(void)
+{
+    size_t i;
+    FILE *fp = fopen("spectre-tables-auto.h", "w");
+    fprintf(fp,
+            "/*\n"
+            " * Autogenerated transition tables for the Spectre tiling.\n"
+            " * Generated by auxiliary/spectre-gen.c.\n"
+            " */\n\n");
+
+    for (i = 0; i < 9; i++) {
+        char buf[64];
+        sprintf(buf, "hexmap_%s.svg", hex_names[i]);
+        Graphics *gr = gr_new(buf, -11, +11, -20, +4.5, 13);
+        lay_out_hexagons(i, gr, fp);
+        gr_free(gr);
+    }
+    for (i = 0; i < 9; i++) {
+        char buf[64];
+        sprintf(buf, "specmap_%s.svg", hex_names[i]);
+        Graphics *gr = gr_new(buf, (i == HEX_S ? -14 : -11.5),
+                              (i == HEX_G ? +10 : 0.5),
+                              -2, +12, 15);
+        lay_out_spectres(i, gr, fp);
+        gr_free(gr);
+    }
+    for (i = 0; i < 9; i++) {
+        char buf[64];
+        sprintf(buf, "basehex_%s.svg", hex_names[i]);
+        Graphics *gr = gr_new(buf, -4, +4, -4.2, +4.5, 15);
+        draw_base_hex(i, gr);
+        gr_free(gr);
+    }
+    for (i = 0; i < 9; i++) {
+        char buf[64];
+        sprintf(buf, "jigsawhex_%s.svg", hex_names[i]);
+        Graphics *gr = gr_new(buf, -4, +4, -4.2, +4.5, 20);
+        gr->jigsaw_mode = true;
+        gr->vertex_blobs = false;
+        gr->number_edges = false;
+        draw_base_hex(i, gr);
+        gr_free(gr);
+    }
+    {
+        Graphics *gr = gr_new("basehex_null.svg", -4, +4, -4.2, +4.5, 20);
+        gr->vertex_blobs = false;
+        draw_base_hex(NO_HEX, gr);
+        gr_free(gr);
+    }
+    {
+        Graphics *gr = gr_new("basespec_null.svg", -7, +6, -14, +1, 15);
+        gr->vertex_blobs = false;
+        draw_one_spectre(gr);
+        gr_free(gr);
+    }
+    {
+        Graphics *gr = gr_new("hexmap_null.svg", -11, +11, -20, +4.5, 10);
+        gr->vertex_blobs = false;
+        gr->number_edges = false;
+        gr->hex_arrows = false;
+        lay_out_hexagons(NO_HEX, gr, NULL);
+        gr_free(gr);
+    }
+    {
+        Graphics *gr = gr_new("specmap_null.svg", -11.5, +10, -2, +12, 15);
+        gr->vertex_blobs = false;
+        gr->number_edges = false;
+        gr->hex_arrows = false;
+        lay_out_spectres(NO_HEX, gr, NULL);
+        gr_free(gr);
+    }
+    for (i = 0; i < 2; i++) {
+        char buf[64];
+        sprintf(buf, "jigsawexpand_%s.svg", hex_names[i]);
+        Graphics *gr = gr_new(buf, -11, +11, -20, +4.5, 10);
+        gr->jigsaw_mode = true;
+        gr->vertex_blobs = false;
+        gr->number_edges = false;
+        lay_out_hexagons(i, gr, fp);
+        gr_free(gr);
+    }
+    make_parent_tables(fp);
+    fclose(fp);
+    return 0;
+}

auxiliary/spectre-help.c

@@ -0,0 +1,417 @@
+/*
+ * Common code between spectre-test and spectre-gen, since both of
+ * them want to output SVG graphics.
+ */
+
+#include <assert.h>
+#include <errno.h>
+#include <math.h>
+#include <stdio.h>
+#include <string.h>
+
+#include "puzzles.h"
+#include "tree234.h"
+#include "spectre-internal.h"
+#include "spectre-tables-extra.h"
+#include "spectre-help.h"
+
+struct HexData {
+    const int *edges;
+};
+
+static const struct HexData hexdata[] = {
+    #define HEXDATA_ENTRY(x) { edges_##x },
+    HEX_LETTERS(HEXDATA_ENTRY)
+    #undef HEXDATA_ENTRY
+};
+
+const char *hex_names[10] = {
+    "G", "D", "J", "L", "X", "P", "S", "F", "Y",
+    "" /* NO_HEX */
+};
+
+Graphics *gr_new(const char *filename, double xmin, double xmax,
+                 double ymin, double ymax, double scale)
+{
+    Graphics *gr = snew(Graphics);
+    if (!strcmp(filename, "-")) {
+        gr->fp = stdout;
+        gr->close_file = false;
+    } else {
+        gr->fp = fopen(filename, "w");
+        if (!gr->fp) {
+            fprintf(stderr, "%s: open: %s\n", filename, strerror(errno));
+            exit(1);
+        }
+        gr->close_file = true;
+    }
+
+    fprintf(gr->fp, "<?xml version=\"1.0\" encoding=\"UTF-8\" "
+            "standalone=\"no\"?>\n");
+    fprintf(gr->fp, "<svg xmlns=\"http://www.w3.org/2000/svg\" "
+            "version=\"1.1\" width=\"%f\" height=\"%f\">\n",
+            (xmax - xmin) * scale, (ymax - ymin) * scale);
+
+    gr->absscale = fabs(scale);
+    gr->xoff = -xmin * scale;
+    gr->xscale = scale;
+    /* invert y axis for SVG top-down coordinate system */
+    gr->yoff = ymax * scale;
+    gr->yscale = -scale;
+
+    /* Defaults, which can be overridden by the caller immediately
+     * after this constructor returns */
+    gr->jigsaw_mode = false;
+    gr->vertex_blobs = true;
+    gr->number_cells = true;
+    gr->four_colour = false;
+    gr->arcs = false;
+    gr->linewidth = 1.5;
+
+    gr->started = false;
+
+    return gr;
+}
+
+void gr_free(Graphics *gr)
+{
+    if (!gr)
+        return;
+    fprintf(gr->fp, "</svg>\n");
+    if (gr->close_file)
+        fclose(gr->fp);
+    sfree(gr);
+}
+
+static void gr_ensure_started(Graphics *gr)
+{
+    if (gr->started)
+        return;
+
+    fprintf(gr->fp, "<style type=\"text/css\">\n");
+    fprintf(gr->fp, "path { fill: none; stroke: black; stroke-width: %f; "
+            "stroke-linejoin: round; stroke-linecap: round; }\n",
+            gr->linewidth);
+    fprintf(gr->fp, "text { fill: black; font-family: Sans; "
+            "text-anchor: middle; text-align: center; }\n");
+    if (gr->four_colour) {
+        fprintf(gr->fp, ".c0 { fill: rgb(255, 178, 178); }\n");
+        fprintf(gr->fp, ".c1 { fill: rgb(255, 255, 178); }\n");
+        fprintf(gr->fp, ".c2 { fill: rgb(178, 255, 178); }\n");
+        fprintf(gr->fp, ".c3 { fill: rgb(153, 153, 255); }\n");
+    } else {
+        fprintf(gr->fp, ".G { fill: rgb(255, 128, 128); }\n");
+        fprintf(gr->fp, ".G1 { fill: rgb(255, 64, 64); }\n");
+        fprintf(gr->fp, ".F { fill: rgb(255, 192, 128); }\n");
+        fprintf(gr->fp, ".Y { fill: rgb(255, 255, 128); }\n");
+        fprintf(gr->fp, ".S { fill: rgb(128, 255, 128); }\n");
+        fprintf(gr->fp, ".D { fill: rgb(128, 255, 255); }\n");
+        fprintf(gr->fp, ".P { fill: rgb(128, 128, 255); }\n");
+        fprintf(gr->fp, ".X { fill: rgb(192, 128, 255); }\n");
+        fprintf(gr->fp, ".J { fill: rgb(255, 128, 255); }\n");
+        fprintf(gr->fp, ".L { fill: rgb(128, 128, 128); }\n");
+        fprintf(gr->fp, ".optional { stroke-dasharray: 5; }\n");
+        fprintf(gr->fp, ".arrow { fill: rgba(0, 0, 0, 0.2); "
+                "stroke: none; }\n");
+    }
+    fprintf(gr->fp, "</style>\n");
+
+    gr->started = true;
+}
+
+/* Logical coordinates in our mathematical space */
+GrCoords gr_logcoords(Point p)
+{
+    double rt3o2 = sqrt(3) / 2;
+    GrCoords r = {
+        p.coeffs[0] + rt3o2 * p.coeffs[1] + 0.5 * p.coeffs[2],
+        p.coeffs[3] + rt3o2 * p.coeffs[2] + 0.5 * p.coeffs[1],
+    };
+    return r;
+}
+
+/* Physical coordinates in the output image */
+GrCoords gr_log2phys(Graphics *gr, GrCoords c)
+{
+    c.x = gr->xoff + gr->xscale * c.x;
+    c.y = gr->yoff + gr->yscale * c.y;
+    return c;
+}
+GrCoords gr_physcoords(Graphics *gr, Point p)
+{
+    return gr_log2phys(gr, gr_logcoords(p));
+}
+
+void gr_draw_text(Graphics *gr, GrCoords logpos, double logheight,
+                  const char *text)
+{
+    GrCoords pos;
+    double height;
+
+    if (!gr)
+        return;
+    gr_ensure_started(gr);
+
+    pos = gr_log2phys(gr, logpos);
+    height = gr->absscale * logheight;
+    fprintf(gr->fp, "<text style=\"font-size: %fpx\" x=\"%f\" y=\"%f\">"
+            "%s</text>\n", height, pos.x, pos.y + 0.35 * height, text);
+}
+
+void gr_draw_path(Graphics *gr, const char *classes, const GrCoords *phys,
+                  size_t n, bool closed)
+{
+    size_t i;
+
+    if (!gr)
+        return;
+    gr_ensure_started(gr);
+
+    fprintf(gr->fp, "<path class=\"%s\" d=\"", classes);
+    for (i = 0; i < n; i++) {
+        GrCoords c = phys[i];
+        if (i == 0)
+            fprintf(gr->fp, "M %f %f", c.x, c.y);
+        else if (gr->arcs)
+            fprintf(gr->fp, "A %f %f 10 0 %zu %f %f",
+                    gr->absscale, gr->absscale, i&1, c.x, c.y);
+        else
+            fprintf(gr->fp, "L %f %f", c.x, c.y);
+    }
+    if (gr->arcs) {
+        /* Explicitly return to the starting point so as to curve the
+         * final edge */
+        fprintf(gr->fp, "A %f %f 10 0 0 %f %f",
+                gr->absscale, gr->absscale, phys[0].x, phys[0].y);
+    }
+    if (closed)
+        fprintf(gr->fp, " z");
+    fprintf(gr->fp, "\"/>\n");
+}
+
+void gr_draw_blob(Graphics *gr, const char *classes, GrCoords log,
+                  double logradius)
+{
+    GrCoords centre;
+
+    if (!gr)
+        return;
+    gr_ensure_started(gr);
+
+    centre = gr_log2phys(gr, log);
+    fprintf(gr->fp, "<circle class=\"%s\" cx=\"%f\" cy=\"%f\" r=\"%f\"/>\n",
+            classes, centre.x, centre.y, gr->absscale * logradius);
+}
+
+void gr_draw_hex(Graphics *gr, unsigned index, Hex htype,
+                 const Point *vertices)
+{
+    size_t i;
+    Point centre;
+
+    if (!gr)
+        return;
+    gr_ensure_started(gr);
+
+    /* Draw the actual hexagon, in its own colour */
+    if (!gr->jigsaw_mode) {
+        GrCoords phys[6];
+        for (i = 0; i < 6; i++)
+            phys[i] = gr_physcoords(gr, vertices[i]);
+        gr_draw_path(gr, (index == 7 && htype == NO_HEX ?
+                          "optional" : hex_names[htype]), phys, 6, true);
+    } else {
+        GrCoords phys[66];
+        size_t pos = 0;
+        const struct HexData *hd = &hexdata[htype];
+
+        for (i = 0; i < 6; i++) {
+            int edge_type = hd->edges[i];
+            int sign = edge_type < 0 ? -1 : +1;
+            int edge_abs = abs(edge_type);
+            int left_sign = (edge_abs & 4) ? sign : edge_type == 0 ? +1 : 0;
+            int mid_sign = (edge_abs & 2) ? sign : 0;
+            int right_sign = (edge_abs & 1) ? sign : edge_type == 0 ? -1 : 0;
+
+            GrCoords start = gr_physcoords(gr, vertices[i]);
+            GrCoords end = gr_physcoords(gr, vertices[(i+1) % 6]);
+            GrCoords x = { (end.x - start.x) / 7, (end.y - start.y) / 7 };
+            GrCoords y = { -x.y, +x.x };
+
+#define addpoint(X, Y) do {                             \
+                GrCoords p = {                          \
+                    start.x + (X) * x.x + (Y) * y.x,    \
+                    start.y + (X) * x.y + (Y) * y.y,    \
+                };                                      \
+                phys[pos++] = p;                        \
+            } while (0)
+
+            if (sign < 0) {
+                int tmp = right_sign;
+                right_sign = left_sign;
+                left_sign = tmp;
+            }
+
+            addpoint(0, 0);
+            if (left_sign) {
+                addpoint(1, 0);
+                addpoint(2, left_sign);
+                addpoint(2, 0);
+            }
+            if (mid_sign) {
+                addpoint(3, 0);
+                addpoint(3, mid_sign);
+                addpoint(4, mid_sign);
+                addpoint(4, 0);
+            }
+            if (right_sign) {
+                addpoint(5, 0);
+                addpoint(5, right_sign);
+                addpoint(6, 0);
+            }
+
+#undef addpoint
+
+        }
+        gr_draw_path(gr, hex_names[htype], phys, pos, true);
+    }
+
+    /* Find the centre of the hex */
+    for (i = 0; i < 4; i++)
+        centre.coeffs[i] = 0;
+    for (i = 0; i < 6; i++)
+        centre = point_add(centre, vertices[i]);
+    for (i = 0; i < 4; i++)
+        centre.coeffs[i] /= 6;
+
+    /* Draw an arrow towards vertex 0 of the hex */
+    if (gr->hex_arrows) {
+        double ext = 0.6;
+        double headlen = 0.3, thick = 0.08, headwid = 0.25;
+        GrCoords top = gr_physcoords(gr, vertices[0]);
+        GrCoords bot = gr_physcoords(gr, vertices[3]);
+        GrCoords mid = gr_physcoords(gr, centre);
+        GrCoords base = { mid.x + ext * (bot.x - mid.x),
+                          mid.y + ext * (bot.y - mid.y) };
+        GrCoords tip = { mid.x + ext * (top.x - mid.x),
+                         mid.y + ext * (top.y - mid.y) };
+        GrCoords len = { tip.x - base.x, tip.y - base.y };
+        GrCoords perp = { -len.y, +len.x };
+        GrCoords basep = { base.x+perp.x*thick, base.y+perp.y*thick };
+        GrCoords basen = { base.x-perp.x*thick, base.y-perp.y*thick };
+        GrCoords hbase = { tip.x-len.x*headlen, tip.y-len.y*headlen };
+        GrCoords headp = { hbase.x+perp.x*thick, hbase.y+perp.y*thick };
+        GrCoords headn = { hbase.x-perp.x*thick, hbase.y-perp.y*thick };
+        GrCoords headP = { hbase.x+perp.x*headwid, hbase.y+perp.y*headwid };
+        GrCoords headN = { hbase.x-perp.x*headwid, hbase.y-perp.y*headwid };
+
+        GrCoords phys[] = {
+            basep, headp, headP, tip, headN, headn, basen
+        };
+
+        gr_draw_path(gr, "arrow", phys, lenof(phys), true);
+    }
+
+    /*
+     * Label the hex with its index and type.
+     */
+    if (gr->number_cells) {
+        char buf[64];
+        if (index == (unsigned)-1) {
+            if (htype == NO_HEX)
+                buf[0] = '\0';
+            else
+                strcpy(buf, hex_names[htype]);
+        } else {
+            if (htype == NO_HEX)
+                sprintf(buf, "%u", index);
+            else
+                sprintf(buf, "%u (%s)", index, hex_names[htype]);
+        }
+        if (buf[0])
+            gr_draw_text(gr, gr_logcoords(centre), 1.2, buf);
+    }
+}
+
+void gr_draw_spectre(Graphics *gr, Hex container, unsigned index,
+                     const Point *vertices)
+{
+    size_t i;
+    GrCoords log[14];
+    GrCoords centre;
+
+    if (!gr)
+        return;
+    gr_ensure_started(gr);
+
+    for (i = 0; i < 14; i++)
+        log[i] = gr_logcoords(vertices[i]);
+
+    /* Draw the actual Spectre */
+    {
+        GrCoords phys[14];
+        char class[16];
+        for (i = 0; i < 14; i++)
+            phys[i] = gr_log2phys(gr, log[i]);
+        if (gr->four_colour) {
+            sprintf(class, "c%u", index);
+        } else if (index == 1 && container == NO_HEX) {
+            sprintf(class, "optional");
+        } else {
+            sprintf(class, "%s%.0u", hex_names[container], index);
+        }
+        gr_draw_path(gr, class, phys, 14, true);
+    }
+
+    /* Pick a point to use as the centre of the Spectre for labelling */
+    centre.x = (log[5].x + log[6].x + log[11].x + log[12].x) / 4;
+    centre.y = (log[5].y + log[6].y + log[11].y + log[12].y) / 4;
+
+    /*
+     * Label the hex with its index and type.
+     */
+    if (gr->number_cells && index != (unsigned)-1) {
+        char buf[64];
+        sprintf(buf, "%u", index);
+        gr_draw_text(gr, centre, 1.2, buf);
+    }
+}
+
+void gr_draw_spectre_from_coords(Graphics *gr, SpectreCoords *sc,
+                                 const Point *vertices)
+{
+    Hex h;
+    unsigned index;
+
+    if (!gr)
+        return;
+    gr_ensure_started(gr);
+
+    if (gr->four_colour) {
+        h = NO_HEX;
+        if (sc->index == 1)
+            index = 3;        /* special colour for odd G1 Spectres */
+        else
+            index = sc->hex_colour;
+    } else if (sc) {
+        h = sc->c[0].type;
+        index = sc->index;
+    } else {
+        h = NO_HEX;
+        index = -1;
+    }
+    gr_draw_spectre(gr, h, index, vertices);
+}
+
+void gr_draw_extra_edge(Graphics *gr, Point a, Point b)
+{
+    GrCoords phys[2];
+
+    if (!gr)
+        return;
+    gr_ensure_started(gr);
+
+    phys[0] = gr_physcoords(gr, a);
+    phys[1] = gr_physcoords(gr, b);
+    gr_draw_path(gr, "extraedge", phys, 2, false);
+}

auxiliary/spectre-help.h

@@ -0,0 +1,51 @@
+/*
+ * Header for spectre-help.c
+ */
+
+/* Dummy value indicating no specific hexagon, used in some diagrams
+ * for the accompanying article. */
+#define NO_HEX (Hex)9
+
+/*
+ * String constants for the hex names, including an extra entry
+ * mapping NO_HEX to the empty string.
+ */
+extern const char *hex_names[10];
+
+typedef struct Graphics {
+    FILE *fp;
+    bool close_file; /* if it's not stdout */
+    bool started;    /* have we written the header yet? */
+    double xoff, xscale, yoff, yscale, absscale, linewidth;
+    bool jigsaw_mode; /* draw protrusions on hex edges */
+    bool vertex_blobs; /* draw blobs marking hex vertices */
+    bool hex_arrows; /* draw arrows orienting each hex */
+    bool number_edges; /* number the edges of everything */
+    bool number_cells; /* number the things themselves */
+    bool four_colour;  /* four-colour Spectres instead of semantically */
+    bool arcs; /* draw Spectre edges as arcs */
+} Graphics;
+
+typedef struct GrCoords {
+    double x, y;
+} GrCoords;
+
+Graphics *gr_new(const char *filename, double xmin, double xmax,
+                 double ymin, double ymax, double scale);
+void gr_free(Graphics *gr);
+GrCoords gr_logcoords(Point p);
+GrCoords gr_log2phys(Graphics *gr, GrCoords c);
+GrCoords gr_physcoords(Graphics *gr, Point p);
+void gr_draw_text(Graphics *gr, GrCoords logpos, double logheight,
+                  const char *text);
+void gr_draw_path(Graphics *gr, const char *classes, const GrCoords *phys,
+                  size_t n, bool closed);
+void gr_draw_blob(Graphics *gr, const char *classes, GrCoords log,
+                  double logradius);
+void gr_draw_hex(Graphics *gr, unsigned index, Hex htype,
+                 const Point *vertices);
+void gr_draw_spectre(Graphics *gr, Hex container, unsigned index,
+                     const Point *vertices);
+void gr_draw_spectre_from_coords(Graphics *gr, SpectreCoords *sc,
+                                 const Point *vertices);
+void gr_draw_extra_edge(Graphics *gr, Point a, Point b);

auxiliary/spectre-tables-extra.h

@@ -0,0 +1,334 @@
+/*
+ * Further data tables used to generate the final transition maps.
+ */
+
+/*
+ * Locations in the plane of the centres of the 8 hexagons in the
+ * expansion of each hex.
+ *
+ * We take the centre-to-centre distance to be 6 units, so that other
+ * locations in the hex tiling (e.g. edge midpoints and vertices) will
+ * still have integer coefficients.
+ *
+ * These locations are represented using the same Point type used for
+ * the whole tiling, but all our angles are 60 degrees, so we don't
+ * ever need the coefficients of d or d^3, only of 1 and d^2.
+ */
+static const Point hex_centres[] = {
+    {{0, 0, 0, 0}}, {{6, 0, 0, 0}},                        /*   0 1 */
+    {{0, 0, -6, 0}}, {{6, 0, -6, 0}},                      /*  2 3  */
+    {{0, 0, -12, 0}}, {{6, 0, -12, 0}}, {{12, 0, -12, 0}}, /* 4 5 6 */
+    {{12, 0, -18, 0}},                                     /*    7  */
+};
+
+/*
+ * Orientations of all the sub-hexes in the expansion of each hex.
+ * Measured anticlockwise (that is, as a power of s) from 0, where 0
+ * means the hex is upright, with its own vertex #0 at the top.
+ */
+
+static const unsigned orientations_G[] = {
+    2, /* HEX_F */
+    1, /* HEX_X */
+    0, /* HEX_G */
+    1, /* HEX_S */
+    4, /* HEX_P */
+    5, /* HEX_D */
+    0, /* HEX_J */
+    /* hex #7 is not present for this tile */
+};
+static const unsigned orientations_D[] = {
+    2, /* HEX_F */
+    1, /* HEX_P */
+    0, /* HEX_G */
+    1, /* HEX_S */
+    4, /* HEX_X */
+    5, /* HEX_D */
+    0, /* HEX_F */
+    5, /* HEX_X */
+};
+static const unsigned orientations_J[] = {
+    2, /* HEX_F */
+    1, /* HEX_P */
+    0, /* HEX_G */
+    1, /* HEX_S */
+    4, /* HEX_Y */
+    5, /* HEX_D */
+    0, /* HEX_F */
+    5, /* HEX_P */
+};
+static const unsigned orientations_L[] = {
+    2, /* HEX_F */
+    1, /* HEX_P */
+    0, /* HEX_G */
+    1, /* HEX_S */
+    4, /* HEX_Y */
+    5, /* HEX_D */
+    0, /* HEX_F */
+    5, /* HEX_X */
+};
+static const unsigned orientations_X[] = {
+    2, /* HEX_F */
+    1, /* HEX_Y */
+    0, /* HEX_G */
+    1, /* HEX_S */
+    4, /* HEX_Y */
+    5, /* HEX_D */
+    0, /* HEX_F */
+    5, /* HEX_P */
+};
+static const unsigned orientations_P[] = {
+    2, /* HEX_F */
+    1, /* HEX_Y */
+    0, /* HEX_G */
+    1, /* HEX_S */
+    4, /* HEX_Y */
+    5, /* HEX_D */
+    0, /* HEX_F */
+    5, /* HEX_X */
+};
+static const unsigned orientations_S[] = {
+    2, /* HEX_L */
+    1, /* HEX_P */
+    0, /* HEX_G */
+    1, /* HEX_S */
+    4, /* HEX_X */
+    5, /* HEX_D */
+    0, /* HEX_F */
+    5, /* HEX_X */
+};
+static const unsigned orientations_F[] = {
+    2, /* HEX_F */
+    1, /* HEX_P */
+    0, /* HEX_G */
+    1, /* HEX_S */
+    4, /* HEX_Y */
+    5, /* HEX_D */
+    0, /* HEX_F */
+    5, /* HEX_Y */
+};
+static const unsigned orientations_Y[] = {
+    2, /* HEX_F */
+    1, /* HEX_Y */
+    0, /* HEX_G */
+    1, /* HEX_S */
+    4, /* HEX_Y */
+    5, /* HEX_D */
+    0, /* HEX_F */
+    5, /* HEX_Y */
+};
+
+/*
+ * For each hex type, indicate the point on the boundary of the
+ * expansion that corresponds to vertex 0 of the superhex. Also,
+ * indicate the initial direction we head in to go round the edge.
+ */
+#define HEX_OUTLINE_START_COMMON {{ -4, 0, -10, 0 }}, {{ +2, 0, +2, 0 }}
+#define HEX_OUTLINE_START_RARE {{ -2, 0, -14, 0 }}, {{ -2, 0, +4, 0 }}
+#define HEX_OUTLINE_START_G HEX_OUTLINE_START_COMMON
+#define HEX_OUTLINE_START_D HEX_OUTLINE_START_RARE
+#define HEX_OUTLINE_START_J HEX_OUTLINE_START_COMMON
+#define HEX_OUTLINE_START_L HEX_OUTLINE_START_COMMON
+#define HEX_OUTLINE_START_X HEX_OUTLINE_START_COMMON
+#define HEX_OUTLINE_START_P HEX_OUTLINE_START_COMMON
+#define HEX_OUTLINE_START_S HEX_OUTLINE_START_RARE
+#define HEX_OUTLINE_START_F HEX_OUTLINE_START_COMMON
+#define HEX_OUTLINE_START_Y HEX_OUTLINE_START_COMMON
+
+/*
+ * Similarly, for each hex type, indicate the point on the boundary of
+ * its Spectre expansion that corresponds to hex vertex 0.
+ *
+ * This time, it's easiest just to indicate which vertex of which
+ * sub-Spectre we take in each case, because the Spectre outlines
+ * don't take predictable turns between the edge expansions, so the
+ * routine consuming this data will have to look things up in its
+ * edgemap anyway.
+ */
+#define SPEC_OUTLINE_START_COMMON 0, 9
+#define SPEC_OUTLINE_START_RARE 0, 8
+#define SPEC_OUTLINE_START_G SPEC_OUTLINE_START_COMMON
+#define SPEC_OUTLINE_START_D SPEC_OUTLINE_START_RARE
+#define SPEC_OUTLINE_START_J SPEC_OUTLINE_START_COMMON
+#define SPEC_OUTLINE_START_L SPEC_OUTLINE_START_COMMON
+#define SPEC_OUTLINE_START_X SPEC_OUTLINE_START_COMMON
+#define SPEC_OUTLINE_START_P SPEC_OUTLINE_START_COMMON
+#define SPEC_OUTLINE_START_S SPEC_OUTLINE_START_RARE
+#define SPEC_OUTLINE_START_F SPEC_OUTLINE_START_COMMON
+#define SPEC_OUTLINE_START_Y SPEC_OUTLINE_START_COMMON
+
+/*
+ * The paper also defines a set of 8 different classes of edges for
+ * the hexagons. (You can imagine these as different shapes of
+ * jigsaw-piece tab, constraining how the hexes can fit together). So
+ * for each hex, we need a list of its edge types.
+ *
+ * Most edge types come in two matching pairs, which the paper labels
+ * with the same lowercase Greek letter and a + or - superscript, e.g.
+ * alpha^+ and alpha^-. The usual rule is that when two edges meet,
+ * they have to be the + and - versions of the same letter. The
+ * exception to this rule is the 'eta' edge, which has no sign: it's
+ * symmetric, so any two eta edges can validly meet.
+ *
+ * We express this here by defining an enumeration in which eta = 0
+ * and all other edge types have positive values, so that integer
+ * negation can be used to indicate the other edge that fits with this
+ * one (and for eta, it doesn't change the value).
+ */
+enum Edge {
+    edge_eta = 0,
+    edge_alpha,
+    edge_beta,
+    edge_gamma,
+    edge_delta,
+    edge_epsilon,
+    edge_zeta,
+    edge_theta,
+};
+
+/*
+ * Edge types for each hex are specified anticlockwise, starting from
+ * the top vertex, so that edge #0 is the top-left diagonal edge, edge
+ * #1 the left-hand vertical edge, etc.
+ */
+static const int edges_G[6] = {
+    -edge_beta, -edge_alpha, +edge_alpha,
+    -edge_gamma, -edge_delta, +edge_beta,
+};
+static const int edges_D[6] = {
+    -edge_zeta, +edge_gamma, +edge_beta,
+    -edge_epsilon, +edge_alpha, -edge_gamma,
+};
+static const int edges_J[6] = {
+    -edge_beta, +edge_gamma, +edge_beta,
+    +edge_theta, +edge_beta, edge_eta,
+};
+static const int edges_L[6] = {
+    -edge_beta, +edge_gamma, +edge_beta,
+    -edge_epsilon, +edge_alpha, -edge_theta,
+};
+static const int edges_X[6] = {
+    -edge_beta, -edge_alpha, +edge_epsilon,
+    +edge_theta, +edge_beta, edge_eta,
+};
+static const int edges_P[6] = {
+    -edge_beta, -edge_alpha, +edge_epsilon,
+    -edge_epsilon, +edge_alpha, -edge_theta,
+};
+static const int edges_S[6] = {
+    +edge_delta, +edge_zeta, +edge_beta,
+    -edge_epsilon, +edge_alpha, -edge_gamma,
+};
+static const int edges_F[6] = {
+    -edge_beta, +edge_gamma, +edge_beta,
+    -edge_epsilon, +edge_epsilon, edge_eta,
+};
+static const int edges_Y[6] = {
+    -edge_beta, -edge_alpha, +edge_epsilon,
+    -edge_epsilon, +edge_epsilon, edge_eta,
+};
+
+/*
+ * Now specify the actual shape of each edge type, in terms of the
+ * angles of turns as you traverse the edge.
+ *
+ * Edges around the outline of a hex expansion are traversed
+ * _clockwise_, because each expansion step flips the handedness of
+ * the whole system.
+ *
+ * Each array has one fewer element than the number of sub-edges in
+ * the edge shape (for the usual reason - n edges in a path have only
+ * n-1 vertices separating them).
+ *
+ * These arrays show the positive version of each edge type. The
+ * negative version is obtained by reversing the order of the turns
+ * and also the sign of each turn.
+ */
+static const int hex_edge_shape_eta[] = { +2, +2, -2, -2 };
+static const int hex_edge_shape_alpha[] = { +2, -2 };
+static const int hex_edge_shape_beta[] = { -2 };
+static const int hex_edge_shape_gamma[] = { +2, -2, -2, +2 };
+static const int hex_edge_shape_delta[] = { -2, +2, -2, +2 };
+static const int hex_edge_shape_epsilon[] = { +2, -2, -2 };
+static const int hex_edge_shape_zeta[] = { -2, +2 };
+static const int hex_edge_shape_theta[] = { +2, +2, -2, -2, +2 };
+
+static const int *const hex_edge_shapes[] = {
+    hex_edge_shape_eta,
+    hex_edge_shape_alpha,
+    hex_edge_shape_beta,
+    hex_edge_shape_gamma,
+    hex_edge_shape_delta,
+    hex_edge_shape_epsilon,
+    hex_edge_shape_zeta,
+    hex_edge_shape_theta,
+};
+static const size_t hex_edge_lengths[] = {
+    lenof(hex_edge_shape_eta) + 1,
+    lenof(hex_edge_shape_alpha) + 1,
+    lenof(hex_edge_shape_beta) + 1,
+    lenof(hex_edge_shape_gamma) + 1,
+    lenof(hex_edge_shape_delta) + 1,
+    lenof(hex_edge_shape_epsilon) + 1,
+    lenof(hex_edge_shape_zeta) + 1,
+    lenof(hex_edge_shape_theta) + 1,
+};
+
+static const int spec_edge_shape_eta[] = { 0 };
+static const int spec_edge_shape_alpha[] = { -2, +3 };
+static const int spec_edge_shape_beta[] = { +3, -2 };
+static const int spec_edge_shape_gamma[] = { +2 };
+static const int spec_edge_shape_delta[] = { +2, +3, +2, -3, +2 };
+static const int spec_edge_shape_epsilon[] = { +3 };
+static const int spec_edge_shape_zeta[] = { -2 };
+/* In expansion to Spectres, a theta edge corresponds to just one
+ * Spectre edge, so its turns array would be completely empty! */
+
+static const int *const spec_edge_shapes[] = {
+    spec_edge_shape_eta,
+    spec_edge_shape_alpha,
+    spec_edge_shape_beta,
+    spec_edge_shape_gamma,
+    spec_edge_shape_delta,
+    spec_edge_shape_epsilon,
+    spec_edge_shape_zeta,
+    NULL, /* theta has no turns */
+};
+static const size_t spec_edge_lengths[] = {
+    lenof(spec_edge_shape_eta) + 1,
+    lenof(spec_edge_shape_alpha) + 1,
+    lenof(spec_edge_shape_beta) + 1,
+    lenof(spec_edge_shape_gamma) + 1,
+    lenof(spec_edge_shape_delta) + 1,
+    lenof(spec_edge_shape_epsilon) + 1,
+    lenof(spec_edge_shape_zeta) + 1,
+    1, /* theta is only one edge long */
+};
+
+/*
+ * Each edge type corresponds to a fixed number of edges of the
+ * hexagon layout in the expansion of each hex, and also to a fixed
+ * number of edges of the Spectre(s) that each hex expands to in the
+ * final step.
+ */
+static const int edgelen_hex[] = {
+    5, /* edge_eta */
+    3, /* edge_alpha */
+    2, /* edge_beta */
+    5, /* edge_gamma */
+    5, /* edge_delta */
+    4, /* edge_epsilon */
+    3, /* edge_zeta */
+    6, /* edge_theta */
+};
+
+static const int edgelen_spectre[] = {
+    2, /* edge_eta */
+    3, /* edge_alpha */
+    3, /* edge_beta */
+    2, /* edge_gamma */
+    6, /* edge_delta */
+    2, /* edge_epsilon */
+    2, /* edge_zeta */
+    1, /* edge_theta */
+};

auxiliary/spectre-test.c

@@ -0,0 +1,534 @@
+/*
+ * Standalone test program for spectre.c.
+ */
+
+#include <assert.h>
+#ifdef NO_TGMATH_H
+#  include <math.h>
+#else
+#  include <tgmath.h>
+#endif
+#include <stdarg.h>
+#include <stdio.h>
+#include <string.h>
+
+#include "puzzles.h"
+#include "spectre-internal.h"
+#include "spectre-tables-manual.h"
+#include "spectre-tables-auto.h"
+#include "spectre-help.h"
+
+static void step_tests(void)
+{
+    SpectreContext ctx[1];
+    random_state *rs;
+    SpectreCoords *sc;
+    unsigned outedge;
+
+    rs = random_new("12345", 5);
+    spectrectx_init_random(ctx, rs);
+
+    /* Simplest possible transition: between the two Spectres making
+     * up a G hex. */
+    sc = spectre_coords_new();
+    spectre_coords_make_space(sc, 1);
+    sc->index = 0;
+    sc->nc = 1;
+    sc->c[0].type = HEX_G;
+    sc->c[0].index = -1;
+    spectrectx_step(ctx, sc, 12, &outedge);
+    assert(outedge == 5);
+    assert(sc->index == 1);
+    assert(sc->nc == 1);
+    assert(sc->c[0].type == HEX_G);
+    assert(sc->c[0].index == -1);
+    spectre_coords_free(sc);
+
+    /* Test the double Spectre transition. Here, within a F superhex,
+     * we attempt to step from the G subhex to the S one, in such a
+     * way that the place where we enter the Spectre corresponding to
+     * the S hex is on its spur of detached edge, causing us to
+     * immediately transition back out of the other side of that spur
+     * and end up in the D subhex instead. */
+    sc = spectre_coords_new();
+    spectre_coords_make_space(sc, 2);
+    sc->index = 1;
+    sc->nc = 2;
+    sc->c[0].type = HEX_G;
+    sc->c[0].index = 2;
+    sc->c[1].type = HEX_F;
+    sc->c[1].index = -1;
+    spectrectx_step(ctx, sc, 1, &outedge);
+    assert(outedge == 6);
+    assert(sc->index == 0);
+    assert(sc->nc == 2);
+    assert(sc->c[0].type == HEX_D);
+    assert(sc->c[0].index == 5);
+    assert(sc->c[1].type == HEX_F);
+    assert(sc->c[1].index == -1);
+    spectre_coords_free(sc);
+
+    /* However, _this_ transition leaves the same G subhex by the same
+     * edge of the hexagon, but further along it, so that we land in
+     * the S Spectre and stay there, without needing a double
+     * transition. */
+    sc = spectre_coords_new();
+    spectre_coords_make_space(sc, 2);
+    sc->index = 1;
+    sc->nc = 2;
+    sc->c[0].type = HEX_G;
+    sc->c[0].index = 2;
+    sc->c[1].type = HEX_F;
+    sc->c[1].index = -1;
+    spectrectx_step(ctx, sc, 13, &outedge);
+    assert(outedge == 4);
+    assert(sc->index == 0);
+    assert(sc->nc == 2);
+    assert(sc->c[0].type == HEX_S);
+    assert(sc->c[0].index == 3);
+    assert(sc->c[1].type == HEX_F);
+    assert(sc->c[1].index == -1);
+    spectre_coords_free(sc);
+
+    /* A couple of randomly generated transition tests that go a long
+     * way up the stack. */
+    sc = spectre_coords_new();
+    spectre_coords_make_space(sc, 7);
+    sc->index = 0;
+    sc->nc = 7;
+    sc->c[0].type = HEX_S;
+    sc->c[0].index = 3;
+    sc->c[1].type = HEX_Y;
+    sc->c[1].index = 7;
+    sc->c[2].type = HEX_Y;
+    sc->c[2].index = 4;
+    sc->c[3].type = HEX_Y;
+    sc->c[3].index = 4;
+    sc->c[4].type = HEX_F;
+    sc->c[4].index = 0;
+    sc->c[5].type = HEX_X;
+    sc->c[5].index = 1;
+    sc->c[6].type = HEX_G;
+    sc->c[6].index = -1;
+    spectrectx_step(ctx, sc, 13, &outedge);
+    assert(outedge == 12);
+    assert(sc->index == 0);
+    assert(sc->nc == 7);
+    assert(sc->c[0].type == HEX_Y);
+    assert(sc->c[0].index == 1);
+    assert(sc->c[1].type == HEX_P);
+    assert(sc->c[1].index == 1);
+    assert(sc->c[2].type == HEX_D);
+    assert(sc->c[2].index == 5);
+    assert(sc->c[3].type == HEX_Y);
+    assert(sc->c[3].index == 4);
+    assert(sc->c[4].type == HEX_X);
+    assert(sc->c[4].index == 7);
+    assert(sc->c[5].type == HEX_S);
+    assert(sc->c[5].index == 3);
+    assert(sc->c[6].type == HEX_G);
+    assert(sc->c[6].index == -1);
+    spectre_coords_free(sc);
+
+    sc = spectre_coords_new();
+    spectre_coords_make_space(sc, 7);
+    sc->index = 0;
+    sc->nc = 7;
+    sc->c[0].type = HEX_Y;
+    sc->c[0].index = 7;
+    sc->c[1].type = HEX_F;
+    sc->c[1].index = 6;
+    sc->c[2].type = HEX_Y;
+    sc->c[2].index = 4;
+    sc->c[3].type = HEX_X;
+    sc->c[3].index = 7;
+    sc->c[4].type = HEX_L;
+    sc->c[4].index = 0;
+    sc->c[5].type = HEX_S;
+    sc->c[5].index = 3;
+    sc->c[6].type = HEX_F;
+    sc->c[6].index = -1;
+    spectrectx_step(ctx, sc, 0, &outedge);
+    assert(outedge == 1);
+    assert(sc->index == 0);
+    assert(sc->nc == 7);
+    assert(sc->c[0].type == HEX_P);
+    assert(sc->c[0].index == 1);
+    assert(sc->c[1].type == HEX_F);
+    assert(sc->c[1].index == 0);
+    assert(sc->c[2].type == HEX_Y);
+    assert(sc->c[2].index == 7);
+    assert(sc->c[3].type == HEX_F);
+    assert(sc->c[3].index == 0);
+    assert(sc->c[4].type == HEX_G);
+    assert(sc->c[4].index == 2);
+    assert(sc->c[5].type == HEX_D);
+    assert(sc->c[5].index == 5);
+    assert(sc->c[6].type == HEX_F);
+    assert(sc->c[6].index == -1);
+    spectre_coords_free(sc);
+
+    spectrectx_cleanup(ctx);
+    random_free(rs);
+}
+
+struct genctx {
+    Graphics *gr;
+    FILE *fp; /* for non-graphical output modes */
+    random_state *rs;
+    Coord xmin, xmax, ymin, ymax;
+};
+
+static void gctx_set_size(
+    struct genctx *gctx, int width, int height, double scale,
+    int *xmin, int *xmax, int *ymin, int *ymax)
+{
+    *xmax = ceil(width/(2*scale));
+    *xmin = -*xmax;
+    *ymax = ceil(height/(2*scale));
+    *ymin = -*ymax;
+
+    /* point_x() and point_y() double their output to avoid having
+     * to use fractions, so double the bounds we'll compare their
+     * results against */
+    gctx->xmin.c1 = *xmin * 2; gctx->xmin.cr3 = 0;
+    gctx->xmax.c1 = *xmax * 2; gctx->xmax.cr3 = 0;
+    gctx->ymin.c1 = *ymin * 2; gctx->ymin.cr3 = 0;
+    gctx->ymax.c1 = *ymax * 2; gctx->ymax.cr3 = 0;
+}
+
+static bool callback(void *vctx, const Spectre *spec)
+{
+    struct genctx *gctx = (struct genctx *)vctx;
+    size_t i;
+
+    for (i = 0; i < 14; i++) {
+        Point p = spec->vertices[i];
+        Coord x = point_x(p), y = point_y(p);
+        if (coord_cmp(x, gctx->xmin) >= 0 && coord_cmp(x, gctx->xmax) <= 0 &&
+            coord_cmp(y, gctx->ymin) >= 0 && coord_cmp(y, gctx->ymax) <= 0)
+            goto ok;
+    }
+    return false;
+
+  ok:
+    gr_draw_spectre_from_coords(gctx->gr, spec->sc, spec->vertices);
+    if (gctx->fp) {
+        /*
+         * Emit calls to a made-up Python 'spectre()' function which
+         * takes the following parameters:
+         *
+         *  - lowest-level hexagon type (one-character string)
+         *  - index of Spectre within hexagon (0 or rarely 1)
+         *  - array of 14 point coordinates. Each is a 2-tuple
+         *    containing x and y. Each of those in turn is a 2-tuple
+         *    containing coordinates of 1 and sqrt(3).
+         */
+        fprintf(gctx->fp, "spectre('%s', %d, [",
+                hex_names[spec->sc->c[0].type], spec->sc->index);
+        for (i = 0; i < 14; i++) {
+            Point p = spec->vertices[i];
+            Coord x = point_x(p), y = point_y(p);
+            fprintf(gctx->fp, "%s((%d,%d),(%d,%d))", i ? ", " : "",
+                    x.c1, x.cr3, y.c1, y.cr3);
+        }
+        fprintf(gctx->fp, "])\n");
+    }
+    return true;
+}
+
+static void generate(struct genctx *gctx)
+{
+    SpectreContext ctx[1];
+    
+    spectrectx_init_random(ctx, gctx->rs);
+    ctx->prototype->hex_colour = random_upto(gctx->rs, 3);
+    ctx->prototype->prev_hex_colour = (ctx->prototype->hex_colour + 1 +
+                                       random_upto(gctx->rs, 2)) % 3;
+    ctx->prototype->incoming_hex_edge = random_upto(gctx->rs, 2);
+
+    spectrectx_generate(ctx, callback, gctx);
+
+    spectrectx_cleanup(ctx);
+}
+
+static inline Point reflected(Point p)
+{
+    /*
+     * This reflection operation is used as a conjugation, so it
+     * doesn't matter _what_ reflection it is, only that it reverses
+     * sense.
+     */
+    Point r;
+    size_t i;
+    for (i = 0; i < 4; i++)
+        r.coeffs[i] = p.coeffs[3-i];
+    return r;
+}
+static void reflect_spectre(Spectre *spec)
+{
+    size_t i;
+    for (i = 0; i < 14; i++)
+        spec->vertices[i] = reflected(spec->vertices[i]);
+}
+
+static void periodic_cheat(struct genctx *gctx)
+{
+    Spectre start, sh, sv;
+    size_t i;
+
+    start.sc = NULL;
+    {
+        Point u = {{ 0, 0, 0, 0 }};
+        Point v = {{ 1, 0, 0, 1 }};
+        v = point_mul(v, point_rot(1));
+        spectre_place(&start, u, v, 0);
+    }
+
+    sh = start;
+    while (callback(gctx, &sh)) {
+        sv = sh;
+        i = 0;
+        do {
+            if (i) {
+                spectre_place(&sv, sv.vertices[6], sv.vertices[7], 0);
+            } else {
+                spectre_place(&sv, reflected(sv.vertices[6]),
+                              reflected(sv.vertices[7]), 0);
+                reflect_spectre(&sv);
+            }
+            i ^= 1;
+        } while (callback(gctx, &sv));
+
+        sv = sh;
+        i = 0;
+        do {
+            if (i) {
+                spectre_place(&sv, sv.vertices[0], sv.vertices[1], 6);
+            } else {
+                spectre_place(&sv, reflected(sv.vertices[0]),
+                              reflected(sv.vertices[1]), 6);
+                reflect_spectre(&sv);
+            }
+            i ^= 1;
+        } while (callback(gctx, &sv));
+
+        spectre_place(&sh, sh.vertices[12], sh.vertices[11], 4);
+    }
+
+    sh = start;
+    do {
+        spectre_place(&sh, sh.vertices[5], sh.vertices[4], 11);
+
+        sv = sh;
+        i = 0;
+        do {
+            if (i) {
+                spectre_place(&sv, sv.vertices[6], sv.vertices[7], 0);
+            } else {
+                spectre_place(&sv, reflected(sv.vertices[6]),
+                              reflected(sv.vertices[7]), 0);
+                reflect_spectre(&sv);
+            }
+            i ^= 1;
+        } while (callback(gctx, &sv));
+
+        sv = sh;
+        i = 0;
+        do {
+            if (i) {
+                spectre_place(&sv, sv.vertices[0], sv.vertices[1], 6);
+            } else {
+                spectre_place(&sv, reflected(sv.vertices[0]),
+                              reflected(sv.vertices[1]), 6);
+                reflect_spectre(&sv);
+            }
+            i ^= 1;
+        } while (callback(gctx, &sv));
+    } while (callback(gctx, &sh));
+}
+
+static void generate_hexes(struct genctx *gctx)
+{
+    SpectreContext ctx[1];
+    spectrectx_init_random(ctx, gctx->rs);
+    SpectreCoords *sc;
+    unsigned orient, outedge, inedge;
+    bool printed_any = false;
+    size_t r = 1, ri = 0, rj = 0;
+
+    Point centre = {{ 0, 0, 0, 0 }};
+    const Point six = {{ 6, 0, 0, 0 }};
+
+    sc = spectre_coords_copy(ctx->prototype);
+    orient = random_upto(gctx->rs, 6);
+
+    while (true) {
+        Point top = {{ -2, 0, 4, 0 }};
+        Point vertices[6];
+        bool print_this = false;
+        size_t i;
+
+        for (i = 0; i < 6; i++) {
+            vertices[i] = point_add(centre, point_mul(
+                                        top, point_rot(2 * (orient + i))));
+            Coord x = point_x(vertices[i]), y = point_y(vertices[i]);
+            if (coord_cmp(x, gctx->xmin) >= 0 &&
+                coord_cmp(x, gctx->xmax) <= 0 &&
+                coord_cmp(y, gctx->ymin) >= 0 &&
+                coord_cmp(y, gctx->ymax) <= 0)
+                print_this = true;
+        }
+
+        if (print_this) {
+            printed_any = true;
+            gr_draw_hex(gctx->gr, -1, sc->c[0].type, vertices);
+        }
+
+        /*
+         * Decide which way to step next. We spiral outwards from a
+         * central hexagon.
+         */
+        outedge = (ri == 0 && rj == 0) ? 5 : ri;
+        if (++rj >= r) {
+            rj = 0;
+            if (++ri >= 6) {
+                ri = 0;
+                if (!printed_any)
+                    break;
+                printed_any = false;
+                ++r;
+            }
+        }
+
+        spectrectx_step_hex(ctx, sc, 0, (outedge + 6 - orient) % 6, &inedge);
+        orient = (outedge + 9 - inedge) % 6;
+
+        centre = point_add(centre, point_mul(six, point_rot(4 + 2 * outedge)));
+    }
+
+    spectre_coords_free(sc);
+    spectrectx_cleanup(ctx);
+}
+
+int main(int argc, char **argv)
+{
+    const char *random_seed = "12345";
+    const char *outfile = "-";
+    bool four_colour = false;
+    enum { TESTS, TILING, CHEAT, HEXES } mode = TILING;
+    enum { SVG, PYTHON } outmode = SVG;
+    double scale = 10, linewidth = 1.5;
+    int width = 1024, height = 768;
+    bool arcs = false;
+
+    while (--argc > 0) {
+        const char *arg = *++argv;
+        if (!strcmp(arg, "--help")) {
+            printf("  usage: spectre-test [FIXME]\n"
+                   "   also: spectre-test --test\n");
+            return 0;
+        } else if (!strcmp(arg, "--test")) {
+            mode = TESTS;
+        } else if (!strcmp(arg, "--hex")) {
+            mode = HEXES;
+        } else if (!strcmp(arg, "--cheat")) {
+            mode = CHEAT;
+        } else if (!strcmp(arg, "--python")) {
+            outmode = PYTHON;
+        } else if (!strcmp(arg, "--arcs")) {
+            arcs = true;
+        } else if (!strncmp(arg, "--seed=", 7)) {
+            random_seed = arg+7;
+        } else if (!strcmp(arg, "--fourcolour")) {
+            four_colour = true;
+        } else if (!strncmp(arg, "--scale=", 8)) {
+            scale = atof(arg+8);
+        } else if (!strncmp(arg, "--width=", 8)) {
+            width = atof(arg+8);
+        } else if (!strncmp(arg, "--height=", 9)) {
+            height = atof(arg+9);
+        } else if (!strncmp(arg, "--linewidth=", 12)) {
+            linewidth = atof(arg+12);
+        } else if (!strcmp(arg, "-o")) {
+            if (--argc <= 0) {
+                fprintf(stderr, "expected argument to '%s'\n", arg);
+                return 1;
+            }
+            outfile = *++argv;
+        } else {
+            fprintf(stderr, "unexpected extra argument '%s'\n", arg);
+            return 1;
+        }
+    }
+
+    switch (mode) {
+      case TESTS: {
+        step_tests();
+        break;
+      }
+
+      case TILING:
+      case CHEAT: {
+        struct genctx gctx[1];
+        bool close_output = false;
+        int xmin, xmax, ymin, ymax;
+
+        gctx_set_size(gctx, width, height, scale, &xmin, &xmax, &ymin, &ymax);
+
+        switch (outmode) {
+          case SVG:
+            gctx->gr = gr_new(outfile, xmin, xmax, ymin, ymax, scale);
+            gctx->gr->number_cells = false;
+            gctx->gr->four_colour = four_colour;
+            gctx->gr->linewidth = linewidth;
+            gctx->gr->arcs = arcs;
+            gctx->fp = NULL;
+            break;
+          case PYTHON:
+            gctx->gr = NULL;
+            if (!strcmp(outfile, "-")) {
+                gctx->fp = stdout;
+            } else {
+                gctx->fp = fopen(outfile, "w");
+                close_output = true;
+            }
+            break;
+        }
+
+        gctx->rs = random_new(random_seed, strlen(random_seed));
+        switch (mode) {
+          case TILING:
+            generate(gctx);
+            break;
+          case CHEAT:
+            periodic_cheat(gctx);
+            break;
+          default: /* shouldn't happen */
+            break;
+        }
+        random_free(gctx->rs);
+        gr_free(gctx->gr);
+        if (close_output)
+            fclose(gctx->fp);
+        break;
+      }
+
+      case HEXES: {
+        struct genctx gctx[1];
+        int xmin, xmax, ymin, ymax;
+
+        gctx_set_size(gctx, width, height, scale, &xmin, &xmax, &ymin, &ymax);
+
+        gctx->gr = gr_new(outfile, xmin, xmax, ymin, ymax, scale);
+        gctx->gr->jigsaw_mode = true;
+        gctx->gr->number_edges = false;
+        gctx->gr->linewidth = linewidth;
+        gctx->rs = random_new(random_seed, strlen(random_seed));
+        generate_hexes(gctx);          /* FIXME: bounds */
+        random_free(gctx->rs);
+        gr_free(gctx->gr);
+        break;
+      }
+    }
+}

grid.c

@@ -24,6 +24,7 @@
 #include "grid.h"
 #include "penrose.h"
 #include "hat.h"
+#include "spectre.h"
 
 /* Debugging options */
 
@@ -3562,6 +3563,316 @@ static grid *grid_new_hats(int width, int height, const char *desc)
     return ctx->g;
 }
 
+#define SPECTRE_TILESIZE 32
+#define SPECTRE_SQUARELEN 7
+#define SPECTRE_UNIT 8
+
+static const char *grid_validate_params_spectres(
+    int width, int height)
+{
+    int l = SPECTRE_UNIT * SPECTRE_SQUARELEN;
+
+    if (width > INT_MAX / l ||                  /* xextent */
+        height > INT_MAX / l ||                 /* yextent */
+        width > (INT_MAX / SPECTRE_SQUARELEN /
+                 SPECTRE_SQUARELEN / height))   /* max_faces */
+        return "Grid must not be unreasonably large";
+    return NULL;
+}
+
+static void grid_size_spectres(int width, int height,
+                               int *tilesize, int *xextent, int *yextent)
+{
+    *tilesize = SPECTRE_TILESIZE;
+    *xextent = width * SPECTRE_UNIT * SPECTRE_SQUARELEN;
+    *yextent = height * SPECTRE_UNIT * SPECTRE_SQUARELEN;
+}
+
+static char *grid_new_desc_spectres(
+    grid_type type, int width, int height, random_state *rs)
+{
+    char *buf;
+    size_t i;
+    struct SpectrePatchParams sp;
+
+    spectre_tiling_randomise(&sp, width * SPECTRE_SQUARELEN,
+                             height * SPECTRE_SQUARELEN, rs);
+
+    buf = snewn(sp.ncoords + 3, char);
+    buf[0] = (sp.orientation < 10 ? '0' + sp.orientation :
+              'A' + sp.orientation - 10);
+    for (i = 0; i < sp.ncoords; i++) {
+        assert(sp.coords[i] < 10);    /* all indices are 1 digit */
+        buf[i+1] = '0' + sp.coords[i];
+    }
+    buf[sp.ncoords+1] = sp.final_hex;
+    buf[sp.ncoords+2] = '\0';
+
+    sfree(sp.coords);
+    return buf;
+}
+
+/* Shared code between validating and reading grid descs.
+ * Always allocates sp->coords, whether or not it returns an error. */
+static const char *grid_desc_to_spectre_params(
+    const char *desc, struct SpectrePatchParams *sp)
+{
+    size_t i;
+
+    if (!*desc)
+        return "empty grid description";
+
+    sp->ncoords = strlen(desc) - 2;
+    sp->coords = snewn(sp->ncoords, unsigned char);
+
+    {
+        char c = desc[0];
+        if (isdigit((unsigned char)c))
+            sp->orientation = c - '0';
+        else if (c == 'A' || c == 'B')
+            sp->orientation = 10 + c - 'A';
+        else
+            return "expected digit or A,B at start of grid description";
+    }
+
+    for (i = 0; i < sp->ncoords; i++) {
+        char c = desc[i+1];
+        if (!isdigit((unsigned char)c))
+            return "expected digit in grid description";
+        sp->coords[i] = c - '0';
+    }
+
+    sp->final_hex = desc[sp->ncoords+1];
+
+    return NULL;
+}
+
+static const char *grid_validate_desc_spectres(
+    grid_type type, int width, int height, const char *desc)
+{
+    struct SpectrePatchParams sp;
+    const char *error = NULL;
+
+    if (!desc)
+        return "Missing grid description string.";
+
+    error = grid_desc_to_spectre_params(desc, &sp);
+    if (!error)
+        error = spectre_tiling_params_invalid(&sp);
+
+    sfree(sp.coords);
+    return error;
+}
+
+struct spectrecontext {
+    grid *g;
+    tree234 *points;
+};
+
+/*
+ * Calculate the nearest integer to n*sqrt(3), via a bitwise algorithm
+ * that avoids floating point.
+ *
+ * (It would probably be OK in practice to use floating point, but I
+ * felt like overengineering it for fun. With FP, there's at least a
+ * theoretical risk of rounding the wrong way, due to the three
+ * successive roundings involved - rounding sqrt(3), rounding its
+ * product with n, and then rounding to the nearest integer. This
+ * approach avoids that: it's exact.)
+ */
+static int mul_root3(int n_signed)
+{
+    unsigned x, r, m;
+    int sign = n_signed < 0 ? -1 : +1;
+    unsigned n = n_signed * sign;
+    unsigned bitpos;
+
+    /*
+     * Method:
+     *
+     * We transform m gradually from zero into n, by multiplying it by
+     * 2 in each step and optionally adding 1, so that it's always
+     * floor(n/2^something).
+     *
+     * At the start of each step, x is the largest integer less than
+     * or equal to m*sqrt(3). We transform m to 2m+bit, and therefore
+     * we must transform x to 2x+something to match. The 'something'
+     * we add to 2x is at most 3. (Worst case is if m sqrt(3) was
+     * equal to x + 1-eps for some tiny eps, and then the incoming bit
+     * of m is 1, so that (2m+1)sqrt(3) = 2x+2+2eps+sqrt(3), i.e.
+     * about 2x + 3.732...)
+     *
+     * To compute this, we also track the residual value r such that
+     * x^2+r = 3m^2.
+     *
+     * The algorithm below is very similar to the usual approach for
+     * taking the square root of an integer in binary. The wrinkle is
+     * that we have an integer multiplier, i.e. we're computing
+     * P*sqrt(Q) (with P=n and Q=3 in this case) rather than just
+     * sqrt(Q). Of course in principle we could just take sqrt(P^2Q),
+     * but we'd need an integer twice the width to hold P^2. Pulling
+     * out P and treating it specially makes overflow less likely.
+     */
+
+    x = r = m = 0;
+
+    for (bitpos = UINT_MAX & ~(UINT_MAX >> 1); bitpos; bitpos >>= 1) {
+        unsigned a, b = (n & bitpos) ? 1 : 0;
+
+        /*
+         * Check invariants. We expect that x^2 + r = 3m^2 (i.e. our
+         * residual term is correct), and also that r < 2x+1 (because
+         * if not, then we could replace x with x+1 and still get a
+         * value that made r non-negative, i.e. x would not be the
+         * _largest_ integer less than m sqrt(3)).
+         */
+        assert(x*x + r == 3*m*m);
+        assert(r < 2*x+1);
+
+        /*
+         * We're going to replace m with 2m+b, and x with 2x+a for
+         * some a we haven't decided on yet.
+         *
+         * The new value of the residual will therefore be
+         *
+         *   3 (2m+b)^2 - (2x+a)^2
+         * = (12m^2 + 12mb + 3b^2) - (4x^2 + 4xa + a^2)
+         * = 4 (3m^2 - x^2) + 12mb + 3b^2 - 4xa - a^2
+         * = 4r + 12mb + 3b^2 - 4xa - a^2          (because r = 3m^2 - x^2)
+         * = 4r + (12m + 3)b - 4xa - a^2           (b is 0 or 1, so b = b^2)
+         */
+        for (a = 0; a < 4; a++) {
+            /* If we made this routine handle square roots of numbers
+             * other than 3 then it would be sensible to make this a
+             * binary search. Here, it hardly seems important. */
+            unsigned pos = 4*r + b*(12*m + 3);
+            unsigned neg = 4*a*x + a*a;
+            if (pos < neg)
+                break;                 /* this value of a is too big */
+        }
+
+        /* The above loop will have terminated with a one too big,
+         * whether that's because we hit the break statement or fell
+         * off the end with a=4. So now decrementing a will give us
+         * the right value to add. */
+        a--;
+
+        r = 4*r + b*(12*m + 3) - (4*a*x + a*a);
+        m = 2*m+b;
+        x = 2*x+a;
+    }
+
+    /*
+     * Finally, round to the nearest integer. At present, x is the
+     * largest integer that is _at most_ m sqrt(3). But we want the
+     * _nearest_ integer, whether that's rounded up or down. So check
+     * whether (x + 1/2) is still less than m sqrt(3), i.e. whether
+     * (x + 1/2)^2 < 3m^2; if it is, then we increment x.
+     *
+     * We have 3m^2 - (x + 1/2)^2 = 3m^2 - x^2 - x - 1/4
+     *                            = r - x - 1/4
+     *
+     * and since r and x are integers, this is greater than 0 if and
+     * only if r > x.
+     *
+     * (There's no need to worry about tie-breaking exact halfway
+     * rounding cases. sqrt(3) is irrational, so none such exist.)
+     */
+    if (r > x)
+        x++;
+
+    /*
+     * Put the sign back on, and convert back from unsigned to int.
+     */
+    if (sign == +1) {
+        return x;
+    } else {
+        /* Be a little careful to avoid compilers deciding I've just
+         * perpetrated signed-integer overflow. This should optimise
+         * down to no actual code. */
+        return INT_MIN + (int)(-x - (unsigned)INT_MIN);
+    }
+}
+
+static void grid_spectres_callback(void *vctx, const int *coords)
+{
+    struct spectrecontext *ctx = (struct spectrecontext *)vctx;
+    size_t i;
+
+    grid_face_add_new(ctx->g, SPECTRE_NVERTICES);
+    for (i = 0; i < SPECTRE_NVERTICES; i++) {
+        grid_dot *d = grid_get_dot(
+            ctx->g, ctx->points,
+            (coords[4*i+0] * SPECTRE_UNIT +
+             mul_root3(coords[4*i+1] * SPECTRE_UNIT)),
+            (coords[4*i+2] * SPECTRE_UNIT +
+             mul_root3(coords[4*i+3] * SPECTRE_UNIT)));
+        grid_face_set_dot(ctx->g, d, i);
+    }
+}
+
+static grid *grid_new_spectres(int width, int height, const char *desc)
+{
+    struct SpectrePatchParams sp;
+    const char *error = NULL;
+    int width2 = width * SPECTRE_SQUARELEN;
+    int height2 = height * SPECTRE_SQUARELEN;
+
+    error = grid_desc_to_spectre_params(desc, &sp);
+    assert(error == NULL && "grid_validate_desc_spectres should have failed");
+
+    /*
+     * Bound on the number of faces: the area of a single face in the
+     * output coordinates is 24 + 24 rt3, which is between 65 and 66.
+     * Every face fits strictly inside the target rectangle, so the
+     * number of faces times a lower bound on their area can't exceed
+     * the area of the rectangle we give to spectre_tiling_generate.
+     */
+    int max_faces = width2 * height2 / 65;
+
+    /*
+     * Bound on number of dots: 14*faces is certainly enough, but
+     * quite wasteful given that _most_ dots are shared between at
+     * least two faces. But to get a better estimate we'd have to
+     * figure out a bound for the number of dots on the perimeter,
+     * which is the number by which the count exceeds 14*faces/2.
+     */
+    int max_dots = 14 * max_faces;
+
+    struct spectrecontext ctx[1];
+
+    ctx->g = grid_empty();
+    ctx->g->tilesize = SPECTRE_TILESIZE;
+    ctx->g->faces = snewn(max_faces, grid_face);
+    ctx->g->dots = snewn(max_dots, grid_dot);
+
+    ctx->points = newtree234(grid_point_cmp_fn);
+
+    spectre_tiling_generate(&sp, width2, height2, grid_spectres_callback, ctx);
+
+    freetree234(ctx->points);
+    sfree(sp.coords);
+
+    grid_trim_vigorously(ctx->g);
+    grid_make_consistent(ctx->g);
+
+    /*
+     * As with the Penrose tiling, we're likely to have different
+     * sized margins due to the lack of a neat grid that this tiling
+     * fits on. So now we know what tiles we're left with, recentre
+     * them.
+     */
+    {
+        int w = width2 * SPECTRE_UNIT, h = height2 * SPECTRE_UNIT;
+        ctx->g->lowest_x -= (w - (ctx->g->highest_x - ctx->g->lowest_x))/2;
+        ctx->g->lowest_y -= (h - (ctx->g->highest_y - ctx->g->lowest_y))/2;
+        ctx->g->highest_x = ctx->g->lowest_x + w;
+        ctx->g->highest_y = ctx->g->lowest_y + h;
+    }
+
+    return ctx->g;
+}
+
 /* ----------- End of grid generators ------------- */
 
 #define FNVAL(upper,lower) &grid_validate_params_ ## lower,
@@ -3588,6 +3899,8 @@ char *grid_new_desc(grid_type type, int width, int height, random_state *rs)
         return grid_new_desc_penrose(type, width, height, rs);
     } else if (type == GRID_HATS) {
         return grid_new_desc_hats(type, width, height, rs);
+    } else if (type == GRID_SPECTRES) {
+        return grid_new_desc_spectres(type, width, height, rs);
     } else if (type == GRID_TRIANGULAR) {
         return dupstr("0"); /* up-to-date version of triangular grid */
     } else {
@@ -3602,6 +3915,8 @@ const char *grid_validate_desc(grid_type type, int width, int height,
         return grid_validate_desc_penrose(type, width, height, desc);
     } else if (type == GRID_HATS) {
         return grid_validate_desc_hats(type, width, height, desc);
+    } else if (type == GRID_SPECTRES) {
+        return grid_validate_desc_spectres(type, width, height, desc);
     } else if (type == GRID_TRIANGULAR) {
         return grid_validate_desc_triangular(type, width, height, desc);
     } else {

grid.h

@@ -111,6 +111,7 @@ typedef struct grid {
   A(PENROSE_P2,penrose_p2_kite) \
   A(PENROSE_P3,penrose_p3_thick) \
   A(HATS,hats) \
+  A(SPECTRES,spectres) \
   /* end of list */
 
 #define ENUM(upper,lower) GRID_ ## upper,

loopy.c

@@ -285,6 +285,7 @@ static void check_caches(const solver_state* sstate);
     A("Kagome",KAGOME,3,3)                                      \
     A("Compass-Dodecagonal",COMPASSDODECAGONAL,2,2)             \
     A("Hats",HATS,6,6)                                          \
+    A("Spectres",SPECTRES,6,6)                                  \
     /* end of list */
 
 #define GRID_NAME(title,type,amin,omin) title,
@@ -557,6 +558,8 @@ static const game_params loopy_presets_more[] = {
     {  3,  3, DIFF_HARD,   LOOPY_GRID_GREATDODECAGONAL },
     {  3,  2, DIFF_HARD,   LOOPY_GRID_GREATGREATDODECAGONAL },
     {  3,  3, DIFF_HARD,   LOOPY_GRID_COMPASSDODECAGONAL },
+    {  6,  6, DIFF_HARD,   LOOPY_GRID_HATS },
+    {  6,  6, DIFF_HARD,   LOOPY_GRID_SPECTRES },
 #else
     { 10, 10, DIFF_HARD,   LOOPY_GRID_HONEYCOMB },
     {  5,  4, DIFF_HARD,   LOOPY_GRID_GREATHEXAGONAL },
@@ -568,6 +571,7 @@ static const game_params loopy_presets_more[] = {
     {  5,  3, DIFF_HARD,   LOOPY_GRID_GREATGREATDODECAGONAL },
     {  5,  4, DIFF_HARD,   LOOPY_GRID_COMPASSDODECAGONAL },
     { 10, 10, DIFF_HARD,   LOOPY_GRID_HATS },
+    { 10, 10, DIFF_HARD,   LOOPY_GRID_SPECTRES },
 #endif
 };
 

spectre-internal.h

@@ -0,0 +1,277 @@
+#include "spectre.h"
+
+/*
+ * List macro of the names for hexagon types, which will be reused all
+ * over the place.
+ *
+ * (I have to call the parameter to this list macro something other
+ * than X, because here, X is also one of the macro arguments!)
+ */
+#define HEX_LETTERS(Z) Z(G) Z(D) Z(J) Z(L) Z(X) Z(P) Z(S) Z(F) Z(Y)
+
+typedef enum Hex {
+    #define HEX_ENUM_DECL(x) HEX_##x,
+    HEX_LETTERS(HEX_ENUM_DECL)
+    #undef HEX_ENUM_DECL
+} Hex;
+
+static inline unsigned num_subhexes(Hex h)
+{
+    return h == HEX_G ? 7 : 8;
+}
+
+static inline unsigned num_spectres(Hex h)
+{
+    return h == HEX_G ? 2 : 1;
+}
+
+/*
+ * Data types used in the lookup tables.
+ */
+struct MapEntry {
+    bool internal;
+    unsigned char hi, lo;
+};
+struct MapEdge {
+    unsigned char startindex, len;
+};
+struct Possibility {
+    unsigned char hi, lo;
+    unsigned long prob;
+};
+
+/*
+ * Coordinate system for tracking Spectres and their hexagonal
+ * metatiles.
+ *
+ * SpectreCoords will store the index of a single Spectre within a
+ * smallest-size hexagon, plus an array of HexCoord each indexing a
+ * hexagon within the expansion of a larger hexagon.
+ *
+ * The last coordinate stored, sc->c[sc->nc-1], will have a hex type
+ * but no index (represented by index==-1). This means "we haven't
+ * decided yet what this level of metatile needs to be". If we need to
+ * refer to this level during the hatctx_step algorithm, we make it up
+ * at random, based on a table of what metatiles each type can
+ * possibly be part of, at what index.
+ */
+typedef struct HexCoord {
+    int index; /* index within that tile, or -1 if not yet known */
+    Hex type;  /* type of this hexagon */
+} HexCoord;
+
+typedef struct SpectreCoords {
+    int index;       /* index of Spectre within the order-0 hexagon */
+    HexCoord *c;
+    size_t nc, csize;
+
+    /* Used by spectre-test to four-colour output tilings, and
+     * maintained unconditionally because it's easier than making it
+     * conditional */
+    unsigned char hex_colour, prev_hex_colour, incoming_hex_edge;
+} SpectreCoords;
+
+SpectreCoords *spectre_coords_new(void);
+void spectre_coords_free(SpectreCoords *hc);
+void spectre_coords_make_space(SpectreCoords *hc, size_t size);
+SpectreCoords *spectre_coords_copy(SpectreCoords *hc_in);
+
+/*
+ * Coordinate system for locating Spectres in the plane.
+ *
+ * The 'Point' structure represents a single point by means of an
+ * integer linear combination of {1, d, d^2, d^3}, where d is the
+ * complex number exp(i pi/6) representing 1/12 of a turn about the
+ * origin.
+ *
+ * The 'Spectre' structure represents an entire Spectre in a tiling,
+ * giving both the locations of all of its vertices and its
+ * combinatorial coordinates. It also contains a linked-list pointer,
+ * used during breadth-first search to generate all the Spectres in an
+ * area.
+ */
+typedef struct Point {
+    int coeffs[4];
+} Point;
+typedef struct Spectre Spectre;
+struct Spectre {
+    Point vertices[14];
+    SpectreCoords *sc;
+    Spectre *next; /* used in breadth-first search */
+};
+
+/* Fill in all the coordinates of a Spectre starting from any single edge */
+void spectre_place(Spectre *spec, Point u, Point v, int index_of_u);
+
+/*
+ * A Point is really a complex number, so we can add, subtract and
+ * multiply them.
+ */
+static inline Point point_add(Point a, Point b)
+{
+    Point r;
+    size_t i;
+    for (i = 0; i < 4; i++)
+        r.coeffs[i] = a.coeffs[i] + b.coeffs[i];
+    return r;
+}
+static inline Point point_sub(Point a, Point b)
+{
+    Point r;
+    size_t i;
+    for (i = 0; i < 4; i++)
+        r.coeffs[i] = a.coeffs[i] - b.coeffs[i];
+    return r;
+}
+static inline Point point_mul_by_d(Point x)
+{
+    Point r;
+    /* Multiply by d by using the identity d^4 - d^2 + 1 = 0, so d^4 = d^2+1 */
+    r.coeffs[0] = -x.coeffs[3];
+    r.coeffs[1] = x.coeffs[0];
+    r.coeffs[2] = x.coeffs[1] + x.coeffs[3];
+    r.coeffs[3] = x.coeffs[2];
+    return r;
+}
+static inline Point point_mul(Point a, Point b)
+{
+    size_t i, j;
+    Point r;
+
+    /* Initialise r to be a, scaled by b's d^3 term */
+    for (j = 0; j < 4; j++)
+        r.coeffs[j] = a.coeffs[j] * b.coeffs[3];
+
+    /* Now iterate r = d*r + (next coefficient down), by Horner's rule */
+    for (i = 3; i-- > 0 ;) {
+        r = point_mul_by_d(r);
+        for (j = 0; j < 4; j++)
+            r.coeffs[j] += a.coeffs[j] * b.coeffs[i];
+    }
+
+    return r;
+}
+static inline bool point_equal(Point a, Point b)
+{
+    size_t i;
+    for (i = 0; i < 4; i++)
+        if (a.coeffs[i] != b.coeffs[i])
+            return false;
+    return true;
+}
+
+/*
+ * Return the Point corresponding to a rotation of s steps around the
+ * origin, i.e. a rotation by 30*s degrees or s*pi/6 radians.
+ */
+static inline Point point_rot(int s)
+{
+    Point r = {{ 1, 0, 0, 0 }};
+    Point dpower = {{ 0, 1, 0, 0 }};
+
+    /* Reduce to a sensible range */
+    s = s % 12;
+    if (s < 0)
+        s += 12;
+
+    while (true) {
+        if (s & 1)
+            r = point_mul(r, dpower);
+        s >>= 1;
+        if (!s)
+            break;
+        dpower = point_mul(dpower, dpower);
+    }
+
+    return r;
+}
+
+/*
+ * SpectreContext is the shared context of a whole run of the
+ * algorithm. Its 'prototype' SpectreCoords object represents the
+ * coordinates of the starting Spectre, and is extended as necessary;
+ * any other SpectreCoord that needs extending will copy the
+ * higher-order values from ctx->prototype as needed, so that once
+ * each choice has been made, it remains consistent.
+ *
+ * When we're inventing a random piece of tiling in the first place,
+ * we append to ctx->prototype by choosing a random (but legal)
+ * higher-level metatile for the current topmost one to turn out to be
+ * part of. When we're replaying a generation whose parameters are
+ * already stored, we don't have a random_state, and we make fixed
+ * decisions if not enough coordinates were provided, as in the
+ * corresponding hat.c system.
+ *
+ * For a normal (non-testing) caller, spectrectx_generate() is the
+ * main useful function. It breadth-first searches a whole area to
+ * generate all the Spectres in it, starting from a (typically
+ * central) one with the coordinates of ctx->prototype. The callback
+ * function processes each Spectre as it's generated, and returns true
+ * or false to indicate whether that Spectre is within the bounds of
+ * the target area (and therefore the search should continue exploring
+ * its neighbours).
+ */
+typedef struct SpectreContext {
+    random_state *rs;
+    bool must_free_rs;
+    Point start_vertices[2]; /* vertices 0,1 of the starting Spectre */
+    int orientation;         /* orientation to put in SpectrePatchParams */
+    SpectreCoords *prototype;
+} SpectreContext;
+
+void spectrectx_init_random(SpectreContext *ctx, random_state *rs);
+void spectrectx_init_from_params(
+    SpectreContext *ctx, const struct SpectrePatchParams *ps);
+void spectrectx_cleanup(SpectreContext *ctx);
+SpectreCoords *spectrectx_initial_coords(SpectreContext *ctx);
+void spectrectx_extend_coords(SpectreContext *ctx, SpectreCoords *hc,
+                              size_t n);
+void spectrectx_step(SpectreContext *ctx, SpectreCoords *sc,
+                     unsigned edge, unsigned *outedge);
+void spectrectx_generate(SpectreContext *ctx,
+                         bool (*callback)(void *cbctx, const Spectre *spec),
+                         void *cbctx);
+
+/* For spectre-test to directly generate a tiling of hexes */
+void spectrectx_step_hex(SpectreContext *ctx, SpectreCoords *sc,
+                         size_t depth, unsigned edge, unsigned *outedge);
+
+/* For extracting the point coordinates */
+typedef struct Coord {
+    int c1, cr3;      /* coefficients of 1 and sqrt(3) respectively */
+} Coord;
+
+static inline Coord point_x(Point p)
+{
+    Coord x = { 2 * p.coeffs[0] + p.coeffs[2], p.coeffs[1] };
+    return x;
+}
+
+static inline Coord point_y(Point p)
+{
+    Coord y = { 2 * p.coeffs[3] + p.coeffs[1], p.coeffs[2] };
+    return y;
+}
+
+static inline int coord_sign(Coord x)
+{
+    if (x.c1 == 0 && x.cr3 == 0)
+        return 0;
+    if (x.c1 >= 0 && x.cr3 >= 0)
+        return +1;
+    if (x.c1 <= 0 && x.cr3 <= 0)
+        return -1;
+
+    if (x.c1 * x.c1 > 3 * x.cr3 * x.cr3)
+        return x.c1 < 0 ? -1 : +1;
+    else
+        return x.cr3 < 0 ? -1 : +1;
+}
+
+static inline int coord_cmp(Coord a, Coord b)
+{
+    Coord diff;
+    diff.c1 = a.c1 - b.c1;
+    diff.cr3 = a.cr3 - b.cr3;
+    return coord_sign(diff);
+}

spectre-tables-manual.h

@@ -0,0 +1,160 @@
+/*
+ * Handwritten data tables for the Spectre tiling.
+ *
+ * This file is used by both the final tiling generator in spectre.c,
+ * and by spectre-gen.c which auto-generates further tables to go with
+ * these.
+ */
+
+/*
+ * We generate the Spectre tiling based on the substitution system of
+ * 9 types of marked hexagon shown in the paper.
+ *
+ * The substitution expands each hexagon into 8 others, except for the
+ * G hex which expands to only seven. The layout, numbered with the
+ * indices we use in the arrays here, is as follows:
+ *
+ *     0 1
+ *    2 3
+ *   4 5 6
+ *      7
+ *
+ * That is: the hexes are oriented with a pair of vertical edges.
+ * Hexes 0 and 1 are horizontally adjacent; 2 and 3 are adjacent on
+ * the next row, with 3 nestling between 0 and 1; 4,5,6 are on the
+ * third row with 5 between 2 and 3; and 7 is by itself on a fourth
+ * row, between 5 and 6. In the expansion of the G hex, #7 is the
+ * missing one, so its indices are still consecutive from 0.
+ *
+ * These arrays list the type of each child hex. The hexes also have
+ * orientations, but those aren't listed here, because only
+ * spectre-gen needs to know them - by the time it's finished
+ * autogenerating transition tables, the orientations are baked into
+ * those and don't need to be consulted separately.
+ */
+
+static const Hex subhexes_G[] = {
+    HEX_F,
+    HEX_X,
+    HEX_G,
+    HEX_S,
+    HEX_P,
+    HEX_D,
+    HEX_J,
+    /* hex #7 is not present for this tile */
+};
+static const Hex subhexes_D[] = {
+    HEX_F,
+    HEX_P,
+    HEX_G,
+    HEX_S,
+    HEX_X,
+    HEX_D,
+    HEX_F,
+    HEX_X,
+};
+static const Hex subhexes_J[] = {
+    HEX_F,
+    HEX_P,
+    HEX_G,
+    HEX_S,
+    HEX_Y,
+    HEX_D,
+    HEX_F,
+    HEX_P,
+};
+static const Hex subhexes_L[] = {
+    HEX_F,
+    HEX_P,
+    HEX_G,
+    HEX_S,
+    HEX_Y,
+    HEX_D,
+    HEX_F,
+    HEX_X,
+};
+static const Hex subhexes_X[] = {
+    HEX_F,
+    HEX_Y,
+    HEX_G,
+    HEX_S,
+    HEX_Y,
+    HEX_D,
+    HEX_F,
+    HEX_P,
+};
+static const Hex subhexes_P[] = {
+    HEX_F,
+    HEX_Y,
+    HEX_G,
+    HEX_S,
+    HEX_Y,
+    HEX_D,
+    HEX_F,
+    HEX_X,
+};
+static const Hex subhexes_S[] = {
+    HEX_L,
+    HEX_P,
+    HEX_G,
+    HEX_S,
+    HEX_X,
+    HEX_D,
+    HEX_F,
+    HEX_X,
+};
+static const Hex subhexes_F[] = {
+    HEX_F,
+    HEX_P,
+    HEX_G,
+    HEX_S,
+    HEX_Y,
+    HEX_D,
+    HEX_F,
+    HEX_Y,
+};
+static const Hex subhexes_Y[] = {
+    HEX_F,
+    HEX_Y,
+    HEX_G,
+    HEX_S,
+    HEX_Y,
+    HEX_D,
+    HEX_F,
+    HEX_Y,
+};
+
+/*
+ * Shape of the Spectre itself.
+ *
+ * Vertex 0 is taken to be at the top of the Spectre's "head"; vertex
+ * 1 is the adjacent vertex, in the direction of the shorter edge of
+ * its "cloak".
+ *
+ * This array indicates how far to turn at each vertex, in 1/12 turns.
+ * All edges are the same length (counting the double-edge as two
+ * edges, which we do).
+ */
+static const int spectre_angles[14] = {
+    -3, -2, 3, -2, -3, 2, -3, 2, -3, -2, 0, -2, 3, -2,
+};
+
+/*
+ * The relative probabilities of the nine hex types, in the limit, as
+ * the expansion process is iterated more and more times. Used to
+ * choose the initial hex coordinates as if the segment of tiling came
+ * from the limiting distribution across the whole plane.
+ *
+ * This is obtained by finding the matrix that says how many hexes of
+ * each type are expanded from each starting hex, and taking the
+ * eigenvector that goes with its limiting eigenvalue.
+ */
+#define PROB_G 10000000                /* 1 */
+#define PROB_D 10000000                /* 1 */
+#define PROB_J  1270167                /* 4 - sqrt(15) */
+#define PROB_L  1270167                /* 4 - sqrt(15) */
+#define PROB_X  7459667                /* 2 sqrt(15) - 7 */
+#define PROB_P  7459667                /* 2 sqrt(15) - 7 */
+#define PROB_S 10000000                /* 1 */
+#define PROB_F 17459667                /* 2 sqrt(15) - 6 */
+#define PROB_Y 13810500                /* 13 - 3 sqrt(15) */

spectre.c

@@ -0,0 +1,594 @@
+/*
+ * Code to generate patches of the aperiodic 'spectre' tiling
+ * discovered in 2023.
+ */
+
+#include <assert.h>
+#include <string.h>
+
+#include "puzzles.h"
+#include "tree234.h"
+
+#include "spectre-internal.h"
+
+#include "spectre-tables-manual.h"
+#include "spectre-tables-auto.h"
+
+static const char *const letters =
+    #define STRINGIFY(x) #x
+    HEX_LETTERS(STRINGIFY)
+    #undef STRINGIFY
+    ;
+
+bool spectre_valid_hex_letter(char letter)
+{
+    return strchr(letters, letter) != NULL;
+}
+
+static Hex hex_from_letter(char letter)
+{
+    char buf[2];
+    buf[0] = letter;
+    buf[1] = '\0';
+    return strcspn(letters, buf);
+}
+
+static Hex hex_to_letter(unsigned char letter)
+{
+    return letters[letter];
+}
+
+struct HexData {
+    const struct MapEntry *hexmap, *hexin, *specmap, *specin;
+    const struct MapEdge *hexedges, *specedges;
+    const Hex *subhexes;
+    const struct Possibility *poss;
+    size_t nposs;
+};
+
+static const struct HexData hexdata[] = {
+    #define HEXDATA_ENTRY(x) { hexmap_##x, hexin_##x, specmap_##x, \
+            specin_##x, hexedges_##x, specedges_##x, subhexes_##x, \
+            poss_##x, lenof(poss_##x) },
+    HEX_LETTERS(HEXDATA_ENTRY)
+    #undef HEXDATA_ENTRY
+};
+
+static const struct Possibility *choose_poss(
+    random_state *rs, const struct Possibility *poss, size_t nposs)
+{
+    /*
+     * If we needed to do this _efficiently_, we'd rewrite all those
+     * tables above as cumulative frequency tables and use binary
+     * search. But this happens about log n times in a grid of area n,
+     * so it hardly matters, and it's easier to keep the tables
+     * legible.
+     */
+    unsigned long limit = 0, value;
+    size_t i;
+
+    for (i = 0; i < nposs; i++)
+        limit += poss[i].prob;
+
+    value = random_upto(rs, limit);
+
+    for (i = 0; i+1 < nposs; i++) {
+        if (value < poss[i].prob)
+            return &poss[i];
+        value -= poss[i].prob;
+    }
+
+    assert(i == nposs - 1);
+    assert(value < poss[i].prob);
+    return &poss[i];
+}
+
+SpectreCoords *spectre_coords_new(void)
+{
+    SpectreCoords *sc = snew(SpectreCoords);
+    sc->nc = sc->csize = 0;
+    sc->c = NULL;
+    return sc;
+}
+
+void spectre_coords_free(SpectreCoords *sc)
+{
+    if (sc) {
+        sfree(sc->c);
+        sfree(sc);
+    }
+}
+
+void spectre_coords_make_space(SpectreCoords *sc, size_t size)
+{
+    if (sc->csize < size) {
+        sc->csize = sc->csize * 5 / 4 + 16;
+        if (sc->csize < size)
+            sc->csize = size;
+        sc->c = sresize(sc->c, sc->csize, HexCoord);
+    }
+}
+
+SpectreCoords *spectre_coords_copy(SpectreCoords *sc_in)
+{
+    SpectreCoords *sc_out = spectre_coords_new();
+    spectre_coords_make_space(sc_out, sc_in->nc);
+    memcpy(sc_out->c, sc_in->c, sc_in->nc * sizeof(*sc_out->c));
+    sc_out->nc = sc_in->nc;
+    sc_out->index = sc_in->index;
+    sc_out->hex_colour = sc_in->hex_colour;
+    sc_out->prev_hex_colour = sc_in->prev_hex_colour;
+    sc_out->incoming_hex_edge = sc_in->incoming_hex_edge;
+    return sc_out;
+}
+
+void spectre_place(Spectre *spec, Point u, Point v, int index_of_u)
+{
+    size_t i;
+    Point disp;
+
+    /* Vector from u to v */
+    disp = point_sub(v, u);
+
+    for (i = 0; i < 14; i++) {
+        spec->vertices[(i + index_of_u) % 14] = u;
+        u = point_add(u, disp);
+        disp = point_mul(disp, point_rot(
+                             spectre_angles[(i + 1 + index_of_u) % 14]));
+    }
+}
+
+static Spectre *spectre_initial(Point u, Point v, int index_of_u,
+                                SpectreCoords *sc)
+{
+    Spectre *spec = snew(Spectre);
+    spectre_place(spec, u, v, index_of_u);
+    spec->sc = spectre_coords_copy(sc);
+    return spec;
+}
+
+static Spectre *spectre_adjacent(
+    SpectreContext *ctx, const Spectre *src_spec, unsigned src_edge)
+{
+    unsigned dst_edge;
+    Spectre *dst_spec = snew(Spectre);
+    dst_spec->sc = spectre_coords_copy(src_spec->sc);
+    spectrectx_step(ctx, dst_spec->sc, src_edge, &dst_edge);
+    spectre_place(dst_spec, src_spec->vertices[(src_edge+1) % 14],
+                  src_spec->vertices[src_edge], dst_edge);
+    return dst_spec;
+}
+
+static int spectre_cmp(void *av, void *bv)
+{
+    Spectre *a = (Spectre *)av, *b = (Spectre *)bv;
+    size_t i, j;
+
+    /* We should only ever need to compare the first two vertices of
+     * any Spectre, because those force the rest */
+    for (i = 0; i < 2; i++) {
+        for (j = 0; j < 4; j++) {
+            int ac = a->vertices[i].coeffs[j], bc = b->vertices[i].coeffs[j];
+            if (ac < bc)
+                return -1;
+            if (ac > bc)
+                return +1;
+        }
+    }
+
+    return 0;
+}
+
+static void spectre_free(Spectre *spec)
+{
+    spectre_coords_free(spec->sc);
+    sfree(spec);
+}
+
+static void spectrectx_start_vertices(SpectreContext *ctx, int orientation)
+{
+    Point minus_sqrt3 = point_add(point_rot(5), point_rot(-5));
+    Point basicedge = point_mul(point_add(point_rot(0), point_rot(-3)),
+                                point_rot(orientation));
+    Point diagonal = point_add(basicedge, point_mul(basicedge, point_rot(-3)));
+    ctx->start_vertices[0] = point_mul(diagonal, minus_sqrt3);
+    ctx->start_vertices[1] = point_add(ctx->start_vertices[0], basicedge);
+    ctx->orientation = orientation;
+}
+
+void spectrectx_init_random(SpectreContext *ctx, random_state *rs)
+{
+    const struct Possibility *poss;
+
+    ctx->rs = rs;
+    ctx->must_free_rs = false;
+    ctx->prototype = spectre_coords_new();
+    spectre_coords_make_space(ctx->prototype, 1);
+    poss = choose_poss(rs, poss_spectre, lenof(poss_spectre));
+    ctx->prototype->index = poss->lo;
+    ctx->prototype->c[0].type = poss->hi;
+    ctx->prototype->c[0].index = -1;
+    ctx->prototype->nc = 1;
+
+    /*
+     * Choose a random orientation for the starting Spectre.
+     *
+     * The obvious thing is to choose the orientation out of all 12
+     * possibilities. But we do it a more complicated way.
+     *
+     * The Spectres in a tiling can be partitioned into two
+     * equivalence classes under the relation 'orientation differs by
+     * a multiple of 1/6 turn'. One class is much more common than the
+     * other class: the 'odd'-orientation Spectres occur rarely (very
+     * like the rare reflected hats in the hats tiling).
+     *
+     * I think it's nicer to arrange that there's a consistent
+     * orientation for the _common_ class of Spectres, so that there
+     * will always be plenty of them in the 'canonical' orientation
+     * with the head upwards. So if the starting Spectre is in the
+     * even class, we pick an even orientation for it, and if it's in
+     * the odd class, we pick an odd orientation.
+     *
+     * An odd-class Spectre is easy to identify from SpectreCoords.
+     * They're precisely the ones expanded from a G hex with index 1,
+     * which means they're the ones that have index 1 _at all_.
+     */
+    spectrectx_start_vertices(ctx, random_upto(rs, 6) * 2 +
+                              ctx->prototype->index);
+
+    /* Initialiise the colouring fields deterministically but unhelpfully.
+     * spectre-test will set these up properly if it wants to */
+    ctx->prototype->hex_colour = 0;
+    ctx->prototype->prev_hex_colour = 0;
+    ctx->prototype->incoming_hex_edge = 0;
+}
+
+void spectrectx_init_from_params(
+    SpectreContext *ctx, const struct SpectrePatchParams *ps)
+{
+    size_t i;
+
+    ctx->rs = NULL;
+    ctx->must_free_rs = false;
+    ctx->prototype = spectre_coords_new();
+    spectre_coords_make_space(ctx->prototype, ps->ncoords);
+
+    ctx->prototype->index = ps->coords[0];
+    for (i = 1; i < ps->ncoords; i++)
+        ctx->prototype->c[i-1].index = ps->coords[i];
+    ctx->prototype->c[ps->ncoords-1].index = -1;
+    ctx->prototype->nc = ps->ncoords;
+
+    ctx->prototype->c[ps->ncoords-1].type = hex_from_letter(ps->final_hex);
+    for (i = ps->ncoords - 1; i-- > 0 ;) {
+        const struct HexData *h = &hexdata[ctx->prototype->c[i+1].type];
+        ctx->prototype->c[i].type = h->subhexes[ctx->prototype->c[i].index];
+    }
+
+    spectrectx_start_vertices(ctx, ps->orientation);
+
+    ctx->prototype->hex_colour = 0;
+    ctx->prototype->prev_hex_colour = 0;
+    ctx->prototype->incoming_hex_edge = 0;
+}
+
+void spectrectx_cleanup(SpectreContext *ctx)
+{
+    if (ctx->must_free_rs)
+        random_free(ctx->rs);
+    spectre_coords_free(ctx->prototype);
+}
+
+SpectreCoords *spectrectx_initial_coords(SpectreContext *ctx)
+{
+    return spectre_coords_copy(ctx->prototype);
+}
+
+/*
+ * Extend sc until it has at least n coordinates in, by copying from
+ * ctx->prototype if needed, and extending ctx->prototype if needed in
+ * order to do that.
+ */
+void spectrectx_extend_coords(SpectreContext *ctx, SpectreCoords *sc, size_t n)
+{
+    if (ctx->prototype->nc < n) {
+        spectre_coords_make_space(ctx->prototype, n);
+        while (ctx->prototype->nc < n) {
+            const struct HexData *h = &hexdata[
+                ctx->prototype->c[ctx->prototype->nc-1].type];
+            const struct Possibility *poss;
+
+            if (!ctx->rs) {
+                /*
+                 * If there's no random_state available, it must be
+                 * because we were given an explicit coordinate string
+                 * and ran off the end of it.
+                 *
+                 * The obvious thing to do here would be to make up an
+                 * answer non-randomly. But in fact there's a danger
+                 * that this leads to endless recursion within a
+                 * single coordinate step, if the hex edge we were
+                 * trying to traverse turns into another copy of
+                 * itself at the higher level. That happened in early
+                 * testing before I put the random_state in at all.
+                 *
+                 * To avoid that risk, in this situation - which
+                 * _shouldn't_ come up at all in sensibly play - we
+                 * make up a random_state, and free it when the
+                 * context goes away.
+                 */
+                ctx->rs = random_new("dummy", 5);
+                ctx->must_free_rs = true;
+            }
+
+            poss = choose_poss(ctx->rs, h->poss, h->nposs);
+            ctx->prototype->c[ctx->prototype->nc-1].index = poss->lo;
+            ctx->prototype->c[ctx->prototype->nc].type = poss->hi;
+            ctx->prototype->c[ctx->prototype->nc].index = -1;
+            ctx->prototype->nc++;
+        }
+    }
+
+    spectre_coords_make_space(sc, n);
+    while (sc->nc < n) {
+        assert(sc->c[sc->nc - 1].index == -1);
+        assert(sc->c[sc->nc - 1].type == ctx->prototype->c[sc->nc - 1].type);
+        sc->c[sc->nc - 1].index = ctx->prototype->c[sc->nc - 1].index;
+        sc->c[sc->nc].index = -1;
+        sc->c[sc->nc].type = ctx->prototype->c[sc->nc].type;
+        sc->nc++;
+    }
+}
+
+void spectrectx_step_hex(SpectreContext *ctx, SpectreCoords *sc,
+                         size_t depth, unsigned edge, unsigned *outedge)
+{
+    const struct HexData *h;
+    const struct MapEntry *m;
+
+    spectrectx_extend_coords(ctx, sc, depth+2);
+
+    assert(0 <= sc->c[depth].index);
+    assert(sc->c[depth].index < num_subhexes(sc->c[depth].type));
+    assert(0 <= edge);
+    assert(edge < 6);
+
+    h = &hexdata[sc->c[depth+1].type];
+    m = &h->hexmap[6 * sc->c[depth].index + edge];
+    if (!m->internal) {
+        unsigned recedge;
+        const struct MapEdge *me;
+        spectrectx_step_hex(ctx, sc, depth+1, m->hi, &recedge);
+        assert(recedge < 6);
+        h = &hexdata[sc->c[depth+1].type];
+        me = &h->hexedges[recedge];
+        assert(m->lo < me->len);
+        m = &h->hexin[me->startindex + me->len - 1 - m->lo];
+        assert(m->internal);
+    }
+    sc->c[depth].index = m->hi;
+    sc->c[depth].type = h->subhexes[sc->c[depth].index];
+    *outedge = m->lo;
+
+    if (depth == 0) {
+        /*
+         * Update the colouring fields to track the colour of the new
+         * hexagon.
+         */
+        unsigned char new_hex_colour;
+
+        if (!((edge ^ sc->incoming_hex_edge) & 1)) {
+            /* We're going out via the same parity of edge we came in
+             * on, so the new hex colour is the same as the previous
+             * one. */
+            new_hex_colour = sc->prev_hex_colour;
+        } else {
+            /* We're going out via the opposite parity of edge, so the
+             * new colour is the one of {0,1,2} that is neither this
+             * _nor_ the previous colour. */
+            new_hex_colour = 0+1+2 - sc->hex_colour - sc->prev_hex_colour;
+        }
+
+        sc->prev_hex_colour = sc->hex_colour;
+        sc->hex_colour = new_hex_colour;
+        sc->incoming_hex_edge = m->lo;
+    }
+}
+
+void spectrectx_step(SpectreContext *ctx, SpectreCoords *sc,
+                     unsigned edge, unsigned *outedge)
+{
+    const struct HexData *h;
+    const struct MapEntry *m;
+
+    assert(0 <= sc->index);
+    assert(sc->index < num_spectres(sc->c[0].type));
+    assert(0 <= edge);
+    assert(edge < 14);
+
+    h = &hexdata[sc->c[0].type];
+    m = &h->specmap[14 * sc->index + edge];
+
+    while (!m->internal) {
+        unsigned recedge;
+        const struct MapEdge *me;
+        spectrectx_step_hex(ctx, sc, 0, m->hi, &recedge);
+        assert(recedge < 6);
+        h = &hexdata[sc->c[0].type];
+        me = &h->specedges[recedge];
+        assert(m->lo < me->len);
+        m = &h->specin[me->startindex + me->len - 1 - m->lo];
+    }
+    sc->index = m->hi;
+    *outedge = m->lo;
+}
+
+void spectrectx_generate(SpectreContext *ctx,
+                         bool (*callback)(void *cbctx, const Spectre *spec),
+                         void *cbctx)
+{
+    tree234 *placed = newtree234(spectre_cmp);
+    Spectre *qhead = NULL, *qtail = NULL;
+
+    {
+        SpectreCoords *sc = spectrectx_initial_coords(ctx);
+        Spectre *spec = spectre_initial(ctx->start_vertices[0],
+                                        ctx->start_vertices[1], 0, sc);
+        spectre_coords_free(sc);
+
+        add234(placed, spec);
+
+        spec->next = NULL;
+
+        if (callback(cbctx, spec))
+            qhead = qtail = spec;
+    }
+
+    while (qhead) {
+        unsigned edge;
+        Spectre *spec = qhead;
+
+        for (edge = 0; edge < 14; edge++) {
+            Spectre *new_spec;
+
+            new_spec = spectre_adjacent(ctx, spec, edge);
+
+            if (find234(placed, new_spec, NULL)) {
+                spectre_free(new_spec);
+                continue;
+            }
+
+            if (!callback(cbctx, new_spec)) {
+                spectre_free(new_spec);
+                continue;
+            }
+
+            add234(placed, new_spec);
+            qtail->next = new_spec;
+            qtail = new_spec;
+            new_spec->next = NULL;
+        }
+
+        qhead = qhead->next;
+    }
+
+    {
+        Spectre *spec;
+        while ((spec = delpos234(placed, 0)) != NULL)
+            spectre_free(spec);
+        freetree234(placed);
+    }
+}
+
+const char *spectre_tiling_params_invalid(
+    const struct SpectrePatchParams *params)
+{
+    size_t i;
+    Hex h;
+
+    if (params->ncoords == 0)
+        return "expected at least one numeric coordinate";
+    if (!spectre_valid_hex_letter(params->final_hex))
+        return "invalid final hexagon type";
+
+    h = hex_from_letter(params->final_hex);
+    for (i = params->ncoords; i-- > 0 ;) {
+        unsigned limit = (i == 0) ? num_spectres(h) : num_subhexes(h);
+        if (params->coords[i] >= limit)
+            return "coordinate out of range";
+
+        if (i > 0)
+            h = hexdata[h].subhexes[params->coords[i]];
+    }
+
+    return NULL;
+}
+
+struct SpectreCallbackContext {
+    int xoff, yoff;
+    Coord xmin, xmax, ymin, ymax;
+
+    spectre_tile_callback_fn external_cb;
+    void *external_cbctx;
+};
+
+static bool spectre_internal_callback(void *vctx, const Spectre *spec)
+{
+    struct SpectreCallbackContext *ctx = (struct SpectreCallbackContext *)vctx;
+    size_t i;
+    int output_coords[4*14];
+
+    for (i = 0; i < 14; i++) {
+        Point p = spec->vertices[i];
+        Coord x = point_x(p), y = point_y(p);
+        if (coord_cmp(x, ctx->xmin) < 0 || coord_cmp(x, ctx->xmax) > 0 ||
+            coord_cmp(y, ctx->ymin) < 0 || coord_cmp(y, ctx->ymax) > 0)
+            return false;
+
+        output_coords[4*i + 0] = ctx->xoff + x.c1;
+        output_coords[4*i + 1] = x.cr3;
+        output_coords[4*i + 2] = ctx->yoff - y.c1;
+        output_coords[4*i + 3] = -y.cr3;
+    }
+
+    if (ctx->external_cb)
+        ctx->external_cb(ctx->external_cbctx, output_coords);
+
+    return true;
+}
+
+static void spectre_set_bounds(struct SpectreCallbackContext *cbctx,
+                               int w, int h)
+{
+    cbctx->xoff = w/2;
+    cbctx->yoff = h/2;
+    cbctx->xmin.c1 = -cbctx->xoff;
+    cbctx->xmax.c1 = -cbctx->xoff + w;
+    cbctx->ymin.c1 = cbctx->yoff - h;
+    cbctx->ymax.c1 = cbctx->yoff;
+    cbctx->xmin.cr3 = 0;
+    cbctx->xmax.cr3 = 0;
+    cbctx->ymin.cr3 = 0;
+    cbctx->ymax.cr3 = 0;
+}
+
+void spectre_tiling_randomise(struct SpectrePatchParams *ps, int w, int h,
+                              random_state *rs)
+{
+    SpectreContext ctx[1];
+    struct SpectreCallbackContext cbctx[1];
+    size_t i;
+
+    spectre_set_bounds(cbctx, w, h);
+    cbctx->external_cb = NULL;
+    cbctx->external_cbctx = NULL;
+
+    spectrectx_init_random(ctx, rs);
+    spectrectx_generate(ctx, spectre_internal_callback, cbctx);
+
+    ps->orientation = ctx->orientation;
+    ps->ncoords = ctx->prototype->nc;
+    ps->coords = snewn(ps->ncoords, unsigned char);
+    ps->coords[0] = ctx->prototype->index;
+    for (i = 1; i < ps->ncoords; i++)
+        ps->coords[i] = ctx->prototype->c[i-1].index;
+    ps->final_hex = hex_to_letter(ctx->prototype->c[ps->ncoords-1].type);
+
+    spectrectx_cleanup(ctx);
+}
+
+void spectre_tiling_generate(
+    const struct SpectrePatchParams *params, int w, int h,
+    spectre_tile_callback_fn external_cb, void *external_cbctx)
+{
+    SpectreContext ctx[1];
+    struct SpectreCallbackContext cbctx[1];
+
+    spectre_set_bounds(cbctx, w, h);
+    cbctx->external_cb = external_cb;
+    cbctx->external_cbctx = external_cbctx;
+
+    spectrectx_init_from_params(ctx, params);
+    spectrectx_generate(ctx, spectre_internal_callback, cbctx);
+    spectrectx_cleanup(ctx);
+}

spectre.h

@@ -0,0 +1,72 @@
+#ifndef PUZZLES_SPECTRE_H
+#define PUZZLES_SPECTRE_H
+
+struct SpectrePatchParams {
+    /*
+     * A patch of Spectre tiling is identified by giving
+     *
+     *  - the coordinates of the central Spectre, using a
+     *    combinatorial coordinate system similar to the Hat tiling in
+     *    hat.h
+     *
+     *  - the orientation of that Spectre, as a number from 0 to 11 (a
+     *    multiple of 30 degrees), with 0 representing the 'head' of
+     *    the Spectre facing upwards, and numbers counting
+     *    anticlockwise.
+     *
+     * Coordinates are a sequence of small non-negative integers. The
+     * valid range for each coordinate depends on the next coordinate,
+     * or on final_hex if it's the last one in the list. The largest
+     * valid range is {0,...,7}.
+     *
+     * 'final_hex' is one of the letters GDJLXPSFY.
+     * spectre_valid_hex_letter() will check that.
+     */
+    int orientation;
+    size_t ncoords;
+    unsigned char *coords;
+    char final_hex;
+};
+
+bool spectre_valid_hex_letter(char c);
+
+/*
+ * Fill in SpectrePatchParams with a randomly selected set of
+ * coordinates, in enough detail to generate a patch of tiling filling
+ * a w x h area. The unit of measurement is 1/(2*sqrt(2)) of a Spectre
+ * edge, i.e. such that a single Spectre edge at 45 degrees would
+ * correspond to the vector (2,2).
+ *
+ * The 'coords' field of the structure will be filled in with a new
+ * dynamically allocated array. Any previous pointer in that field
+ * will be overwritten.
+ */
+void spectre_tiling_randomise(struct SpectrePatchParams *params, int w, int h,
+                              random_state *rs);
+
+/*
+ * Validate a SpectrePatchParams to ensure it contains no illegal
+ * coordinates. Returns NULL if it's acceptable, or an error string if
+ * not.
+ */
+const char *spectre_tiling_params_invalid(
+    const struct SpectrePatchParams *params);
+
+/*
+ * Generate the actual set of Spectre tiles from a SpectrePatchParams,
+ * passing each one to a callback. The callback receives the vertices
+ * of each point, in the form of an array of 4*14 integers. Each
+ * vertex is represented by four consecutive integers in this array,
+ * with the first two giving the x coordinate and the last two the y
+ * coordinate. Each pair of integers a,b represent a single coordinate
+ * whose value is a + b*sqrt(3). The unit of measurement is as above.
+ */
+typedef void (*spectre_tile_callback_fn)(void *ctx, const int *coords);
+
+#define SPECTRE_NVERTICES 14
+
+void spectre_tiling_generate(
+    const struct SpectrePatchParams *params, int w, int h,
+    spectre_tile_callback_fn cb, void *cbctx);
+
+#endif