Collatz 予想

キミならどう書く 2.0 - ROUND 2 - が発表されました。前回の問題「100までの整数から素数を列挙せよ」にもPerl で答えたので、今回も Perl で自分なりに挑戦してみます。(^_^)
ちょっと長いですが、以下問題の引用です。

Collatz予想の収束までのステップ数の最大値を求めよ.
Collatz予想とは,1以上の自然数 n に対して,次の関数 f(n) が必ず1を返すものとする.

f(n) = 1 , n = 1 のとき
= f (n/2) , n が偶数のとき
= f (3*n + 1) , n が奇数のとき

たとえば,f(3)が 1 になるまでの経緯は,

f(3) -> f(10) -> f(5) -> f(16) -> f(8) -> f(4) -> f(2) -> f(1) -> 1

となり, f(n) の呼出は8回である.これをステップ数と呼ぶことにする.
ここで f(n) が 1 になるまでのステップ数を g(n) とする.つまり,g(3) = 8 である.
このとき

h(n) = k, 1 ≦ k ≦ n ∧ g(k) = max (g(1),g(2),…,g(n))

について h(100) を求めよ.

また,余力のある者は大きな n についても h(n) を求めよ.

なかなか難しそうですね。早く解くには工夫が必要そうです。


というわけで、今回も他の方の回答は見ずに、自分なりに書いてみました。Perl というよりも、C言語っぽい書き方になってしまいましたが……。(^_^;)

use strict;

sub collatz($) {
  my $n = shift;
  my $step = 1;

  while ($n > 1) {
    $n = ($n & 1) ? (3 * $n + 1) : ($n >> 1);
    $step++;
  }

  return $step;
}

print collatz(3), "\n";
print collatz(100), "\n";

実行結果は以下の通りです。

8
26

つまり、問題への答え(h(100) = collatz(100))は「26」となります。追記:問題を読み違えていました。(^_^;) おっちょこちょいですいません……。正しい回答はこのエントリの最後に追記いたしました。


ちなみに、途中の様子を出力させると、このようになりました。

f(3)->f(10)->f(5)->f(16)->f(8)->f(4)->f(2)->f(1)->1
8
f(100)->f(50)->f(25)->f(76)->f(38)->f(19)->f(58)->f(29)->f(88)->f(44)->f(22)->f(
11)->f(34)->f(17)->f(52)->f(26)->f(13)->f(40)->f(20)->f(10)->f(5)->f(16)->f(8)->
f(4)->f(2)->f(1)->1
26

さらに、大きな数でもほぼ一瞬で解く事ができました。

print collatz(1000000), "\n";
print collatz(100000000), "\n";
153
108


以下、スクリプトの解説です。といっても、難しいのはこの式ぐらいでしょうか。

$n = ($n & 1) ? (3 * $n + 1) : ($n >> 1);

「A ? B : C」は三項演算子と呼ばれる「?:」を利用しています。「A」が真の時は「B」を、偽の時は「C」を返す演算子です。本当は、Perl演算子の優先順位から言えば、カッコは必要ないのですが、少しでもわかりやすくする為につけています。
「$n & 1」の「&」はビット演算子 AND ですね。この式は $n の1ビット目が 1 ならば真を、0 ならば偽を返すわけですが、実はこの「1ビット目」を見るだけで「偶数か奇数か」が判断できてしまうのです。詳しくは長くなるので省略しますが、2進数について調べてみると良いでしょう。
コンピュータはビット単位での計算がとても得意なので、この方法は 2 で割って余りを出すよりも断然速く偶数か奇数かを求める事ができます。
「3 * $n + 1」は、問題にある通り n が奇数の時の式ですね。
そして「$n >> 1」ですが、これは $n が正の数の時に限り「$n / 2」と同じ意味になります。つまり、2 で割っているのです。「>>」は「ビットシフト演算子」と呼ばれる演算子で、この場合「$n を2進数で表して各桁を1つ右にずらす」という操作になります。これもやはり2進数の性質に関係しています。
やはりビット単位での計算なので、普通に 2 で除算するよりも速くなります。


ビット単位での計算は、普段は目に見えませんが、コンピュータの中では頻繁に行なわれています。なぜなら、コンピュータ内部では全てのデータが2進数(0と1)で表現されているからです。
プログラマが、2進数の性質について詳しく知っていれば、それだけ効率よくコンピュータに指示を与える事ができ、プログラムの高速化につながります。ただし、こういった「コンピュータが理解しやすい書き方」は、プログラマにとっては「読みづらい」という事も、覚えておかなければいけませんね。アセンブラ言語と、Perlのような高級言語、どちらが人間にとって読みやすいかは明白です。
「高速化」をとるか「可読性」(読みやすさ)をとるか、このバランスの取り方もプログラマの腕の見せ所ですね。(^_^)

追記(2006-07-13)

問題を読み間違えていました。とんでもない馬鹿な間違いを……。ですが、失敗は成功の母という事で、少し変えて再度挑戦します。
つまり、問題の

h(n) = k, 1 ≦ k ≦ n ∧ g(k) = max (g(1),g(2),…,g(n))

の部分を読み落として、g(100) を求めていたわけです。そりゃいくらなんでも簡単すぎますね……。おっちょこちょいです。(^_^;)
というわけで、再度回答を。

use strict;

sub g($) {
  my $n = shift;
  my $step = 1;

  while ($n > 1) {
    $n = ($n & 1) ? (3 * $n + 1) : ($n >> 1);
    $step++;
  }

  return $step;
}

sub h($) {
  my $n = shift;
  my $max = 0;

  for my $i (1..$n) {
    my $step = g($i);
    $max = $step if $step > $max;
  }

  return $max;
}

print h(100), "\n";

実行結果は以下の通りです。

119

ああ、恥ずかしい。(^_^;)
でも、良い経験になりました。教訓:「問題はよく読みましょう」。

追記の追記

あのDan Kogaiさんから言及を頂きました。ありがとうございます。(^_^)
g の中で $n を変化させていく過程で、$n が大きい数になってしまうとオーバーフロー(桁あふれ)するようです。奇数の時は3倍しているので、そのせいですね。use bignum をすれば対策できますが、Danさんによると「低速ゆえお薦めしない」との事。勉強になります。
そもそも上のコードでは h(1000000) を求めるのに、とてつもなく時間がかかってしまうので、さらなる高速化が必要です。うーん……。