Loading [Contrib]/a11y/accessibility-menu.js
$\newcommand{\O}{\mathrm{O}}$ My Algorithm : kopricky アルゴリズムライブラリ

kopricky アルゴリズムライブラリ

Miller Rabbin

コードについての説明(個人的メモ)

ミラーラビン法は素数判定を高速に行うモンテカルロアルゴリズムである. つまり間違える可能性はあるが高い可能性で正解を返すアルゴリズムである. $2^{32}$ および $2^{64}$以下などの数については経験的に $3$ 個および $7$ 個のある要素について調べれば十分であることが分かっているためそういう意味では決定的アルゴリズムと考えてもよい.
素数ならばフェルマーの小定理を満たすという必要条件をチェックしていくフェルマーテストでは素数でないのに必ず素数と判定してしまうようなもの(カーマイケル数)が存在するなどの問題があったためこちらの方がより正確な素数判定アルゴリズムと言えるかもしれない.
実装上の注意として long long 型の数($2^{32}$ 以上の数) について素数判定を行う場合は __int128 型の剰余算を途中で行う必要があるが, これがびっくりするくらい遅い(たぶん long long 型の剰余算の $10$ 倍近く遅い).
ただ計算量的には 1 つの数の判定に $\O (|numset| \ast \log (VAL))$ ($VAL$ は判定する数) time とかなのでたくさんの数の素数判定を行わないなら特に気にする必要のない注意である.
numset の値については こちらのサイト が参考になる.

(関数)
check$(n)$ : $n$ が素数であるかどうか(true / false)を返す

時間計算量: $\O (|numset| \ast \log (VAL))$ ($VAL$ は判定する数)

コード

  1. const unsigned int numset[] = {2,7,61}; // < 4,759,123,141 ≒ 4×10^9 までは決定的
  2. // const unsigned long long numset[] = {2, 325, 9375, 28178, 450775, 9780504, 1795265022ULL}; // 少なくとも 2^64 までは決定的
  3.  
  4. // int 型の素数判定の場合
  5.  
  6. unsigned int mod_pow(unsigned int x, unsigned int k, unsigned int mod){
  7.     unsigned int res = 1;
  8.     while(k){
  9.         if(k & 1) res = (unsigned long long)res * x % mod;
  10.         x = (unsigned long long)x * x % mod;
  11.         k >>= 1;
  12.     }
  13.     return res;
  14. }
  15.  
  16. bool check(unsigned int n){
  17.     if(n < 2 || ((n % 6 != 1) && (n % 6 != 5))) return (n == 2) || (n == 3);
  18.     unsigned int d = n - 1, s = 0;
  19.     while(d % 2 == 0){
  20.         d /= 2;
  21.         s++;
  22.     }
  23.     for(unsigned int a : numset){
  24.         if(a == n) return true;
  25.         unsigned int res = mod_pow(a, d, n);
  26.         if(res == 1) continue;
  27.         bool ok = true;
  28.         for(unsigned int r = 0; r < s; r++){
  29.             if(res == n-1){
  30.                 ok = false;
  31.                 break;
  32.             }
  33.             res = (unsigned long long)res * res % n;
  34.         }
  35.         if(ok) return false;
  36.     }
  37.     return true;
  38. }
  39.  
  40. // long long型の素数判定の場合
  41.  
  42. // unsigned long long mod_mul(unsigned long long a, unsigned long long b, unsigned long long mod) {
  43. // 	long long ret = a * b - mod * (unsigned long long)((long double)(a) * (long double)(b) / (long double)(mod));
  44. // 	return ret + mod * (ret < 0) - mod * (ret >= (ll)mod);
  45. // }
  46.  
  47. // unsigned long long mod_pow(unsigned long long x, unsigned long long k, unsigned long long mod){
  48. //     unsigned long long res = 1;
  49. //     while(k){
  50. //         if(k & 1) res = mod_mul(res, x, mod);
  51. //         x = mod_mul(x, x, mod);
  52. //         k >>= 1;
  53. //     }
  54. //     return res;
  55. // }
  56.  
  57. // bool check(unsigned long long n){
  58. //     if(n < 2 || ((n % 6 != 1) && (n % 6 != 5))) return (n == 2) || (n == 3);
  59. //     unsigned long long d = n - 1, s = 0;
  60. //     while(d % 2 == 0){
  61. //         d /= 2;
  62. //         s++;
  63. //     }
  64. //     for(unsigned long long a : numset){
  65. //         if(a % n == 0) return true;
  66. //         unsigned long long res = mod_pow(a, d, n);
  67. //         if(res == 1) continue;
  68. //         bool ok = true;
  69. //         for(unsigned int r = 0; r < s; r++){
  70. //             if(res == n-1){
  71. //                 ok = false;
  72. //                 break;
  73. //             }
  74. //             res = mod_mul(res, res, n);
  75. //         }
  76. //         if(ok) return false;
  77. //     }
  78. //     return true;
  79. // }

verify 用の問題

yukicoder : ミラー・ラビン素数判定法のテスト 提出コード