素数判定法5 繰り返し二乗法の続きになります。

前回はフェルマーテストの変数をlong long型で宣言していました。

long long型でも巨大な数字になると溢れてしまいます。

巨大な素数に対応するためにフェルマーテストをgmpライブラリーを使ってプログラミングしました。

えらい大変やった。脳みそ使うわ。

fermatgmp.c

#include<stdio.h>
#include<gmp.h>

int main()
{
mpz_t m, n, s, a;
mpz_t one, tmp1, tmp2, tmp3;

printf("数字を入力してください\n");
mpz_init(n);
gmp_scanf("%Zd", n);

mpz_init(m);
mpz_sub_ui(m, n, 1); //m = n - 1;

mpz_init(s);
mpz_set_str(s, "1", 10); // s = 1;

mpz_init(a);
mpz_set_str(a, "2", 10); // a = 2;

mpz_init(one);
mpz_set_str(one, "1", 10);

mpz_init(tmp1);
mpz_init(tmp2);
mpz_init(tmp3);

while(mpz_sgn(m)) { //while(m != 0) {
mpz_and(tmp1, m, one); //if (m & 1) {
if (mpz_sgn(tmp1)) {
mpz_mul(tmp2, s, a); //s = (s * a) % n;
mpz_mod(s, tmp2, n);

}
printf("m = ");
mpz_out_str(stdout, 10, m);
printf(" s = ");
mpz_out_str(stdout, 10, s);
printf(" a = ");
mpz_out_str(stdout, 10, a);
printf("\n");

//printf("m = %lld s = %lld a = %lld\n", m, s, a);

mpz_mul(tmp3, a, a); //a = (a * a) % n;
mpz_mod(a, tmp3, n);

mpz_div_2exp(m, m, 1); //m = m>>1;
}

mpz_out_str(stdout, 10, n);
if(!mpz_cmp_ui(s, 1)) {
printf("は素数です\n");
} else {
printf("は素数じゃないよ\n");
}

mpz_clear(m);
mpz_clear(n);
mpz_clear(s);
mpz_clear(a);
mpz_clear(tmp1);
mpz_clear(tmp2);
mpz_clear(tmp3);
return 0;

}

$ gcc fermatgmp.c -o fermatgmp -lgmp
$ ./fermatgmp
数字を入力してください
26588814640432479998443
26588814640432479998443は素数です

私のMac BookProだと一瞬で計算が終わります。
対応するC言語をコメントで残しています。
検証用に計算過程を表示させています。
邪魔であればコメントアウトしてくださぁい。

フェルマーテストの場合、素数を誤判定する可能性があります。

それで改良されたのがミラーラビン法です。

ミラーラビン法を使うことで誤判定の確率が減ります。

ミラーラビン法へ続く予定。いや続かないかも。多分続きます。