並列素因数分解 －Elixir見習いの実習ー

はじめに

Elixirに習熟するために並列計算による素因数分解に取り組みました。考えたことなどを書き留めます。

アルゴリズム

素朴な方法です。初等数学で習う試し割の方法を区間区分する方法です。

自然数ｎを素因数分解するのに、2,3,5,7と順に試し割していく方法です。最初だけ2で試し割しますが、2以後は素数は奇数しかありませんので、以後は奇数だけで試し割をします。逐次計算ですとnの平方根まで試し割をすれば素因数が得られます。

これを区間を分割する方法を考えました。例えば10^11を素因数分解する場合を考えます。平方根は100,000です。素因数分解する場合には次のように10区間に分割します。

2..10,000
10,001..20,000
...
90,001..100,000

それぞれにつきプロセスを生成しその区間内で試し割をします。そうしますと、最も小さな数の区間ではそのままでいいのですが、より大きな数の区間では合成数が含まれてしまいます。それらの合成数は下位の区間の素因数で割りきれますから、篩を使って合成数を除去します。

10のプロセスからは順不同で因数が返ってきます。これを結合し昇順にソートします。そしてもっとも小さい素因数で残りの要素が割りきれるかどうか？で合成数を除去し、素因数だけを残します。

素因数分解する数が小さな場合には区間区分すると却って非効率です。10,000より小さい場合には逐次計算の方法に拠っています。

コード

normal/1 比較用の逐次計算による素因数分解
factorization/1 並列計算による素因数分解
factors1/5　試し割区間を10分割し10のプロセスを生成する。
factors2/2　10のプロセスから因数の結果をリストで受け取る。
compress/1　[2,2,2,5,5] のような素因数リストをべき乗の[[2,3][5,2]]に変換する。
sieve/1　　　リストの先頭の要素でそれ以後の数のうち割りきれるものを除去する。
sieve1/1　　sieve/1 の補助関数
adjust/2 その数自身が素数だった場合に調整します。
part　　　　ワーカープロセス。区間などのデータをメッセージで受け取りfactor/4で計算してメッ　　　　　　セージで返信する。
factor/4　　区間内の因数を求める。
time/1　　　実行時間計測用のマクロ

defmodule M do
  def normal(n) do
    W.factor([],n,2,trunc(:math.sqrt(n)))
    |> compress
  end

  def factorization(n) when n < 10000 do
    W.factor([],n,2,trunc(:math.sqrt(n)))
    |> compress
  end
  def factorization(n) do
    limit = trunc(:math.sqrt(n))
    factors1(10,n,2,div(limit,10)+1,limit)
    factors2(10,[])
    |> Enum.sort |> sieve |> adjust(n)|> compress |> Enum.reverse
  end

  def compress(ls) do
    Enum.chunk_by(ls,fn(n) -> n end)
    |> Enum.map(fn(x) -> [hd(x),length(x)] end)
    |> Enum.reverse
  end

  def factors1(1,n,start,_,limit) do
    pid = spawn(W,:part,[])
    send pid, {self(),{n,start,limit}}
  end
  def factors1(c,n,start,add,limit) do
    stop = start + add
    pid = spawn(W,:part,[])
    send pid, {self(),{n,start,stop}}
    factors1(c-1,n,stop+1,add,limit)
  end

  def factors2(0,ls) do ls end
  def factors2(n,ls) do
    receive do
      {:answer, msg } ->
        factors2(n-1,msg++ls)
    end
  end

  def sieve([]) do [] end
  def sieve([l|ls]) do
    [l] ++ sieve(sieve1(ls,l,[]))
  end

  def sieve1([],_,res) do Enum.reverse(res) end
  def sieve1([l|ls],x,res) do
    if rem(l,x) == 0 && x != l do
      sieve1(ls,x,res)
    else
      sieve1(ls,x,[l]++res)
    end
  end

  def adjust(l,x) do
    if Enum.member?(l,x) do
      [x]
    else
      l
    end
  end

end

defmodule W do
  def part do
    receive do
      {sender,{n,start,stop}} -> send sender,{:answer, W.factor([],n,start,stop) }
    end
  end

  def factor(p,n,x,limit) do
    cond do
      n == 1 -> p
      x > limit -> [n|p]
      rem(n,x) == 0 -> factor([x|p],div(n,x),x,limit)
      x == 2 -> factor(p,n,3,limit)
      rem(x,2) == 0 -> factor(p,n,x+1,limit)
      true -> factor(p,n,x+2,limit)
    end
  end

end

defmodule T do
  defmacro time(exp) do
    quote do
    {time, dict} = :timer.tc(fn() -> unquote(exp) end)
    IO.inspect "time: #{time} micro second"
    IO.inspect "-------------"
    dict
    end
  end
end

実行結果

でたらめな数

でたらめな数を与えてみました。

iex(50)> T.time(M.factorization(18237192731987123))
"time: 609000 micro second"
"-------------"
[[67, 1], [211, 1], [227, 1], [5682963577, 1]]
iex(51)> T.time(M.normal(18237192731987123))
"time: 2422004 micro second"
"-------------"
[[67, 1], [211, 1], [227, 1], [5682963577, 1]]
iex(52)>

およそ4倍、高速になっています。マシンはicore5です。

大きな素数

大きな素数を与えてみました。　

iex(12)> T.time(M.factorization(1111111111111111111))
"time: 13672000 micro second"
"-------------"
[[1111111111111111111, 1]]
iex(13)> T.time(M.normal(1111111111111111111))
"time: 51469000 micro second"
"-------------"
[[1111111111111111111, 1]]
iex(14)>

およそ4倍程度です。

大きな素因数を含まない数

小さな素因数で構成される数を与えてみました。

iex(14)> T.time(M.factorization(13123123121232))
"time: 32000 micro second"
"-------------"
[[2, 4], [3, 2], [17, 1], [8311, 1], [645019, 2]]
iex(15)> T.time(M.normal(13123123121232))
"time: 0 micro second"
"-------------"
[[2, 4], [3, 2], [17, 1], [8311, 1], [645019, 1]]
iex(16)>

並列動作の方が遅くなっています。

暗算でもできそうな数

暗算でもできそうな簡単な数だと、おもいっきり遅いです。

iex(8)>  T.time(M.factorization(10000000000000000))
"time: 453000 micro second"
"-------------"
[[2, 16], [5, 16]]
iex(9)>  T.time(M.normal(10000000000000000))
"time: 0 micro second"
"-------------"
[[2, 16], [5, 16]]
iex(10)>

大きな素因数を２つもつ数

iex(4)> T.time(M.factorization(18446744073709551609))
"time: 54187000 micro second"
"-------------"
[[3, 3], [818923289, 1], [2502845209, 1]]
iex(5)> T.time(M.normal(18446744073709551609))
"time: 68750000 micro second"
"-------------"
[[3, 2], [818923289, 1], [2502845209, 1]]
iex(6)>

あまり差がつきませんでした。がっかり。

考察

大きな素因数を持つ数、あるいはその数が素数である場合には並列動作の方が高速でした。しかし、大きな素因数を持たない場合にはかえって並列動作の方が遅くなっています。これは大きな数の区間で生ずる合成数の因数を篩で除去するコストの方が大きいためだろうと思います。

終わりに

素朴な方法で並列動作による素因数分解を試してみました。ググってみますと本格的な並列素因数分解には高等数学を使うもののようです。http://www.ntt.co.jp/news2010/1001/100108a.html
アルゴリズムの上での間違いなど助言ありましたら、コメントをお願いいたします。