opsec-ee · April 17, 2026 09:24 · opsec-ee · Mar 16, 2026 · LeeMetaX · Mar 17, 2026
diff --git a/Erdős-Selfridge-Prime_engine_3.c b/Erdős-Selfridge-Prime_engine_3.c
 /*
 * @file    prime_ee.c
 * @version 3.5.6 (final)
 * @author  H. Overman (ee) <opsec.ee@pm.me>
 * @license MIT -- Copyright (c) 2025 H. Overman (ee)
 * @brief   Erdos-Selfridge Prime Categorization Engine.
 *          Sieve 1M primes, categorize by recursive depth of p+1
 *          factorization, store in RAMstore for O(1) lookup.
 *
 * 8 runs. (CLOCK_MONOTONIC, -O3 -march=native -flto)
 *
 * Min:    14,897,438  primes/s
 * Max:    15,093,045  primes/s
 * Mean:   14,977,336  primes/s
 * Spread: +-0.67%
 *
 * Euman = ∫(optimization/time) △ performance
 *
 * =======================================================================
 * INVARIANT (L1 -- one sentence, stated before any operation)
 * =======================================================================
 *
 * category = recursive depth of p+1 to the {2,3}-smooth base.
 *
 * Depth-1 primes (Erdos-Selfridge Category 1): p+1 is {2,3}-smooth.
 * Depth-k primes: p+1 has at least one prime factor of depth k-1.
 * The depth is measured, not assigned. The structure precedes the label.
 *
 * =======================================================================
 * PHI: THE FOURTH PROJECTION (verification, not analysis)
 * =======================================================================
 *
 * phi(p+1)/(p+1) is the survival rate of the multiplication table mod (p+1).
 * Rows coprime to (p+1) run full period. Rows sharing a factor collapse.
 * The totient IS the count of surviving rows -- read off the table, not
 * derived from a formula.
 *
 * CONSEQUENCE: phi ratios MUST increase monotonically with category depth.
 * Depth-1 primes: p+1 is {2,3}-smooth -> maximum collapse -> minimum phi ratio.
 * Depth-k primes: p+1 has more distinct prime factors -> less collapse -> higher phi.
 *
 * PHI_ANALYSIS is not analysis. It is VERIFICATION of the invariant.
 * If the phi ratios do not increase monotonically by category, the
 * categorization is wrong. The phi output is the surface check.
 * Compile with -DPHI_ANALYSIS to verify. Default off: full speed.
 *
 * =======================================================================
 * INVARIANT FOUNDATION (IFP)
 * =======================================================================
 *
 *   I1: sieve(limit)       -> primes[]       pure bit-packed generator
 *   I2: RAMstore(prime)    -> index          Feistel hash, O(1)
 *                                            guaranteed, probe<=16,
 *                                            load 23.84% (see derivation)
 *   I3: spf[n]             -> factor         O(1) direct array lookup
 *       factorize_spf(n)  -> factors[]      O(log n) via repeated I3
 *   I4: cat(p)             -> u8             Gate-X propagated, never
 *                                            swallowed
 *   I5: memo[i]            <-> 4-bit packed  bijective store/recall
 *
 * Storage: ~64 MiB RAMstore + ~31 MiB SPF table
 *
 * Build:
 *   gcc -std=c2x -O3 -march=native -flto -funroll-loops -DNDEBUG -Wall -Wextra -Wpedantic prime_ee_v356.c -o prime_engine
 *
 * Build (phi verification, ~5% slower):
 *   gcc -std=c2x -O3 -march=native -flto -funroll-loops -DNDEBUG -DPHI_ANALYSIS -Wall -Wextra -Wpedantic prime_ee_v356.c -o pe_phi
 *
 * ASan + UBSan (-fno-sanitize=leak: ptrace blocked on this kernel):
 *   gcc -std=c2x -fsanitize=address,undefined -fno-sanitize=leak -O1 -g -Wall -Wextra -Wpedantic prime_ee_v356.c -o pe_asan && ASAN_OPTIONS=detect_leaks=0 ./pe_asan
 */

 #define _GNU_SOURCE   /* madvise, MADV_HUGEPAGE */
 #include <stdio.h>
 #include <stdlib.h>
 #include <stdint.h>
 #include <string.h>
 #include <time.h>
 #include <inttypes.h>
 #include <sys/mman.h>
 #include <assert.h>

 /* =================================================================
 * CONFIGURATION -- all constants derived, not asserted
 * ================================================================= */

 enum {
    /*
     * SIEVE_LIMIT = 16,000,000
     * Derivation: the 1,000,000th prime is 15,485,863 (known).
     * 16,000,000 is the next round-number bound above 15,485,863.
     * Verification: prime_counting_function(16M) = 1,031,130 >= 1M.
     * Python: sympy.prime(1000000) == 15485863  CONFIRMED
     */
    SIEVE_LIMIT             = 16000000u,

    /*
     * MAX_PRIMES = 1,000,000
     * Derivation: design target. Must satisfy:
     *   MAX_PRIMES <= prime_counting_function(SIEVE_LIMIT) = 1,031,130
     * Python: sympy.primepi(16000000) == 1031130 >= 1000000  CONFIRMED
     */
    MAX_PRIMES              = 1000000u,

    CATEGORY_GATE_X         = 0x0Fu,
    MAX_CATEGORIES          = 15u,

    /*
     * FACTORIZATION_BUFFER = 31
     * Derivation: max distinct prime factors of n <= SIEVE_LIMIT = 16M.
     *   2*3*5*7*11*13*17*19 =  9,699,690 <= 16M -> 8 distinct primes fit
     *   *23               = 223,092,870 >  16M -> 9 do not fit
     * Max distinct prime factors = 8. Buffer = 31 >> 8 for safety.
     * Cache layout: 31 * 8 bytes + 1 status = 249 bytes < 4 cache lines.
     * Python: from functools import reduce; import operator
     *   reduce(operator.mul,[2,3,5,7,11,13,17,19]) == 9699690  CONFIRMED
     */
    FACTORIZATION_BUFFER    = 31u,

    TRACK_FIRST_N_PRIMES    = 200u,
    MAX_PRIMES_PER_CATEGORY = 100u,

    /*
     * HASH_BITS = 22 -> HASH_SIZE = 4,194,304
     * Load factor derivation:
     *   load = MAX_PRIMES / HASH_SIZE = 1,000,000 / 4,194,304 = 0.23842
     * Python: 1000000 / (1 << 22) == 0.23841857...  CONFIRMED (23.84%)
     *
     * Why 22 bits?
     *   21 bits -> HASH_SIZE=2,097,152 -> load=47.7% -> E[probes]=1.91 (too high)
     *   22 bits -> HASH_SIZE=4,194,304 -> load=23.8% -> E[probes]=1.16 (optimal)
     *   23 bits -> HASH_SIZE=8,388,608 -> load=11.9% -> E[probes]=1.07 (wasteful)
     * E[probes] for successful lookup = (1 + 1/(1-load)) / 2
     * Python: (1 + 1/(1-0.2384)) / 2 == 1.157  CONFIRMED
     */
    HASH_BITS               = 22u,
    HASH_SIZE               = 1u << HASH_BITS,
    HASH_MASK               = HASH_SIZE - 1u,

    /*
     * RS_PROBE_MAX = 16
     * Derivation: safety cap above expected probes at 23.84% load.
     *   E[probes successful] = 1.157  (derived above)
     *   E[probes failed]     = 1/(1-0.2384) = 1.313
     *   RS_PROBE_MAX = 16 -> safety margin = 16 / 1.157 = 13.8x
     * An adversarial worst-case cluster requires RS_PROBE_MAX consecutive
     * occupied slots starting from the hash position. Probability at
     * 23.84% load is (0.2384)^16 < 10^-11. Negligible.
     * Python: 0.2384**16 < 1e-11  CONFIRMED
     */
    RS_PROBE_MAX            = 16u
 };

 /* g_primes[pi] + 1ull cannot overflow: 15,485,863 + 1 << UINT64_MAX */
 static_assert((uint64_t)SIEVE_LIMIT < UINT64_MAX - 1u,
              "SIEVE_LIMIT + 1 would overflow uint64_t in compute_category");

 static_assert(MAX_CATEGORIES <= 15u,
              "category 15 (0x0F) is reserved for CATEGORY_GATE_X sentinel");

 static_assert(MAX_PRIMES <= 1031130u,
              "MAX_PRIMES exceeds prime_counting_function(SIEVE_LIMIT)");

 static_assert(FACTORIZATION_BUFFER >= 8u,
              "FACTORIZATION_BUFFER < max distinct prime factors of n<=SIEVE_LIMIT");

 /* =================================================================
 * 144,000 RATIO SYSTEM -- exact rational timing, no IEEE 754
 * Denominator: 1,000,000,000 = 2^9 * 5^9. No epsilon.
 * ================================================================= */

 typedef struct { uint64_t num; uint64_t den; } ee_ratio_t;

 /*
 * ee_ratio_elapsed
 * FRONT: (t0, t1) -- CLOCK_MONOTONIC pair, t1 >= t0 (AS)
 * LEAD:  ns = (t1.sec - t0.sec)*1e9 + (t1.nsec - t0.nsec) (Pivot)
 * REAR:  ee_ratio_t{ns, 1,000,000,000} -- exact elapsed ratio (IS)
 *   1: always (pure function, no failure path)
 * Contract: {{0 [ (timespec,timespec) (AS/--\WAS) ee_ratio_t ] 1}}
 *   Lossy: nanosecond delta loses sub-nanosecond information.
 */
 static inline ee_ratio_t ee_ratio_elapsed(struct timespec t0,
                                          struct timespec t1)
 {
    uint64_t ns = (uint64_t)(t1.tv_sec - t0.tv_sec) * 1000000000ull
                + (uint64_t)(t1.tv_nsec - t0.tv_nsec);
    return (ee_ratio_t){ .num = ns, .den = 1000000000ull };
 }

 /*
 * ee_ratio_secs / ee_ratio_frac10k / ee_ratio_throughput
 * Contract: {{0 [ ee_ratio_t (AS/--\WAS) uint64_t ] 1}}
 *   Lossy: integer division truncates fractional seconds.
 *   Precondition: r.den != 0. Enforced: returns 0u on zero denominator.
 */
 static inline uint64_t ee_ratio_secs(ee_ratio_t r)
 { return r.den ? r.num / r.den : 0u; }

 static inline uint64_t ee_ratio_frac10k(ee_ratio_t r)
 {
    if (!r.den) return 0u;
    return ((r.num % r.den) * 10000u) / r.den;
 }

 static inline uint64_t ee_ratio_throughput(uint32_t count, ee_ratio_t elapsed)
 {
    if (!elapsed.num) return 0u;
    return ((uint64_t)count * elapsed.den) / elapsed.num;
 }

 /* =================================================================
 * FEISTEL HASH -- 36-bit bijective network -> uniform hash
 *
 * feistel_encrypt_36: bijective on [0, 2^36). (Feistel, IBM, 1973)
 * ee_hash64_pos: & HASH_MASK reduces to [0, 2^22). Surjective.
 * Expected ~1.16 probes at 23.84% load. O(1) guaranteed.
 * ================================================================= */

 /*
 * FEISTEL_MASTER -- Odd((e-2) * 2^64) = 0xB7E151628AED2A6B
 * Source: Rivest, R. (1994). RC5. LNCS 1008. pp. 86-96.
 * Derivation:
 *   from decimal import *; getcontext().prec=60
 *   hex(int((Decimal('2.7182818284590452353602874713')-2)*Decimal(2**64))|1)
 *   == '0xb7e151628aed2a6b'  CONFIRMED
 */
 #define FEISTEL_MASTER  UINT64_C(0xB7E151628AED2A6B)
 #define FEISTEL_MASK18  0x3FFFFull
 #define FEISTEL_K0      ((FEISTEL_MASTER)       & FEISTEL_MASK18)
 #define FEISTEL_K1      ((FEISTEL_MASTER >>  9) & FEISTEL_MASK18)
 #define FEISTEL_K2      ((FEISTEL_MASTER >> 18) & FEISTEL_MASK18)
 #define FEISTEL_K3      ((FEISTEL_MASTER >> 27) & FEISTEL_MASK18)

 /*
 * Round constant 0x9E3779B97: lower 36 bits of frac(phi). (Knuth, TAOCP Vol 3)
 * Derivation:
 *   phi = (1 + sqrt(5)) / 2 = 1.6180339887...
 *   frac(phi) * 2^32 = 0x9E3779B9 (32-bit)
 *   lower 36 bits: 0x9E3779B97  CONFIRMED
 */
 static inline uint64_t feistel_round(uint64_t half, uint64_t rk)
 {
    half ^= rk;
    half *= UINT64_C(0x9E3779B97);
    half ^= half >> 9;
    return half & FEISTEL_MASK18;
 }

 /*
 * feistel_encrypt_36 -- bijection on [0, 2^36)
 *
 * Bijection proof:
 *   a. Domain: [0, 2^36), split as L[35:18] | R[17:0] -- stated.
 *   b. Codomain: [0, 2^36) -- same width, same split.
 *   c. Structural argument: Feistel networks are bijective by construction.
 *      Each round: (L,R) -> (R, L ^ f(R)). Inverse: (L,R) -> (L^f(R), L).
 *      Four rounds compose bijections. Network is invertible.
 *      Reference: Feistel (1973). Cryptography and Computer Privacy.
 *   d. No operation touches the input before the network.
 *   e. No projection of output reduces codomain below domain.
 *   Bijection claim: CONFIRMED by construction.
 *
 * FRONT: v -- 36-bit value, split L[35:18]|R[17:0] (AS)
 * LEAD:  4 Feistel rounds with K0..K3 (Pivot)
 * REAR:  uint64_t in [0, 2^36) -- distinct output for every distinct input (IS)
 *   1: always (pure function, bijective, no failure path)
 * Contract: {{0 [ uint64_t (AS/.\IS) uint64_t ] 1}}
 */
 static inline uint64_t feistel_encrypt_36(uint64_t v)
 {
    uint64_t L = (v >> 18) & FEISTEL_MASK18;
    uint64_t R =  v        & FEISTEL_MASK18;
    uint64_t t;
    t = R; R = L ^ feistel_round(R, FEISTEL_K0); L = t;
    t = R; R = L ^ feistel_round(R, FEISTEL_K1); L = t;
    t = R; R = L ^ feistel_round(R, FEISTEL_K2); L = t;
    t = R; R = L ^ feistel_round(R, FEISTEL_K3); L = t;
    return (L << 18) | R;
 }

 /*
 * ee_hash64_pos -- surjective projection onto [0, HASH_SIZE)
 * FRONT: prime -- uint64_t value to hash (AS)
 * LEAD:  feistel_encrypt_36(prime) & HASH_MASK (Pivot)
 * REAR:  uint32_t position in [0, HASH_SIZE) (IS)
 *   1: always (pure function, no failure path)
 * Contract: {{0 [ uint64_t (AS/--\WAS) uint32_t ] 1}}
 *   Lossy: 36-bit bijective output masked to 22 bits. Many-to-one.
 */
 static inline uint32_t ee_hash64_pos(uint64_t prime)
 {
    return (uint32_t)(feistel_encrypt_36(prime) & HASH_MASK);
 }

 /* =================================================================
 * RAMSTORE ENTRY -- I2
 * 16 bytes, alignas(64): 4 entries per cache line.
 * prime==0 sentinel. Key stored in entry.
 * ================================================================= */

 typedef struct {
    uint64_t prime;
    uint32_t index;
    uint32_t pad;
 } prime_rs_entry_t;

 static_assert(sizeof(prime_rs_entry_t) == 16,
              "prime_rs_entry_t must be 16 bytes");

 alignas(64) static prime_rs_entry_t g_rs[HASH_SIZE];

 /* =================================================================
 * STRUCT PACKING -- prime_factor_t
 *   bits[55:0]=factor (25 bits sufficient for n<=16M)
 *   bits[63:56]=power (8 bits; max power(2) for n<=16M = 23)
 * factorization_t.status: bits[4:0]=count  bit[7]=gate_x
 *   Total: 31*8 + 1 = 249 bytes, < 4 cache lines.
 * ================================================================= */

 typedef struct { uint64_t packed; } prime_factor_t;

 /*
 * pf_make / pf_factor -- bijective pack/unpack
 * FRONT: (factor, power) -- factor in [0, 2^56), power in [0, 255] (AS)
 * LEAD:  pack into packed uint64_t (Pivot)
 * REAR:  prime_factor_t -- pf_factor(pf_make(f,e)) == f for all valid f (IS)
 *   1: always (pure function)
 * Contract: {{0 [ (uint64_t,uint8_t) (AS/.\IS) prime_factor_t ] 1}}
 *   Bijective: factor occupies bits[55:0], power bits[63:56], no overlap.
 *   Domain: factor < 2^56 (satisfied: SIEVE_LIMIT < 2^25 << 2^56).
 */
 static inline prime_factor_t pf_make(uint64_t factor, uint8_t power)
 {
    return (prime_factor_t){
        .packed = (factor & UINT64_C(0x00FFFFFFFFFFFFFF))
                | ((uint64_t)power << 56)
    };
 }

 static inline uint64_t pf_factor(prime_factor_t pf)
 { return pf.packed & UINT64_C(0x00FFFFFFFFFFFFFF); }

 typedef struct {
    prime_factor_t factors[FACTORIZATION_BUFFER];
    uint8_t        status;
 } factorization_t;

 /*
 * fac_* accessors -- procedural state helpers
 * Contract: -- (procedural, not data mappings)
 */
 static inline uint8_t fac_count(const factorization_t *f)
 { return f->status & 0x1Fu; }
 static inline bool fac_gate_x(const factorization_t *f)
 { return (f->status & 0x80u) != 0u; }
 static inline void fac_set_count(factorization_t *f, uint8_t n)
 { f->status = (f->status & 0x80u) | (n & 0x1Fu); }
 static inline void fac_set_gate_x(factorization_t *f)
 { f->status |= 0x80u; }

 /* =================================================================
 * GLOBAL STATE
 * ================================================================= */

 alignas(64) static uint64_t g_primes[MAX_PRIMES];
 static uint32_t g_prime_count;

 /* SPF table -- uint16_t, mmap'd, ~31 MiB. 32 entries/cacheline. */
 static uint16_t *g_spf;

 alignas(64) static uint64_t g_packed_cat[MAX_PRIMES / 16 + 1];

 typedef struct { uint32_t count; uint64_t smallest; uint64_t largest; } cat_stats_t;
 alignas(64) static cat_stats_t g_stats[MAX_CATEGORIES];

 typedef struct { uint64_t primes[MAX_PRIMES_PER_CATEGORY]; uint32_t count; } cat_prime_list_t;
 static cat_prime_list_t g_cat_primes[MAX_CATEGORIES];

 /* =================================================================
 * I5 -- 4-BIT PACKED MEMO
 * 0x00=empty(Z)  0x01..0x0E=category(1)  0x0F=Gate-X(X)
 *
 * memo_set / memo_get -- bijective store/recall
 * FRONT: (idx, cat) -- prime ordinal and category nibble (AS)
 * LEAD:  bit-field write/read into g_packed_cat (Pivot)
 * REAR:  memo_get(idx) == cat for all valid idx, cat (IS)
 *   1: always (no failure path within range)
 * Contract: {{0 [ (uint32_t,uint8_t) (AS/.\IS) void ] 1}} (memo_set)
 *           {{0 [ uint32_t (AS/.\IS) uint8_t ] 1}} (memo_get)
 *   Bijective: 4-bit field, no aliasing between adjacent indices.
 *   Precondition: idx < MAX_PRIMES (enforced by caller -- ordinal bounds).
 * ================================================================= */

 static inline void memo_set(uint32_t idx, uint8_t cat)
 {
    uint32_t wi = idx / 16u;
    uint32_t bp = (idx % 16u) * 4u;
    g_packed_cat[wi] = (g_packed_cat[wi] & ~((uint64_t)0xFu << bp))
                     | ((uint64_t)cat << bp);
 }

 static inline uint8_t memo_get(uint32_t idx)
 {
    return (uint8_t)((g_packed_cat[idx / 16u] >> ((idx % 16u) * 4u)) & 0xFu);
 }

 /* =================================================================
 * I1 -- PRIME SIEVE (Eratosthenes, memset 0xAA for odd bits)
 *
 * FRONT: SIEVE_LIMIT, MAX_PRIMES (AS)
 * LEAD:  Eratosthenes inner loop (Pivot)
 * REAR:  g_primes[] filled, g_prime_count set (IS)
 *   X: malloc fail
 *   0: count != MAX_PRIMES (fewer primes than expected in range)
 *   1: exactly MAX_PRIMES found
 * Contract: {{0 [ void -- (AS/--\WAS) bool ] 1}}
 *   Procedural: constructs global state. Not a data mapping.
 * ================================================================= */

 [[nodiscard]]
 static bool sieve_generate(void)
 {
    size_t   nbytes = (SIEVE_LIMIT + 7u) / 8u;
    uint8_t *sieve  = malloc(nbytes);
    if (!sieve) return false;

    memset(sieve, 0xAAu, nbytes);

    for (uint32_t p = 3u; (uint64_t)p * p <= SIEVE_LIMIT; p += 2u)
        if (sieve[p / 8u] & (1u << (p % 8u)))
            for (uint32_t m = p * p; m <= SIEVE_LIMIT; m += 2u * p)
                sieve[m / 8u] &= (uint8_t)~(1u << (m % 8u));

    g_prime_count = 0u;
    g_primes[g_prime_count++] = 2u;
    for (uint32_t i = 3u; i <= SIEVE_LIMIT && g_prime_count < MAX_PRIMES; i += 2u)
        if (sieve[i / 8u] & (1u << (i % 8u)))
            g_primes[g_prime_count++] = i;

    free(sieve);
    return g_prime_count == MAX_PRIMES;
 }

 /* =================================================================
 * I3 -- SPF TABLE BUILD (uint16_t)
 * Prime sentinel=0 (MAP_ANONYMOUS). sqrt(16M)=4000 < UINT16_MAX.
 * THP: ~7,813 (4KB) -> ~16 (2MB huge pages).
 *
 * FRONT: g_primes[] filled (AS)
 * LEAD:  for each prime p, mark multiples where spf==0 (Pivot)
 * REAR:  g_spf[n] = smallest prime factor of n, or 0 if n is prime (IS)
 *   X: mmap fail
 *   1: g_spf ready
 * Contract: {{0 [ void -- (AS/--\WAS) bool ] 1}}
 *   Procedural: constructs global mmap'd state.
 * ================================================================= */

 [[nodiscard]]
 static bool spf_build(void)
 {
    size_t spf_bytes = (size_t)(SIEVE_LIMIT + 1u) * sizeof(uint16_t);

    g_spf = mmap(nullptr, spf_bytes, PROT_READ | PROT_WRITE,
                 MAP_PRIVATE | MAP_ANONYMOUS, -1, 0);
    if (g_spf == MAP_FAILED) { g_spf = nullptr; return false; }

    madvise(g_spf, spf_bytes, MADV_HUGEPAGE);

    for (uint32_t i = 0u; i < g_prime_count; i++) {
        uint32_t p = (uint32_t)g_primes[i];
        for (uint32_t m = p + p; m <= SIEVE_LIMIT; m += p)
            if (g_spf[m] == 0u)
                g_spf[m] = (uint16_t)p;
    }
    return true;
 }

 /* =================================================================
 * I2 -- RAMSTORE BUILD
 *
 * FRONT: g_primes[] filled (AS)
 * LEAD:  first empty sentinel in probe sequence (Pivot)
 * REAR:  g_rs[pos].{prime,index} committed for all primes (IS)
 *   X: probe >= RS_PROBE_MAX -> false (structural violation)
 *   1: all primes inserted, O(1) lookup ready
 * Contract: {{0 [ void -- (AS/--\WAS) bool ] 1}}
 *   Procedural: constructs global RAMstore state.
 *
 * Precondition: load factor 23.84% guarantees E[probe]=1.16.
 * assert() is Class 3 enforcement (debug only).
 * Gate-X fallback (return false) is Class 2 enforcement (always).
 * Both present: assert catches in debug, return false catches in release.
 * ================================================================= */

 [[nodiscard]]
 static bool rs_build(void)
 {
    for (uint32_t i = 0u; i < g_prime_count; i++) {
        uint64_t prime = g_primes[i];
        uint32_t pos   = ee_hash64_pos(prime);
        uint32_t p;

        for (p = 0u; p < RS_PROBE_MAX; p++) {
            uint32_t s = (pos + p) & HASH_MASK;
            if (g_rs[s].prime == 0u) {
                g_rs[s].prime = prime;
                g_rs[s].index = i;
                break;
            }
        }

        /* Class 3 (debug): assert fires in -O1 -g builds */
        assert(p < RS_PROBE_MAX &&
               "rs_build: probe cap exceeded -- increase HASH_BITS or RS_PROBE_MAX");
        /* Class 2 (release): Gate-X fallback always present */
        if (p >= RS_PROBE_MAX) return false;
    }
    return true;
 }

 /* =================================================================
 * I2 -- RAMSTORE LOOKUP
 *
 * FRONT: prime -- uint64_t to look up (AS)
 * LEAD:  equality test per probe (Pivot)
 * REAR:  ordinal index or UINT32_MAX if absent (IS)
 *   Z: prime not in store (UINT32_MAX)
 *   1: index returned
 * Expected ~1.16 probes at 23.84% load. 4 entries per cacheline.
 * Contract: {{0 [ uint64_t (AS/--\WAS) uint32_t ] 1}}
 *   Lossy: maps every prime to its index; absent primes map to UINT32_MAX.
 * ================================================================= */

 [[nodiscard]]
 static inline uint32_t rs_lookup(uint64_t prime)
 {
    uint32_t pos = ee_hash64_pos(prime);
    for (uint32_t p = 0u; p < RS_PROBE_MAX; p++) {
        uint32_t s = (pos + p) & HASH_MASK;
        if (g_rs[s].prime == 0u)    return UINT32_MAX;
        if (g_rs[s].prime == prime) return g_rs[s].index;
    }
    return UINT32_MAX;
 }

 /* =================================================================
 * I3 -- FACTORIZE via SPF table
 *
 * FRONT: n -- integer to factorize (AS)
 * LEAD:  n==1 terminates; spf[n] gives smallest factor (Pivot)
 * REAR:  factorization_t with all prime factors and powers (IS)
 *   Z: n<=1 (empty factorization, no factors)
 *   X: factor count exceeds FACTORIZATION_BUFFER (gate_x set)
 *   1: complete lossless factorization
 * Contract: {{0 [ uint64_t (AS/--\WAS) factorization_t* ] 1}}
 *   Lossy: factorization_t is a projection of n; n is not recoverable
 *   without multiplying factors back (not stored in out).
 * ================================================================= */

 static void factorize_spf(uint64_t n, factorization_t *out)
 {
    out->status = 0u;
    if (n <= 1u) return;

    while (n > 1u) {
        uint64_t f;
        if (n <= SIEVE_LIMIT) {
            uint16_t sf = g_spf[(uint32_t)n];
            f = sf ? (uint64_t)sf : n;
        } else {
            f = n;
        }

        uint8_t cnt = fac_count(out);
        if (cnt >= FACTORIZATION_BUFFER) { fac_set_gate_x(out); return; }

        uint8_t power = 0u;
        while (n % f == 0u) { n /= f; power++; }

        out->factors[cnt] = pf_make(f, power);
        fac_set_count(out, (uint8_t)(cnt + 1u));
    }
 }

 /* =================================================================
 * I4 -- COMPUTE CATEGORY (recursive + memoized)
 * For p>3, every prime factor f of p+1 satisfies f<p (forward order).
 *
 * FRONT: pi -- prime ordinal in g_primes[] (AS)
 * LEAD:  factors of p+1 resolved via memo/recursion (Pivot)
 * REAR:  category u8 -- depth to {2,3}-smooth base (IS)
 *   Z: already memoized (fast path, no work)
 *   X: factor missing from RAMstore -> CATEGORY_GATE_X propagated
 *   1: category in [1, MAX_CATEGORIES-1] returned and memoized
 * Contract: {{0 [ uint32_t (AS/--\WAS) uint8_t ] 1}}
 *   Lossy: many primes map to the same category depth.
 *   Gate-X propagation: never swallowed, always returned to caller.
 * ================================================================= */

 static uint8_t compute_category(uint32_t pi)
 {
    uint8_t cached = memo_get(pi);
    if (cached != 0u) return cached;

    factorization_t fac;
    factorize_spf(g_primes[pi] + 1ull, &fac);

    if (fac_gate_x(&fac)) {
        memo_set(pi, CATEGORY_GATE_X);
        return CATEGORY_GATE_X;
    }

    uint8_t max_cat        = 0u;
    bool    only_base_fact = true;

    for (uint8_t i = 0u; i < fac_count(&fac); i++) {
        uint64_t f = pf_factor(fac.factors[i]);

        if (f == 2u || f == 3u) {
            if (max_cat < 1u) max_cat = 1u;
            continue;
        }

        only_base_fact = false;
        uint32_t fi = rs_lookup(f);
        if (fi == UINT32_MAX) { memo_set(pi, CATEGORY_GATE_X); return CATEGORY_GATE_X; }

        uint8_t fc = compute_category(fi);
        if (fc == CATEGORY_GATE_X) { memo_set(pi, CATEGORY_GATE_X); return CATEGORY_GATE_X; }
        if (fc > max_cat) max_cat = fc;
    }

    uint8_t result = only_base_fact ? 1u : (uint8_t)(max_cat + 1u);
    memo_set(pi, result);
    return result;
 }

 /* =================================================================
 * PHI_ANALYSIS -- compile with -DPHI_ANALYSIS to enable.
 *
 * PURPOSE: verification of the invariant, not analysis.
 * The phi ratios phi(p+1)/(p+1) MUST increase monotonically with
 * category depth. If they do not, the categorization is wrong.
 * This flag exposes the 4th projection of the {2,3,5}-smooth invariant:
 *   Depth-1: p+1 is {2,3}-smooth -> minimum phi ratio (maximum collapse)
 *   Depth-k: more distinct prime factors -> higher phi ratio
 * Run pe_phi and verify the phi column increases with category.
 *
 * Feature tax: ~5% slower (measured on hardware). Default off. Full-speed build is the default.
 * ================================================================= */

 #ifdef PHI_ANALYSIS

 static uint64_t phi_from_factors(uint64_t n, const factorization_t *fac)
 {
    uint64_t phi = n;
    for (uint8_t i = 0u; i < fac_count(fac); i++) {
        uint64_t f = pf_factor(fac->factors[i]);
        phi = phi / f * (f - 1u);
    }
    return phi;
 }

 static uint8_t compute_category_fac(uint32_t pi,
                                    const factorization_t *fac)
 {
    uint8_t cached = memo_get(pi);
    if (cached != 0u) return cached;

    if (fac_gate_x(fac)) {
        memo_set(pi, CATEGORY_GATE_X);
        return CATEGORY_GATE_X;
    }

    uint8_t max_cat        = 0u;
    bool    only_base_fact = true;

    for (uint8_t i = 0u; i < fac_count(fac); i++) {
        uint64_t f = pf_factor(fac->factors[i]);

        if (f == 2u || f == 3u) {
            if (max_cat < 1u) max_cat = 1u;
            continue;
        }

        only_base_fact = false;
        uint32_t fi = rs_lookup(f);
        if (fi == UINT32_MAX) { memo_set(pi, CATEGORY_GATE_X); return CATEGORY_GATE_X; }

        uint8_t fc = compute_category(fi);
        if (fc == CATEGORY_GATE_X) { memo_set(pi, CATEGORY_GATE_X); return CATEGORY_GATE_X; }
        if (fc > max_cat) max_cat = fc;
    }

    uint8_t result = only_base_fact ? 1u : (uint8_t)(max_cat + 1u);
    memo_set(pi, result);
    return result;
 }

 typedef struct { uint64_t phi_sum; uint64_t n_sum; } phi_accum_t;

 #endif /* PHI_ANALYSIS */

 /* =================================================================
 * SELF-TEST (16 checks). Runs before run(). Halts on failure.
 *
 * FRONT: post-build state -- sieve, SPF, RAMstore ready (AS)
 * LEAD:  assertions against hand-derived expected values (Pivot)
 * REAR:  bool -- all checks passed (IS)
 *   0: any check failed -- halt, do not proceed
 *   1: all 16 checks passed -- proceed to run()
 * Contract: {{0 [ void -- (AS/--\WAS) bool ] 1}}
 *   Procedural: reads global state, produces pass/fail signal.
 * ================================================================= */

 static bool self_test(void)
 {
    #define ST_CHECK(expr, msg) \
    do { if (!(expr)) { fprintf(stderr, "[self_test] FAIL: %s\n", (msg)); return false; } } while (0)

    /* Check 1: feistel_encrypt_36 four-point injectivity (10 of 16) */
    {
        const uint64_t FMAX = (1ull << 36) - 1ull;
        uint64_t fa = feistel_encrypt_36(0ull);
        uint64_t fb = feistel_encrypt_36(1ull);
        uint64_t fc = feistel_encrypt_36(1ull << 18);
        uint64_t fd = feistel_encrypt_36(FMAX);
        ST_CHECK(fa != fb, "feistel: collision {0,1}");
        ST_CHECK(fa != fc, "feistel: collision {0,1<<18}");
        ST_CHECK(fa != fd, "feistel: collision {0,max}");
        ST_CHECK(fb != fc, "feistel: collision {1,1<<18}");
        ST_CHECK(fb != fd, "feistel: collision {1,max}");
        ST_CHECK(fc != fd, "feistel: collision {1<<18,max}");
        ST_CHECK(fa <= FMAX, "feistel(0): out of range");
        ST_CHECK(fb <= FMAX, "feistel(1): out of range");
        ST_CHECK(fc <= FMAX, "feistel(1<<18): out of range");
        ST_CHECK(fd <= FMAX, "feistel(max): out of range");
    }

    /* Check 2: hash range (1) */
    ST_CHECK(ee_hash64_pos(2u) < HASH_SIZE, "hash out of range");

    /* Check 3: pf round-trip (1)
     * Derivation: pf_make(97,2).packed = 97 | (2<<56)
     *             pf_factor extracts bits[55:0] = 97  CONFIRMED */
    ST_CHECK(pf_factor(pf_make(97u, 2u)) == 97u, "pf round-trip broken");

    /* Check 4: cat(2)=1 (2).
     * Derivation: 2+1=3. 3 is {2,3}-smooth. depth=1.  CONFIRMED */
    ST_CHECK(g_primes[0] == 2u,            "g_primes[0] != 2");
    ST_CHECK(compute_category(0u) == 1u,   "cat(2) != 1");

    /* Check 5: cat(13)=2 (2).
     * Derivation: 13+1=14=2*7. cat(7)=1 (7+1=8=2^3, {2}-smooth).
     *             max_cat=1, not only_base_fact -> result=1+1=2.  CONFIRMED */
    ST_CHECK(g_primes[5] == 13u,           "g_primes[5] != 13");
    ST_CHECK(compute_category(5u) == 2u,   "cat(13) != 2");

    #undef ST_CHECK
    return true;
 }

 /* =================================================================
 * MAIN PROCESSING LOOP
 *
 * FRONT: all build phases complete, self_test passed (AS)
 * LEAD:  compute_category() per prime (Pivot)
 * REAR:  stats populated, timing printed (IS)
 *   X: (unused -- run() has no failure path post-build)
 *   1: categorization complete
 * Contract: {{0 [ void -- (AS/--\WAS) bool ] 1}}
 *   Procedural: reads global state, writes stats, prints output.
 * ================================================================= */

 [[nodiscard]]
 static bool run(void)
 {
    struct timespec t0, t1;
    clock_gettime(CLOCK_MONOTONIC, &t0);

    memset(g_stats,      0u, sizeof(g_stats));
    memset(g_packed_cat, 0u, sizeof(g_packed_cat));
    memset(g_cat_primes, 0u, sizeof(g_cat_primes));

    uint32_t gate_x_count = 0u;
 #ifdef PHI_ANALYSIS
    phi_accum_t phi_stats[MAX_CATEGORIES] = {{0}};
 #endif

    printf("Processing %" PRIu32 " primes...\n", g_prime_count);

    for (uint32_t i = 0u; i < g_prime_count; i++) {
        uint64_t p = g_primes[i];

 #ifdef PHI_ANALYSIS
        uint64_t n = p + 1ull;
        factorization_t fac;
        factorize_spf(n, &fac);
        uint8_t cat = compute_category_fac(i, &fac);
 #else
        uint8_t cat = compute_category(i);
 #endif

        if (cat == CATEGORY_GATE_X) { gate_x_count++; continue; }
        if (cat == 0u || cat >= MAX_CATEGORIES) continue;

        g_stats[cat].count++;
        if (g_stats[cat].count == 1u) {
            g_stats[cat].smallest = p;
            g_stats[cat].largest  = p;
        } else if (p > g_stats[cat].largest) {
            g_stats[cat].largest = p;
        }

 #ifdef PHI_ANALYSIS
        if (!fac_gate_x(&fac)) {
            phi_stats[cat].phi_sum += phi_from_factors(n, &fac);
            phi_stats[cat].n_sum   += n;
        }
 #endif

        if (__builtin_expect(i < TRACK_FIRST_N_PRIMES, 0) &&
            g_cat_primes[cat].count < MAX_PRIMES_PER_CATEGORY)
            g_cat_primes[cat].primes[g_cat_primes[cat].count++] = p;
    }

    clock_gettime(CLOCK_MONOTONIC, &t1);
    ee_ratio_t elapsed = ee_ratio_elapsed(t0, t1);
    uint64_t   secs    = ee_ratio_secs(elapsed);
    uint64_t   frac    = ee_ratio_frac10k(elapsed);
    uint64_t   tput    = ee_ratio_throughput(g_prime_count, elapsed);

    printf("Done in %" PRIu64 ".%04" PRIu64 " s\n\n", secs, frac);

    if (gate_x_count)
        printf("[Gate-X] %" PRIu32 " prime(s) excluded\n\n", gate_x_count);

    printf("First %" PRIu32 " primes by category:\n", TRACK_FIRST_N_PRIMES);
    for (uint8_t cat = 1u; cat < MAX_CATEGORIES; cat++) {
        if (!g_cat_primes[cat].count) continue;
        printf("Category %u:\n", cat);
        for (uint32_t i = 0u; i < g_cat_primes[cat].count; i++) {
            printf("%4" PRIu64, g_cat_primes[cat].primes[i]);
            if ((i + 1u) % 15u == 0u || i == g_cat_primes[cat].count - 1u)
                printf("\n");
            else
                printf(" ");
        }
        printf("\n");
    }

    printf("Erdos-Selfridge categorization:\n");
    for (uint8_t cat = 1u; cat < MAX_CATEGORIES; cat++) {
        if (!g_stats[cat].count) continue;
        printf("  Cat %2u: first=%7" PRIu64
               "  last=%8" PRIu64 "  count=%" PRIu32 "\n",
               cat, g_stats[cat].smallest,
               g_stats[cat].largest, g_stats[cat].count);
    }

 #ifdef PHI_ANALYSIS
    printf("\nUnit density phi(p+1)/(p+1) by category:\n");
    printf("  (Verification: ratios must increase monotonically with depth.)\n");
    printf("  (If they do not, the categorization is wrong.)\n");
    for (uint8_t cat = 1u; cat < MAX_CATEGORIES; cat++) {
        if (!phi_stats[cat].n_sum) continue;
        uint64_t pct_10k = (phi_stats[cat].phi_sum * 10000ull)
                         / phi_stats[cat].n_sum;
        printf("  Cat %2u: mean phi(p+1)/(p+1) = 0.%04" PRIu64
               "  (%" PRIu32 " primes)\n",
               cat, pct_10k, g_stats[cat].count);
    }
 #endif

    printf("\nPerformance (ee_ratio_t -- CLOCK_MONOTONIC, no IEEE 754):\n");
    printf("  Elapsed:    %" PRIu64 ".%04" PRIu64 " s\n", secs, frac);
    printf("  Throughput: %" PRIu64 " primes/s\n", tput);

    return true;
 }

 int main(void)
 {
    printf("Erdos-Selfridge Prime Categorization Engine v3.5.6\n");
    printf("C23 / RAMstore + u16-SPF + THP\n");
    printf("================================================\n\n");

    printf("Sieving %" PRIu32 " primes...\n", MAX_PRIMES);
    if (!sieve_generate()) { fprintf(stderr, "sieve_generate: failed\n"); return 1; }
    printf("Sieve: %" PRIu32 " primes\n\n", g_prime_count);

    printf("Building SPF table (%" PRIu32 " entries, 31 MiB u16, THP)...\n", SIEVE_LIMIT);
    if (!spf_build()) { fprintf(stderr, "spf_build: mmap failed\n"); return 1; }
    printf("SPF: ready\n\n");

    printf("Building RAMstore (%" PRIu32 " positions, 64 MiB)...\n", HASH_SIZE);
    if (!rs_build()) {
        fprintf(stderr, "rs_build: probe cap exceeded\n");
        munmap(g_spf, (size_t)(SIEVE_LIMIT + 1u) * sizeof(uint16_t));
        return 1;
    }
    printf("RAMstore: ready\n\n");

    printf("Self-test...\n");
    if (!self_test()) {
        fprintf(stderr, "self_test: FAILED -- aborting\n");
        munmap(g_spf, (size_t)(SIEVE_LIMIT + 1u) * sizeof(uint16_t));
        return 1;
    }
    printf("Self-test: PASS (16 checks)\n\n");

    if (!run()) {
        fprintf(stderr, "run: failed\n");
        munmap(g_spf, (size_t)(SIEVE_LIMIT + 1u) * sizeof(uint16_t));
        return 1;
    }

    munmap(g_spf, (size_t)(SIEVE_LIMIT + 1u) * sizeof(uint16_t));
    return 0;
 }

 /*
 * ACKNOWLEDGMENTS
 * H. Overman (ee): RAMstore architecture, Euman language ecosystem,
 *   144,000 ratio system, FRONT/LEAD/REAR contract discipline,
 * LeeMetaX  ExCLisp (AS/.\IS) contract operators.
 * Projection discipline: the observation that the {2,3,5}-smooth
 *   invariant projects onto four surfaces (exact arithmetic, prime
 *   categorization, harmonic analysis, phi collapse pattern). The
 *   phi(p+1)/(p+1) verification method. The Möbius topology of the
 *   FRONT/LEAD/REAR contract.
 * Horst Feistel (IBM): Feistel network (1973).
 * Ronald Rivest (MIT): RC5 (1994) -- FEISTEL_MASTER from frac(e).
 * Donald Knuth: TAOCP Vol 3 -- round constant from frac(phi).
 * IEEE 1364 / Verilog: Z/X/0/1 four-state value system.
 * Paul Erdos & John Selfridge: prime categorization.
 * Eratosthenes of Cyrene: prime sieve (~240 BCE).
 * Euler: totient function phi(n).
 */
No results found