裏 RjpWiki

Julia ときどき R, Python によるコンピュータプログラム,コンピュータ・サイエンス,統計学

小難しい方法でフィボナッチ数列を求めるメリットは?

2015年05月01日 | ブログラミング

フィボナッチ数列の各項には,行列を用いて以下のような関係がある。

この行列を活用して,数同士の足し算や掛け算の回数を 200 回以内に抑えてフィボナッチ数列の第 1000 項の下 5 桁を求めよ。
解答としては「フィボナッチ数列の第 1000 項の下 5 桁」と「計算中の数同士の足し算や掛け算の回数」を答えよ。
ただし、数を下 5 桁にする計算は足し算や掛け算として含めない。

行列を t とすると,t の n 乗 の 2 行目の 1 列目の要素がフィボナッチ数列の n 項目になる。

> a = t = matrix(c(1,1,1,0), 2)
> t
     [,1] [,2]
[1,]    1    1
[2,]    1    0

> for (i in 2:10) a = a %*% t
> a
     [,1] [,2]
[1,]   89   55
[2,]   55   34

> library(gmp)
> fibnum(10)
Big Integer ('bigz') :
[1] 55

1000 項目は t の 1000 乗を求めればよいが,t の 512, 256, 128, 64, 32, 8 乗 の積を求めるのが効率がよい。

なお,フィボナッチ数列の末尾 5 桁を求めればよいので,行列の積を求めるごとに,行列の要素を 100000 で割った余りを演算結果として残せばよい。

t = matrix(c(1,1,1,0), 2)
t2 = t %*% t %% 100000
t4 = t2 %*% t2 %% 100000
t8 = t4 %*% t4 %% 100000
t16 = t8 %*% t8 %% 100000
t32 = t16 %*% t16 %% 100000
t64 = t32 %*% t32 %% 100000
t128 = t64 %*% t64 %% 100000
t256 = t128 %*% t128 %% 100000
t512 = t256 %*% t256 %% 100000
a = t512 %*% t256 %% 100000
a = a %*% t128 %% 100000
a = a %*% t64 %% 100000
a = a %*% t32 %% 100000
a = a %*% t8 %% 100000

> a
      [,1]  [,2]
[1,]  3501 28875
[2,] 28875 74626

1000 項目の末尾 5 桁は a[2,1] で 28875 である。

単純に 3,4,5,..., 1000 項を求めるプログラムは以下のようになる。足し算は 998 回だ。

b = integer(1000)
b[1] = b[2] = 1
for (i in 3:1000) b[i] = (b[i-2]+b[i-1]) %% 100000
b[1000]

> b[1000]
[1] 28875

同じ答えが得られることが確かめられた。

2×2 行列の積は,結果の各要素ごとに 2 回の掛け算と,1 回の足し算が必要である。上記プログラムは 14 回の行列積を行う。
よって,全ての計算において 14 * (2 * 4) = 112 回の掛け算と,14 * (1 * 4) = 56 回の足し算が行われる。

998 回の足し算と
112 回の掛け算と 56 回の足し算
どちらがコストが高いか。

> fib1 = function() {
+ t = matrix(c(1,1,1,0), 2)
+ t2 = t %*% t %% 100000
+ t4 = t2 %*% t2 %% 100000
+ t8 = t4 %*% t4 %% 100000
+ t16 = t8 %*% t8 %% 100000
+ t32 = t16 %*% t16 %% 100000
+ t64 = t32 %*% t32 %% 100000
+ t128 = t64 %*% t64 %% 100000
+ t256 = t128 %*% t128 %% 100000
+ t512 = t256 %*% t256 %% 100000
+ a = t512 %*% t256 %% 100000
+ a = a %*% t128 %% 100000
+ a = a %*% t64 %% 100000
+ a = a %*% t32 %% 100000
+ a = a %*% t8 %% 100000
+ }
>
> fib2 = function() {
+ b = integer(1000)
+ b[1] = b[2] = 1
+ for (i in 3:1000) b[i] = (b[i-2]+b[i-1]) %% 100000
+ b[1000]
+ }
>
> library(rbenchmark)
> benchmark(fib1, fib2, replications=100000)
  test replications elapsed relative user.self sys.self user.child sys.child
1 fib1       100000   0.417    1.185     0.364    0.040          0         0
2 fib2       100000   0.352    1.000     0.348    0.007          0         0

ほとんど差はない,わざわざ複雑な計算方法を採る意味がないということになった。(最初から分かっている?)

コメント
  • X
  • Facebookでシェアする
  • はてなブックマークに追加する
  • LINEでシェアする

PVアクセスランキング にほんブログ村

PVアクセスランキング にほんブログ村