blob: b2ddc744a73edebfe57acdef67295547ec3a3805 (
plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
|
// <auto-generated />
//
// Mono.Math.Prime.PrimalityTests.cs - Test for primality
//
// Authors:
// Ben Maurer
//
// Copyright (c) 2003 Ben Maurer. All rights reserved
//
using System;
using System.Security.Cryptography;
namespace Mono.Math.Prime {
//[CLSCompliant(false)]
internal delegate bool PrimalityTest (BigInteger bi, ConfidenceFactor confidence);
//[CLSCompliant(false)]
internal sealed class PrimalityTests {
#region SPP Test
private static int GetSPPRounds (BigInteger bi, ConfidenceFactor confidence)
{
int bc = bi.bitCount();
int Rounds;
// Data from HAC, 4.49
if (bc <= 100 ) Rounds = 27;
else if (bc <= 150 ) Rounds = 18;
else if (bc <= 200 ) Rounds = 15;
else if (bc <= 250 ) Rounds = 12;
else if (bc <= 300 ) Rounds = 9;
else if (bc <= 350 ) Rounds = 8;
else if (bc <= 400 ) Rounds = 7;
else if (bc <= 500 ) Rounds = 6;
else if (bc <= 600 ) Rounds = 5;
else if (bc <= 800 ) Rounds = 4;
else if (bc <= 1250) Rounds = 3;
else Rounds = 2;
switch (confidence) {
case ConfidenceFactor.ExtraLow:
Rounds >>= 2;
return Rounds != 0 ? Rounds : 1;
case ConfidenceFactor.Low:
Rounds >>= 1;
return Rounds != 0 ? Rounds : 1;
case ConfidenceFactor.Medium:
return Rounds;
case ConfidenceFactor.High:
return Rounds <<= 1;
case ConfidenceFactor.ExtraHigh:
return Rounds <<= 2;
case ConfidenceFactor.Provable:
throw new Exception ("The Rabin-Miller test can not be executed in a way such that its results are provable");
default:
throw new ArgumentOutOfRangeException ("confidence");
}
}
/// <summary>
/// Probabilistic prime test based on Rabin-Miller's test
/// </summary>
/// <param name="bi" type="BigInteger.BigInteger">
/// <para>
/// The number to test.
/// </para>
/// </param>
/// <param name="confidence" type="int">
/// <para>
/// The number of chosen bases. The test has at least a
/// 1/4^confidence chance of falsely returning True.
/// </para>
/// </param>
/// <returns>
/// <para>
/// True if "this" is a strong pseudoprime to randomly chosen bases.
/// </para>
/// <para>
/// False if "this" is definitely NOT prime.
/// </para>
/// </returns>
public static bool RabinMillerTest (BigInteger bi, ConfidenceFactor confidence)
{
int Rounds = GetSPPRounds (bi, confidence);
// calculate values of s and t
BigInteger p_sub1 = bi - 1;
int s = p_sub1.LowestSetBit ();
BigInteger t = p_sub1 >> s;
int bits = bi.bitCount ();
BigInteger a = null;
RandomNumberGenerator rng = RandomNumberGenerator.Create ();
BigInteger.ModulusRing mr = new BigInteger.ModulusRing (bi);
for (int round = 0; round < Rounds; round++) {
while (true) { // generate a < n
a = BigInteger.genRandom (bits, rng);
// make sure "a" is not 0
if (a > 1 && a < bi)
break;
}
if (a.gcd (bi) != 1) return false;
BigInteger b = mr.Pow (a, t);
if (b == 1) continue; // a^t mod p = 1
bool result = false;
for (int j = 0; j < s; j++) {
if (b == p_sub1) { // a^((2^j)*t) mod p = p-1 for some 0 <= j <= s-1
result = true;
break;
}
b = (b * b) % bi;
}
if (result == false)
return false;
}
return true;
}
public static bool SmallPrimeSppTest (BigInteger bi, ConfidenceFactor confidence)
{
int Rounds = GetSPPRounds (bi, confidence);
// calculate values of s and t
BigInteger p_sub1 = bi - 1;
int s = p_sub1.LowestSetBit ();
BigInteger t = p_sub1 >> s;
BigInteger.ModulusRing mr = new BigInteger.ModulusRing (bi);
for (int round = 0; round < Rounds; round++) {
BigInteger b = mr.Pow (BigInteger.smallPrimes [round], t);
if (b == 1) continue; // a^t mod p = 1
bool result = false;
for (int j = 0; j < s; j++) {
if (b == p_sub1) { // a^((2^j)*t) mod p = p-1 for some 0 <= j <= s-1
result = true;
break;
}
b = (b * b) % bi;
}
if (result == false)
return false;
}
return true;
}
#endregion
// TODO: Implement the Lucus test
// TODO: Implement other new primality tests
// TODO: Implement primality proving
}
}
|