179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giropackage org.bouncycastle.math;
279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giroimport java.math.BigInteger;
479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giroimport java.security.SecureRandom;
579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giroimport org.bouncycastle.crypto.Digest;
779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giroimport org.bouncycastle.util.Arrays;
879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giroimport org.bouncycastle.util.BigIntegers;
979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
1079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro/**
1179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro * Utility methods for generating primes and testing for primality.
1279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro */
1379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giropublic abstract class Primes
1479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro{
1579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    public static final int SMALL_FACTOR_LIMIT = 211;
1679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
1779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static final BigInteger ONE = BigInteger.valueOf(1);
1879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static final BigInteger TWO = BigInteger.valueOf(2);
1979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static final BigInteger THREE = BigInteger.valueOf(3);
2079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
2179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    /**
2279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * Used to return the output from the
2379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * {@linkplain Primes#enhancedMRProbablePrimeTest(BigInteger, SecureRandom, int) Enhanced
2479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * Miller-Rabin Probabilistic Primality Test}
2579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     */
2679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    public static class MROutput
2779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
2879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private static MROutput probablyPrime()
2979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
3079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return new MROutput(false, null);
3179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
3279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
3379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private static MROutput provablyCompositeWithFactor(BigInteger factor)
3479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
3579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return new MROutput(true, factor);
3679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
3779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
3879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private static MROutput provablyCompositeNotPrimePower()
3979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
4079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return new MROutput(true, null);
4179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
4279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
4379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private boolean provablyComposite;
4479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private BigInteger factor;
4579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
4679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private MROutput(boolean provablyComposite, BigInteger factor)
4779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
4879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            this.provablyComposite = provablyComposite;
4979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            this.factor = factor;
5079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
5179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
5279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        public BigInteger getFactor()
5379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
5479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return factor;
5579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
5679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
5779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        public boolean isProvablyComposite()
5879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
5979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return provablyComposite;
6079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
6179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
6279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        public boolean isNotPrimePower()
6379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
6479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return provablyComposite && factor == null;
6579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
6679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
6779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
6879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    /**
6979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * Used to return the output from the
7079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * {@linkplain Primes#generateSTRandomPrime(Digest, int, byte[]) Shawe-Taylor Random_Prime
7179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * Routine}
7279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     */
7379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    public static class STOutput
7479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
7579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private BigInteger prime;
7679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private byte[] primeSeed;
7779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private int primeGenCounter;
7879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
7979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        private STOutput(BigInteger prime, byte[] primeSeed, int primeGenCounter)
8079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
8179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            this.prime = prime;
8279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            this.primeSeed = primeSeed;
8379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            this.primeGenCounter = primeGenCounter;
8479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
8579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
8679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        public BigInteger getPrime()
8779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
8879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return prime;
8979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
9079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
9179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        public byte[] getPrimeSeed()
9279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
9379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return primeSeed;
9479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
9579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
9679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        public int getPrimeGenCounter()
9779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
9879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return primeGenCounter;
9979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
10079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
10179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
10279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    /**
10379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * FIPS 186-4 C.6 Shawe-Taylor Random_Prime Routine
10479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
10579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * Construct a provable prime number using a hash function.
10679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
10779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param hash
10879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the {@link Digest} instance to use (as "Hash()"). Cannot be null.
10979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param length
11079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the length (in bits) of the prime to be generated. Must be at least 2.
11179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param inputSeed
11279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the seed to be used for the generation of the requested prime. Cannot be null or
11379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            empty.
11479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @return an {@link STOutput} instance containing the requested prime.
11579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     */
11679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    public static STOutput generateSTRandomPrime(Digest hash, int length, byte[] inputSeed)
11779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
11879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (hash == null)
11979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
12079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("'hash' cannot be null");
12179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
12279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (length < 2)
12379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
12479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("'length' must be >= 2");
12579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
12679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (inputSeed == null || inputSeed.length == 0)
12779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
12879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("'inputSeed' cannot be null or empty");
12979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
13079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
13179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        return implSTRandomPrime(hash, length, Arrays.clone(inputSeed));
13279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
13379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
13479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    /**
13579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * FIPS 186-4 C.3.2 Enhanced Miller-Rabin Probabilistic Primality Test
13679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
13779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * Run several iterations of the Miller-Rabin algorithm with randomly-chosen bases. This is an
13879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * alternative to {@link #isMRProbablePrime(BigInteger, SecureRandom, int)} that provides more
13979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * information about a composite candidate, which may be useful when generating or validating
14079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * RSA moduli.
14179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
14279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param candidate
14379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the {@link BigInteger} instance to test for primality.
14479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param random
14579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the source of randomness to use to choose bases.
14679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param iterations
14779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the number of randomly-chosen bases to perform the test for.
14879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @return an {@link MROutput} instance that can be further queried for details.
14979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     */
15079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    public static MROutput enhancedMRProbablePrimeTest(BigInteger candidate, SecureRandom random, int iterations)
15179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
15279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        checkCandidate(candidate, "candidate");
15379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
15479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (random == null)
15579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
15679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("'random' cannot be null");
15779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
15879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (iterations < 1)
15979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
16079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("'iterations' must be > 0");
16179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
16279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
16379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (candidate.bitLength() == 2)
16479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
16579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return MROutput.probablyPrime();
16679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
16779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (!candidate.testBit(0))
16879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
16979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return MROutput.provablyCompositeWithFactor(TWO);
17079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
17179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
17279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger w = candidate;
17379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger wSubOne = candidate.subtract(ONE);
17479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger wSubTwo = candidate.subtract(TWO);
17579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
17679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int a = wSubOne.getLowestSetBit();
17779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger m = wSubOne.shiftRight(a);
17879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
17979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        for (int i = 0; i < iterations; ++i)
18079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
18179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            BigInteger b = BigIntegers.createRandomInRange(TWO, wSubTwo, random);
18279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            BigInteger g = b.gcd(w);
18379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
18479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (g.compareTo(ONE) > 0)
18579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
18679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                return MROutput.provablyCompositeWithFactor(g);
18779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
18879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
18979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            BigInteger z = b.modPow(m, w);
19079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
19179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (z.equals(ONE) || z.equals(wSubOne))
19279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
19379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                continue;
19479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
19579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
19679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            boolean primeToBase = false;
19779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
19879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            BigInteger x = z;
19979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            for (int j = 1; j < a; ++j)
20079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
20179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                z = z.modPow(TWO, w);
20279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
20379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                if (z.equals(wSubOne))
20479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                {
20579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    primeToBase = true;
20679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    break;
20779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                }
20879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
20979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                if (z.equals(ONE))
21079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                {
21179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    break;
21279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                }
21379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
21479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                x = z;
21579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
21679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
21779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (!primeToBase)
21879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
21979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                if (!z.equals(ONE))
22079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                {
22179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    x = z;
22279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    z = z.modPow(TWO, w);
22379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
22479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    if (!z.equals(ONE))
22579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    {
22679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                        x = z;
22779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    }
22879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                }
22979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
23079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                g = x.subtract(ONE).gcd(w);
23179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
23279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                if (g.compareTo(ONE) > 0)
23379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                {
23479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    return MROutput.provablyCompositeWithFactor(g);
23579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                }
23679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
23779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                return MROutput.provablyCompositeNotPrimePower();
23879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
23979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
24079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
24179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        return MROutput.probablyPrime();
24279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
24379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
24479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    /**
24579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * A fast check for small divisors, up to some implementation-specific limit.
24679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
24779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param candidate
24879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the {@link BigInteger} instance to test for division by small factors.
24979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
25079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @return <code>true</code> if the candidate is found to have any small factors,
25179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *         <code>false</code> otherwise.
25279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     */
25379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    public static boolean hasAnySmallFactors(BigInteger candidate)
25479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
25579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        checkCandidate(candidate, "candidate");
25679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
25779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        return implHasAnySmallFactors(candidate);
25879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
25979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
26079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    /**
26179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * FIPS 186-4 C.3.1 Miller-Rabin Probabilistic Primality Test
26279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
26379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * Run several iterations of the Miller-Rabin algorithm with randomly-chosen bases.
26479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
26579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param candidate
26679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the {@link BigInteger} instance to test for primality.
26779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param random
26879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the source of randomness to use to choose bases.
26979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param iterations
27079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the number of randomly-chosen bases to perform the test for.
27179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @return <code>false</code> if any witness to compositeness is found amongst the chosen bases
27279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *         (so <code>candidate</code> is definitely NOT prime), or else <code>true</code>
27379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *         (indicating primality with some probability dependent on the number of iterations
27479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *         that were performed).
27579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     */
27679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    public static boolean isMRProbablePrime(BigInteger candidate, SecureRandom random, int iterations)
27779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
27879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        checkCandidate(candidate, "candidate");
27979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
28079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (random == null)
28179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
28279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("'random' cannot be null");
28379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
28479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (iterations < 1)
28579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
28679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("'iterations' must be > 0");
28779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
28879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
28979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (candidate.bitLength() == 2)
29079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
29179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
29279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
29379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (!candidate.testBit(0))
29479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
29579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return false;
29679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
29779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
29879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger w = candidate;
29979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger wSubOne = candidate.subtract(ONE);
30079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger wSubTwo = candidate.subtract(TWO);
30179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
30279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int a = wSubOne.getLowestSetBit();
30379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger m = wSubOne.shiftRight(a);
30479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
30579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        for (int i = 0; i < iterations; ++i)
30679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
30779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            BigInteger b = BigIntegers.createRandomInRange(TWO, wSubTwo, random);
30879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
30979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (!implMRProbablePrimeToBase(w, wSubOne, m, a, b))
31079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
31179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                return false;
31279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
31379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
31479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
31579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        return true;
31679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
31779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
31879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    /**
31979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * FIPS 186-4 C.3.1 Miller-Rabin Probabilistic Primality Test (to a fixed base).
32079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
32179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * Run a single iteration of the Miller-Rabin algorithm against the specified base.
32279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *
32379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param candidate
32479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the {@link BigInteger} instance to test for primality.
32579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @param base
32679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *            the base value to use for this iteration.
32779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     * @return <code>false</code> if the specified base is a witness to compositeness (so
32879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     *         <code>candidate</code> is definitely NOT prime), or else <code>true</code>.
32979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro     */
33079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    public static boolean isMRProbablePrimeToBase(BigInteger candidate, BigInteger base)
33179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
33279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        checkCandidate(candidate, "candidate");
33379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        checkCandidate(base, "base");
33479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
33579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (base.compareTo(candidate.subtract(ONE)) >= 0)
33679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
33779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("'base' must be < ('candidate' - 1)");
33879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
33979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
34079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (candidate.bitLength() == 2)
34179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
34279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
34379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
34479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
34579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger w = candidate;
34679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger wSubOne = candidate.subtract(ONE);
34779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
34879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int a = wSubOne.getLowestSetBit();
34979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger m = wSubOne.shiftRight(a);
35079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
35179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        return implMRProbablePrimeToBase(w, wSubOne, m, a, base);
35279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
35379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
35479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static void checkCandidate(BigInteger n, String name)
35579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
35679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (n == null || n.signum() < 1 || n.bitLength() < 2)
35779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
35879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("'" + name + "' must be non-null and >= 2");
35979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
36079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
36179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
36279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static boolean implHasAnySmallFactors(BigInteger x)
36379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
36479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        /*
36579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         * Bundle trial divisors into ~32-bit moduli then use fast tests on the ~32-bit remainders.
36679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         */
36779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int m = 2 * 3 * 5 * 7 * 11 * 13 * 17 * 19 * 23;
36879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int r = x.mod(BigInteger.valueOf(m)).intValue();
36979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 2) == 0 || (r % 3) == 0 || (r % 5) == 0 || (r % 7) == 0 || (r % 11) == 0 || (r % 13) == 0
37079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            || (r % 17) == 0 || (r % 19) == 0 || (r % 23) == 0)
37179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
37279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
37379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
37479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
37579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        m = 29 * 31 * 37 * 41 * 43;
37679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        r = x.mod(BigInteger.valueOf(m)).intValue();
37779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 29) == 0 || (r % 31) == 0 || (r % 37) == 0 || (r % 41) == 0 || (r % 43) == 0)
37879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
37979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
38079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
38179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
38279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        m = 47 * 53 * 59 * 61 * 67;
38379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        r = x.mod(BigInteger.valueOf(m)).intValue();
38479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 47) == 0 || (r % 53) == 0 || (r % 59) == 0 || (r % 61) == 0 || (r % 67) == 0)
38579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
38679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
38779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
38879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
38979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        m = 71 * 73 * 79 * 83;
39079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        r = x.mod(BigInteger.valueOf(m)).intValue();
39179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 71) == 0 || (r % 73) == 0 || (r % 79) == 0 || (r % 83) == 0)
39279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
39379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
39479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
39579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
39679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        m = 89 * 97 * 101 * 103;
39779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        r = x.mod(BigInteger.valueOf(m)).intValue();
39879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 89) == 0 || (r % 97) == 0 || (r % 101) == 0 || (r % 103) == 0)
39979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
40079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
40179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
40279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
40379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        m = 107 * 109 * 113 * 127;
40479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        r = x.mod(BigInteger.valueOf(m)).intValue();
40579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 107) == 0 || (r % 109) == 0 || (r % 113) == 0 || (r % 127) == 0)
40679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
40779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
40879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
40979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
41079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        m = 131 * 137 * 139 * 149;
41179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        r = x.mod(BigInteger.valueOf(m)).intValue();
41279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 131) == 0 || (r % 137) == 0 || (r % 139) == 0 || (r % 149) == 0)
41379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
41479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
41579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
41679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
41779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        m = 151 * 157 * 163 * 167;
41879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        r = x.mod(BigInteger.valueOf(m)).intValue();
41979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 151) == 0 || (r % 157) == 0 || (r % 163) == 0 || (r % 167) == 0)
42079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
42179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
42279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
42379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
42479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        m = 173 * 179 * 181 * 191;
42579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        r = x.mod(BigInteger.valueOf(m)).intValue();
42679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 173) == 0 || (r % 179) == 0 || (r % 181) == 0 || (r % 191) == 0)
42779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
42879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
42979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
43079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
43179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        m = 193 * 197 * 199 * 211;
43279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        r = x.mod(BigInteger.valueOf(m)).intValue();
43379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((r % 193) == 0 || (r % 197) == 0 || (r % 199) == 0 || (r % 211) == 0)
43479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
43579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
43679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
43779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
43879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        /*
43979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         * NOTE: Unit tests depend on SMALL_FACTOR_LIMIT matching the
44079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         * highest small factor tested here.
44179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         */
44279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        return false;
44379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
44479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
44579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static boolean implMRProbablePrimeToBase(BigInteger w, BigInteger wSubOne, BigInteger m, int a, BigInteger b)
44679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
44779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger z = b.modPow(m, w);
44879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
44979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (z.equals(ONE) || z.equals(wSubOne))
45079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
45179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return true;
45279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
45379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
45479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        boolean result = false;
45579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
45679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        for (int j = 1; j < a; ++j)
45779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
45879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            z = z.modPow(TWO, w);
45979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
46079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (z.equals(wSubOne))
46179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
46279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                result = true;
46379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                break;
46479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
46579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
46679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (z.equals(ONE))
46779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
46879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                return false;
46979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
47079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
47179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
47279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        return result;
47379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
47479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
47579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static STOutput implSTRandomPrime(Digest d, int length, byte[] primeSeed)
47679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
47779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int dLen = d.getDigestSize();
47879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
47979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (length < 33)
48079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
48179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            int primeGenCounter = 0;
48279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
48379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            byte[] c0 = new byte[dLen];
48479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            byte[] c1 = new byte[dLen];
48579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
48679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            for (;;)
48779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
48879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                hash(d, primeSeed, c0, 0);
48979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                inc(primeSeed, 1);
49079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
49179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                hash(d, primeSeed, c1, 0);
49279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                inc(primeSeed, 1);
49379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
49479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                int c = extract32(c0) ^ extract32(c1);
49579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                c &= (-1 >>> (32 - length));
49679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                c |= (1 << (length - 1)) | 1;
49779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
49879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                ++primeGenCounter;
49979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
50079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                long c64 = c & 0xFFFFFFFFL;
50179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                if (isPrime32(c64))
50279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                {
50379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    return new STOutput(BigInteger.valueOf(c64), primeSeed, primeGenCounter);
50479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                }
50579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
50679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                if (primeGenCounter > (4 * length))
50779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                {
50879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    throw new IllegalStateException("Too many iterations in Shawe-Taylor Random_Prime Routine");
50979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                }
51079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
51179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
51279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
51379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        STOutput rec = implSTRandomPrime(d, (length + 3) / 2, primeSeed);
51479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
51579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger c0 = rec.getPrime();
51679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        primeSeed = rec.getPrimeSeed();
51779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int primeGenCounter = rec.getPrimeGenCounter();
51879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
51979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int outlen = 8 * dLen;
52079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int iterations = (length - 1) / outlen;
52179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
52279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int oldCounter = primeGenCounter;
52379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
52479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger x = hashGen(d, primeSeed, iterations + 1);
52579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        x = x.mod(ONE.shiftLeft(length - 1)).setBit(length - 1);
52679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
52779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger c0x2 = c0.shiftLeft(1);
52879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger tx2 = x.subtract(ONE).divide(c0x2).add(ONE).shiftLeft(1);
52979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int dt = 0;
53079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
53179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        BigInteger c = tx2.multiply(c0).add(ONE);
53279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
53379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        /*
53479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         * TODO Since the candidate primes are generated by constant steps ('c0x2'), sieving could
53579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         * be used here in place of the 'hasAnySmallFactors' approach.
53679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         */
53779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        for (;;)
53879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
53979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (c.bitLength() > length)
54079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
54179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                tx2 = ONE.shiftLeft(length - 1).subtract(ONE).divide(c0x2).add(ONE).shiftLeft(1);
54279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                c = tx2.multiply(c0).add(ONE);
54379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
54479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
54579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            ++primeGenCounter;
54679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
54779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            /*
54879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro             * This is an optimization of the original algorithm, using trial division to screen out
54979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro             * many non-primes quickly.
55079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro             *
55179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro             * NOTE: 'primeSeed' is still incremented as if we performed the full check!
55279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro             */
55379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (!implHasAnySmallFactors(c))
55479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
55579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                BigInteger a = hashGen(d, primeSeed, iterations + 1);
55679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                a = a.mod(c.subtract(THREE)).add(TWO);
55779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
55879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                tx2 = tx2.add(BigInteger.valueOf(dt));
55979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                dt = 0;
56079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
56179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                BigInteger z = a.modPow(tx2, c);
56279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
56379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                if (c.gcd(z.subtract(ONE)).equals(ONE) && z.modPow(c0, c).equals(ONE))
56479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                {
56579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    return new STOutput(c, primeSeed, primeGenCounter);
56679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                }
56779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
56879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            else
56979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
57079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                inc(primeSeed, iterations + 1);
57179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
57279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
57379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (primeGenCounter >= ((4 * length) + oldCounter))
57479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
57579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                throw new IllegalStateException("Too many iterations in Shawe-Taylor Random_Prime Routine");
57679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
57779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
57879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            dt += 2;
57979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            c = c.add(c0x2);
58079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
58179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
58279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
58379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static int extract32(byte[] bs)
58479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
58579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int result = 0;
58679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
58779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int count = Math.min(4, bs.length);
58879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        for (int i = 0; i < count; ++i)
58979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
59079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            int b = bs[bs.length - (i + 1)] & 0xFF;
59179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            result |= (b << (8 * i));
59279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
59379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
59479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        return result;
59579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
59679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
59779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static void hash(Digest d, byte[] input, byte[] output, int outPos)
59879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
59979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        d.update(input, 0, input.length);
60079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        d.doFinal(output, outPos);
60179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
60279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
60379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static BigInteger hashGen(Digest d, byte[] seed, int count)
60479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
60579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int dLen = d.getDigestSize();
60679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int pos = count * dLen;
60779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        byte[] buf = new byte[pos];
60879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        for (int i = 0; i < count; ++i)
60979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
61079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            pos -= dLen;
61179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            hash(d, seed, buf, pos);
61279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            inc(seed, 1);
61379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
61479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        return new BigInteger(1, buf);
61579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
61679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
61779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static void inc(byte[] seed, int c)
61879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
61979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        int pos = seed.length;
62079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        while (c > 0 && --pos >= 0)
62179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
62279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            c += (seed[pos] & 0xFF);
62379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            seed[pos] = (byte)c;
62479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            c >>>= 8;
62579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
62679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
62779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
62879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    private static boolean isPrime32(long x)
62979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    {
63079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (x >>> 32 != 0L)
63179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
63279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            throw new IllegalArgumentException("Size limit exceeded");
63379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
63479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
63579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        /*
63679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         * Use wheel factorization with 2, 3, 5 to select trial divisors.
63779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro         */
63879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
63979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if (x <= 5L)
64079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
64179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return x == 2L || x == 3L || x == 5L;
64279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
64379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
64479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        if ((x & 1L) == 0L || (x % 3L) == 0L || (x % 5L) == 0L)
64579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
64679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            return false;
64779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
64879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
64979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        long[] ds = new long[]{ 1L, 7L, 11L, 13L, 17L, 19L, 23L, 29L };
65079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        long base = 0L;
65179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        for (int pos = 1;; pos = 0)
65279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        {
65379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            /*
65479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro             * Trial division by wheel-selected divisors
65579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro             */
65679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            while (pos < ds.length)
65779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
65879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                long d = base + ds[pos];
65979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                if (x % d == 0L)
66079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                {
66179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                    return x < 30L;
66279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                }
66379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                ++pos;
66479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
66579d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
66679d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            base += 30L;
66779d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro
66879d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            if (base * base >= x)
66979d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            {
67079d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro                return true;
67179d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro            }
67279d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro        }
67379d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro    }
67479d3bf2425a53baab7feb744dad710b6c15533c9Sergio Giro}