Lucas-Lehmerの判定法のgmpxxによるプログラム例

Lucas-Lehmer の判定法で Mersenne 数が素数かどうかを判定する、多倍長計算ライブラリ GMP (GNU Multi-Precision Library) を利用した C++ プログラムの例。

2012年4月19日更新 (プログラム修正)。

サンプルプログラム

#include <iostream>
#include <vector>
#include <gmpxx.h>

using namespace std;

// Sieve of Eratosthenes
// list of prime numbers at most n
void sieve( vector<unsigned>& primes, unsigned n )
{
  if (n < 2) { return; }

  // 2 は素数
  primes.push_back(2);
  if (n < 3) { return; }

  // 奇素数を列挙する
  const unsigned lastIdx = (n - 3) / 2;
  vector<char> flag( lastIdx + 1, 1 );  // すべてに 1 をセット
  for ( unsigned i = 0; i <= lastIdx; i++ ) {
    if ( flag[i] ) {
      unsigned p = 2 * i + 3;  // 添字 i は 2*i+3 に対応
      primes.push_back(p);
      if ( p <= n / p ) {
        // 残った最小の数を残し、その倍数をすべて取り除く
        for ( unsigned j = i + p; j <= lastIdx; j += p ) {
          flag[j] = 0;
        }
      }
    }
  }
}

// Lucas-Lehmer Test
bool llTest( unsigned p )
{ 
  // M2 = 2^2 - 1 = 3 は素数
  if ( p == 2 ) { return true; }

  // Lucas-Lehmer Test
  mpz_class mp;
  mpz_mul_2exp( mp.get_mpz_t(), mpz_class(1).get_mpz_t(), p );
  mp -= 1;  // mp = 2^p - 1
  mpz_class x = 4;
  for ( unsigned i = 2; i <= p - 1; i++ ) { 
    x = (x * x - 2) % mp; 
  }

  return ( x == 0 );
}

int main() 
{
  const unsigned n = 100;
  vector<unsigned> primes;

  // n 以下の素数を求める
  sieve( primes, n );

  vector<unsigned>::const_iterator it = primes.begin();
  while ( it != primes.end() ) {
    int p = *it;
    if ( llTest(p) ) { 
      cout << p << " ";
    }
    ++it;
  }

  return 0;
}

実行結果:

2 3 5 7 13 17 19 31 61 89

関連記事

Lucas-Lehmerの判定法
Windows XPにおけるGMPのセットアップ作業のメモ 2012年3月版

【theme : プログラミング
【genre : コンピュータ

プロフィール

よしいず

Author:よしいず
MATHEMATICS.PDFというウェブサイトを運営しています。

管理の都合上、トラックバックとコメントはオフにしてあります。ブログ経験者なら分かっていただけると思いますが、スパム(アダルトやその他の宣伝)ばかりなのが現実です。

リンクは自由です。当サイトの記事に対する間違いの指摘・意見・感想などを述べた記事からのリンクは歓迎です。ただし、ブログ記事アップ直後はミスが多く、頻繁に修正します。場合によっては削除する可能性もあります。その際、何も断りもなく修正・削除しますがご了承ください。内容を参考にする場合には投稿後一週間ほど様子を見てからにしてください(笑)。

記事の間違いを指摘するときは、その具体的箇所、理由(仕様に反するなど)・根拠(参考にした文献など)、代替案(同じ結果を得るための正しいやり方)も教えてください。そうしないと、(指摘される側および第三者はその時点では無知の状態なので、)どこが間違いなのか分かりませんし、本当に間違っているのかどうかが判断・検証できません。実際、間違いだと指摘されたことが結局は正しかったというケースもありますので。

このブログのタイトル一覧

リンク
月別アーカイブ
カテゴリ
最新記事
検索フォーム
RSSリンクの表示