素数の判定 ================================================================================ つぎのアルゴリズムで、素数を判定してみました。また、実行時間も測定しています。 1. 1~1億までの範囲の素数表引き 2. 試し割りアルゴリズム 3. フェルマーの因数分解アルゴリズム 4. 正規表現を利用 参考までに、 5. BigIntegerのisProbablePrime()メソッドによる判定 それぞれのかかった時間(ミリ秒)を記します。時間は、それぞれ、3~5回実行したものの 平均です。判定する範囲と個数は違いますが、参考まで。 なお、時間の確認はJUnitの結果です。 --------------- --------- ----------------------- 1~999 99990001~99999999の奇数 --------------- --------- ----------------------- 1. 素数表 0 0.008 2. 試し割り 0 0.360 3. フェルマー 0.008 538,223.238 4. 正規表現 0.312 - 5. BigInteger 0.357 0.143 --------------- --------- ----------------------- フェルマーの因数分解アルゴリズムは、かなり遅いです。これは、計算のたびに平方根を 計算しているからです。因数がsqrt(n)に近い場合は効率はよいですが、そうでない場合 は実用に耐えそうにありません。 正規表現は考えは面白いですが、1並び数R(n)を作るのが大変です。実際に、R(99990001) は、「java.lang.OutOfMemoryError: Java heap space」で動きません。 上記の実行時間の大半はR(n)を作っている時間でしょう。 BigIntegerのisProbablePrime()メソッドでは、99990001~99999999の奇数の方が速いので 驚きです。 ■1. 1~1億までの範囲の素数表引き 1億までの素数は、 兵庫県立明石北高等学校のホームページ http://www.hyogo-c.ed.jp/~meihoku-hs/club/astronomy-p.html からダウンロードしたファイルをを使いました。 また、表引きは2分検索として、ArraysクラスのbinarySearch()メソッドを使いました。 ■2. 試し割りアルゴリズム 3以上の奇数の場合に、3,5,7, ..., sqrt(n) までで割って、割り切れれば合成数、いず れの数値でも割り切れなければ素数になります。 ■3. フェルマーの因数分解アルゴリズム 3以上の奇数nとします。 n = x2 - y2 = (x - y)(x + y)、x >= y >= 0 となるような非負整数 x、y を見つける。見つかれば、(x - y)、(x + y) は n の因数と なる。 n が平方数の場合は、n = x2 として、y = 0 となる。因数は、x。 n が平方数でない、すなわち、y > 0 のとき、 x = sqrt(n + y2) > sqrt(n) なので、x = [sqrt(n)] から探索を始めればよい。 また、探索の終わりは、nの因数を見つけることなので、nの半分まででよい。 すなわち、n は奇数であるので、(n + 1) / 2 になるまで因数が見つからなければ探索終 了してよい。 ■4. 正規表現を利用 ネットで調べていたら、「正規表現で素数判定」 http://d.hatena.ne.jp/ytakano/20100721/1279736214 というサイトを見つけました。 指定された正整数n(int型の最大値以下)が素数がどうかを判断するのに、nの1並び数 R(n) = (10n - 1) / 9 の文字列表現が、次の正規表現とマッチすれば合成数である。 ^1?$|^(11+?)\1+$ なかなか、面白そうなので、試してみました。ただ、このままでもよいのですが、試し割 りアルゴリズムでも、フェルマーの因数分解アルゴリズムでも、3以上の奇数に限定して いるので、正規表現でも3以上の奇数限定にしました。正規表現はつぎのようになります。 ^(1(11)+?)\1+$ □メタ記号 ----- ---------------------------------------------------- 記号 意味 ----- ---------------------------------------------------- ^ 行の先頭 $ 行の末尾 X+? X、1 回以上(最短一致数量子) (X) X、前方参照を行う正規表現グループ \n マッチした n 番目の前方参照を行う正規表現グループ X+ X、1 回以上(最長一致数量子) ----- ---------------------------------------------------- □動作のロジック 1(11)+? なので、"111" と最短一致でマッチするかを確認する。これは3以上の奇数であ れば必ずマッチする。 (1(11)+?)\1 は、"111" と最短一致でマッチしたときの1番目のグループなので、これは 文字列 "111" を意味する。 3の場合は、"111" なので、最短一致でマッチした残り "" と、 (111)+ がマッチするか を確認し、 もうマッチしないので、素数となる。 5の場合は、"11111" なので、最短一致でマッチした残り "11" と、 (111)+ がマッチす るかを確認し、マッチしない。 つぎに、1(11)+? なので、"11111" と最短一致でマッチするかを確認する。 マッチするので、残り "" が (11111)+ とマッチするか確認し、これもマッチしないので 素数となる。 9の場合は、"111111111" なので、最短一致でマッチした残り "111111" と、 (111)+ がマッチするかを確認し、これはマッチする。したがって、合成数となる。 以下、3の倍数は、(111)+ でマッチするので、合成数になる。 25の場合は、最短一致でマッチした残り "11...1"(22個) と、(111)+ がマッチするかを 確認し、これはマッチしない。 つぎに、1(11)+? なので、"11111" と最短一致でマッチするかを確認する。 マッチするので、残り "11...1"(20個) が (11111)+ とマッチするか確認し、これはマッ チする。 したがって、合成数となる。 以下、5の倍数は、(11111)+ でマッチするので、合成数になる。 いくつかの数値についてつぎにまとめる。 --- --------- -------------------------- ------ ------ 数 最短一致 グループのマッチ 結果 判定 --- --------- -------------------------- ------ ------ 3 111 "" と (111)+ false 素数 --- --------- -------------------------- ------ ------ 5 111 "11" と (111)+ false 11111 "" と (11111)+ false 素数 --- --------- -------------------------- ------ ------ 7 111 "1111" と (111)+ false 11111 "11" と (11111)+ false 1111111 "" と (1111111)+ false 素数 --- --------- -------------------------- ------ ------ 9 111 "111111" と (111)+ true 合成数 --- --------- -------------------------- ------ ------ 25 111 "111"(22個) と (111)+ false 11111 "111"(20個) と (11111)+ true 合成数 --- --------- -------------------------- ------ ------ 77 111 "111"(74個) と (111)+ false 11111 "111"(72個) と (11111)+ false 1111111 "111"(70個) と (1111111)+ true 合成数 --- --------- -------------------------- ------ ------ ■5. BigIntegerのisProbablePrime()メソッドによる判定 ソースコードをみると、つぎの2つの方法を使っている。 Miller-Rabin 法 Lucas-Lehmer 法 isProbablePrime()メソッドは、名前のとおり、素数である可能性が高い場合にtrue、必 ず合成数である場合にfalseが返るので、今回のテストでも、素数なのに合成数と判定さ れることがあった。 このメソッドは、呼び出し側が許容しない 確率の尺度を引数で渡せる。今回は、4を渡し たので、trueで、素数である確率は 1 - 1/24 = 1 - 1/16 = 0.9375 である。 ■プログラム --- import java.io.File; import java.io.FileNotFoundException; import java.math.BigInteger; import java.util.Arrays; import java.util.Scanner; import java.util.regex.Pattern; public class Factorization02 { /** * べき乗演算子の文字列表現 */ public static final String POWER_OPERATOR = "**"; /** * 正規表現を使って素数を判断するパターン(3以上の奇数に対応) */ public static Pattern PRIME_NUMBER_REGEX_ODD = Pattern.compile("^(1(11)+?)\\1+$"); /** * 100000000以下の素数の配列 */ private static int[] PRIME_NUMBERS_TABLE = new int[5761455]; static { read("prime.txt"); } /** * 単純なアルゴリズムで、指定された正整数が素数がどうかを判断する。 * @param n 正整数 * @return 素数の場合true */ public boolean isPrimeNumberBySimple(long n) { if (n < 1) { throw new IllegalArgumentException("n = " + n); } if (n == 1) { return false; } else if (n == 2) { return true; } else if ((n % 2) == 0) { return false; } for (int i = 3; i <= (int) Math.sqrt(n); i += 2) { if (n % i == 0) { return false; } } return true; } /** * フェルマーの素因数分解アルゴリズムで、指定された正整数が素数がどうかを判断する。 * @param n 正整数 * @return 素数の場合true */ public boolean isPrimeNumberByFermat(long n) { if (n < 1) { throw new IllegalArgumentException("n = " + n); } if (n == 1) { return false; } else if (n == 2) { return true; } else if ((n % 2) == 0) { return false; } long x = (long) Math.sqrt(n); if (n == x * x) { return false; // n は合成数 } for (x++; x < (n + 1) / 2; x++) { double dy = Math.sqrt(x * x - n); long y = (long) dy; if ((y * y) == (x * x - n) ) { // y が整数の場合は、合成数 return false; // n = (x - y)(x + y) } } return true; // n は素数 } /** * 正規表現で、指定された正整数n(int型の最大値以下)が素数がどうかを判断する。 * 正整数が3以上の奇数であれば、nの1並び数(R(n) = (10n - 1) / 9)の文字列表現が * 次の正規表現とマッチすれば合成数である。 * ^(1(11)+?)\1+$ * @param n int型の最大値以下の正整数 * @return 素数の場合true */ public boolean isPrimeByRegex(long n) { if (n < 1 || n > Integer.MAX_VALUE) { throw new IllegalArgumentException("n = " + n); } if (n == 1) { return false; } else if (n == 2) { return true; } else if ((n % 2) == 0) { return false; } // 1並び数を作成する StringBuilder nn = new StringBuilder((int) n); nn.append("1"); for (long i = 0; i < n / 2; i++) { nn.append("11"); } return !PRIME_NUMBER_REGEX_ODD.matcher(nn.toString()).matches(); } /** * 素数表を使って、指定された正整数(100000000以下)が素数がどうかを判断する。 * @param n 100000000以下の正整数 * @return 素数の場合true */ public boolean isPrimeByTable(long n) { if (n < 1 || n > 100000000) { throw new IllegalArgumentException("n = " + n); } return Arrays.binarySearch(PRIME_NUMBERS_TABLE, (int) n) >= 0; } /** * BigIntegerで、指定された正整数が素数がどうかを判断する。 * @param n 正整数 * @return 素数の可能性がある場合true、必ず合成数である場合false */ public boolean isProbablePrimeByBigInteger(long n, int certainty) { if (n < 1) { throw new IllegalArgumentException("n = " + n); } return new BigInteger(Long.toString(n)).isProbablePrime(certainty); // n は素数 } private static int read(String fileName) { int i = 0; try { Scanner in = new Scanner(new File(fileName)); for (i = 0; i < PRIME_NUMBERS_TABLE.length && in.hasNextInt(); i++) { PRIME_NUMBERS_TABLE[i] = in.nextInt(); } } catch (FileNotFoundException e) { System.out.println(e); } return i; } } --- ■テスト用のプログラム --- import java.util.Arrays; import junit.framework.TestCase; public class Factorization02Test extends TestCase { private static final int[] PRIME_NUMBERS_1_999 = { 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571, 577, 587, 593, 599, 601, 607, 613, 617, 619, 631, 641, 643, 647, 653, 659, 661, 673, 677, 683, 691, 701, 709, 719, 727, 733, 739, 743, 751, 757, 761, 769, 773, 787, 797, 809, 811, 821, 823, 827, 829, 839, 853, 857, 859, 863, 877, 881, 883, 887, 907, 911, 919, 929, 937, 941, 947, 953, 967, 971, 977, 983, 991, 997, }; private static final int[] PRIME_NUMBERS_99990001_99999999 = { 99990001, 99990029, 99990043, 99990049, 99990053, 99990083, 99990091, 99990113, 99990139, 99990161, 99990181, 99990239, 99990263, 99990281, 99990287, 99990299, 99990307, 99990329, 99990337, 99990389, 99990419, 99990461, 99990467, 99990487, 99990491, 99990493, 99990511, 99990533, 99990547, 99990551, 99990571, 99990577, 99990581, 99990587, 99990589, 99990641, 99990647, 99990679, 99990703, 99990731, 99990739, 99990757, 99990799, 99990871, 99990887, 99990893, 99990911, 99990929, 99990967, 99990971, 99990973, 99991009, 99991033, 99991049, 99991063, 99991081, 99991103, 99991109, 99991207, 99991237, 99991249, 99991271, 99991273, 99991277, 99991303, 99991343, 99991363, 99991379, 99991387, 99991447, 99991481, 99991499, 99991517, 99991519, 99991559, 99991561, 99991583, 99991609, 99991667, 99991673, 99991687, 99991721, 99991733, 99991753, 99991757, 99991769, 99991789, 99991819, 99991831, 99991877, 99991891, 99991897, 99991901, 99991933, 99991937, 99991951, 99991981, 99991987, 99991999, 99992027, 99992051, 99992059, 99992089, 99992099, 99992119, 99992141, 99992149, 99992161, 99992173, 99992209, 99992213, 99992219, 99992239, 99992251, 99992303, 99992309, 99992311, 99992363, 99992381, 99992401, 99992413, 99992429, 99992437, 99992447, 99992449, 99992461, 99992467, 99992471, 99992483, 99992509, 99992527, 99992533, 99992539, 99992617, 99992671, 99992681, 99992747, 99992771, 99992773, 99992777, 99992791, 99992797, 99992803, 99992807, 99992813, 99992831, 99992843, 99992869, 99992897, 99992903, 99992911, 99992917, 99992939, 99992941, 99992947, 99992993, 99993001, 99993053, 99993071, 99993097, 99993109, 99993119, 99993161, 99993163, 99993169, 99993193, 99993199, 99993211, 99993253, 99993263, 99993281, 99993287, 99993317, 99993329, 99993331, 99993343, 99993349, 99993359, 99993407, 99993419, 99993427, 99993433, 99993449, 99993461, 99993463, 99993469, 99993497, 99993511, 99993527, 99993583, 99993601, 99993611, 99993671, 99993679, 99993743, 99993767, 99993781, 99993787, 99993799, 99993821, 99993833, 99993841, 99993889, 99993893, 99993911, 99993917, 99993923, 99993947, 99993991, 99993997, 99994009, 99994021, 99994043, 99994061, 99994073, 99994079, 99994123, 99994129, 99994157, 99994177, 99994183, 99994201, 99994211, 99994229, 99994231, 99994253, 99994267, 99994273, 99994277, 99994309, 99994319, 99994327, 99994361, 99994373, 99994379, 99994471, 99994481, 99994553, 99994567, 99994571, 99994579, 99994607, 99994621, 99994633, 99994637, 99994669, 99994673, 99994679, 99994691, 99994703, 99994711, 99994723, 99994751, 99994759, 99994771, 99994801, 99994819, 99994849, 99994861, 99994877, 99994883, 99994907, 99994919, 99994963, 99995011, 99995039, 99995041, 99995059, 99995087, 99995111, 99995113, 99995171, 99995177, 99995197, 99995209, 99995213, 99995227, 99995249, 99995257, 99995261, 99995267, 99995279, 99995281, 99995297, 99995303, 99995317, 99995407, 99995417, 99995453, 99995507, 99995527, 99995549, 99995551, 99995563, 99995569, 99995573, 99995591, 99995603, 99995629, 99995641, 99995647, 99995663, 99995669, 99995671, 99995681, 99995699, 99995711, 99995723, 99995747, 99995767, 99995771, 99995783, 99995807, 99995813, 99995839, 99995851, 99995869, 99995891, 99995893, 99995899, 99995923, 99995933, 99995939, 99995957, 99995983, 99995999, 99996011, 99996059, 99996077, 99996079, 99996097, 99996121, 99996131, 99996133, 99996139, 99996151, 99996179, 99996181, 99996191, 99996229, 99996241, 99996251, 99996257, 99996271, 99996287, 99996293, 99996317, 99996343, 99996353, 99996389, 99996409, 99996419, 99996473, 99996509, 99996511, 99996517, 99996521, 99996551, 99996577, 99996583, 99996619, 99996647, 99996653, 99996667, 99996671, 99996679, 99996691, 99996697, 99996719, 99996749, 99996769, 99996779, 99996797, 99996803, 99996821, 99996833, 99996863, 99996899, 99996931, 99996947, 99996971, 99996983, 99997033, 99997063, 99997067, 99997069, 99997081, 99997127, 99997129, 99997151, 99997159, 99997187, 99997189, 99997243, 99997267, 99997357, 99997367, 99997369, 99997393, 99997411, 99997427, 99997441, 99997459, 99997493, 99997511, 99997523, 99997531, 99997543, 99997571, 99997609, 99997619, 99997649, 99997669, 99997679, 99997687, 99997691, 99997759, 99997769, 99997783, 99997789, 99997837, 99997841, 99997867, 99997871, 99997873, 99997879, 99997897, 99997901, 99997907, 99997913, 99997949, 99997957, 99997973, 99997981, 99998021, 99998051, 99998071, 99998099, 99998111, 99998147, 99998191, 99998207, 99998231, 99998251, 99998257, 99998291, 99998299, 99998323, 99998357, 99998363, 99998377, 99998413, 99998417, 99998429, 99998441, 99998447, 99998449, 99998473, 99998491, 99998497, 99998513, 99998543, 99998551, 99998557, 99998603, 99998609, 99998611, 99998663, 99998749, 99998783, 99998797, 99998819, 99998851, 99998869, 99998891, 99998909, 99998929, 99998933, 99998953, 99998971, 99998999, 99999043, 99999073, 99999077, 99999079, 99999089, 99999103, 99999113, 99999131, 99999157, 99999167, 99999187, 99999217, 99999247, 99999257, 99999259, 99999307, 99999323, 99999329, 99999343, 99999353, 99999373, 99999401, 99999437, 99999439, 99999481, 99999509, 99999517, 99999539, 99999541, 99999547, 99999551, 99999563, 99999587, 99999589, 99999611, 99999617, 99999623, 99999643, 99999677, 99999703, 99999721, 99999773, 99999787, 99999821, 99999827, 99999839, 99999847, 99999931, 99999941, 99999959, 99999971, 99999989, }; public void testSetUp() throws Exception { new Factorization02(); } public void testIsPrimeNumberByTable1() { Factorization02 f = new Factorization02(); for (int i = 1; i < 1000; i++) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_1_999, i) >= 0, f.isPrimeByTable(i)); } } public void testIsPrimeNumberByTable2() { Factorization02 f = new Factorization02(); for (int i = 99990001; i < 100000000; i += 2) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_99990001_99999999, i) >= 0, f.isPrimeByTable(i)); } } public void testIsPrimeNumberSimple1() { Factorization02 f = new Factorization02(); for (int i = 1; i < 1000; i++) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_1_999, i) >= 0, f.isPrimeNumberBySimple(i)); } } public void testIsPrimeNumberSimple2() { Factorization02 f = new Factorization02(); for (int i = 99990001; i < 100000000; i += 2) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_99990001_99999999, i) >= 0, f.isPrimeNumberBySimple(i)); } } public void testIsPrimeByRegex1() { Factorization02 f = new Factorization02(); for (int i = 1; i < 1000; i++) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_1_999, i) >= 0, f.isPrimeByRegex(i)); } } public void testIsPrimeByRegex2() { Factorization02 f = new Factorization02(); for (int i = 99990001; i < 100000000; i += 2) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_99990001_99999999, i) >= 0, f.isPrimeByRegex(i)); System.out.print(i + " "); System.out.flush(); } /* * java.lang.OutOfMemoryError: Java heap space */ } public void testIsPrimeNumberFermat1() { Factorization02 f = new Factorization02(); for (int i = 1; i < 1000; i++) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_1_999, i) >= 0, f.isPrimeNumberByFermat(i)); } } public void testIsPrimeNumberFermat2() { Factorization02 f = new Factorization02(); for (int i = 99990001; i < 100000000; i += 2) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_99990001_99999999, i) >= 0, f.isPrimeNumberByFermat(i)); } } public void testIsProbablePrimeByBigInteger1() { Factorization02 f = new Factorization02(); for (int i = 1; i < 1000; i++) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_1_999, i) >= 0, f.isProbablePrimeByBigInteger(i, 4)); } } public void testIsProbablePrimeByBigInteger2() { Factorization02 f = new Factorization02(); for (int i = 99990001; i < 100000000; i += 2) { assertEquals(Arrays.binarySearch(PRIME_NUMBERS_99990001_99999999, i) >= 0, f.isProbablePrimeByBigInteger(i, 4)); } } } ---