2 * (c) Lambros Lambrou 2008
4 * Code for working with general grids, which can be any planar graph
5 * with faces, edges and vertices (dots). Includes generators for a few
6 * types of grid, including square, hexagonal, triangular and others.
22 /* Debugging options */
28 /* ----------------------------------------------------------------------
29 * Deallocate or dereference a grid
31 void grid_free(grid *g)
36 if (g->refcount == 0) {
38 for (i = 0; i < g->num_faces; i++) {
39 sfree(g->faces[i].dots);
40 sfree(g->faces[i].edges);
42 for (i = 0; i < g->num_dots; i++) {
43 sfree(g->dots[i].faces);
44 sfree(g->dots[i].edges);
53 /* Used by the other grid generators. Create a brand new grid with nothing
54 * initialised (all lists are NULL) */
55 static grid *grid_empty(void)
61 g->num_faces = g->num_edges = g->num_dots = 0;
63 g->lowest_x = g->lowest_y = g->highest_x = g->highest_y = 0;
67 /* Helper function to calculate perpendicular distance from
68 * a point P to a line AB. A and B mustn't be equal here.
70 * Well-known formula for area A of a triangle:
72 * 2A = determinant of matrix | px ax bx |
75 * Also well-known: 2A = base * height
76 * = perpendicular distance * line-length.
78 * Combining gives: distance = determinant / line-length(a,b)
80 static double point_line_distance(long px, long py,
84 long det = ax*by - bx*ay + bx*py - px*by + px*ay - ax*py;
87 len = sqrt(SQ(ax - bx) + SQ(ay - by));
91 /* Determine nearest edge to where the user clicked.
92 * (x, y) is the clicked location, converted to grid coordinates.
93 * Returns the nearest edge, or NULL if no edge is reasonably
96 * Just judging edges by perpendicular distance is not quite right -
97 * the edge might be "off to one side". So we insist that the triangle
98 * with (x,y) has acute angles at the edge's dots.
105 * | edge2 is OK, but edge1 is not, even though
106 * | edge1 is perpendicularly closer to (x,y)
110 grid_edge *grid_nearest_edge(grid *g, int x, int y)
112 grid_edge *best_edge;
113 double best_distance = 0;
118 for (i = 0; i < g->num_edges; i++) {
119 grid_edge *e = &g->edges[i];
120 long e2; /* squared length of edge */
121 long a2, b2; /* squared lengths of other sides */
124 /* See if edge e is eligible - the triangle must have acute angles
125 * at the edge's dots.
126 * Pythagoras formula h^2 = a^2 + b^2 detects right-angles,
127 * so detect acute angles by testing for h^2 < a^2 + b^2 */
128 e2 = SQ((long)e->dot1->x - (long)e->dot2->x) + SQ((long)e->dot1->y - (long)e->dot2->y);
129 a2 = SQ((long)e->dot1->x - (long)x) + SQ((long)e->dot1->y - (long)y);
130 b2 = SQ((long)e->dot2->x - (long)x) + SQ((long)e->dot2->y - (long)y);
131 if (a2 >= e2 + b2) continue;
132 if (b2 >= e2 + a2) continue;
134 /* e is eligible so far. Now check the edge is reasonably close
135 * to where the user clicked. Don't want to toggle an edge if the
136 * click was way off the grid.
137 * There is room for experimentation here. We could check the
138 * perpendicular distance is within a certain fraction of the length
139 * of the edge. That amounts to testing a rectangular region around
141 * Alternatively, we could check that the angle at the point is obtuse.
142 * That would amount to testing a circular region with the edge as
144 dist = point_line_distance((long)x, (long)y,
145 (long)e->dot1->x, (long)e->dot1->y,
146 (long)e->dot2->x, (long)e->dot2->y);
147 /* Is dist more than half edge length ? */
148 if (4 * SQ(dist) > e2)
151 if (best_edge == NULL || dist < best_distance) {
153 best_distance = dist;
159 /* ----------------------------------------------------------------------
169 #define FACE_COLOUR "red"
170 #define EDGE_COLOUR "blue"
171 #define DOT_COLOUR "black"
173 static void grid_output_svg(FILE *fp, grid *g, int which)
178 <?xml version=\"1.0\" encoding=\"UTF-8\" standalone=\"no\"?>\n\
179 <!DOCTYPE svg PUBLIC \"-//W3C//DTD SVG 20010904//EN\"\n\
180 \"http://www.w3.org/TR/2001/REC-SVG-20010904/DTD/svg10.dtd\">\n\
182 <svg xmlns=\"http://www.w3.org/2000/svg\"\n\
183 xmlns:xlink=\"http://www.w3.org/1999/xlink\">\n\n");
185 if (which & SVG_FACES) {
186 fprintf(fp, "<g>\n");
187 for (i = 0; i < g->num_faces; i++) {
188 grid_face *f = g->faces + i;
189 fprintf(fp, "<polygon points=\"");
190 for (j = 0; j < f->order; j++) {
191 grid_dot *d = f->dots[j];
192 fprintf(fp, "%s%d,%d", (j == 0) ? "" : " ",
195 fprintf(fp, "\" style=\"fill: %s; fill-opacity: 0.2; stroke: %s\" />\n",
196 FACE_COLOUR, FACE_COLOUR);
198 fprintf(fp, "</g>\n");
200 if (which & SVG_EDGES) {
201 fprintf(fp, "<g>\n");
202 for (i = 0; i < g->num_edges; i++) {
203 grid_edge *e = g->edges + i;
204 grid_dot *d1 = e->dot1, *d2 = e->dot2;
206 fprintf(fp, "<line x1=\"%d\" y1=\"%d\" x2=\"%d\" y2=\"%d\" "
207 "style=\"stroke: %s\" />\n",
208 d1->x, d1->y, d2->x, d2->y, EDGE_COLOUR);
210 fprintf(fp, "</g>\n");
213 if (which & SVG_DOTS) {
214 fprintf(fp, "<g>\n");
215 for (i = 0; i < g->num_dots; i++) {
216 grid_dot *d = g->dots + i;
217 fprintf(fp, "<ellipse cx=\"%d\" cy=\"%d\" rx=\"%d\" ry=\"%d\" fill=\"%s\" />",
218 d->x, d->y, g->tilesize/20, g->tilesize/20, DOT_COLOUR);
220 fprintf(fp, "</g>\n");
223 fprintf(fp, "</svg>\n");
230 static void grid_try_svg(grid *g, int which)
232 char *svg = getenv("PUZZLES_SVG_GRID");
234 FILE *svgf = fopen(svg, "w");
236 grid_output_svg(svgf, g, which);
239 fprintf(stderr, "Unable to open file `%s': %s", svg, strerror(errno));
245 /* Show the basic grid information, before doing grid_make_consistent */
246 static void grid_debug_basic(grid *g)
248 /* TODO: Maybe we should generate an SVG image of the dots and lines
249 * of the grid here, before grid_make_consistent.
250 * Would help with debugging grid generation. */
253 printf("--- Basic Grid Data ---\n");
254 for (i = 0; i < g->num_faces; i++) {
255 grid_face *f = g->faces + i;
256 printf("Face %d: dots[", i);
258 for (j = 0; j < f->order; j++) {
259 grid_dot *d = f->dots[j];
260 printf("%s%d", j ? "," : "", (int)(d - g->dots));
266 grid_try_svg(g, SVG_FACES);
270 /* Show the derived grid information, computed by grid_make_consistent */
271 static void grid_debug_derived(grid *g)
276 printf("--- Derived Grid Data ---\n");
277 for (i = 0; i < g->num_edges; i++) {
278 grid_edge *e = g->edges + i;
279 printf("Edge %d: dots[%d,%d] faces[%d,%d]\n",
280 i, (int)(e->dot1 - g->dots), (int)(e->dot2 - g->dots),
281 e->face1 ? (int)(e->face1 - g->faces) : -1,
282 e->face2 ? (int)(e->face2 - g->faces) : -1);
285 for (i = 0; i < g->num_faces; i++) {
286 grid_face *f = g->faces + i;
288 printf("Face %d: faces[", i);
289 for (j = 0; j < f->order; j++) {
290 grid_edge *e = f->edges[j];
291 grid_face *f2 = (e->face1 == f) ? e->face2 : e->face1;
292 printf("%s%d", j ? "," : "", f2 ? (int)(f2 - g->faces) : -1);
297 for (i = 0; i < g->num_dots; i++) {
298 grid_dot *d = g->dots + i;
300 printf("Dot %d: dots[", i);
301 for (j = 0; j < d->order; j++) {
302 grid_edge *e = d->edges[j];
303 grid_dot *d2 = (e->dot1 == d) ? e->dot2 : e->dot1;
304 printf("%s%d", j ? "," : "", (int)(d2 - g->dots));
307 for (j = 0; j < d->order; j++) {
308 grid_face *f = d->faces[j];
309 printf("%s%d", j ? "," : "", f ? (int)(f - g->faces) : -1);
315 grid_try_svg(g, SVG_DOTS | SVG_EDGES | SVG_FACES);
319 /* Helper function for building incomplete-edges list in
320 * grid_make_consistent() */
321 static int grid_edge_bydots_cmpfn(void *v1, void *v2)
327 /* Pointer subtraction is valid here, because all dots point into the
328 * same dot-list (g->dots).
329 * Edges are not "normalised" - the 2 dots could be stored in any order,
330 * so we need to take this into account when comparing edges. */
332 /* Compare first dots */
333 da = (a->dot1 < a->dot2) ? a->dot1 : a->dot2;
334 db = (b->dot1 < b->dot2) ? b->dot1 : b->dot2;
337 /* Compare last dots */
338 da = (a->dot1 < a->dot2) ? a->dot2 : a->dot1;
339 db = (b->dot1 < b->dot2) ? b->dot2 : b->dot1;
347 * 'Vigorously trim' a grid, by which I mean deleting any isolated or
348 * uninteresting faces. By which, in turn, I mean: ensure that the
349 * grid is composed solely of faces adjacent to at least one
350 * 'landlocked' dot (i.e. one not in contact with the infinite
351 * exterior face), and that all those dots are in a single connected
354 * This function operates on, and returns, a grid satisfying the
355 * preconditions to grid_make_consistent() rather than the
356 * postconditions. (So call it first.)
358 static void grid_trim_vigorously(grid *g)
360 int *dotpairs, *faces, *dots;
362 int i, j, k, size, newfaces, newdots;
365 * First construct a matrix in which each ordered pair of dots is
366 * mapped to the index of the face in which those dots occur in
369 dotpairs = snewn(g->num_dots * g->num_dots, int);
370 for (i = 0; i < g->num_dots; i++)
371 for (j = 0; j < g->num_dots; j++)
372 dotpairs[i*g->num_dots+j] = -1;
373 for (i = 0; i < g->num_faces; i++) {
374 grid_face *f = g->faces + i;
375 int dot0 = f->dots[f->order-1] - g->dots;
376 for (j = 0; j < f->order; j++) {
377 int dot1 = f->dots[j] - g->dots;
378 dotpairs[dot0 * g->num_dots + dot1] = i;
384 * Now we can identify landlocked dots: they're the ones all of
385 * whose edges have a mirror-image counterpart in this matrix.
387 dots = snewn(g->num_dots, int);
388 for (i = 0; i < g->num_dots; i++) {
390 for (j = 0; j < g->num_dots; j++) {
391 if ((dotpairs[i*g->num_dots+j] >= 0) ^
392 (dotpairs[j*g->num_dots+i] >= 0))
393 dots[i] = FALSE; /* non-duplicated edge: coastal dot */
398 * Now identify connected pairs of landlocked dots, and form a dsf
401 dsf = snew_dsf(g->num_dots);
402 for (i = 0; i < g->num_dots; i++)
403 for (j = 0; j < i; j++)
404 if (dots[i] && dots[j] &&
405 dotpairs[i*g->num_dots+j] >= 0 &&
406 dotpairs[j*g->num_dots+i] >= 0)
407 dsf_merge(dsf, i, j);
410 * Now look for the largest component.
414 for (i = 0; i < g->num_dots; i++) {
416 if (dots[i] && dsf_canonify(dsf, i) == i &&
417 (newsize = dsf_size(dsf, i)) > size) {
424 * Work out which faces we're going to keep (precisely those with
425 * at least one dot in the same connected component as j) and
426 * which dots (those required by any face we're keeping).
428 * At this point we reuse the 'dots' array to indicate the dots
429 * we're keeping, rather than the ones that are landlocked.
431 faces = snewn(g->num_faces, int);
432 for (i = 0; i < g->num_faces; i++)
434 for (i = 0; i < g->num_dots; i++)
436 for (i = 0; i < g->num_faces; i++) {
437 grid_face *f = g->faces + i;
439 for (k = 0; k < f->order; k++)
440 if (dsf_canonify(dsf, f->dots[k] - g->dots) == j)
444 for (k = 0; k < f->order; k++)
445 dots[f->dots[k]-g->dots] = TRUE;
450 * Work out the new indices of those faces and dots, when we
451 * compact the arrays containing them.
453 for (i = newfaces = 0; i < g->num_faces; i++)
454 faces[i] = (faces[i] ? newfaces++ : -1);
455 for (i = newdots = 0; i < g->num_dots; i++)
456 dots[i] = (dots[i] ? newdots++ : -1);
459 * Free the dynamically allocated 'dots' pointer lists in faces
460 * we're going to discard.
462 for (i = 0; i < g->num_faces; i++)
464 sfree(g->faces[i].dots);
467 * Go through and compact the arrays.
469 for (i = 0; i < g->num_dots; i++)
471 grid_dot *dnew = g->dots + dots[i], *dold = g->dots + i;
472 *dnew = *dold; /* structure copy */
474 for (i = 0; i < g->num_faces; i++)
476 grid_face *fnew = g->faces + faces[i], *fold = g->faces + i;
477 *fnew = *fold; /* structure copy */
478 for (j = 0; j < fnew->order; j++) {
480 * Reindex the dots in this face.
482 k = fnew->dots[j] - g->dots;
483 fnew->dots[j] = g->dots + dots[k];
486 g->num_faces = newfaces;
487 g->num_dots = newdots;
495 /* Input: grid has its dots and faces initialised:
496 * - dots have (optionally) x and y coordinates, but no edges or faces
497 * (pointers are NULL).
498 * - edges not initialised at all
499 * - faces initialised and know which dots they have (but no edges yet). The
500 * dots around each face are assumed to be clockwise.
502 * Output: grid is complete and valid with all relationships defined.
504 static void grid_make_consistent(grid *g)
507 tree234 *incomplete_edges;
508 grid_edge *next_new_edge; /* Where new edge will go into g->edges */
512 /* ====== Stage 1 ======
516 /* We know how many dots and faces there are, so we can find the exact
517 * number of edges from Euler's polyhedral formula: F + V = E + 2 .
518 * We use "-1", not "-2" here, because Euler's formula includes the
519 * infinite face, which we don't count. */
520 g->num_edges = g->num_faces + g->num_dots - 1;
521 g->edges = snewn(g->num_edges, grid_edge);
522 next_new_edge = g->edges;
524 /* Iterate over faces, and over each face's dots, generating edges as we
525 * go. As we find each new edge, we can immediately fill in the edge's
526 * dots, but only one of the edge's faces. Later on in the iteration, we
527 * will find the same edge again (unless it's on the border), but we will
528 * know the other face.
529 * For efficiency, maintain a list of the incomplete edges, sorted by
531 incomplete_edges = newtree234(grid_edge_bydots_cmpfn);
532 for (i = 0; i < g->num_faces; i++) {
533 grid_face *f = g->faces + i;
535 for (j = 0; j < f->order; j++) {
536 grid_edge e; /* fake edge for searching */
537 grid_edge *edge_found;
542 e.dot2 = f->dots[j2];
543 /* Use del234 instead of find234, because we always want to
544 * remove the edge if found */
545 edge_found = del234(incomplete_edges, &e);
547 /* This edge already added, so fill out missing face.
548 * Edge is already removed from incomplete_edges. */
549 edge_found->face2 = f;
551 assert(next_new_edge - g->edges < g->num_edges);
552 next_new_edge->dot1 = e.dot1;
553 next_new_edge->dot2 = e.dot2;
554 next_new_edge->face1 = f;
555 next_new_edge->face2 = NULL; /* potentially infinite face */
556 add234(incomplete_edges, next_new_edge);
561 freetree234(incomplete_edges);
563 /* ====== Stage 2 ======
564 * For each face, build its edge list.
567 /* Allocate space for each edge list. Can do this, because each face's
568 * edge-list is the same size as its dot-list. */
569 for (i = 0; i < g->num_faces; i++) {
570 grid_face *f = g->faces + i;
572 f->edges = snewn(f->order, grid_edge*);
573 /* Preload with NULLs, to help detect potential bugs. */
574 for (j = 0; j < f->order; j++)
578 /* Iterate over each edge, and over both its faces. Add this edge to
579 * the face's edge-list, after finding where it should go in the
581 for (i = 0; i < g->num_edges; i++) {
582 grid_edge *e = g->edges + i;
584 for (j = 0; j < 2; j++) {
585 grid_face *f = j ? e->face2 : e->face1;
587 if (f == NULL) continue;
588 /* Find one of the dots around the face */
589 for (k = 0; k < f->order; k++) {
590 if (f->dots[k] == e->dot1)
591 break; /* found dot1 */
593 assert(k != f->order); /* Must find the dot around this face */
595 /* Labelling scheme: as we walk clockwise around the face,
596 * starting at dot0 (f->dots[0]), we hit:
597 * (dot0), edge0, dot1, edge1, dot2,...
607 * Therefore, edgeK joins dotK and dot{K+1}
610 /* Around this face, either the next dot or the previous dot
611 * must be e->dot2. Otherwise the edge is wrong. */
615 if (f->dots[k2] == e->dot2) {
616 /* dot1(k) and dot2(k2) go clockwise around this face, so add
617 * this edge at position k (see diagram). */
618 assert(f->edges[k] == NULL);
622 /* Try previous dot */
626 if (f->dots[k2] == e->dot2) {
627 /* dot1(k) and dot2(k2) go anticlockwise around this face. */
628 assert(f->edges[k2] == NULL);
632 assert(!"Grid broken: bad edge-face relationship");
636 /* ====== Stage 3 ======
637 * For each dot, build its edge-list and face-list.
640 /* We don't know how many edges/faces go around each dot, so we can't
641 * allocate the right space for these lists. Pre-compute the sizes by
642 * iterating over each edge and recording a tally against each dot. */
643 for (i = 0; i < g->num_dots; i++) {
644 g->dots[i].order = 0;
646 for (i = 0; i < g->num_edges; i++) {
647 grid_edge *e = g->edges + i;
651 /* Now we have the sizes, pre-allocate the edge and face lists. */
652 for (i = 0; i < g->num_dots; i++) {
653 grid_dot *d = g->dots + i;
655 assert(d->order >= 2); /* sanity check */
656 d->edges = snewn(d->order, grid_edge*);
657 d->faces = snewn(d->order, grid_face*);
658 for (j = 0; j < d->order; j++) {
663 /* For each dot, need to find a face that touches it, so we can seed
664 * the edge-face-edge-face process around each dot. */
665 for (i = 0; i < g->num_faces; i++) {
666 grid_face *f = g->faces + i;
668 for (j = 0; j < f->order; j++) {
669 grid_dot *d = f->dots[j];
673 /* Each dot now has a face in its first slot. Generate the remaining
674 * faces and edges around the dot, by searching both clockwise and
675 * anticlockwise from the first face. Need to do both directions,
676 * because of the possibility of hitting the infinite face, which
677 * blocks progress. But there's only one such face, so we will
678 * succeed in finding every edge and face this way. */
679 for (i = 0; i < g->num_dots; i++) {
680 grid_dot *d = g->dots + i;
681 int current_face1 = 0; /* ascends clockwise */
682 int current_face2 = 0; /* descends anticlockwise */
684 /* Labelling scheme: as we walk clockwise around the dot, starting
685 * at face0 (d->faces[0]), we hit:
686 * (face0), edge0, face1, edge1, face2,...
698 * So, for example, face1 should be joined to edge0 and edge1,
699 * and those edges should appear in an anticlockwise sense around
700 * that face (see diagram). */
702 /* clockwise search */
704 grid_face *f = d->faces[current_face1];
708 /* find dot around this face */
709 for (j = 0; j < f->order; j++) {
713 assert(j != f->order); /* must find dot */
715 /* Around f, required edge is anticlockwise from the dot. See
716 * the other labelling scheme higher up, for why we subtract 1
722 d->edges[current_face1] = e; /* set edge */
724 if (current_face1 == d->order)
728 d->faces[current_face1] =
729 (e->face1 == f) ? e->face2 : e->face1;
730 if (d->faces[current_face1] == NULL)
731 break; /* cannot progress beyond infinite face */
734 /* If the clockwise search made it all the way round, don't need to
735 * bother with the anticlockwise search. */
736 if (current_face1 == d->order)
737 continue; /* this dot is complete, move on to next dot */
739 /* anticlockwise search */
741 grid_face *f = d->faces[current_face2];
745 /* find dot around this face */
746 for (j = 0; j < f->order; j++) {
750 assert(j != f->order); /* must find dot */
752 /* Around f, required edge is clockwise from the dot. */
756 if (current_face2 == -1)
757 current_face2 = d->order - 1;
758 d->edges[current_face2] = e; /* set edge */
761 if (current_face2 == current_face1)
763 d->faces[current_face2] =
764 (e->face1 == f) ? e->face2 : e->face1;
765 /* There's only 1 infinite face, so we must get all the way
766 * to current_face1 before we hit it. */
767 assert(d->faces[current_face2]);
771 /* ====== Stage 4 ======
772 * Compute other grid settings
775 /* Bounding rectangle */
776 for (i = 0; i < g->num_dots; i++) {
777 grid_dot *d = g->dots + i;
779 g->lowest_x = g->highest_x = d->x;
780 g->lowest_y = g->highest_y = d->y;
782 g->lowest_x = min(g->lowest_x, d->x);
783 g->highest_x = max(g->highest_x, d->x);
784 g->lowest_y = min(g->lowest_y, d->y);
785 g->highest_y = max(g->highest_y, d->y);
789 grid_debug_derived(g);
792 /* Helpers for making grid-generation easier. These functions are only
793 * intended for use during grid generation. */
795 /* Comparison function for the (tree234) sorted dot list */
796 static int grid_point_cmp_fn(void *v1, void *v2)
801 return p2->y - p1->y;
803 return p2->x - p1->x;
805 /* Add a new face to the grid, with its dot list allocated.
806 * Assumes there's enough space allocated for the new face in grid->faces */
807 static void grid_face_add_new(grid *g, int face_size)
810 grid_face *new_face = g->faces + g->num_faces;
811 new_face->order = face_size;
812 new_face->dots = snewn(face_size, grid_dot*);
813 for (i = 0; i < face_size; i++)
814 new_face->dots[i] = NULL;
815 new_face->edges = NULL;
816 new_face->has_incentre = FALSE;
819 /* Assumes dot list has enough space */
820 static grid_dot *grid_dot_add_new(grid *g, int x, int y)
822 grid_dot *new_dot = g->dots + g->num_dots;
824 new_dot->edges = NULL;
825 new_dot->faces = NULL;
831 /* Retrieve a dot with these (x,y) coordinates. Either return an existing dot
832 * in the dot_list, or add a new dot to the grid (and the dot_list) and
834 * Assumes g->dots has enough capacity allocated */
835 static grid_dot *grid_get_dot(grid *g, tree234 *dot_list, int x, int y)
844 ret = find234(dot_list, &test, NULL);
848 ret = grid_dot_add_new(g, x, y);
849 add234(dot_list, ret);
853 /* Sets the last face of the grid to include this dot, at this position
854 * around the face. Assumes num_faces is at least 1 (a new face has
855 * previously been added, with the required number of dots allocated) */
856 static void grid_face_set_dot(grid *g, grid_dot *d, int position)
858 grid_face *last_face = g->faces + g->num_faces - 1;
859 last_face->dots[position] = d;
863 * Helper routines for grid_find_incentre.
865 static int solve_2x2_matrix(double mx[4], double vin[2], double vout[2])
869 det = (mx[0]*mx[3] - mx[1]*mx[2]);
873 inv[0] = mx[3] / det;
874 inv[1] = -mx[1] / det;
875 inv[2] = -mx[2] / det;
876 inv[3] = mx[0] / det;
878 vout[0] = inv[0]*vin[0] + inv[1]*vin[1];
879 vout[1] = inv[2]*vin[0] + inv[3]*vin[1];
883 static int solve_3x3_matrix(double mx[9], double vin[3], double vout[3])
888 det = (mx[0]*mx[4]*mx[8] + mx[1]*mx[5]*mx[6] + mx[2]*mx[3]*mx[7] -
889 mx[0]*mx[5]*mx[7] - mx[1]*mx[3]*mx[8] - mx[2]*mx[4]*mx[6]);
893 inv[0] = (mx[4]*mx[8] - mx[5]*mx[7]) / det;
894 inv[1] = (mx[2]*mx[7] - mx[1]*mx[8]) / det;
895 inv[2] = (mx[1]*mx[5] - mx[2]*mx[4]) / det;
896 inv[3] = (mx[5]*mx[6] - mx[3]*mx[8]) / det;
897 inv[4] = (mx[0]*mx[8] - mx[2]*mx[6]) / det;
898 inv[5] = (mx[2]*mx[3] - mx[0]*mx[5]) / det;
899 inv[6] = (mx[3]*mx[7] - mx[4]*mx[6]) / det;
900 inv[7] = (mx[1]*mx[6] - mx[0]*mx[7]) / det;
901 inv[8] = (mx[0]*mx[4] - mx[1]*mx[3]) / det;
903 vout[0] = inv[0]*vin[0] + inv[1]*vin[1] + inv[2]*vin[2];
904 vout[1] = inv[3]*vin[0] + inv[4]*vin[1] + inv[5]*vin[2];
905 vout[2] = inv[6]*vin[0] + inv[7]*vin[1] + inv[8]*vin[2];
910 void grid_find_incentre(grid_face *f)
912 double xbest, ybest, bestdist;
914 grid_dot *edgedot1[3], *edgedot2[3];
922 * Find the point in the polygon with the maximum distance to any
925 * Such a point must exist which is in contact with at least three
926 * edges and/or vertices. (Proof: if it's only in contact with two
927 * edges and/or vertices, it can't even be at a _local_ maximum -
928 * any such circle can always be expanded in some direction.) So
929 * we iterate through all 3-subsets of the combined set of edges
930 * and vertices; for each subset we generate one or two candidate
931 * points that might be the incentre, and then we vet each one to
932 * see if it's inside the polygon and what its maximum radius is.
934 * (There's one case which this algorithm will get noticeably
935 * wrong, and that's when a continuum of equally good answers
936 * exists due to parallel edges. Consider a long thin rectangle,
937 * for instance, or a parallelogram. This algorithm will pick a
938 * point near one end, and choose the end arbitrarily; obviously a
939 * nicer point to choose would be in the centre. To fix this I
940 * would have to introduce a special-case system which detected
941 * parallel edges in advance, set aside all candidate points
942 * generated using both edges in a parallel pair, and generated
943 * some additional candidate points half way between them. Also,
944 * of course, I'd have to cope with rounding error making such a
945 * point look worse than one of its endpoints. So I haven't done
946 * this for the moment, and will cross it if necessary when I come
949 * We don't actually iterate literally over _edges_, in the sense
950 * of grid_edge structures. Instead, we fill in edgedot1[] and
951 * edgedot2[] with a pair of dots adjacent in the face's list of
952 * vertices. This ensures that we get the edges in consistent
953 * orientation, which we could not do from the grid structure
954 * alone. (A moment's consideration of an order-3 vertex should
955 * make it clear that if a notional arrow was written on each
956 * edge, _at least one_ of the three faces bordering that vertex
957 * would have to have the two arrows tip-to-tip or tail-to-tail
958 * rather than tip-to-tail.)
964 for (i = 0; i+2 < 2*f->order; i++) {
966 edgedot1[nedges] = f->dots[i];
967 edgedot2[nedges++] = f->dots[(i+1)%f->order];
969 dots[ndots++] = f->dots[i - f->order];
971 for (j = i+1; j+1 < 2*f->order; j++) {
973 edgedot1[nedges] = f->dots[j];
974 edgedot2[nedges++] = f->dots[(j+1)%f->order];
976 dots[ndots++] = f->dots[j - f->order];
978 for (k = j+1; k < 2*f->order; k++) {
979 double cx[2], cy[2]; /* candidate positions */
980 int cn = 0; /* number of candidates */
983 edgedot1[nedges] = f->dots[k];
984 edgedot2[nedges++] = f->dots[(k+1)%f->order];
986 dots[ndots++] = f->dots[k - f->order];
989 * Find a point, or pair of points, equidistant from
990 * all the specified edges and/or vertices.
994 * Three edges. This is a linear matrix equation:
995 * each row of the matrix represents the fact that
996 * the point (x,y) we seek is at distance r from
997 * that edge, and we solve three of those
998 * simultaneously to obtain x,y,r. (We ignore r.)
1000 double matrix[9], vector[3], vector2[3];
1003 for (m = 0; m < 3; m++) {
1004 int x1 = edgedot1[m]->x, x2 = edgedot2[m]->x;
1005 int y1 = edgedot1[m]->y, y2 = edgedot2[m]->y;
1006 int dx = x2-x1, dy = y2-y1;
1009 * ((x,y) - (x1,y1)) . (dy,-dx) = r |(dx,dy)|
1011 * => x dy - y dx - r |(dx,dy)| = (x1 dy - y1 dx)
1014 matrix[3*m+1] = -dx;
1015 matrix[3*m+2] = -sqrt((double)dx*dx+(double)dy*dy);
1016 vector[m] = (double)x1*dy - (double)y1*dx;
1019 if (solve_3x3_matrix(matrix, vector, vector2)) {
1020 cx[cn] = vector2[0];
1021 cy[cn] = vector2[1];
1024 } else if (nedges == 2) {
1026 * Two edges and a dot. This will end up in a
1027 * quadratic equation.
1029 * First, look at the two edges. Having our point
1030 * be some distance r from both of them gives rise
1031 * to a pair of linear equations in x,y,r of the
1034 * (x-x1) dy - (y-y1) dx = r sqrt(dx^2+dy^2)
1036 * We eliminate r between those equations to give
1037 * us a single linear equation in x,y describing
1038 * the locus of points equidistant from both lines
1039 * - i.e. the angle bisector.
1041 * We then choose one of x,y to be a parameter t,
1042 * and derive linear formulae for x,y,r in terms
1043 * of t. This enables us to write down the
1044 * circular equation (x-xd)^2+(y-yd)^2=r^2 as a
1045 * quadratic in t; solving that and substituting
1046 * in for x,y gives us two candidate points.
1048 double eqs[2][4]; /* a,b,c,d : ax+by+cr=d */
1049 double eq[3]; /* a,b,c: ax+by=c */
1050 double xt[2], yt[2], rt[2]; /* a,b: {x,y,r}=at+b */
1051 double q[3]; /* a,b,c: at^2+bt+c=0 */
1054 /* Find equations of the two input lines. */
1055 for (m = 0; m < 2; m++) {
1056 int x1 = edgedot1[m]->x, x2 = edgedot2[m]->x;
1057 int y1 = edgedot1[m]->y, y2 = edgedot2[m]->y;
1058 int dx = x2-x1, dy = y2-y1;
1062 eqs[m][2] = -sqrt(dx*dx+dy*dy);
1063 eqs[m][3] = x1*dy - y1*dx;
1066 /* Derive the angle bisector by eliminating r. */
1067 eq[0] = eqs[0][0]*eqs[1][2] - eqs[1][0]*eqs[0][2];
1068 eq[1] = eqs[0][1]*eqs[1][2] - eqs[1][1]*eqs[0][2];
1069 eq[2] = eqs[0][3]*eqs[1][2] - eqs[1][3]*eqs[0][2];
1071 /* Parametrise x and y in terms of some t. */
1072 if (fabs(eq[0]) < fabs(eq[1])) {
1073 /* Parameter is x. */
1074 xt[0] = 1; xt[1] = 0;
1075 yt[0] = -eq[0]/eq[1]; yt[1] = eq[2]/eq[1];
1077 /* Parameter is y. */
1078 yt[0] = 1; yt[1] = 0;
1079 xt[0] = -eq[1]/eq[0]; xt[1] = eq[2]/eq[0];
1082 /* Find a linear representation of r using eqs[0]. */
1083 rt[0] = -(eqs[0][0]*xt[0] + eqs[0][1]*yt[0])/eqs[0][2];
1084 rt[1] = (eqs[0][3] - eqs[0][0]*xt[1] -
1085 eqs[0][1]*yt[1])/eqs[0][2];
1087 /* Construct the quadratic equation. */
1088 q[0] = -rt[0]*rt[0];
1089 q[1] = -2*rt[0]*rt[1];
1090 q[2] = -rt[1]*rt[1];
1091 q[0] += xt[0]*xt[0];
1092 q[1] += 2*xt[0]*(xt[1]-dots[0]->x);
1093 q[2] += (xt[1]-dots[0]->x)*(xt[1]-dots[0]->x);
1094 q[0] += yt[0]*yt[0];
1095 q[1] += 2*yt[0]*(yt[1]-dots[0]->y);
1096 q[2] += (yt[1]-dots[0]->y)*(yt[1]-dots[0]->y);
1099 disc = q[1]*q[1] - 4*q[0]*q[2];
1105 t = (-q[1] + disc) / (2*q[0]);
1106 cx[cn] = xt[0]*t + xt[1];
1107 cy[cn] = yt[0]*t + yt[1];
1110 t = (-q[1] - disc) / (2*q[0]);
1111 cx[cn] = xt[0]*t + xt[1];
1112 cy[cn] = yt[0]*t + yt[1];
1115 } else if (nedges == 1) {
1117 * Two dots and an edge. This one's another
1118 * quadratic equation.
1120 * The point we want must lie on the perpendicular
1121 * bisector of the two dots; that much is obvious.
1122 * So we can construct a parametrisation of that
1123 * bisecting line, giving linear formulae for x,y
1124 * in terms of t. We can also express the distance
1125 * from the edge as such a linear formula.
1127 * Then we set that equal to the radius of the
1128 * circle passing through the two points, which is
1129 * a Pythagoras exercise; that gives rise to a
1130 * quadratic in t, which we solve.
1132 double xt[2], yt[2], rt[2]; /* a,b: {x,y,r}=at+b */
1133 double q[3]; /* a,b,c: at^2+bt+c=0 */
1137 /* Find parametric formulae for x,y. */
1139 int x1 = dots[0]->x, x2 = dots[1]->x;
1140 int y1 = dots[0]->y, y2 = dots[1]->y;
1141 int dx = x2-x1, dy = y2-y1;
1142 double d = sqrt((double)dx*dx + (double)dy*dy);
1144 xt[1] = (x1+x2)/2.0;
1145 yt[1] = (y1+y2)/2.0;
1146 /* It's convenient if we have t at standard scale. */
1150 /* Also note down half the separation between
1151 * the dots, for use in computing the circle radius. */
1155 /* Find a parametric formula for r. */
1157 int x1 = edgedot1[0]->x, x2 = edgedot2[0]->x;
1158 int y1 = edgedot1[0]->y, y2 = edgedot2[0]->y;
1159 int dx = x2-x1, dy = y2-y1;
1160 double d = sqrt((double)dx*dx + (double)dy*dy);
1161 rt[0] = (xt[0]*dy - yt[0]*dx) / d;
1162 rt[1] = ((xt[1]-x1)*dy - (yt[1]-y1)*dx) / d;
1165 /* Construct the quadratic equation. */
1167 q[1] = 2*rt[0]*rt[1];
1170 q[2] -= halfsep*halfsep;
1173 disc = q[1]*q[1] - 4*q[0]*q[2];
1179 t = (-q[1] + disc) / (2*q[0]);
1180 cx[cn] = xt[0]*t + xt[1];
1181 cy[cn] = yt[0]*t + yt[1];
1184 t = (-q[1] - disc) / (2*q[0]);
1185 cx[cn] = xt[0]*t + xt[1];
1186 cy[cn] = yt[0]*t + yt[1];
1189 } else if (nedges == 0) {
1191 * Three dots. This is another linear matrix
1192 * equation, this time with each row of the matrix
1193 * representing the perpendicular bisector between
1194 * two of the points. Of course we only need two
1195 * such lines to find their intersection, so we
1196 * need only solve a 2x2 matrix equation.
1199 double matrix[4], vector[2], vector2[2];
1202 for (m = 0; m < 2; m++) {
1203 int x1 = dots[m]->x, x2 = dots[m+1]->x;
1204 int y1 = dots[m]->y, y2 = dots[m+1]->y;
1205 int dx = x2-x1, dy = y2-y1;
1208 * ((x,y) - (x1,y1)) . (dx,dy) = 1/2 |(dx,dy)|^2
1210 * => 2x dx + 2y dy = dx^2+dy^2 + (2 x1 dx + 2 y1 dy)
1212 matrix[2*m+0] = 2*dx;
1213 matrix[2*m+1] = 2*dy;
1214 vector[m] = ((double)dx*dx + (double)dy*dy +
1215 2.0*x1*dx + 2.0*y1*dy);
1218 if (solve_2x2_matrix(matrix, vector, vector2)) {
1219 cx[cn] = vector2[0];
1220 cy[cn] = vector2[1];
1226 * Now go through our candidate points and see if any
1227 * of them are better than what we've got so far.
1229 for (m = 0; m < cn; m++) {
1230 double x = cx[m], y = cy[m];
1233 * First, disqualify the point if it's not inside
1234 * the polygon, which we work out by counting the
1235 * edges to the right of the point. (For
1236 * tiebreaking purposes when edges start or end on
1237 * our y-coordinate or go right through it, we
1238 * consider our point to be offset by a small
1239 * _positive_ epsilon in both the x- and
1243 for (e = 0; e < f->order; e++) {
1244 int xs = f->edges[e]->dot1->x;
1245 int xe = f->edges[e]->dot2->x;
1246 int ys = f->edges[e]->dot1->y;
1247 int ye = f->edges[e]->dot2->y;
1248 if ((y >= ys && y < ye) || (y >= ye && y < ys)) {
1250 * The line goes past our y-position. Now we need
1251 * to know if its x-coordinate when it does so is
1254 * The x-coordinate in question is mathematically
1255 * (y - ys) * (xe - xs) / (ye - ys), and we want
1256 * to know whether (x - xs) >= that. Of course we
1257 * avoid the division, so we can work in integers;
1258 * to do this we must multiply both sides of the
1259 * inequality by ye - ys, which means we must
1260 * first check that's not negative.
1262 int num = xe - xs, denom = ye - ys;
1267 if ((x - xs) * denom >= (y - ys) * num)
1274 double mindist = HUGE_VAL;
1277 double mindist = DBL_MAX;
1279 #error No way to get maximum floating-point number.
1285 * This point is inside the polygon, so now we check
1286 * its minimum distance to every edge and corner.
1287 * First the corners ...
1289 for (d = 0; d < f->order; d++) {
1290 int xp = f->dots[d]->x;
1291 int yp = f->dots[d]->y;
1292 double dx = x - xp, dy = y - yp;
1293 double dist = dx*dx + dy*dy;
1299 * ... and now also check the perpendicular distance
1300 * to every edge, if the perpendicular lies between
1301 * the edge's endpoints.
1303 for (e = 0; e < f->order; e++) {
1304 int xs = f->edges[e]->dot1->x;
1305 int xe = f->edges[e]->dot2->x;
1306 int ys = f->edges[e]->dot1->y;
1307 int ye = f->edges[e]->dot2->y;
1310 * If s and e are our endpoints, and p our
1311 * candidate circle centre, the foot of a
1312 * perpendicular from p to the line se lies
1313 * between s and e if and only if (p-s).(e-s) lies
1314 * strictly between 0 and (e-s).(e-s).
1316 int edx = xe - xs, edy = ye - ys;
1317 double pdx = x - xs, pdy = y - ys;
1318 double pde = pdx * edx + pdy * edy;
1319 long ede = (long)edx * edx + (long)edy * edy;
1320 if (0 < pde && pde < ede) {
1322 * Yes, the nearest point on this edge is
1323 * closer than either endpoint, so we must
1324 * take it into account by measuring the
1325 * perpendicular distance to the edge and
1326 * checking its square against mindist.
1329 double pdre = pdx * edy - pdy * edx;
1330 double sqlen = pdre * pdre / ede;
1332 if (mindist > sqlen)
1338 * Right. Now we know the biggest circle around this
1339 * point, so we can check it against bestdist.
1341 if (bestdist < mindist) {
1365 assert(bestdist > 0);
1367 f->has_incentre = TRUE;
1368 f->ix = xbest + 0.5; /* round to nearest */
1369 f->iy = ybest + 0.5;
1372 /* ------ Generate various types of grid ------ */
1374 /* General method is to generate faces, by calculating their dot coordinates.
1375 * As new faces are added, we keep track of all the dots so we can tell when
1376 * a new face reuses an existing dot. For example, two squares touching at an
1377 * edge would generate six unique dots: four dots from the first face, then
1378 * two additional dots for the second face, because we detect the other two
1379 * dots have already been taken up. This list is stored in a tree234
1380 * called "points". No extra memory-allocation needed here - we store the
1381 * actual grid_dot* pointers, which all point into the g->dots list.
1382 * For this reason, we have to calculate coordinates in such a way as to
1383 * eliminate any rounding errors, so we can detect when a dot on one
1384 * face precisely lands on a dot of a different face. No floating-point
1388 #define SQUARE_TILESIZE 20
1390 static void grid_size_square(int width, int height,
1391 int *tilesize, int *xextent, int *yextent)
1393 int a = SQUARE_TILESIZE;
1396 *xextent = width * a;
1397 *yextent = height * a;
1400 static grid *grid_new_square(int width, int height, const char *desc)
1404 int a = SQUARE_TILESIZE;
1406 /* Upper bounds - don't have to be exact */
1407 int max_faces = width * height;
1408 int max_dots = (width + 1) * (height + 1);
1412 grid *g = grid_empty();
1414 g->faces = snewn(max_faces, grid_face);
1415 g->dots = snewn(max_dots, grid_dot);
1417 points = newtree234(grid_point_cmp_fn);
1419 /* generate square faces */
1420 for (y = 0; y < height; y++) {
1421 for (x = 0; x < width; x++) {
1427 grid_face_add_new(g, 4);
1428 d = grid_get_dot(g, points, px, py);
1429 grid_face_set_dot(g, d, 0);
1430 d = grid_get_dot(g, points, px + a, py);
1431 grid_face_set_dot(g, d, 1);
1432 d = grid_get_dot(g, points, px + a, py + a);
1433 grid_face_set_dot(g, d, 2);
1434 d = grid_get_dot(g, points, px, py + a);
1435 grid_face_set_dot(g, d, 3);
1439 freetree234(points);
1440 assert(g->num_faces <= max_faces);
1441 assert(g->num_dots <= max_dots);
1443 grid_make_consistent(g);
1447 #define HONEY_TILESIZE 45
1448 /* Vector for side of hexagon - ratio is close to sqrt(3) */
1452 static void grid_size_honeycomb(int width, int height,
1453 int *tilesize, int *xextent, int *yextent)
1458 *tilesize = HONEY_TILESIZE;
1459 *xextent = (3 * a * (width-1)) + 4*a;
1460 *yextent = (2 * b * (height-1)) + 3*b;
1463 static grid *grid_new_honeycomb(int width, int height, const char *desc)
1469 /* Upper bounds - don't have to be exact */
1470 int max_faces = width * height;
1471 int max_dots = 2 * (width + 1) * (height + 1);
1475 grid *g = grid_empty();
1476 g->tilesize = HONEY_TILESIZE;
1477 g->faces = snewn(max_faces, grid_face);
1478 g->dots = snewn(max_dots, grid_dot);
1480 points = newtree234(grid_point_cmp_fn);
1482 /* generate hexagonal faces */
1483 for (y = 0; y < height; y++) {
1484 for (x = 0; x < width; x++) {
1491 grid_face_add_new(g, 6);
1493 d = grid_get_dot(g, points, cx - a, cy - b);
1494 grid_face_set_dot(g, d, 0);
1495 d = grid_get_dot(g, points, cx + a, cy - b);
1496 grid_face_set_dot(g, d, 1);
1497 d = grid_get_dot(g, points, cx + 2*a, cy);
1498 grid_face_set_dot(g, d, 2);
1499 d = grid_get_dot(g, points, cx + a, cy + b);
1500 grid_face_set_dot(g, d, 3);
1501 d = grid_get_dot(g, points, cx - a, cy + b);
1502 grid_face_set_dot(g, d, 4);
1503 d = grid_get_dot(g, points, cx - 2*a, cy);
1504 grid_face_set_dot(g, d, 5);
1508 freetree234(points);
1509 assert(g->num_faces <= max_faces);
1510 assert(g->num_dots <= max_dots);
1512 grid_make_consistent(g);
1516 #define TRIANGLE_TILESIZE 18
1517 /* Vector for side of triangle - ratio is close to sqrt(3) */
1518 #define TRIANGLE_VEC_X 15
1519 #define TRIANGLE_VEC_Y 26
1521 static void grid_size_triangular(int width, int height,
1522 int *tilesize, int *xextent, int *yextent)
1524 int vec_x = TRIANGLE_VEC_X;
1525 int vec_y = TRIANGLE_VEC_Y;
1527 *tilesize = TRIANGLE_TILESIZE;
1528 *xextent = (width+1) * 2 * vec_x;
1529 *yextent = height * vec_y;
1532 static char *grid_validate_desc_triangular(grid_type type, int width,
1533 int height, const char *desc)
1536 * Triangular grids: an absent description is valid (indicating
1537 * the old-style approach which had 'ears', i.e. triangles
1538 * connected to only one other face, at some grid corners), and so
1539 * is a description reading just "0" (indicating the new-style
1540 * approach in which those ears are trimmed off). Anything else is
1544 if (!desc || !strcmp(desc, "0"))
1547 return "Unrecognised grid description.";
1550 /* Doesn't use the previous method of generation, it pre-dates it!
1551 * A triangular grid is just about simple enough to do by "brute force" */
1552 static grid *grid_new_triangular(int width, int height, const char *desc)
1555 int version = (desc == NULL ? -1 : atoi(desc));
1557 /* Vector for side of triangle - ratio is close to sqrt(3) */
1558 int vec_x = TRIANGLE_VEC_X;
1559 int vec_y = TRIANGLE_VEC_Y;
1563 /* convenient alias */
1566 grid *g = grid_empty();
1567 g->tilesize = TRIANGLE_TILESIZE;
1569 if (version == -1) {
1571 * Old-style triangular grid generation, preserved as-is for
1572 * backwards compatibility with old game ids, in which it's
1573 * just a little asymmetric and there are 'ears' (faces linked
1574 * to only one other face) at two grid corners.
1576 * Example old-style game ids, which should still work, and in
1577 * which you should see the ears in the TL/BR corners on the
1578 * first one and in the TL/BL corners on the second:
1580 * 5x5t1:2c2a1a2a201a1a1a1112a1a2b1211f0b21a2a2a0a
1581 * 5x6t1:a022a212h1a1d1a12c2b11a012b1a20d1a0a12e
1584 g->num_faces = width * height * 2;
1585 g->num_dots = (width + 1) * (height + 1);
1586 g->faces = snewn(g->num_faces, grid_face);
1587 g->dots = snewn(g->num_dots, grid_dot);
1591 for (y = 0; y <= height; y++) {
1592 for (x = 0; x <= width; x++) {
1593 grid_dot *d = g->dots + index;
1594 /* odd rows are offset to the right */
1598 d->x = x * 2 * vec_x + ((y % 2) ? vec_x : 0);
1604 /* generate faces */
1606 for (y = 0; y < height; y++) {
1607 for (x = 0; x < width; x++) {
1608 /* initialise two faces for this (x,y) */
1609 grid_face *f1 = g->faces + index;
1610 grid_face *f2 = f1 + 1;
1613 f1->dots = snewn(f1->order, grid_dot*);
1614 f1->has_incentre = FALSE;
1617 f2->dots = snewn(f2->order, grid_dot*);
1618 f2->has_incentre = FALSE;
1620 /* face descriptions depend on whether the row-number is
1623 f1->dots[0] = g->dots + y * w + x;
1624 f1->dots[1] = g->dots + (y + 1) * w + x + 1;
1625 f1->dots[2] = g->dots + (y + 1) * w + x;
1626 f2->dots[0] = g->dots + y * w + x;
1627 f2->dots[1] = g->dots + y * w + x + 1;
1628 f2->dots[2] = g->dots + (y + 1) * w + x + 1;
1630 f1->dots[0] = g->dots + y * w + x;
1631 f1->dots[1] = g->dots + y * w + x + 1;
1632 f1->dots[2] = g->dots + (y + 1) * w + x;
1633 f2->dots[0] = g->dots + y * w + x + 1;
1634 f2->dots[1] = g->dots + (y + 1) * w + x + 1;
1635 f2->dots[2] = g->dots + (y + 1) * w + x;
1642 * New-style approach, in which there are never any 'ears',
1643 * and if height is even then the grid is nicely 4-way
1646 * Example new-style grids:
1648 * 5x5t1:0_21120b11a1a01a1a00c1a0b211021c1h1a2a1a0a
1649 * 5x6t1:0_a1212c22c2a02a2f22a0c12a110d0e1c0c0a101121a1
1651 tree234 *points = newtree234(grid_point_cmp_fn);
1652 /* Upper bounds - don't have to be exact */
1653 int max_faces = height * (2*width+1);
1654 int max_dots = (height+1) * (width+1) * 4;
1656 g->faces = snewn(max_faces, grid_face);
1657 g->dots = snewn(max_dots, grid_dot);
1659 for (y = 0; y < height; y++) {
1661 * Each row contains (width+1) triangles one way up, and
1662 * (width) triangles the other way up. Which way up is
1663 * which varies with parity of y. Also, the dots around
1664 * each face will flip direction with parity of y, so we
1665 * set up n1 and n2 to cope with that easily.
1668 y0 = y1 = y * vec_y;
1677 for (x = 0; x <= width; x++) {
1678 int x0 = 2*x * vec_x, x1 = x0 + vec_x, x2 = x1 + vec_x;
1681 * If the grid has odd height, then we skip the first
1682 * and last triangles on this row, otherwise they'll
1685 if (height % 2 == 1 && y == height-1 && (x == 0 || x == width))
1688 grid_face_add_new(g, 3);
1689 grid_face_set_dot(g, grid_get_dot(g, points, x0, y0), 0);
1690 grid_face_set_dot(g, grid_get_dot(g, points, x1, y1), n1);
1691 grid_face_set_dot(g, grid_get_dot(g, points, x2, y0), n2);
1694 for (x = 0; x < width; x++) {
1695 int x0 = (2*x+1) * vec_x, x1 = x0 + vec_x, x2 = x1 + vec_x;
1697 grid_face_add_new(g, 3);
1698 grid_face_set_dot(g, grid_get_dot(g, points, x0, y1), 0);
1699 grid_face_set_dot(g, grid_get_dot(g, points, x1, y0), n2);
1700 grid_face_set_dot(g, grid_get_dot(g, points, x2, y1), n1);
1704 freetree234(points);
1705 assert(g->num_faces <= max_faces);
1706 assert(g->num_dots <= max_dots);
1709 grid_make_consistent(g);
1713 #define SNUBSQUARE_TILESIZE 18
1714 /* Vector for side of triangle - ratio is close to sqrt(3) */
1715 #define SNUBSQUARE_A 15
1716 #define SNUBSQUARE_B 26
1718 static void grid_size_snubsquare(int width, int height,
1719 int *tilesize, int *xextent, int *yextent)
1721 int a = SNUBSQUARE_A;
1722 int b = SNUBSQUARE_B;
1724 *tilesize = SNUBSQUARE_TILESIZE;
1725 *xextent = (a+b) * (width-1) + a + b;
1726 *yextent = (a+b) * (height-1) + a + b;
1729 static grid *grid_new_snubsquare(int width, int height, const char *desc)
1732 int a = SNUBSQUARE_A;
1733 int b = SNUBSQUARE_B;
1735 /* Upper bounds - don't have to be exact */
1736 int max_faces = 3 * width * height;
1737 int max_dots = 2 * (width + 1) * (height + 1);
1741 grid *g = grid_empty();
1742 g->tilesize = SNUBSQUARE_TILESIZE;
1743 g->faces = snewn(max_faces, grid_face);
1744 g->dots = snewn(max_dots, grid_dot);
1746 points = newtree234(grid_point_cmp_fn);
1748 for (y = 0; y < height; y++) {
1749 for (x = 0; x < width; x++) {
1752 int px = (a + b) * x;
1753 int py = (a + b) * y;
1755 /* generate square faces */
1756 grid_face_add_new(g, 4);
1758 d = grid_get_dot(g, points, px + a, py);
1759 grid_face_set_dot(g, d, 0);
1760 d = grid_get_dot(g, points, px + a + b, py + a);
1761 grid_face_set_dot(g, d, 1);
1762 d = grid_get_dot(g, points, px + b, py + a + b);
1763 grid_face_set_dot(g, d, 2);
1764 d = grid_get_dot(g, points, px, py + b);
1765 grid_face_set_dot(g, d, 3);
1767 d = grid_get_dot(g, points, px + b, py);
1768 grid_face_set_dot(g, d, 0);
1769 d = grid_get_dot(g, points, px + a + b, py + b);
1770 grid_face_set_dot(g, d, 1);
1771 d = grid_get_dot(g, points, px + a, py + a + b);
1772 grid_face_set_dot(g, d, 2);
1773 d = grid_get_dot(g, points, px, py + a);
1774 grid_face_set_dot(g, d, 3);
1777 /* generate up/down triangles */
1779 grid_face_add_new(g, 3);
1781 d = grid_get_dot(g, points, px + a, py);
1782 grid_face_set_dot(g, d, 0);
1783 d = grid_get_dot(g, points, px, py + b);
1784 grid_face_set_dot(g, d, 1);
1785 d = grid_get_dot(g, points, px - a, py);
1786 grid_face_set_dot(g, d, 2);
1788 d = grid_get_dot(g, points, px, py + a);
1789 grid_face_set_dot(g, d, 0);
1790 d = grid_get_dot(g, points, px + a, py + a + b);
1791 grid_face_set_dot(g, d, 1);
1792 d = grid_get_dot(g, points, px - a, py + a + b);
1793 grid_face_set_dot(g, d, 2);
1797 /* generate left/right triangles */
1799 grid_face_add_new(g, 3);
1801 d = grid_get_dot(g, points, px + a, py);
1802 grid_face_set_dot(g, d, 0);
1803 d = grid_get_dot(g, points, px + a + b, py - a);
1804 grid_face_set_dot(g, d, 1);
1805 d = grid_get_dot(g, points, px + a + b, py + a);
1806 grid_face_set_dot(g, d, 2);
1808 d = grid_get_dot(g, points, px, py - a);
1809 grid_face_set_dot(g, d, 0);
1810 d = grid_get_dot(g, points, px + b, py);
1811 grid_face_set_dot(g, d, 1);
1812 d = grid_get_dot(g, points, px, py + a);
1813 grid_face_set_dot(g, d, 2);
1819 freetree234(points);
1820 assert(g->num_faces <= max_faces);
1821 assert(g->num_dots <= max_dots);
1823 grid_make_consistent(g);
1827 #define CAIRO_TILESIZE 40
1828 /* Vector for side of pentagon - ratio is close to (4+sqrt(7))/3 */
1832 static void grid_size_cairo(int width, int height,
1833 int *tilesize, int *xextent, int *yextent)
1835 int b = CAIRO_B; /* a unused in determining grid size. */
1837 *tilesize = CAIRO_TILESIZE;
1838 *xextent = 2*b*(width-1) + 2*b;
1839 *yextent = 2*b*(height-1) + 2*b;
1842 static grid *grid_new_cairo(int width, int height, const char *desc)
1848 /* Upper bounds - don't have to be exact */
1849 int max_faces = 2 * width * height;
1850 int max_dots = 3 * (width + 1) * (height + 1);
1854 grid *g = grid_empty();
1855 g->tilesize = CAIRO_TILESIZE;
1856 g->faces = snewn(max_faces, grid_face);
1857 g->dots = snewn(max_dots, grid_dot);
1859 points = newtree234(grid_point_cmp_fn);
1861 for (y = 0; y < height; y++) {
1862 for (x = 0; x < width; x++) {
1868 /* horizontal pentagons */
1870 grid_face_add_new(g, 5);
1872 d = grid_get_dot(g, points, px + a, py - b);
1873 grid_face_set_dot(g, d, 0);
1874 d = grid_get_dot(g, points, px + 2*b - a, py - b);
1875 grid_face_set_dot(g, d, 1);
1876 d = grid_get_dot(g, points, px + 2*b, py);
1877 grid_face_set_dot(g, d, 2);
1878 d = grid_get_dot(g, points, px + b, py + a);
1879 grid_face_set_dot(g, d, 3);
1880 d = grid_get_dot(g, points, px, py);
1881 grid_face_set_dot(g, d, 4);
1883 d = grid_get_dot(g, points, px, py);
1884 grid_face_set_dot(g, d, 0);
1885 d = grid_get_dot(g, points, px + b, py - a);
1886 grid_face_set_dot(g, d, 1);
1887 d = grid_get_dot(g, points, px + 2*b, py);
1888 grid_face_set_dot(g, d, 2);
1889 d = grid_get_dot(g, points, px + 2*b - a, py + b);
1890 grid_face_set_dot(g, d, 3);
1891 d = grid_get_dot(g, points, px + a, py + b);
1892 grid_face_set_dot(g, d, 4);
1895 /* vertical pentagons */
1897 grid_face_add_new(g, 5);
1899 d = grid_get_dot(g, points, px, py);
1900 grid_face_set_dot(g, d, 0);
1901 d = grid_get_dot(g, points, px + b, py + a);
1902 grid_face_set_dot(g, d, 1);
1903 d = grid_get_dot(g, points, px + b, py + 2*b - a);
1904 grid_face_set_dot(g, d, 2);
1905 d = grid_get_dot(g, points, px, py + 2*b);
1906 grid_face_set_dot(g, d, 3);
1907 d = grid_get_dot(g, points, px - a, py + b);
1908 grid_face_set_dot(g, d, 4);
1910 d = grid_get_dot(g, points, px, py);
1911 grid_face_set_dot(g, d, 0);
1912 d = grid_get_dot(g, points, px + a, py + b);
1913 grid_face_set_dot(g, d, 1);
1914 d = grid_get_dot(g, points, px, py + 2*b);
1915 grid_face_set_dot(g, d, 2);
1916 d = grid_get_dot(g, points, px - b, py + 2*b - a);
1917 grid_face_set_dot(g, d, 3);
1918 d = grid_get_dot(g, points, px - b, py + a);
1919 grid_face_set_dot(g, d, 4);
1925 freetree234(points);
1926 assert(g->num_faces <= max_faces);
1927 assert(g->num_dots <= max_dots);
1929 grid_make_consistent(g);
1933 #define GREATHEX_TILESIZE 18
1934 /* Vector for side of triangle - ratio is close to sqrt(3) */
1935 #define GREATHEX_A 15
1936 #define GREATHEX_B 26
1938 static void grid_size_greathexagonal(int width, int height,
1939 int *tilesize, int *xextent, int *yextent)
1944 *tilesize = GREATHEX_TILESIZE;
1945 *xextent = (3*a + b) * (width-1) + 4*a;
1946 *yextent = (2*a + 2*b) * (height-1) + 3*b + a;
1949 static grid *grid_new_greathexagonal(int width, int height, const char *desc)
1955 /* Upper bounds - don't have to be exact */
1956 int max_faces = 6 * (width + 1) * (height + 1);
1957 int max_dots = 6 * width * height;
1961 grid *g = grid_empty();
1962 g->tilesize = GREATHEX_TILESIZE;
1963 g->faces = snewn(max_faces, grid_face);
1964 g->dots = snewn(max_dots, grid_dot);
1966 points = newtree234(grid_point_cmp_fn);
1968 for (y = 0; y < height; y++) {
1969 for (x = 0; x < width; x++) {
1971 /* centre of hexagon */
1972 int px = (3*a + b) * x;
1973 int py = (2*a + 2*b) * y;
1978 grid_face_add_new(g, 6);
1979 d = grid_get_dot(g, points, px - a, py - b);
1980 grid_face_set_dot(g, d, 0);
1981 d = grid_get_dot(g, points, px + a, py - b);
1982 grid_face_set_dot(g, d, 1);
1983 d = grid_get_dot(g, points, px + 2*a, py);
1984 grid_face_set_dot(g, d, 2);
1985 d = grid_get_dot(g, points, px + a, py + b);
1986 grid_face_set_dot(g, d, 3);
1987 d = grid_get_dot(g, points, px - a, py + b);
1988 grid_face_set_dot(g, d, 4);
1989 d = grid_get_dot(g, points, px - 2*a, py);
1990 grid_face_set_dot(g, d, 5);
1992 /* square below hexagon */
1993 if (y < height - 1) {
1994 grid_face_add_new(g, 4);
1995 d = grid_get_dot(g, points, px - a, py + b);
1996 grid_face_set_dot(g, d, 0);
1997 d = grid_get_dot(g, points, px + a, py + b);
1998 grid_face_set_dot(g, d, 1);
1999 d = grid_get_dot(g, points, px + a, py + 2*a + b);
2000 grid_face_set_dot(g, d, 2);
2001 d = grid_get_dot(g, points, px - a, py + 2*a + b);
2002 grid_face_set_dot(g, d, 3);
2005 /* square below right */
2006 if ((x < width - 1) && (((x % 2) == 0) || (y < height - 1))) {
2007 grid_face_add_new(g, 4);
2008 d = grid_get_dot(g, points, px + 2*a, py);
2009 grid_face_set_dot(g, d, 0);
2010 d = grid_get_dot(g, points, px + 2*a + b, py + a);
2011 grid_face_set_dot(g, d, 1);
2012 d = grid_get_dot(g, points, px + a + b, py + a + b);
2013 grid_face_set_dot(g, d, 2);
2014 d = grid_get_dot(g, points, px + a, py + b);
2015 grid_face_set_dot(g, d, 3);
2018 /* square below left */
2019 if ((x > 0) && (((x % 2) == 0) || (y < height - 1))) {
2020 grid_face_add_new(g, 4);
2021 d = grid_get_dot(g, points, px - 2*a, py);
2022 grid_face_set_dot(g, d, 0);
2023 d = grid_get_dot(g, points, px - a, py + b);
2024 grid_face_set_dot(g, d, 1);
2025 d = grid_get_dot(g, points, px - a - b, py + a + b);
2026 grid_face_set_dot(g, d, 2);
2027 d = grid_get_dot(g, points, px - 2*a - b, py + a);
2028 grid_face_set_dot(g, d, 3);
2031 /* Triangle below right */
2032 if ((x < width - 1) && (y < height - 1)) {
2033 grid_face_add_new(g, 3);
2034 d = grid_get_dot(g, points, px + a, py + b);
2035 grid_face_set_dot(g, d, 0);
2036 d = grid_get_dot(g, points, px + a + b, py + a + b);
2037 grid_face_set_dot(g, d, 1);
2038 d = grid_get_dot(g, points, px + a, py + 2*a + b);
2039 grid_face_set_dot(g, d, 2);
2042 /* Triangle below left */
2043 if ((x > 0) && (y < height - 1)) {
2044 grid_face_add_new(g, 3);
2045 d = grid_get_dot(g, points, px - a, py + b);
2046 grid_face_set_dot(g, d, 0);
2047 d = grid_get_dot(g, points, px - a, py + 2*a + b);
2048 grid_face_set_dot(g, d, 1);
2049 d = grid_get_dot(g, points, px - a - b, py + a + b);
2050 grid_face_set_dot(g, d, 2);
2055 freetree234(points);
2056 assert(g->num_faces <= max_faces);
2057 assert(g->num_dots <= max_dots);
2059 grid_make_consistent(g);
2063 #define OCTAGONAL_TILESIZE 40
2064 /* b/a approx sqrt(2) */
2065 #define OCTAGONAL_A 29
2066 #define OCTAGONAL_B 41
2068 static void grid_size_octagonal(int width, int height,
2069 int *tilesize, int *xextent, int *yextent)
2071 int a = OCTAGONAL_A;
2072 int b = OCTAGONAL_B;
2074 *tilesize = OCTAGONAL_TILESIZE;
2075 *xextent = (2*a + b) * width;
2076 *yextent = (2*a + b) * height;
2079 static grid *grid_new_octagonal(int width, int height, const char *desc)
2082 int a = OCTAGONAL_A;
2083 int b = OCTAGONAL_B;
2085 /* Upper bounds - don't have to be exact */
2086 int max_faces = 2 * width * height;
2087 int max_dots = 4 * (width + 1) * (height + 1);
2091 grid *g = grid_empty();
2092 g->tilesize = OCTAGONAL_TILESIZE;
2093 g->faces = snewn(max_faces, grid_face);
2094 g->dots = snewn(max_dots, grid_dot);
2096 points = newtree234(grid_point_cmp_fn);
2098 for (y = 0; y < height; y++) {
2099 for (x = 0; x < width; x++) {
2102 int px = (2*a + b) * x;
2103 int py = (2*a + b) * y;
2105 grid_face_add_new(g, 8);
2106 d = grid_get_dot(g, points, px + a, py);
2107 grid_face_set_dot(g, d, 0);
2108 d = grid_get_dot(g, points, px + a + b, py);
2109 grid_face_set_dot(g, d, 1);
2110 d = grid_get_dot(g, points, px + 2*a + b, py + a);
2111 grid_face_set_dot(g, d, 2);
2112 d = grid_get_dot(g, points, px + 2*a + b, py + a + b);
2113 grid_face_set_dot(g, d, 3);
2114 d = grid_get_dot(g, points, px + a + b, py + 2*a + b);
2115 grid_face_set_dot(g, d, 4);
2116 d = grid_get_dot(g, points, px + a, py + 2*a + b);
2117 grid_face_set_dot(g, d, 5);
2118 d = grid_get_dot(g, points, px, py + a + b);
2119 grid_face_set_dot(g, d, 6);
2120 d = grid_get_dot(g, points, px, py + a);
2121 grid_face_set_dot(g, d, 7);
2124 if ((x > 0) && (y > 0)) {
2125 grid_face_add_new(g, 4);
2126 d = grid_get_dot(g, points, px, py - a);
2127 grid_face_set_dot(g, d, 0);
2128 d = grid_get_dot(g, points, px + a, py);
2129 grid_face_set_dot(g, d, 1);
2130 d = grid_get_dot(g, points, px, py + a);
2131 grid_face_set_dot(g, d, 2);
2132 d = grid_get_dot(g, points, px - a, py);
2133 grid_face_set_dot(g, d, 3);
2138 freetree234(points);
2139 assert(g->num_faces <= max_faces);
2140 assert(g->num_dots <= max_dots);
2142 grid_make_consistent(g);
2146 #define KITE_TILESIZE 40
2147 /* b/a approx sqrt(3) */
2151 static void grid_size_kites(int width, int height,
2152 int *tilesize, int *xextent, int *yextent)
2157 *tilesize = KITE_TILESIZE;
2158 *xextent = 4*b * width + 2*b;
2159 *yextent = 6*a * (height-1) + 8*a;
2162 static grid *grid_new_kites(int width, int height, const char *desc)
2168 /* Upper bounds - don't have to be exact */
2169 int max_faces = 6 * width * height;
2170 int max_dots = 6 * (width + 1) * (height + 1);
2174 grid *g = grid_empty();
2175 g->tilesize = KITE_TILESIZE;
2176 g->faces = snewn(max_faces, grid_face);
2177 g->dots = snewn(max_dots, grid_dot);
2179 points = newtree234(grid_point_cmp_fn);
2181 for (y = 0; y < height; y++) {
2182 for (x = 0; x < width; x++) {
2184 /* position of order-6 dot */
2190 /* kite pointing up-left */
2191 grid_face_add_new(g, 4);
2192 d = grid_get_dot(g, points, px, py);
2193 grid_face_set_dot(g, d, 0);
2194 d = grid_get_dot(g, points, px + 2*b, py);
2195 grid_face_set_dot(g, d, 1);
2196 d = grid_get_dot(g, points, px + 2*b, py + 2*a);
2197 grid_face_set_dot(g, d, 2);
2198 d = grid_get_dot(g, points, px + b, py + 3*a);
2199 grid_face_set_dot(g, d, 3);
2201 /* kite pointing up */
2202 grid_face_add_new(g, 4);
2203 d = grid_get_dot(g, points, px, py);
2204 grid_face_set_dot(g, d, 0);
2205 d = grid_get_dot(g, points, px + b, py + 3*a);
2206 grid_face_set_dot(g, d, 1);
2207 d = grid_get_dot(g, points, px, py + 4*a);
2208 grid_face_set_dot(g, d, 2);
2209 d = grid_get_dot(g, points, px - b, py + 3*a);
2210 grid_face_set_dot(g, d, 3);
2212 /* kite pointing up-right */
2213 grid_face_add_new(g, 4);
2214 d = grid_get_dot(g, points, px, py);
2215 grid_face_set_dot(g, d, 0);
2216 d = grid_get_dot(g, points, px - b, py + 3*a);
2217 grid_face_set_dot(g, d, 1);
2218 d = grid_get_dot(g, points, px - 2*b, py + 2*a);
2219 grid_face_set_dot(g, d, 2);
2220 d = grid_get_dot(g, points, px - 2*b, py);
2221 grid_face_set_dot(g, d, 3);
2223 /* kite pointing down-right */
2224 grid_face_add_new(g, 4);
2225 d = grid_get_dot(g, points, px, py);
2226 grid_face_set_dot(g, d, 0);
2227 d = grid_get_dot(g, points, px - 2*b, py);
2228 grid_face_set_dot(g, d, 1);
2229 d = grid_get_dot(g, points, px - 2*b, py - 2*a);
2230 grid_face_set_dot(g, d, 2);
2231 d = grid_get_dot(g, points, px - b, py - 3*a);
2232 grid_face_set_dot(g, d, 3);
2234 /* kite pointing down */
2235 grid_face_add_new(g, 4);
2236 d = grid_get_dot(g, points, px, py);
2237 grid_face_set_dot(g, d, 0);
2238 d = grid_get_dot(g, points, px - b, py - 3*a);
2239 grid_face_set_dot(g, d, 1);
2240 d = grid_get_dot(g, points, px, py - 4*a);
2241 grid_face_set_dot(g, d, 2);
2242 d = grid_get_dot(g, points, px + b, py - 3*a);
2243 grid_face_set_dot(g, d, 3);
2245 /* kite pointing down-left */
2246 grid_face_add_new(g, 4);
2247 d = grid_get_dot(g, points, px, py);
2248 grid_face_set_dot(g, d, 0);
2249 d = grid_get_dot(g, points, px + b, py - 3*a);
2250 grid_face_set_dot(g, d, 1);
2251 d = grid_get_dot(g, points, px + 2*b, py - 2*a);
2252 grid_face_set_dot(g, d, 2);
2253 d = grid_get_dot(g, points, px + 2*b, py);
2254 grid_face_set_dot(g, d, 3);
2258 freetree234(points);
2259 assert(g->num_faces <= max_faces);
2260 assert(g->num_dots <= max_dots);
2262 grid_make_consistent(g);
2266 #define FLORET_TILESIZE 150
2267 /* -py/px is close to tan(30 - atan(sqrt(3)/9))
2268 * using py=26 makes everything lean to the left, rather than right
2270 #define FLORET_PX 75
2271 #define FLORET_PY -26
2273 static void grid_size_floret(int width, int height,
2274 int *tilesize, int *xextent, int *yextent)
2276 int px = FLORET_PX, py = FLORET_PY; /* |( 75, -26)| = 79.43 */
2277 int qx = 4*px/5, qy = -py*2; /* |( 60, 52)| = 79.40 */
2279 /* rx unused in determining grid size. */
2281 *tilesize = FLORET_TILESIZE;
2282 *xextent = (6*px+3*qx)/2 * (width-1) + 4*qx + 2*px;
2283 *yextent = (5*qy-4*py) * (height-1) + 4*qy + 2*ry;
2286 static grid *grid_new_floret(int width, int height, const char *desc)
2289 /* Vectors for sides; weird numbers needed to keep puzzle aligned with window
2290 * -py/px is close to tan(30 - atan(sqrt(3)/9))
2291 * using py=26 makes everything lean to the left, rather than right
2293 int px = FLORET_PX, py = FLORET_PY; /* |( 75, -26)| = 79.43 */
2294 int qx = 4*px/5, qy = -py*2; /* |( 60, 52)| = 79.40 */
2295 int rx = qx-px, ry = qy-py; /* |(-15, 78)| = 79.38 */
2297 /* Upper bounds - don't have to be exact */
2298 int max_faces = 6 * width * height;
2299 int max_dots = 9 * (width + 1) * (height + 1);
2303 grid *g = grid_empty();
2304 g->tilesize = FLORET_TILESIZE;
2305 g->faces = snewn(max_faces, grid_face);
2306 g->dots = snewn(max_dots, grid_dot);
2308 points = newtree234(grid_point_cmp_fn);
2310 /* generate pentagonal faces */
2311 for (y = 0; y < height; y++) {
2312 for (x = 0; x < width; x++) {
2315 int cx = (6*px+3*qx)/2 * x;
2316 int cy = (4*py-5*qy) * y;
2318 cy -= (4*py-5*qy)/2;
2319 else if (y && y == height-1)
2320 continue; /* make better looking grids? try 3x3 for instance */
2322 grid_face_add_new(g, 5);
2323 d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
2324 d = grid_get_dot(g, points, cx+2*rx , cy+2*ry ); grid_face_set_dot(g, d, 1);
2325 d = grid_get_dot(g, points, cx+2*rx+qx, cy+2*ry+qy); grid_face_set_dot(g, d, 2);
2326 d = grid_get_dot(g, points, cx+2*qx+rx, cy+2*qy+ry); grid_face_set_dot(g, d, 3);
2327 d = grid_get_dot(g, points, cx+2*qx , cy+2*qy ); grid_face_set_dot(g, d, 4);
2329 grid_face_add_new(g, 5);
2330 d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
2331 d = grid_get_dot(g, points, cx+2*qx , cy+2*qy ); grid_face_set_dot(g, d, 1);
2332 d = grid_get_dot(g, points, cx+2*qx+px, cy+2*qy+py); grid_face_set_dot(g, d, 2);
2333 d = grid_get_dot(g, points, cx+2*px+qx, cy+2*py+qy); grid_face_set_dot(g, d, 3);
2334 d = grid_get_dot(g, points, cx+2*px , cy+2*py ); grid_face_set_dot(g, d, 4);
2336 grid_face_add_new(g, 5);
2337 d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
2338 d = grid_get_dot(g, points, cx+2*px , cy+2*py ); grid_face_set_dot(g, d, 1);
2339 d = grid_get_dot(g, points, cx+2*px-rx, cy+2*py-ry); grid_face_set_dot(g, d, 2);
2340 d = grid_get_dot(g, points, cx-2*rx+px, cy-2*ry+py); grid_face_set_dot(g, d, 3);
2341 d = grid_get_dot(g, points, cx-2*rx , cy-2*ry ); grid_face_set_dot(g, d, 4);
2343 grid_face_add_new(g, 5);
2344 d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
2345 d = grid_get_dot(g, points, cx-2*rx , cy-2*ry ); grid_face_set_dot(g, d, 1);
2346 d = grid_get_dot(g, points, cx-2*rx-qx, cy-2*ry-qy); grid_face_set_dot(g, d, 2);
2347 d = grid_get_dot(g, points, cx-2*qx-rx, cy-2*qy-ry); grid_face_set_dot(g, d, 3);
2348 d = grid_get_dot(g, points, cx-2*qx , cy-2*qy ); grid_face_set_dot(g, d, 4);
2350 grid_face_add_new(g, 5);
2351 d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
2352 d = grid_get_dot(g, points, cx-2*qx , cy-2*qy ); grid_face_set_dot(g, d, 1);
2353 d = grid_get_dot(g, points, cx-2*qx-px, cy-2*qy-py); grid_face_set_dot(g, d, 2);
2354 d = grid_get_dot(g, points, cx-2*px-qx, cy-2*py-qy); grid_face_set_dot(g, d, 3);
2355 d = grid_get_dot(g, points, cx-2*px , cy-2*py ); grid_face_set_dot(g, d, 4);
2357 grid_face_add_new(g, 5);
2358 d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
2359 d = grid_get_dot(g, points, cx-2*px , cy-2*py ); grid_face_set_dot(g, d, 1);
2360 d = grid_get_dot(g, points, cx-2*px+rx, cy-2*py+ry); grid_face_set_dot(g, d, 2);
2361 d = grid_get_dot(g, points, cx+2*rx-px, cy+2*ry-py); grid_face_set_dot(g, d, 3);
2362 d = grid_get_dot(g, points, cx+2*rx , cy+2*ry ); grid_face_set_dot(g, d, 4);
2366 freetree234(points);
2367 assert(g->num_faces <= max_faces);
2368 assert(g->num_dots <= max_dots);
2370 grid_make_consistent(g);
2374 /* DODEC_* are used for dodecagonal and great-dodecagonal grids. */
2375 #define DODEC_TILESIZE 26
2376 /* Vector for side of triangle - ratio is close to sqrt(3) */
2380 static void grid_size_dodecagonal(int width, int height,
2381 int *tilesize, int *xextent, int *yextent)
2386 *tilesize = DODEC_TILESIZE;
2387 *xextent = (4*a + 2*b) * (width-1) + 3*(2*a + b);
2388 *yextent = (3*a + 2*b) * (height-1) + 2*(2*a + b);
2391 static grid *grid_new_dodecagonal(int width, int height, const char *desc)
2397 /* Upper bounds - don't have to be exact */
2398 int max_faces = 3 * width * height;
2399 int max_dots = 14 * width * height;
2403 grid *g = grid_empty();
2404 g->tilesize = DODEC_TILESIZE;
2405 g->faces = snewn(max_faces, grid_face);
2406 g->dots = snewn(max_dots, grid_dot);
2408 points = newtree234(grid_point_cmp_fn);
2410 for (y = 0; y < height; y++) {
2411 for (x = 0; x < width; x++) {
2413 /* centre of dodecagon */
2414 int px = (4*a + 2*b) * x;
2415 int py = (3*a + 2*b) * y;
2420 grid_face_add_new(g, 12);
2421 d = grid_get_dot(g, points, px + ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 0);
2422 d = grid_get_dot(g, points, px + ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 1);
2423 d = grid_get_dot(g, points, px + (2*a + b), py - ( a )); grid_face_set_dot(g, d, 2);
2424 d = grid_get_dot(g, points, px + (2*a + b), py + ( a )); grid_face_set_dot(g, d, 3);
2425 d = grid_get_dot(g, points, px + ( a + b), py + ( a + b)); grid_face_set_dot(g, d, 4);
2426 d = grid_get_dot(g, points, px + ( a ), py + (2*a + b)); grid_face_set_dot(g, d, 5);
2427 d = grid_get_dot(g, points, px - ( a ), py + (2*a + b)); grid_face_set_dot(g, d, 6);
2428 d = grid_get_dot(g, points, px - ( a + b), py + ( a + b)); grid_face_set_dot(g, d, 7);
2429 d = grid_get_dot(g, points, px - (2*a + b), py + ( a )); grid_face_set_dot(g, d, 8);
2430 d = grid_get_dot(g, points, px - (2*a + b), py - ( a )); grid_face_set_dot(g, d, 9);
2431 d = grid_get_dot(g, points, px - ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 10);
2432 d = grid_get_dot(g, points, px - ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 11);
2434 /* triangle below dodecagon */
2435 if ((y < height - 1 && (x < width - 1 || !(y % 2)) && (x > 0 || (y % 2)))) {
2436 grid_face_add_new(g, 3);
2437 d = grid_get_dot(g, points, px + a, py + (2*a + b)); grid_face_set_dot(g, d, 0);
2438 d = grid_get_dot(g, points, px , py + (2*a + 2*b)); grid_face_set_dot(g, d, 1);
2439 d = grid_get_dot(g, points, px - a, py + (2*a + b)); grid_face_set_dot(g, d, 2);
2442 /* triangle above dodecagon */
2443 if ((y && (x < width - 1 || !(y % 2)) && (x > 0 || (y % 2)))) {
2444 grid_face_add_new(g, 3);
2445 d = grid_get_dot(g, points, px - a, py - (2*a + b)); grid_face_set_dot(g, d, 0);
2446 d = grid_get_dot(g, points, px , py - (2*a + 2*b)); grid_face_set_dot(g, d, 1);
2447 d = grid_get_dot(g, points, px + a, py - (2*a + b)); grid_face_set_dot(g, d, 2);
2452 freetree234(points);
2453 assert(g->num_faces <= max_faces);
2454 assert(g->num_dots <= max_dots);
2456 grid_make_consistent(g);
2460 static void grid_size_greatdodecagonal(int width, int height,
2461 int *tilesize, int *xextent, int *yextent)
2466 *tilesize = DODEC_TILESIZE;
2467 *xextent = (6*a + 2*b) * (width-1) + 2*(2*a + b) + 3*a + b;
2468 *yextent = (3*a + 3*b) * (height-1) + 2*(2*a + b);
2471 static grid *grid_new_greatdodecagonal(int width, int height, const char *desc)
2474 /* Vector for side of triangle - ratio is close to sqrt(3) */
2478 /* Upper bounds - don't have to be exact */
2479 int max_faces = 30 * width * height;
2480 int max_dots = 200 * width * height;
2484 grid *g = grid_empty();
2485 g->tilesize = DODEC_TILESIZE;
2486 g->faces = snewn(max_faces, grid_face);
2487 g->dots = snewn(max_dots, grid_dot);
2489 points = newtree234(grid_point_cmp_fn);
2491 for (y = 0; y < height; y++) {
2492 for (x = 0; x < width; x++) {
2494 /* centre of dodecagon */
2495 int px = (6*a + 2*b) * x;
2496 int py = (3*a + 3*b) * y;
2501 grid_face_add_new(g, 12);
2502 d = grid_get_dot(g, points, px + ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 0);
2503 d = grid_get_dot(g, points, px + ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 1);
2504 d = grid_get_dot(g, points, px + (2*a + b), py - ( a )); grid_face_set_dot(g, d, 2);
2505 d = grid_get_dot(g, points, px + (2*a + b), py + ( a )); grid_face_set_dot(g, d, 3);
2506 d = grid_get_dot(g, points, px + ( a + b), py + ( a + b)); grid_face_set_dot(g, d, 4);
2507 d = grid_get_dot(g, points, px + ( a ), py + (2*a + b)); grid_face_set_dot(g, d, 5);
2508 d = grid_get_dot(g, points, px - ( a ), py + (2*a + b)); grid_face_set_dot(g, d, 6);
2509 d = grid_get_dot(g, points, px - ( a + b), py + ( a + b)); grid_face_set_dot(g, d, 7);
2510 d = grid_get_dot(g, points, px - (2*a + b), py + ( a )); grid_face_set_dot(g, d, 8);
2511 d = grid_get_dot(g, points, px - (2*a + b), py - ( a )); grid_face_set_dot(g, d, 9);
2512 d = grid_get_dot(g, points, px - ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 10);
2513 d = grid_get_dot(g, points, px - ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 11);
2515 /* hexagon below dodecagon */
2516 if (y < height - 1 && (x < width - 1 || !(y % 2)) && (x > 0 || (y % 2))) {
2517 grid_face_add_new(g, 6);
2518 d = grid_get_dot(g, points, px + a, py + (2*a + b)); grid_face_set_dot(g, d, 0);
2519 d = grid_get_dot(g, points, px + 2*a, py + (2*a + 2*b)); grid_face_set_dot(g, d, 1);
2520 d = grid_get_dot(g, points, px + a, py + (2*a + 3*b)); grid_face_set_dot(g, d, 2);
2521 d = grid_get_dot(g, points, px - a, py + (2*a + 3*b)); grid_face_set_dot(g, d, 3);
2522 d = grid_get_dot(g, points, px - 2*a, py + (2*a + 2*b)); grid_face_set_dot(g, d, 4);
2523 d = grid_get_dot(g, points, px - a, py + (2*a + b)); grid_face_set_dot(g, d, 5);
2526 /* hexagon above dodecagon */
2527 if (y && (x < width - 1 || !(y % 2)) && (x > 0 || (y % 2))) {
2528 grid_face_add_new(g, 6);
2529 d = grid_get_dot(g, points, px - a, py - (2*a + b)); grid_face_set_dot(g, d, 0);
2530 d = grid_get_dot(g, points, px - 2*a, py - (2*a + 2*b)); grid_face_set_dot(g, d, 1);
2531 d = grid_get_dot(g, points, px - a, py - (2*a + 3*b)); grid_face_set_dot(g, d, 2);
2532 d = grid_get_dot(g, points, px + a, py - (2*a + 3*b)); grid_face_set_dot(g, d, 3);
2533 d = grid_get_dot(g, points, px + 2*a, py - (2*a + 2*b)); grid_face_set_dot(g, d, 4);
2534 d = grid_get_dot(g, points, px + a, py - (2*a + b)); grid_face_set_dot(g, d, 5);
2537 /* square on right of dodecagon */
2538 if (x < width - 1) {
2539 grid_face_add_new(g, 4);
2540 d = grid_get_dot(g, points, px + 2*a + b, py - a); grid_face_set_dot(g, d, 0);
2541 d = grid_get_dot(g, points, px + 4*a + b, py - a); grid_face_set_dot(g, d, 1);
2542 d = grid_get_dot(g, points, px + 4*a + b, py + a); grid_face_set_dot(g, d, 2);
2543 d = grid_get_dot(g, points, px + 2*a + b, py + a); grid_face_set_dot(g, d, 3);
2546 /* square on top right of dodecagon */
2547 if (y && (x < width - 1 || !(y % 2))) {
2548 grid_face_add_new(g, 4);
2549 d = grid_get_dot(g, points, px + ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 0);
2550 d = grid_get_dot(g, points, px + (2*a ), py - (2*a + 2*b)); grid_face_set_dot(g, d, 1);
2551 d = grid_get_dot(g, points, px + (2*a + b), py - ( a + 2*b)); grid_face_set_dot(g, d, 2);
2552 d = grid_get_dot(g, points, px + ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 3);
2555 /* square on top left of dodecagon */
2556 if (y && (x || (y % 2))) {
2557 grid_face_add_new(g, 4);
2558 d = grid_get_dot(g, points, px - ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 0);
2559 d = grid_get_dot(g, points, px - (2*a + b), py - ( a + 2*b)); grid_face_set_dot(g, d, 1);
2560 d = grid_get_dot(g, points, px - (2*a ), py - (2*a + 2*b)); grid_face_set_dot(g, d, 2);
2561 d = grid_get_dot(g, points, px - ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 3);
2566 freetree234(points);
2567 assert(g->num_faces <= max_faces);
2568 assert(g->num_dots <= max_dots);
2570 grid_make_consistent(g);
2574 typedef struct setface_ctx
2576 int xmin, xmax, ymin, ymax;
2582 static double round_int_nearest_away(double r)
2584 return (r > 0.0) ? floor(r + 0.5) : ceil(r - 0.5);
2587 static int set_faces(penrose_state *state, vector *vs, int n, int depth)
2589 setface_ctx *sf_ctx = (setface_ctx *)state->ctx;
2593 if (depth < state->max_depth) return 0;
2594 #ifdef DEBUG_PENROSE
2595 if (n != 4) return 0; /* triangles are sent as debugging. */
2598 for (i = 0; i < n; i++) {
2599 double tx = v_x(vs, i), ty = v_y(vs, i);
2601 xs[i] = (int)round_int_nearest_away(tx);
2602 ys[i] = (int)round_int_nearest_away(ty);
2604 if (xs[i] < sf_ctx->xmin || xs[i] > sf_ctx->xmax) return 0;
2605 if (ys[i] < sf_ctx->ymin || ys[i] > sf_ctx->ymax) return 0;
2608 grid_face_add_new(sf_ctx->g, n);
2609 debug(("penrose: new face l=%f gen=%d...",
2610 penrose_side_length(state->start_size, depth), depth));
2611 for (i = 0; i < n; i++) {
2612 grid_dot *d = grid_get_dot(sf_ctx->g, sf_ctx->points,
2614 grid_face_set_dot(sf_ctx->g, d, i);
2615 debug((" ... dot 0x%x (%d,%d) (was %2.2f,%2.2f)",
2616 d, d->x, d->y, v_x(vs, i), v_y(vs, i)));
2622 #define PENROSE_TILESIZE 100
2624 static void grid_size_penrose(int width, int height,
2625 int *tilesize, int *xextent, int *yextent)
2627 int l = PENROSE_TILESIZE;
2630 *xextent = l * width;
2631 *yextent = l * height;
2634 static grid *grid_new_penrose(int width, int height, int which, const char *desc); /* forward reference */
2636 static char *grid_new_desc_penrose(grid_type type, int width, int height, random_state *rs)
2638 int tilesize = PENROSE_TILESIZE, startsz, depth, xoff, yoff, aoff;
2639 double outer_radius;
2642 int which = (type == GRID_PENROSE_P2 ? PENROSE_P2 : PENROSE_P3);
2646 /* We want to produce a random bit of penrose tiling, so we
2647 * calculate a random offset (within the patch that penrose.c
2648 * calculates for us) and an angle (multiple of 36) to rotate
2651 penrose_calculate_size(which, tilesize, width, height,
2652 &outer_radius, &startsz, &depth);
2654 /* Calculate radius of (circumcircle of) patch, subtract from
2655 * radius calculated. */
2656 inner_radius = (int)(outer_radius - sqrt(width*width + height*height));
2658 /* Pick a random offset (the easy way: choose within outer
2659 * square, discarding while it's outside the circle) */
2661 xoff = random_upto(rs, 2*inner_radius) - inner_radius;
2662 yoff = random_upto(rs, 2*inner_radius) - inner_radius;
2663 } while (sqrt(xoff*xoff+yoff*yoff) > inner_radius);
2665 aoff = random_upto(rs, 360/36) * 36;
2667 debug(("grid_desc: ts %d, %dx%d patch, orad %2.2f irad %d",
2668 tilesize, width, height, outer_radius, inner_radius));
2669 debug((" -> xoff %d yoff %d aoff %d", xoff, yoff, aoff));
2671 sprintf(gd, "G%d,%d,%d", xoff, yoff, aoff);
2674 * Now test-generate our grid, to make sure it actually
2675 * produces something.
2677 g = grid_new_penrose(width, height, which, gd);
2682 /* If not, go back to the top of this while loop and try again
2683 * with a different random offset. */
2689 static char *grid_validate_desc_penrose(grid_type type, int width, int height,
2692 int tilesize = PENROSE_TILESIZE, startsz, depth, xoff, yoff, aoff, inner_radius;
2693 double outer_radius;
2694 int which = (type == GRID_PENROSE_P2 ? PENROSE_P2 : PENROSE_P3);
2698 return "Missing grid description string.";
2700 penrose_calculate_size(which, tilesize, width, height,
2701 &outer_radius, &startsz, &depth);
2702 inner_radius = (int)(outer_radius - sqrt(width*width + height*height));
2704 if (sscanf(desc, "G%d,%d,%d", &xoff, &yoff, &aoff) != 3)
2705 return "Invalid format grid description string.";
2707 if (sqrt(xoff*xoff + yoff*yoff) > inner_radius)
2708 return "Patch offset out of bounds.";
2709 if ((aoff % 36) != 0 || aoff < 0 || aoff >= 360)
2710 return "Angle offset out of bounds.";
2713 * Test-generate to ensure these parameters don't end us up with
2716 g = grid_new_penrose(width, height, which, desc);
2718 return "Patch coordinates do not identify a usable grid fragment";
2725 * We're asked for a grid of a particular size, and we generate enough
2726 * of the tiling so we can be sure to have enough random grid from which
2730 static grid *grid_new_penrose(int width, int height, int which, const char *desc)
2732 int max_faces, max_dots, tilesize = PENROSE_TILESIZE;
2733 int xsz, ysz, xoff, yoff, aoff;
2742 penrose_calculate_size(which, tilesize, width, height,
2743 &rradius, &ps.start_size, &ps.max_depth);
2745 debug(("penrose: w%d h%d, tile size %d, start size %d, depth %d",
2746 width, height, tilesize, ps.start_size, ps.max_depth));
2748 ps.new_tile = set_faces;
2751 max_faces = (width*3) * (height*3); /* somewhat paranoid... */
2752 max_dots = max_faces * 4; /* ditto... */
2755 g->tilesize = tilesize;
2756 g->faces = snewn(max_faces, grid_face);
2757 g->dots = snewn(max_dots, grid_dot);
2759 points = newtree234(grid_point_cmp_fn);
2761 memset(&sf_ctx, 0, sizeof(sf_ctx));
2763 sf_ctx.points = points;
2766 if (sscanf(desc, "G%d,%d,%d", &xoff, &yoff, &aoff) != 3)
2767 assert(!"Invalid grid description.");
2769 xoff = yoff = aoff = 0;
2772 xsz = width * tilesize;
2773 ysz = height * tilesize;
2775 sf_ctx.xmin = xoff - xsz/2;
2776 sf_ctx.xmax = xoff + xsz/2;
2777 sf_ctx.ymin = yoff - ysz/2;
2778 sf_ctx.ymax = yoff + ysz/2;
2780 debug(("penrose: centre (%f, %f) xsz %f ysz %f",
2781 0.0, 0.0, xsz, ysz));
2782 debug(("penrose: x range (%f --> %f), y range (%f --> %f)",
2783 sf_ctx.xmin, sf_ctx.xmax, sf_ctx.ymin, sf_ctx.ymax));
2785 penrose(&ps, which, aoff);
2787 freetree234(points);
2788 assert(g->num_faces <= max_faces);
2789 assert(g->num_dots <= max_dots);
2791 debug(("penrose: %d faces total (equivalent to %d wide by %d high)",
2792 g->num_faces, g->num_faces/height, g->num_faces/width));
2795 * Return NULL if we ended up with an empty grid, either because
2796 * the initial generation was over too small a rectangle to
2797 * encompass any face or because grid_trim_vigorously ended up
2798 * removing absolutely everything.
2800 if (g->num_faces == 0 || g->num_dots == 0) {
2804 grid_trim_vigorously(g);
2805 if (g->num_faces == 0 || g->num_dots == 0) {
2810 grid_make_consistent(g);
2813 * Centre the grid in its originally promised rectangle.
2815 g->lowest_x -= ((sf_ctx.xmax - sf_ctx.xmin) -
2816 (g->highest_x - g->lowest_x)) / 2;
2817 g->highest_x = g->lowest_x + (sf_ctx.xmax - sf_ctx.xmin);
2818 g->lowest_y -= ((sf_ctx.ymax - sf_ctx.ymin) -
2819 (g->highest_y - g->lowest_y)) / 2;
2820 g->highest_y = g->lowest_y + (sf_ctx.ymax - sf_ctx.ymin);
2825 static void grid_size_penrose_p2_kite(int width, int height,
2826 int *tilesize, int *xextent, int *yextent)
2828 grid_size_penrose(width, height, tilesize, xextent, yextent);
2831 static void grid_size_penrose_p3_thick(int width, int height,
2832 int *tilesize, int *xextent, int *yextent)
2834 grid_size_penrose(width, height, tilesize, xextent, yextent);
2837 static grid *grid_new_penrose_p2_kite(int width, int height, const char *desc)
2839 return grid_new_penrose(width, height, PENROSE_P2, desc);
2842 static grid *grid_new_penrose_p3_thick(int width, int height, const char *desc)
2844 return grid_new_penrose(width, height, PENROSE_P3, desc);
2847 /* ----------- End of grid generators ------------- */
2849 #define FNNEW(upper,lower) &grid_new_ ## lower,
2850 #define FNSZ(upper,lower) &grid_size_ ## lower,
2852 static grid *(*(grid_news[]))(int, int, const char*) = { GRIDGEN_LIST(FNNEW) };
2853 static void(*(grid_sizes[]))(int, int, int*, int*, int*) = { GRIDGEN_LIST(FNSZ) };
2855 char *grid_new_desc(grid_type type, int width, int height, random_state *rs)
2857 if (type == GRID_PENROSE_P2 || type == GRID_PENROSE_P3) {
2858 return grid_new_desc_penrose(type, width, height, rs);
2859 } else if (type == GRID_TRIANGULAR) {
2860 return dupstr("0"); /* up-to-date version of triangular grid */
2866 char *grid_validate_desc(grid_type type, int width, int height,
2869 if (type == GRID_PENROSE_P2 || type == GRID_PENROSE_P3) {
2870 return grid_validate_desc_penrose(type, width, height, desc);
2871 } else if (type == GRID_TRIANGULAR) {
2872 return grid_validate_desc_triangular(type, width, height, desc);
2875 return "Grid description strings not used with this grid type";
2880 grid *grid_new(grid_type type, int width, int height, const char *desc)
2882 char *err = grid_validate_desc(type, width, height, desc);
2883 if (err) assert(!"Invalid grid description.");
2885 return grid_news[type](width, height, desc);
2888 void grid_compute_size(grid_type type, int width, int height,
2889 int *tilesize, int *xextent, int *yextent)
2891 grid_sizes[type](width, height, tilesize, xextent, yextent);
2894 /* ----------- End of grid helpers ------------- */
2896 /* vim: set shiftwidth=4 tabstop=8: */