꽤나 큰 N이 소수 인지를 판별할 때 쓰는 방법입니다.
2부터 √N까지 나누어보는 방법 은 상당히 느리므로 밀러-라빈 소수판별법 을 기초로 합니다.
밀러-라빈 소수판별법은 입력된 수가 홀수일 때 깔끔하게 작동하니 2를 제외한 짝수는 미리 걸러주도록 합니다.
또, N이 정말 랜덤한 큰 수 라는 가정 하에 2, 3, 5...의 배수가 아닐 확률은 크지 않습니다.
예를 들어, 2의 배수가 아닐 확률은 12, 3의 배수가 아닐 확률은 23, 5의 배수가 아닐 확률은 45...인데, 이를 4041까지 다 곱해보면 19110297600131710070791=0.145093670..이 되는데, 이미 이 과정에서 소수 판별 대상이 대략 17로 줄어든 것입니다.
또, 페르마의 소정리 와 파이썬의 파이썬 pow 함수 를 이용해 pow(2, N - 1, N) != 1인지 확인하여 (2^{N - 1} \not\equiv 1 \pmod{N}인지 확인) 기타 합성수들을 한번에 걸러줍니다.
이미 이 과정에서 거의 대부분의 합성수들은 날아가고 남은 판별 대상은 진짜 소수들과 유사소수인 Poulet Number 밖에 남지 않습니다.
Poulet Number는 생각보다 많지 않기 때문에, 이 과정에서 대부분의 합성수는 날아갑니다.
풀렛 수들은 밀러-라빈 소수판별법으로 날려주되, 그 전에 발전된 프로트의 정리 를 사용해줍니다.
N이 프로트수인지 판별해주고, d \not\equiv 0 \pmod{3}이면 3^{\frac{N - 1}{2}} \equiv -1 \pmod{N}인지만 판단해도 충분합니다.
| import math | |
| mPL = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41] | |
| def MRTest(a, s, d, N): | |
| b = pow(a, d, N) # b0 구하기 | |
| if b == 1 or b == N - 1: return True | |
| for i in range(0, s - 1): | |
| b = pow(b, 2, N) # b1, b2, b3...b_{s-1} 구하기 | |
| if b == N - 1: return True | |
| return False # b1, b2, b3...b_{s-1} 까지 다 구했는데 전부 아니면 확실한 합성수 | |
| def isPrime(N): # 해설 : https://ii-eugene-ii.github.io/Post/Post01000/00010.html | |
| if N in mPL: return True | |
| for i in mPL: | |
| if N % i == 0: return False | |
| if N < 1849: return N != 1 # 43 미만의 소수는 이미 전처리 해주었으므로... | |
| if pow(2, N - 1, N) != 1: return False # 페르마의 소정리 (빠르게 확실한 합성수 정리) | |
| s, d = 0, N - 1 | |
| while d & 1 == 0: | |
| s += 1 | |
| d >>= 1 | |
| if d < (1 << s): # 발전된 프로트의 정리 | |
| if d % 3 != 0: return pow(3, (N - 1) // 2, N) == N - 1 | |
| elif d % 5 != 0 and N % 5 != 4: return pow(5, (N - 1) // 2, N) == N - 1 | |
| if N < 3215031751: # int형 범위 내 / 모두 참이여야 True 반환 | |
| if N < 2047: return MRTest(2, s, d, N) | |
| elif N < 1373653: return MRTest(2, s, d, N) and MRTest(3, s, d, N) | |
| elif N < 25326001: return MRTest(2, s, d, N) and MRTest(3, s, d, N) and MRTest(5, s, d, N) | |
| else: return MRTest(2, s, d, N) and MRTest(3, s, d, N) and MRTest(5, s, d, N) and MRTest(7, s, d, N) | |
| if N < 3317044064679887385961981: # unsigned long long int형 범위 내 | |
| for i in mPL: | |
| if not MRTest(i, s, d, N): return False | |
| return True | |
| ln2N = int(2 * (math.log(N) ** 2)) + 1 | |
| for i in range(2, ln2N): | |
| if not MRTest(i, s, d, N): return False | |
| return True |