素数

エラトステネスのふるい

n以下の素数をすべて見つける方法に、「エラトステネスの篩(ふるい)」というものがあります。何も難しいことではなく、ある素数の倍数を消していくだけのことです。たとえば、100以下の素数でしたら、2から始めて、3,5・・・の倍数を消していけばよいのです。難しいのはそれを抽象的に考えた場合ですが、高校数学ではそこまでは扱いません。

さて、エラトステネスの篩で素数をリストアップするプログラムを作ります。VBAなどでもできます。そのあと、「リスト」処理ができるCindyScriptならではのもっとエレガントな方法.を紹介します。

まずは、基本形。2から始めて√nまで、倍数を消すことを「繰り返す」方法です。

はじめに、2からnまでの自然数を用意します。それを「リスト」にします。「2からnまでの自然数をリストアップする」と考えてよいでしょう。CindyScriptでは、

suu=2..n

で2からnまでの自然数のリストが作れます。ピリオドが2個です。具体的に使う場合は、nに目的の数を代入します。作った結果をprintln()で表示してみましょう。

エラトステネスのふるいでは、この中から「倍数」を消していきます。その手順を具体的に考えてみましょう。

(1) 2の倍数を消す

(2) 3の倍数を消す

(3) 4は、(1)で消えているので、次の5の倍数を消す

(4) つぎの6は(1)で消えているので、次の7の倍数を消す

これを続けて、10までやればおわり。実際には10は消えていますが。

これでよいですか? いいえ、もっと丁寧に考える必要があります。(1)はどのようにしますか?

(1-1) リストからひとつ取り出す。

(1-2) それが2の倍数かどうか調べる。

(1-3) 2の倍数であればそれをリストから消す。

という手続きが必要ですね。これでよいですか? まだです。(1)の「ひとつ」というのは何番目ですか?

(1-1) リストからk番目の要素をひとつ取り出す。

(1-2) それが2の倍数かどうか調べる。

(1-3) 2の倍数であればそれをリストから消す。

これを繰り返します。「繰り返す」ということは、(1-1)のkを 1,2,3 ・・ と変えていくことになります。ところがここで問題が生じます。(1-3)でリストから消すと、次の要素が繰り上がります。消さなければ繰り上がりません。となると「次の数」はどれにすればよいでしょう。

具体的に考えてみましょう。

2,3,4,5,6,7・・

で、2番目の要素3を取り出して2の倍数かどうか調べます。2の倍数ではないので消しません。したがって、リストはそのままです。

次に3番目の要素4を取り出して2の倍数かどうか調べます。2の倍数ですので消します。すると、リストは次のようになります。

2,3,5,6,7・・

すると、次に取り出すのは4番目の要素ではなく、やっぱり3番目の要素です。

ややこしくなってしまいました。・・・・・ あきらめずにがんばりましょう。

対策を2つ考えます。

その1 k番目を消したら次もk番目、そうでなければ次は(k+1)番目とする。「もし〜ならば〜」なので、条件判断の if を使えばよい。

その2 直接消してしまわないで、もうひとつリストを用意して、そちらに消したかどうかの情報を書く。消したら0,残っていれば1にすればよい。

では、やってみましょう。まず「その1」から。

--------------------------------------------------------------------------

n=100;

suu=2..n;

i=1;

a=suu_i;

while(a<sqrt(n),

 k=1;

 while(i+k<=length(suu),

    b=suu_(i+k);

    if(mod(b,a)==0,

        suu=remove(suu,b),

        k=k+1;

    );

 );

 println(suu+","+length(suu));

 i=i+1;

 a=suu_i;

);

--------------------------------------------------------------------------

結構ややこしいですが、動作を追ってみてください。

a=suu_i; は、suu のリストから i 番目の要素を取り出して、aに代入します。

while(a<sqrt(n),・・・ は、aの値が√nより小さい間、・・・を繰り返します。

if(mod(b,a)==0 の mod(b,a) はbをaで割った余りを求める関数です。==0 ということは、余りが0,つまり割り切れるということです。

println(suu+","+length(suu)); は、残っている数のリストと、その個数を表示します。

つぎに「その2」です。

--------------------------------------------------------------------------

n=100;

suu=1..n;

flag=apply(1..n,#=1);

flag_1=0;

repeat(sqrt(n),i,

 if(flag_i!=0,

   a=suu_i;

   repeat(n-i,k,

     b=suu_(i+k);

     if(mod(b,a)==0,

        flag_(i+k)=0;

     );

   );

 );

);

repeat(n,i,

  if(flag_i==1,

     print(suu_i+",");

  );

);

--------------------------------------------------------------------------

flag=apply(1..n,#=1);

flag_1=0;

の2行で、第1要素だけ0,あとは1のflagという名前のリストを作ります。

repeat の繰り返しが3つあります。

2番目がエラトステネスの篩です。少し無駄がありますが、aからあとの数を全部調べて、倍数ならflagリストの該当要素を0にしています。

4番目は結果の表示です。flagリストのうち、1のままのものが消去されずに残った素数です。

リスト処理の関数を使って素数のリストを作る

エラトステネスの篩では、while または repeat の繰り返しを使って、素数だけを残していきます。しかし、CindyScriptにはリストを処理する関数が用意されており、これを使うとプログラムが驚くほど簡単になります。

約数のリストを作る

まず、ある数の約数のリストを作ることを考えます。自然数nの約数は、nを割り切ることのできる数です。したがって、あるnに対し、1からnまでの数について、nを割りきる数を選び出します。リストの中からある条件をみたすものを選び出すには select という関数を使います。次の1行で、nの約数のリストが yaku に作られます。

yaku=select(1..n,mod(n,#)==0);

どういう仕掛けになっているのかを見ていきましょう。

select(list,boolexpr) という関数は、list のリストから、boolexpr の条件を満たすものを選び出す関数です。今の場合、1..n が対象となるリストになります。すなわち、1..n によって、1からnまでの自然数を要素とするリストが内部的に作られ、そのリストを対象とすることになります。boolexpr の条件は mod(n,#)==0 ですが、#は「実行変数」と呼ばれるもので、やはり内部的に作られ、ここではリストの要素を指すことになります。select()関数は、リストの最初から順に要素を#へ代入し、条件を満たすかどうかを判定していきます。その結果、1からnまでの自然数のうち、nの約数が選ばれて、yakuに代入されます。yakuはリストです。

さて、多様なnの値に対応するために、これを関数として定義します。CindyScriptでは、先ほどのコードと組み合わせて、次のように定義できます。

yakusuu(n):=select(1..n,mod(n,#)==0);

さきほどと違うのは、イコール「=」がコロン・イコール「:=」に変わっていることです。これで関数の定義ができるのです。呼び出し方は

yaku=yakusuu(72)

です。これでリストyakuに72の約数が入ります。

素数の判定

1とn以外に約数を持たない数nが素数です。CindyScriptのリスト処理では、この定義をそのまま使います。すなわち、約数のリストを作っておいて、その個数が2(1とnだけ)ならば素数であると判定できるわけです。たとえば、191が素数かどうかを判定してみましょう。

さらに、これを応用して、100までの素数をリストアップすることもできます。使うのは、やはりselect()です。つまり、1から100までの自然数のリストを作っておき、そこから、素数であるものをselectすればいいのです。

このように、リストを使うと繰り返し処理が簡単にできます。

素因数分解

与えられた数を素因数分解します。2から始めて、順に素数で割っていきます。割り切れたら、その商をさらに同じ素数で割ります。割り切れなかったら次の素数で割ります。これを繰り返していくのですが、どこまで繰り返せばよいでしょう。nまで? いいえ、n以下の素数のうち一番大きな数まででよいですね。実際にはそれまでには商が1になりますので、そこで終わりになります。(nが素数の場合は最後までやることになります)そこで、次のように考えましょう。

(1) n以下の素数のリストを作る

(2) 2から始めて上の手続きを行なう。そのとき、ある素数で何回割り切れたか記録する。

(3) 商が1になったら終わり、結果を表示する。

実際にスクリプトを書いてみると、(3)の結果の表示が結構面倒です。CindyScriptではTeX形式での数式表現ができますので、指数を使って結果を表示することができます。

次のスクリプトがその例です。

--------------------------------------------------------------------------

yakusuu(n):=select(1..n,mod(n,#)==0);    // 関数の定義

isprime(n):=length(yakusuu(n))==2;

prime(n):=select(1..n,isprime(#));

// n に素因数分解したい数を代入して実行

n=360;

// ーーーーーーーーー

plist=prime(n);              // n以下の素数のリスト

pn=length(plist);           // n以下の素数の数

kosu=apply(1..pn,#=0);    // 割り切れた素数の個数を記録するリスト

m=n;            // m は最初がn そのあとは素数で割った商

i=1;

while(m>1,

  p=plist_i;                     // 素数のリストから一つ取り出して

  if(mod(m,p)==0,           // 割り切れたら

    m=div(m,p);               // mは商にして

    kosu_i=kosu_i+1,       // 割り切れた素数をカウント

    i=i+1;                        // 割り切れなかったら、次の素数に進むために番号を増やす

  );

);

disp="";                          // このあとは結果の表示をするスクリプト Texの式を作る

repeat(pn,i,

  if(kosu_i!=0,

    if(length(disp)>0,disp=disp+"・");

    disp=disp+"$"+text(plist_i);

    if(kosu_i>1,disp=disp+"^"+text(kosu_i));

    disp=disp+"$";

  );

);

drawtext((0,1),n+"="+disp,size->20);     // 描画面に結果を表示

--------------------------------------------------------------------------

戻る