// For each cap C (one per orbit, from capsK.txt) and each t = 1..tmax, compute
// minExt(C,t) = min over t-point sets X disjoint from C of #lines(C u X), by branch and bound.
// Usage: ./ext capsfile tmax [cutoff]   (cutoff: stop refining below; default = search exact min)
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <vector>
#include <algorithm>
typedef unsigned __int128 u128;
static int THIRD[81][81];
static inline int enc(int a, int b, int c, int d) { return a * 27 + b * 9 + c * 3 + d; }
static inline bool bit(u128 m, int i) { return (m >> i) & 1; }
struct Search {
    u128 S; int pts[81]; int n;
    int best; int t;
    int cand[81]; int nc;          // candidate outside points, sorted by secant count w.r.t. C
    int secC[81];
    // marginal secants of x w.r.t. current set S
    int marg(int x) { int c = 0; for (int i = 0; i < n; i++) if (bit(S, THIRD[x][pts[i]])) c++; return c / 2; }
    void rec(int idx, int added, int lines) {
        if (lines >= best) return;
        if (added == t) { best = lines; return; }
        int need = t - added;
        // bound: lines + sum of the 'need' smallest secC among remaining candidates (secC <= marginal)
        // candidates are sorted by secC so the smallest remaining are cand[idx..idx+need)
        if (nc - idx < need) return;
        int lb = lines; for (int k = 0; k < need; k++) lb += secC[cand[idx + k]];
        if (lb >= best) return;
        for (int i = idx; i <= nc - need; i++) {
            int x = cand[i]; int m = marg(x);
            if (lines + m >= best) continue;
            S |= (u128)1 << x; pts[n++] = x;
            rec(i + 1, added + 1, lines + m);
            n--; S &= ~((u128)1 << x);
        }
    }
    int run(u128 C, int tt, int cutoff) {
        S = C; n = 0; for (int i = 0; i < 81; i++) if (bit(C, i)) pts[n++] = i;
        t = tt; nc = 0;
        for (int x = 0; x < 81; x++) if (!bit(C, x)) { secC[x] = marg(x); cand[nc++] = x; }
        std::sort(cand, cand + nc, [&](int a, int b) { return secC[a] < secC[b]; });
        best = cutoff;  // we only care whether min < cutoff; report exact min if below
        rec(0, 0, 0);
        return best;
    }
};
int main(int argc, char** argv) {
    for (int x = 0; x < 81; x++) for (int y = 0; y < 81; y++) {
        int xc[4] = {x / 27, (x / 9) % 3, (x / 3) % 3, x % 3}, yc[4] = {y / 27, (y / 9) % 3, (y / 3) % 3, y % 3};
        THIRD[x][y] = enc((6 - xc[0] - yc[0]) % 3, (6 - xc[1] - yc[1]) % 3, (6 - xc[2] - yc[2]) % 3, (6 - xc[3] - yc[3]) % 3);
    }
    FILE* f = fopen(argv[1], "r"); int tmax = atoi(argv[2]); int cutoff = argc > 3 ? atoi(argv[3]) : 1 << 20;
    std::vector<u128> caps; char line[4096];
    while (fgets(line, sizeof line, f)) { u128 m = 0; char* p = strtok(line, " \n"); while (p) { m |= (u128)1 << atoi(p); p = strtok(nullptr, " \n"); } if (m) caps.push_back(m); }
    int k = 0; for (int i = 0; i < 81; i++) if (bit(caps[0], i)) k++;
    printf("%zu caps of size %d from %s\n", caps.size(), k, argv[1]);
    Search Sch;
    std::vector<int> overall(tmax + 1, 1 << 20);
    for (size_t ci = 0; ci < caps.size(); ci++) {
        // secant distribution
        int hist[32] = {0}; Sch.S = caps[ci]; Sch.n = 0; for (int i = 0; i < 81; i++) if (bit(caps[ci], i)) Sch.pts[Sch.n++] = i;
        for (int x = 0; x < 81; x++) if (!bit(caps[ci], x)) hist[Sch.marg(x)]++;
        printf("cap %zu secant-dist:", ci); for (int s = 0; s < 32; s++) if (hist[s]) printf(" %d:%d", s, hist[s]);
        printf(" | minExt(t):");
        for (int t = 1; t <= tmax; t++) { int r = Sch.run(caps[ci], t, cutoff); overall[t] = std::min(overall[t], r); printf(" t=%d:%s%d", t, r >= cutoff ? ">=" : "", r); }
        printf("\n"); fflush(stdout);
    }
    printf("OVERALL min over caps (>= means cutoff reached):");
    for (int t = 1; t <= tmax; t++) printf(" t=%d:%s%d", t, overall[t] >= cutoff ? ">=" : "", overall[t]);
    printf("\n");
}