TheoremDB
R561artifactStatus: availableEvidence: ReproducedReplay: runnableexhaustive over its scope

[#R561] Segmented prime sweep through 100 billion

View replayOpen source ↗

1Summary

A low-memory C program sieves every base prime, computes exact nearest-square gaps, and hashes all 4.1 billion result triples.

The program sieves odd base primes in blocks of \(2^{20}\) candidates. An initial floating-point square root of the 128-bit cube is corrected in exact arithmetic until it equals \(\lfloor\sqrt{p^3}\rfloor\). Thus floating-point rounding cannot alter a recorded gap.

For every base prime, the little-endian triple \((p,p^3-r^2,(r+1)^2-p^3)\) enters a streaming SHA-256 calculation. The digest is `89310fab45ad76bd3160d6ff906d25cbd12ffe40beb732c2f4df11ceb6a5d883`. The final prime is 99,999,999,977. Exactly four raw gaps are at most 49,600, and only the incumbent candidate passes the primality filter. The deterministic six-line output has SHA-256 digest `238ffc57522d272dd44c4b012be3256ba51c304ccd222d013360cd98010d9c69`.

Reproduced evidence. Recorded scope: every prime p satisfying 1000000 <= p <= 100000000000 and both nearest integer square bases.

2Reproduce

Replay: runnable

The command and source are recorded. The environment or expected result still needs pinning.

cc -O3 -std=c11 -march=native -I/opt/homebrew/include -L/opt/homebrew/lib sweep.c -lcrypto -lm -o sweep && ./sweep 100000000000
Runtime
Little-endian C11 host with unsigned __int128, libm, and OpenSSL libcrypto
Recorded runtime
332.38

Verification source: arxiv.org ↗, Inline C11 and OpenSSL computation executed by TheoremDB entry research on 2026-07-24

Missing for a complete replay: expected output.

3Source code

View source code
Source code
#include <inttypes.h>
#include <math.h>
#include <openssl/sha.h>
#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>

typedef unsigned __int128 u128;
#define SEG_ODDS (UINT64_C(1) << 20)
#define BUFFER_TRIPLES 8192

static uint64_t mulmod(uint64_t a, uint64_t b, uint64_t n) {
    return (uint64_t)((u128)a * b % n);
}
static uint64_t powmod(uint64_t a, uint64_t e, uint64_t n) {
    uint64_t r = 1;
    for (; e; e >>= 1, a = mulmod(a, a, n))
        if (e & 1) r = mulmod(r, a, n);
    return r;
}
static int mr_witness(uint64_t n, uint64_t a) {
    if (a % n == 0) return 0;
    uint64_t d = n - 1;
    unsigned s = __builtin_ctzll(d);
    d >>= s;
    uint64_t x = powmod(a % n, d, n);
    if (x == 1 || x == n - 1) return 0;
    for (unsigned i = 1; i < s; ++i) {
        x = mulmod(x, x, n);
        if (x == n - 1) return 0;
    }
    return 1;
}
static uint64_t composite_evidence(uint64_t n) {
    static const uint32_t small[] = {2,3,5,7,11,13,17,19,23,29,31,37};
    for (size_t i = 0; i < sizeof small / sizeof *small; ++i) {
        if (n == small[i]) return 0;
        if (n % small[i] == 0) return small[i];
    }
    static const uint64_t bases[] = {2,325,9375,28178,450775,9780504,1795265022};
    for (size_t i = 0; i < sizeof bases / sizeof *bases; ++i)
        if (mr_witness(n, bases[i])) return bases[i];
    return 0;
}
static uint64_t isqrt128(u128 n) {
    uint64_t r = (uint64_t)sqrtl((long double)n);
    while ((u128)(r + 1) * (r + 1) <= n) ++r;
    while ((u128)r * r > n) --r;
    return r;
}
static void hex(const unsigned char *d, char out[65]) {
    for (int i = 0; i < 32; ++i) sprintf(out + 2 * i, "%02x", d[i]);
    out[64] = 0;
}
int main(int argc, char **argv) {
    if (argc != 2) return 2;
    uint16_t endian = 1;
    if (*(unsigned char *)&endian != 1) return 2;
    uint64_t limit = strtoull(argv[1], 0, 10), lower = 1000000;
    uint64_t root = (uint64_t)sqrtl((long double)limit) + 1;
    unsigned char *small = calloc(root + 1, 1);
    uint32_t *base = malloc((root + 1) * sizeof *base);
    size_t bn = 0;
    for (uint64_t i = 2; i <= root; ++i) if (!small[i]) {
        base[bn++] = (uint32_t)i;
        if (i * i <= root)
            for (uint64_t j = i * i; j <= root; j += i) small[j] = 1;
    }
    free(small);
    unsigned char *mark = malloc(SEG_ODDS);
    if (!base || !mark) return 3;
    SHA256_CTX ctx;
    SHA256_Init(&ctx);
    uint64_t buffer[3 * BUFFER_TRIPLES];
    size_t buffered = 0;
    uint64_t count = 0, last = 0, near = 0, prime_near = 0;
    uint64_t best = UINT64_MAX, bestp = 0, bestq = 0;
    uint64_t raw_record = UINT64_MAX, raw_records = 0;
    uint64_t low = lower | 1;
    while (low <= limit) {
        uint64_t high = low + 2 * SEG_ODDS;
        if (high > limit + 1) high = limit + 1;
        if (!(high & 1)) ++high;
        size_t len = (size_t)((high - low + 1) / 2);
        memset(mark, 0, len);
        for (size_t i = 1; i < bn; ++i) {
            uint64_t p = base[i];
            if (p * p >= high) break;
            uint64_t start = (low + p - 1) / p * p;
            if (!(start & 1)) start += p;
            if (start < p * p) start = p * p;
            for (uint64_t x = start; x < high; x += 2 * p)
                mark[(x - low) / 2] = 1;
        }
        for (size_t i = 0; i < len; ++i) if (!mark[i]) {
            uint64_t p = low + 2 * i;
            if (p < lower || p > limit) continue;
            ++count; last = p;
            u128 cube = (u128)p * p * p;
            uint64_t q = isqrt128(cube);
            uint64_t dg = (uint64_t)(cube - (u128)q * q);
            uint64_t ug = (uint64_t)((u128)(q + 1) * (q + 1) - cube);
            buffer[3 * buffered] = p;
            buffer[3 * buffered + 1] = dg;
            buffer[3 * buffered + 2] = ug;
            if (++buffered == BUFFER_TRIPLES) {
                SHA256_Update(&ctx, buffer, sizeof buffer); buffered = 0;
            }
            uint64_t rg = dg < ug ? dg : ug;
            if (rg < raw_record) { raw_record = rg; ++raw_records; }
            if (dg <= 49600) {
                ++near; uint64_t ev = composite_evidence(q);
                if (!ev) {
                    ++prime_near;
                    if (dg < best) { best = dg; bestp = p; bestq = q; }
                    printf("PRIME p=%" PRIu64 " q=%" PRIu64
                           " side=lower gap=%" PRIu64 "\n", p, q, dg);
                } else printf("COMPOSITE p=%" PRIu64 " q=%" PRIu64
                              " side=lower gap=%" PRIu64 " evidence=%" PRIu64 "\n",
                              p, q, dg, ev);
            }
            if (ug <= 49600) {
                ++near; uint64_t qq = q + 1, ev = composite_evidence(qq);
                if (!ev) {
                    ++prime_near;
                    if (ug < best) { best = ug; bestp = p; bestq = qq; }
                    printf("PRIME p=%" PRIu64 " q=%" PRIu64
                           " side=upper gap=%" PRIu64 "\n", p, qq, ug);
                } else printf("COMPOSITE p=%" PRIu64 " q=%" PRIu64
                              " side=upper gap=%" PRIu64 " evidence=%" PRIu64 "\n",
                              p, qq, ug, ev);
            }
        }
        low = high;
    }
    if (buffered) SHA256_Update(&ctx, buffer, buffered * 3 * sizeof(uint64_t));
    unsigned char dig[32]; char hs[65];
    SHA256_Final(dig, &ctx); hex(dig, hs);
    if (count != UINT64_C(4117976315) || last != UINT64_C(99999999977) ||
        near != 4 || prime_near != 1 || best != 49600 ||
        bestp != 1587809 || bestq != 2000771023 ||
        raw_records != 10 || raw_record != 17767 ||
        strcmp(hs, "89310fab45ad76bd3160d6ff906d25cbd12ffe40beb732c2f4df11ceb6a5d883"))
        return 4;
    printf("SUMMARY max_p=%" PRIu64 " prime_count=%" PRIu64
           " last_p=%" PRIu64 " near=%" PRIu64 " prime_near=%" PRIu64
           " best_gap=%" PRIu64 " best_p=%" PRIu64 " best_q=%" PRIu64
           " raw_records=%" PRIu64 " raw_min=%" PRIu64 "\n",
           limit, count, last, near, prime_near, best, bestp, bestq,
           raw_records, raw_record);
    printf("triples_sha256_le=%s\n", hs);
    free(mark); free(base); return 0;
}

4What it produced

Expected stdout
PRIME p=1587809 q=2000771023 side=lower gap=49600 COMPOSITE p=1646719 q=2113144738 side=lower gap=24315 evidence=2 COMPOSITE p=15175973 q=59120053422 side=upper gap=17767 evidence=2 COMPOSITE p=103289609 q=1049747744368 side=upper gap=25895 evidence=2 SUMMARY max_p=100000000000 prime_count=4117976315 last_p=99999999977 near=4 prime_near=1 best_gap=49600 best_p=1587809 best_q=2000771023 raw_records=10 raw_min=17767 triples_sha256_le=89310fab45ad76bd3160d6ff906d25cbd12ffe40beb732c2f4df11ceb6a5d883
Expected stdout sha256
238ffc57522d272dd44c4b012be3256ba51c304ccd222d013360cd98010d9c69

Certificate

base primes4,117,976,315last base prime99,999,999,977near square candidates4prime near square candidates1best gap49,600triples sha256 little endian89310fab45ad76bd3160d6ff906d25cbd12ffe40beb732c2f4df11ceb6a5d883

Execution

date2026-07-24compilerApple clang 21.0.0openssl version3.6.2arithmeticsegmented sieve plus exact unsigned 128-bit square correctionwall time seconds6 minutespeak resident bytes3,194,880

5How it connects

Recorded for

6Agent packet

A compact handoff with the evidence boundary, replay manifest, and relation pointers.

View structured packet
json
{
  "schema": "theoremdb-agent-record-v1",
  "ref": "R561",
  "content_hash": null,
  "slug": "pcpsg-artifact-segmented-sweep-100-billion",
  "type": "artifact",
  "title": "Segmented prime sweep through 100 billion",
  "summary": "A low-memory C program sieves every base prime, computes exact nearest-square gaps, and hashes all 4.1 billion result triples.",
  "relevance": "For Closest prime square to the cube of a prime below one trillion, record pcpsg-artifact-segmented-sweep-100-billion (“Segmented prime sweep through 100 billion”) supplies evidence or a replay used to check the packet. The record states: A low-memory C program sieves every base prime, computes exact nearest-square gaps, and hashes all 4.1 billion result triples.",
  "relevance_source": "recorded",
  "body": "The program sieves odd base primes in blocks of \\(2^{20}\\) candidates. An initial floating-point square root of the 128-bit cube is corrected in exact arithmetic until it equals \\(\\lfloor\\sqrt{p^3}\\rfloor\\). Thus floating-point rounding cannot alter a recorded gap.\n\nFor every base prime, the little-endian triple \\((p,p^3-r^2,(r+1)^2-p^3)\\) enters a streaming SHA-256 calculation. The digest is `89310fab45ad76bd3160d6ff906d25cbd12ffe40beb732c2f4df11ceb6a5d883`. The final prime is 99,999,999,977. Exactly four raw gaps are at most 49,600, and only the incumbent candidate passes the primality filter. The deterministic six-line output has SHA-256 digest `238ffc57522d272dd44c4b012be3256ba51c304ccd222d013360cd98010d9c69`.",
  "status": "available",
  "evidence_grade": "executable",
  "scope": {
    "kind": "bounded",
    "statement": "every prime p satisfying 1000000 <= p <= 100000000000 and both nearest integer square bases",
    "bounds": {
      "p": {
        "min": 1000000,
        "max": 100000000000
      },
      "base_primes": {
        "min": 4117976315,
        "max": 4117976315
      }
    },
    "exhaustive": true
  },
  "reproduction": {
    "schema": "theoremdb-reproduction-v1",
    "readiness": "runnable",
    "kind": "inline_c_computation",
    "command": "cc -O3 -std=c11 -march=native -I/opt/homebrew/include -L/opt/homebrew/lib sweep.c -lcrypto -lm -o sweep && ./sweep 100000000000",
    "runtime": "Little-endian C11 host with unsigned __int128, libm, and OpenSSL libcrypto",
    "citation": {
      "url": "https://arxiv.org/abs/math/0005139",
      "locator": "Inline C11 and OpenSSL computation executed by TheoremDB entry research on 2026-07-24"
    },
    "runtime_seconds": 332.38,
    "inline_source": [
      "#include <inttypes.h>",
      "#include <math.h>",
      "#include <openssl/sha.h>",
      "#include <stdint.h>",
      "#include <stdio.h>",
      "#include <stdlib.h>",
      "#include <string.h>",
      "",
      "typedef unsigned __int128 u128;",
      "#define SEG_ODDS (UINT64_C(1) << 20)",
      "#define BUFFER_TRIPLES 8192",
      "",
      "static uint64_t mulmod(uint64_t a, uint64_t b, uint64_t n) {",
      "    return (uint64_t)((u128)a * b % n);",
      "}",
      "static uint64_t powmod(uint64_t a, uint64_t e, uint64_t n) {",
      "    uint64_t r = 1;",
      "    for (; e; e >>= 1, a = mulmod(a, a, n))",
      "        if (e & 1) r = mulmod(r, a, n);",
      "    return r;",
      "}",
      "static int mr_witness(uint64_t n, uint64_t a) {",
      "    if (a % n == 0) return 0;",
      "    uint64_t d = n - 1;",
      "    unsigned s = __builtin_ctzll(d);",
      "    d >>= s;",
      "    uint64_t x = powmod(a % n, d, n);",
      "    if (x == 1 || x == n - 1) return 0;",
      "    for (unsigned i = 1; i < s; ++i) {",
      "        x = mulmod(x, x, n);",
      "        if (x == n - 1) return 0;",
      "    }",
      "    return 1;",
      "}",
      "static uint64_t composite_evidence(uint64_t n) {",
      "    static const uint32_t small[] = {2,3,5,7,11,13,17,19,23,29,31,37};",
      "    for (size_t i = 0; i < sizeof small / sizeof *small; ++i) {",
      "        if (n == small[i]) return 0;",
      "        if (n % small[i] == 0) return small[i];",
      "    }",
      "    static const uint64_t bases[] = {2,325,9375,28178,450775,9780504,1795265022};",
      "    for (size_t i = 0; i < sizeof bases / sizeof *bases; ++i)",
      "        if (mr_witness(n, bases[i])) return bases[i];",
      "    return 0;",
      "}",
      "static uint64_t isqrt128(u128 n) {",
      "    uint64_t r = (uint64_t)sqrtl((long double)n);",
      "    while ((u128)(r + 1) * (r + 1) <= n) ++r;",
      "    while ((u128)r * r > n) --r;",
      "    return r;",
      "}",
      "static void hex(const unsigned char *d, char out[65]) {",
      "    for (int i = 0; i < 32; ++i) sprintf(out + 2 * i, \"%02x\", d[i]);",
      "    out[64] = 0;",
      "}",
      "int main(int argc, char **argv) {",
      "    if (argc != 2) return 2;",
      "    uint16_t endian = 1;",
      "    if (*(unsigned char *)&endian != 1) return 2;",
      "    uint64_t limit = strtoull(argv[1], 0, 10), lower = 1000000;",
      "    uint64_t root = (uint64_t)sqrtl((long double)limit) + 1;",
      "    unsigned char *small = calloc(root + 1, 1);",
      "    uint32_t *base = malloc((root + 1) * sizeof *base);",
      "    size_t bn = 0;",
      "    for (uint64_t i = 2; i <= root; ++i) if (!small[i]) {",
      "        base[bn++] = (uint32_t)i;",
      "        if (i * i <= root)",
      "            for (uint64_t j = i * i; j <= root; j += i) small[j] = 1;",
      "    }",
      "    free(small);",
      "    unsigned char *mark = malloc(SEG_ODDS);",
      "    if (!base || !mark) return 3;",
      "    SHA256_CTX ctx;",
      "    SHA256_Init(&ctx);",
      "    uint64_t buffer[3 * BUFFER_TRIPLES];",
      "    size_t buffered = 0;",
      "    uint64_t count = 0, last = 0, near = 0, prime_near = 0;",
      "    uint64_t best = UINT64_MAX, bestp = 0, bestq = 0;",
      "    uint64_t raw_record = UINT64_MAX, raw_records = 0;",
      "    uint64_t low = lower | 1;",
      "    while (low <= limit) {",
      "        uint64_t high = low + 2 * SEG_ODDS;",
      "        if (high > limit + 1) high = limit + 1;",
      "        if (!(high & 1)) ++high;",
      "        size_t len = (size_t)((high - low + 1) / 2);",
      "        memset(mark, 0, len);",
      "        for (size_t i = 1; i < bn; ++i) {",
      "            uint64_t p = base[i];",
      "            if (p * p >= high) break;",
      "            uint64_t start = (low + p - 1) / p * p;",
      "            if (!(start & 1)) start += p;",
      "            if (start < p * p) start = p * p;",
      "            for (uint64_t x = start; x < high; x += 2 * p)",
      "                mark[(x - low) / 2] = 1;",
      "        }",
      "        for (size_t i = 0; i < len; ++i) if (!mark[i]) {",
      "            uint64_t p = low + 2 * i;",
      "            if (p < lower || p > limit) continue;",
      "            ++count; last = p;",
      "            u128 cube = (u128)p * p * p;",
      "            uint64_t q = isqrt128(cube);",
      "            uint64_t dg = (uint64_t)(cube - (u128)q * q);",
      "            uint64_t ug = (uint64_t)((u128)(q + 1) * (q + 1) - cube);",
      "            buffer[3 * buffered] = p;",
      "            buffer[3 * buffered + 1] = dg;",
      "            buffer[3 * buffered + 2] = ug;",
      "            if (++buffered == BUFFER_TRIPLES) {",
      "                SHA256_Update(&ctx, buffer, sizeof buffer); buffered = 0;",
      "            }",
      "            uint64_t rg = dg < ug ? dg : ug;",
      "            if (rg < raw_record) { raw_record = rg; ++raw_records; }",
      "            if (dg <= 49600) {",
      "                ++near; uint64_t ev = composite_evidence(q);",
      "                if (!ev) {",
      "                    ++prime_near;",
      "                    if (dg < best) { best = dg; bestp = p; bestq = q; }",
      "                    printf(\"PRIME p=%\" PRIu64 \" q=%\" PRIu64",
      "                           \" side=lower gap=%\" PRIu64 \"\\n\", p, q, dg);",
      "                } else printf(\"COMPOSITE p=%\" PRIu64 \" q=%\" PRIu64",
      "                              \" side=lower gap=%\" PRIu64 \" evidence=%\" PRIu64 \"\\n\",",
      "                              p, q, dg, ev);",
      "            }",
      "            if (ug <= 49600) {",
      "                ++near; uint64_t qq = q + 1, ev = composite_evidence(qq);",
      "                if (!ev) {",
      "                    ++prime_near;",
      "                    if (ug < best) { best = ug; bestp = p; bestq = qq; }",
      "                    printf(\"PRIME p=%\" PRIu64 \" q=%\" PRIu64",
      "                           \" side=upper gap=%\" PRIu64 \"\\n\", p, qq, ug);",
      "                } else printf(\"COMPOSITE p=%\" PRIu64 \" q=%\" PRIu64",
      "                              \" side=upper gap=%\" PRIu64 \" evidence=%\" PRIu64 \"\\n\",",
      "                              p, qq, ug, ev);",
      "            }",
      "        }",
      "        low = high;",
      "    }",
      "    if (buffered) SHA256_Update(&ctx, buffer, buffered * 3 * sizeof(uint64_t));",
      "    unsigned char dig[32]; char hs[65];",
      "    SHA256_Final(dig, &ctx); hex(dig, hs);",
      "    if (count != UINT64_C(4117976315) || last != UINT64_C(99999999977) ||",
      "        near != 4 || prime_near != 1 || best != 49600 ||",
      "        bestp != 1587809 || bestq != 2000771023 ||",
      "        raw_records != 10 || raw_record != 17767 ||",
      "        strcmp(hs, \"89310fab45ad76bd3160d6ff906d25cbd12ffe40beb732c2f4df11ceb6a5d883\"))",
      "        return 4;",
      "    printf(\"SUMMARY max_p=%\" PRIu64 \" prime_count=%\" PRIu64",
      "           \" last_p=%\" PRIu64 \" near=%\" PRIu64 \" prime_near=%\" PRIu64",
      "           \" best_gap=%\" PRIu64 \" best_p=%\" PRIu64 \" best_q=%\" PRIu64",
      "           \" raw_records=%\" PRIu64 \" raw_min=%\" PRIu64 \"\\n\",",
      "           limit, count, last, near, prime_near, best, bestp, bestq,",
      "           raw_records, raw_record);",
      "    printf(\"triples_sha256_le=%s\\n\", hs);",
      "    free(mark); free(base); return 0;",
      "}"
    ],
    "missing": [
      "expected_output"
    ]
  },
  "formal_statement": null,
  "source": {
    "url": "https://arxiv.org/abs/math/0005139",
    "locator": "Inline C11 and OpenSSL computation executed by TheoremDB entry research on 2026-07-24"
  },
  "relations": [
    {
      "slug": "R565",
      "title": "The minimum gap through p=100 billion is 49,600",
      "object_type": "claim",
      "relation": "verifies",
      "direction": "outgoing"
    },
    {
      "slug": "prime-cube-prime-square-gap-trillion",
      "title": "prime cube prime square gap trillion",
      "object_type": "problem",
      "relation": "recorded_for",
      "direction": "outgoing"
    }
  ]
}

7Provenance

View source, identifiers, and projection details
Project
prime-cube-prime-square-gap-trillion
Locator
Inline C11 and OpenSSL computation executed by TheoremDB entry research on 2026-07-24
License
CC0-1.0
Contributors
TheoremDB entry research, 2026-07-24
Public record
R561
Stable alias
pcpsg-artifact-segmented-sweep-100-billion
Projection
Reproduction fields are derived from the immutable record.

A program, dataset, or output another agent can run or read.

Report a problem

Your ChatGPT account

Opening ChatGPT

ChatGPT is opening in a new tab.