コンピュータ「オイラと数学しよう」 ――オイラープロジェクトの問題にチャレンジ(その2)無償ソフトで技術計算しよう【プログラミング応用編】(3)(2/3 ページ)

» 2014年12月12日 10時00分 公開
[伊藤孝宏,MONOist]

問題6

Problem 6  Sum square difference

The sum of the squares of the first ten natural numbers is, 12 + 22 + ... + 102 = 385

The square of the sum of the first ten natural numbers is, (1 + 2 + ... + 10)2 = 552 = 3025

Hence the difference between the sum of the squares of the first ten natural numbers and the square of the sum is 3025-385 = 2640.

Find the difference between the sum of the squares of the first one hundred natural numbers and the square of the sum.

 意味は、「1〜10までの自然数の二乗和は385になります。一方、和の二乗は3025になります。従って、和の二乗と二乗和との差は2640となります。1〜100までの自然数の和の二乗と二乗和との差を求めなさい」。

 1〜10までの配列nを作り、和の二乗は(sum(n))^2となり、二乗和はドット演算子を使い、sum(n.^2)で求まります。答えは下記のようになります。これもサービス問題ですね。

--> n=1:10;(sum(n))^2-sum(n.^2)
ans = 2640
--> n=1:100;(sum(n))^2-sum(n.^2)
ans = 25164150

問題7

Problem 7 10001st prime

By listing the first six prime numbers: 2, 3, 5, 7, 11, and 13, we can see that the 6th prime is 13.

What is the 10 001st prime number?

 意味は、「最初から6番目までの素数は、2、3、5、7、11、13となり、6番目の素数は13です。では、10001番目の素数を求めなさい」。

 エラトステネスの篩を使います。今回の場合、n番目の素数と、あらかじめ素数の数が分かっているので、ps=[2,3,5,7,11,zeros(1,n-5)]として、素数の配列psを用意します。既に分かっている素数2、3、5、7、11の後は0がnまで続く配列で、0の部分に素数が格納されていきます。nnを11から増加させ、psで割り切れるかどうかで素数であるかを調べます。

 ここで、偶数はあきらかに素数にはならないので、nnは11から2ずつ増加させます。素数かどうかの判別は、find(mod(nn,ps)==0)として、割り切れる、つまりmodが0となる位置を調べます。mod(nn,0)はnnとなるので、まだ素数が格納されていない部分でmodが0となることはありません。

 modが0となる要素が見つからない場合、空の配列が返るので、「isempty」コマンドで検出します。isemptyコマンドは、引数が空配列だと1を返し、それ以外では0を返します。modに0がない、すなわち、素数の場合は、psに追加し、cntを1増加させます。cntがnになるまで、以上の操作を繰り返すと、n個の素数が得られ、目的の素数は素数の配列psのn番目の要素になります。

function p7=ex331(n)
    ps=[2,3,5,7,11,zeros(1,n-5)];
    cnt=5;nn=11;
    while(cnt~=n)
        if(isempty(find(mod(nn,ps)==0)))
            cnt=cnt+1;ps(cnt)=nn;
        end
        nn=nn+2;
    end
    p7=ps(n);
ex331.m

>>「ex331.m」ダウンロード

 例で試してみると、下記のように、正しいことが分かります。

--> ex331(6)
ans = 13

 本番で試してみると、残念ながら、FreeMat自体の処理速度が遅いせいなのか、答えが得られませんでした。そこで、また掟(おきて)破りのようで申し訳ないのですが、MATLABでex331.mを動作させてみました。答えは以下が得られました。

>> ex331(10001)
ans = 104743

Copyright © ITmedia, Inc. All Rights Reserved.