Skip to content

Commit c20d822

Browse files
authored
Merge pull request kodecocodes#262 from scha00/master
Miller-Rabin primality test: probabilistic version
2 parents 10abdbc + a2fa24d commit c20d822

File tree

7 files changed

+313
-0
lines changed

7 files changed

+313
-0
lines changed
Loading
Lines changed: 24 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,24 @@
1+
//: Playground - noun: a place where people can play
2+
3+
// Real primes
4+
mrPrimalityTest(5)
5+
mrPrimalityTest(439)
6+
mrPrimalityTest(1201)
7+
mrPrimalityTest(143477)
8+
mrPrimalityTest(1299869)
9+
mrPrimalityTest(15487361)
10+
mrPrimalityTest(179426363)
11+
mrPrimalityTest(32416187747)
12+
13+
// Fake primes
14+
mrPrimalityTest(15)
15+
mrPrimalityTest(435)
16+
mrPrimalityTest(1207)
17+
mrPrimalityTest(143473)
18+
mrPrimalityTest(1291869)
19+
mrPrimalityTest(15487161)
20+
mrPrimalityTest(178426363)
21+
mrPrimalityTest(32415187747)
22+
23+
// With iteration
24+
mrPrimalityTest(32416190071, iteration: 10)
Lines changed: 111 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,111 @@
1+
//
2+
// MRPrimality.swift
3+
//
4+
//
5+
// Created by Sahn Cha on 2016. 10. 18..
6+
//
7+
//
8+
9+
import Foundation
10+
11+
/*
12+
Miller-Rabin Primality Test.
13+
One of the most used algorithms for primality testing.
14+
Since this is a probabilistic algorithm, it needs to be executed multiple times to increase accuracy.
15+
16+
Inputs:
17+
n: UInt64 { target integer to be tested for primality }
18+
k: Int { optional. number of iterations }
19+
20+
Outputs:
21+
true { probably prime }
22+
false { composite }
23+
*/
24+
public func mrPrimalityTest(_ n: UInt64, iteration k: Int = 1) -> Bool {
25+
guard n > 2 && n % 2 == 1 else { return false }
26+
27+
var d = n - 1
28+
var s = 0
29+
30+
while d % 2 == 0 {
31+
d /= 2
32+
s += 1
33+
}
34+
35+
let range = UInt64.max - UInt64.max % (n - 2)
36+
var r: UInt64 = 0
37+
repeat {
38+
arc4random_buf(&r, MemoryLayout.size(ofValue: r))
39+
} while r >= range
40+
41+
r = r % (n - 2) + 2
42+
43+
for _ in 1 ... k {
44+
var x = powmod64(r, d, n)
45+
if x == 1 || x == n - 1 { continue }
46+
47+
if s == 1 { s = 2 }
48+
49+
for _ in 1 ... s - 1 {
50+
x = powmod64(x, 2, n)
51+
if x == 1 { return false }
52+
if x == n - 1 { break }
53+
}
54+
55+
if x != n - 1 { return false }
56+
}
57+
58+
return true
59+
}
60+
61+
/*
62+
Calculates (base^exp) mod m.
63+
64+
Inputs:
65+
base: UInt64 { base }
66+
exp: UInt64 { exponent }
67+
m: UInt64 { modulus }
68+
69+
Outputs:
70+
the result
71+
*/
72+
private func powmod64(_ base: UInt64, _ exp: UInt64, _ m: UInt64) -> UInt64 {
73+
if m == 1 { return 0 }
74+
75+
var result: UInt64 = 1
76+
var b = base % m
77+
var e = exp
78+
79+
while e > 0 {
80+
if e % 2 == 1 { result = mulmod64(result, b, m) }
81+
b = mulmod64(b, b, m)
82+
e >>= 1
83+
}
84+
85+
return result
86+
}
87+
88+
/*
89+
Calculates (first * second) mod m, hopefully without overflow. :]
90+
91+
Inputs:
92+
first: UInt64 { first integer }
93+
second: UInt64 { second integer }
94+
m: UInt64 { modulus }
95+
96+
Outputs:
97+
the result
98+
*/
99+
private func mulmod64(_ first: UInt64, _ second: UInt64, _ m: UInt64) -> UInt64 {
100+
var result: UInt64 = 0
101+
var a = first
102+
var b = second
103+
104+
while a != 0 {
105+
if a % 2 == 1 { result = ((result % m) + (b % m)) % m } // This may overflow if 'm' is a 64bit number && both 'result' and 'b' are very close to but smaller than 'm'.
106+
a >>= 1
107+
b = (b << 1) % m
108+
}
109+
110+
return result
111+
}
Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,4 @@
1+
<?xml version="1.0" encoding="UTF-8" standalone="yes"?>
2+
<playground version='5.0' target-platform='macos'>
3+
<timeline fileName='timeline.xctimeline'/>
4+
</playground>
Lines changed: 111 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,111 @@
1+
//
2+
// MRPrimality.swift
3+
//
4+
//
5+
// Created by Sahn Cha on 2016. 10. 18..
6+
//
7+
//
8+
9+
import Foundation
10+
11+
/*
12+
Miller-Rabin Primality Test.
13+
One of the most used algorithms for primality testing.
14+
Since this is a probabilistic algorithm, it needs to be executed multiple times to increase accuracy.
15+
16+
Inputs:
17+
n: UInt64 { target integer to be tested for primality }
18+
k: Int { optional. number of iterations }
19+
20+
Outputs:
21+
true { probably prime }
22+
false { composite }
23+
*/
24+
public func mrPrimalityTest(_ n: UInt64, iteration k: Int = 1) -> Bool {
25+
guard n > 2 && n % 2 == 1 else { return false }
26+
27+
var d = n - 1
28+
var s = 0
29+
30+
while d % 2 == 0 {
31+
d /= 2
32+
s += 1
33+
}
34+
35+
let range = UInt64.max - UInt64.max % (n - 2)
36+
var r: UInt64 = 0
37+
repeat {
38+
arc4random_buf(&r, MemoryLayout.size(ofValue: r))
39+
} while r >= range
40+
41+
r = r % (n - 2) + 2
42+
43+
for _ in 1 ... k {
44+
var x = powmod64(r, d, n)
45+
if x == 1 || x == n - 1 { continue }
46+
47+
if s == 1 { s = 2 }
48+
49+
for _ in 1 ... s - 1 {
50+
x = powmod64(x, 2, n)
51+
if x == 1 { return false }
52+
if x == n - 1 { break }
53+
}
54+
55+
if x != n - 1 { return false }
56+
}
57+
58+
return true
59+
}
60+
61+
/*
62+
Calculates (base^exp) mod m.
63+
64+
Inputs:
65+
base: UInt64 { base }
66+
exp: UInt64 { exponent }
67+
m: UInt64 { modulus }
68+
69+
Outputs:
70+
the result
71+
*/
72+
private func powmod64(_ base: UInt64, _ exp: UInt64, _ m: UInt64) -> UInt64 {
73+
if m == 1 { return 0 }
74+
75+
var result: UInt64 = 1
76+
var b = base % m
77+
var e = exp
78+
79+
while e > 0 {
80+
if e % 2 == 1 { result = mulmod64(result, b, m) }
81+
b = mulmod64(b, b, m)
82+
e >>= 1
83+
}
84+
85+
return result
86+
}
87+
88+
/*
89+
Calculates (first * second) mod m, hopefully without overflow. :]
90+
91+
Inputs:
92+
first: UInt64 { first integer }
93+
second: UInt64 { second integer }
94+
m: UInt64 { modulus }
95+
96+
Outputs:
97+
the result
98+
*/
99+
private func mulmod64(_ first: UInt64, _ second: UInt64, _ m: UInt64) -> UInt64 {
100+
var result: UInt64 = 0
101+
var a = first
102+
var b = second
103+
104+
while a != 0 {
105+
if a % 2 == 1 { result = ((result % m) + (b % m)) % m } // This may overflow if 'm' is a 64bit number && both 'result' and 'b' are very close to but smaller than 'm'.
106+
a >>= 1
107+
b = (b << 1) % m
108+
}
109+
110+
return result
111+
}
Lines changed: 43 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,43 @@
1+
# Miller-Rabin Primality Test
2+
3+
In 1976, Gray Miller introduced an algorithm, through his ph.d thesis[1], which determines a primality of the given number. The original algorithm was deterministic under the Extended Reimann Hypothesis, which is yet to be proven. After four years, Michael O. Rabin improved the algorithm[2] by using probabilistic approach and it no longer assumes the unproven hypothesis.
4+
5+
## Probabilistic
6+
7+
The result of the test is simply a boolean. However, `true` does not implicate _the number is prime_. In fact, it means _the number is **probably** prime_. But `false` does mean _the number is composite_.
8+
9+
In order to increase the accuracy of the test, it needs to be iterated few times. If it returns `true` in every single iteration, then we can say with confidence that _the number is pro......bably prime_.
10+
11+
## Algorithm
12+
13+
Let `n` be the given number, and write `n-1` as `2^s·d`, where `d` is odd. And choose a random number `a` within the range from `2` to `n - 1`.
14+
15+
Now make a sequence, in modulo `n`, as following:
16+
17+
> a^d, a^(2·d), a^(4·d), ... , a^((2^(s-1))·d), a^((2^s)·d) = a^(n-1)
18+
19+
And we say the number `n` passes the test, _probably prime_, if 1) `a^d` is congruence to `1` in modulo `n`, or 2) `a^((2^k)·d)` is congruence to `-1` for some `k = 1, 2, ..., s-1`.
20+
21+
### Pseudo Code
22+
23+
The following pseudo code is excerpted from Wikipedia[3]:
24+
25+
![Image of Pseudocode](./Images/img_pseudo.png)
26+
27+
## Usage
28+
29+
```swift
30+
mrPrimalityTest(7) // test if 7 is prime. (default iteration = 1)
31+
mrPrimalityTest(7, iteration: 10) // test if 7 is prime && iterate 10 times.
32+
```
33+
34+
## Reference
35+
1. G. L. Miller, "Riemann's Hypothesis and Tests for Primality". _J. Comput. System Sci._ 13 (1976), 300-317.
36+
2. M. O. Rabin, "Probabilistic algorithm for testing primality". _Journal of Number Theory._ 12 (1980), 128-138.
37+
3. Miller–Rabin primality test - Wikipedia, the free encyclopedia
38+
39+
_Written for Swift Algorithm Club by **Sahn Cha**, @scha00_
40+
41+
[1]: https://cs.uwaterloo.ca/research/tr/1975/CS-75-27.pdf
42+
[2]: http://www.sciencedirect.com/science/article/pii/0022314X80900840
43+
[3]: https://en.wikipedia.org/wiki/Miller–Rabin_primality_test

swift-algorithm-club.xcworkspace/contents.xcworkspacedata

Lines changed: 20 additions & 0 deletions
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

0 commit comments

Comments
 (0)