Search on the blog

2011年7月31日日曜日

カーマイケル数の列挙

カーマイケル数を列挙してみる。カーマイケル数といえば、フェルマーテストの偽陽性を引き起こすクセモノ。あらかじめカーマイケル数を列挙しておけば、フェルマーテストと併用することで、正確な素数判定ができる。

いろいろな関数を使うことになった。
  • 最大公約数
  • 最小公倍数
  • 素因数分解
  • 繰り返し二乗法
など基本的な数学系処理の練習にいいかも。。

さて、問題のカーマイケル数の部分だが、カーマイケル関数を使って、


となる最小のmを求める。このmをλ(n)と表記する。ここで、aとnはcoprimeである。カーマイケル数は、


となるnなので、(n-1)がλ(n)で割り切れるとき、nはカーマイケル数の候補となる。さらに、ここから素数を除く必要がある。素数の場合は、φ(n) = n-1となるので、それを除けばいい。逆に、カーマイケル数は、合成数なので、totient functionはn-1より小さくなるため、必ずφ(n) != n-1となり、この除去で偽陰性となることはない。
map<int, int>factorize(int x) {
map<int, int> ret;

for (int i = 2; i*i <= x; i++) {
while (x % i == 0) {
ret[i]++;
x /= i;
}
}

if (x != 1)
ret[x] = 1;

return ret;
}

int mpow(int x, int p) {
int ret = 1;

while (p) {
if (p & 1)
ret = ret * x;
x = x * x;
p >>= 1;
}

return ret;
}

int gcd(int a, int b) {
if (b)
return gcd(b, a % b);
return a;
}

int lcm(int a, int b) {
return (long long int) a * b / gcd(a, b);
}

int carmichaelFunction(int x) {
map<int, int> fact = factorize(x);

int ret = 1;
for (map<int, int>::iterator itr = fact.begin(); itr != fact.end(); itr++) {
int p = itr->first;
int k = itr->second;

int lambda = (p > 2 || k < 3) ? mpow(p, k-1) * (p-1) : mpow(p, k-2);
ret = lcm(ret, lambda);
}

return ret;
}

int eulerFunction(int x) {
map<int, int> fact = factorize(x);

double ret = x;
for (map<int, int>::iterator itr = fact.begin(); itr != fact.end(); itr++)
ret *= (1 - 1./itr->first);

return ret + 0.5;
}

int main() {

for (int i = 2; i < 100000; i++) {
int lambda = carmichaelFunction(i);
int phi = eulerFunction(i);

if ((i-1) % lambda == 0 && phi != i-1) {
printf("%d\n", i);
}
}

return 0;
}

カーマイケル数は、561, 1105, 1729, 2465, 2821, 6601, 8911, 10585, 15841, 29341, 41041, 46657, 52633, 62745, 63973, 75361...と続く。

0 件のコメント:

コメントを投稿