首页 > 代码库 > Miller-Rabin素数测试
Miller-Rabin素数测试
费尔马小定理:如果p是一个素数,且0<a<p,则a^(p-1)%p=1.利用费尔马小定理,对于给定的整数n,可以设计素数判定算法,
通过计算d=a^(n-1)%n来判断n的素性,当d!=1时,n肯定不是素数,当d=1时,n 很可能是素数.
二次探测定理:如果
利用二次探测定理,可以再利用费尔马小定理计算a^(n-1)%n的过程 中增加对整数n的二次探测,
一旦发现违背二次探测条件,即得出n不是素数的结论.
先贴出long long 范围的代码：
1 /* 2 long long 范围的Miller-Rabin素数测试 3 */ 4 5 #include <iostream> 6 #include <stdio.h> 7 #include <string.h> 8 #include <stdlib.h> 9 #include <time.h> 10 #define ll long long 11 using namespace std; 12 13 //计算a*b%mod 14 ll mod_mul(ll a, ll b, ll n){ 15 ll cnt = 0LL; 16 while(b){ 17 if(b&1LL) cnt = (cnt+a)%n; 18 a=(a+a)%n; 19 b >>= 1LL; 20 } 21 return cnt; 22 } 23 //计算a^b%mod 24 ll mod_exp(ll a, ll b, ll n){ 25 ll res = 1LL; 26 while(b){ 27 if(b&1LL) res = mod_mul(res,a,n); 28 a = mod_mul(a,a,n); 29 b >>= 1LL; 30 } 31 return res; 32 } 33 //Miller-Rabin测试，测试n是否为素数 34 bool Miller_Rabin(ll n, int respat){ 35 if(n==2LL || n==3LL || n==5LL || n==7LL || n==11LL) return true; 36 if(n==1 || !(n%2) || !(n%3) || !(n%5) || !(n%7) || !(n%11)) return false; 37 38 int k = 0; 39 ll d = n-1; //要求x^u%n 不为1一定不是素数 40 41 while(!(d&1LL)){ 42 k++; d >>= 1LL; 43 } 44 srand((ll)time(0)); 45 for(int i = 0; i < respat; i ++) { 46 ll a = rand()%(n-2)+2; //在[2,n)中取整数 47 ll x = mod_exp(a,d,n); 48 ll y = 0LL; 49 for(int j = 0; j < k; j ++){ 50 y = mod_mul(x,x,n); 51 if(1LL==y && 1LL!=x && n-1LL!=x)return false; //二次探测定理，这里如果y = 1则x必须等于 1 52 //或n-1,否则可以判断不是素数 53 x = y; 54 } 55 if(1LL != y) return false; //费马小定理 56 } 57 return true; 58 } 59 int main(){ 60 int n; 61 cin>>n; 62 if(Miller_Rabin(n,5))cout << "Yes\n"; 63 else cout << "No\n"; 64 return 0; 65 }
当数字达到10^30时，就要用到大数了，贴下模板。
1 #include <cstdio> 2 #include <cstring> 3 #include <cstdlib> 4 5 #define MAXL 4 6 #define M10 1000000000 7 #define Z10 9 8 9 const int zero[MAXL - 1] = {0}; 10 11 struct bnum{ 12 int data[MAXL]; // 断成每截9个长度 13 14 // 读取字符串并转存 15 void read(){ 16 memset(data, 0, sizeof(data)); 17 char buf[32]; 18 scanf("%s", buf); 19 int len = (int)strlen(buf); 20 int i = 0, k; 21 while (len >= Z10){ 22 for (k = len - Z10; k < len; ++k){ 23 data[i] = data[i] * 10 + buf[k] - ‘0‘; 24 } 25 ++i; 26 len -= Z10; 27 } 28 if (len > 0){ 29 for (k = 0; k < len; ++k){ 30 data[i] = data[i] * 10 + buf[k] - ‘0‘; 31 } 32 } 33 } 34 35 bool operator == (const bnum &x) { 36 return memcmp(data, x.data, sizeof(data)) == 0; 37 } 38 39 bnum & operator = (const int x){ 40 memset(data, 0, sizeof(data)); 41 data[0] = x; 42 return *this; 43 } 44 45 bnum operator + (const bnum &x){ 46 int i, carry = 0; 47 bnum ans; 48 for (i = 0; i < MAXL; ++i){ 49 ans.data[i] = data[i] + x.data[i] + carry; 50 carry = ans.data[i] / M10; 51 ans.data[i] %= M10; 52 } 53 return ans; 54 } 55 56 bnum operator - (const bnum &x){ 57 int i, carry = 0; 58 bnum ans; 59 for (i = 0; i < MAXL; ++i){ 60 ans.data[i] = data[i] - x.data[i] - carry; 61 if (ans.data[i] < 0){ 62 ans.data[i] += M10; 63 carry = 1; 64 } 65 else{ 66 carry = 0; 67 } 68 } 69 return ans; 70 } 71 72 // assume *this < x * 2 73 bnum operator % (const bnum &x){ 74 int i; 75 for (i = MAXL - 1; i >= 0; --i){ 76 if (data[i] < x.data[i]){ 77 return *this; 78 } 79 else if (data[i] > x.data[i]){ 80 break; 81 } 82 } 83 return ((*this) - x); 84 } 85 86 bnum & div2(){ 87 int i, carry = 0, tmp; 88 for (i = MAXL - 1; i >= 0; --i){ 89 tmp = data[i] & 1; 90 data[i] = (data[i] + carry) >> 1; 91 carry = tmp * M10; 92 } 93 return *this; 94 } 95 96 bool is_odd(){ 97 return (data[0] & 1) == 1; 98 } 99 100 bool is_zero(){ 101 for (int i = 0; i < MAXL; ++i){ 102 if (data[i]){ 103 return false; 104 } 105 } 106 return true; 107 } 108 }; 109 110 void mulmod(bnum &a0, bnum &b0, bnum &p, bnum &ans){ 111 bnum tmp = a0, b = b0; 112 ans = 0; 113 while (!b.is_zero()){ 114 if (b.is_odd()){ 115 ans = (ans + tmp) % p; 116 } 117 tmp = (tmp + tmp) % p; 118 b.div2(); 119 } 120 } 121 122 void powmod(bnum &a0, bnum &b0, bnum &p, bnum &ans){ 123 bnum tmp = a0, b = b0; 124 ans = 1; 125 while (!b.is_zero()){ 126 if (b.is_odd()){ 127 mulmod(ans, tmp, p, ans); 128 } 129 mulmod(tmp, tmp, p, tmp); 130 b.div2(); 131 } 132 } 133 134 bool MillerRabinTest(bnum &p, int iter){ 135 int i, small = 0, j, d = 0; 136 for (i = 1; i < MAXL; ++i){ 137 if (p.data[i]){ 138 break; 139 } 140 } 141 if (i == MAXL){ 142 // small integer test 143 if (p.data[0] < 2){ 144 return false; 145 } 146 if (p.data[0] == 2){ 147 return true; 148 } 149 small = 1; 150 } 151 if (!p.is_odd()){ 152 return false; // even number 153 } 154 bnum a, s, m, one, pd1; 155 one = 1; 156 s = pd1 = p - one; 157 while (!s.is_odd()){ 158 s.div2(); 159 ++d; 160 } 161 162 for (i = 0; i < iter; ++i){ 163 a = rand(); 164 if (small){ 165 a.data[0] = a.data[0] % (p.data[0] - 1) + 1; 166 } 167 else{ 168 a.data[1] = a.data[0] / M10; 169 a.data[0] %= M10; 170 } 171 if (a == one){ 172 continue; 173 } 174 175 powmod(a, s, p, m); 176 177 for (j = 0; j < d && !(m == one) && !(m == pd1); ++j){ 178 mulmod(m, m, p, m); 179 } 180 if (!(m == pd1) && j > 0){ 181 return false; 182 } 183 } 184 return true; 185 } 186 187 int main(){ 188 bnum x; 189 190 x.read(); 191 puts(MillerRabinTest(x, 5) ? "Yes" : "No"); 192 193 return 0; 194 }
Miller-Rabin素数测试
声明：以上内容来自用户投稿及互联网公开渠道收集整理发布，本网站不拥有所有权，未作人工编辑处理，也不承担相关法律责任，若内容有误或涉及侵权可进行投诉： 投诉/举报 工作人员会在5个工作日内联系你，一经查实，本站将立刻删除涉嫌侵权内容。