1 /*
2 ==============================================================================
3
4 This file is part of the JUCE library.
5 Copyright (c) 2020 - Raw Material Software Limited
6
7 JUCE is an open source library subject to commercial or open-source
8 licensing.
9
10 By using JUCE, you agree to the terms of both the JUCE 6 End-User License
11 Agreement and JUCE Privacy Policy (both effective as of the 16th June 2020).
12
13 End User License Agreement: www.juce.com/juce-6-licence
14 Privacy Policy: www.juce.com/juce-privacy-policy
15
16 Or: You may also use this code under the terms of the GPL v3 (see
17 www.gnu.org/licenses).
18
19 JUCE IS PROVIDED "AS IS" WITHOUT ANY WARRANTY, AND ALL WARRANTIES, WHETHER
20 EXPRESSED OR IMPLIED, INCLUDING MERCHANTABILITY AND FITNESS FOR PURPOSE, ARE
21 DISCLAIMED.
22
23 ==============================================================================
24 */
25
26 namespace juce
27 {
28
29 namespace PrimesHelpers
30 {
createSmallSieve(const int numBits,BigInteger & result)31 static void createSmallSieve (const int numBits, BigInteger& result)
32 {
33 result.setBit (numBits);
34 result.clearBit (numBits); // to enlarge the array
35
36 result.setBit (0);
37 int n = 2;
38
39 do
40 {
41 for (int i = n + n; i < numBits; i += n)
42 result.setBit (i);
43
44 n = result.findNextClearBit (n + 1);
45 }
46 while (n <= (numBits >> 1));
47 }
48
bigSieve(const BigInteger & base,const int numBits,BigInteger & result,const BigInteger & smallSieve,const int smallSieveSize)49 static void bigSieve (const BigInteger& base, const int numBits, BigInteger& result,
50 const BigInteger& smallSieve, const int smallSieveSize)
51 {
52 jassert (! base[0]); // must be even!
53
54 result.setBit (numBits);
55 result.clearBit (numBits); // to enlarge the array
56
57 int index = smallSieve.findNextClearBit (0);
58
59 do
60 {
61 const unsigned int prime = ((unsigned int) index << 1) + 1;
62
63 BigInteger r (base), remainder;
64 r.divideBy (prime, remainder);
65
66 unsigned int i = prime - remainder.getBitRangeAsInt (0, 32);
67
68 if (r.isZero())
69 i += prime;
70
71 if ((i & 1) == 0)
72 i += prime;
73
74 i = (i - 1) >> 1;
75
76 while (i < (unsigned int) numBits)
77 {
78 result.setBit ((int) i);
79 i += prime;
80 }
81
82 index = smallSieve.findNextClearBit (index + 1);
83 }
84 while (index < smallSieveSize);
85 }
86
findCandidate(const BigInteger & base,const BigInteger & sieve,const int numBits,BigInteger & result,const int certainty)87 static bool findCandidate (const BigInteger& base, const BigInteger& sieve,
88 const int numBits, BigInteger& result, const int certainty)
89 {
90 for (int i = 0; i < numBits; ++i)
91 {
92 if (! sieve[i])
93 {
94 result = base + (unsigned int) ((i << 1) + 1);
95
96 if (Primes::isProbablyPrime (result, certainty))
97 return true;
98 }
99 }
100
101 return false;
102 }
103
passesMillerRabin(const BigInteger & n,int iterations)104 static bool passesMillerRabin (const BigInteger& n, int iterations)
105 {
106 const BigInteger one (1), two (2);
107 const BigInteger nMinusOne (n - one);
108
109 BigInteger d (nMinusOne);
110 const int s = d.findNextSetBit (0);
111 d >>= s;
112
113 BigInteger smallPrimes;
114 int numBitsInSmallPrimes = 0;
115
116 for (;;)
117 {
118 numBitsInSmallPrimes += 256;
119 createSmallSieve (numBitsInSmallPrimes, smallPrimes);
120
121 const int numPrimesFound = numBitsInSmallPrimes - smallPrimes.countNumberOfSetBits();
122
123 if (numPrimesFound > iterations + 1)
124 break;
125 }
126
127 int smallPrime = 2;
128
129 while (--iterations >= 0)
130 {
131 smallPrime = smallPrimes.findNextClearBit (smallPrime + 1);
132
133 BigInteger r (smallPrime);
134 r.exponentModulo (d, n);
135
136 if (r != one && r != nMinusOne)
137 {
138 for (int j = 0; j < s; ++j)
139 {
140 r.exponentModulo (two, n);
141
142 if (r == nMinusOne)
143 break;
144 }
145
146 if (r != nMinusOne)
147 return false;
148 }
149 }
150
151 return true;
152 }
153 }
154
155 //==============================================================================
createProbablePrime(const int bitLength,const int certainty,const int * randomSeeds,int numRandomSeeds)156 BigInteger Primes::createProbablePrime (const int bitLength,
157 const int certainty,
158 const int* randomSeeds,
159 int numRandomSeeds)
160 {
161 using namespace PrimesHelpers;
162 int defaultSeeds [16];
163
164 if (numRandomSeeds <= 0)
165 {
166 randomSeeds = defaultSeeds;
167 numRandomSeeds = numElementsInArray (defaultSeeds);
168 Random r1, r2;
169
170 for (int j = 10; --j >= 0;)
171 {
172 r1.setSeedRandomly();
173
174 for (int i = numRandomSeeds; --i >= 0;)
175 defaultSeeds[i] ^= r1.nextInt() ^ r2.nextInt();
176 }
177 }
178
179 BigInteger smallSieve;
180 const int smallSieveSize = 15000;
181 createSmallSieve (smallSieveSize, smallSieve);
182
183 BigInteger p;
184
185 for (int i = numRandomSeeds; --i >= 0;)
186 {
187 BigInteger p2;
188
189 Random r (randomSeeds[i]);
190 r.fillBitsRandomly (p2, 0, bitLength);
191
192 p ^= p2;
193 }
194
195 p.setBit (bitLength - 1);
196 p.clearBit (0);
197
198 const int searchLen = jmax (1024, (bitLength / 20) * 64);
199
200 while (p.getHighestBit() < bitLength)
201 {
202 p += 2 * searchLen;
203
204 BigInteger sieve;
205 bigSieve (p, searchLen, sieve,
206 smallSieve, smallSieveSize);
207
208 BigInteger candidate;
209
210 if (findCandidate (p, sieve, searchLen, candidate, certainty))
211 return candidate;
212 }
213
214 jassertfalse;
215 return BigInteger();
216 }
217
isProbablyPrime(const BigInteger & number,const int certainty)218 bool Primes::isProbablyPrime (const BigInteger& number, const int certainty)
219 {
220 using namespace PrimesHelpers;
221
222 if (! number[0])
223 return false;
224
225 if (number.getHighestBit() <= 10)
226 {
227 const unsigned int num = number.getBitRangeAsInt (0, 10);
228
229 for (unsigned int i = num / 2; --i > 1;)
230 if (num % i == 0)
231 return false;
232
233 return true;
234 }
235 else
236 {
237 if (number.findGreatestCommonDivisor (2 * 3 * 5 * 7 * 11 * 13 * 17 * 19 * 23) != 1)
238 return false;
239
240 return passesMillerRabin (number, certainty);
241 }
242 }
243
244 } // namespace juce
245