今日の料理

f:id:nisemono_san:20160510104739j:plain

コラッツの問題

コラッツの予問題というとピンとこないかもしれないが、別名「3n+1問題」と言われたら、何処かで聞いたことがある人もいるんじゃないんだろうか。内容自体は簡単で、「任意の整数が与えられたとき、その整数が偶数ならば2で割り、奇数なら3倍に1を足す。そうしてできた数を先のルールに従って、同じように数を作る。こりを繰り返し、1になったら終了する」というもので、これが停止するかどうか、という話である。。

どんなアルゴリズムであっても、それが停止することが保証されていなければ、実際に使うときに非常に困ったことになる。このコラッツの問題が直接なんの役に立つのかは置いておくとして、このような単純なアルゴリズムでさえ、停止するかどうかの証明は難しいらしいのだが、この説明を聞いたら、直感的には「これは何となく止まりそうだ」と思うんじゃないかと思う。実際、確かめられているところによると、 ${2^{60}}$ までは停止することが確認されている。

で、先の問題は非常に面白い挙動をすることが確認されている。例えば、上記にあわせて「7」を変化させた場合、下のようなステップを取る。

    7 -> 22 -> 11 -> 34 -> 17 
-> 52 -> 26 -> 13 -> 40 -> 20 
-> 10 ->  5 -> 16 ->  8 ->  4 
-> 2 -> 1

これくらいであるならば、だいたいそんなものだろうと思うだろう。しかし、このアルゴリズムの面白いところは、不定期にステップ数が爆発する初期値が存在していることだ。例えば、「27」だと、次のような変節を辿る。

   27 -> 82 -> 41 -> 124 -> 62
-> 31 -> 94 -> 47 -> 142 -> 71 
-> 214 -> 107 -> 322 -> 161 -> 484
-> 242 -> 121 -> 364 -> 182 -> 91
-> 274 -> 137 -> 412 -> 206 -> 103
-> 310 -> 155 -> 466 -> 233 -> 700
-> 350 -> 175 -> 526 -> 263 -> 790
-> 395 -> 1186 -> 593 -> 1780 -> 890
-> 445 -> 1336 -> 668 -> 334 -> 167
-> 502 -> 251 -> 754 -> 377 -> 1132
-> 566 -> 283 -> 850 -> 425 -> 1276
-> 638 -> 319 -> 958 -> 479 -> 1438
-> 719 -> 2158 -> 1079 -> 3238 -> 1619
-> 4858 -> 2429 -> 7288 -> 3644 -> 1822
-> 911 -> 2734 -> 1367 -> 4102 -> 2051
-> 6154 -> 3077 -> 9232 -> 4616 -> 2308
-> 1154 -> 577 -> 1732 -> 866 -> 433 
-> 1300 -> 650 -> 325 -> 976 -> 488 
-> 244 -> 122 -> 61 -> 184 -> 92 
-> 46 -> 23 -> 70 -> 35 -> 106
 -> 53 -> 160 -> 80 -> 40 -> 20 
-> 10 -> 5 -> 16 -> 8 -> 4 -> 2 -> 1

だいたい、総ステップ数111と尋常では無い数となる。このように、コラッツの問題が生みだす数列は、面白い問題を孕んでいると言える。

なぜか最初にHaskellで書く

このあたりのステップ数を調べることに興味が出てきたので、そういうスクリプトを書こうと思って、最初におもいついたのがHaskellだった。というのは、自分の経験上、割と暴力的なコードを書いたとしても、割となんとかしてくれるのがHaskellだからだ。テストコードとしては、下のような感じ。

import Data.List
import Control.Parallel.Strategies

collatz :: Int -> (Int, Int)
collatz n = collatz_prosess n 0 n

collatz_prosess :: Int -> Int -> Int -> (Int, Int)
collatz_prosess n t f
  | n == 1    = (f, t)
  | even n    = collatz_prosess even_n (t + 1) f
  | otherwise = collatz_prosess odd_n  (t + 1) f
    where
      odd_n = n * 3 + 1
      even_n  = n `div` 2

testrange :: Int
testrange = 10000000

testdata :: [(Int, Int)]
testdata = parMap rpar collatz [1..testrange]

collatz_times_compare :: (Int, Int) -> (Int, Int) -> Ordering
collatz_times_compare (_, x) (_, y) = compare x y

collatz_times_rank :: [(Int, Int)]
collatz_times_rank = sortBy collatz_times_compare $
                     filter (\(_, i) -> i > 450)
                     $ testdata

main :: IO ()
main = do
  {- forM と mapM_ の違いを抑える -}
  mapM_ print collatz_times_rank
  return ()

コードを書いていて、もうすこし高速な処理をしたいという欲がでたので、並列処理をしようとしたような断片が残っているのだけれども、それはおいおい解説できればと思う。

コラッツの問題自体は、ステップ数は不明だけれども、しかし数自体は御互いに依存している関係ではないため、非常に並列処理を練習するのには向いている。なので、機会があれば、このHaskellのコードも、もうすこし高速化できるようにしたいと思う。

そこでGolangの出番

じゃあ並列処理に向いているプログラミング言語っていったら何かな、と言われて思いうかんだのがGolangだった。そういえば、最近Golangを採用するプロダクトも増えているしなー、最近書いてないし、申しわけないな、いというよこしまな思いで、Golangに書きなおしたりしていた。

全体のソース自体はgistのコードを読んでもらうとして、ここでは個人的に「なるほどなー」と思ったポイントをメモしておく。Golangが得意な人だと、いまいち要領のえないコードだと思うので、そこは微笑ましく指摘して頂ければと思う。

並列化の諸事情(特にsyncに関して)

Golangにおいては、並列化をするアプローチとして、二つの方法がある。それは、syncライブラリにあるWaitGroupがある。ただ、これを素朴に使うとすると、罠があって:

type CollatzResults []CollatzResult

// ...

    var cs CollatzResults  
    var wg sync.WaitGroup
    for i := 1; context.Number > i; i++ {
        wg.Add(1)
        go func(j int) {
            defer wg.Done()
            var c CollatzResult = collatz(j, !*context.Simple)
            if c.Step > *context.FilterMax {
                cs = append(cs, c)
            }
        }(i)
    }
    wg.Wait()

みたいに、なんらかの配列に結果をストックしていきたいよねー、みたいに考えていると、実行するたびに結果が変わって大変なことになる。別段、WaitGroupで処理漏れがおきているわけではなく、実行のコンテキストによって、csの配列の値が違うことが原因になる。要するに、この関数が実行されたcsの配列の値が実行されているため、順番によって変動するということのようなのである。

だったらチャネルを使えばいいじゃない、というのも一つの提案ではあるのだが、そこまで大層なものではないし、今回の場合、大量に関数を発行するので、チャネルを圧迫して破滅する。問題はある処理が終了するまで、配列に追加する処理に限っては、処理を一旦待ちうけるようにするといいはずだ。

そのときに使われるライブラリというのが、どうやらMutexというもので、こいつを使うことによって、うまく待ちうけて処理をすることが可能になる。

   var wg sync.WaitGroup
    var mt sync.Mutex

    for i := 1; context.Number > i; i++ {
        wg.Add(1)
        go func(j int) {
            defer wg.Done()
            var c CollatzResult = collatz(j, !*context.Simple)

            if c.Step > *context.FilterMax {
                mt.Lock()
                cs = append(cs, c)
                mt.Unlock()
            }
        }(i)
    }
    wg.Wait()

ちょっとした並列処理の関数を作るときはクロージャで

わざわざ言う必要もないのかもしれないけれど、上のゴルーチンを走らせている部分において、一度iの値を渡している。なんでこんなことをわざわざしているのかというと、関数の内部に対して、実行されたコンテキストを保存するためだ。

実験してみるとわかることだが、このようにクロージャを作らない場合、関数を発行しきったあとのiを参照することになってしまうため、意図しない結果が出てきてしまう。なので、こういうクロージャを作っておく必要がある。

小分けしてチャネルに投げこむ

とはいえ、ここまでやっても、正直速くなった気がしなかったので、チャネルで小分けに投げこんでみることにする。

func collatzGenerate(context CommandType) {
    const IncrementNumber = 10000
    var cs CollatzResults
    var current_n int
    var next_n int = IncrementNumber
    var ch chan CollatzResult = make(chan CollatzResult)

    gocollatz := func(j int) {
        ch <- collatz(j, !*context.Simple)
    }

    for next_n < context.Number {
        for i := 1 + current_n; next_n > i; i++ {
            go gocollatz(i)
        }

        for i := 1 + current_n; next_n > i; i++ {
            c := <-ch
            if c.Step > *context.FilterMax {
                cs = append(cs, c)
            }
        }
        current_n += IncrementNumber
        next_n += IncrementNumber
    }

    var reminder_ap int = current_n

    for current_n < context.Number {
        go gocollatz(current_n)
        current_n++
    }

    for reminder_ap < context.Number {
        c := <-ch
        if c.Step > *context.FilterMax {
z           cs = append(cs, c)
        }
        reminder_ap++
    }
}

こういう感じで小出しに処理すると上手くいって、timeではかると、3分以上かかっていたものが、2分を切るようになっていた。とはいえ、直感的にはチャネルのほうがコストが高そうな気がするので、なんでちょっと速くなっているのか不安なので、機会があれば調べたいところ。

flagライブラリ便利

で、ついでなので、オプションとかを渡せると便利だよね、ということなのでflagライブラリを使ったりしていた。

type CommandType struct {
    Simple    *bool
    Sort      *bool
    Increment *bool
    FilterMax *int
    Number    int
}

func initialize() (n CommandType) {
    // flag
    n.Simple = flag.Bool(
        "s", false, "Simple Output. Number and Step only")
    n.Increment = flag.Bool(
        "i", false, "Increment Output")
    n.FilterMax = flag.Int(
        "fmax", 0, "Under Step, not show")
    n.Sort = flag.Bool(
        "sort", false, "Sort by Step")

    // Get Command Line
    flag.Parse()
    f := flag.Args()

    if flag.NArg() == 0 {
        n.Number = 21
    } else {
        n.Number, _ = strconv.Atoi(f[0])
    }

    if n.Number < 1 {
        fmt.Println("under zero number cannot end")
        os.Exit(1)
    }
    return
}

微妙な英語は許して頂くとして、単純なオブションを渡したければflag.Boolを、なんらかの値と一緒に渡したい(上記なら-fmax=1000みたいな感じ)ならflag.Intを使う。また、flag.Parse()のときに、実際に渡されたオプションの解析がおこなわれる。ちなみに、オプションとは関係のない引数は、flag.Args()に格納されるようになっているので便利だ。

このライブラリについてはこのあたりが便利だった。

ソートするためのsortライブラリ

あともう一つのポイントとしては、今回構造体をソートする必要があって、そのあたりをどういう風に実現するといいのかしら、と調べていたらsortというライブラリがあるらしい。具体的には、下のようなインターフェイスがあれば、sort.Sortが使えるようになる。

type CollatzResults []CollatzResult

func (cs CollatzResults) Len() int {
    return len(cs)
}

func (cs CollatzResults) Less(i, j int) bool {
    return cs[i].Step < cs[j].Step
}

func (cs CollatzResults) Swap(i, j int) {
    cs[i], cs[j] = cs[j], cs[i]
}

このあたりについてはここに詳しい。