M.Hiroi's Home Page

Functional Programming

お気楽 Standard ML of New Jersey 入門

[ PrevPage | SML/NJ | NextPage ]

関数型電卓プログラムの改良 (付録D)

関数型電卓プログラム fcalc の使用例として、Yet Another SML/NJ Problems からいくつか問題を選んで fcalc で解いてみました。

●問題1

自然数 n を素因数分解する関数 factorization(n) を定義してください。返り値はリスト ((p . q) ...) で、(p . q) は pq を表します。

Calc> factorization(6);
((2 . 1) (3 . 1))
Calc> factorization(12345678);
((2 . 1) (3 . 2) (47 . 1) (14593 . 1))
Calc> factorization(123456789);
((3 . 2) (3607 . 1) (3803 . 1))
Calc> factorization(1234567890);
((2 . 1) (3 . 2) (5 . 1) (3607 . 1) (3803 . 1))
Calc> factorization(1111111111);
((11 . 1) (41 . 1) (271 . 1) (9091 . 1))

解答

●問題2

自然数 n の約数の個数を求める関数 divisor_num を定義してください。

Calc> divisor_num(6);
4
Calc> divisor_num(12345678);
24
Calc> divisor_num(123456789);
12
Calc> divisor_num(1234567890);
48
Calc> divisor_num(1111111111);
16

解答

●問題3

自然数 n の約数の合計値を求める関数 divisor_sum を定義してください。

Calc> divisor_sum(6);
12
Calc> divisor_sum(12345678);
27319968
Calc> divisor_sum(123456789);
178422816
Calc> divisor_sum(1234567890);
3211610688
Calc> divisor_sum(1111111111);
1246404096

解答

●問題4

自然数 n の約数をリストに格納して返す関数 divisor を定義してください。

Calc> divisor(6);
(1 2 3 6)
Calc> 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)
Calc> divisor(123456789);
(1 3 9 3607 10821 32463 3803 11409 34227 13717421 41152263 123456789)
Calc> divisor(1234567890);
(1 2 3 6 9 18 5 10 15 30 45 90 3607 7214 10821 21642 32463 64926 18035 36070 
54105 108210 162315 324630 3803 7606 11409 22818 34227 68454 19015 38030 57045 
114090 171135 342270 13717421 27434842 41152263 82304526 123456789 246913578 
68587105 137174210 205761315 411522630 617283945 1234567890)
Calc> divisor(1111111111);
(1 11 41 451 271 2981 11111 122221 9091 100001 372731 4100041 2463661 27100271
 101010101 1111111111)

解答

●問題5

完全数 - Wikipedia によると、『完全数(かんぜんすう,perfect number)とは、その数自身を除く約数の和が、その数自身と等しい自然数のことである。』 とのことです。自然数 n 以下の完全数を求める関数 perfect_number を定義してください。

Calc> perfect_number(10000);
6 28 496 8128

解答

●問題6

友愛数 - Wikipedia によると、『友愛数(ゆうあいすう)とは、異なる2つの自然数の組で、自分自身を除いた約数の和が、互いに他方と等しくなるような数をいう。』 とのことです。自然数 n 以下の友愛数を求める関数 yuuai_number を定義してください。

Calc> yuuai_number(100000);
(220 284)
(1184 1210)
(2620 2924)
(5020 5564)
(6232 6368)
(10744 10856)
(12285 14595)
(17296 18416)
(63020 76084)
(66928 66992)
(67095 71145)
(69615 87633)
(79750 88730)

解答

●問題7

リスト xs から重複を許して n 個の要素を選ぶ順列を求める関数 repeat_perm(f, n, xs) を定義してください。

Calc> repeat_perm(print, 3, list(1,2,3));
(1 1 1)(1 1 2)(1 1 3)(1 2 1)(1 2 2)(1 2 3)(1 3 1)(1 3 2)(1 3 3)(2 1 1)(2 1 2)
(2 1 3)(2 2 1)(2 2 2)(2 2 3)(2 3 1)(2 3 2)(2 3 3)(3 1 1)(3 1 2)(3 1 3)(3 2 1)
(3 2 2)(3 2 3)(3 3 1)(3 3 2)(3 3 3)

解答

●問題8

リスト xs から重複を許して n 個の要素を選ぶ組み合わせを求める関数 repeat_comb(f, n, xs) を定義してください。

Calc> repeat_comb(print, 3, list(1,2,3,4));
(1 1 1)(1 1 2)(1 1 3)(1 1 4)(1 2 2)(1 2 3)(1 2 4)(1 3 3)(1 3 4)(1 4 4)
(2 2 2)(2 2 3)(2 2 4)(2 3 3)(2 3 4)(2 4 4)(3 3 3)(3 3 4)(3 4 4)(4 4 4)

解答

●問題9

m 個の整数 1, 2, ..., m の順列を考えます。このとき、i 番目の要素が整数 i ではない順列を「完全順列」といいます。1 から m までの整数値で完全順列を生成する高階関数 derangement(f, m) を定義してください。

Calc> derangement(fn(x) print(x), print("\n") end, 3);
(2 3 1)
(3 1 2)

Calc> derangement(fn(x) print(x), print("\n") end, 4);
(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)

解答

●問題10

完全順列の総数を「モンモール数 (Montmort number)」といいます。モンモール数は次の漸化式で求めることができます。

\(\begin{array}{l} A_1 = 0 \\ A_2 = 1 \\ A_n = (n - 1) \times (A_{n-1} + A_{n-2}) \quad (n \geq 3) \end{array}\)

モンモール数を求める関数 montmort_number を定義してください。

Calc> montmort_number(1);
0
Calc> montmort_number(2);
1
Calc> montmort_number(3);
2
Calc> montmort_number(4);
9
Calc> montmort_number(5);
44
Calc> montmort_number(10);
1334961
Calc> montmort_number(20);
895014631192902121
Calc> montmort_number(30);
97581073836835777732377428235481

解答

●問題11

バランスの取れた n 対のカッコ列を生成する高階関数 kakko(f, n) を定義してください。カッコ列は ( と ) からなる列のことで、バランスが取れているカッコ列は、右カッコで閉じることができる、つまり右カッコに対応する左カッコがある状態のことをいいます。たとえば n = 1 の場合、( ) はバランスの取れたカッコ列ですが、) ( はバランスが取れていません。

Calc> kakko(print_kakko, 3);
((()))
(()())
(())()
()(())
()()()
0
Calc> kakko(print_kakko, 4);
(((())))
((()()))
((())())
((()))()
(()(()))
(()()())
(()())()
(())(())
(())()()
()((()))
()(()())
()(())()
()()(())
()()()()
0

解答

●問題12

バランスの取れた n 対のカッコ列の総数を求める関数 kakko_num n を定義してください。

Calc> kakko_num(1);
1
Calc> kakko_num(2);
2
Calc> kakko_num(3);
5
Calc> kakko_num(4);
14
Calc> kakko_num(5);
42
Calc> kakko_num(10);
16796
Calc> kakko_num(50);
1978261657756160653623774456
Calc> kakko_num(100);
896519947090131496687170070074100632420837521538745909320

解答

●問題13

整数 n を 1 以上の自然数の和で表すことを考えます。これを「整数の分割」といいます。整数を分割するとき、同じ自然数を何回使ってもかまいませんが、並べる順序が違うだけのものは同じ分割とします。簡単な例を示します。

n = 6
6 分割 : 1 + 1 + 1 + 1 + 1 + 1
5 分割 : 2 + 1 + 1 + 1 + 1
4 分割 : 3 + 1 + 1 + 1
         2 + 2 + 1 + 1
3 分割 : 4 + 1 + 1
         3 + 2 + 1
         2 + 2 + 2
2 分割 : 5 + 1
         4 + 2
         3 + 3
1 分割 : 6

6 の場合、分割の仕方は 11 通りあります。この数を「分割数」といいます。整数 n の分割の仕方をすべて求める高階関数 partition_of_integer(f, n) を定義してください。

Calc> partition_of_integer(fn(x) print(x), print("\n") end, 5);
(5)
(4 1)
(3 2)
(3 1 1)
(2 2 1)
(2 1 1 1)
(1 1 1 1 1)

Calc> partition_of_integer(fn(x) print(x), print("\n") end, 6);
(6)
(5 1)
(4 2)
(4 1 1)
(3 3)
(3 2 1)
(3 1 1 1)
(2 2 2)
(2 2 1 1)
(2 1 1 1 1)
(1 1 1 1 1 1)

解答

●問題14

自然数 n の分割数を求める関数 partition_number を定義してください。

Calc> partition_number(1);
1
Calc> partition_number(2);
2
Calc> partition_number(3);
3
Calc> partition_number(4);
5
Calc> partition_number(5);
7
Calc> partition_number(10);
42
Calc> partition_number(20);
627
Calc> partition_number(50);
204226

解答14

●問題15

リストで表した集合 ls を分割することを考えます。たとえば、集合 (1 2 3) は次のように分割することができます。

1 分割 : ((1 2 3))
2 分割 : ((1 2) (3)), ((1 3) (2)), ((1) (2 3))
3 分割 ; ((1) (2) (3))

このように、分割した集合 xs は元の集合 ls の部分集合になります。分割した部分集合の積は空集合になり、分割した部分集合のすべての和を求めると元の集合になります。

ls の分割の仕方をすべて求める高階関数 parititon_of_set(f, xs) を定義してください。

Calc> partition_of_set(fn(x) print(x), print("\n") end, list(1,2,3));
((1 2 3))
((1) (2 3))
((1 2) (3))
((1 3) (2))
((1) (2) (3))

Calc> partition_of_set(fn(x) print(x), print("\n") end, list(1,2,3,4));
((1 2 3 4))
((1) (2 3 4))
((1 2) (3 4))
((1 3 4) (2))
((1) (2) (3 4))
((1 2 3) (4))
((1 4) (2 3))
((1) (2 3) (4))
((1 2 4) (3))
((1 3) (2 4))
((1) (2 4) (3))
((1 2) (3) (4))
((1 3) (2) (4))
((1 4) (2) (3))
((1) (2) (3) (4))

解答

●問題16

集合を分割する方法の総数を「ベル数 (Bell Number)」といい、次の漸化式で求めることができます。

\(\begin{array}{l} B(0) = B(1) = 1 \\ B(n+1) = \displaystyle \sum_{k=0}^n {}_n \mathrm{C}_k B(k) \quad (n \geq 1) \end{array}\)

ベル数を求める関数 bell_number(n) を定義してください。

Calc> bell_number(1);
1
Calc> bell_number(2);
2
Calc> bell_number(3);
5
Calc> bell_number(4);
15
Calc> bell_number(5);
52
Calc> bell_number(10);
115975
Calc> bell_number(20);
51724158235372
Calc> bell_number(30);
846749014511809332450147
Calc> bell_number(40);
157450588391204931289324344702531067
Calc> bell_number(50);
185724268771078270438257767181908917499221852770

解答

●問題17

k 個の要素をもつ集合 ls を要素数が等しい m 個の部分集合に分割することを考えます。部分集合の要素数 n は k / m になります。分割の仕方をすべて求める高階関数 group_partition(f, n, m, ls) を定義してください。

Calc> group_partition(fn(x) print(x), print("\n") end, 2, 2, iota(1, 4));
((1 2) (3 4))
((1 4) (2 3))
((1 3) (2 4))
0
Calc> group_partition(fn(x) print(x), print("\n") end, 2, 3, iota(1, 6));
((1 2) (3 4) (5 6))
((1 4) (2 3) (5 6))
((1 3) (2 4) (5 6))
((1 2) (3 6) (4 5))
((1 6) (2 3) (4 5))
((1 3) (2 6) (4 5))
((1 2) (3 5) (4 6))
((1 5) (2 3) (4 6))
((1 3) (2 5) (4 6))
((1 6) (2 5) (3 4))
((1 5) (2 6) (3 4))
((1 6) (2 4) (3 5))
((1 4) (2 6) (3 5))
((1 5) (2 4) (3 6))
((1 4) (2 5) (3 6))
0

解答

●問題18

集合を group_partition で分割するとき、その仕方の総数を求める関数 group_partition_number(n, m) を定義してください。引数 n は部分集合の要素数、m は部分集合の個数です。

Calc> group_partition_number(2, 2);
3
Calc> group_partition_number(2, 3);
15
Calc> group_partition_number(3, 3);
280
Calc> group_partition_number(3, 5);
1401400

解答

●問題19

整数 n 未満の素数をすべて求める関数 sieve(f, n) を作ってください。

Calc> sieve(fn(x) print(x), print(" ") end, 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
Calc> sieve(fn(x) print(x), print(" ") end, 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

解答

●問題20

「ラテン方陣」は数独の枠の条件を無くした方陣です。ラテン方陣の定義を 参考文献 [1] より引用します。

『ラテン方陣を一般的にいうなら、n 行 n 列の正方形の枡に n 種類の記号を n 個ずつ配列して、各行各列に記号の重複のないものを n 次のラテン方陣というのです。』

このラテン方陣をパズルに応用したものが数独というわけです。

簡単な例を示しましょう。3 次のラテン方陣は次に示す 12 通りになります。

 1 2 3    1 2 3    1 3 2    1 3 2    2 1 3    2 1 3 
 2 3 1    3 1 2    2 1 3    3 2 1    1 3 2    3 2 1 
 3 1 2    2 3 1    3 2 1    2 1 3    3 2 1    1 3 2 
 標準形

 2 3 1    2 3 1    3 1 2    3 1 2    3 2 1    3 2 1 
 1 2 3    3 1 2    1 2 3    2 3 1    1 3 2    2 1 3 
 3 1 2    1 2 3    2 3 1    1 2 3    2 1 3    1 3 2 


               図 : 3 次のラテン方陣

この中で、最初の行と列の要素を昇順に並べたものを「標準形」といいます。3 次のラテン方陣の場合、標準形は 1 種類しかありません。ラテン方陣は任意の行を交換する、または任意の列を交換してもラテン方陣になります。3 次のラテン方陣の場合、標準形から行または列を交換することで、残りの 11 種類のラテン方陣を生成することができます。

4 次の標準形ラテン方陣をすべて求めてください。

解答

-- 参考文献 --------
[1] 大村平, 『数理パズルのはなし』, 日科技連出版社, 1998

●解答1

リスト : 素因数分解

def factorization(n)
  let rec
    a, iter = nil,
              fn(m, c)
                if n % m != 0 then
                  if c > 0 then a = cons(cons(m, c), a) end
                else
                  n = n / m,
                  iter(m, c + 1)
                end
              end
  in
    # 2 で割る
    iter(2, 0),
    # 奇数で割る
    let i = 3 in
      while n > 1 and n >= i * i do
        iter(i, 0),
        i = i + 2
      end
    end,
    nreverse(if n > 1 then cons(cons(n, 1), a) else a end)
  end
end

素因数分解は素数 2, 3, 5, ... で順番に割り算していけばいいのですが、いちいち素数を求めるのは大変なので、2 と 3 以上の奇数列で割り算していきます。局所関数 iter は n を m で割り算します。このとき、m で割り切れる回数を求めます。iter は m と割り算した回数をコンスセルに格納し、それを累積変数 a のリストに追加します。

最初に、iter を呼び出して n を 2 で割り算します。それから、while ループで奇数列を生成します。変数 i は 3 で初期化します。n が 1 になる、または \(\sqrt n \lt i\) になったら繰り返しを終了します。奇数列には素数ではないものがありますが、その前に小さな素数で素因数分解されているので、n がその値で割り切れることはありません。

●解答2

n の素因数分解ができると、約数の個数を求めるのは簡単です。\(n = p^a \times q^b \times r^c\) とすると、約数の個数は \((a + 1) \times (b + 1) \times (c + 1)\) になります。たとえば、12 は \(2^2 \times 3^1\) になるので、約数の個数は 3 * 2 = 6 になります。実際、12 の約数は 1, 2, 3, 4, 6, 12 の 6 個です。

プログラムは次のようになります。

リスト : 約数の個数

def divisor_num(n)
  foldl(fn(x, a) a * (cdr(x) + 1) end, 1, factorization(n))
end

divisor_num は foldl を使って 1 + cdr(x) を a に掛け算していくだけです。

●解答3

n の素因数分解ができると、約数の合計値を求めるのは簡単です。n の素因数分解が \(p^a\) だった場合、その約数の合計値は次の式で求めることができます。

\( \sigma(p, a) = p^a + p^{a-1} + \cdots + p^2 + p + 1 \)

たとえば、8 の素因数分解は \(2^3\) になり、素数の合計値は 8 + 4 + 2 + 1 = 15 になります。

\(p^a\) の約数の合計値を \(\sigma(p, a)\) で表すことにします。\(n = p^a \times q^b \times r^c\) の場合、n の約数の合計値は \(\sigma(p, a) \times \sigma(q, b) \times \sigma(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)
  let rec
    div_sum = fn(p, n, a)
      if n == 0 then
        a + 1
      else
        div_sum(p, n - 1, a + expt(p, n))
      end
    end
  in
    foldl(fn(x, a) a * div_sum(car(x), cdr(x), 0) end, 1, factorization(n))
  end
end

局所関数 div_sum は \(\sigma(p, n)\) を計算します。あとは foldl で div_sum の返り値を累積変数 a に掛け算していくだけです。

●解答4

p が素数の場合、\(p^a\) の約数は次のように簡単に求めることができます。

\( p^{a}, \ p^{a-1}, \ \cdots, \ p^{2}, \ p, \ 1 \)

n の素因数分解が \(p^a \times q^b\) だったとすると、その約数は次のようになります。

\(\begin{array}{l} (p^{a}, \ p^{a-1}, \ \cdots, \ p^{2}, \ p, \ 1) \times q^b, \\ (p^{a}, \ p^{a-1}, \ \cdots, \ p^{2}, \ p, \ 1) \times q^{b-1}, \\ \quad \quad \cdots \cdots \\ (p^{a}, \ p^{a-1}, \ \cdots, \ p^{2}, \ p, \ 1) \times q^2, \\ (p^{a}, \ p^{a-1}, \ \cdots, \ p^{2}, \ p, \ 1) \times q^1, \\ (p^{a}, \ p^{a-1}, \ \cdots, \ p^{2}, \ p, \ 1) \times 1 \end{array}\)

たとえば、12 の約数は \(2^4\) = (1, 2, 4) と 3 = (1, 3) から、(1, 2, 4) * 1 と (1, 2, 4) * 3 のすべての要素 (1, 2, 4, 3, 6, 12) になります。

プログラムは次のようになります。

リスト : 約数をすべて求める

def divisor(n)
  let rec
    div_sub = fn(p, n, a)
      if n == 0 then
        cons(1, a)
      else
        div_sub(p, n - 1, cons(expt(p, n), a))
      end
    end
  in
    let xs = factorization(n) in
      foldl(fn(x, a)
              map(fn(x) car(x) * cdr(x) end,
                  product(div_sub(car(x), cdr(x), nil), a))
            end,
            div_sub(caar(xs), cdar(xs), nil),
            cdr(xs))
    end
  end
end

局所関数 div_sub は pn の約数をリストに格納して返します。最初に、foldl で素因数分解した結果を順番に取り出し、(p . n) を div_sub でリストに変換します。あとは、その結果と累積変数 a のリストから product で直積集合を生成し、map で 2 つの要素を掛け算していくだけです。

●解答5

リスト : 完全数

def perfect_number(n)
  let x = 2 in
    while x <= n do
      if divisor_sum(x) - x == x then
        print(x), print(" ")
      end,
      x = x + 1
    end,
    nil
  end
end

完全数を求める perfect_number は簡単です。x の約数の合計値を divisor_sum で求め、その値から x を引いた値が x と等しければ完全数です。print で x を表示します。

●解答6

リスト : 友愛数

def yuuai_number(n)
  let x = 2 in
    while x <= n do
      let m = divisor_sum(x) - x in
        if x < m and x == divisor_sum(m) - m then
          print(list(x, m)), print("\n")
        end,
        x = x + 1
      end
    end,
    nil
  end
end

友愛数を求める yuuai_number も簡単です。divisor_sum で x の約数の合計値を求め、その値から x を引いた値を変数 m にセットします。m の約数の合計値から m を引いた値が x と等しければ、x と m は友愛数です。print で x と m を表示します。同じ組を表示しないようにするため、x < m を条件に入れています。

●解答7

リスト : 重複順列

def repeat_perm(f, n, xs)
  let rec
    iter = fn(n, xs, a)
      if n == 0 then
        f(reverse(a))
      else
        foreach(fn(x) iter(n - 1, xs, cons(x, a)) end, xs)
      end
    end
  in
    iter(n, xs, nil)
  end
end

重複順列は簡単にプログラムできます。選んだ要素を取り除く必要がないので、局所変数 iter を再帰呼び出しするとき、リスト xs をそのまま渡すだけです。

●解答8

リスト : 重複組み合わせ

def repeat_comb(f, n, xs)
  let rec
    iter = fn(n, xs, a)
      if n == 0 then
        f(reverse(a))
      else
        if single(xs) then
          f(append(reverse(a), makelist(n, car(xs))))
        else
          iter(n - 1, xs, cons(car(xs), a)),
          iter(n, cdr(xs), a)
        end
      end
    end
  in
    if n > length(xs) then nil else iter(n, xs, nil) end
  end
end

重複組み合わせを求める repeat_comb も簡単です。局所変数 iter で、リスト xs に要素が一つしかない場合は、その要素を n 個選びます。makelist で car(xs) を n 個格納したリストを生成します。最後の節では、先頭の要素を選んだあと、それを取り除かないで xs から n - 1 個の要素を選びます。

●解答9

リスト : 完全順列

def derangement(f, m)
  let rec
    perm = fn(n, xs, a)
      if null(xs) then
        f(reverse(a))
      else
        foreach(fn(x)
                  if n != x then
                    perm(n + 1, remove(fn(y) y == x end, xs), cons(x, a))
                  end
                end,
                xs)
      end
    end
  in
    perm(1, iota(1, m), nil)
  end
end

derangement は簡単です。実際の処理は局所関数 perm で行います。iota で 1 から m までの数値を格納したリストを生成し、それを引数 xs に渡します。引数 n が順番を表します。foreach の匿名関数の中で、数字 x が n と等しくない場合、その数字を選択することできます。等しい場合は選択しません。xs が空リストになったら、reverse で a を反転して関数 f を評価します。これで完全順列を生成することができます。

●解答10

リスト : 完全順列の総数

def montmort_number(n)
  if n == 1 then
    0
  else
    if n == 2 then
      1
    else
      (n - 1) * (montmort_number(n - 1) + montmort_number(n - 2))
    end
  end
end

# 別解
def montmort_number1(n)
  let rec
    iter = fn(i, a, b)
      if i == n then
        a
      else
        iter(i + 1, b, (i + 1) * (a + b))
      end
    end
  in
    iter(1, 0, 1)
  end
end

関数 montmort_number は公式をそのままプログラムしただけです。二重再帰になっているので、実行速度はとても遅くなります。これを繰り返しに変換すると別解のようになります。考え方はフィボナッチ数列と同じです。累積変数 a に i 番目の値を、b に i + 1 番目の値を保存しておきます。すると、i + 2 番目の値は (i + 1) * (a + b) で計算することができます。あとは、b の値を a に、新しい値を b にセットして処理を繰り返すだけです。

●解答11

リスト : カッコ列の生成

def print_kakko(xs)
  foreach(fn(x) print(x) end, xs),
  print("\n")
end

def kakko(f, m)
  let rec
    kakko_sub = fn(x, y, a)
      if x == m and y == m then
        f(reverse(a))
      else
        if x < m then kakko_sub(x + 1, y, cons("(", a)) end,
        if y < x then kakko_sub(x, y + 1, cons(")", a)) end
      end
    end
  in
    kakko_sub(0, 0, nil)
  end
end

カッコ列の生成は簡単です。局所関数 kakko_sub の引数 x が左カッコの個数、引数 y が右カッコの個数を表します。引数 a は累積変数で、文字 "(", ")" を格納したリストです。

バランスの取れたカッコ列の場合、x, y, m には y <= x <= m の関係が成り立ちます。x == m and y == m の場合、カッコ列がひとつ完成しました。リスト a を反転して、引数の関数 f を呼び出します。そうでなければ、kakko_sub を再帰呼び出しします。x < m であれば左カッコを追加し、y < x であれば右カッコを追加します。これでカッコ列を生成することができます。

●解答12

カタラン数 - Wikipedia によると、カッコ列の総数は「カタラン数 (Catalan number)」になるとのことです。カタラン数は次に示す公式で求めることができます。

\( \mathrm{C}_n = \dfrac{(2n)!}{(n+1)!n!} \)

これをそのままプログラムしてもいいのですが、それではちょっと面白くないので別な方法でプログラムを作ってみましょう。カタラン数は次に示す経路図において、A から B までの最短距離の道順を求めるとき、対角線を超えないものの総数に一致します。


              図 : 道順の総数の求め方

A からある地点にいたる最短距離の道順の総数は、左隣と真下の地点の値を足したものになります。一番下の地点は 1 で、対角線を越えた地点は 0 になります。あとは下から順番に足し算していけば、A から B までの道順の総数を求めることができます。上図の場合はカラタン数 C4 に相当し、その値は 14 となります。

プログラムはベクタを使うと簡単です。次の図を見てください。

0 : [1, 1, 1, 1, 1]

1 : [1, 1, 1, 1, 1]

2 : [1, 1, 1+1=2, 2+1=3, 3+1=4]
 => [1, 1, 2, 3, 4]

3 : [1, 1, 2, 3+2=5, 5+4=9]
 => [1, 1, 2, 5, 9]

4 : [1, 1, 2, 5, 5+9=14]
 => [1, 1, 2, 5, 14]

上図は Cn (n = 4) を求める場合です。大きさが n + 1, 要素の値が 1 のベクタを用意します。n = 0, 1 の場合は n 番目の要素をそのまま返します。n が 2 よりも大きい場合、変数 i を 2 に初期化して、i - 1 番目以降の要素の累積和を求めます。

たとえば i = 2 の場合、2 番目の要素は 1 番目の要素と自分自身を加算した値 2 になります。3 番目の要素は 2 番目の要素と自分自身を足した値 3 になり、4 番目の要素は 3 + 1 = 4 になります。次に i を +1 して同じことを繰り返します。3 番目の要素は 2 + 3 = 5 になり、4 番目の要素は 5 + 4 = 9 になります。i = 4 のとき、4 番目の要素は 5 + 9 = 14 となり、C4 の値を求めることができました。

プログラムは次のようになります。

リスト : カッコ列の総数

def kakko_num(n)
  let
    tbl = makeVector(n + 1, 1),
    i = 2
  in
    while i <= n do
      let j = i in
        while j <= n do
          tbl[j] = tbl[j] + tbl[j - 1],
          j = j + 1
        end
      end,
      i = i + 1
    end,
    tbl[n]
  end
end

説明したことをそのままプログラムしただけなので、とくに難しいところはないと思います。

●解答13

リスト : 整数の分割

def partition_of_integer(f, n)
  let rec
    part_int = fn(n, k, a)
      if n <= 1 then
        if n == 0 then f(reverse(a)) else f(reverse(cons(1, a))) end
      else
        if k == 1 then
          f(reverse(append(makelist(n, 1), a)))
        else
          if n - k >= 0 then
            part_int(n - k, k, cons(k, a))
          end,
          part_int(n, k - 1, a)
        end
      end
    end
  in
    part_int(n, n, nil)
  end
end

整数の分割は局所関数 part_int(n, k, a) で行います。part_int は整数 n を k 以下の整数で分割する関数です。a は結果を格納する累積変数です。n が 0 の場合、分割は終了したので関数 f を呼び出します。n が 1 の場合、分割の方法は 1 しかありません。a に 1 を追加してから関数 f を呼び出します。k が 1 の場合は n を 1 で分割するだけです。makelist で n 個の 1 を格納したリストを生成し、それを累積変数 a に連結してから関数 f を呼び出します。

それ以外の場合は、n から k を選ぶ方法と n から k - 1 を選ぶ方法に分割します。n - k が 0 以上の場合、n から k を選ぶことができます。あとは n - k から k を選ぶ方法を試すため、part_int を再帰呼び出しします。k を選ばずに k - 1 を選ぶ場合も、part_int を再帰呼び出しするだけです。最後に part_int(n, n, nil) を呼び出せば、整数 n の分割をすべて求めることができます。

●解答14

リスト : 分割数

def partition_number(n)
  let rec
    part_num = fn(n, k)
      if n < 0 or k < 1 then
        0
      else
        if n <= 1 or k == 1 then
          1
        else
          part_num(n - k, k) + part_num(n, k - 1)
        end
      end
    end
  in
    part_num(n, n)
  end
end

基本的な考え方は partition_of_integer と同じです。局所関数 part_num(n, k) は整数 n を k 以下の整数で分割する場合の総数を求めます。n < 0 または k < 1 の場合は整数を分割できないので 0 を返します。n <= 1 または k == 1 の場合、分割の方法は 1 通りしかないので 1 を返します。あとは part_num(n - k, k) と part_num(n, k - 1) の返り値を足し算するだけです。なお、このプログラムはナイーブな実装なため、実行速度はとても遅いです。ご注意くださいませ。

ところで、分割数は「動的計画法」を使うと高速に求めることができます。ご参考までにプログラムを示します。

リスト : 分割数 (動的計画法)

def partition_number1(n)
  let
    a = makeVector(n + 1, 1),
    k = 2
  in
    while k <= n do
      let m = k in
        while m <= n do
          a[m] = a[m] + a[m - k],
          m = m + 1
        end
      end,
      k = k + 1
    end,
    a[n]
  end
end

●解答15

集合を分割するアルゴリズムは簡単です。たとえば、n -1 個の要素 x1, ..., xn-1 を持つ集合を分割したところ、i 個の部分集合 S1, ..., Si が生成されたとしましょう。ここに、n 番目の要素 xn を追加すると、要素が n 個の集合を分割することができます。

新しい要素を追加する場合は次に示す手順で行います。

  1. 部分集合 Sk (k = 1 から i まで) に要素 xn を追加する
  2. 新しい部分集合 Si+1 (要素が xn だけの集合) を生成する

簡単な例を示しましょう。次の図を見てください。

部分集合を格納するリストを用意します。最初、部分集合は空集合なので空リストに初期化します。次に、要素 1 を追加します。部分集合は空リストなので、手順 1 は適用できません。手順 2 を適用して新しい部分集合 (1) を追加します。

次に要素 2 を追加します。((1)) に 手順 1 を適用すると、部分集合 (1) に要素を追加して ((1 2)) になります。手順 2 を適用すると、新しい部分集合 (2) を追加して ((1) (2)) になります。最後に 3 を追加します。((1 2)) に手順 1 を適用すると ((1 2 3)) に、手順 2 を適用すると ((1 2) (3)) になります。((1) (2)) に手順 1 を適用すると ((1 3) (2)) と ((1) (2 3)) になり、手順 2 を適用すると ((1) (2) (3)) になります。

このように、簡単な方法で集合を分割することができます。実際にプログラムを作る場合、上図を木と考えて、深さ優先で木をたどると簡単です。次のリストを見てください。

リスト : 集合の分割

def partition_of_set(f, xs)
  let rec
    part_set = fn(xs, a)
      if null(xs) then
        f(a)
      else
        foreach(fn(y)
                  part_set(cdr(xs),
                           cons(cons(car(xs), y),
                                remove(fn(x) equal(x, y) end, a)))
                end,
                a),
        part_set(cdr(xs), cons(list(car(xs)), a))
      end
    end
  in
    part_set(reverse(xs), nil)
  end
end

実際の処理は局所関数 part_set で行います。最初の if の then 節が再帰呼び出しの停止条件です。else 節の foreach で、部分集合 y に xs の先頭要素を追加し、それを a から y を取り除いたリストに追加します。最後に、新しい部分集合を a に追加します。ただし、このままでは要素の並び方が逆順になるので、part_set を呼び出す前に reverse でリスト xs を反転しています。これで集合の分割をすべて求めることができます。

●解答16

リスト : ベル数

# 添え字付き畳み込み
def fold_with_index(f, a, xs)
  let rec
    iter = fn(xs, i, a)
      if null(xs) then
        a
      else
        iter(cdr(xs), i + 1, f(car(xs), i, a))
      end
    end
  in
    iter(xs, 0, a)
  end
end

# ベル数
def bell_number(n)
  let rec
    iter = fn(i, bs)
      if i == n then
        car(bs)
      else
        iter(i + 1, 
             cons(fold_with_index(fn(x, k, a) comb(i, k) * x + a end, 0, bs), bs))
      end
    end
  in
    iter(0, list(1))
  end
end

bell_number は公式をそのままプログラムするだけです。累積変数 bs にベル数を逆順で格納します。\({}_n \mathrm{C}_k\) は関数 comb で求めます。\({}_n \mathrm{C}_k \times B(k)\) の総和は関数 fold_with_index で計算します。fold_with_index は添字を関数に渡して畳み込みを行います。匿名関数の引数 x がリストの要素、k が添字、a が累積変数です。bs は逆順になっていますが、二項係数は \({}_n \mathrm{C}_i\) と \({}_n \mathrm{C}_{n-i}\) の値が同じになるので、そのまま計算しても大丈夫です。もちろん、reverse で bs を逆順にしてから計算してもかまいません。

●解答17

リスト : 集合のグループ分け

def group_partition(f, n, m, xs)
  let rec
    group_part = fn(xs, a)
      if null(xs) then
        f(a)
      else
        foreach(
          fn(y)
            if length(y) < n then
              group_part(cdr(xs), cons(cons(car(xs), y),
                                       remove(fn(x) equal(x, y) end, a)))
            end
          end,
          a),
        if length(a) < m then
          group_part(cdr(xs), cons(list(car(xs)), a))
        end
      end
    end
  in
    group_part(reverse(xs), nil)
  end
end

group_partition は partition_of_set を改造するだけで簡単に作成することができます。生成する部分集合の大きさを n に、部分集合の個数を m に制限するだけです。部分集合に要素を追加する場合、length(y) が n 未満であることをチェックします。新しい部分集合を追加する場合、length(a) が m 未満であることをチェックします。これで集合をグループに分けることができます。

●解答18

グループ分けの総数は次の式で求めることができます。

\(\begin{array}{l} k = n \times m \\ \dfrac{{}_k \mathrm{C}_n \times {}_{k-n} \mathrm{C}_n \times {}_{k-2n} \mathrm{C}_n \times \cdots \times {}_{2n} \mathrm{C}_n \times {}_n \mathrm{C}_n}{m!} \end{array}\)

たとえば、n = 3, m = 5 の場合は次のようになります。

\( \dfrac{{}_{15} \mathrm{C}_3 \times {}_{12} \mathrm{C}_3 \times {}_9 \mathrm{C}_3 \times {}_6 \mathrm{C}_3 \times {}_3 \mathrm{C}_3}{5!} = 1401400 \)

これをそのままプログラムすると次のようになります。

リスト : グループ分けの総数

def group_partition_number(n, m)
  let rec
    iter = fn(k, a)
      if k == 0 then
        a / fact(m)
      else
        iter(k - n, a * comb(k, n))
      end
    end
  in
    iter(n * m, 1)
  end
end

階乗は関数 fact で、組み合わせの個数は関数 comb で計算します。要素の個数を変数 k にセットし、累積変数 a に comb(k, n) を乗算します。あとは k から n を減算し、k が 0 でなければ処理を繰り返すだけです。最後に a / fact(m) を計算して返します。

●解答19

素数を求める基本的な考え方は簡単です。最初に、2 から n までの整数列を生成します。先頭の 2 は素数なので、この整数列から 2 で割り切れる整数を取り除き除きます。2 で割り切れる整数が取り除かれたので、残った要素の先頭が素数になります。先頭要素は 3 になるので、今度は 3 で割り切れる整数を取り除けばいいのです。このように、素数を見つけたらそれで割り切れる整数を取り除いていくアルゴリズムを「エラトステネスの篩 (ふるい)」といいます。

プログラムは次のようになります。

リスト : 要素を求める

# 要素の削除 (破壊的修正)
def delete(f, xs)
  let rec
    top, cp = cons(nil, xs), top
  in
    while not null(cdr(cp)) do
      if f(cadr(cp)) then
        setCdr(cp, cddr(cp))
      else
        cp = cdr(cp)
      end
    end,
    cdr(top)
  end
end

# n 未満の素数をすべて求める
def sieve(f, n)
  let rec
    iter = fn(xs)
      let x = car(xs) in
        if x * x > n then
          foreach(f, xs)
        else
          f(x),
          iter(delete(fn(y) y % x == 0 end, cdr(xs)))
        end
      end
    end
  in
    iter(iota(2, n - 1))
  end
end

sieve の処理は局所関数 iter で行います。iota で 2 から n - 1 までの整数列を生成し、それを引数 xs に渡します。iter の中で、リスト xs の先頭要素を取り出して変数 x にセットします。この x が素数になります。f(x) を呼び出したあと、x で割り切れる要素を関数 delete で取り除きます。delete はリストを破壊的に修正します。remove を使うよりも実行速度は少しだけ速くなります。x が √n より大きい場合、xs をふるいにかける必要はありません。foreach で残りの要素に関数 f を適用します。

ご参考までにベクタを使ったプログラムを示します。

リスト : 素数を求める

def sieve1(f, n)
  let
    p = makeVector(n / 2 + 1, nil),  # 奇数を表す配列
    i = 3,                           # 整数を表す
    j = 1,                           # j = i / 2 (配列 p の位置)
    k = n / 2 + 1
  in
    f(2),
    while i * i <= n do
      if null(p[j]) then             # p[j] が nil ならは i は素数
        f(i),
        let m = j + i in             # i の倍数を削除
          while m < k do
            p[m] = 1,
            m = m + i
          end
        end
      end,
      i = i + 2,
      j = j + 1
    end,
    while i < n do
      if null(p[j]) then f(i) end,
      i = i + 2,
      j = j + 1
    end
  end
end

実行速度は sieve よりも sieve1 のほうが速くなります。

●解答20

リスト : 標準形ラテン方陣を求める

def check_latina(n, x, xs)
  find(fn(y) x == y end, map(fn(ys) nth(ys, n - 1) end, xs))
end

def latina(f, size)
  let rec
    solve = fn(n, xs, a, b)
      if null(xs) then
        if size - 1 == length(b) then
          f(reverse(cons(reverse(a), b)))
        else
          let 
            m = length(b) + 2
          in
            solve(2, remove(fn(x) x == m end, iota(1, size)), list(m), cons(reverse(a), b))
          end
        end
      else
        foreach(
          fn(x)
            if not check_latina(n, x, b) then
              solve(n + 1, remove(fn(z) z == x end, xs), cons(x, a), b)
            end
          end,
          xs
        )
      end
    end
  in
    solve(1, iota(1, size), nil, list(iota(1, size)))
  end
end

実際の処理は局所関数 solve で行います。基本的な考え方は完全順列とほぼ同じで、累積変数 a に順列を格納し、完成した順列を累積変数 b に格納します。引数 xs が空リストの場合、順列がひとつ完成しました。b の要素数をチェックして、size - 1 と等しければラテン方陣ができました。reverse(a) を b に追加して関数 f を評価します。そうでなければ solve を再帰呼び出しします。このとき、先頭要素は b の要素数 + 2 になることに注意してください。

順列を生成する場合、関数 check_latina を呼び出して数字 x を選択できるかチェックします。map で xs に格納されたリストの n - 1 番目の要素を nth で取り出します。リストは 0 から数えることに注意してください。そして、x が map の返り値に含まれているか関数 find でチェックします。x が含まれていれば、x を選択することはできません。そうでなければ x を選択します。

実行結果は次のようになります。

Calc> latina(fn(xs) foreach(fn(x) print(x), print("\n") end, xs), print("\n") end, 3);
(1 2 3)
(2 3 1)
(3 1 2)

Calc> latina(fn(xs) foreach(fn(x) print(x), print("\n") end, xs), print("\n") end, 4);
(1 2 3 4)
(2 1 4 3)
(3 4 1 2)
(4 3 2 1)

(1 2 3 4)
(2 1 4 3)
(3 4 2 1)
(4 3 1 2)

(1 2 3 4)
(2 3 4 1)
(3 4 1 2)
(4 1 2 3)

(1 2 3 4)
(2 4 1 3)
(3 1 4 2)
(4 3 2 1)

ちなみに、標準形ラテン方陣の総数は次のようになります。

I4 = 4
I5 = 56
I6 = 9408
I7 = 16942080

高次の標準形ラテン方陣の総数は、簡単に求めることができない非常にハードな問題だといわれています。


初版 2012 年 10 月 7 日
改訂 2021 年 6 月 5 日

Copyright (C) 2012-2021 Makoto Hiroi
All rights reserved.

[ PrevPage | SML/NJ | NextPage ]