prime judge(HOJ1356)

Given a positive integer, your job is writing a program to determine whether it is a prime number or not.

Input

There are several tests. Each test consists of a positive integer n(no more than 2^31) on a single line. Input is terminated by end of file.

Output

For each integer in the input, print a "YES" if it is a prime number, otherwise print a "NO" instead.

Sample Input

4
5
6

Sample Output

NO
YES
NO

　　

　　题目大意就是给出一个不超过2^31的数，来判断它是否为素数，对于此题的规模，一般的素性检测显然不行，要用到Miller rabin, 这个算法主要是基于费尔马小定理，如果 n 为素数，那么对于小于n的数a有a^(n-1) = 1(modn) ('='在这里就代表同余符号)。显然这是一个必要条件，然而只要满足这个条件就基本上是一个素数了，称为‘伪素数’，正确率为1-(1/4)^m，m用不同的基测试的次数，所以多测试几次就可以保证结果的正确了。而初学者在实现Miller rabin算法的时候或多或少都会遇到一些问题， 这里通过此题来给出Miller rabin的一个简单的实现方法。

#include <stdio.h>
#include <stdlib.h>

typedef long long llong;
int Miller_rabbin( llong );
void power( llong, llong, llong );
llong ans;

int main()
{
 int m, times, sel;
 llong n;
 while( scanf( "%lld", &n ) != EOF )
 {
  times = 3, sel = 1; 　//times为每个数进行测试的次数
  if( n <= 1 ) sel = 0;　//sel = 1表示此数为素数
  else if( n == 2 || n == 3 ) sel = 1;
  else
  {
   while( times-- && sel )
    if( !Miller_rabbin( n ) )
     sel = 0;
  }
  if( sel )
   printf( "YES/n" );
  else
   printf( "NO/n" );
 }
 return 0;
}

int Miller_rabbin( llong n )　//核心Miller rabin
{
 llong a;
 a = rand()%(n-1) + 1;
 power( a, n-1, n );
 if( ans % n != 1 ) return 0;
 else return 1;
}

void power( llong a, llong b, llong n ) //快速求幂
{
 if( b == 1 )
 { ans = a; return; }
 power( a, b/2, n );
 if( b % 2 == 0 )
  ans = ans * ans % n;
 else
  ans = ( ans * ans % n ) * a % n;
}