ミラーラビン素数テスト
4977 ワード
カミェル数の存在により,フェルマの小さな定理は1つの数が素数であるか否かを判断できない.
フェルマの小さい定理:pを素数にして、aは任意の整数でそしてaです!三0(mod p)、
a^(p-1)三1(mod p)
//========================================================
カミェル数:それは合数で、1<=a<=nのとき、すべてa^n 3 a(mod n)があります.
//========================================================
カミシェル数のコセット判別法:nを合数とすると、nはカミシェル数当であり、奇数である場合のみ、nを除去する素数pは以下の2つの条件を満たす.
1)p^2不整除n
2)p-1整除n-1
//========================================================
かなり大きな素数を判断するには,合数のラビン・ミラーを用いて定理をテストすることが望ましい.
合数のラビン・ミラー試験の定理:nを奇素数とする、n-1=2^k*qと記す、qは奇数であり、nで整除されないあるaに対して、下記の2つの条件が成立すれば、nは合数である.
a) a^q !三1(mod n);
b)すべてのi=0,1,2,...,k-1, a^((2^i)*q) !三-1(mod n);
//========================================================
ここでは、合数のラビン・ミラーテストの定理aの値を示します.
if n < 1,373,653, it is enough to test a = 2 and 3.
if n < 9,080,191, it is enough to test a = 31 and 73.
if n < 4,759,123,141, it is enough to test a = 2, 7, and 61.
if n < 2,152,302,898,747, it is enough to test a = 2, 3, 5, 7, and 11.
フェルマの小さい定理:pを素数にして、aは任意の整数でそしてaです!三0(mod p)、
a^(p-1)三1(mod p)
//========================================================
カミェル数:それは合数で、1<=a<=nのとき、すべてa^n 3 a(mod n)があります.
//========================================================
カミシェル数のコセット判別法:nを合数とすると、nはカミシェル数当であり、奇数である場合のみ、nを除去する素数pは以下の2つの条件を満たす.
1)p^2不整除n
2)p-1整除n-1
//========================================================
かなり大きな素数を判断するには,合数のラビン・ミラーを用いて定理をテストすることが望ましい.
合数のラビン・ミラー試験の定理:nを奇素数とする、n-1=2^k*qと記す、qは奇数であり、nで整除されないあるaに対して、下記の2つの条件が成立すれば、nは合数である.
a) a^q !三1(mod n);
b)すべてのi=0,1,2,...,k-1, a^((2^i)*q) !三-1(mod n);
//========================================================
ここでは、合数のラビン・ミラーテストの定理aの値を示します.
if n < 1,373,653, it is enough to test a = 2 and 3.
if n < 9,080,191, it is enough to test a = 31 and 73.
if n < 4,759,123,141, it is enough to test a = 2, 7, and 61.
if n < 2,152,302,898,747, it is enough to test a = 2, 3, 5, 7, and 11.
// montgomery (n ^ p) % m, power
unsigned __int64 montgomery(unsigned __int64 n, unsigned __int64 p, unsigned __int64 m)
{
unsigned __int64 r = n % m;
unsigned __int64 tmp = 1;
while (p > 1)
{
if ((p & 1)!=0)
{
tmp = (tmp * r) % m;
}
r = (r * r) % m;
p >>= 1;
}
return (r * tmp) % m;
}
// true:n , false:n
bool R_M_Help(unsigned __int64 a, unsigned __int64 k, unsigned __int64 q, unsigned __int64 n)
{
if ( 1 != montgomery( a, q, n ) )
{
int e = 1;
for ( int i = 0; i < k; ++i )
{
if ( n - 1 == montgomery( a, q * e, n ) )
return false;
e <<= 1;
}
return true;
}
return false;
}
// - true:n , false:n
bool R_M( unsigned __int64 n )
{
if( n < 2 )
throw 0;
if ( n == 2 || n == 3 )
{
return false;
}
if( (n & 1) == 0 )
return true;
// k q, n = 2^k * q + 1;
unsigned __int64 k = 0, q = n - 1;
while( 0 == ( q & 1 ) )
{
q >>= 1;
k++;
}
/*if n < 1,373,653, it is enough to test a = 2 and 3.
if n < 9,080,191, it is enough to test a = 31 and 73.
if n < 4,759,123,141, it is enough to test a = 2, 7, and 61.
if n < 2,152,302,898,747, it is enough to test a = 2, 3, 5, 7, and 11.*/
if( n < 1373653 )
{
if( R_M_Help(2, k, q, n )
|| R_M_Help(3, k, q, n ) )
return true;
}
else if( n < 9080191 )
{
if( R_M_Help(31, k, q, n )
|| R_M_Help(73, k, q, n ) )
return true;
}
else if( n < 4759123141 )
{
if( R_M_Help(2, k, q, n )
|| R_M_Help(3, k, q, n )
|| R_M_Help(5, k, q, n )
|| R_M_Help(11, k, q, n ) )
return true;
}
else if( n < 2152302898747 )
{
if( R_M_Help(2, k, q, n )
|| R_M_Help(3, k, q, n )
|| R_M_Help(5, k, q, n )
|| R_M_Help(7, k, q, n )
|| R_M_Help(11, k, q, n ) )
return true;
}
else
{
if( R_M_Help(2, k, q, n )
|| R_M_Help(3, k, q, n )
|| R_M_Help(5, k, q, n )
|| R_M_Help(7, k, q, n )
|| R_M_Help(11, k, q, n )
|| R_M_Help(31, k, q, n )
|| R_M_Help(61, k, q, n )
|| R_M_Help(73, k, q, n ) )
return true;
}
return false;
}#include<stdio.h>
#include<string.h>
#include<math.h>
#include<time.h>
#include<stdlib.h>
long long prime[]={2,3,5,7,11,13,17,19,23,29};
#define TIME 10 //Miller
long long gcd(long long a,long long b){return b==0?a:gcd(b,a%b);}
long long mod_mult(long long a, long long b, long long n) // (a*b) mod n
{
long long s=0; a=a%n;
while(b){
if (b&1){
s += a;if(s>=n) s-= n;
}
a=a<<1;if(a>=n)a-=n; b=b>>1;
}
return s;
}
long long mod_exp(long long a, long long b, long long n) // (a^b) mod n
{
long long d=1; a=a%n;
while(b>=1) {
if(b&1)d=mod_mult(d,a,n);
a=mod_mult(a,a,n); b=b>>1;
}
return d;
}
bool Wintess(long long a, long long n) // a n Miller
{
long long m, x, y;
int i, j = 0; m = n-1;
while (!(m&1)) // (n-1)=m*(2^j) j m,j=0 m=n-1, 2 m
{
m=m>>1; j++;
}
x = mod_exp(a, m, n);
for (i = 1; i <= j; i++) {
y = mod_mult(x, x, n);
if ((y == 1) && (x != 1) && (x != n - 1)) //
return true; // true ,n
x = y;
}
if (y!=1) return true;
return false;
}
bool miller_rabin(int times, long long n) // n times Miller
{
long long a;
int i;
if (n == 2) return true;
if (n<2||(n&1)==0) return false;
for (i = 1; i <= times; i++) {
a =rand() % (n - 2) + 2;
if (Wintess(a, n)) return false;
}
return true;
}