The primality test uses Fermat's Little Theorem, but I thought it was not 100% perfect. (Or perhaps it is for left-truncatable primes?)
I didn't know about the third parameter in pow(), either. Quite interesting.
I didn't know about the third parameter in pow(), either. Quite interesting.
I should go looking ...
252601: 41 61 101
294409: 37 73 109
399001: 31 61 211
410041: 41 73 137
488881: 37 73 181
512461: 31 61 271
1152271: 43 127 211
1193221: 31 61 631
1461241: 37 73 541
1857241: 31 181 331I could write a more comprehensive test - didn't seem worth it for this case, but it's an interesting exercise. It's easy to extend this test so the list of primes used goes as far as the smallest factor in each of these factorisations, say, 109. Then these would all be correctly identified as composite. It would also be easy to combine it with Pollard Rho or PPP tests, to extend the range before we need to do something more serious.
Anyway, here are the first 30 with factorisations
252601 : 41 61 101
294409 : 37 73 109
399001 : 31 61 211
410041 : 41 73 137
488881 : 37 73 181
512461 : 31 61 271
1152271 : 43 127 211
1193221 : 31 61 631
1461241 : 37 73 541
1615681 : 23 199 353
1857241 : 31 181 331
1909001 : 41 101 461
2508013 : 53 79 599
3057601 : 43 211 337
3581761 : 29 113 1093
3828001 : 101 151 251
4335241 : 53 157 521
5049001 : 31 271 601
5148001 : 41 241 521
5444489 : 29 197 953
5481451 : 31 151 1171
5968873 : 43 127 1093
6189121 : 61 241 421
6733693 : 109 163 379
6868261 : 43 211 757
7519441 : 41 241 761
8134561 : 37 109 2017
8927101 : 31 37 43 181
9439201 : 61 271 571
9494101 : 23 61 67 101The primality test in the post (Fermat's test) is interesting to analyze though. BTW how did you generate your list of false positives? I didn't mention it earlier, but as a shortcut I used the list from https://oeis.org/A001567 (the base-2 pseudoprimes up to 10^12), and found that the function in the post has 2690 of these false positives less than 10^12. I wrote out only the first 10 above (after piping through `factor`). Of these 2690 false positives, these have the 5 largest smallest prime factors:
681719179681: 5743 8527 13921
802204143001: 5981 8971 14951
747734199031: 6791 8731 12611
949631589089: 6917 10193 13469
251994268081: 354961 709921
I think compared to the Fermat test used in the post, the Miller-Rabin test is almost as easy to write (only a few lines longer), and it comes with stronger guarantees. Not only is it a test used in the real world (using it with just a few random bases gives a high degree of confidence), but it also comes with guarantees when a deterministic set of primes (like your idea) is used: for example with the first 13 primes (2 to 41), the first false positive is over 3.3 × 10^24 (https://en.wikipedia.org/w/index.php?title=Miller%E2%80%93Ra...). See also https://miller-rabin.appspot.com/
def is_probable_prime(n, num_random_trials=5):
"""
An implementation of the Miller-Rabin primality test.
If this function returns False, then n is definitely not prime.
If this function returns True, then n is probably prime:
the probability of n being composite is less than 1/4^(num_random_trials), and
nonzero only if n >= 3317044064679887385961981 > 3.3*10^24 > 1.3*2^81 (by https://arxiv.org/abs/1509.00864).
"""
assert n == abs(int(n)), 'Need a nonnegative integer, got %s' % n
if n < 2: return False
if n % 2 == 0: return n == 2
(s, d) = (0, n - 1)
while d % 2 == 0: (s, d) = (s + 1, d // 2)
assert n - 1 == (2**s) * d and d % 2 == 1
small_primes = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41]
for i in range(len(small_primes) + num_random_trials):
a = small_primes[i] if i < len(small_primes) else random.randrange(small_primes[-1], n - 1)
if a >= n: break
if not (pow(a, d, n) == 1 or any(pow(a, 2**j * d, n) == n - 1 for j in range(s))):
return False
return True> Yes the full program for left-truncatable primes happens to be fine ...
Yes - although it's a little lucky, but easy enough to run an extra confirmation.
> how did you generate your list of false positives?
Nothing clever - I literally wrote a slow trial division program, then ran through all integers up to some limit - 10^8 I think - and when the "is_prime" routine said it was prime I called the slow version to see if it was right. Machines are fast.
> the function in the post has 2690 of these false positives less than 10^12. ... Of these 2690 false positives, these have the 5 largest smallest prime factors:
> 681719179681: 5743 8527 13921
> 802204143001: 5981 8971 14951
> 747734199031: 6791 8731 12611
> 949631589089: 6917 10193 13469
> 251994268081: 354961 709921
> And indeed the `small_primes` used by the test needs to be extended up to the prime 6917, if we want to get no false positives at least for numbers up to 10^12.
We can augment this with the Perrin Psuedo-Prime test, but Miller-Rabin is still probably the best way to go to do it "properly". If I re-wrote the material on the Perrin tests[0] I'd probably use M-R.
So we're pulling in the same direction, although making different trade-offs. And it's fun to play and see what happens.
[0] http://www.solipsys.co.uk/new/FindingPerrinPseudoPrimes_Part...