001: /*
002: * Licensed to the Apache Software Foundation (ASF) under one or more
003: * contributor license agreements. See the NOTICE file distributed with
004: * this work for additional information regarding copyright ownership.
005: * The ASF licenses this file to You under the Apache License, Version 2.0
006: * (the "License"); you may not use this file except in compliance with
007: * the License. You may obtain a copy of the License at
008: *
009: * http://www.apache.org/licenses/LICENSE-2.0
010: *
011: * Unless required by applicable law or agreed to in writing, software
012: * distributed under the License is distributed on an "AS IS" BASIS,
013: * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014: * See the License for the specific language governing permissions and
015: * limitations under the License.
016: */
017:
018: package java.math;
019:
020: import java.util.Arrays;
021: import java.util.Random;
022:
023: /**
024: * Provides primality probabilistic methods.
025: * @author Intel Middleware Product Division
026: * @author Instituto Tecnologico de Cordoba
027: */
028: class Primality {
029:
030: /** Just to denote that this class can't be instantiated. */
031: private Primality() {
032: }
033:
034: /* Private Fields */
035:
036: /** All prime numbers with bit length lesser than 10 bits. */
037: private static final int primes[] = { 2, 3, 5, 7, 11, 13, 17, 19,
038: 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83,
039: 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149,
040: 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211,
041: 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277,
042: 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353,
043: 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431,
044: 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499,
045: 503, 509, 521, 523, 541, 547, 557, 563, 569, 571, 577, 587,
046: 593, 599, 601, 607, 613, 617, 619, 631, 641, 643, 647, 653,
047: 659, 661, 673, 677, 683, 691, 701, 709, 719, 727, 733, 739,
048: 743, 751, 757, 761, 769, 773, 787, 797, 809, 811, 821, 823,
049: 827, 829, 839, 853, 857, 859, 863, 877, 881, 883, 887, 907,
050: 911, 919, 929, 937, 941, 947, 953, 967, 971, 977, 983, 991,
051: 997, 1009, 1013, 1019, 1021 };
052:
053: /** All {@code BigInteger} prime numbers with bit length lesser than 8 bits. */
054: private static final BigInteger BIprimes[] = new BigInteger[primes.length];
055:
056: /**
057: * It encodes how many iterations of Miller-Rabin test are need to get an
058: * error bound not greater than {@code 2<sup>(-100)</sup>}. For example:
059: * for a {@code 1000}-bit number we need {@code 4} iterations, since
060: * {@code BITS[3] < 1000 <= BITS[4]}.
061: */
062: private static final int[] BITS = { 0, 0, 1854, 1233, 927, 747,
063: 627, 543, 480, 431, 393, 361, 335, 314, 295, 279, 265, 253,
064: 242, 232, 223, 216, 181, 169, 158, 150, 145, 140, 136, 132,
065: 127, 123, 119, 114, 110, 105, 101, 96, 92, 87, 83, 78, 73,
066: 69, 64, 59, 54, 49, 44, 38, 32, 26, 1 };
067:
068: /**
069: * It encodes how many i-bit primes there are in the table for
070: * {@code i=2,...,10}. For example {@code offsetPrimes[6]} says that from
071: * index {@code 11} exists {@code 7} consecutive {@code 6}-bit prime
072: * numbers in the array.
073: */
074: private static final int[][] offsetPrimes = { null, null, { 0, 2 },
075: { 2, 2 }, { 4, 2 }, { 6, 5 }, { 11, 7 }, { 18, 13 },
076: { 31, 23 }, { 54, 43 }, { 97, 75 } };
077:
078: static {// To initialize the dual table of BigInteger primes
079: for (int i = 0; i < primes.length; i++) {
080: BIprimes[i] = BigInteger.valueOf(primes[i]);
081: }
082: }
083:
084: /* Package Methods */
085:
086: /**
087: * It uses the sieve of Eratosthenes to discard several composite numbers in
088: * some appropriate range (at the moment {@code [this, this + 1024]}). After
089: * this process it applies the Miller-Rabin test to the numbers that were
090: * not discarded in the sieve.
091: *
092: * @see BigInteger#nextProbablePrime()
093: * @see #millerRabin(BigInteger, int)
094: */
095: static BigInteger nextProbablePrime(BigInteger n) {
096: // PRE: n >= 0
097: int i, j;
098: int certainty;
099: int gapSize = 1024; // for searching of the next probable prime number
100: int modules[] = new int[primes.length];
101: boolean isDivisible[] = new boolean[gapSize];
102: BigInteger startPoint;
103: BigInteger probPrime;
104: // If n < "last prime of table" searches next prime in the table
105: if ((n.numberLength == 1) && (n.digits[0] >= 0)
106: && (n.digits[0] < primes[primes.length - 1])) {
107: for (i = 0; n.digits[0] >= primes[i]; i++) {
108: ;
109: }
110: return BIprimes[i];
111: }
112: /*
113: * Creates a "N" enough big to hold the next probable prime Note that: N <
114: * "next prime" < 2*N
115: */
116: startPoint = new BigInteger(1, n.numberLength,
117: new int[n.numberLength + 1]);
118: System.arraycopy(n.digits, 0, startPoint.digits, 0,
119: n.numberLength);
120: // To fix N to the "next odd number"
121: if (n.testBit(0)) {
122: Elementary.inplaceAdd(startPoint, 2);
123: } else {
124: startPoint.digits[0] |= 1;
125: }
126: // To set the improved certainly of Miller-Rabin
127: j = startPoint.bitLength();
128: for (certainty = 2; j < BITS[certainty]; certainty++) {
129: ;
130: }
131: // To calculate modules: N mod p1, N mod p2, ... for first primes.
132: for (i = 0; i < primes.length; i++) {
133: modules[i] = Division.remainder(startPoint, primes[i])
134: - gapSize;
135: }
136: while (true) {
137: // At this point, all numbers in the gap are initialized as
138: // probably primes
139: Arrays.fill(isDivisible, false);
140: // To discard multiples of first primes
141: for (i = 0; i < primes.length; i++) {
142: modules[i] = (modules[i] + gapSize) % primes[i];
143: j = (modules[i] == 0) ? 0 : (primes[i] - modules[i]);
144: for (; j < gapSize; j += primes[i]) {
145: isDivisible[j] = true;
146: }
147: }
148: // To execute Miller-Rabin for non-divisible numbers by all first
149: // primes
150: for (j = 0; j < gapSize; j++) {
151: if (!isDivisible[j]) {
152: probPrime = startPoint.copy();
153: Elementary.inplaceAdd(probPrime, j);
154:
155: if (millerRabin(probPrime, certainty)) {
156: return probPrime;
157: }
158: }
159: }
160: Elementary.inplaceAdd(startPoint, gapSize);
161: }
162: }
163:
164: /**
165: * A random number is generated until a probable prime number is found.
166: *
167: * @see BigInteger#BigInteger(int,int,Random)
168: * @see BigInteger#probablePrime(int,Random)
169: * @see #isProbablePrime(BigInteger, int)
170: */
171: static BigInteger consBigInteger(int bitLength, int certainty,
172: Random rnd) {
173: // PRE: bitLength >= 2;
174: // For small numbers get a random prime from the prime table
175: if (bitLength <= 10) {
176: int rp[] = offsetPrimes[bitLength];
177: return BIprimes[rp[0] + rnd.nextInt(rp[1])];
178: }
179: int shiftCount = (-bitLength) & 31;
180: int last = (bitLength + 31) >> 5;
181: BigInteger n = new BigInteger(1, last, new int[last]);
182:
183: last--;
184: do {// To fill the array with random integers
185: for (int i = 0; i < n.numberLength; i++) {
186: n.digits[i] = rnd.nextInt();
187: }
188: // To fix to the correct bitLength
189: n.digits[last] |= 0x80000000;
190: n.digits[last] >>>= shiftCount;
191: // To create an odd number
192: n.digits[0] |= 1;
193: } while (!isProbablePrime(n, certainty));
194: return n;
195: }
196:
197: /**
198: * @see BigInteger#isProbablePrime(int)
199: * @see #millerRabin(BigInteger, int)
200: * @ar.org.fitc.ref Optimizations: "A. Menezes - Handbook of applied
201: * Cryptography, Chapter 4".
202: */
203: static boolean isProbablePrime(BigInteger n, int certainty) {
204: // PRE: n >= 0;
205: if ((certainty <= 0)
206: || ((n.numberLength == 1) && (n.digits[0] == 2))) {
207: return true;
208: }
209: // To discard all even numbers
210: if (!n.testBit(0)) {
211: return false;
212: }
213: // To check if 'n' exists in the table (it fit in 10 bits)
214: if ((n.numberLength == 1) && ((n.digits[0] & 0XFFFFFC00) == 0)) {
215: return (Arrays.binarySearch(primes, n.digits[0]) >= 0);
216: }
217: // To check if 'n' is divisible by some prime of the table
218: for (int i = 1; i < primes.length; i++) {
219: if (Division.remainderArrayByInt(n.digits, n.numberLength,
220: primes[i]) == 0) {
221: return false;
222: }
223: }
224: // To set the number of iterations necessary for Miller-Rabin test
225: int i;
226: int bitLength = n.bitLength();
227:
228: for (i = 2; bitLength < BITS[i]; i++) {
229: ;
230: }
231: certainty = Math.min(i, 1 + ((certainty - 1) >> 1));
232:
233: return millerRabin(n, certainty);
234: }
235:
236: /* Private Methods */
237:
238: /**
239: * The Miller-Rabin primality test.
240: *
241: * @param n the input number to be tested.
242: * @param t the number of trials.
243: * @return {@code false} if the number is definitely compose, otherwise
244: * {@code true} with probability {@code 1 - 4<sup>(-t)</sup>}.
245: * @ar.org.fitc.ref "D. Knuth, The Art of Computer Programming Vo.2, Section
246: * 4.5.4., Algorithm P"
247: */
248: private static boolean millerRabin(BigInteger n, int t) {
249: // PRE: n >= 0, t >= 0
250: BigInteger x; // x := UNIFORM{2...n-1}
251: BigInteger y; // y := x^(q * 2^j) mod n
252: BigInteger n_minus_1 = n.subtract(BigInteger.ONE); // n-1
253: int bitLength = n_minus_1.bitLength(); // ~ log2(n-1)
254: // (q,k) such that: n-1 = q * 2^k and q is odd
255: int k = n_minus_1.getLowestSetBit();
256: BigInteger q = n_minus_1.shiftRight(k);
257: Random rnd = new Random();
258:
259: for (int i = 0; i < t; i++) {
260: // To generate a witness 'x', first it use the primes of table
261: if (i < primes.length) {
262: x = BIprimes[i];
263: } else {/*
264: * It generates random witness only if it's necesssary. Note
265: * that all methods would call Miller-Rabin with t <= 50 so
266: * this part is only to do more robust the algorithm
267: */
268: do {
269: x = new BigInteger(bitLength, rnd);
270: } while ((x.compareTo(n) >= BigInteger.EQUALS)
271: || (x.sign == 0) || x.isOne());
272: }
273: y = x.modPow(q, n);
274: if (y.isOne() || y.equals(n_minus_1)) {
275: continue;
276: }
277: for (int j = 1; j < k; j++) {
278: if (y.equals(n_minus_1)) {
279: continue;
280: }
281: y = y.multiply(y).mod(n);
282: if (y.isOne()) {
283: return false;
284: }
285: }
286: if (!y.equals(n_minus_1)) {
287: return false;
288: }
289: }
290: return true;
291: }
292:
293: }
|