Miller–Rabin primality test: Difference between revisions
Content added Content deleted
Thundergnat (talk | contribs) m (syntax highlighting fixup automation) |
(New post.) |
||
Line 903: | Line 903: | ||
} |
} |
||
}</syntaxhighlight> |
}</syntaxhighlight> |
||
=={{header|C++}}== |
|||
This test is deterministic for n < 3'474'749'660'383 |
|||
<syntaxhighlight lang="c++"> |
|||
#include <algorithm> |
|||
#include <cmath> |
|||
#include <cstdint> |
|||
#include <iostream> |
|||
#include <vector> |
|||
std::vector<uint32_t> small_primes{2, 3}; |
|||
uint64_t add_modulus(const uint64_t& a, const uint64_t& b, const uint64_t& modulus) { |
|||
uint64_t am = ( a < modulus ) ? a : a % modulus; |
|||
if ( b == 0 ) { |
|||
return am; |
|||
} |
|||
uint64_t bm = ( b < modulus ) ? b : b % modulus; |
|||
uint64_t b_from_m = modulus - bm; |
|||
if ( am >= b_from_m ) { |
|||
return am - b_from_m; |
|||
} |
|||
return am + bm; |
|||
} |
|||
uint64_t multiply_modulus(const uint64_t& a, const uint64_t& b, const uint64_t& modulus) { |
|||
uint64_t am = ( a < modulus ) ? a : a % modulus; |
|||
uint64_t bm = ( b < modulus ) ? b : b % modulus; |
|||
if ( bm > am ) { |
|||
std::swap(am, bm); |
|||
} |
|||
uint64_t result; |
|||
while ( bm > 0 ) { |
|||
if ( ( bm & 1) == 1 ) { |
|||
result = add_modulus(result, am, modulus); |
|||
} |
|||
am = ( am << 1 ) - ( am >= ( modulus - am ) ? modulus : 0 ); |
|||
bm >>= 1; |
|||
} |
|||
return result; |
|||
} |
|||
uint64_t exponentiation_modulus(const uint64_t& base, const uint64_t& exponent, const uint64_t& modulus) { |
|||
uint64_t b = base; |
|||
uint64_t e = exponent; |
|||
uint64_t result = 1; |
|||
while ( e > 0 ) { |
|||
if ( ( e & 1 ) == 1 ) { |
|||
result = multiply_modulus(result, b, modulus); |
|||
} |
|||
e >>= 1; |
|||
b = multiply_modulus(b, b, modulus); |
|||
} |
|||
return result; |
|||
} |
|||
bool is_composite(const uint32_t& a, const uint64_t& d, const uint64_t& n, const uint32_t& s) { |
|||
if ( exponentiation_modulus(a, d, n) == 1 ) { |
|||
return false; |
|||
} |
|||
for ( uint64_t i = 0; i < s; ++i ) { |
|||
if ( exponentiation_modulus(a, pow(2, i) * d, n) == n - 1 ) { |
|||
return false; |
|||
} |
|||
} |
|||
return true; |
|||
} |
|||
bool composite_test(const std::vector<uint32_t>& primes, const uint64_t& d, const uint64_t& n, const uint32_t& s) { |
|||
for ( const uint32_t& prime : primes ) { |
|||
if ( is_composite(prime, d, n, s) ) { |
|||
return true; |
|||
} |
|||
} |
|||
return false; |
|||
} |
|||
bool is_prime(const uint64_t& n) { |
|||
if ( n == 0 || n == 1 ) { |
|||
return false; |
|||
} |
|||
if ( std::find(small_primes.begin(), small_primes.end(), n) != small_primes.end() ) { |
|||
return true; |
|||
} |
|||
if ( std::any_of(small_primes.begin(), small_primes.end(), [n](uint32_t p) { return n % p == 0; }) ) { |
|||
return false; |
|||
} |
|||
uint64_t d = n - 1; |
|||
uint32_t s = 0; |
|||
while ( ! d % 2 ) { |
|||
d >>= 1; |
|||
s++; |
|||
} |
|||
if ( n < 1'373'653 ) { |
|||
return composite_test({ 2, 3 }, d, n, s); |
|||
} |
|||
if ( n < 25'326'001 ) { |
|||
return composite_test({ 2, 3, 5 }, d, n, s); |
|||
} |
|||
if ( n < 118'670'087'467 ) { |
|||
if ( n == 3'215'031'751 ) { |
|||
return false; |
|||
} |
|||
return composite_test({ 2, 3, 5, 7 }, d, n, s); |
|||
} |
|||
if ( n < 2'152'302'898'747 ) { |
|||
return composite_test({ 2, 3, 5, 7, 11 }, d, n, s); |
|||
} |
|||
if ( n < 3'474'749'660'383 ) { |
|||
return composite_test({ 2, 3, 5, 7, 11, 13 }, d, n, s); |
|||
} |
|||
if ( n < 341'550'071'728'321 ) { |
|||
return composite_test({ 2, 3, 5, 7, 11, 13, 17 }, d, n, s); |
|||
} |
|||
const std::vector<uint32_t> test_primes(small_primes.begin(), small_primes.begin() + 16); |
|||
return composite_test(test_primes, d, n, s); |
|||
} |
|||
void create_small_primes() { |
|||
for ( uint32_t i = 5; i < 1'000; i += 2 ) { |
|||
if ( is_prime(i) ) { |
|||
small_primes.emplace_back(i); |
|||
} |
|||
} |
|||
} |
|||
int main() { |
|||
create_small_primes(); |
|||
for ( const uint64_t number : { 1'234'567'890'123'456'733, 1'234'567'890'123'456'737 } ) { |
|||
std::cout << "is_prime(" << number << ") = " << std::boolalpha << is_prime(number) << std::endl; |
|||
} |
|||
} |
|||
</syntaxhighlight> |
|||
{{ out }} |
|||
<pre> |
|||
is_prime(1234567890123456733) = false |
|||
is_prime(1234567890123456737) = true |
|||
</pre> |
|||
=={{header|Clojure}}== |
=={{header|Clojure}}== |