リスト : 九九表 def ninetynine for x in 1 .. 9 for y in 1 .. 9 printf "%2d ", x * y end print "\n" end end def ninetynine1 (1..9).each {|x| (1..9).each {|y| printf "%2d ", x * y } print "\n" } end def ninetynine2 (1..9).map {|x| (1..9).map {|y| x * y}} end
九九表は for 文による二重ループで作成することができます。
irb> ninetynine 1 2 3 4 5 6 7 8 9 2 4 6 8 10 12 14 16 18 3 6 9 12 15 18 21 24 27 4 8 12 16 20 24 28 32 36 5 10 15 20 25 30 35 40 45 6 12 18 24 30 36 42 48 54 7 14 21 28 35 42 49 56 63 8 16 24 32 40 48 56 64 72 9 18 27 36 45 54 63 72 81 => 1..9
ninetynine1() は for 文のかわりにメソッド each() を使う方法です。配列に格納する場合は map() を使うと簡単です。
irb> ninetynine2 => [[1, 2, 3, 4, 5, 6, 7, 8, 9], [2, 4, 6, 8, 10, 12, 14, 16, 18], [3, 6, 9, 12, 15, 18, 21, 24, 27], [4, 8, 12, 16, 20, 24, 28, 32, 36], [5, 10, 15, 20, 25, 30, 35, 40, 45], [6, 12, 18, 24, 30, 36, 42, 48, 54], [7, 14, 21, 28, 35, 42, 49, 56, 63], [8, 16, 24, 32, 40, 48, 56, 64, 72], [9, 18, 27, 36, 45, 54, 63, 72, 81]]
整数 n から m までの数列の和、二乗の和、三乗の和を求めることを考えます。この場合、高階関数 sum_of() を作ると簡単です。
リスト : 整数列の和 def sum_of(n, m, &func) (n..m).reduce(0) {|a, x| a + func.call(x)} end
irb> sum_of(1, 100){|x| x} => 5050 irb> sum_of(1, 100){|x| x * x} => 338350 irb> sum_of(1, 100){|x| x * x * x} => 25502500 irb> sum_of(1, 100, &:itself) => 5050
sum_of() のブロックに引数をそのまま返す関数、二乗する関数、三乗する関数を渡すだけです。メソッド object.itself はレシーバーの object をそのまま返します。関数型言語では、引数をそのまま返す関数を「恒等関数 (identity function)」といいます。
ところで、次に示す数列の和を求める公式を使うと、もっとスマートに値を求めることができます。
興味のある方はプログラムを作ってみてください。
下記に示すデータの度数分布表と累積度数表を求めるプログラムを作ります。
リスト : 身長のデータ Height = [ 148.7, 149.5, 133.7, 157.9, 154.2, 147.8, 154.6, 159.1, 148.2, 153.1, 138.2, 138.7, 143.5, 153.2, 150.2, 157.3, 145.1, 157.2, 152.3, 148.3, 152.0, 146.0, 151.5, 139.4, 158.8, 147.6, 144.0, 145.8, 155.4, 155.5, 153.6, 138.5, 147.1, 149.6, 160.9, 148.9, 157.5, 155.1, 138.9, 153.0, 153.9, 150.9, 144.4, 160.3, 153.4, 163.0, 150.9, 153.3, 146.6, 153.3, 152.3, 153.3, 142.8, 149.0, 149.4, 156.5, 141.7, 146.2, 151.0, 156.5, 150.8, 141.0, 149.0, 163.2, 144.1, 147.1, 167.9, 155.3, 142.9, 148.7, 164.8, 154.1, 150.4, 154.2, 161.4, 155.0, 146.8, 154.2, 152.7, 149.7, 151.5, 154.5, 156.8, 150.3, 143.2, 149.5, 145.6, 140.4, 136.5, 146.9, 158.9, 144.4, 148.1, 155.5, 152.4, 153.3, 142.3, 155.3, 153.1, 152.3 ]
階級 度数 累積度数 ------------------------ 130 - 135 1 1 135 - 140 6 7 140 - 145 12 19 145 - 150 25 44 150 - 155 32 76 155 - 160 17 93 160 - 165 6 99 165 - 170 1 100
階級はデータの範囲を表します。この表では x cm 以上 y cm 未満を x - y で表しています。度数はその階級に出現したデータの個数です。度数を示してある表のことを「度数分布表」といいます。累積度数はその階級までの度数を全部加えたものです。累積度数を示してある表を「累積度数分布表」といいます。
リスト : 度数分布表と累積度数表 def make_frequency m = 8 freq = [0] * m # 度数分布表 cum = [0] * m # 累積度数表 low = 130.0 z = 5.0 # 度数分布表の作成 for x in Height for j in 0 ... m if x < low + z * (j + 1) freq[j] += 1 break end end end # 累積度数表の作成 cum[0] = freq[0] for i in 1 ... m cum[i] = cum[i - 1] + freq[i] end # 表示 for i in 0 ... m printf "%.1f - %.1f | ", low + z * i, low + z * (i + 1) printf "%3d %3d\n", freq[i], cum[i] end end
配列 freq が度数分布表、cum が累積度数表、変数 low が階級の下限値、z が階級の幅を表します。freq は 0 で初期化します。度数分布表の作成は簡単です。最初の for ループで height の要素を取り出します。次の for ループで階級を求めます。変数 j が階級を表し、その上限値は low + z * (j + 1) で求めることができます。height[i] がこの値よりも小さい場合、その要素は階級 j であることがわかります。freq[j] の値を +1 して、for ループを脱出します。
累積度数表の作成も簡単です。cum[0] を freq[0] で初期化します。あとは、度数分布表の値 freq[i] と累積度数表の値 cum[i - 1] を足し算していくだけです。最後に、for ループで 2 つの表を出力します。
実行結果は次のようになります。
irb> make_frequency 130.0 - 135.0 | 1 1 135.0 - 140.0 | 6 7 140.0 - 145.0 | 12 19 145.0 - 150.0 | 25 44 150.0 - 155.0 | 32 76 155.0 - 160.0 | 17 93 160.0 - 165.0 | 6 99 165.0 - 170.0 | 1 100 => 0...8
総和、最小値、最大値は reduce() を使うと簡単です。
リスト : 総和、最大値、最小値、平均値 def sum(a) a.reduce {|m, x| m + x} end def max(a) a.reduce {|m, x| m < x ? x : m} end def min(a) a.reduce {|m, x| m > x ? x : m} end def average(a) sum(a).fdiv(a.size) end
irb> a = [5,6,4,7,3,8,2,9,1,0,10] => [5, 6, 4, 7, 3, 8, 2, 9, 1, 0, 10] irb> sum(a) => 55 irb> min(a) => 0 irb> a.min => 0 irb> max(a) => 10 irb> a.max => 10 irb> average(a) => 5.5
Ruby の Enumerable モジュールには最大値と最小値を求めるメソッド max(), min() があります。C言語のような三項演算子 ?, : もあります。
ここでは xn の n を整数とします。一番簡単な方法は x を n 回乗算することです。
リスト : 累乗 def power0(x, n) a = 1 n.times {|_| a *= x} a end def power1(x, n) (1..n).reduce(1) {|a, _| a * x} end
power0() は単純な繰り返しで、power1() は reduce() を使って累乗を求めます。
irb> power0 2, 32 => 4294967296 irb> power0 2, 64 => 18446744073709551616 irb> power0 2, 128 => 340282366920938463463374607431768211456 irb> power1 2, 32 => 4294967296 irb> power1 2, 64 => 18446744073709551616 irb> power1 2, 128 => 340282366920938463463374607431768211456
式を変形するともっと少ない回数で求めることができます。
x ** 4 = (x ** 2) ** 2 -> 2 回 x ** 8 = (x ** 4) ** 2 -> 3 回 x ** 16 = (x ** 8) ** 2 -> 4 回 一般化すると x ** n = (x ** (n / 2)) ** 2; (n は偶数) x ** n = ((x ** (n / 2)) ** 2) * x; (n は奇数)
階乗の計算では n を n - 1 の計算に置き換えていきますが、累乗の場合は n を n / 2 に置き換えていくことができます。n が半分になっていくので、ひとつずつ n を減らすよりも減少の度合いは大きくなります。その分だけ計算回数が少なくなるわけです。
これを再帰呼び出しでプログラムすると次のようになります。
リスト : 累乗 (2) def power2(x, n) return 1 if n == 0 v = power2(x, n / 2) v *= v v *= x if n % 2 == 1 v end
irb> power2 2, 128 => 340282366920938463463374607431768211456
最初の if 文が再帰呼び出しの停止条件です。次に、n を 2 で割った値で power2() を再帰呼び出しします。この返り値を 2 乗して、n が奇数ならば x をさらに掛け算します。最後に計算結果 v を返します。
ハノイの塔は、棒に刺さっている大きさが異なる複数の円盤を、次の規則に従ってほかの棒に移動させるパズルです。
ハノイの塔は、再帰を使えば簡単に解ける問題です。
リスト : ハノイの塔 def hanoi(n, from, to, via) hanoi(n - 1, from, via, to) if n > 1 printf("disk %d : %c -> %c\n", n, from, to) hanoi(n - 1, via, to, from) if n > 1 end
irb> hanoi 3, ?a, ?b, ?c disk 1 : a -> b disk 2 : a -> c disk 1 : b -> c disk 3 : a -> b disk 1 : c -> a disk 2 : c -> b disk 1 : a -> b => nil
n 番目のフィボナッチ数 F(n) は次の数式で定義されます。
フィボナッチ数も再帰呼び出しを使えば簡単にプログラムできます。
リスト : フィボナッチ数列 def fibo(n) return n if n < 2 fibo(n - 1) + fibo(n - 2) end
irb> 20.times {|n| puts fibo(n)} 0 1 1 2 3 5 8 13 21 34 55 89 144 233 377 610 987 1597 2584 4181 => 20
関数 fibo() は自分自身を 2 回呼び出しています。これを「二重再帰」といいます。fibo() は同じ値を何回も求めているため、効率はとても悪いのです。この場合、二重再帰を「末尾再帰」に変換すると高速化することができます。
リスト : フィボナッチ数列 (2) def fibo1(n, a = 0, b = 1) return a if n == 0 fibo1(n - 1, b, a + b) end
irb> fibo(30) => 832040 irb> fibo1(30) => 832040 irb> fibo1(50) => 12586269025 irb> fibo1(100) => 354224848179261915075
累算変数 a と b の使い方がポイントです。現在のフィボナッチ数を変数 a に、ひとつ先の値を変数 b に格納しておきます。あとは a と b を足し算して、新しいフィボナッチ数を計算すればいいわけです。
フィボナッチ数列を生成するジェネレータはクラス Enumerator を使うと簡単です。
リスト : フィボナッチ数列 (3) Fibos = Enumerator.new {|y| a, b = 0, 1 loop { y << a a, b = b, a + b } }
irb> Fibos => #<Enumerator: ...> irb> 20.times {|_| puts Fibos.next} 0 1 1 2 3 5 8 13 21 34 55 89 144 233 377 610 987 1597 2584 4181 => 20 irb> fibos.first(20) => [0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181] irb> fibos.lazy.drop(20).first(20) => [6765, 10946, 17711, 28657, 46368, 75025, 121393, 196418, 317811, 514229, 832040, 1346269, 2178309, 3524578, 5702887, 9227465, 14930352, 24157817, 39088169, 63245986]
Enumerator で作成したジェネレータはメソッド next() で値を順番に取得することができます。また、メソッド first() や take() で先頭から n 個の要素を取り出すことができます。
ただし、このプログラムは無限数列になっているので、先頭から n 個の要素を取り除くメソッド drop() を使うと、要素を n 個取り除いたあとの配列 (無限個) を生成しようとするため、irb がフリーズしてしまいます。この場合、メソッド Enumerable#lazy を使って遅延リストに変換すると動作します。
組み合わせの数 \({}_n \mathrm{C}_r\) を求めるには、次の公式を使えば簡単です。
皆さんお馴染みの公式ですね。ところが、整数値の範囲が限られているプログラミング言語では、この公式を使うと乗算で「桁あふれ」を起こす恐れがあります。Ruby は多倍長演算をサポートしているので、桁あふれを心配する必要はありません。
この公式をそのままプログラムすることもできますが、次の式を使うともっと簡単にプログラムできます。
この式は \({}_n \mathrm{C}_r\) と \({}_n \mathrm{C}_{r-1}\) の関係を表しています。あとは、再帰定義を使って簡単にプログラムできます。
リスト : 組み合わせの数 def comb(n, r) return 1 if n == r || r == 0 comb(n, r - 1) * (n - r + 1) / r end
irb> comb(10, 5) => 252 irb> comb(30, 15) => 155117520 irb> comb(50, 25) => 126410606437752
下図に「パスカルの三角形」を示します。
1 0C0 / \ / \ 1 1 1C0 1C1 / \ / \ / \ / \ 1 2 1 2C0 2C1 2C2 / \ / \ / \ / \ / \ / \ 1 3 3 1 3C0 3C1 3C2 3C3 / \ / \ / \ / \ / \ / \ / \ / \ 1 4 6 4 1 4C0 4C1 4C2 4C3 4C4 図 : パスカルの三角形
パスカルの三角形は、左側の図のように両側がすべて 1 で、内側の数はその左上と右上の和になっています。これは式 \((a + b)^n\) を展開したときの各項の係数を表しています。そして、その値は右側の図のように組み合わせの数 \({}_n \mathrm{C}_r\) に対応しています。
きれいな三角形にはなりませんが、簡単なプログラムを示します。
リスト : パスカルの三角形 def pascal(x) for n in 0 .. x for r in 0 .. n print comb(n, r), " " end print "\n" end end
irb> pascal(16) 1 1 1 1 2 1 1 3 3 1 1 4 6 4 1 1 5 10 10 5 1 1 6 15 20 15 6 1 1 7 21 35 35 21 7 1 1 8 28 56 70 56 28 8 1 1 9 36 84 126 126 84 36 9 1 1 10 45 120 210 252 210 120 45 10 1 1 11 55 165 330 462 462 330 165 55 11 1 1 12 66 220 495 792 924 792 495 220 66 12 1 1 13 78 286 715 1287 1716 1716 1287 715 286 78 13 1 1 14 91 364 1001 2002 3003 3432 3003 2002 1001 364 91 14 1 1 15 105 455 1365 3003 5005 6435 6435 5005 3003 1365 455 105 15 1 1 16 120 560 1820 4368 8008 11440 12870 11440 8008 4368 1820 560 120 16 1 => 0..16
for ループで二重ループを構成しています。最初の for ループで変数 n の値を 0 から x まで +1 ずつ増やし、次の for ループで変数 r の値を 0 から n まで +1 ずつ増やします。あとは comb() で nCr の値を計算するだけです。
comb() を使わないでプログラムを作ることもできます。
リスト : パスカルの三角形 (2) def pascal1(x) a = [1] (x + 1).times {|_| print a, "\n" a.push 0 a = a.zip(a.reverse).map {|xs| xs[0] + xs[1]} } end
irb> pascal1(16) [1] [1, 1] [1, 2, 1] [1, 3, 3, 1] [1, 4, 6, 4, 1] [1, 5, 10, 10, 5, 1] [1, 6, 15, 20, 15, 6, 1] [1, 7, 21, 35, 35, 21, 7, 1] [1, 8, 28, 56, 70, 56, 28, 8, 1] [1, 9, 36, 84, 126, 126, 84, 36, 9, 1] [1, 10, 45, 120, 210, 252, 210, 120, 45, 10, 1] [1, 11, 55, 165, 330, 462, 462, 330, 165, 55, 11, 1] [1, 12, 66, 220, 495, 792, 924, 792, 495, 220, 66, 12, 1] [1, 13, 78, 286, 715, 1287, 1716, 1716, 1287, 715, 286, 78, 13, 1] [1, 14, 91, 364, 1001, 2002, 3003, 3432, 3003, 2002, 1001, 364, 91, 14, 1] [1, 15, 105, 455, 1365, 3003, 5005, 6435, 6435, 5005, 3003, 1365, 455, 105, 15, 1] [1, 16, 120, 560, 1820, 4368, 8008, 11440, 12870, 11440, 8008, 4368, 1820, 560, 120, 16, 1] => 17
前段の値を配列 a に格納します。a に 0 を追加して、それを反転した配列の要素と足し算すれば、次段の値を求めることができます。もちろん、配列を一つで済ますことも可能です。興味のある方は挑戦してみてください。
たとえば、1 から 5 までの数字の中から 3 個を選ぶ組み合わせは次のようになります。
(1, 2, 3), (1, 2, 4), (1, 2, 5), (1, 3, 4), (1, 3, 5), (1, 4, 5), (2, 3, 4), (2, 3, 5), (2, 4, 5), (3, 4, 5)
最初に 1 を選択した場合、次は (2, 3, 4, 5) の中から 2 個を選べばいいですね。2 番目に 2 を選択したら、次は (3, 4, 5) の中から 1 個を選べばいいわけです。これで、(1, 2, 3), (1, 2, 4), (1, 2, 5) が生成されます。(2, 3, 4, 5) の中から 2 個選ぶとき、2 を選ばない場合があります。この場合は (3, 4, 5) の中から 2 個を選べばいいわけです。ここで 3 を選ぶと (1, 3, 4), (1, 3, 5) が生成できます。同様に、3 を除いた (4, 5) の中から 2 個を選ぶと (1, 4, 5) を生成することができます。
これで 1 を含む組み合わせを生成したので、次は 1 を含まない組み合わせ、つまり (2, 3, 4, 5) から 3 個を選ぶ組み合わせを生成すればいいわけです。けっきょく、この処理の考え方は次に示す組み合わせの公式と同じです。
Ruby でプログラムを作ると次のようになります。
リスト : 組み合わせの生成 def combination(n, xs, ys = [], &func) if n == 0 func.call ys elsif xs.size == n func.call (ys + xs) else combination(n - 1, xs[1..-1], ys + [xs[0]], &func) combination(n, xs[1..-1], ys, &func) end end
combination() は引数 xs の配列から n 個を選ぶ組み合わせを生成し、それをブロック func に渡して実行します。選んだ要素は第 3 引数 ys に格納します。n が 0 になったら組み合わせを一つ生成できたので、func.call(ys) を実行します。
次の節で、配列 xs の長さが n と等しくなったならば、xs の要素を全て選択します。ys + xs で配列を連結してブロック func を呼び出します。あとは combination() を再帰呼び出しするだけです。最初の呼び出しは先頭の要素を選択する場合です。先頭要素 xs[0] を ys に追加して、xs[1..-1] の中から n - 1 個を選びます。最後の呼び出しが先頭の要素を選ばない場合です。xs[1..-1] の中から n 個を選びます。
irb> combination(3, [1,2,3,4,5]){|xs| print xs, "\n"} [1, 2, 3] [1, 2, 4] [1, 2, 5] [1, 3, 4] [1, 3, 5] [1, 4, 5] [2, 3, 4] [2, 3, 5] [2, 4, 5] [3, 4, 5] => nil
なお、Ruby には組み合わせを生成するメソッド combination() が用意されています。
irb> [1,2,3,4,5].combination(3) {|xs| print xs, "\n"} [1, 2, 3] [1, 2, 4] [1, 2, 5] [1, 3, 4] [1, 3, 5] [1, 4, 5] [2, 3, 4] [2, 3, 5] [2, 4, 5] [3, 4, 5] => [1, 2, 3, 4, 5] irb> c = [1,2,3,4,5].combination(3) => #<Enumerator: ...> irb> c.next => [1, 2, 3] irb> c.next => [1, 2, 4] irb> c.next => [1, 2, 5]
メソッド combination() はブロックを省略すると、Enumerator のオブジェクトを返します。Ruby の高階関数はブロックを省略すると Enumerator のオブジェクトを返すものが多くあります。
配列 xs の中から n 個の要素を選択する順列を生成します。プログラムは次のようになります。
リスト : 順列の生成 def permutation(n, xs, ys = [], &func) if n == 0 func.call ys else xs.each {|x| next if ys.member? x permutation(n - 1, xs, ys + [x], &func) } end end
irb> permutation(4, [1,2,3,4]) {|xs| print xs, "\n"} [1, 2, 3, 4] [1, 2, 4, 3] [1, 3, 2, 4] [1, 3, 4, 2] [1, 4, 2, 3] [1, 4, 3, 2] [2, 1, 3, 4] [2, 1, 4, 3] [2, 3, 1, 4] [2, 3, 4, 1] [2, 4, 1, 3] [2, 4, 3, 1] [3, 1, 2, 4] [3, 1, 4, 2] [3, 2, 1, 4] [3, 2, 4, 1] [3, 4, 1, 2] [3, 4, 2, 1] [4, 1, 2, 3] [4, 1, 3, 2] [4, 2, 1, 3] [4, 2, 3, 1] [4, 3, 1, 2] [4, 3, 2, 1] => [1, 2, 3, 4]
なお、Ruby には組み合わせを生成するメソッド permutation() が用意されています。
irb> [1,2,3,4].permutation(3) {|xs| print xs, "\n"} [1, 2, 3] [1, 2, 4] [1, 3, 2] [1, 3, 4] [1, 4, 2] [1, 4, 3] [2, 1, 3] [2, 1, 4] [2, 3, 1] [2, 3, 4] [2, 4, 1] [2, 4, 3] [3, 1, 2] [3, 1, 4] [3, 2, 1] [3, 2, 4] [3, 4, 1] [3, 4, 2] [4, 1, 2] [4, 1, 3] [4, 2, 1] [4, 2, 3] [4, 3, 1] [4, 3, 2] => [1, 2, 3, 4] irb> ps = [1,2,3,4].permutation(3) => #<Enumerator: ...> irb> ps.next => [1, 2, 3] irb> ps.next => [1, 2, 4] irb> ps.next => [1, 3, 2] irb> ps.next => [1, 3, 4]
N 個の整数 1, 2, ..., N の順列を考えます。このとき、i 番目 (先頭が 1 番目) の要素が整数 i ではない順列を完全順列 (derangement) といいます。
リスト : 完全順列 def derangement(n, xs = [], &func) if xs.size == n func.call xs else (1..n).each {|x| next if (xs.size + 1 == x) || (xs.member? x) derangement(n, xs + [x], &func) } end end
1 から n までの数字を n 個選ぶ順列を生成する処理で、数字 x が xs.size + 1 と等しい場合は数字 x を選択しません。xs.size が n と等しい場合は n 個の数字を選んだので func.call(xs) を実行します。これで完全順列を生成することができます。
irb> derangement(3) {|xs| print xs, "\n"} [2, 3, 1] [3, 1, 2] => 1..3 irb> derangement(4) {|xs| print xs, "\n"} [2, 1, 4, 3] [2, 3, 4, 1] [2, 4, 1, 3] [3, 1, 4, 2] [3, 4, 1, 2] [3, 4, 2, 1] [4, 1, 2, 3] [4, 3, 1, 2] [4, 3, 2, 1] => 1..4
整数 n 以下の素数を求めるいちばん簡単な方法は、奇数 3, 5, 7, 9, ... をそれまでに見つけた素数で割ってみることです。この方法は Ruby 入門講座でプログラムを作ったことがあります。ここでは別の方法を紹介しましょう。
最初に、2 から n までの整数列を生成します。先頭の 2 は素数なので、この整数列から 2 で割り切れる整数を取り除き除きます。2 で割り切れる整数が取り除かれたので、残った要素の先頭が素数になります。先頭要素は 3 になるので、今度は 3 で割り切れる整数を取り除けばいいのです。このように、素数を見つけたらそれで割り切れる整数を取り除いていくアルゴリズムを「エラトステネスの篩 (ふるい)」といいます。
プログラムは次のようになります。
リスト : エラトステネスの篩 def sieve(n) primes = [2] table = Array.new(n / 2) { true } x = 3 while x * x <= n y = (x - 3) / 2 if table[y] primes.push x (y + x).step(table.size, x) {|i| table[i] = false } end x += 2 end while x <= n primes.push x if table[(x - 3) / 2] x += 2 end primes end
配列 primes に素数を格納し、table で奇数列 (3, 5, 7, ... ) を表します。true で素数を表し、素数でない場合は false に書き換えます。配列の生成は Array.new() で行うことができます。
Array.new(size) {|index| ... }
引数 size は配列の大きさです。ブロックの引数 index には配列の添字 (0 から size - 1 まで) が渡され、ブロックの返り値が配列の要素になります。ブロックを省略した場合は nil に初期化されます。
奇数を変数 x とし、それに対応する table の添字を変数 y とすると、変数 x は 3, 5, 7, 9, ... に、それに対応する変数 y は 0, 1, 2, 3, ... になります。この場合、x の倍数に対応する y の値は y + x, y + x * 2, y + x * 3, ... になります。たとえば、3, 5, 7 の倍数は次のようになります。
x : 3 5 7 9 11 13 15 17 19 21 23 25 y : 0 1 2 3 4 5 6 7 8 9 10 11 --+------------------------------------- 3 | O 0 O 0 5 | 0 0 0 7 | 0 0
プログラムは簡単です。最初の while ループで、x を √n まで +2 ずつ増やして素数かチェックします。table の添字 y は (x - 3) / 2 で求めることができます。table[y] が true ならば x は素数です。x の倍数を table から削除します。step() はブロックを繰り返し実行するメソッドです。
num.step(limit, step) {|n| ... }
step() は num からはじめ step を足しながら limit を越える前までブロックを繰り返します。ブロックの引数 n には数値が渡されます。times() と違って、step() は整数だけではなく浮動小数点でも動作します。
最後に、while ループで √n よりも大きい素数を求めます。
それでは実行してみましょう。
irb> sieve(100) => [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97] irb> sieve(500) => [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499]
ところで、エラトステネスの篩は次のようにプログラムすることもできます。
リスト : エラトステネスの篩 (2) def sieve1(n) ps = [2] xs = Array.new((n + 1) / 2 - 1) {|x| x * 2 + 3} while true p = xs[0] break if p * p > n ps.push p xs.delete_if {|x| x % p == 0} end ps + xs end
xs に奇数を格納した配列をセットします。あとは while 文の中で、xs から素数で割り切れるものを delete_if() で削除していくだけです。ただし、配列を破壊的に修正するので、sieve() よりも時間がかかると思います。興味のある方は試してみてください。
なお、Rubyには素数を扱うライブラリ prime が用意されています。
irb> require 'prime' => true irb> Prime.each(100) {|n| print n, " "} 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 => nil irb> Prime.first(20) => [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71]
詳しい説明は Ruby のリファレンス class Prime をお読みください。
差が 2 である素数の組を「双子素数 (twin prime)」といいます。あらかじめ素数表が用意されていれば、Ruby のメソッド each_cons() を使うと双子素数は簡単に求めることができます。
obj.each_cons(n) {|xs| ... } obj.each_slice(n) {|xs| ... }
どちらのメソッドも obj から n 個ずつ要素を取り出し、それらを配列に格納してブロックの引数 xs に渡します。たとえば、[1, 2, 3, 4] から要素を 2 つ取り出す場合、each_slice() は [1, 2], [3, 4] のように先頭から 2 つずつ順番に取り出しますが、each_cons() は [1, 2], [2, 3], [3, 4] のように、隣り合った要素を順番に取り出します。
簡単な使用例を示しましょう。
irb> (1..10).each_cons(2) {|a| print a} [1, 2][2, 3][3, 4][4, 5][5, 6][6, 7][7, 8][8, 9][9, 10]=> nil irb> (1..10).each_cons(3) {|a| print a} [1, 2, 3][2, 3, 4][3, 4, 5][4, 5, 6][5, 6, 7][6, 7, 8][7, 8, 9][8, 9, 10]=> nil irb> (1..10).each_slice(2) {|a| print a} [1, 2][3, 4][5, 6][7, 8][9, 10]=> nil irb> (1..10).each_slice(3) {|a| print a} [1, 2, 3][4, 5, 6][7, 8, 9][10]=> nil
双子素数は次のようになります。
irb> sieve(1000).each_cons(2) {|xs| print xs if xs[1] - xs[0] == 2} [3, 5][5, 7][11, 13][17, 19][29, 31][41, 43][59, 61][71, 73][101, 103][107, 109] [137, 139][149, 151][179, 181][191, 193][197, 199][227, 229][239, 241][269, 271] [281, 283][311, 313][347, 349][419, 421][431, 433][461, 463][521, 523][569, 571] [599, 601][617, 619][641, 643][659, 661][809, 811][821, 823][827, 829][857, 859] [881, 883]=> nil
もちろん、Ruby の Prime を使っても簡単に求めることができます。
irb> Prime.each(1000).each_cons(2) {|xs| print xs if xs[1] - xs[0] == 2} ・・・省略・・・
素因数分解はエラトステネスの篩と同じ考え方で行うことができます。プログラムは次のようになります。
リスト : 素因数分解 def factor_sub(n, m) c = 0 while n % m == 0 c += 1 n /= m end return c, n end def factorization(n) buff = [] c, m = factor_sub(n, 2) buff.push [2, c] if c > 0 x = 3 while x * x <= m c, m = factor_sub(m, x) buff.push [x, c] if c > 0 x += 2 end buff.push [m, 1] if m > 1 buff end
最初に 2 で割り算します。割り算するときは、その数で割り切れるあいだは割り算を続けることに注意してください。たとえば、27 を素因数分解すると 3 * 3 * 3 になりますが、3 を一回だけしか割り算しないと、結果は 3 * 9 のように素数ではない数が含まれてしまいます。この処理を関数 factor_sub() で行っています。
あとは、factor_sub() の返り値をチェックして、割り算した回数 c が 0 よりも大きければ、素数と c を配列に格納して、それを buff に追加します。これを変数 m が x * x 以下のあいだ繰り返します。最後に m が 1 よりも大きければ、(m, 1) を buff に追加して、return で buff を返します。
それでは実行してみましょう。
irb> factorization(24) => [[2, 3], [3, 1]] irb> factorization(12345678) => [[2, 1], [3, 2], [47, 1], [14593, 1]] irb> factorization(123456789) => [[3, 2], [3607, 1], [3803, 1]] irb> factorization(1234567890) => [[2, 1], [3, 2], [5, 1], [3607, 1], [3803, 1]] irb> factorization(1111111111) => [[11, 1], [41, 1], [271, 1], [9091, 1]]
なお、これはとても単純なアルゴリズムなので、大きな整数の素因数分解には適していません。巨大な合成数の素因数分解はとても難しい問題です。興味のある方は素因数分解について調べてみてください。
また、Ruby の prime には素因数分解を行うメソッド prime_division() が用意されています。
irb> Prime.prime_division(24) => [[2, 3], [3, 1]] irb> Prime.prime_division(123456789) => [[3, 2], [3607, 1], [3803, 1]]
n の素因数分解ができると、約数の個数を求めるのは簡単です。\(n = p^a \times q^b \times r^c\) とすると、約数の個数は (a + 1) * (b + 1) * (c + 1) になります。たとえば、12 は \(2^2 \times 3^1\) になるので、約数の個数は 3 * 2 = 6 になります。実際、12 の約数は 1, 2, 3, 4, 6, 12 の 6 個です。
プログラムは reduce() を使うと簡単です。
リスト : 約数の個数 def divisor_num(n) factorization(n).reduce(1) {|m, xs| m *= xs[1] + 1} end
irb> divisor_num(24) => 8 irb> divisor_num(12345678) => 24 irb> divisor_num(123456789) => 12 irb> divisor_num(1234567890) => 48 irb> divisor_num(1111111111) => 16
n の素因数分解ができると、約数の合計値を求めるのは簡単です。素因数分解した結果が \(p^a\) の場合、その約数は \(1, p, p^2, \cdots, p^a\) となり、初項が 1 で公比が p の等比数列になります。約数の合計値を \(s(p, a)\) で表すことにすると、等比数列の和の公式から \(s(p, a)\) は以下の式になります。
たとえば、8 の素因数分解は \(2^3\) になり、素数の合計値は 8 + 4 + 2 + 1 = 15 になります。整数 n の素因数分解が \(p^a \times q^b \times r^c\) だとすると、\(s(p, a) \times s(q, b) \times s(r, c)\) は約数の合計値になります。実際に式を展開すると次のようになります。
このように、各項は約数と一対一に対応しているので、その和は約数の合計値になるわけです。たとえば、12 は \(2^2 \times 3\) に素因数分解できますが、その合計値は (4 + 2 + 1) * (3 + 1) = 28 となります。12 の約数は 1, 2, 3, 4, 6, 12 なので、その合計値は確かに 28 になります。
プログラムは次のようになります。。
リスト : 約数の和 def divisor_sum(n) factorization(n).reduce(1) {|m, xs| m * (xs[0] ** (xs[1] + 1) - 1) / (xs[0] - 1)} end
irb> divisor_sum(24) => 60 irb> divisor_sum(12345678) => 27319968 irb> divisor_sum(123456789) => 178422816 irb> divisor_sum(1234567890) => 3211610688 irb> divisor_sum(1111111111) => 1246404096
p が素数の場合、\(p^a\) の約数は次のように簡単に求めることができます。
n の素因数分解が \(p^a \times q^b\) だったとすると、その約数は次のようになります。
たとえば、12 の約数は \(2^2\) = (1, 2, 4) と 3 = (1, 3) から、(1, 2, 4) * 1 と (1, 2, 4) * 3 のすべての要素 (1, 2, 4, 3, 6, 12) になります。この処理はメソッド product() を使うと簡単です。product() は 2 つの配列から直積集合を求めます、簡単な例を示しましょう。
irb> [1,2,3].product([4,5,6]) => [[1, 4], [1, 5], [1, 6], [2, 4], [2, 5], [2, 6], [3, 4], [3, 5], [3, 6]] irb> [1,2,3].product([4,5,6]).map {|x, y| x * y} => [4, 5, 6, 8, 10, 12, 12, 15, 18]
プログラムは次のようになります。
リスト : 約数 def divisor_sub(p, q) (0..q).reduce([]) {|xs, i| xs + [p ** i]} end def divisor(n) factorization(n).reduce([1]) {|ys, xs| ys.product(divisor_sub(*xs)).map {|x, y| x * y} }.sort() end
関数 divisor_sub() は pq の約数を配列に格納して返します。引数 n を factorization() で素因数分解して、その結果を reduce() に渡します。引数 ys が累積変数、xs は p, q を格納した配列です。xs を divisor_sub() に渡して配列に変換し、product() で直積集合を生成し、それから map() で要素 x, y を掛け算します。最後に sort() で配列をソートします。
簡単な実行例を示します。
irb> divisor(24) => [1, 2, 3, 4, 6, 8, 12, 24] irb> divisor(12345678) => [1, 2, 3, 6, 9, 18, 47, 94, 141, 282, 423, 846, 14593, 29186, 43779, 87558, 131337, 262674, 685871, 1371742, 2057613, 4115226, 6172839, 12345678] irb> divisor(123456789) => [1, 3, 9, 3607, 3803, 10821, 11409, 32463, 34227, 13717421, 41152263, 123456789] irb> divisor(1234567890) => [1, 2, 3, 5, 6, 9, 10, 15, 18, 30, 45, 90, 3607, 3803, 7214, 7606, 10821, 11409, 18035, 19015, 21642, 22818, 32463, 34227, 36070, 38030, 54105, 57045, 64926, 68454, 108210, 114090, 162315, 171135, 324630, 342270, 13717421, 27434842, 41152263, 68587105, 82304526, 123456789, 137174210, 205761315, 246913578, 411522630, 617283945, 1234567890] irb> divisor(1111111111) => [1, 11, 41, 271, 451, 2981, 9091, 11111, 100001, 122221, 372731, 2463661, 4100041, 27100271, 101010101, 1111111111]
配列 xs のべき集合を求める関数 power_set(xs) を作ります。たとえば、配列 [1, 2, 3] のべき集合は [], [1], [2], [3], [1, 2], [1, 3], [2, 3], [1, 2, 3] になります。
プログラムは次のようになります。
リスト : べき集合 def power_set(xs) if xs == [] return [[]] else ys = power_set(xs[1..-1]) ys + ys.map {|zs| [xs[0]] + zs} end end
irb> power_set([1,2,3]) => [[], [3], [2], [2, 3], [1], [1, 3], [1, 2], [1, 2, 3]] irb> power_set([1,2,3,4]) => [[], [4], [3], [3, 4], [2], [2, 4], [2, 3], [2, 3, 4], [1], [1, 4], [1, 3], [1, 3, 4], [1, 2], [1, 2, 4], [1, 2, 3], [1, 2, 3, 4]]
xs が空の場合は [] を格納した配列を返します。そうでなければ power_set() を再帰呼び出しして xs[1..-1] のべき集合 ys を求め、ys に先頭要素 xs[0] を追加します。そして、その集合と ys を演算子 + で連結します。
ところで、べき集合は再帰呼び出しを使わなくても、繰り返しでも簡単にプログラムすることができます。配列 xs のべき集合の要素は 2 ** xs.size 個あります。これらの要素を 0 から 2 ** xs.size - 1 まで整数に対応させ、ビットが 1 の位置の要素を選んでいくとべき集合を求めることができます。
Ruby の場合、整数 n に角カッコをつけると、その位置にあるビットの値を求めることができます。
irb> 0xf0[0] => 0 irb> 0xf0[4] => 1 irb> 0xf0[7] => 1
プログラムは次のようになります。
リスト : べき集合 (2) def power_set1(xs, &func) n = 2 ** xs.size for x in 0 ... n ys = [] for i in 0 ... n ys.push xs[i] if x[i] == 1 end func.call(ys) end end
irb> power_set1([1,2,3]) {|xs| print xs, " "} [] [1] [2] [1, 2] [3] [1, 3] [2, 3] [1, 2, 3] => 0...8
1 3 5 B───D───F 0 /│ │ A │ │ \│ │ C───E───G 2 4 6 図 : 経路の探索
スタート (A) からゴール (F) までの経路を求めるプログラムです。アルゴリズムの詳しい説明は拙作のページ「Algorithms with Python: 集合、グラフ、経路の探索」をお読みくださいませ。
リスト : 経路の探索 # 隣接リスト $adjacent = [ [1, 2], # A [0, 2, 3], # B [0, 1, 4], # C [1, 4, 5], # D [2, 3, 6], # E [3], # F [4]] # G # 経路の表示 def print_path(path) path.each do |x| print x, " " end print "\n" end # 深さ優先探索 def depth_search(goal, path) if goal == path[-1] print_path(path) else for n in $adjacent[path[-1]] unless path.member?(n) path.push n depth_search(goal, path) path.pop end end end end # 幅優先探索 def breadth_search(start, goal) queue = [[start]] while queue.size > 0 path = queue.shift if goal == path[-1] print_path(path) else for n in $adjacent[path[-1]] unless path.member?(n) queue.push (path.dup.push n) end end end end end # 反復深化 def id_search(limit, goal, path) if path.size == limit print_path(path) if goal == path[-1] else for n in $adjacent[path[-1]] unless path.member?(n) path.push n id_search(limit, goal, path) path.pop end end end end def keiro print "-- 深さ優先探索 --\n" depth_search(5, [0]) print "-- 幅優先探索 --\n" breadth_search(0, 5) print "-- 反復深化 --\n" for x in 1 .. 7 print x, " moves\n" id_search(x, 5, [0]) end end
irb> keiro -- 深さ優先探索 -- 0 1 2 4 3 5 0 1 3 5 0 2 1 3 5 0 2 4 3 5 -- 幅優先探索 -- 0 1 3 5 0 2 1 3 5 0 2 4 3 5 0 1 2 4 3 5 -- 反復深化 -- 1 moves 2 moves 3 moves 4 moves 0 1 3 5 5 moves 0 2 1 3 5 0 2 4 3 5 6 moves 0 1 2 4 3 5 7 moves => 1..7
リスト : 簡単なプログラム # coding: utf-8 # 九九表 def ninetynine for x in 1 .. 9 for y in 1 .. 9 printf "%2d ", x * y end print "\n" end end def ninetynine1 (1..9).each {|x| (1..9).each {|y| printf "%2d ", x * y } print "\n" } end def ninetynine2 (1..9).map {|x| (1..9).map {|y| x * y}} end # 整数列の和 def sum_of(n, m, &func) (n..m).reduce(0) {|a, x| a + func.call(x)} end # 度数分布表 # 身長 Height = [ 148.7, 149.5, 133.7, 157.9, 154.2, 147.8, 154.6, 159.1, 148.2, 153.1, 138.2, 138.7, 143.5, 153.2, 150.2, 157.3, 145.1, 157.2, 152.3, 148.3, 152.0, 146.0, 151.5, 139.4, 158.8, 147.6, 144.0, 145.8, 155.4, 155.5, 153.6, 138.5, 147.1, 149.6, 160.9, 148.9, 157.5, 155.1, 138.9, 153.0, 153.9, 150.9, 144.4, 160.3, 153.4, 163.0, 150.9, 153.3, 146.6, 153.3, 152.3, 153.3, 142.8, 149.0, 149.4, 156.5, 141.7, 146.2, 151.0, 156.5, 150.8, 141.0, 149.0, 163.2, 144.1, 147.1, 167.9, 155.3, 142.9, 148.7, 164.8, 154.1, 150.4, 154.2, 161.4, 155.0, 146.8, 154.2, 152.7, 149.7, 151.5, 154.5, 156.8, 150.3, 143.2, 149.5, 145.6, 140.4, 136.5, 146.9, 158.9, 144.4, 148.1, 155.5, 152.4, 153.3, 142.3, 155.3, 153.1, 152.3 ] def make_frequency m = 8 freq = [0] * m # 度数分布表 cum = [0] * m # 累積度数表 low = 130.0 z = 5.0 # 度数分布表の作成 for x in Height for j in 0 ... m if x < low + z * (j + 1) freq[j] += 1 break end end end # 累積度数表の作成 cum[0] = freq[0] for i in 1 ... m cum[i] = cum[i - 1] + freq[i] end # 表示 for i in 0 ... m printf "%.1f - %.1f | ", low + z * i, low + z * (i + 1) printf "%3d %3d\n", freq[i], cum[i] end end # 総和、最大値、最小値、平均値 def sum(a) a.reduce {|m, x| m + x} end def max(a) a.reduce {|m, x| m < x ? x : m} end def min(a) a.reduce {|m, x| m > x ? x : m} end def average(a) sum(a).fdiv(a.size) end # 累乗 def power0(x, n) a = 1 n.times {|_| a *= x} a end def power1(x, n) (1..n).reduce(1) {|a, _| a * x} end def power2(x, n) return 1 if n == 0 v = power2(x, n / 2) v *= v v *= x if n % 2 == 1 v end # ハノイの塔 def hanoi(n, from, to, via) hanoi(n - 1, from, via, to) if n > 1 printf("disk %d : %c -> %c\n", n, from, to) hanoi(n - 1, via, to, from) if n > 1 end # フィボナッチ数列 def fibo(n) return n if n < 2 fibo(n - 1) + fibo(n - 2) end def fibo1(n, a = 0, b = 1) return a if n == 0 fibo1(n - 1, b, a + b) end Fibos = Enumerator.new {|y| a, b = 0, 1 loop { y << a a, b = b, a + b } } # 組み合わせの数 def comb(n, r) return 1 if n == r || r == 0 comb(n, r - 1) * (n - r + 1) / r end # パスカルの三角形 def pascal(x) for n in 0 .. x for r in 0 .. n print comb(n, r), " " end print "\n" end end def pascal1(x) a = [1] (x + 1).times {|_| print a, "\n" a.push 0 a = a.zip(a.reverse).map {|xs| xs[0] + xs[1]} } end # 組み合わせの生成 def combination(n, xs, ys = [], &func) if n == 0 func.call ys elsif xs.size == n func.call (ys + xs) else combination(n - 1, xs[1..-1], ys + [xs[0]], &func) combination(n, xs[1..-1], ys, &func) end end # 順列の生成 def permutation(n, xs, ys = [], &func) if n == 0 func.call ys else xs.each {|x| next if ys.member? x permutation(n - 1, xs, ys + [x], &func) } end end # 完全順列 def derangement(n, xs = [], &func) if xs.size == n func.call xs else (1..n).each {|x| next if (xs.size + 1 == x) || (xs.member? x) derangement(n, xs + [x], &func) } end end # エラトステネスの篩 def sieve(n) primes = [2] table = Array.new(n / 2) { true } x = 3 while x * x <= n y = (x - 3) / 2 if table[y] primes.push x (y + x).step(table.size, x) {|i| table[i] = false } end x += 2 end while x <= n primes.push x if table[(x - 3) / 2] x += 2 end primes end def sieve1(n) ps = [2] xs = Array.new((n + 1) / 2 - 1) {|x| x * 2 + 3} while true p = xs[0] break if p * p > n ps.push p xs.delete_if {|x| x % p == 0} end ps + xs end # 素因数分解 def factor_sub(n, m) c = 0 while n % m == 0 c += 1 n /= m end return c, n end def factorization(n) buff = [] c, m = factor_sub(n, 2) buff.push [2, c] if c > 0 x = 3 while x * x <= m c, m = factor_sub(m, x) buff.push [x, c] if c > 0 x += 2 end buff.push [m, 1] if m > 1 buff end # 約数の個数 def divisor_num(n) factorization(n).reduce(1) {|m, xs| m *= xs[1] + 1} end # 約数の和 def divisor_sum(n) factorization(n).reduce(1) {|m, xs| m * (0..xs[1]).reduce(0) {|a, x| a + xs[0] ** x}} end # 約数 def divisor_sub(p, q) (0..q).reduce([]) {|xs, i| xs + [p ** i]} end def divisor(n) factorization(n).reduce([1]) {|ys, xs| ys.product(divisor_sub(*xs)).map {|x, y| x * y} }.sort() end # べき集合 def power_set(xs) if xs == [] return [[]] else ys = power_set(xs[1..-1]) ys + ys.map {|zs| [xs[0]] + zs} end end def power_set1(xs, &func) n = 2 ** xs.size for x in 0 ... n ys = [] for i in 0 ... n ys.push xs[i] if x[i] == 1 end func.call(ys) end end # 経路の探索 # 隣接リスト $adjacent = [ [1, 2], # A [0, 2, 3], # B [0, 1, 4], # C [1, 4, 5], # D [2, 3, 6], # E [3], # F [4]] # G # 経路の表示 def print_path(path) path.each do |x| print x, " " end print "\n" end # 深さ優先探索 def depth_search(goal, path) if goal == path[-1] print_path(path) else for n in $adjacent[path[-1]] unless path.member?(n) path.push n depth_search(goal, path) path.pop end end end end # 幅優先探索 def breadth_search(start, goal) queue = [[start]] while queue.size > 0 path = queue.shift if goal == path[-1] print_path(path) else for n in $adjacent[path[-1]] unless path.member?(n) queue.push (path.dup.push n) end end end end end # 反復深化 def id_search(limit, goal, path) if path.size == limit print_path(path) if goal == path[-1] else for n in $adjacent[path[-1]] unless path.member?(n) path.push n id_search(limit, goal, path) path.pop end end end end def keiro print "-- 深さ優先探索 --\n" depth_search(5, [0]) print "-- 幅優先探索 --\n" breadth_search(0, 5) print "-- 反復深化 --\n" for x in 1 .. 7 print x, " moves\n" id_search(x, 5, [0]) end end