Yuki Nakata's Blog

One color just reflects another

カテゴリ:プログラミング > C言語

素数判定プログラム2 C言語GMPライブラリ の続きです。

GMPライブラリを使って素数判定プログラムを書いていきます。

gmpには素数を判定する関数が用意されています。

int mpz_probab_prime_p(mpz_t n, int reps);

戻り値で素数かどうかを判定します。

2が帰ってきたときは絶対素数 (definitely prime)
1のときは素数っぽい、(probably prime)
0のときは合成数つまり素数ではない、です。

repsはミラーラビン法のパラメータですが、通常25を指定しておけばいいそうです。

実行してみてください。一瞬で計算が終わりますよ。

$ gcc prime4.c -o prime4 -lgmp

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

int main(void)
{
mpz_t number;
int answer = 0;

  mpz_init(number);

mpz_set_str(number, "940461086986004843694934910131104208392131088686023657900173332902199657733778583489", 10);
 
answer =  mpz_probab_prime_p(number, 25); 

mpz_out_str(stdout, 10, number);
if (answer) {
printf(" は素数だよ\n");
} else {
printf(" は素数じゃないよ\n");
}
return 0;
}

素数判定プログラム C言語の続きです。

今度はGMPライブラリを使って素数判定プログラムを書いていきます。

GMPプログラミング初めてです。

ぶっつけ本番でも大丈夫でしょう。

変数をそれぞれ保存しないといけないのでプログラムが汚くなってしまいました。

C言語風のコメントも書いておきました。見づらくなったかも。w

GMPライブラリーもマスター出来ました。

//$ gcc prime3.c -o prime3 -lgmp

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

int isPrime(mpz_t number)
{
mpz_t i, four, twomod, threemod, sqrt, tmpmod, tmpmod2, iplus2;

mpz_init(i);
mpz_set_str(i, "5", 10); //10進数 5で初期化 i = 5;

mpz_init(sqrt);
mpz_sqrt(sqrt, number); 

mpz_init(four);
mpz_set_str(four, "4", 10);//10進数 4で初期化 four = 4;
mpz_init(twomod);
mpz_set_str(twomod, "0", 10);//10進数 0で初期化 twomod = 0;

mpz_init(threemod);
mpz_set_str(threemod, "0", 10);//10進数 0で初期化 threemod = 0;

mpz_init(tmpmod);
mpz_init(tmpmod2);
mpz_init(iplus2);

mpz_mod_ui(twomod, number, 2); //twomod = number % 2;
mpz_mod_ui(threemod, number, 3); //threemod = number % 3;

if (mpz_cmp_ui(number, 1) > 0 && 0 < mpz_cmp(four, number)) {
return 1;
} else if (mpz_get_ui(twomod) == 0 || mpz_get_ui(threemod) == 0) {
return 0;
} else {
for (i; 0 <= mpz_cmp(sqrt, i); mpz_add_ui(i, i, 6)) {
mpz_mod(tmpmod, number, i); //tmpmod = number % i;
mpz_add_ui(iplus2, i, 2); //iplus2 = i + 2;
mpz_mod(tmpmod2, number, iplus2); //tmpmod2 = number % iplus2;
if (mpz_get_ui(tmpmod) == 0 || mpz_get_ui(tmpmod2) == 0) {
return 0;
}
}
}

mpz_clear(i);
mpz_clear(sqrt);
mpz_clear(four);
mpz_clear(twomod);
mpz_clear(threemod);
mpz_clear(tmpmod);
mpz_clear(tmpmod2);
mpz_clear(iplus2);
return 1;
}

int main(void)
{
mpz_t number;
int answer = 0;

  mpz_init(number);

mpz_set_str(number, "150094772735952593", 10);

answer = isPrime(number);

mpz_out_str(stdout, 10, number);
if (answer) {
printf(" は素数だよ\n");
} else {
printf(" は素数じゃないよ\n");
}
mpz_clear(number);
return 0;
}

150094772735952593 は18桁の素数です。

ここに判定する数字を入力してください。
 
素数サンプル 楽勝レベル
31389448217
282437925089
1853028778786433
5559077746424707
150094772735952593

PCだと厳しい 一気に時間がかかります。終わらないかもね。
8862938260389989451257
26588814640432479998443

実はGMPライブラリには素数を判定する関数が用意されています。

それを呼べば一発で素数かどうかを判定することができます。

次回に続きます。 素数判定プログラム3 C言語GMPライブラリ

参考サイト
多倍長整数演算ライブラリ (GNU MP) 

素数一覧 1000万 

素数表 

素数を判定するプログラムを書いていきます。

単純にC言語で書いたのがコレ。

//$ gcc prime1.c -o prime1

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

int isPrime(unsigned long number)
{
int i;
for (i = 2; i < number; i++) {
if (number % i == 0) {
return 0;
}
}
return 1;
}

int main(void)
{
unsigned long number = 9999991;
int answer = 0;
answer = isPrime(number);
if (answer) {
printf("%lu は素数だよ\n", number);
} else {
printf("%lu は素数じゃないよ\n", number);
}
return 0;
}

numberのところに判定する数字を入力してください。

3箇所アルゴリズムを改良して早くしていきます。 

1. まず2や3で割り切れるなら素数ではありません。
 候補ではない数字を最初に除外します。これをエラトステネスの篩といいます。

2. 判定する数字を最後までチェックする必要はありません。
 判定する数字を2つの積で表します。
 たとえば100を2つの積で表してみましょう。
 2 * 50     2で素数かどうかチェック済み
    4 * 25 4で素数かどうかチェック済み
    5 * 20 5で素数かどうかチェック済み
   10 * 10 10で素数かどうかチェック済み

 つまりわざわざ100まで調べなくても良いのです。100のルートを取り10を上限とします。
 素数かどうかを調べるには10以下まで調べれば良いです。
  
3. 素数の数字に注目する。
 5以上100までの素数を書いてみます。
 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97
 素数というのは6の倍数の隣の数字になります。確認してください。
 6k ± 1で表せます。

 この3つを考慮してプログラムを書くと次になります。

//$ gcc prime2.c -o prime2

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

int isPrime(unsigned long number)
{
unsigned int i;
if (number <= 3 && number > 1) {
return 1;
} else if (number % 2 == 0 || number % 3 == 0) {
return 0;
} else {
for (i = 5; i * i <= number; i += 6) {
if (number % i == 0 || number % (i + 2) == 0) {
return 0;
}
}
}
return 1;
}

int main(void)
{
unsigned long number = 9999991;
int answer = 0;
answer = isPrime(number);
if (answer) {
printf("%lu は素数だよ\n", number);
} else {
printf("%lu は素数じゃないよ\n", number);
}
return 0;
}


これでだいぶ早くなったはずです。
for分がややこしいかもしれませんが、一度紙に書いてトレースして確認してみてください。
 
次回はGMPライブラリーを使って素数判定プログラムに改良していきます。

素数判定プログラム C言語 GMPライブラリー

↑このページのトップヘ