既約分数

id:Nabetaniさんの日記*1にコメントをした。
せっかくなので、自分の日記にも少し書いておこう。

題材は、結城浩さんが2003年の暮れに公開されていた問題*2。(この問題は、ちょっとしたドラマ仕立てで書かれている。でもそこに登場する出題者の少女は、いったい何者なんだかわからずじまいだ。)

問題:正の整数Nが与えられているとき、以下の条件を満たす既約分数p/qを「すべて」求めるアルゴリズムを示してください。条件は:

* p, qは整数(pは0以上で、qは1以上N以下).
* gcd(p, q) = 1 (pとqの最大公約数は1).
* 0 <= p/q <= 1.

整数論では、faray series という数列が知られているらしい。
その数列を使うと効率的なアルゴリズムが作れるそうだ。
が、それはここでは置いておく。

私はHaskell初心者なので、Haskellで表現することで短くわかりやすく書けることが嬉しい。だから、それを試してみたいのだ。

irreducibles q | 1 == q = [(0,1), (1,1)]
               | 1 <  q = irreducibles (q - 1) ++ [(p,q)|p <- [1..q], gcd p q == 1]

これだけだ。
ただし、分数については、ここでは真分数を(分子,分母)の形で表すことにした。

irreducibles q という関数は、 1 == q の場合には、 [(0,1), (1,1)] と定義する。(つまり、0/1と1/1の2つの分数を要素とするリストだ。)

そうでなく、 1 < q の場合には、 irreducibles (q - 1) ++ [(p,q)|p <- [1..q], gcd p q == 1] と定義する。(つまり、qよりひとつだけ小さい分母(q-1)についてのこの関数の結果があれば、それにさらに追加のリストとして、qを分母とするリストを連結したものが、分母qについてのこの関数の結果だ。)

例を挙げる。

irreducibles 1 = [(0,1), (1,1)]
irreducibles 2 = [(0,1), (1,1), (1,2)]
irreducibles 3 = [(0,1), (1,1), (1,2), (1,3), (2,3)]
irreducibles 4 = [(0,1), (1,1), (1,2), (1,3), (2,3), (1,4), (3,4)]
irreducibles 5 = [(0,1), (1,1), (1,2), (1,3), (2,3), (1,4), (3,4), (1,5), (2,5), (3,5), (4,5)]

書き換えると、

irreducibles 1 = [(0,1), (1,1)]
irreducibles 2 = irreducibles 1 ++ [(1,2)]
irreducibles 3 = irreducibles 2 ++ [(1,3), (2,3)]
irreducibles 4 = irreducibles 3 ++ [(1,4), (3,4)]
irreducibles 5 = irreducibles 4 ++ [(1,5), (2,5), (3,5), (4,5)]

1 < q の場合には、再帰的に関数 irreducibles を適用することで帰納的に定義される。

この中の、

[(p,q)|p <- [1..q], gcd p q == 1]

という部分は、評価結果が (p,q) つまり分数を表す要素からなるリストになる。
pには、[1..q]つまり(たとえば[1..5]なら[1,2,3,4,5]というような)1からqまでの連続する整数からなるリストの各要素を順番に適用していく。
しかも、各分数は gcd p q が 1 に等しい、という条件を満たすべきであることを指定しているので、そうでない場合は分数(p,q)は含まれない。

次のように書いたほうが、あるいは多少わかりやすいかもしれない。

irreducibleFractions n = filter irreducible (fractions n)
  where irreducible (p,q) | p /= 0 = gcd p q == 1
                          | p == 0 = q == 1
        fractions q | 1 == q = [(0,1), (1,1)]
                    | 1 <  q = fractions (q - 1) ++ [(p,q)|p <- [1..q]]

ちょっと長くなったが、それでもまだかなり短い。

このコードは、記号名を日本語で書いたとすると、こんな感じになる。

既約分数たち n = ろ過器 既約な (分数たち n)
  ここで 既約な (p,q) | p /= 0 = 最大公約数 p q == 1
                      | p == 0 = q == 1
         分数たち q | 1 == q = [(0,1), (1,1)]
                    | 1 <  q = 分数たち (q - 1) ++ [(p,q)|p <- [1..q]]

このコードでは、既約でないのも既約なのも含めて、指定された分母またはそれ以下の分母を持つ分数すべて(fractions関数)をまず考える。
そのすべての分数のリストに対して、さらに既約なものだけという制約を課した結果が既約分数のリストなのだというわけだ。
(あまりにも愚直なアルゴリズムではある。その分、間違いが入り込む余地は小さい。)

filterはHaskellの標準関数だ。
filterに関数irreducibleを渡し、さらにそれに続いてリストを渡せば、そのリストから条件に適う要素だけが選び取られる。
具体的には、リストの各要素について、最初に渡した関数irreducibleを適用した結果がTrue(真)になるものだけが選び取られる。

テスト実行する場合は、showFractionを定義して、map関数を使って各要素に適用したほうが表示結果がわかりやすい。

showFraction (p, q) = show p ++ "/" ++ show q

例として、N=12を実行する場合。

Main> map showFraction (irreducibles 12)
["0/1","1/1","1/2","1/3","2/3","1/4","3/4","1/5","2/5","3/5","4/5","1/6","5/6","1/7","2/7","3/7","4/7","5/7","6/7","1/8","3/8","5/8","7/8","1/9","2/9","4/9","5/9","7/9","8/9","1/10","3/10","7/10","9/10","1/11","2/11","3/11","4/11","5/11","6/11","7/11","8/11","9/11","10/11","1/12","5/12","7/12","11/12"]

同じく、例として、N=12を実行する場合。

Main> map showFraction (irreducibleFractions 12)
["0/1","1/1","1/2","1/3","2/3","1/4","3/4","1/5","2/5","3/5","4/5","1/6","5/6","1/7","2/7","3/7","4/7","5/7","6/7","1/8","3/8","5/8","7/8","1/9","2/9","4/9","5/9","7/9","8/9","1/10","3/10","7/10","9/10","1/11","2/11","3/11","4/11","5/11","6/11","7/11","8/11","9/11","10/11","1/12","5/12","7/12","11/12"]
Main>

当然ながら、irreducibles 12とirreducibleFractions 12の結果は同じになる。

なんだかくどい文章になった。
いつかHaskellを忘れた頃に読み直してもわかるかな。