API tools faq
paste
Guest User

ntqfind.c

a guest
Nov 29th, 2017
218
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C 47.24 KB | None | 0 0
raw download clone embed print report
  1. /* ntqfind v0.2
  2. ** A spaceship search program by Matthias Merzenich, with modifications by Aidan Pierce
  3. ** Based on code by David Eppstein, Paul Tooke, and "zdr"
  4. ** This is a first attempt at combining gfind and zfind.
  5. ** The code is currently a mess, but I hope to improve it over time.
  6. */
  7.  
  8. #include <stdio.h>
  9. #include <stdlib.h>
  10. #include <string.h>
  11. #include <stdint.h>
  12. #include <omp.h>
  13. #include "nttable.c"
  14.  
  15. #define BANNER "ntqfind v0.2 by Matthias Merzenich and Aidan Pierce, 29 November 2017"
  16. #define FILEVERSION ((unsigned long) 2017070501)  //yyyymmddnn, same as last qfind release (differs from the above)
  17.  
  18. #define CHUNK_SIZE 64
  19.  
  20. #define QBITS 23
  21.  
  22. #define HASHBITS 21
  23. #define DEFAULT_DEPTHLIMIT (qBits-3)
  24. #define SEARCHLIMIT 50
  25.  
  26. #define P_RULE 0
  27. #define P_WIDTH 1
  28. #define P_PERIOD 2
  29. #define P_OFFSET 3
  30. #define P_SYMMETRY 4
  31. #define P_NUM_SHIPS 5
  32. #define P_REORDER 6
  33. #define P_CHECKPOINT 7
  34. #define P_BASEBITS 8
  35. #define P_QBITS 9
  36. #define P_HASHBITS 10
  37. #define P_DEPTHLIMIT 11
  38. #define P_NUMTHREADS 12
  39. #define P_INIT_ROWS 13
  40. #define P_MINDEEP 14
  41.  
  42. #define NUM_PARAMS 15
  43.  
  44. #define SYM_ASYM 1
  45. #define SYM_ODD 2
  46. #define SYM_EVEN 3
  47. #define SYM_GUTTER 4
  48.  
  49. int bc[8] = {0, 1, 1, 2, 1, 2, 2, 3};
  50.  
  51. int params[NUM_PARAMS];
  52. int width;
  53. int deepeningAmount;
  54. int lastdeep;
  55. int nRowsInState;
  56. int phase;
  57.  
  58. int period;
  59. #define MAXPERIOD 30
  60.  
  61. int fwdOff[MAXPERIOD], backOff[MAXPERIOD], doubleOff[MAXPERIOD], tripleOff[MAXPERIOD];
  62.  
  63. int offset;
  64. unsigned long rule;
  65.  
  66. int aborting;
  67. int nFound;
  68.  
  69. enum mode {
  70.    asymmetric,          /* basic orthogonal or diagonal pattern */
  71.    odd, even,           /* orthogonal with bilateral symmetry */
  72.    gutter,              /* orthogonal bilateral symmetry with empty column in middle */
  73. } mode;
  74.  
  75. /* the big data structures */
  76. #define qBits params[P_QBITS]
  77. #define QSIZE (1<<qBits)
  78.  
  79. #define hashBits params[P_HASHBITS]
  80. #define HASHSIZE (1<<hashBits)
  81. #define HASHMASK (HASHSIZE - 1)
  82.  
  83. typedef uint32_t node;
  84. typedef uint16_t row;
  85.  
  86. row * rows;
  87. node * base;
  88. node * hash;
  89.  
  90. /*
  91. ** Representation of vertices.
  92. **
  93. ** Each vertex is represented by an entry in the rows[] array.
  94. ** That entry consists of the bits of the last row of the vertex's pattern
  95. ** concatenated with a number, the "offset" of the vertex.
  96. ** The parent of vertex i is formed by adding the offset of i
  97. ** with the number in base[i/BASEFACTOR].
  98. **
  99. ** If node i has offset -1, the entry is considered to be empty and is skipped.
  100. ** This is to deal with the case that base[i/BASEFACTOR] is already set to
  101. ** something too far away to deal with via the offset.
  102. **
  103. ** qIsEmpty() checks the size of the queue
  104. ** enqueue(n,r) adds a new node with last row r and parent n
  105. ** dequeue() returns the index of the next unprocessed node
  106. ** pop() undoes the previous enqueue operation
  107. ** resetQ() clears the queue
  108. **
  109. ** ROW(b) returns the last row of b
  110. ** PARENT(b) returns the index of the parent of b
  111. */
  112.  
  113. #define MAXWIDTH (10)
  114.  
  115. #define ROWBITS ((1<<width)-1)
  116. #define BASEBITS (params[P_BASEBITS])
  117. #define BASEFACTOR (1<<BASEBITS)
  118. #define MAXOFFSET ((((row) -1) >> width) - 1)
  119.  
  120. #define ROW(i) (rows[i] & ROWBITS)
  121. #define OFFSET(i) (rows[i] >> width)
  122. #define EMPTY(i) (rows[i] == (row)-1)
  123. #define MAKEEMPTY(i) rows[i] = (row)-1
  124. #define PARENT(i) (base[(i)>>BASEBITS]+OFFSET(i))
  125. #define FIRSTBASE(i) (((i) & ((1<<BASEBITS) - 1)) == 0)
  126.  
  127. #define MINDEEP  ((params[P_MINDEEP]>0) ? params[P_MINDEEP] : period)
  128.  
  129. void printRow(row theRow){
  130.    int i;
  131.    for(i = width - 1; i >= 0; i--) printf("%c",(theRow & 1 << i ? 'o' : '.'));
  132.    printf("\n");
  133. }
  134.  
  135. /* =========================================== */
  136. /*  Lookup Tables to determine successor rows  */
  137. /* =========================================== */
  138.  
  139. uint32_t *gInd;
  140. uint32_t *gcount;
  141. uint16_t *gRows;
  142.  
  143. //uint32_t *pInd, *pRemain;
  144. //uint16_t *pRows;
  145.  
  146. unsigned char *causesBirth;
  147.  
  148. int evolveBit(int row1, int row2, int row3, int bshift){
  149.    return nttable[(((row2>>bshift) & 2)<<7) | (((row1>>bshift) & 2)<<6)
  150.                 | (((row1>>bshift) & 4)<<4) | (((row2>>bshift) & 4)<<3)
  151.                 | (((row3>>bshift) & 7)<<2) | (((row2>>bshift) & 1)<<1)
  152.                 |  ((row1>>bshift) & 1)<<0];
  153. }
  154.  
  155. int evolveRow(int row1, int row2, int row3){
  156.    int row4;
  157.    int row1_s,row2_s,row3_s;
  158.    int j,s = 0;
  159.    if(params[P_SYMMETRY] == SYM_ODD) s = 1;
  160.    if(evolveBit(row1, row2, row3, width - 1)) return -1;
  161.    if(params[P_SYMMETRY] == SYM_ASYM && evolveBit(row1 << 2, row2 << 2, row3 << 2, 0)) return -1;
  162.    if(params[P_SYMMETRY] == SYM_ODD || params[P_SYMMETRY] == SYM_EVEN){
  163.       row1_s = (row1 << 1) + ((row1 >> s) & 1);
  164.       row2_s = (row2 << 1) + ((row2 >> s) & 1);
  165.       row3_s = (row3 << 1) + ((row3 >> s) & 1);
  166.    }
  167.    else{
  168.       row1_s = (row1 << 1);
  169.       row2_s = (row2 << 1);
  170.       row3_s = (row3 << 1);
  171.    }
  172.    row4 = evolveBit(row1_s, row2_s, row3_s, 0);
  173.    for(j = 1; j < width; j++)row4 += evolveBit(row1, row2, row3, j - 1) << j;
  174.    return row4;
  175. }
  176.  
  177. void sortRows(uint32_t rowSet){
  178.    uint32_t totalRows = gInd[rowSet + 1] - gInd[rowSet];
  179.    uint16_t *theRow = &(gRows[gInd[rowSet]]);
  180.    uint32_t i;
  181.    int64_t j;
  182.    uint16_t t;
  183.    for(i = 1; i < totalRows; ++i){
  184.       t = theRow[i];
  185.       j = i - 1;
  186.       while(j >= 0 && gcount[theRow[j]] < gcount[t]){
  187.          theRow[j+1] = theRow[j];
  188.          --j;
  189.       }
  190.       theRow[j+1] = t;
  191.    }
  192. }
  193.  
  194. void makeTables(){
  195.    printf("\nBuilding lookup tables... ");
  196.  
  197.    causesBirth = malloc((long long)sizeof(*causesBirth) << width);
  198.  
  199.    gInd = malloc(((long long)sizeof(*gInd) << (width * 3)) + sizeof(*gInd));
  200.    gcount = malloc((long long)sizeof(*gcount) * (1 << width));
  201.    uint32_t i;
  202.    int row1,row2,row3,row4;
  203.    long int rows123,rows124;
  204.    uint32_t numValid = 0;
  205.  
  206.    for(row1 = 0; row1 < 1 << width; row1++) causesBirth[row1] = (evolveRow(row1,0,0) ? 1 : 0);
  207.  
  208.    for(i = 0; i < 1 << width; ++i) gcount[i] = 0;
  209.    for(i = 0; i < ((1 << (3 * width)) + 1); i++)gInd[i] = 0;
  210.    rows123 = -1;     //represents row1, row2, and row3 stacked vertically
  211.    for(row1 = 0; row1 < 1 << width; row1++)for(row2 = 0; row2 < 1 << width; row2++)for(row3 = 0; row3 < 1 << width; row3++){
  212.       rows123++;
  213.       row4 = evolveRow(row1,row2,row3);
  214.       if(row4 < 0) continue;
  215.       ++gcount[row4];
  216.       gInd[rows123 - row3 + row4]++;
  217.       numValid++;
  218.    }
  219.    gRows = malloc(2 * numValid);
  220.    for(rows124 = 1; rows124 < 1 << (3 * width); rows124++) gInd[rows124] += gInd[rows124 - 1];
  221.    gInd[1 << (3 * width)] = gInd[(1 << (3 * width)) - 1];  //extra needed for last set to calculate number
  222.    rows123 = -1;
  223.    for(row1 = 0; row1 < 1 << width; row1++)for(row2 = 0; row2 < 1 << width; row2++)for(row3 = 0; row3 < 1 << width; row3++){
  224.       rows123++;
  225.       row4 = evolveRow(row1,row2,row3);
  226.       if(row4 < 0) continue;
  227.       rows124 = rows123 - row3 + row4;
  228.       gInd[rows124]--;
  229.       gRows[gInd[rows124]] = (uint16_t)row3;
  230.    }
  231.    printf("Lookup tables built.\n");
  232.  
  233.    gcount[0] = UINT32_MAX;
  234.    if(params[P_REORDER]){
  235.       printf("Sorting lookup table..... ");
  236.       for(rows124 = 0; rows124 < 1 << (3 * width); ++rows124){
  237.          sortRows(rows124);
  238.       }
  239.       printf("Lookup table sorted.\n");
  240.    }
  241.    free(gcount);
  242. }
  243.  
  244. void makePhases(){
  245.    int i;
  246.    for (i = 0; i < period; i++) backOff[i] = -1;
  247.    i = 0;
  248.    for (;;) {
  249.       int j = offset;
  250.       while (backOff[(i+j)%period] >= 0 && j < period) j++;
  251.       if (j == period) {
  252.          backOff[i] = period-i;
  253.          break;
  254.       }
  255.       backOff[i] = j;
  256.       i = (i+j)%period;
  257.    }
  258.    for (i = 0; i < period; i++)
  259.       fwdOff[(i+backOff[i])%period] = backOff[i];
  260.    for (i = 0; i < period; i++) {
  261.       int j = (i - fwdOff[i]);
  262.       if (j < 0) j += period;
  263.       doubleOff[i] = fwdOff[i] + fwdOff[j];
  264.    }
  265.    for (i = 0; i <  period; i++){
  266.       int j = (i - fwdOff[i]);
  267.       if (j < 0) j += period;
  268.       tripleOff[i] = fwdOff[i] + doubleOff[j];
  269.    }
  270. }
  271.  
  272. /* ====================================================== */
  273. /*  Hash table for detecting equivalent partial patterns  */
  274. /* ====================================================== */
  275.  
  276. void resetHash() { if (hash != 0) memset(hash,0,4*HASHSIZE); }
  277.  
  278. int hashPhase = 0;
  279.  
  280. static inline long hashFunction(node b, row r) {
  281.    long h = r;
  282.    int i;
  283.    for (i = 0; i < nRowsInState; i++) {
  284.       h = (h * 269) + ROW(b);
  285.       b = PARENT(b);
  286.    }
  287.    h += (h>>16)*269;
  288.    h += (h>>8)*269;
  289.    return h & HASHMASK;
  290. }
  291.  
  292. /* test if q+r is same as p */
  293. static inline int same(node p, node q, row r) {
  294.    int i;
  295.    for (i = 0; i < nRowsInState; i++) {
  296.       if (p >= QSIZE || q >= QSIZE || EMPTY(p) || EMPTY(q)) return 0;   /* sanity check */
  297.       if (ROW(p) != r) return 0;
  298.       p = PARENT(p);
  299.       r = ROW(q);
  300.       q = PARENT(q);
  301.    }
  302.    return 1;
  303. }
  304.  
  305. /* test if we've seen this child before */
  306. //void success(node);
  307. static inline int isVisited(node b, row r) {
  308.    if (same(0,b,r)) return 1;
  309.    if (hash != 0) {
  310.       int hashVal = hashFunction(b,r);
  311.       node hashNode = hash[hashVal];
  312.       if (hashNode == 0) return 0;
  313.       if (same(hashNode,b,r)){
  314.          return 1;
  315.       }
  316.    }
  317.    return 0;
  318. }
  319.  
  320. /* set node (NOT child) to be visited */
  321. static inline void setVisited(node b) {
  322.    if (hash != 0) hash[ hashFunction(PARENT(b),ROW(b)) ] = b;
  323. }
  324.  
  325.  
  326.  
  327. /* ============================================== */
  328. /*  Output patterns found by successful searches  */
  329. /* ============================================== */
  330.  
  331. #define MAXRLELINEWIDTH 63
  332. int RLEcount = 0;
  333. int RLElineWidth = 0;
  334. char RLEchar;
  335.  
  336. void sendRLE(char c) {
  337.   if (RLEcount > 0 && c != RLEchar) {
  338.       if (RLElineWidth++ >= MAXRLELINEWIDTH) {
  339.          if (RLEchar != '\n') putchar('\n');
  340.          RLElineWidth = 0;
  341.       }
  342.     if (RLEcount == 1) putchar(RLEchar);
  343.     else {
  344.        printf("%d%c", RLEcount, RLEchar);
  345.        RLElineWidth ++;
  346.        if (RLEcount > 9) RLElineWidth++;
  347.     }
  348.     RLEcount = 0;
  349.     if (RLEchar == '\n') RLElineWidth = 0;
  350.   }
  351.   if (c != '\0') {
  352.     RLEcount++;
  353.     RLEchar = c;
  354.   } else RLElineWidth = 0;
  355. }
  356.  
  357. int outputParity = 0;
  358.  
  359. void putRow(unsigned long rr, unsigned long r, int shift) {
  360.    while (r | rr) {
  361.       if (shift == 0)
  362.          sendRLE(r & 1 ? 'o' : 'b');
  363.       else shift--;
  364.       r >>= 1;
  365.       if (rr & 1) r |= (1<<31);
  366.       rr >>= 1;
  367.    }
  368.    sendRLE('$');
  369. }
  370.  
  371. /*void printRule() {
  372.   int i;
  373.   printf("B");
  374.   for (i = 0; i < 9; i++) if (rule & (1 << (9+i))) putchar(i+'0');
  375.   printf("/S");
  376.   for (i = 0; i < 9; i++) if (rule & (1 << i)) putchar(i+'0');
  377. }*/
  378.  
  379. void printRule() {
  380.    int i;
  381.    printf("B");
  382.    for(i = 0; i < 9; i++) if(rule & (1 << i)) putchar(i+'0');
  383.    printf("/S");
  384.    for(i = 0; i < 9; i++) if(rule & (1 << (i+9))) putchar(i+'0');
  385. }
  386.  
  387. /*void printRule() {
  388.    int i;
  389.    printf("B");
  390.    for(i = 0; i < 9; i++){
  391.       if(rule & (1 << i)) printf("%d",i);
  392.    }
  393.    printf("/S");
  394.    for(i = 9; i < 18; i++){
  395.       if(rule & (1 << i)) printf("%d",i - 9);
  396.    }
  397. }*/
  398.  
  399.  
  400. int modeWidth() {
  401.    switch(mode) {
  402.       case asymmetric:
  403.          return width;
  404.       case odd:
  405.          return 2*width-1;
  406.       case even:
  407.          return 2*width;
  408.       case gutter:
  409.          return 2*width+1;
  410.    }
  411.    return 0;
  412. }
  413.  
  414. /* PATCH - Keep static array for output to prevent memory leak */
  415. //int sxsAllocRows =  0;
  416. //unsigned long *  sxsAllocData;
  417. //unsigned long *  sxsAllocData2;
  418.  
  419. /* PATCH - Avoid Intel shift bug */
  420. static inline unsigned long
  421. safeShift(unsigned long r, int i)
  422. {
  423.     unsigned long rr = r;
  424.     while (i>16) { rr >>= 16; i-=16;}
  425.     return (rr>>i);
  426. }
  427.  
  428.  
  429. void success(node b, row *pRows, int nodeRow, uint32_t lastRow){
  430.    node c;
  431.    int nrows = 0;
  432.    int skewAmount = 0;
  433.    int swidth;
  434.    int p, i, j, margin;
  435.    //row *srows, *ssrows, *drows, *ddrows;
  436.    unsigned long *srows, *ssrows, *drows, *ddrows;
  437.    static unsigned long *oldsrows = 0, *oldssrows = 0;
  438.    static unsigned long *olddrows = 0, *oldddrows = 0;
  439.    static int oldnrows = 0;
  440.  
  441.    uint32_t currRow = lastRow;
  442.    int nDeepRows = 0;
  443.    int nodeDiff;
  444.  
  445.    /* check if output disabled while deepening */
  446.    /*if (perdidor) {
  447.       perdidor = 2;
  448.       return;
  449.    }*/
  450.  
  451.    if(pRows != NULL){
  452.       while(pRows[currRow] == 0){
  453.          if(currRow == 0){
  454.             printf("Success called on search root!\n");
  455.             aborting = 1;
  456.             return;
  457.          }
  458.          currRow--;
  459.       }
  460.       nDeepRows = (currRow / period) - 1;
  461.       nodeDiff = nodeRow - period - (currRow%period);
  462.       nodeRow -= nodeDiff;
  463.  
  464.       for(j = 0; j < nodeDiff; j++){
  465.          b = PARENT(b);
  466.       }
  467.       currRow = currRow - period + 1;
  468.       nrows = nDeepRows;
  469.      
  470.    }
  471.    else{
  472.       /* shift until we find nonzero row.
  473.          then shift PERIOD-1 more times to get leading edge of glider */
  474.       while (ROW(b) == 0) {
  475.          b = PARENT(b);
  476.          if (b == 0) {
  477.             printf("Success called on search root!\n");
  478.             aborting = 1;
  479.             return;
  480.          }
  481.       }
  482.    }
  483.    if(nrows < 0) nrows = 0;
  484.    
  485.    for (p = 0; p < period-1; p++) b = PARENT(b);
  486.    if (b == 0) {
  487.       printf("Success called on search root!\n");
  488.       aborting = 1;
  489.       return;
  490.    }
  491.    
  492.    /* count rows */
  493.    c = b;
  494.    while (c != 0) {
  495.       for (p = 0; p < period; p++)
  496.          c = PARENT(c);
  497.       nrows++;
  498.    }
  499.    
  500.    /* build data structure of rows so we can reduce width etc */
  501.    srows = malloc((nrows+MAXWIDTH+1) * sizeof(unsigned long));
  502.    ssrows = malloc((nrows+MAXWIDTH+1) * sizeof(unsigned long));
  503.    drows = srows; ddrows = ssrows; /* save orig ptr for free() */
  504.    for (i = 0; i <= nrows+MAXWIDTH; i++) srows[i]=ssrows[i]=0;
  505.    for (i = nrows - 1; i >= 0; i--) {
  506.       row r;
  507.       if(nDeepRows > 0){
  508.          r = pRows[currRow];
  509.          currRow -= period;
  510.          nDeepRows--;
  511.       }
  512.       else{
  513.          r = ROW(b);
  514.          for (p = 0; p < period; p++) {
  515.             b = PARENT(b);
  516.          }
  517.       }
  518.       //row rx;
  519.       switch(mode) {
  520.          case asymmetric:
  521.             srows[i] = r;
  522.             break;
  523.  
  524.          case odd:
  525.             srows[i] = r << (MAXWIDTH - 1);
  526.             ssrows[i] = r >> (32 - (MAXWIDTH - 1));
  527.             for (j = 1; j < MAXWIDTH; j++)
  528.                if (r & (1<<j))
  529.                   srows[i] |= 1 << (MAXWIDTH - 1 - j);
  530.             break;
  531.  
  532.          case even:
  533.             srows[i] = r << MAXWIDTH;
  534.             ssrows[i] = r >> (32 - MAXWIDTH);
  535.             for (j = 0; j < MAXWIDTH; j++)
  536.                if (r & (1<<j))
  537.                   srows[i] |= 1 << (MAXWIDTH - 1 - j);
  538.             break;
  539.  
  540.          case gutter:
  541.             srows[i] = r << (MAXWIDTH + 1);
  542.             ssrows[i] = r >> (32 - (MAXWIDTH + 1));
  543.             for (j = 0; j < MAXWIDTH; j++)
  544.                if (r & (1<<j))
  545.                   srows[i+skewAmount] |= 1 << (MAXWIDTH - 1 - j);
  546.             break;
  547.  
  548.          default:
  549.             printf("Unexpected mode in success!\n");
  550.             aborting = 1;
  551.             return;
  552.       }
  553.    }
  554.    
  555.    /* normalize nrows to only include blank rows */
  556.    nrows += MAXWIDTH;
  557.    while (srows[nrows-1] == 0 && ssrows[nrows-1] == 0 && nrows>0) nrows--;
  558.    while (srows[0] == 0 && ssrows[0] == 0 && nrows>0) {
  559.       srows++;
  560.       ssrows++;
  561.       nrows--;
  562.    }
  563.    
  564.    /* sanity check: are all rows empty? */
  565.    int allEmpty = 1;
  566.    for(i = 0; i < nrows; i++){
  567.       if(srows[i]){
  568.          allEmpty = 0;
  569.          break;
  570.       }
  571.    }
  572.    if(allEmpty) return;
  573.    
  574.    /* make at least one row have nonzero first bit */
  575.    i = 0;
  576.    while ((srows[i] & 1) == 0) {
  577.       for (i = 0; (srows[i] & 1) == 0 && i < nrows; i++) { }
  578.       if (i == nrows) {
  579.          for (i = 0; i < nrows; i++) {
  580.             srows[i] >>= 1;
  581.             if (ssrows[i] & 1) srows[i] |= (1<<31);
  582.             ssrows[i] >>= 1;
  583.          }
  584.          i = 0;
  585.       }
  586.    }
  587.    
  588.    swidth = 0;
  589.    for (i = 0; i < nrows; i++)
  590.       while (safeShift(ssrows[i],swidth))
  591.          swidth++;
  592.    if (swidth) swidth += 32;
  593.    for (i = 0; i < nrows; i++)
  594.       while (safeShift(srows[i],swidth))
  595.          swidth++;
  596.    
  597.      
  598.    /* compute margin on other side of width */
  599.    margin = 0;
  600.  
  601.    /* make sure we didn't just output the exact same pattern (happens a lot for puffer) */
  602.    if (nrows == oldnrows) {
  603.       int different = 0;
  604.       for (i = 0; i < nrows && !different; i++)
  605.          different = (srows[i] != oldsrows[i] || ssrows[i] != oldssrows[i]);
  606.       if (!different) {
  607.          free(drows);
  608.          free(ddrows);
  609.          return;
  610.       }
  611.    }
  612.    if (olddrows != 0) free(olddrows);
  613.    if (oldddrows != 0) free(oldddrows);
  614.    oldsrows = srows;
  615.    oldssrows = ssrows;
  616.    olddrows = drows;
  617.    oldddrows = ddrows;
  618.    oldnrows = nrows;
  619.  
  620.    /* output it all */
  621.    printf("\nx = %d, y = %d, rule = ", swidth - margin, nrows);
  622.    printRule();
  623.    putchar('\n');
  624.  
  625.    while (nrows-- > 0) {
  626.       if (margin > nrows) putRow(ssrows[nrows], srows[nrows], margin - nrows);
  627.       else putRow(ssrows[nrows], srows[nrows], 0);
  628.    }
  629.    RLEchar = '!';
  630.    sendRLE('\0');
  631.    printf("\n\n");
  632.    fflush(stdout);
  633.    //if (++nFound >= findLimit) aborting = 1;
  634. }
  635.  
  636.  
  637. /* Is this a node at which we can stop? */
  638. int terminal(node n){
  639.    int p;
  640.  
  641.    for (p = 0; p < period; p++) {   /* last row in each phase must be zero */
  642.       if (ROW(n) != 0) return 0;
  643.       n = PARENT(n);
  644.    }
  645.  
  646.    for (p = 0; p < period; p++) {
  647.       if(causesBirth[ROW(n)]) return 0;
  648.       n = PARENT(n);
  649.    }
  650.    return 1;
  651. }
  652.  
  653.  
  654.  
  655.  
  656. /* ================================================ */
  657. /*  Queue of partial patterns still to be examined  */
  658. /* ================================================ */
  659.  
  660. /* SU patch */
  661. node qHead,qTail;
  662.  
  663. /* PATCH queue dimensions required during save/restore */
  664. node qStart; /* index of first node in queue */
  665. node qEnd;   /* index of first unused node after end of queue */
  666.  
  667. /* Maintain phase of queue nodes. After dequeue(), the global variable phase
  668.    gives the phase of the dequeued item.  If the queue is compacted, this information
  669.    needs to be reinitialized by a call to rephase(), after which phase will not be
  670.    valid until the next call to dequeue().  Variable nextRephase points to the next
  671.    node for which dequeue will need to increment the phase. Phase is not maintained
  672.    when treating queue as a stack (using pop()) -- caller must do it in that case.
  673.    It's ok to change phase since we maintain a separate copy in queuePhase. */
  674.  
  675. int queuePhase = 0;
  676. node nextRephase = 0;
  677. void rephase() {
  678.    node x, y;
  679.    while (qHead < qTail && EMPTY(qHead)) qHead++;   /* skip past empty queue cells */
  680.    x = qHead;   /* find next item in queue */
  681.    queuePhase = period - 1;
  682.    while (x != 0) {
  683.       x = PARENT(x);
  684.       queuePhase++;
  685.    }
  686.    queuePhase %= period;
  687.  
  688.    /* now walk forward through queue finding breakpoints between each generation
  689.       invariants: y is always the first in its generation */
  690.    x = 0; y = 0;
  691.    while (y <= qHead) {
  692.       ++x;
  693.       if (x >= qTail || (!EMPTY(x) && PARENT(x) >= y)) y = x;
  694.    }
  695.    nextRephase = y;
  696. }
  697.  
  698. /* phase of an item on the queue */
  699. int peekPhase(node i) {
  700.    return (i < nextRephase? queuePhase : (queuePhase+1)%period);
  701. }
  702.  
  703. /* Test queue status */
  704. static inline int qIsEmpty() {
  705.    while (qHead < qTail && EMPTY(qHead)) qHead++;
  706.    return (qTail == qHead);
  707. }
  708.  
  709. void qFull() {
  710.     if (aborting != 2) {
  711.       printf("Exceeded %d node limit, search aborted\n", QSIZE);
  712.       fflush(stdout);
  713.       aborting = 2;
  714.    }
  715. }
  716.  
  717. static inline void enqueue(node b, row r) {
  718.    node i = qTail++;
  719.    if (i >= QSIZE) qFull();
  720.    else if (FIRSTBASE(i)) {
  721.       base[i>>BASEBITS] = b;
  722.       rows[i] = r;
  723.    } else {
  724.       long o = b - base[i>>BASEBITS];
  725.       if (o < 0 || o >(long) MAXOFFSET) {   /* offset out of range */
  726.          while (!FIRSTBASE(i)) {
  727.             rows[i] = -1;
  728.             i = qTail++;
  729.             if (i >= QSIZE) qFull();
  730.          }
  731.          base[i>>BASEBITS] = b;
  732.          rows[i] = r;
  733.       } else rows[i] = (o << width) + r;
  734.    }
  735. }
  736.  
  737. static inline node dequeue() {
  738.    while (qHead < qTail && EMPTY(qHead)) qHead++;
  739.    if (qHead >= nextRephase) {
  740.       queuePhase = (queuePhase+1)%period;
  741.       nextRephase = qTail;
  742.    }
  743.    phase = queuePhase;
  744.    return qHead++;
  745. }
  746.  
  747. static inline void pop() {
  748.    qTail--;
  749.    while (qTail > qHead && EMPTY(qTail-1)) qTail--;
  750. }
  751.  
  752. void resetQ() { qHead = qTail = 0; }
  753.  
  754. static inline int qTop() { return qTail - 1; }
  755.  
  756.  
  757. /* ================================= */
  758. /* PATCH08 - dump state              */
  759. /* ================================= */
  760.  
  761. int dumpNum = 1;
  762. char dumpFile[12];
  763. #define DUMPROOT "dump"
  764. int dumpFlag = 0; /* Dump status flags, possible values follow */
  765. #define DUMPPENDING (1)
  766. #define DUMPFAILURE (2)
  767. #define DUMPSUCCESS (3)
  768.  
  769. int dumpandexit = 0;
  770.  
  771. FILE * openDumpFile()
  772. {
  773.     FILE * fp;
  774.  
  775.     while (dumpNum < 10000)
  776.     {
  777.         sprintf(dumpFile,"%s%04d",DUMPROOT,dumpNum++);
  778.         if ((fp=fopen(dumpFile,"r")))
  779.             fclose(fp);
  780.         else
  781.             return fopen(dumpFile,"w");
  782.     }
  783.     return (FILE *) 0;
  784. }
  785.  
  786. void dumpState()
  787. {
  788.     FILE * fp;
  789.     int i;
  790.     dumpFlag = DUMPFAILURE;
  791.     if (!(fp = openDumpFile())) return;
  792.     fprintf(fp,"%lu\n",FILEVERSION);
  793.     for (i = 0; i < NUM_PARAMS; i++)
  794.         fprintf(fp,"%d\n",params[i]);
  795.     fprintf(fp,"%d\n",width);
  796.     //fprintf(fp,"%d\n",widthReduction);
  797.     fprintf(fp,"%d\n",period);
  798.     fprintf(fp,"%d\n",offset);
  799.     fprintf(fp,"%d\n",mode);
  800.     //fprintf(fp,"%d\n",diagonal);
  801.     //fprintf(fp,"%d\n",nFound);
  802.     fprintf(fp,"%d\n",lastdeep);
  803.  
  804.     fprintf(fp,"%u\n",qHead-qStart);
  805.     fprintf(fp,"%u\n",qEnd-qStart);
  806.     for (i = qStart; i < qEnd; i++)
  807.         fprintf(fp,"%u\n",rows[i]);
  808.     fclose(fp);
  809.     dumpFlag = DUMPSUCCESS;
  810. }
  811.  
  812.  
  813.  
  814.  
  815.  
  816. /* ================================= */
  817. /*  Compaction of nearly full queue  */
  818. /* ================================= */
  819.  
  820. void putnum(long n) {
  821.    char suffix;
  822.    if (n >= 1000000) {
  823.       n /= 100000;
  824.       suffix = 'M';
  825.    } else if (n >= 1000) {
  826.       n /= 100;
  827.       suffix = 'k';
  828.    } else {
  829.       printf("%ld", n);
  830.       return;
  831.    }
  832.  
  833.    if (n >= 100) printf("%ld", n/10);
  834.    else printf("%ld.%ld", n/10, n%10);
  835.    putchar(suffix);
  836. }
  837.  
  838. long currentDepth() {
  839.    long i;
  840.    node x;
  841.    x = qTail - 1;
  842.    i = 1;
  843.    while (x != 0) {
  844.       x = PARENT(x);
  845.       i++;
  846.    }
  847.    return i;
  848. }
  849.  
  850.  
  851.  
  852. /*
  853.    doCompact() now has two parts.  The first part compresses
  854.    the queue.  The second part consists of the last loop which
  855.    converts parent bits to back parent pointers. The search
  856.    state may be saved in between.  The queue dimensions, which
  857.    were previously saved in local variables are saved in globals.
  858. */
  859.  
  860. void doCompactPart1()
  861. {
  862.    node x,y;
  863.    qEnd = qTail;
  864.    
  865.    /* make a pass backwards from the end finding unused nodes at or before qHead */
  866.    x = qTail - 1;
  867.    y = qHead - 1;
  868.    while (y > 0) {
  869.       /* invariants: everything after y is still active.
  870.                      everything after x points to something after y.
  871.                      x is nonempty and points to y or something before y.
  872.                      so, if x doesnt point to y, y must be unused and can be removed. */
  873.       if (!EMPTY(y)) {
  874.          if (y > PARENT(x)) rows[y] = -1;
  875.          else while (EMPTY(x) || PARENT(x) == y) x--;
  876.       }
  877.       y--;
  878.    }
  879.    
  880.    /* make a pass forwards converting parent pointers to offset from prev parent ptr.
  881.       note that after unused nodes are eliminated, all these offsets are zero or one. */
  882.    y = 0;
  883.    for (x = 0; x < qTail; x++) if (!EMPTY(x)) {
  884.       if (PARENT(x) == y) rows[x] = ROW(x);
  885.       else {
  886.          y = PARENT(x);
  887.          rows[x] = (1<<width) + ROW(x);
  888.       }
  889.    }
  890.    
  891.    /*
  892.       Make a pass backwards compacting gaps.
  893.    
  894.       For most times we run this, it could be combined with the next phase, but
  895.       every once in a while the process of repacking the remaining items causes them
  896.       to use *MORE* space than they did before they were repacked (because of the need
  897.       to leave empty space when OFFSET gets too big) and without this phase the repacked
  898.       stuff overlaps the not-yet-repacked stuff causing major badness.
  899.      
  900.       For this phase, y points to the current item to be repacked, and x points
  901.       to the next free place to pack an item.
  902.     */
  903.    x = y = qTail-1;
  904.    for (;;) {
  905.       if (qHead == y) qHead = x;
  906.       if (!EMPTY(y)) {
  907.          rows[x] = rows[y];
  908.          x--;
  909.       }
  910.       if (y-- == 0) break;   /* circumlocution for while (y >= 0) because x is unsigned */
  911.    }
  912.     qStart = ++x;    /* mark start of queue */
  913. }
  914.  
  915. void doCompactPart2()
  916. {
  917.     node x,y;
  918.  
  919.    /*
  920.       Make a pass forwards converting parent bits back to parent pointers.
  921.      
  922.       For this phase, x points to the current item to be repacked, and y points
  923.       to the parent of the previously repacked item.
  924.       After the previous pass, x is initialized to first nonempty item,
  925.       and all items after x are nonempty.
  926.     */
  927.    qTail = 0; y = 0;
  928.    resetHash();
  929.    for (x = qStart; x < qEnd; x++) {
  930.       if (OFFSET(x)) {   /* skip forward to next parent */
  931.          y++;
  932.          while (EMPTY(y)) y++;
  933.       }
  934.       enqueue(y,ROW(x));
  935.       if (aborting) return;
  936.       if (qHead == x) qHead = qTail - 1;
  937.       setVisited(qTail - 1);
  938.    }
  939.    rephase();
  940. }
  941.  
  942. void doCompact()
  943. {
  944.    /* make sure we still have something left in the queue */
  945.    if (qIsEmpty()) {
  946.       qTail = qHead = 0;   /* nothing left, make an extremely compact queue */
  947.       return;
  948.    }
  949.     /* First loop of part 1 requires qTail-1 to be non-empty. Make it so */
  950.     while(EMPTY(qTail-1))
  951.         qTail--;
  952.  
  953.     doCompactPart1();
  954.     if (dumpFlag == DUMPPENDING) dumpState();
  955.     doCompactPart2();
  956. }
  957.  
  958.  
  959.  
  960. /* =================================== */
  961. /* PATCH08 - preview partial results   */
  962. /* =================================== */
  963.  
  964. static void preview(int allPhases)
  965. {
  966.     node i,j,k;
  967.     int ph;
  968.  
  969.     for (i = qHead; (i<qTail) && EMPTY(i); i++);
  970.     for (j = qTail-1; (j>i) && EMPTY(j); j--);
  971.     if (j<i) return;
  972. /*
  973.     for (ph = 0; ph < period; ph++)
  974.     {
  975.         i=PARENT(i);
  976.         j=PARENT(j);
  977.     }
  978. */
  979.     while (j>=i && !aborting)
  980.     {
  981.         if (!EMPTY(j))
  982.         {
  983.             success(j, NULL, 0, 0);
  984.             //success(j);
  985.             if (allPhases == 0)
  986.             {
  987.                 k=j;
  988.                 for (ph = 1; ph < period; ph++)
  989.                 {
  990.                     k=PARENT(k);
  991.                     success(k, NULL, 0, 0);
  992.                     //success(k);
  993.                 }
  994.             }
  995.         }
  996.         j--;
  997.     }
  998. }
  999.  
  1000.  
  1001. /* ================================= */
  1002. /* PATCH08 - resume from saved state */
  1003. /* ================================= */
  1004.  
  1005. char * loadFile;
  1006.  
  1007. void loadFail()
  1008. {
  1009.     printf("Load from file %s failed\n",loadFile);
  1010.     exit(1);
  1011. }
  1012.  
  1013. signed int loadInt(FILE *fp)
  1014. {
  1015.     signed int v;
  1016.     if (fscanf(fp,"%d\n",&v) != 1) loadFail();
  1017.     return v;
  1018. }
  1019.  
  1020. unsigned int loadUInt(FILE *fp)
  1021. {
  1022.     unsigned int v;
  1023.     if (fscanf(fp,"%u\n",&v) != 1) loadFail();
  1024.     return v;
  1025. }
  1026.  
  1027. void loadState(char * cmd, char * file)
  1028. {
  1029.     FILE * fp;
  1030.     int i;
  1031.  
  1032.     loadFile = file;
  1033.     fp = fopen(loadFile, "r");
  1034.     if (!fp) loadFail();
  1035.     if (loadUInt(fp) != FILEVERSION)
  1036.     {
  1037.         printf("Incompatible file version\n");
  1038.         exit(1);
  1039.     }
  1040.  
  1041.     /* Load parameters and set stuff that can be derived from them */
  1042.     for (i = 0; i < NUM_PARAMS; i++)
  1043.         params[i] = loadInt(fp);
  1044.    // rule = (params[P_BIRTHS] << 9) + params[P_SURVIVES];
  1045.    rule = params[P_RULE];
  1046.  
  1047.     /* Load / initialise globals */
  1048.     width          = loadInt(fp);
  1049.     //widthReduction = loadInt(fp);
  1050.     period         = loadInt(fp);
  1051.     offset         = loadInt(fp);
  1052.     mode           = loadInt(fp);
  1053.     //diagonal       = loadInt(fp);
  1054.     //nFound         = loadInt(fp);
  1055.     lastdeep       = loadInt(fp);
  1056.    deepeningAmount = period; /* Currently redundant, since it's recalculated */
  1057.    //perdidor        = 0;
  1058.    aborting        = 0;
  1059.    nRowsInState = period+period;   /* how many rows needed to compute successor graph? */
  1060.  
  1061.     /* Allocate space for the data structures */
  1062.    base = malloc((QSIZE>>BASEBITS)*sizeof(node));
  1063.    rows = malloc(QSIZE*sizeof(row));
  1064.    if (base == 0 || rows == 0) {
  1065.       printf("Unable to allocate BFS queue!\n");
  1066.       exit(0);
  1067.    }
  1068.    
  1069.    if (hashBits == 0) hash = 0;
  1070.    else {
  1071.       hash = malloc(HASHSIZE*sizeof(node));
  1072.       if (hash == 0) printf("Unable to allocate hash table, duplicate elimination disabled\n");
  1073.    }
  1074.  
  1075.     /* Various bits that need doing before we load the BFS queue */
  1076.    // makeReversals();
  1077.    //makeUnterm();
  1078.    // set4x();
  1079.    
  1080.     /* Load up BFS queue and complete compaction */
  1081.     qHead  = loadUInt(fp);
  1082.     qEnd   = loadUInt(fp);
  1083.     qStart = QSIZE - qEnd;
  1084.     qEnd   = QSIZE;
  1085.     qHead += qStart;
  1086.     if (qStart > QSIZE || qStart < QSIZE/16)
  1087.     {
  1088.       printf("BFS queue is too small for saved state\n");
  1089.       exit(0);
  1090.    }
  1091.    for (i = qStart; i < qEnd; i++)
  1092.         rows[i] = (row) loadUInt(fp);
  1093.     fclose(fp);
  1094. /*
  1095.    printf("qHead:  %d qStart: %d qEnd: %d\n",qHead,qStart,qEnd);
  1096.    printf("rows[0]: %d\n",rows[qStart]);
  1097.    printf("rows[1]: %d\n",rows[qStart+1]);
  1098.    printf("rows[2]: %d\n",rows[qStart+2]);
  1099.    fflush(stdout);
  1100.    exit(0);
  1101. */
  1102.     doCompactPart2();
  1103.  
  1104.  
  1105.     /* Let the user know that we got this far */
  1106.    printf("State successfully loaded from file %s",loadFile);
  1107.    
  1108.    if(!strcmp(cmd,"p") || !strcmp(cmd,"P")){
  1109.       preview(1);
  1110.       exit(0);
  1111.    }
  1112.    
  1113.    fflush(stdout);
  1114.    
  1115.    omp_set_num_threads(params[P_NUMTHREADS]);
  1116. }
  1117.  
  1118.  
  1119.  
  1120.  
  1121.  
  1122. int lookAhead(row *pRows, int a, int pPhase){
  1123.    uint32_t ri11, ri12, ri13, ri22, ri23;  //indices: first number represents vertical offset, second number represents generational offset
  1124.    uint32_t rowSet11, rowSet12, rowSet13, rowSet22, rowSet23, rowSet33;
  1125.    uint32_t riStart11, riStart12, riStart13, riStart22, riStart23;
  1126.    uint32_t numRows11, numRows12, numRows13, numRows22, numRows23;
  1127.    uint32_t row11, row12, row13, row22, row23;
  1128.    
  1129.    
  1130.    rowSet11 = (pRows[a - params[P_PERIOD] - fwdOff[pPhase]] << (2 * params[P_WIDTH]))
  1131.              +(pRows[a - fwdOff[pPhase]] << params[P_WIDTH])
  1132.              + pRows[a];
  1133.    riStart11 = gInd[rowSet11];
  1134.    numRows11 = gInd[rowSet11 + 1] - riStart11;
  1135.    
  1136.    if(!numRows11) return 0;
  1137.  
  1138.    rowSet12 = (pRows[a - params[P_PERIOD] - doubleOff[pPhase]] << (2 * params[P_WIDTH]))
  1139.              +(pRows[a - doubleOff[pPhase]] << params[P_WIDTH])
  1140.              + pRows[a - fwdOff[pPhase]];
  1141.    riStart12 = gInd[rowSet12];
  1142.    numRows12 = gInd[rowSet12 + 1] - riStart12;
  1143.  
  1144.    if(tripleOff[pPhase] >= params[P_PERIOD]){
  1145.       //riStart13 = pInd[a + params[P_PERIOD] - tripleOff[pPhase]] + pRemain[a + params[P_PERIOD] - tripleOff[pPhase]];
  1146.       numRows13 = 1;
  1147.    }
  1148.    else{
  1149.       rowSet13 = (pRows[a - params[P_PERIOD] - tripleOff[pPhase]] << (2 * params[P_WIDTH]))
  1150.                 +(pRows[a - tripleOff[pPhase]] << params[P_WIDTH])
  1151.                 + pRows[a - doubleOff[pPhase]];
  1152.       riStart13 = gInd[rowSet13];
  1153.       numRows13 = gInd[rowSet13 + 1] - riStart13;
  1154.    }
  1155.  
  1156.    for(ri11 = 0; ri11 < numRows11; ++ri11){
  1157.       row11 = gRows[ri11 + riStart11];
  1158.       for(ri12 = 0; ri12 < numRows12; ++ri12){
  1159.          row12 = gRows[ri12 + riStart12];
  1160.          rowSet22 = (pRows[a - doubleOff[pPhase]] << (2 * params[P_WIDTH]))
  1161.                    +(row12 << params[P_WIDTH])
  1162.                    + row11;
  1163.          riStart22 = gInd[rowSet22];
  1164.          numRows22 = gInd[rowSet22 + 1] - riStart22;
  1165.          if(!numRows22) continue;
  1166.  
  1167.          for(ri13 = 0; ri13 < numRows13; ++ri13){
  1168.             if(tripleOff[pPhase] >= params[P_PERIOD]){
  1169.                row13 = pRows[a + params[P_PERIOD] - tripleOff[pPhase]];
  1170.             }
  1171.             else{
  1172.                row13 = gRows[ri13 + riStart13];
  1173.             }
  1174.             rowSet23 = (pRows[a - tripleOff[pPhase]] << (2 * params[P_WIDTH]))
  1175.                       +(row13 << params[P_WIDTH])
  1176.                       + row12;
  1177.             riStart23 = gInd[rowSet23];
  1178.             numRows23 = gInd[rowSet23 + 1] - riStart23;
  1179.             if(!numRows23) continue;
  1180.  
  1181.             for(ri22 = 0; ri22 < numRows22; ++ri22){
  1182.                row22 = gRows[ri22 + riStart22];
  1183.                for(ri23 = 0; ri23 < numRows23; ++ri23){
  1184.                   row23 = gRows[ri23 + riStart23];
  1185.                   rowSet33 = (row13 << (2 * params[P_WIDTH]))
  1186.                             +(row23 << params[P_WIDTH])
  1187.                             + row22;
  1188.                   if(gInd[rowSet33] != gInd[rowSet33 + 1]) return 1;
  1189.                }
  1190.             }
  1191.          }
  1192.       }
  1193.    }
  1194.    return 0;
  1195. }
  1196.  
  1197.  
  1198.  
  1199. void process(node theNode)
  1200. {
  1201.    long long int i;
  1202.    int firstRow = 0;
  1203.    uint32_t numRows;
  1204.    uint32_t newRowSet;
  1205.    node x = theNode;
  1206.    int pPhase = peekPhase(x);
  1207.    row pRows[3*MAXPERIOD];
  1208.    int currRow = 2*period + pPhase + 1;
  1209.    for(i = currRow - 1; i >= 0; --i){
  1210.       pRows[i] = ROW(x);
  1211.       x = PARENT(x);
  1212.    }
  1213.  
  1214.    ++pPhase;
  1215.    if(pPhase == period) pPhase = 0;
  1216.  
  1217.    newRowSet = (pRows[currRow - 2 * period] << (2 * params[P_WIDTH]))
  1218.               +(pRows[currRow - period] << params[P_WIDTH])
  1219.               + pRows[currRow - period + backOff[pPhase]];
  1220.    numRows = gInd[newRowSet + 1] - gInd[newRowSet];
  1221.  
  1222.  
  1223.    if(theNode == 0){
  1224.       firstRow = 1;
  1225.    }
  1226.  
  1227.    for(i = firstRow; i < numRows; ++i){
  1228.       pRows[currRow] = gRows[gInd[newRowSet] + i];
  1229.       if (!isVisited(theNode, pRows[currRow]) && lookAhead(pRows, currRow, pPhase)){
  1230.          enqueue(theNode, pRows[currRow]);
  1231.          if (terminal(qTail-1)) success(qTail-1, NULL, 0, 0);
  1232.          setVisited(qTail - 1);
  1233.       }
  1234.    }
  1235. }
  1236.  
  1237. int depthFirst(node theNode, long howDeep, uint32_t *pInd, uint32_t *pRemain, row *pRows){
  1238. //   uint32_t *pInd;
  1239. //   uint32_t *pRemain;
  1240.  
  1241.    int pPhase;
  1242. //   row *pRows;
  1243.  
  1244.    uint32_t newRowSet;
  1245.    pPhase = peekPhase(theNode);
  1246.  
  1247.  
  1248.  
  1249.  
  1250. //   pInd = malloc((long long)sizeof(*pInd) * (howDeep + 3 * params[P_PERIOD]));
  1251. //   pRemain = malloc((long long)sizeof(*pRemain) * (howDeep + 3 * params[P_PERIOD]));
  1252. //   pRows = malloc((long long)sizeof(*pRows) * (howDeep + 3 * params[P_PERIOD]));
  1253.  
  1254.    node x = theNode;
  1255.    uint32_t startRow = 2*period + pPhase + 1;
  1256.    uint32_t currRow = startRow;
  1257.    
  1258.    int i;
  1259.    for(i = currRow - 1; i >= 0; --i){
  1260.       pRows[i] = ROW(x);
  1261.       x = PARENT(x);
  1262.    }
  1263.  
  1264.    ++pPhase;
  1265.    if(pPhase == period) pPhase = 0;
  1266.  
  1267.    newRowSet = (pRows[currRow - 2 * period] << (2 * width))
  1268.               |(pRows[currRow - period] << width)
  1269.               | pRows[currRow - period + backOff[pPhase]];
  1270.  
  1271.    pRemain[currRow] = gInd[newRowSet + 1] - gInd[newRowSet];
  1272.    pInd[currRow] = gInd[newRowSet + 1];
  1273.  
  1274.    
  1275.    for(;;){
  1276.  
  1277.       if(!pRemain[currRow]){
  1278.          --currRow;
  1279.          if(pPhase == 0) pPhase = period;
  1280.          --pPhase;
  1281.          if(currRow < startRow){
  1282. //            free(pInd);
  1283. //            free(pRemain);
  1284. //            free(pRows);
  1285.             return 0;
  1286.          }
  1287.          continue;
  1288.       }
  1289.       //--pRemain[currRow];
  1290.  
  1291.       pRows[currRow] = gRows[pInd[currRow] - pRemain[currRow]];
  1292.       --pRemain[currRow];
  1293.       if(!lookAhead(pRows, currRow, pPhase)) continue;
  1294.  
  1295.       ++currRow;
  1296.       ++pPhase;
  1297.       if(pPhase == period) pPhase = 0;
  1298.       if(currRow > startRow + howDeep){
  1299. //         free(pInd);
  1300. //         free(pRemain);
  1301. //         free(pRows);
  1302.          for(i = 1; i <= period; ++i){
  1303.             if(pRows[currRow - i]) return 1;
  1304.          }
  1305.          currRow -= period;
  1306.          for(i = 1; i<= period; ++i){
  1307.             if(causesBirth[pRows[currRow - i]]) return 1;
  1308.          }
  1309.          success(theNode, pRows, startRow - 1, currRow + period - 1);
  1310.          
  1311.          return 1;
  1312.          
  1313.       }
  1314.  
  1315.       newRowSet = (pRows[currRow - 2 * period] << (2 * params[P_WIDTH]))
  1316.                  +(pRows[currRow - period] << params[P_WIDTH])
  1317.                  + pRows[currRow - period + backOff[pPhase]];
  1318.       pRemain[currRow] = gInd[newRowSet + 1] - gInd[newRowSet];
  1319.       pInd[currRow] = gInd[newRowSet + 1];
  1320.  
  1321.    }
  1322. }
  1323.  
  1324.  
  1325. static void deepen(){
  1326.    node i;
  1327.    //node j;
  1328.  
  1329. //   if (findLimit > 1) perdidor = 1;   /* disable success if want mult pattern output */
  1330.  
  1331.    /* compute amount to deepen, apply reduction if too deep */
  1332. #ifdef PATCH07
  1333.    timeStamp();
  1334. #endif
  1335.    printf("Queue full");
  1336.    i = currentDepth();
  1337.    if (i >= lastdeep) deepeningAmount = MINDEEP;
  1338.    else deepeningAmount = lastdeep + MINDEEP - i;   /* go at least MINDEEP deeper */
  1339.  
  1340.    lastdeep = i + deepeningAmount;
  1341.  
  1342.    /* start report of what's happening */
  1343.    printf(", depth %ld, deepening %d, ", (long int) i, deepeningAmount);
  1344.    putnum(qTail - qHead);
  1345.    printf("/");
  1346.    putnum(qTail);
  1347.    fflush(stdout);
  1348.  
  1349.  
  1350.    /* go through queue, deepening each one */
  1351.    
  1352.    #pragma omp parallel
  1353.    {
  1354.    node j;
  1355.    uint32_t *pInd;
  1356.    uint32_t *pRemain;
  1357.    row *pRows;
  1358.    
  1359.    
  1360.    pInd = malloc((long long)sizeof(*pInd) * (deepeningAmount + 4 * params[P_PERIOD]));
  1361.    pRemain = malloc((long long)sizeof(*pRemain) * (deepeningAmount + 4 * params[P_PERIOD]));
  1362.    pRows = malloc((long long)sizeof(*pRows) * (deepeningAmount + 4 * params[P_PERIOD]));
  1363.    
  1364.    #pragma omp for schedule(dynamic, CHUNK_SIZE)
  1365.    for (j = qHead; j < qTail; j++) {
  1366.       if (!EMPTY(j) && !depthFirst(j, deepeningAmount, pInd, pRemain, pRows))
  1367.          MAKEEMPTY(j);
  1368.    }
  1369.    free(pInd);
  1370.    free(pRemain);
  1371.    free(pRows);
  1372.    }
  1373.    
  1374.    if (deepeningAmount > period) deepeningAmount--; /* allow for gradual depth reduction */
  1375.    
  1376.    /* before reporting new queue size, shrink tree back down */
  1377.    printf(" -> ");
  1378.    fflush(stdout);
  1379. //   hash = savedHash;
  1380.    //perdidor = 0;
  1381.  
  1382.    /* signal time for dump */
  1383.    if (params[P_CHECKPOINT]) dumpFlag = DUMPPENDING;
  1384.  
  1385.    doCompact();
  1386.  
  1387.    /* now finish report */
  1388.    putnum(qTail - qHead);
  1389.    printf("/");
  1390.    putnum(qTail);
  1391.    printf("\n");
  1392.  
  1393.    /* Report successful/unsuccessful dump */
  1394.    if (dumpFlag == DUMPSUCCESS)
  1395.    {
  1396. #ifdef PATCH07
  1397.       timeStamp();
  1398. #endif
  1399.        printf("State dumped to %s\n",dumpFile);
  1400.        /*analyse();
  1401.        if (chainWidth)
  1402.            printf("[%d/%d]\n",chainDepth,chainWidth+1);
  1403.        else
  1404.            printf("[%d/-]\n",chainDepth);*/
  1405.    }
  1406.    else if (dumpFlag == DUMPFAILURE)
  1407.    {
  1408. #ifdef PATCH07
  1409.       timeStamp();
  1410. #endif
  1411.       printf("State dump unsuccessful\n");
  1412.    }
  1413.  
  1414.    
  1415.    fflush(stdout);
  1416. }
  1417.  
  1418.  
  1419. static void breadthFirst()
  1420. {
  1421.    while (!aborting && !qIsEmpty()) {
  1422.       if (qTail - qHead >= (1<<params[P_DEPTHLIMIT]) || qTail >= QSIZE - QSIZE/16 ||
  1423.           qTail >= QSIZE - (deepeningAmount << 2)) deepen();
  1424.       else process(dequeue());
  1425.    }
  1426. }
  1427.  
  1428.  
  1429. int gcd(int a, int b) {
  1430.    if (a > b) return gcd(b,a);
  1431.    else if (a == 0) return b;
  1432.    else return gcd(b-a,a);
  1433. }
  1434.  
  1435. void echoParams(){
  1436.    printf("Rule: ");
  1437.    printRule();
  1438.    printf("\n");
  1439.    printf("Period: %d\n",params[P_PERIOD]);
  1440.    printf("Offset: %d\n",params[P_OFFSET]);
  1441.    printf("Width:  %d\n",params[P_WIDTH]);
  1442.    if(params[P_SYMMETRY] == SYM_ASYM) printf("Symmetry: asymmetric\n");
  1443.    else if(params[P_SYMMETRY] == SYM_ODD) printf("Symmetry: odd\n");
  1444.    else if(params[P_SYMMETRY] == SYM_EVEN) printf("Symmetry: even\n");
  1445.    else if(params[P_SYMMETRY] == SYM_GUTTER) printf("Symmetry: gutter\n");
  1446.    if(params[P_CHECKPOINT]) printf("Dump state after queue compaction\n");
  1447.    if(!params[P_REORDER]) printf("Use naive search order.\n");
  1448.    printf("Queue size: 2^%d\n",params[P_QBITS]);
  1449.    printf("Hash table size: 2^%d\n",params[P_HASHBITS]);
  1450.    printf("Minimum deepening increment: %d\n",MINDEEP);
  1451.    printf("Number of threads: %d\n",params[P_NUMTHREADS]);
  1452. }
  1453.  
  1454.  
  1455. void usage(){
  1456.    printf("%s\n",BANNER);
  1457.    printf("\n");
  1458.    printf("Usage: \"qfind options\"\n");
  1459.    printf("  e.g. \"qfind B3/S23 p3 k1 w6 v\" searches Life (rule B3/S23) for\n");
  1460.    printf("  c/3 orthogonal spaceships with even bilateral symmetry and a\n");
  1461.    printf("  search width of 6 (full width 12).\n");
  1462.    printf("\n");
  1463.    printf("Available options:\n");
  1464.    printf("  bNN/sNN searches for spaceships in the specified rule (default: b3/s23)\n");
  1465.    printf("\n");
  1466.    printf("  pNN  searches for spaceships with period NN\n");
  1467.    printf("  kNN  searches for spaceships that travel NN cells every period\n");
  1468.    printf("  wNN  searches for spaceships with search width NN\n");
  1469.    printf("       (full width depends on symmetry type)\n");
  1470.    printf("\n");
  1471.    printf("  tNN  runs search using NN threads during deepening step (default: 1)\n");
  1472.    printf("  hNN  sets the hash table size to 2^NN (default: %d)\n",HASHBITS);
  1473.    printf("       Use h0 to disable duplicate elimination.\n");
  1474.    printf("  qNN  sets the BFS queue size to 2^NN (default: %d)\n",QBITS);
  1475.    printf("  iNN  groups 2^NN queue entries to an index node (default: 4)\n");
  1476.    printf("\n");
  1477.    printf("  mNN  sets minimum deepening increment to NN (default: period)\n");
  1478.    //printf("  lNN  terminates the search if it reaches a depth of NN (default: %d)\n",DEFAULT_DEPTH_LIMIT);
  1479.    //printf("  mNN  disallows spaceships longer than a depth of NN\n");
  1480.    //printf("       (the spaceship length is approximately depth/period)\n");
  1481.    //printf("  fNN  disallows spaceships that do not have the full period by a depth of NN\n");
  1482.    //printf("  tNN  disallows full-period rows of width greater than NN\n");
  1483.    //printf("  sNN  terminates the search if NN spaceships are found (default: 1)\n");
  1484.    printf("\n");
  1485.    printf("  d    dumps the search state after each queue compaction\n");
  1486.    //printf("  j    dumps the state at start of search\n");
  1487.    printf("\n");
  1488.    printf("  a    searches for asymmetric spaceships\n");
  1489.    printf("  u    searches for odd bilaterally symmetric spaceships\n");
  1490.    printf("  v    searches for even bilaterally symmetric spaceships\n");
  1491.    printf("  g    searches for symmetric spaceships with gutters (empty center column)\n");
  1492.    printf("\n");
  1493.    printf("  o    uses naive search order (not recommended)\n");
  1494.    printf("\n");
  1495.    printf("  e FF uses rows in the file FF as the initial rows for the search\n");
  1496.    printf("       (use the companion Golly python script to easily generate the\n");
  1497.    printf("       initial row file)\n");
  1498.    printf("\n");
  1499.    printf("\"qfind command file\" reloads the state from the specified file\n");
  1500.    printf("and performs the command. Available commands: \n");
  1501.    printf("  s    resumes search from the loaded state\n");
  1502.    printf("  p    previews partial results\n");
  1503. }
  1504.  
  1505. void loadInitRows(char * file){
  1506.    FILE * fp;
  1507.    int i,j;
  1508.    char rowStr[MAXWIDTH];
  1509.    row theRow = 0;
  1510.    
  1511.    loadFile = file;
  1512.    fp = fopen(loadFile, "r");
  1513.    if (!fp) loadFail();
  1514.    
  1515.    printf("Starting search from rows in %s:\n",loadFile);
  1516.    
  1517.    for(i = 0; i < 2 * period; i++){
  1518.       fscanf(fp,"%s",rowStr);
  1519.       for(j = 0; j < width; j++){
  1520.          theRow |= ((rowStr[width - j - 1] == '.') ? 0:1) << j;
  1521.       }
  1522.       printRow(theRow);
  1523.       enqueue(dequeue(),theRow);
  1524.       theRow = 0;
  1525.    }
  1526.    fclose(fp);
  1527. }
  1528.  
  1529. int main(int argc, char *argv[]){
  1530.    printf("%s\n",BANNER);
  1531.    
  1532.    params[P_RULE] = 6152;         //first 9 bits represent births; next 9 bits represent survivals
  1533.    params[P_WIDTH] = 0;
  1534.    params[P_PERIOD] = 0;
  1535.    params[P_OFFSET] = 0;
  1536.    params[P_SYMMETRY] = 0;
  1537.    params[P_INIT_ROWS] = 0;
  1538.    //params[P_FULL_PERIOD] = 0;
  1539.    params[P_NUM_SHIPS] = 1;
  1540.    //params[P_FULL_WIDTH] = 0;
  1541.    params[P_REORDER] = 1;
  1542.    params[P_BASEBITS] = 4;
  1543.    params[P_QBITS] = QBITS;
  1544.    params[P_HASHBITS] = HASHBITS;
  1545.    params[P_NUMTHREADS] = 1;
  1546.    params[P_INIT_ROWS] = 0;
  1547.    params[P_MINDEEP] = 0;
  1548.  
  1549.    int loadDumpFlag = 0;
  1550.  
  1551.  
  1552.    //int dumpandexit = 0;
  1553.    int skipNext = 0;
  1554.    //int div1,div2;
  1555.    int s;
  1556.    if(argc == 2 && !strcmp(argv[1],"c")){
  1557.       usage();
  1558.       return 0;
  1559.    }
  1560.    if(argc == 3 && (!strcmp(argv[1],"s") || !strcmp(argv[1],"S") || !strcmp(argv[1],"p") || !strcmp(argv[1],"P"))) loadDumpFlag = 1;
  1561.    else{
  1562.       for(s = 1; s < argc; s++){    //read input parameters
  1563.          if(skipNext){
  1564.             skipNext = 0;
  1565.             continue;
  1566.          }
  1567.          switch(argv[s][0]){
  1568.             case 'b': case 'B':     //read rule
  1569.                params[P_RULE] = 0;
  1570.                int sshift = 0;
  1571.                int i;
  1572.                for(i = 1; i < 100; i++){
  1573.                   int rnum = argv[s][i];
  1574.                   if(!rnum)break;
  1575.                   if(rnum == 's' || rnum == 'S')sshift = 9;
  1576.                   if(rnum >= '0' && rnum <= '8')params[P_RULE] += 1 << (sshift + rnum - '0');
  1577.                }
  1578.             break;
  1579.             case 'w': case 'W': sscanf(&argv[s][1], "%d", &params[P_WIDTH]); break;
  1580.             case 'p': case 'P': sscanf(&argv[s][1], "%d", &params[P_PERIOD]); break;
  1581.             case 'k': case 'K': sscanf(&argv[s][1], "%d", &params[P_OFFSET]); break;
  1582.             case 'u': case 'U': params[P_SYMMETRY] = SYM_ODD; mode = odd; break;
  1583.             case 'v': case 'V': params[P_SYMMETRY] = SYM_EVEN; mode = even; break;
  1584.             case 'a': case 'A': params[P_SYMMETRY] = SYM_ASYM; mode = asymmetric; break;
  1585.             case 'g': case 'G': params[P_SYMMETRY] = SYM_GUTTER; mode = gutter; break;
  1586.             case 'd': case 'D': params[P_CHECKPOINT] = 1; break;
  1587.             //case 'j': case 'J': dumpandexit = 1; break;
  1588.             case 'e': case 'E': params[P_INIT_ROWS] = s + 1; skipNext = 1; break;
  1589.             case 'm': case 'M': sscanf(&argv[s][1], "%d", &params[P_MINDEEP]); break;
  1590.             //case 'f': case 'F': sscanf(&argv[s][1], "%d", &params[P_FULL_PERIOD]); break;
  1591.             case 's': case 'S': sscanf(&argv[s][1], "%d", &params[P_NUM_SHIPS]); break;
  1592.             case 't': case 'T': sscanf(&argv[s][1], "%d", &params[P_NUMTHREADS]); break;
  1593.             case 'o': case 'O': params[P_REORDER] = 0; break;
  1594.             case 'q': case 'Q': sscanf(&argv[s][1], "%d", &params[P_QBITS]); break;
  1595.             case 'h': case 'H': sscanf(&argv[s][1], "%d", &params[P_HASHBITS]); break;
  1596.             case 'i': case 'I': sscanf(&argv[s][1], "%d", &params[P_BASEBITS]); break;
  1597.          }
  1598.       }
  1599.    }
  1600.  
  1601.    
  1602.    
  1603.    
  1604.    if(loadDumpFlag) loadState(argv[1],argv[2]);
  1605.    else{
  1606.      
  1607.       width = params[P_WIDTH];
  1608.       period = params[P_PERIOD];
  1609.       offset = params[P_OFFSET];
  1610.       deepeningAmount = period;
  1611.       lastdeep = 0;
  1612.       hashPhase = (gcd(period,offset)>1);
  1613.    
  1614.       rule = params[P_RULE];
  1615.    
  1616.       nRowsInState = period+period;
  1617.    
  1618.       params[P_DEPTHLIMIT] = DEFAULT_DEPTHLIMIT;
  1619.    
  1620.       //if(params[P_SYMMETRY] == SYM_ASYM){
  1621.       //   mode = asymmetric;
  1622.       //}
  1623.    
  1624.       base = malloc((QSIZE>>BASEBITS)*sizeof(node));
  1625.       rows = malloc(QSIZE*sizeof(row));
  1626.       if (base == 0 || rows == 0) {
  1627.          printf("Unable to allocate BFS queue!\n");
  1628.          exit(0);
  1629.       }
  1630.    
  1631.       if (hashBits == 0) hash = 0;
  1632.       else {
  1633.          hash = malloc(HASHSIZE*sizeof(node));
  1634.          if (hash == 0) printf("Unable to allocate hash table, duplicate elimination disabled\n");
  1635.       }
  1636.      
  1637.       omp_set_num_threads(params[P_NUMTHREADS]);
  1638.      
  1639.       enqueue(0,0);
  1640.      
  1641.       if(params[P_INIT_ROWS]) loadInitRows(argv[params[P_INIT_ROWS]]);
  1642.    }
  1643.    
  1644.    echoParams();
  1645.    
  1646.    makePhases();
  1647.    makeTables();
  1648.  
  1649.    //enqueue(0,0);
  1650.    rephase();
  1651.  
  1652.  
  1653.    /*printRule();
  1654.    printf("\nperiod: %d\n",period);
  1655.    printf("width: %d\n",width);
  1656.    printf("offset: %d\n",offset);
  1657.    printf("Symmetry: %d\n\n",params[P_SYMMETRY]);*/
  1658.    //printf("Hash bits: %d\n\n",hashBits);
  1659.  
  1660.  
  1661.    printf("Starting search\n");
  1662.    
  1663.    breadthFirst();
  1664.  
  1665.    printf("Search complete.\n");
  1666.  
  1667.    return 0;
  1668. }
Add Comment
Please, Sign In to add comment
Public Pastes
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!