Project Euler – 問題3

13195の素因数は5,7,13,29である。

600851475143の一番大きな素因数を求めよ

せっかくHaskellを使っているので、まずは問題をシンプルかつ直感的にソース コードにしたようなナイーブなやり方を考える方針とする。もしそれが実行速度が遅すぎて答えが求まらないような場合は実行効率も考えたプログラムを考えることにする。

最近の高級言語はレジスタに入らないような大きな数も自動的にうまく扱ってくれるので問題の数字が大きいことは問題ない。これはこういう演習的な問題を解く上での労力の削減に大きく影響すると思う。HaskellだとInteger型で大きさを制限されない整数を扱うことが出来る。ただし、Haskellには強力な型推論があり、実際には型を明示しなくても勝手にうまく扱ってくれる。

素因数とは素数である約数なので、乱暴にいえば素数で順番に割っていき、割り切れた一番大きな素数を求めれば答えということになる。ここで素数をエラトステネスの篩を使って求めることにするが、これがまたHaskell的には得意分 野であると思う。よくある定義は以下のとおり。

sieve (x:xs) = x : sieve [y | y <- xs, y `mod` x /= 0]
primes = sieve $ 2 : [3, 5 ..]

そして問題を上記の乱暴なやり方でプログラムを書くと次のようになる。まず 目的の数より小さい素数のリストをtakeWhileでとりだす。次に、その中から目 的の数の約数をfilterで取り出す。ここで\x -> n `mod` x == 0というのは filterの条件をλ表記の無名の関数で定義している。引数がxで目的の数nをxで 割ったあまりが0 かどうか比較する関数という意味だ。これで目的の数を割り 切ることができる素数を取り出すことが出来る。最後にlastでリストの一番最 後の要素を取り出すことで一番大きな値を求めることが出来る。

この無名関数を用いたところを、Haskellの場合0と比較する関数とあまりを求 める関数を合成するという方法で書くこともできる。個人的にはそういう表記 にはなれていないこともあって、ここでは無名関数で書いた。ちなみに関数合 成を利用してよけいな変数を明示しないスタイルをポイントフリースタイルと いうそうだ。関数合成の演算子は.で、ここでは(== 0) . mod nと書くことで目 的の数nを引数で割ったあまりが0かどうかを比較する関数を表現できる。

euler003 = last $ filter (\x -> n `mod` x == 0) $ takeWhile (<= n) primes
    where n = 600851475143

一見これで問題解決のようだが、実際には上のプログラムは遅すぎて使い物に ならない。もし凄く速い実行環境を持っている場合は上のプログラムでもいい のだが、ここではもう少し高速なバージョンを考えることにする。

上記のプログラムで実行に時間がかかっているのは素数を求める部分である。 エラトステネスの篩で求める場合、求まった全ての素数について、それ以降の 数が倍数かどうかを試すことになるので、目的の数についてn^2程度のオーダー の時間がかかると考えられる。上記のやりかたでは600851475143という非常に 大きな数までの素数を求めなくてはいけないので計算が終了しない。

実はというか当然ではあるが素因数を求めるのに全ての素数で割ってみるとい うやり方はソースコードは単純になるが好ましいものではない。もっと普通の 求め方は、例えば13195の素因数を求める場合、まず5で割れるとわかるので5で 割ると2639になる。次に2639が7で割れるので7で割ると377になる。377が13で 割れるので13で割ると29になる。最大の素因数は29 とわかる。このとき13195 の最大の素因数を求めるのに29までの素数しか必要としない。

上で説明した素因数分解のやり方をそのままプログラムにしたのが関数 primeFactorizeである。引数に素因数分解したい数と素数のリストを渡すと素 因数のリストを返す。ここで素数のリストは無限リストであるがHaskellでは必 要なときに必要な分しか計算されない。

この関数では小さい素数から順に割れるかどうか試していくので、対象になる 数の平方根以上で割りきれることはない。Haskellでは|の後ろに条件を書くこ とで、条件ごとに関数の処理を分けることが出来るが、1行目そのチェックであ り、このときnは素因数なので返り値のリストに追加する。Haskellでは変数の 型を変更することが出来ないので、Floating型の値を扱う関数sqrtと整数とを 同時に扱うにはfromInteger関数で明示的に変換する必要がある。

2行目は素数で割り切れたときで、対象となる数を素数で割った値で次のチェッ クを行う。その素数は素因数なのでリストに追加する。同じ素数で複数回割り 切れることもあるので素数のリストは変更しない。3行目は割り切れなかったと きで、その素数では割り切れないので取り除いて次の素数でチェックする。

Haskellでは:を用いてp:psのようにリストを表記でき、pが最初の要素、psが2 番目以降の要素のリストとなる。関数の引数に(p:ps)とい書いた場合はパター ンマッチとなり引数に渡された値が、pにリストの先頭の要素、psにリストの残 りの要素という風にバインドされる。パターンにマッチしない引数(この場合 では空リスト)の場合はこの関数定義は使用されない。ここでは素数のリスト は無限リストであるが、呼び出しにマッチするパターンがなかった場合はエラー となる。

最後にprimeFactorize関数を使って600851475143 の素因数のリストを求め、そ の一番最後の値(最大の値)が問題の答えとなる。

primeFactorize n (p:ps)
    | fromInteger (p)
      > sqrt (fromInteger n) = [n]
    | n `mod` p == 0         = p : primeFactorize (n `div` p) (p:ps)
    | otherwise              = primeFactorize n ps

euler003 = last $ primeFactorize 600851475143 primes