F# Programming / F# Junk Scripts

遅延シーケンス

F# のシーケンス seq<'a> は遅延シーケンス (遅延リスト) として動作しますが、lazy 式を使えば私達でも遅延シーケンスをプログラムすることができます。今回は F# のお勉強ということで、あえて遅延シーケンスを実装してみましょう。

●仕様

module Lseq
  exception Empty_lseq
  type 'a lseq =
    private | Nil
            | Cons of 'a * Lazy<'a lseq>
  val empty: 'a lseq
  val isEmpty: _arg1: 'a lseq -> bool
  val cons: x: 'a -> s: Lazy<'a lseq> -> 'a lseq
  val head: _arg1: 'a lseq -> 'a
  val tail: _arg1: 'a lseq -> 'a lseq
  val unfold: gen: ('a -> ('b * 'a) option) -> state: 'a -> 'b lseq
  val initInfinite: proc: (int -> 'a) -> 'a lseq
  val item: n: int -> s: 'a lseq -> 'a
  val take: n: int -> s: 'a lseq -> 'a lseq
  val skip: n: int -> s: 'a lseq -> 'a lseq
  val ofList: _arg1: 'a list -> 'a lseq
  val toList: s: 'a lseq -> 'a list
  val append: s1: 'a lseq -> s2: 'a lseq -> 'a lseq
  val appendLazy: s1: 'a lseq -> s2: Lazy<'a lseq> -> 'a lseq
  val interleave: s1: 'a lseq -> s2: 'a lseq -> 'a lseq
  val concat: ss: 'a lseq lseq -> 'a lseq
  val map: proc: ('a -> 'b) -> s: 'a lseq -> 'b lseq
  val map2: proc: ('a -> 'b -> 'c) -> s1: 'a lseq -> s2: 'b lseq -> 'c lseq
  val concatMap: proc: ('a -> 'b lseq) -> s: 'a lseq -> 'b lseq
  val filter: pred: ('a -> bool) -> s: 'a lseq -> 'a lseq
  val fold: proc: ('a -> 'b -> 'a) -> a: 'a -> s: 'b lseq -> 'a
  val foldBack: proc: ('a -> 'b -> 'b) -> s: 'a lseq -> a: 'b -> 'b
  val iter: proc: ('a -> unit) -> s: 'a lseq -> unit
  val contains: x: 'a -> s: 'a lseq -> bool when 'a: equality
  val exists: pred: ('a -> bool) -> s: 'a lseq -> bool
  val forall: pred: ('a -> bool) -> s: 'a lseq -> bool
  val tryFind: pred: ('a -> bool) -> s: 'a lseq -> 'a option
  val tryFindIndex: pred: ('a -> bool) -> s: 'a lseq -> int option
  val takeWhile: pred: ('a -> bool) -> s: 'a lseq -> 'a lseq
  val skipWhile: pred: ('a -> bool) -> s: 'a lseq -> 'a lseq
  val union: s1: 'a lseq -> s2: 'a lseq -> 'a lseq when 'a: comparison
  val intersect: s1: 'a lseq -> s2: 'a lseq -> 'a lseq when 'a: comparison

●基本的な操作

> Lseq.empty;;
val it: 'a Lseq.lseq

> Lseq.empty |> Lseq.isEmpty;;
val it: bool = true

> let s0 = Lseq.cons 1 (lazy (Lseq.cons 2 (lazy (Lseq.cons 3 (lazy Lseq.empty)))));;
val s0: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.head s0;;
val it: int = 1

> Lseq.tail s0 |> Lseq.head;;
val it: int = 2

> Lseq.tail s0 |> Lseq.tail |> Lseq.head;;
val it: int = 3

> Lseq.tail s0 |> Lseq.tail |> Lseq.tail;;
val it: int Lseq.lseq = Nil

> let s1 = Lseq.initInfinite (fun x -> x + 1);;
val s1: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.isEmpty s1;;
val it: bool = false

> Lseq.head s1;;
val it: int = 1

> Lseq.tail s1 |> Lseq.head;;
val it: int = 2

> Lseq.tail s1 |> Lseq.tail |> Lseq.head;;
val it: int = 3

> Lseq.item 10 s1;;
val it: int = 11

> Lseq.item 100 s1;;
val it: int = 101

> let a = Lseq.ofList [1..10];;
val a: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.toList a;;
val it: int list = [1; 2; 3; 4; 5; 6; 7; 8; 9; 10]

●フィボナッチ数列

> let fibo = Lseq.unfold (fun (a, b) -> Some (a, (b, b + a))) (0I, 1I);;
val fibo: System.Numerics.BigInteger Lseq.lseq =
  Cons (0, Value is not created.)

> Lseq.take 10 fibo |> Lseq.iter (fun x -> printfn "%A" x);;
0
1
1
2
3
5
8
13
21
34
val it: unit = ()

> Lseq.skip 50 fibo |> Lseq.take 10 |> Lseq.iter (fun x -> printfn "%A" x);;
12586269025
20365011074
32951280099
53316291173
86267571272
139583862445
225851433717
365435296162
591286729879
956722026041
val it: unit = ()

> Lseq.item 100 fibo;;
val it: System.Numerics.BigInteger =
  354224848179261915075 {IsEven = false;
                         IsOne = false;
                         IsPowerOfTwo = false;
                         IsZero = false;
                         Sign = 1;}

●シーケンスの連結

> let s2 = Lseq.append (Lseq.ofList [1;2;3;4]) (Lseq.ofList [5;6;7;8]);;
val s2: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.take 8 s2 |> Lseq.toList;;
val it: int list = [1; 2; 3; 4; 5; 6; 7; 8]

> Lseq.interleave (Lseq.ofList [1;2;3;4]) (Lseq.ofList [5;6;7;8]) |> Lseq.take 8 |> Lseq.toList;;
val it: int list = [1; 5; 2; 6; 3; 7; 4; 8]

> let rec ones = Lseq.cons 1 (lazy ones);;
val ones: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.take 10 ones |> Lseq.toList;;
val it: int list = [1; 1; 1; 1; 1; 1; 1; 1; 1; 1]

> let rec twos = Lseq.cons 2 (lazy twos);;
val twos: int Lseq.lseq = Cons (2, Value is not created.)

> Lseq.take 10 twos |> Lseq.toList;;
val it: int list = [2; 2; 2; 2; 2; 2; 2; 2; 2; 2]

> Lseq.interleave ones twos |> Lseq.take 10 |> Lseq.toList;;
val it: int list = [1; 2; 1; 2; 1; 2; 1; 2; 1; 2]

> let s3 = Lseq.ofList [Lseq.ofList [1;2;3]; Lseq.ofList [4;5;6]; Lseq.ofList [7;8;9]] |> Lseq.concat;;
val s3: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.take 9 s3 |> Lseq.toList;;
val it: int list = [1; 2; 3; 4; 5; 6; 7; 8; 9]

●高階関数

> let ints = Lseq.initInfinite (fun x -> x + 1);;
val ints: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.take 10 ints |> Lseq.toList;;
val it: int list = [1; 2; 3; 4; 5; 6; 7; 8; 9; 10]

> let squares = Lseq.map (fun x -> x * x) ints;;
val squares: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.take 10 squares |> Lseq.toList;;
val it: int list = [1; 4; 9; 16; 25; 36; 49; 64; 81; 100]

> let addseq s1 s2 = Lseq.map2 (+) s1 s2;;
val addseq: s1: int Lseq.lseq -> s2: int Lseq.lseq -> int Lseq.lseq

> let evens = addseq ints ints;;
val evens: int Lseq.lseq = Cons (2, Value is not created.)

> Lseq.take 10 evens |> Lseq.toList;;
val it: int list = [2; 4; 6; 8; 10; 12; 14; 16; 18; 20]

> let s4 = Lseq.concatMap (fun x -> Lseq.ofList([1..x])) ints;;
val s4: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.take 21 s4 |> Lseq.toList;;
val it: int list =
  [1; 1; 2; 1; 2; 3; 1; 2; 3; 4; 1; 2; 3; 4; 5; 1; 2; 3; 4; 5; 6]

> Lseq.filter (fun x -> x % 2 = 0) ints |> Lseq.take 10 |> Lseq.toList;;
val it: int list = [2; 4; 6; 8; 10; 12; 14; 16; 18; 20]

> Lseq.filter (fun x -> x % 2 <> 0) ints |> Lseq.take 10 |> Lseq.toList;;
val it: int list = [1; 3; 5; 7; 9; 11; 13; 15; 17; 19]

> Lseq.take 100 ints |> Lseq.fold (+) 0;;
val it: int = 5050

> Lseq.foldBack (+) (Lseq.take 100 ints) 0;;
val it: int = 5050

●フィボナッチ数列 (2)

> #nowarn "40";;
> let addseq (s1: bigint Lseq.lseq) (s2: bigint Lseq.lseq) = Lseq.map2 (+) s1 s2;;
val addseq:
  s1: bigint Lseq.lseq -> s2: bigint Lseq.lseq
    -> System.Numerics.BigInteger Lseq.lseq

> let rec fibs = Lseq.cons 0I (lazy (Lseq.cons 1I (lazy (addseq (Lseq.tail fibs) fibs))));;
val fibs: System.Numerics.BigInteger Lseq.lseq =
   Cons (0, Value is not created.)

> Lseq.take 10 fibs |> Lseq.iter (fun x -> printfn "%A" x);;
0
1
1
2
3
5
8
13
21
34
val it: unit = ()

> Lseq.skip 50 fibs |> Lseq.take 10 |> Lseq.iter (fun x -> printfn "%A" x);;
12586269025
20365011074
32951280099
53316291173
86267571272
139583862445
225851433717
365435296162
591286729879
956722026041
val it: unit = ()

> for i = 100 to 109 do Lseq.item i fibs |> printfn "%A";;
354224848179261915075
573147844013817084101
927372692193078999176
1500520536206896083277
2427893228399975082453
3928413764606871165730
6356306993006846248183
10284720757613717413913
16641027750620563662096
26925748508234281076009
val it: unit = ()

●検査と探索

> let randgen = System.Random();;
val randgen: System.Random

> let randseq = Lseq.initInfinite (fun _ -> randgen.Next(0x7fffffff));;
val randseq: int Lseq.lseq = Cons (497462975, Value is not created.)

> Lseq.take 10 randseq |> Lseq.iter (fun x -> printfn "%d" x);;
497462975
1244616573
1111397384
2092954616
569463218
1534196420
1407322786
317110337
400384826
793399835
val it: unit = ()

> Lseq.take 10 randseq |> Lseq.contains 400384826;;
val it: bool = true

> Lseq.take 10 randseq |> Lseq.contains 400384827;;
val it: bool = false

> Lseq.take 10 randseq |> Lseq.exists (fun x -> x % 2 = 0);;
val it: bool = true

> Lseq.take 10 randseq |> Lseq.exists (fun x -> x > 2092954616);;
val it: bool = false

> Lseq.take 10 randseq |> Lseq.forall (fun x -> x <= 2092954616);;
val it: bool = true

> Lseq.take 10 randseq |> Lseq.forall (fun x -> x >= 2092954616);;
val it: bool = false

> Lseq.take 10 randseq |> Lseq.tryFind (fun x -> x % 2 = 0);;
val it: int option = Some 1111397384

> Lseq.take 10 randseq |> Lseq.tryFindIndex (fun x -> x % 2 = 0);;
val it: int option = Some 2

> Lseq.take 10 randseq |> Lseq.tryFind (fun x -> x = 0);;
val it: int option = None

> Lseq.take 10 randseq |> Lseq.tryFindIndex (fun x -> x = 0);;
val it: int option = None

> Lseq.takeWhile (fun x -> x < 16) ints |> Lseq.toList;;
val it: int list = [1; 2; 3; 4; 5; 6; 7; 8; 9; 10; 11; 12; 13; 14; 15]

> Lseq.skipWhile (fun x -> x < 16) ints |> Lseq.take 10 |> Lseq.toList;;
val it: int list = [16; 17; 18; 19; 20; 21; 22; 23; 24; 25]

●集合演算

> let s1 = Lseq.initInfinite (fun x -> (x + 1) * (x + 2) / 2);;   // 三角数
val s1: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.take 10 s1 |> Lseq.toList;;
val it: int list = [1; 3; 6; 10; 15; 21; 28; 36; 45; 55]

> let s2 = Lseq.initInfinite (fun x -> (x + 1) * (x + 1));;       // 四角数 (平方数)
val s2: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.take 10 s2 |> Lseq.toList;;
val it: int list = [1; 4; 9; 16; 25; 36; 49; 64; 81; 100]

> Lseq.union s1 s2 |> Lseq.take 20 |> Lseq.toList;;
val it: int list =
  [1; 3; 4; 6; 9; 10; 15; 16; 21; 25; 28; 36; 45; 49; 55; 64; 66; 78; 81; 91]

> Lseq.intersect s1 s2 |> Lseq.take 5 |> Lseq.toList;;             // 平方三角数
val it: int list = [1; 36; 1225; 41616; 1413721]

●ハミングの問題

ここで union を使うと簡単に解ける問題を紹介しましょう。

[ハミングの問題]

7 以上の素数で割り切れない正の整数を小さい順に N 個求めよ

参考文献 : 奥村晴彦,『Ｃ言語による最新アルゴリズム事典』, 技術評論社, 1991 (361 ページより引用)

7 以上の素数で割り切れない正の整数は、素因子が 2, 3, 5 しかない自然数のことで、これを「ハミング数 (Hamming Numbers)」といいます。ハミング数は素因数分解したとき、2ⁱ * 3^j * 5^k (i, j, k >= 0) の形式になります。たとえば、100 以下のハミング数は次のようになります。

1, 2, 3, 4, 5, 6, 8, 9, 10, 12, 15, 16, 18, 20, 24, 25, 27, 30, 32, 36, 40, 45, 48, 50, 
54, 60, 64, 72, 75, 80, 81, 90, 96, 100

遅延シーケンスを使うと「ハミングの問題」は簡単に解くことができます。小さい順にハミング数を出力する遅延シーケンスを hs としましょう。hs は 1 から始まるので次のように定義できます。

let rec hs = cons 1 (lazy (...))

最初の要素は 1 なので、それに 2, 3, 5 を掛け算した値 (2, 3, 5) もハミング数になります。この値は次の式で生成することができます。

map (fun x -> x * 2) hs
map (fun x -> x * 3) hs
map (fun x -> x * 5) hs

あとは、これらの遅延ストリームを union でひとつにまとめて、小さい順に出力すればいいわけです。

プログラムと実行結果を示します。

リスト : ハミングの問題

let rec hs = 
  Lseq.cons 1
            (lazy (Lseq.union
                    (Lseq.map (fun x -> x * 2) hs)
                    (Lseq.union
                      (Lseq.map (fun x -> x * 3) hs)
                      (Lseq.map (fun x -> x * 5) hs))))

> let rec hs =
-   Lseq.cons 1
-             (lazy (Lseq.union
-                     (Lseq.map (fun x -> x * 2) hs)
-                     (Lseq.union
-                       (Lseq.map (fun x -> x * 3) hs)
-                       (Lseq.map (fun x -> x * 5) hs))));;
val hs: int Lseq.lseq = Cons (1, Value is not created.)

> Lseq.take 100 hs |> Lseq.toList;;
val it: int list =
  [1; 2; 3; 4; 5; 6; 8; 9; 10; 12; 15; 16; 18; 20; 24; 25; 27; 30; 32; 36; 40;
   45; 48; 50; 54; 60; 64; 72; 75; 80; 81; 90; 96; 100; 108; 120; 125; 128;
   135; 144; 150; 160; 162; 180; 192; 200; 216; 225; 240; 243; 250; 256; 270;
   288; 300; 320; 324; 360; 375; 384; 400; 405; 432; 450; 480; 486; 500; 512;
   540; 576; 600; 625; 640; 648; 675; 720; 729; 750; 768; 800; 810; 864; 900;
   960; 972; 1000; 1024; 1080; 1125; 1152; 1200; 1215; 1250; 1280; 1296; 1350;
   1440; 1458; 1500; 1536]

●素数列

リスト : エラトステネスの篩

let rec sieve (s: int Lseq.lseq) =
  let x = Lseq.head s
  Lseq.cons x (lazy (sieve (Lseq.filter (fun n -> n % x <> 0) (Lseq.tail s))))

> let rec sieve (s: int Lseq.lseq) =
-   let x = Lseq.head s
-   Lseq.cons x (lazy (sieve (Lseq.filter (fun n -> n % x <> 0) (Lseq.tail s))));;
val sieve: s: int Lseq.lseq -> int Lseq.lseq

> Lseq.initInfinite (fun x -> x + 2) |> sieve |> Lseq.take 25 |> Lseq.toList;;
val it: int list =
  [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]

リスト : 素数列の生成 (高速化)

let rec primesFrom (n: int) =
  if isPrime n then Lseq.cons n (lazy (primesFrom (n + 2)))
  else primesFrom (n + 2)
and isPrime n =
  let rec iter s =
    let p = Lseq.head s
    if p * p > n then true
    else if n % p = 0 then false
    else Lseq.tail s |> iter
  Lseq.tail primes |> iter
and primes = Lseq.cons 2 (lazy (Lseq.cons 3 (lazy (Lseq.cons 5 (lazy (primesFrom 7))))))

> let rec primesFrom (n: int) =
-   if isPrime n then Lseq.cons n (lazy (primesFrom (n + 2)))
-   else primesFrom (n + 2)
- and isPrime n =
-   let rec iter s =
-     let p = Lseq.head s
-     if p * p > n then true
-     else if n % p = 0 then false
-     else Lseq.tail s |> iter
-   Lseq.tail primes |> iter
- and primes = Lseq.cons 2 (lazy (Lseq.cons 3 (lazy (Lseq.cons 5 (lazy (primesFrom 7))))));;
val primesFrom: n: int -> int Lseq.lseq
val isPrime: n: int -> bool
val primes: int Lseq.lseq = Cons (2, Value is not created.)

> Lseq.take 25 primes |> Lseq.toList;;
val it: int list =
  [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]

> Lseq.item 99 primes;;
val it: int = 541

> Lseq.item 500 primes;;
val it: int = 3581

●双子素数

差が 2 である素数の組を「双子素数 (twin prime)」といいます。素数列 primes を使うと双子素数は簡単に求めることができます。

> let twin = Lseq.filter
-                (fun (x, y) -> y - x = 2)
-                (Lseq.map2 (fun x y -> (x, y)) primes (Lseq.tail primes));;
val twin: (int * int) Lseq.lseq = Cons ((3, 5), Value is not created.)

> Lseq.take 50 twin |> Lseq.toList;;
val it: (int * int) list =
  [(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); (1019, 1021); (1031, 1033);
   (1049, 1051); (1061, 1063); (1091, 1093); (1151, 1153); (1229, 1231);
   (1277, 1279); (1289, 1291); (1301, 1303); (1319, 1321); (1427, 1429);
   (1451, 1453); (1481, 1483); (1487, 1489)]

ところで、双子素数 - Wikipedia によると、『双子素数は無数に存在するかという問題、いわゆる「双子素数の予想」や「双子素数の問題」は、いまだに数学上の未解決問題である。無数に存在するだろう、とは、多くの数論学者が予想している。』とのことです。

●素因数分解

素数列 primes を使った単純な素因数分解です。

リスト : 素因数分解

let factorization n =
  let rec factor_sub n m c =
    if n % m <> 0 then (c, n)
    else factor_sub (n / m) m (c + 1)

  let rec iter ps n a =
    let i = Lseq.head ps
    if n = 1 then List.rev a
    else if n < i * i then List.rev ((n, 1)::a)
    else
      let (c, m) = factor_sub n i 0
      if c = 0 then iter (Lseq.tail ps) n a
      else iter (Lseq.tail ps) m ((i, c)::a)

  iter primes n []

> factorization 12345678;;
val it: (int * int) list = [(2, 1); (3, 2); (47, 1); (14593, 1)]

> factorization 123456789;;
val it: (int * int) list = [(3, 2); (3607, 1); (3803, 1)]

> factorization 1234567890;;
val it: (int * int) list = [(2, 1); (3, 2); (5, 1); (3607, 1); (3803, 1)]

> factorization 1111111111;;
val it: (int * int) list = [(11, 1); (41, 1); (271, 1); (9091, 1)]

●参考 URL

●プログラムリスト

//
// lseq.fsx : 遅延シーケンス
//
//            Copyright (C) 2022 Makoto Hiroi
//
module Lseq

// 例外
exception Empty_lseq

// データ型の定義
type 'a lseq = private Nil | Cons of 'a * Lazy<'a lseq>

// 終端
let empty = Nil

// 終端の判定
let isEmpty = function
  Nil -> true
| _ -> false

// コンストラクタ
let cons x s = Cons(x, s)

// 先頭要素を取り出す
let head = function
  Nil -> raise Empty_lseq
| Cons (x, _) -> x

// 先頭要素を取り除く
let tail = function
  Nil -> raise Empty_lseq
| Cons (_, x) -> x.Force()

// 解きほぐし
let rec unfold gen state =
  match gen state with
    None -> Nil
  | Some (x, nextState) -> Cons (x, lazy (unfold gen nextState))

// 無限シーケンス
let rec initInfinite proc =
  unfold (fun x -> Some (proc x, x + 1)) 0

// n 番目の要素を求める
let rec item n s =
  if isEmpty s then raise Empty_lseq
  else if n = 0 then head s
  else item (n - 1) (tail s)

// 先頭から n 個の要素を取り出す
let rec take n s =
  if isEmpty s || n = 0 then Nil
  else Cons(head s, lazy (take (n - 1) (tail s)))

// 先頭から n 個の要素を取り除く
let rec skip n s =
  if isEmpty s || n = 0 then s
  else skip (n - 1) (tail s)

// list -> lseq
let rec ofList = function
  [] -> Nil
| x::xs -> Cons (x, lazy (ofList xs))

// lseq -> list
let rec toList s =
  if isEmpty s then []
  else head s :: toList (tail s)

// 連結
let rec append s1 s2 =
  if isEmpty s1 then s2
  else Cons (head s1, lazy (append (tail s1) s2))

let rec appendLazy<'a> s1 (s2: Lazy<'a lseq>) =
  if isEmpty s1 then s2.Force()
  else Cons (head s1, lazy (appendLazy (tail s1) s2))

let rec interleave s1 s2 =
  if isEmpty s1 then s2
  else Cons (head s1, lazy (interleave s2 (tail s1)))

let rec concat ss =
  if isEmpty ss then Nil
  else appendLazy (head ss) (lazy (concat (tail ss)))

// 高階関数
let rec map proc s =
  if isEmpty s then Nil
  else Cons (proc (head s), lazy (map proc (tail s)))

let rec map2 proc s1 s2 =
  if isEmpty s1 || isEmpty s2 then Nil
  else Cons (proc (head s1) (head s2), lazy (map2 proc (tail s1) (tail s2)))

let rec concatMap proc s =
  if isEmpty s then Nil
  else appendLazy (proc (head s)) (lazy (concatMap proc (tail s)))

let rec filter pred s =
  if isEmpty s then Nil
  else if pred (head s) then Cons(head s, lazy (filter pred (tail s)))
  else filter pred (tail s)

let rec fold proc a s =
  if isEmpty s then a
  else fold proc (proc a (head s)) (tail s)

let rec foldBack proc s a =
  if isEmpty s then a
  else proc (head s) (foldBack proc (tail s) a)

let rec iter proc s =
  if isEmpty s then ()
  else (
    proc (head s)
    iter proc (tail s)
  )

// 検査と探索
let rec contains x s =
  if isEmpty s then false
  else if (head s) = x then true
  else contains x (tail s)

let rec exists pred s =
  if isEmpty s then false
  else if pred (head s) then true
  else exists pred (tail s)

let rec forall pred s =
  if isEmpty s then true
  else if not (pred (head s)) then false
  else forall pred (tail s)

let rec tryFind pred s =
  if isEmpty s then None
  else if pred (head s) then Some (head s)
  else tryFind pred (tail s)

let tryFindIndex pred s =
  let rec iter i s =
    if isEmpty s then None
    else if pred (head s) then Some i
    else iter (i + 1) (tail s)
  iter 0 s

let rec takeWhile pred s =
  if isEmpty s || not (pred (head s)) then Nil
  else Cons (head s, lazy (takeWhile pred (tail s)))

let rec skipWhile pred s =
  if isEmpty s || not (pred (head s)) then s
  else skipWhile pred (tail s)

// 和集合
let rec union s1 s2 =
  if isEmpty s1 then s2
  else if isEmpty s2 then s1
  else
    let x = head s1
    let y = head s2
    if x = y then Cons (x, lazy (union (tail s1) (tail s2)))
    else if x < y then Cons (x, lazy (union (tail s1) s2))
    else Cons (y, lazy (union s1 (tail s2)))

// 積集合
let rec intersect s1 s2 =
  if isEmpty s1 || isEmpty s2 then Nil
  else
    let x = head s1
    let y = head s2
    if x = y then Cons (x, lazy (intersect (tail s1) (tail s2)))
    else if x < y then intersect (tail s1) s2
    else intersect s1 (tail s2)

スライドパズル

●NO-OFF パズル

[問題]

問題 A, B から GOAL までの最短手順を求めてください。

スライドパズル NO-OFF は、問題 A の "ON-OFF" を GOAL のように "NO-OFF" にチェンジするパズルです。NO-OFF は芦ヶ原伸之氏が考案されたパズルです。このパズルは局面の総数が 540 通りしかないにもかかわらず、手数がけっこうかかる面白いパズルです。詳しい説明は拙作のページ Puzzle DE Programming NO-OFF パズルをお読みくださいませ。

●プログラムリスト

//
// no_off.fsx : NO-OFF puzzle
//
//              Copyright (C) 2022 Makoto Hiroi
//

// キュー
#load "queue.fs"

// 駒
type P =
  | S  = 0
  | L1 = 1
  | L2 = 2
  | N  = 3
  | F  = 4
  | O  = 5

// 局面
type state =
  | Nil
  | State of int * P[] * state

  member this.Space =
    match this with
      Nil -> failwith "empty state"
    | State (s, _, _) -> s

  member this.Board =
    match this with
      Nil -> failwith "empty state"
    | State (_, b, _) -> b

// 隣接リスト
// 0 1 2 3
// 4 5 6 7
let adjacent = [|
  [1; 4]     // 0
  [0; 2; 5]  // 1
  [1; 3; 6]  // 2
  [2; 7]     // 3
  [0; 5]     // 4
  [1; 4; 6]  // 5
  [2; 5; 7]  // 6
  [3; 6]  |] // 7

// 問題 A : 44 move
let qa = [|
  P.L1; P.L2; P.O; P.N;
  P.O;  P.F;  P.F; P.S  |]

// 問題 B: 56 move
let qb = [|
  P.N; P.O; P.L1; P.L2;
  P.F; P.O; P.F;  P.S  |]

// ゴール
let goal = [|
  P.L1; P.L2; P.N; P.O;
  P.O;  P.F;  P.F; P.S  |]

// 駒を文字列に変換
let pToStr = function
  P.L1 -> "L"
| P.L2 -> "L"
| P.N  -> "N"
| P.F  -> "F"
| P.O  -> "O"
| _    -> "_"

// 盤面の表示
let print_board (board: P[]) =
  List.iter (fun x -> board.[x] |> pToStr |> printf "%s ") [0..3]
  printf "\n"
  List.iter (fun x -> board.[x] |> pToStr |> printf "%s ") [4..7]
  printf "\n\n"

// 手順の表示
let rec print_answer = function
  Nil -> ()
| State (_, b, prev) -> print_answer prev; print_board b

// 駒の移動
let move_piece x = function
  Nil -> failwith "empty state"
| State (s, b, _) as st ->
    let p = b.[x]
    let newb = Array.copy b
    if p = P.L1 then
      if s > x then None
      else (
        // 電球を左へ動かす
        newb.[s] <- newb.[x]
        newb.[x] <- newb.[x + 1]
        newb.[x + 1] <- P.S
        Some (State(x + 1, newb, st))
      )
    else if p = P.L2 then
      if s > 3 then None
      else (
        // 電球を右へ動かす
        newb.[s] <- newb.[x]
        newb.[x] <- newb.[x - 1];
        newb.[x - 1] <- P.S
        Some (State(x - 1, newb, st))
      )
    else (
      newb.[s] <- newb.[x]
      newb.[x] <- P.S
      Some (State(x, newb, st))
    )

// 幅優先探索
let bfs start goal =
  let rec iter (que, table) =
    if Queue.is_empty que then ()
    else
      let st: state = Queue.top que
      if st.Board = goal then
        print_answer st
      else
        List.fold
          (fun (q, tbl) x ->
            match move_piece x st with
              None -> (q, tbl)
            | Some newst ->
                if List.contains newst.Board tbl then (q, tbl)
                else (Queue.enqueue newst q, newst.Board::tbl))
          (Queue.dequeue que, table)
          adjacent.[st.Space] |> iter
  let st = State (Array.findIndex (fun x -> x = P.S) start, start, Nil)
  (Queue.enqueue st Queue.create, [start]) |> iter

●実行結果

> open No_off;;
> bfs qa goal;;
L L O N
O F F _

L L O _
O F F N

・・省略・・

L L N _
O F F O

L L N O
O F F _

val it: unit = ()

> bfs qb goal;;
N O L L
F O F _

N O L L
F O _ F

・・省略・・

L L N _
O F F O

L L N O
O F F _

val it: unit = ()

●問題 A の解答

L L が電球を表し、_ が空き場所を表します。

  (0)        (1)        (2)        (3)        (4)        (5)        (6)        (7)
  L L O N    L L O _    L L _ O    _ L L O    O L L O    O L L O    O L L O    O L L O 
  O F F _    O F F N    O F F N    O F F N    _ F F N    F _ F N    F F _ N    F F N _ 

  (8)        (9)        (10)       (11)       (12)       (13)       (14)       (15)
  O L L _    O _ L L    _ O L L    F O L L    F O L L    F _ L L    F L L _    F L L O 
  F F N O    F F N O    F F N O    _ F N O    F _ N O    F O N O    F O N O    F O N _ 

  (16)       (17)       (18)       (19)       (20)       (21)       (22)       (23)
  F L L O    F L L O    F L L O    _ L L O    L L _ O    L L O _    L L O N    L L O N 
  F O _ N    F _ O N    _ F O N    F F O N    F F O N    F F O N    F F O _    F F _ O 

  (24)       (25)       (26)       (27)       (28)       (29)       (30)       (31)
  L L _ N    _ L L N    F L L N    F L L N    F L L N    F L L N    F L L _    F _ L L 
  F F O O    F F O O    _ F O O    F _ O O    F O _ O    F O O _    F O O N    F O O N 

  (32)       (33)       (34)       (35)       (36)       (37)       (38)       (39)
  F O L L    F O L L    _ O L L    O _ L L    O L L _    O L L N    O L L N    O L L N 
  F _ O N    _ F O N    F F O N    F F O N    F F O N    F F O _    F F _ O    F _ F O 

  (40)       (41)       (42)       (43)       (44)
  O L L N    _ L L N    L L _ N    L L N _    L L N O    
  _ F F O    O F F O    O F F O    O F F O    O F F _

●問題 B の解答

L L が電球を表し、_ が空き場所を表します。

  (0)        (1)        (2)        (3)        (4)        (5)        (6)        (7)
  N O L L    N O L L    N O L L    N _ L L    N L L _    N L L F    N L L F    N L L F 
  F O F _    F O _ F    F _ O F    F O O F    F O O F    F O O _    F O _ O    F _ O O 

  (8)        (9)        (10)       (11)       (12)       (13)       (14)       (15)
  N L L F    _ L L F    L L _ F    L L F _    L L F O    L L F O    L L _ O    _ L L O 
  _ F O O    N F O O    N F O O    N F O O    N F O _    N F _ O    N F F O    N F F O 

  (16)       (17)       (18)       (19)       (20)       (21)       (22)       (23)
  N L L O    N L L O    N L L O    N L L O    N L L _    N _ L L    _ N L L    F N L L 
  _ F F O    F _ F O    F F _ O    F F O _    F F O O    F F O O    F F O O    _ F O O 

  (24)       (25)       (26)       (27)       (28)       (29)       (30)       (31)
  F N L L    F _ L L    F L L _    F L L O    F L L O    F L L O    F L L O    _ L L O 
  F _ O O    F N O O    F N O O    F N O _    F N _ O    F _ N O    _ F N O    F F N O 

  (32)       (33)       (34)       (35)       (36)       (37)       (38)       (39)
  L L _ O    L L N O    L L N O    L L N _    L L _ N    _ L L N    F L L N    F L L N 
  F F N O    F F _ O    F F O _    F F O O    F F O O    F F O O    _ F O O    F _ O O 

  (40)       (41)       (42)       (43)       (44)       (45)       (46)       (47)
  F L L N    F L L N    F L L _    F _ L L    F O L L    F O L L    _ O L L    O _ L L 
  F O _ O    F O O _    F O O N    F O O N    F _ O N    _ F O N    F F O N    F F O N 

  (48)       (49)       (50)       (51)       (52)       (53)       (54)       (55)
  O L L _    O L L N    O L L N    O L L N    O L L N    _ L L N    L L _ N    L L N _ 
  F F O N    F F O _    F F _ O    F _ F O    _ F F O    O F F O    O F F O    O F F O 

  (56)
  L L N O 
  O F F _

●8 パズル

皆さんお馴染みの「15 パズル」をひとまわり小さくした「8 パズル」を解くプログラムです。解法の詳しい説明は、以下の拙作のページをお読みくださいませ。

Ｃ言語入門: スライドパズル
OCaml 入門: パズルの解法 (2)
Scala 入門: パズルの解法 (3)
Common Lisp 入門: 幅優先探索とスライドパズル
Julia 入門: Puzzle DE Julia!! スライドパズル編

●プログラムリスト

//
// eight.ml : 8 Puzzle
//
//            Copyright (C) 2022 Makoto Hiroi
//
type hset<'a> = System.Collections.Generic.HashSet<'a>
#load "queue.fs"

// 隣接リスト
let adjacent = [|
  [1; 3]
  [0; 2; 4]
  [1; 5]
  [0; 4; 6]
  [1; 3; 5; 7]
  [2; 4; 8]
  [3; 7]
  [4; 6; 8]
  [5; 7] |]

// 局面
[<CustomEquality; NoComparison>]
type state =
  | Nil
  | State of int * int[] * state

  // 空き場所
  member this.Space =
    match this with
      Nil -> failwith "empty state"
    | State (sp, _, _) -> sp

  // 盤面
  member this.Board =
    match this with
      Nil -> failwith "empty state"
    | State (_, b, _) -> b

  // 駒の移動
  member this.move_piece (x) =
    match this with
      Nil -> failwith "empty state"
    | State(sp, b, _) -> let newb = Array.copy b
                         newb.[sp] <- newb.[x]
                         newb.[x] <- 0
                         State(x, newb, this)

  // 等値演算子
  override this.Equals (other: obj) =
    let that = other :?> state
    match (this, that) with
      (Nil, Nil) -> true
    | (State (_, b1, _), State (_, b2, _)) -> b1 = b2
    | (_, _) -> false

  // ハッシュ関数
  override this.GetHashCode () =
    match this with
      Nil -> 0
    | State(_, b, _) -> Array.fold (fun a x -> a * 10 + x) 0 b

// 手順の表示
let rec print_answer = function
  Nil -> ()
| State (_, board, prev) ->
    print_answer prev;
    printfn "%A" board

// 幅優先探索

// 同一局面のチェックに Set を使用
let bfs start goal =
  let start_st = State(Array.findIndex (fun x -> x = 0) start, start, Nil)
  let goal_st  = State(Array.findIndex (fun x -> x = 0) goal,  goal,  Nil)
  let rec iter (que, table) =
    if Queue.is_empty que then ()
    else
      let st = Queue.top que
      if st = goal_st then print_answer st
      else
        List.fold
          (fun (q, tbl: Set<int[]>) x ->
            let newst = st.move_piece(x)
            if tbl.Contains(newst.Board) then (q, tbl)
            else (Queue.enqueue newst q, tbl.Add(newst.Board)))
          (Queue.dequeue que, table)
          adjacent.[st.Space] |> iter
  iter (Queue.enqueue start_st Queue.create, Set [start])

// 同一局面のチェックに HashSet を使用
let bfs1 start goal =
  let start_st = State(Array.findIndex (fun x -> x = 0) start, start, Nil)
  let goal_st  = State(Array.findIndex (fun x -> x = 0) goal,  goal,  Nil)
  let tbl = hset<state>()
  let rec iter que =
    if Queue.is_empty que then ()
    else
      let st = Queue.top que
      if st = goal_st then print_answer st
      else
        List.fold
          (fun q x ->
            let newst = st.move_piece(x)
            if tbl.Add(newst) then
              Queue.enqueue newst q
            else q)
          (Queue.dequeue que)
          adjacent.[st.Space] |> iter
  tbl.Add(start_st) |> ignore
  iter (Queue.enqueue start_st Queue.create)

// 最長手数の局面を求める
let bfs2 start =
  let start_st = State(Array.findIndex (fun x -> x = 0) start, start, Nil)
  let tbl = hset<state>()
  let rec iter i xs =
    let ys = List.fold (fun a (st: state) ->
                         List.fold
                           (fun b x ->
                             let newst = st.move_piece(x)
                             if tbl.Add(newst) then (newst::b) else b)
                           a
                           adjacent.[st.Space])
                       []
                       xs
    if List.isEmpty ys then (
      printfn "max = %d" i
      List.iter (fun (st: state) -> printfn "%A" st.Board) xs
    ) else iter (i + 1) ys

  tbl.Add(start_st) |> ignore
  iter 0 [start_st]

//
// 反復深化＋下限値枝刈り法
//

// 移動距離表
let distance_tbl = [|
  [|0; 0; 0; 0; 0; 0; 0; 0; 0|]
  [|0; 1; 2; 1; 2; 3; 2; 3; 4|]
  [|1; 0; 1; 2; 1; 2; 3; 2; 3|]
  [|2; 1; 0; 3; 2; 1; 4; 3; 2|]
  [|1; 2; 3; 0; 1; 2; 1; 2; 3|]
  [|2; 1; 2; 1; 0; 1; 2; 1; 2|]
  [|3; 2; 1; 2; 1; 0; 3; 2; 1|]
  [|2; 3; 4; 1; 2; 3; 0; 1; 2|]
  [|3; 2; 3; 2; 1; 2; 1; 0; 1|]  |]

// 移動距離を求める
let get_distance (board: int[]) =
  Array.fold (fun a n -> a + distance_tbl.[board.[n]].[n]) 0 board

let ids (board: int[]) goal =
  let mutable count = 0
  let rec dfs n limit space move_list lower =
    if n = limit then
      if board = goal then (
          count <- count + 1
          printfn "%A" (List.tail (List.rev move_list))
      ) else ()
    else
      List.iter
        (fun x ->
          let p = board.[x]
          if p <> List.head move_list then
            let new_lower = lower - distance_tbl.[p].[x] + distance_tbl.[p].[space]
            if new_lower + n <= limit then (
              // 駒を動かす
              board.[space] <- p
              board.[x] <- 0
              // 再帰呼び出し
              dfs (n+1) limit x (p::move_list) new_lower
              // 元に戻す
              board.[space] <- 0
              board.[x] <- p
            ) else ()
          else ())
        adjacent.[space]
  let lower = get_distance board
  let mutable i = lower
  while (i <= 31 && count = 0) do
    printfn "----- %d -----" i
    dfs 0 i (Array.findIndex (fun x -> x = 0) board) [-1] lower
    i <- i + 1

let solver  start = bfs  start [|1;2;3;4;5;6;7;8;0|]
let solver1 start = bfs1 start [|1;2;3;4;5;6;7;8;0|]
let solver2 start = ids  start [|1;2;3;4;5;6;7;8;0|]
let solverMax ()  = bfs2 [|1;2;3;4;5;6;7;8;0|]

●実行結果

> #time;;

--> 今すぐタイミング オン

> Eight.solver [|8; 6; 7; 2; 5; 4; 3; 0; 1|];;
[|8; 6; 7; 2; 5; 4; 3; 0; 1|]
[|8; 6; 7; 2; 0; 4; 3; 5; 1|]
[|8; 0; 7; 2; 6; 4; 3; 5; 1|]

  ・・・省略・・・

[|1; 2; 3; 4; 5; 6; 0; 7; 8|]
[|1; 2; 3; 4; 5; 6; 7; 0; 8|]
[|1; 2; 3; 4; 5; 6; 7; 8; 0|]
リアル: 00:00:10.841、CPU: 00:00:10.890、GC 全般0: 1123, 全般1: 15, 全般2: 1
val it: unit = ()

> Eight.solver1 [|8; 6; 7; 2; 5; 4; 3; 0; 1|];;
[|8; 6; 7; 2; 5; 4; 3; 0; 1|]
[|8; 6; 7; 2; 0; 4; 3; 5; 1|]
[|8; 0; 7; 2; 6; 4; 3; 5; 1|]

  ・・・省略・・・

[|1; 2; 3; 4; 5; 6; 0; 7; 8|]
[|1; 2; 3; 4; 5; 6; 7; 0; 8|]
[|1; 2; 3; 4; 5; 6; 7; 8; 0|]
リアル: 00:00:00.513、CPU: 00:00:00.640、GC 全般0: 16, 全般1: 6, 全般2: 1
val it: unit = ()

> Eight.solverMax();;
max = 31
[|8; 6; 7; 2; 5; 4; 3; 0; 1|]
[|6; 4; 7; 8; 5; 0; 3; 2; 1|]
リアル: 00:00:00.339、CPU: 00:00:00.369、GC 全般0: 10, 全般1: 3, 全般2: 0
val it: unit = ()

> Eight.solver2 [|8;6;7;2;5;4;3;0;1|];;
----- 21 -----
----- 22 -----
----- 23 -----
----- 24 -----
----- 25 -----
----- 26 -----
----- 27 -----
----- 28 -----
----- 29 -----
----- 30 -----
----- 31 -----
[5; 6; 8; 2; 3; 5; 1; 4; 7; 8; 6; 3; 5; 1; 4; 7; 8; 6; 3; 5; 1; 4; 7; 8; 6; 3; 2;
 1; 4; 7; 8]
[5; 6; 7; 4; 6; 2; 3; 5; 1; 6; 2; 3; 8; 7; 4; 2; 3; 1; 5; 8; 7; 4; 1; 5; 8; 7; 4;
 1; 2; 3; 6]

  ・・・省略・・・

[1; 4; 5; 2; 3; 1; 4; 5; 7; 6; 2; 3; 8; 2; 3; 8; 1; 4; 8; 7; 5; 8; 7; 5; 6; 3; 2;
 1; 4; 7; 8]
[1; 4; 5; 2; 3; 1; 4; 5; 7; 6; 2; 3; 8; 2; 3; 8; 1; 4; 5; 7; 8; 5; 7; 8; 6; 3; 2;
 1; 4; 7; 8]
リアル: 00:00:00.100、CPU: 00:00:00.090、GC 全般0: 7, 全般1: 1, 全般2: 0
val it: unit = ()