2023W-EFFPROG

Magic Hexagon.
Log | Files | Refs | README

magichex.c (12839B)


      1 #include <stdio.h>
      2 #include <stdlib.h>
      3 #include <assert.h>
      4 #include <string.h>
      5 #include <stdbool.h>
      6 
      7 #define PLACEHOLDER_ENTRY_ID -1 // Placeholder entries are not actually part of the hexagon but used for padding.
      8 
      9 /* constraint variable; if lo==hi, this is the variable's value */
     10 struct hexagonEntry {
     11   long id; /* variable id; id == PLACEHOLDER_ENTRY_ID if the variable is not part of the hexagon */
     12   long lo; /* lower bound */
     13   long hi; /* upper bound */
     14 };
     15 
     16 typedef struct hexagonEntry HexagonEntry;
     17 
     18 /* representation of a hexagon of order n: (2n-1)^2 square array
     19    for a hexagon of order 2:
     20      A B
     21     C D E
     22      F G
     23    the representation is:
     24     A B .
     25     C D E
     26     . F G
     27    The . slots have lo>hi.
     28    The . slots have id == PLACEHOLDER_ENTRY_ID
     29 
     30    The variable array is organized as a single-dimension array with accesses
     31    hexagon[y*r+x]
     32    This allows to access the diagonal with stride 2*order
     33 
     34    Variable names n, r, H, S according to German Wikipedia Article
     35    Instead of "i", the deviation variable is called "d" (d=0 means
     36    that the sum=0; to have the lowest value 1, d=2)
     37    
     38    n is the order (number of elements of a side of the hexagon).
     39    r = 2n-1 (length of the middle row/diagonals)
     40    H = 3n^2-3n+1 (number of variables)
     41    M = dH (sum of each row/diagonal)
     42    lowest  value = dr - (H-1)/2
     43    highest value = dr + (H-1)/2
     44 */
     45 static long n;
     46 static long number_hex_entries;
     47 static long d;
     48 static unsigned long r;
     49 static unsigned long H;
     50 static unsigned long M;
     51 
     52 unsigned long solutions = 0; /* counter of solutions */
     53 unsigned long leafs = 0; /* counter of leaf nodes visited in the search tree */
     54 
     55 long min(long a, long b)
     56 {
     57   return (a<b)?a:b;
     58 }
     59 
     60 long max(long a, long b)
     61 {
     62   return (a>b)?a:b;
     63 }
     64 
     65 
     66 /* unless otherwise noted, the following functions return
     67 
     68    0 if there is no solution (i.e., the action eliminates all values
     69    from a variable),
     70 
     71    1 if there was a change 
     72 
     73    2 if there was no change 
     74 */
     75 typedef enum RESULTTYPE { NOSOLUTION = 0, CHANGED = 1, NOCHANGE = 2 } CHANGE_IDENTIFIER;  
     76 
     77 
     78 CHANGE_IDENTIFIER sethi(HexagonEntry *hexagonEntry, long x) {
     79   assert(hexagonEntry->id != PLACEHOLDER_ENTRY_ID);
     80   if (x < hexagonEntry->hi) {
     81     hexagonEntry->hi = x;
     82     if (hexagonEntry->lo > hexagonEntry->hi)
     83       return NOSOLUTION;
     84     else
     85       return CHANGED;
     86   }
     87   return NOCHANGE;
     88 }
     89 
     90 CHANGE_IDENTIFIER setlo(HexagonEntry *hexagonEntry, long x) {
     91   assert(hexagonEntry->id != PLACEHOLDER_ENTRY_ID);
     92   if (x > hexagonEntry->lo) {
     93     hexagonEntry->lo = x;
     94     if (hexagonEntry->lo > hexagonEntry->hi)
     95       return NOSOLUTION;
     96     else
     97       return CHANGED;
     98   }
     99   return NOCHANGE;
    100 }
    101 
    102 /* returns NOSOLUTION if there is no solution, CHANGED if one of the variables has changed */
    103 CHANGE_IDENTIFIER lessthan(HexagonEntry *v1, HexagonEntry *v2)
    104 {
    105   assert(v1->id != PLACEHOLDER_ENTRY_ID);
    106   assert(v2->id != PLACEHOLDER_ENTRY_ID);
    107   CHANGE_IDENTIFIER f = sethi(v1, v2->hi-1);
    108   if (f == NOSOLUTION || f == CHANGED)
    109     return f;
    110   return (setlo(v2, v1->lo+1));
    111 }
    112 
    113 CHANGE_IDENTIFIER sum(HexagonEntry hexagon[], unsigned long nv, unsigned long stride,
    114         HexagonEntry *hexagonStart, HexagonEntry *hexagonEnd)
    115 {
    116   register unsigned long i;
    117   long hi = M;
    118   long lo = M;
    119   HexagonEntry *hexagonEntry_p;
    120 #if 0
    121   printf("sum(hexagonStart+%ld, %lu, %lu, %ld, hexagonStart, hexagonStart+%ld)   ",hexagon-hexagonStart,nv,stride,sum,hexagonEnd-hexagonStart); fflush(stdout);
    122   for (i=0, hexagonEntry_p=hexagon; i<nv; i++, hexagonEntry_p+=stride) {
    123     assert(hexagonEntry_p>=hexagonStart);
    124     assert(hexagonEntry_p<hexagonEnd);
    125     assert(hexagonEntry_p->id != PLACEHOLDER_ENTRY_ID);
    126     printf("hexagonEntry%ld ",hexagonEntry_p->id);
    127   }
    128   printf("\n");
    129 #endif
    130   for (i=0, hexagonEntry_p=hexagon; i<nv; i++, hexagonEntry_p+=stride) {
    131     assert(hexagonEntry_p>=hexagonStart);
    132     assert(hexagonEntry_p<hexagonEnd);
    133     assert(hexagonEntry_p->id != PLACEHOLDER_ENTRY_ID);
    134     hi -= hexagonEntry_p->lo;
    135     lo -= hexagonEntry_p->hi;
    136   }
    137 
    138   if(hi < lo){
    139     return NOSOLUTION;
    140   }
    141   /* hi is the upper bound of sum-sum(hexagon), lo the lower bound */
    142   for (i=0, hexagonEntry_p=hexagon; i<nv; i++, hexagonEntry_p+=stride) {
    143     CHANGE_IDENTIFIER f = sethi(hexagonEntry_p,hi+hexagonEntry_p->lo); /* readd hexagonEntry_p->lo to get an upper bound of hexagonEntry_p */
    144     assert(hexagonEntry_p>=hexagonStart);
    145     assert(hexagonEntry_p<hexagonEnd);
    146     assert(hexagonEntry_p->id != PLACEHOLDER_ENTRY_ID);
    147     if (f == NOSOLUTION || f == CHANGED)
    148       return f;
    149     f = setlo(hexagonEntry_p,lo+hexagonEntry_p->hi); /* likewise, readd hexagonEntry_p->hi */
    150     if (f == NOSOLUTION || f == CHANGED)
    151       return f;
    152   }
    153   return NOCHANGE;
    154 }
    155     
    156 /* reduce the ranges of the variables as much as possible (with the
    157    constraints we use);  returns 1 if all variables still have a
    158    non-empty range left, 0 if one has an empty range */
    159 bool solve(HexagonEntry hexagon[])
    160 {
    161   long o = d*r - (H-1)/2; /* offset in occupation array */
    162   unsigned long occupation[H]; /* if hexagon[i] has value x, occupation[x-o]==i, 
    163                                   if no hexagon[*] has value x, occupation[x-o]==H */
    164   unsigned long corners[] = {0, n-1, (n-1)*r+0, (n-1)*r+r-1, (r-1)*r+n-1, (r-1)*r+r-1};
    165   register unsigned long i;
    166   //printf("(re)start\n");
    167   /* deal with the alldifferent constraint */
    168   for (i=0; i<H; i++)
    169     occupation[i] = number_hex_entries;
    170  restart:
    171   for (i=0; i<number_hex_entries; i++) {
    172     HexagonEntry *hexagonEntry = &hexagon[i];
    173     if (hexagonEntry->lo == hexagonEntry->hi && occupation[hexagonEntry->lo-o] != i) {
    174       if (occupation[hexagonEntry->lo-o] < number_hex_entries)
    175         return false; /* another variable has the same value */
    176       occupation[hexagonEntry->lo-o] = i; /* occupy hexagonEntry->lo */
    177       //goto restart;
    178     }
    179   }
    180   bool k = false;
    181   /* now propagate the alldifferent results to the bounds */
    182   for (i=0; i<number_hex_entries; i++) {
    183     HexagonEntry *hexagonEntry = &hexagon[i];
    184     if (hexagonEntry->lo < hexagonEntry->hi) {
    185       if (occupation[hexagonEntry->lo-o] < number_hex_entries) {
    186         hexagonEntry->lo++;
    187         k = true;
    188       }
    189       if (occupation[hexagonEntry->hi-o] < number_hex_entries) {
    190         hexagonEntry->hi--;
    191         k = true;
    192       }
    193     }
    194   }
    195   if(k){
    196     goto restart;
    197   }
    198   /* the < constraints; all other corners are smaller than the first
    199      one (eliminate rotational symmetry) */
    200   for (i=1; i<sizeof(corners)/sizeof(corners[0]); i++) {
    201     CHANGE_IDENTIFIER f = lessthan(&hexagon[corners[0]],&hexagon[corners[i]]);
    202     if (f==NOSOLUTION) return false;
    203     if (f==CHANGED) goto restart;
    204   }
    205   /* eliminate the mirror symmetry between the corners to the right
    206      and left of the first corner */
    207   {
    208     CHANGE_IDENTIFIER f = lessthan(&hexagon[corners[2]],&hexagon[corners[1]]); 
    209     if (f==NOSOLUTION) return false;
    210     if (f==CHANGED) goto restart;
    211   }
    212   /* sum constraints: each line and diagonal sums up to M */
    213   /* line sum constraints */
    214   for (i=0; i<r; i++) {
    215     CHANGE_IDENTIFIER f;
    216     /* line */
    217     f = sum(hexagon+r*i+max(0,i+1-n), min(i+n,r+n-i-1), 1, hexagon, hexagon+number_hex_entries);
    218     if (f==NOSOLUTION) return false;
    219     if (f==CHANGED) k = true;
    220     /* column (diagonal down-left in the hexagon) */
    221     f = sum(hexagon+i+max(0,i+1-n)*r, min(i+n,r+n-i-1), r, hexagon, hexagon+number_hex_entries);
    222     if (f==NOSOLUTION) return false;
    223     if (f==CHANGED) k = true;
    224     /* diagonal (down-right) */
    225     f = sum(hexagon-n+1+i+max(0,n-i-1)*(r+1), min(i+n,r+n-i-1), r+1, hexagon, hexagon+number_hex_entries);
    226     if (f==NOSOLUTION) return false;
    227     if (f==CHANGED) k = true;
    228   }
    229   if(k){
    230     goto restart;
    231   }
    232   return true;  // all done
    233 }
    234 
    235 void printhexagon(HexagonEntry hexagon[])
    236 {
    237   register unsigned long i,j;
    238 
    239   for (i=0; i<r; i++) {
    240     unsigned long l=0;
    241     unsigned long h=r;
    242     if (i+1>n)
    243       l = i+1-n;
    244     if (i+1<n)
    245       h = n+i;
    246     for (j=h-l; j<r; j++)
    247       printf("    ");
    248     for (j=l; j<h; j++) {
    249       assert(i<r);
    250       assert(j<r);
    251       HexagonEntry *hexagonEntry=&hexagon[i*r+j];
    252       assert(hexagonEntry->lo <= hexagonEntry->hi);
    253 #if 0
    254       printf("%6ld  ",hexagonEntry->id);
    255 #else
    256       if (hexagonEntry->lo < hexagonEntry->hi)
    257         printf("%4ld-%-3ld",hexagonEntry->lo,hexagonEntry->hi);
    258       else
    259         printf("%6ld  ",hexagonEntry->lo);
    260 #endif
    261     }
    262     printf("\n");
    263   }
    264 }
    265 
    266 /* make a spiral permutation of the hexagon variables */
    267 unsigned long *makeSpiralPermutation()
    268 {
    269   unsigned long level;
    270   unsigned long j;
    271   unsigned long i = 0;
    272   unsigned long index = 0;
    273 
    274   unsigned long *spiral = calloc(H,sizeof(unsigned long));
    275   unsigned long steps[] = {1, r+1, r, -1, -r-1, -r};
    276 
    277   for (level = n-1; level != 0 ;level--) {
    278     for(j = 0; j < 6 * level; j++) {
    279       spiral[i++] = index;
    280       index += steps[j/level];
    281     }
    282     index += r+1;
    283   }
    284 
    285   spiral[i] = index;
    286   return spiral;
    287 }
    288 
    289 /* make a spiral permutation of the hexagon variables with the corners
    290    in the first six entries */
    291 unsigned long *makeCornerSpiralPermutation()
    292 {
    293   unsigned long *spiral = makeSpiralPermutation();
    294   long i;
    295   long j;
    296   long k;
    297   
    298   if (n > 2) {
    299     unsigned long corners[6];
    300     for (i=0; i<6; i++) {
    301       corners[i] = spiral[i*(n-1)];
    302     }
    303     for (i=4*(n-1)+1, j=i+1; i>0; i-=(n-1),j-=(n-2)) {
    304       for (k=n-3; k>=0; k--) {
    305         spiral[j+k] = spiral[i+k];
    306       } 
    307     }
    308     for (i=0; i<6; i++) {
    309       spiral[i] = corners[i];
    310     }
    311   }
    312 
    313   return spiral;
    314 }
    315 
    316 /* assign values to hexagon[index] and all later variables in hexagon such that
    317    the constraints hold */
    318 void labeling(HexagonEntry hexagon[], unsigned long index, unsigned long *order)
    319 {
    320   unsigned long entryIndex = order[index];
    321 
    322   HexagonEntry *hexagonEntry = &hexagon[entryIndex];
    323 
    324   if (index >= H) {
    325     // All hexagonEntries fully constrained, print solution.
    326     printhexagon(hexagon);
    327     solutions++;
    328     leafs++;
    329     printf("leafs visited: %lu\n\n",leafs);
    330     return;
    331   }
    332 
    333   if (hexagonEntry->id == PLACEHOLDER_ENTRY_ID)
    334     // Do not process placeholder entries, continue to next hexagon index.
    335     return labeling(hexagon,index+1,order); 
    336   if(hexagonEntry->lo == hexagonEntry->hi)
    337     // hexagonEntry already fully constrained, continue to next hexagon index.
    338     return labeling(hexagon,index+1,order);
    339 
    340 
    341 
    342   HexagonEntry newHexagon[r*r];
    343   // bisect the range [lo,hi] to new Hexagons and recurse on them.
    344   long midpoint = (hexagonEntry->lo + hexagonEntry->hi)/2;
    345   if(hexagonEntry->lo + hexagonEntry->hi < 0 && (hexagonEntry->lo + hexagonEntry->hi) % 2 != 0){
    346     midpoint--;  // nudge midpoint for negative values to avoid endless recursion
    347   }
    348   // Initalize newHexagon to hexagon to explore lower halve of original range [lo,hi]
    349   memmove(newHexagon,hexagon,number_hex_entries*sizeof(HexagonEntry));
    350   newHexagon[entryIndex].lo = hexagonEntry->lo;
    351   newHexagon[entryIndex].hi = midpoint;
    352   if (solve(newHexagon)) {
    353     labeling(newHexagon,index,order);
    354   }
    355   else
    356     leafs++;
    357   // Initalize newHexagon to hexagon to explore upper halve of original range [lo,hi]
    358   memmove(newHexagon,hexagon,number_hex_entries*sizeof(HexagonEntry));
    359   newHexagon[entryIndex].lo = midpoint+1;
    360   newHexagon[entryIndex].hi = hexagonEntry->hi;
    361   if (solve(newHexagon)) {
    362     labeling(newHexagon,index,order); 
    363     }
    364   else
    365     leafs++;
    366 }
    367 
    368 HexagonEntry *makehexagon()
    369 {
    370   register unsigned long i,j;
    371 
    372   
    373   HexagonEntry *hexagon = calloc(number_hex_entries,sizeof(HexagonEntry));
    374   unsigned long id = 0;
    375   for (i=0; i<number_hex_entries; i++) {
    376     hexagon[i].id = PLACEHOLDER_ENTRY_ID;
    377     hexagon[i].lo = 1;
    378     hexagon[i].hi = 0;
    379   }
    380   for (i=0; i<r; i++) {
    381     unsigned long l=0;
    382     unsigned long h=r;
    383     if (i+1>n)
    384       l = i+1-n;
    385     if (i+1<n)
    386       h = n+i;
    387     for (j=l; j<h; j++) {
    388       assert(i<r);
    389       assert(j<r);
    390       assert(hexagon[i*r+j].lo>hexagon[i*r+j].hi);
    391       hexagon[i*r+j].id = id++;
    392       hexagon[i*r+j].lo = d*r - (H-1)/2;
    393       hexagon[i*r+j].hi = d*r + (H-1)/2;
    394     }
    395   }
    396   return hexagon;
    397 }
    398 
    399 int main(int argc, char *argv[])
    400 {
    401   register unsigned long i;
    402   register unsigned long j=0;
    403 
    404   if (argc < 3) {
    405     fprintf(stderr, "Usage: %s <order> <deviation> <value> ... <value>\n", argv[0]);
    406     exit(1);
    407   }
    408   n = strtoul(argv[1],NULL,10);
    409   if (n<1) {
    410     fprintf(stderr, "order must be >=1\n");
    411     exit(1);
    412   }
    413   d = strtol(argv[2],NULL,10);
    414   r = 2*n-1;
    415   number_hex_entries = r*r;
    416   H = 3*n*n-3*n+1;
    417   M = d*H;
    418 
    419   HexagonEntry *hexagon = makehexagon();
    420   for (i=3; i<argc; i++) {
    421     while (hexagon[j].id == PLACEHOLDER_ENTRY_ID)
    422       j++; // skip this entry as it is a placeholder
    423     hexagon[j].lo = hexagon[j].hi = strtol(argv[i],NULL,10);
    424     j++;
    425   }
    426   labeling(hexagon,0,makeCornerSpiralPermutation());
    427   printf("%lu solution(s), %lu leafs visited\n",solutions, leafs);
    428   //(void)solve(n, d, hexagon);
    429   //printhexagon(n, hexagon);
    430   return 0;
    431 }