Started work on k=2 tables; further split h48 gendata

4 files changed, 283 insertions, 208 deletions

diff --git a/src/nissy.c b/src/nissy.c
index 9df775c..fac7a20 100644
--- a/src/nissy.c
+++ b/src/nissy.c
@@ -164,22 +164,23 @@ nissy_gendata(
 )
 {
         int64_t ret;
-        uint8_t maxdepth, h, i, j;
+        uint8_t i;
+        gendata_h48_arg_t arg;
 
+        arg.buf = data;
         if (!strcmp(solver, "h48")) {
-                /*  options are in the form "h;maxdepth" */
+                /*  options are in the form "h;k;maxdepth" */
+                arg.h = atoi(options);
                 for (i = 0; options[i] != ';'; i++) ;
-                for (j = i; options[j]; j++) ;
-                h = atoi(options);
-                if (h != 0) {
-                        LOG("Temporarily only h=0 is supported\n");
-                        ret = -1;
-                } else {
-                        maxdepth = atoi(&options[i+1]);
-                        ret = gendata_h48h0k4(data, maxdepth);
-                }
+                arg.k = atoi(&options[i+1]);
+                for (i = i+1; options[i] != ';'; i++) ;
+                arg.maxdepth = atoi(&options[i+1]);
+                ret = gendata_h48(&arg);
         } else if (!strcmp(solver, "h48stats")) {
-                ret = gendata_h48h0k4(data, 20);
+                arg.h = 0;
+                arg.k = 4;
+                arg.maxdepth = 20;
+                ret = gendata_h48(&arg);
         } else {
                 LOG("gendata: implemented only for h48 solver\n");
                 ret = -1;
diff --git a/src/solvers/h48/gendata_cocsep.h b/src/solvers/h48/gendata_cocsep.h
new file mode 100644
index 0000000..13ce68b
--- /dev/null
+++ b/src/solvers/h48/gendata_cocsep.h
@@ -0,0 +1,169 @@
+#define COCSEP_CLASSES        ((size_t)3393)
+#define COCSEP_TABLESIZE      ((size_t)_3p7 << (size_t)7)
+#define COCSEP_VISITEDSIZE    ((COCSEP_TABLESIZE + (size_t)7) / (size_t)8)
+#define COCSEP_FULLSIZE       ((size_t)4 * (COCSEP_TABLESIZE + (size_t)12))
+
+#define VISITED_IND(i)        ((uint32_t)(i) / UINT32_C(8))
+#define VISITED_MASK(i)       (UINT32_C(1) << ((uint32_t)(i) % UINT32_C(8)))
+
+#define CBOUND_MASK           UINT32_C(0xFF)
+#define CBOUND(x)             ((x) & CBOUND_MASK)
+
+typedef struct {
+        cube_t cube;
+        uint8_t depth;
+        uint8_t maxdepth;
+        uint16_t *n;
+        uint32_t *buf32;
+        uint8_t *visited;
+        uint64_t *selfsim;
+        cube_t *rep;
+} dfsarg_cocsep_t;
+
+_static_inline bool get_visited(const uint8_t *, int64_t);
+_static_inline void set_visited(uint8_t *, int64_t);
+
+_static size_t gendata_cocsep(void *, uint64_t *, cube_t *);
+_static uint32_t gendata_cocsep_dfs(dfsarg_cocsep_t *);
+
+_static_inline int8_t get_h48_cdata(cube_t, uint32_t *, uint32_t *);
+
+/*
+Each element of the cocsep table is a uint32_t used as follows:
+  - Lowest 8-bit block: pruning value
+  - Second-lowest 8-bit block: "ttrep" (transformation to representative)
+  - Top 16-bit block: symcoord value
+After the data as described above, more auxiliary information is appended:
+  - A uint32_t representing the number of symmetry classes
+  - A uint32_t representing the highest value of the pruning table
+  - One uint32_t for each "line" of the pruning table, representing the number
+    of positions having that pruning value.
+*/
+_static size_t
+gendata_cocsep(void *buf, uint64_t *selfsim, cube_t *rep)
+{
+        uint32_t *buf32, *info, cc;
+        uint16_t n;
+        uint8_t i, j, visited[COCSEP_VISITEDSIZE];
+        dfsarg_cocsep_t arg;
+
+        if (buf == NULL)
+                goto gendata_cocsep_return_size;
+
+        buf32 = (uint32_t *)buf;
+        info = buf32 + COCSEP_TABLESIZE;
+        memset(buf32, 0xFF, sizeof(uint32_t) * COCSEP_TABLESIZE);
+        if (selfsim != NULL)
+                memset(selfsim, 0, sizeof(uint64_t) * COCSEP_CLASSES);
+
+        arg = (dfsarg_cocsep_t) {
+                .cube = solved,
+                .n = &n,
+                .buf32 = buf32,
+                .visited = visited,
+                .selfsim = selfsim,
+                .rep = rep
+        };
+        for (i = 0, n = 0, cc = 0; i < 10; i++) {
+                LOG("cocsep: generating depth %" PRIu8 "\n", i);
+                memset(visited, 0, COCSEP_VISITEDSIZE);
+                arg.depth = 0;
+                arg.maxdepth = i;
+                cc = gendata_cocsep_dfs(&arg);
+                info[i+2] = cc;
+                LOG("found %" PRIu32 "\n", cc);
+        }
+
+        info[0] = (uint32_t)n;
+        info[1] = 9; /* Known max pruning value */
+        DBG_ASSERT(n == COCSEP_CLASSES, 0,
+            "cocsep: computed %" PRIu16 " symmetry classes, "
+            "expected %zu\n", n, COCSEP_CLASSES);
+
+        LOG("cocsep data computed\n");
+        LOG("Symmetry classes: %" PRIu32 "\n", info[0]);
+        LOG("Maximum pruning value: %" PRIu32 "\n", info[1]);
+        LOG("Pruning value distribution:\n");
+        for (j = 0; j < 10; j++)
+                LOG("%" PRIu8 ":\t%" PRIu32 "\n", j, info[j+2]);
+
+gendata_cocsep_return_size:
+        return COCSEP_FULLSIZE;
+}
+
+_static uint32_t
+gendata_cocsep_dfs(dfsarg_cocsep_t *arg)
+{
+        uint8_t m;
+        uint32_t cc, class, ttrep, depth, olddepth, tinv;
+        uint64_t t;
+        int64_t i, j;
+        cube_t d;
+        dfsarg_cocsep_t nextarg;
+
+        i = coord_cocsep(arg->cube);
+        olddepth = (uint8_t)(arg->buf32[i] & 0xFF);
+        if (olddepth < arg->depth || get_visited(arg->visited, i))
+                return 0;
+        set_visited(arg->visited, i);
+
+        if (arg->depth == arg->maxdepth) {
+                if ((arg->buf32[i] & 0xFF) != 0xFF)
+                        return 0;
+
+                if (arg->rep != NULL)
+                        arg->rep[*arg->n] = arg->cube;
+                for (t = 0, cc = 0; t < 48; t++) {
+                        d = transform_corners(arg->cube, t);
+                        j = coord_cocsep(d);
+                        if (i == j && arg->selfsim != NULL)
+                                arg->selfsim[*arg->n] |= UINT64_C(1) << t;
+                        if (COCLASS(arg->buf32[j]) != UINT32_C(0xFFFF))
+                                continue;
+                        set_visited(arg->visited, j);
+                        tinv = inverse_trans(t);
+                        olddepth = arg->buf32[j] & 0xFF;
+                        cc += olddepth == 0xFF;
+
+                        class = (uint32_t)(*arg->n) << UINT32_C(16);
+                        ttrep = (uint32_t)tinv << UINT32_C(8);
+                        depth = (uint32_t)arg->depth;
+                        arg->buf32[j] = class | ttrep | depth;
+                }
+                (*arg->n)++;
+
+                return cc;
+        }
+
+        memcpy(&nextarg, arg, sizeof(dfsarg_cocsep_t));
+        nextarg.depth++;
+        for (m = 0, cc = 0; m < 18; m++) {
+                nextarg.cube = move(arg->cube, m);
+                cc += gendata_cocsep_dfs(&nextarg);
+        }
+
+        return cc;
+}
+
+_static_inline bool
+get_visited(const uint8_t *a, int64_t i)
+{
+        return a[VISITED_IND(i)] & VISITED_MASK(i);
+}
+
+_static_inline void
+set_visited(uint8_t *a, int64_t i)
+{
+        a[VISITED_IND(i)] |= VISITED_MASK(i);
+}
+
+_static_inline int8_t
+get_h48_cdata(cube_t cube, uint32_t *cocsepdata, uint32_t *cdata)
+{
+        int64_t coord;
+
+        coord = coord_cocsep(cube);
+        *cdata = cocsepdata[coord];
+
+        return CBOUND(*cdata);
+}
diff --git a/src/solvers/h48/gendata.h b/src/solvers/h48/gendata_full.h
index d6bf85b..be234d9 100644
--- a/src/solvers/h48/gendata.h
+++ b/src/solvers/h48/gendata_full.h
@@ -1,8 +1,3 @@
-#define COCSEP_CLASSES        ((size_t)3393)
-#define COCSEP_TABLESIZE      ((size_t)_3p7 << (size_t)7)
-#define COCSEP_VISITEDSIZE    ((COCSEP_TABLESIZE + (size_t)7) / (size_t)8)
-#define COCSEP_FULLSIZE       ((size_t)4 * (COCSEP_TABLESIZE + (size_t)12))
-
 #define ESEP_NOEO             (COCSEP_CLASSES * (size_t)_12c4 * (size_t)_8c4)
 #define ESEP_MAX(h)           (ESEP_NOEO << (size_t)(h))
 #define ESEP_TABLESIZE(h, k)  (ESEP_MAX((h)) / ((size_t)8 / (size_t)(k)))
@@ -10,11 +5,6 @@
 #define ESEP_IND(i)           ((uint32_t)(i) / UINT32_C(8))
 #define ESEP_SHIFT(i)         (UINT32_C(4) * ((uint32_t)(i) % UINT32_C(8)))
 #define ESEP_MASK(i)          ((_bit_u32(4) - (uint32_t)(1)) << ESEP_SHIFT(i))
-#define VISITED_IND(i)        ((uint32_t)(i) / UINT32_C(8))
-#define VISITED_MASK(i)       (UINT32_C(1) << ((uint32_t)(i) % UINT32_C(8)))
-
-#define CBOUND_MASK           UINT32_C(0xFF)
-#define CBOUND(x)             ((x) & CBOUND_MASK)
 
 #define MAXLEN 20
 
@@ -38,25 +28,19 @@ _t by _ttrep).
         }
 
 typedef struct {
-        cube_t cube;
-        uint8_t depth;
+        uint8_t h;
+        uint8_t k;
         uint8_t maxdepth;
-        uint16_t *n;
-        uint32_t *buf32;
-        uint8_t *visited;
-        uint64_t *selfsim;
-        cube_t *rep;
-} dfsarg_cocsep_t;
+        void * buf;
+} gendata_h48_arg_t;
 
-/* TODO keep or not? */
 typedef struct {
-        cube_t cube;
-        int8_t nmoves;
-        int8_t depth;
-        uint8_t moves[MAXLEN];
-        uint32_t *cocsepdata;
-        h48map_t *visited;
-} dfsarg_genh48set_t;
+        uint8_t maxdepth;
+        const uint32_t *cocsepdata;
+        const cube_t *crep;
+        const uint64_t *selfsim;
+        h48map_t *map;
+} gendata_h48short_arg_t;
 
 typedef struct {
         uint8_t depth;
@@ -67,148 +51,22 @@ typedef struct {
         cube_t *crep;
 } bfsarg_esep_t;
 
-_static_inline bool get_visited(const uint8_t *, int64_t);
-_static_inline void set_visited(uint8_t *, int64_t);
 _static_inline uint8_t get_esep_pval(const uint32_t *, int64_t);
 _static_inline void set_esep_pval(uint32_t *, int64_t, uint8_t);
 
-_static size_t gendata_cocsep(void *, uint64_t *, cube_t *);
-_static uint32_t gendata_cocsep_dfs(dfsarg_cocsep_t *);
-_static uint64_t gen_h48short(
-    uint8_t, const uint32_t *, const cube_t *, const uint64_t *, h48map_t *);
+_static uint64_t gen_h48short(gendata_h48short_arg_t *);
+_static size_t gendata_h48(gendata_h48_arg_t *);
 _static size_t gendata_h48h0k4(void *, uint8_t);
 _static int64_t gendata_h48h0k4_bfs(bfsarg_esep_t *);
 _static int64_t gendata_h48h0k4_bfs_fromdone(bfsarg_esep_t *);
 _static int64_t gendata_h48h0k4_bfs_fromnew(bfsarg_esep_t *);
+_static size_t gendata_h48k2(void *, uint8_t, uint8_t);
 
-_static_inline int8_t get_h48_cdata(cube_t, uint32_t *, uint32_t *);
 _static_inline int8_t get_h48_bound(cube_t, uint32_t, uint8_t, uint32_t *);
 
-/*
-Each element of the cocsep table is a uint32_t used as follows:
-  - Lowest 8-bit block: pruning value
-  - Second-lowest 8-bit block: "ttrep" (transformation to representative)
-  - Top 16-bit block: symcoord value
-After the data as described above, more auxiliary information is appended:
-  - A uint32_t representing the number of symmetry classes
-  - A uint32_t representing the highest value of the pruning table
-  - One uint32_t for each "line" of the pruning table, representing the number
-    of positions having that pruning value.
-*/
-_static size_t
-gendata_cocsep(void *buf, uint64_t *selfsim, cube_t *rep)
-{
-        uint32_t *buf32, *info, cc;
-        uint16_t n;
-        uint8_t i, j, visited[COCSEP_VISITEDSIZE];
-        dfsarg_cocsep_t arg;
-
-        if (buf == NULL)
-                goto gendata_cocsep_return_size;
-
-        buf32 = (uint32_t *)buf;
-        info = buf32 + COCSEP_TABLESIZE;
-        memset(buf32, 0xFF, sizeof(uint32_t) * COCSEP_TABLESIZE);
-        if (selfsim != NULL)
-                memset(selfsim, 0, sizeof(uint64_t) * COCSEP_CLASSES);
-
-        arg = (dfsarg_cocsep_t) {
-                .cube = solved,
-                .n = &n,
-                .buf32 = buf32,
-                .visited = visited,
-                .selfsim = selfsim,
-                .rep = rep
-        };
-        for (i = 0, n = 0, cc = 0; i < 10; i++) {
-                LOG("cocsep: generating depth %" PRIu8 "\n", i);
-                memset(visited, 0, COCSEP_VISITEDSIZE);
-                arg.depth = 0;
-                arg.maxdepth = i;
-                cc = gendata_cocsep_dfs(&arg);
-                info[i+2] = cc;
-                LOG("found %" PRIu32 "\n", cc);
-        }
-
-        info[0] = (uint32_t)n;
-        info[1] = 9; /* Known max pruning value */
-        DBG_ASSERT(n == COCSEP_CLASSES, 0,
-            "cocsep: computed %" PRIu16 " symmetry classes, "
-            "expected %zu\n", n, COCSEP_CLASSES);
-
-        LOG("cocsep data computed\n");
-        LOG("Symmetry classes: %" PRIu32 "\n", info[0]);
-        LOG("Maximum pruning value: %" PRIu32 "\n", info[1]);
-        LOG("Pruning value distribution:\n");
-        for (j = 0; j < 10; j++)
-                LOG("%" PRIu8 ":\t%" PRIu32 "\n", j, info[j+2]);
-
-gendata_cocsep_return_size:
-        return COCSEP_FULLSIZE;
-}
-
-_static uint32_t
-gendata_cocsep_dfs(dfsarg_cocsep_t *arg)
-{
-        uint8_t m;
-        uint32_t cc, class, ttrep, depth, olddepth, tinv;
-        uint64_t t;
-        int64_t i, j;
-        cube_t d;
-        dfsarg_cocsep_t nextarg;
-
-        i = coord_cocsep(arg->cube);
-        olddepth = (uint8_t)(arg->buf32[i] & 0xFF);
-        if (olddepth < arg->depth || get_visited(arg->visited, i))
-                return 0;
-        set_visited(arg->visited, i);
-
-        if (arg->depth == arg->maxdepth) {
-                if ((arg->buf32[i] & 0xFF) != 0xFF)
-                        return 0;
-
-                if (arg->rep != NULL)
-                        arg->rep[*arg->n] = arg->cube;
-                for (t = 0, cc = 0; t < 48; t++) {
-                        d = transform_corners(arg->cube, t);
-                        j = coord_cocsep(d);
-                        if (i == j && arg->selfsim != NULL)
-                                arg->selfsim[*arg->n] |= UINT64_C(1) << t;
-                        if (COCLASS(arg->buf32[j]) != UINT32_C(0xFFFF))
-                                continue;
-                        set_visited(arg->visited, j);
-                        tinv = inverse_trans(t);
-                        olddepth = arg->buf32[j] & 0xFF;
-                        cc += olddepth == 0xFF;
-
-                        class = (uint32_t)(*arg->n) << UINT32_C(16);
-                        ttrep = (uint32_t)tinv << UINT32_C(8);
-                        depth = (uint32_t)arg->depth;
-                        arg->buf32[j] = class | ttrep | depth;
-                }
-                (*arg->n)++;
-
-                return cc;
-        }
-
-        memcpy(&nextarg, arg, sizeof(dfsarg_cocsep_t));
-        nextarg.depth++;
-        for (m = 0, cc = 0; m < 18; m++) {
-                nextarg.cube = move(arg->cube, m);
-                cc += gendata_cocsep_dfs(&nextarg);
-        }
-
-        return cc;
-}
-
 _static uint64_t
-gen_h48short(
-        uint8_t n,
-        const uint32_t *cocsepdata,
-        const cube_t *crep,
-        const uint64_t *selfsim,
-        h48map_t *map
-) {
+gen_h48short(gendata_h48short_arg_t *arg)
+{
         uint8_t i, m;
         int64_t coord;
         uint64_t j, oldn;
@@ -216,33 +74,48 @@ gen_h48short(
         cube_t cube, d;
 
         cube = solvedcube();
-        coord = coord_h48(cube, cocsepdata, 11);
-        h48map_insertmin(map, coord, 0);
+        coord = coord_h48(cube, arg->cocsepdata, 11);
+        h48map_insertmin(arg->map, coord, 0);
         oldn = 0;
-        LOG("Short h48: generating depth 0\nfound %" PRIu8 "\n", map->n-oldn);
-        for (i = 0; i < n; i++) {
-                LOG("Short h48: generating depth %" PRIu8 "\n", i+1);
+        LOG("Short h48: depth 0\nfound %" PRIu8 "\n", arg->map->n-oldn);
+        for (i = 0; i < arg->maxdepth; i++) {
+                LOG("Short h48: depth %" PRIu8 "\n", i+1);
                 j = 0;
-                oldn = map->n;
-                for (kv = h48map_nextkvpair(map, &j);
-                     j != map->capacity;
-                     kv = h48map_nextkvpair(map, &j)
+                oldn = arg->map->n;
+                for (kv = h48map_nextkvpair(arg->map, &j);
+                     j != arg->map->capacity;
+                     kv = h48map_nextkvpair(arg->map, &j)
                 ) {
                         if (kv.val != i)
                                 continue;
-                        cube = invcoord_h48(kv.key, crep, 11);
+                        cube = invcoord_h48(kv.key, arg->crep, 11);
                         for (m = 0; m < 18; m++) {
                                 d = move(cube, m);
-                                _foreach_h48sim(d, cocsepdata, selfsim, 11,
-                                        coord = coord_h48(d, cocsepdata, 11);
-                                        h48map_insertmin(map, coord, i+1);
+                                _foreach_h48sim(
+                                        d, arg->cocsepdata, arg->selfsim, 11,
+                                        coord = coord_h48(d, arg->cocsepdata, 11);
+                                        h48map_insertmin(arg->map, coord, i+1);
                                 )
                         }
                 }
-                LOG("found %" PRIu8 "\n", map->n-oldn);
+                LOG("found %" PRIu8 "\n", arg->map->n-oldn);
         }
 
-        return map->n;
+        return arg->map->n;
+}
+
+/* Generic function that dispatches to the data generators */
+_static size_t
+gendata_h48(gendata_h48_arg_t *arg)
+{
+        if (arg->h == 0 && arg->k == 4) {
+                return gendata_h48h0k4(arg->buf, arg->maxdepth);
+        } else if (arg->k == 2) {
+                return gendata_h48k2(arg->buf, arg->h, arg->maxdepth);
+        }
+
+        LOG("Cannot generate data for h = %" PRIu8 " and k = %" PRIu8
+            " (not implemented yet)\n", arg->h, arg->k);
 }
 
 /*
@@ -259,14 +132,13 @@ gendata_h48h0k4(void *buf, uint8_t maxdepth)
         cube_t crep[COCSEP_CLASSES];
         size_t cocsepsize, infosize;
 
+        if (buf == NULL)
+                goto gendata_h48h0k4_return_size;
+
         /* TODO: move info at start of tables (all tables!) */
-        infosize = 4 * maxdepth;
         cocsepsize = gendata_cocsep(buf, selfsim, crep);
         infosize = 88;
 
-        if (buf == NULL)
-                goto gendata_h48h0k4_return_size;
-
         esep_max = (int64_t)ESEP_MAX(0);
         cocsepdata = (uint32_t *)buf;
         buf32 = cocsepdata + cocsepsize / 4;
@@ -379,16 +251,59 @@ gendata_h48h0k4_bfs_fromnew(bfsarg_esep_t *arg)
         return cc;
 }
 
-_static_inline bool
-get_visited(const uint8_t *a, int64_t i)
+_static size_t
+gendata_h48k2(void *buf, uint8_t h, uint8_t maxdepth)
 {
-        return a[VISITED_IND(i)] & VISITED_MASK(i);
-}
+        static uint64_t capacity = 10000019; /* First prime after 1e8 */
+        static uint64_t randomizer = 10000079; /* Second prime after 1e8 */
+        static uint8_t base[] = {
+                [0]  = 8,
+                [1]  = 8,
+                [2]  = 8,
+                [3]  = 8,
+                [4]  = 9,
+                [5]  = 9,
+                [6]  = 9,
+                [7]  = 9,
+                [8]  = 10,
+                [9]  = 10,
+                [10] = 10,
+                [11] = 10
+        };
 
-_static_inline void
-set_visited(uint8_t *a, int64_t i)
-{
-        a[VISITED_IND(i)] |= VISITED_MASK(i);
+        uint64_t nshort;
+        uint32_t *buf32, *info, *cocsepdata;
+        h48map_t depth8cubes;
+        gendata_h48short_arg_t shortarg;
+        uint64_t selfsim[COCSEP_CLASSES];
+        cube_t crep[COCSEP_CLASSES];
+        size_t cocsepsize, infosize;
+
+        DBG_ASSERT(base[h] == 8, 0, "Only implemented for h <= 3 (base 8)\n");
+
+        if (buf == NULL)
+                goto gendata_h48k2_return_size;
+
+        cocsepdata = (uint32_t *)buf;
+        cocsepsize = gendata_cocsep(buf, selfsim, crep);
+        infosize = 88;
+
+        h48map_create(&depth8cubes, capacity, randomizer);
+        shortarg = (gendata_h48short_arg_t) {
+                .maxdepth = 8,
+                .cocsepdata = cocsepdata,
+                .crep = crep,
+                .selfsim = selfsim,
+                .map = &depth8cubes
+        };
+
+        nshort = gen_h48short(&shortarg);
+        LOG("%" PRIu64 "\n", nshort);
+
+        h48map_destroy(&depth8cubes);
+
+gendata_h48k2_return_size:
+        return cocsepsize + ESEP_TABLESIZE(h, 2) + infosize;
 }
 
 _static_inline uint8_t
@@ -405,17 +320,6 @@ set_esep_pval(uint32_t *buf32, int64_t i, uint8_t val)
 }
 
 _static_inline int8_t
-get_h48_cdata(cube_t cube, uint32_t *cocsepdata, uint32_t *cdata)
-{
-        int64_t coord;
-
-        coord = coord_cocsep(cube);
-        *cdata = cocsepdata[coord];
-
-        return CBOUND(*cdata);
-}
-
-_static_inline int8_t
 get_h48_bound(cube_t cube, uint32_t cdata, uint8_t h, uint32_t *h48data)
 {
         int64_t coord;
diff --git a/src/solvers/h48/h48.h b/src/solvers/h48/h48.h
index d33ad63..d9e8af3 100644
--- a/src/solvers/h48/h48.h
+++ b/src/solvers/h48/h48.h
@@ -1,4 +1,5 @@
 #include "coordinate.h"
 #include "map.h"
-#include "gendata.h"
+#include "gendata_cocsep.h"
+#include "gendata_full.h"
 #include "solve.h"