[ACM]伪素数

Time Limit: 1000MS Memory Limit: 65536KB
Total Submissions: 81 Accepted: 9

Description:

Fermat's theorem states that for any prime number p and for any integer a > 1, a^p == a (mod p). That is, if we raise a to the pth power and divide by p, the remainder is a. Some (but not very many) non-prime values of p, known as base-a pseudoprimes, have this property for some a. (And some, known as Carmichael Numbers, are base-a pseudoprimes for all a.)
Given 2 < p ≤ 1,000,000,000 and 1 < a < p, determine whether or not p is a base-a pseudoprime.

Input:

Input contains several test cases followed by a line containing "0 0". Each test case consists of a line containing p and a.

Output:

For each test case, output "yes" if p is a base-a pseudoprime; otherwise output "no".

Sample Input:

3 2
10 3
341 2
341 3
1105 2
1105 3
0 0

Sample Output:

no
no
yes
no
yes
yes

Hint:

Source:

23 September, 2007 - Waterloo local contest

解题思路：

这题主要是考察素数的判定方法。首先要确保不是素数，然后还满足a^p==a mod p。

#include <iostream>
#include <cstring>
#include <cstdlib>
#include <cstdio>
using namespace std;
typedef unsigned __int64 u64;
const int MAX = 100;
const int MAXN = 30;
u64 len, dig, limit;
u64 mod(u64 a, u64 b, u64 n) {
     if(! a) return 0;
     return (((a & dig) * b) % n + (mod(a >> len, b, n) << len) % n) % n;
}
u64 by(u64 a, u64 b, u64 n) {
      u64 p;
      p = 8, len = 61;
      while(p < n) {
           p <<= 4;
           len -= 4;
      }
      dig = ((limit / p) << 1) - 1;
      return mod(a, b, n);
}
u64 random() {
      u64 a;
      a = rand();
      a *= rand();
      a *= rand();
      a *= rand();
      return a;
}
u64 square_multiply(u64 x, u64 c, u64 n) {
      u64 z = 1;
      while(c) {
           if(c & 1) z = by(z, x, n);
           x = by(x, x, n);
           c >>= 1;
      }
      return z;
}
bool Miller_Rabin(u64 n) {
      if(n < 2) return false;
      if(n == 2) return true;
      if(! (n & 1)) return false;
      u64 i, j, k, m, a;
      m = n - 1;
      k = 0;
      while(m % 2 == 0) {
           m >>= 1; k ++;
      }
      for(i = 0; i < MAX; i ++) {
           a = square_multiply(random() % (n - 1) + 1, m, n);
           if(a == 1) continue;
           for(j = 0; j < k; j ++) {
                if(a == n - 1) break;
                a = by(a, a, n);
           }
           if(j < k) continue;
           return false;
      }
      return true;
}
int main() {
      u64 a, p;
      while(scanf("%I64d %I64d", &p, &a) && (p || a)) {
           if(Miller_Rabin(p)) {
                 printf("no ");
           }else if(square_multiply(a, p, p) == a) {
                 printf("yes ");
           }else {
                 printf("no ");
           }
      }
      return(0);
}

上面的代码中用到了Miller_Rabin素性检测算法，虽然是概率算法，但是在实际中还是够用的，如果用传统的确定性素数检测算法必然导致超时。