Collatz問題 #23

C

以前に、ある範囲内のそれぞれの値が何回の操作で1に至るかについて、値と操作回数の散布図を作成した。
今回は初期値からどのような経路を通って1に到達するのかをグラフではなく値と操作回数の平面上でプロットしてみる。

collatz13.c
#include <stdio.h>
#include <gmp.h>

void collatz_path(mpz_t x);

void collatz_path(mpz_t x)
{
    mpz_t n, m;
    mpz_init_set(n, x);
    mpz_init_set_ui(m, 1UL);
    while (1) {
        gmp_printf("%Zu %Zu ", n, m);
        if (mpz_odd_p(n)) {
            mpz_mul_ui(n, n, 3UL);
            mpz_add_ui(n, n, 1UL);
            gmp_printf("1\n%Zu %Zu 1\n", n, m);
        } else {
            mpz_tdiv_q_2exp(n, n, 1);
            gmp_printf("0\n%Zu %Zu 0\n", n, m);
        }
        if (mpz_cmp_ui(n, 1UL) <= 0) break;
        gmp_printf("\n");
        mpz_add_ui(m, m, 1UL);
    }
    mpz_clear(n);
}

int main(int argc, char *argv[])
{
    mpz_t x;
    if (argc == 2) {
        if (mpz_init_set_str(x, argv[1], 10) == 0 && mpz_sgn(x) > 0) {
            if (mpz_cmp_ui(x, 1UL) == 0) {
                gmp_fprintf(stderr, "warning: init value is already 1\n");
            } else {
                collatz_path(x);
            }
        } else {
            gmp_fprintf(stderr, "error: init value \"%s\" is invalid\n", argv[1]);
        }
        mpz_clear(x); /* variable x is initialized even if mpz_init_set_str causes error */
    } else {
        gmp_fprintf(stderr, "usage: %s <init value>\n", argv[0][0] != '\0' ? argv[0] : "collatz13");
    }
    return 0;
}

プログラムは引数に初期値を与えると、一回の操作ごとに操作前後の値と初期値からの操作回数とを出力する。
図示にgnuplotを利用するための仕掛けとして出力データの形式は冗長性のあるものとなっている。

$ ./collatz13 6 > p6.dat
$ cat p6.dat
6 1 0
3 1 0

3 2 1
10 2 1

10 3 0
5 3 0

5 4 1
16 4 1

16 5 0
8 5 0

8 6 0
4 6 0

4 7 0
2 7 0

2 8 0
1 8 0

1回の操作に関するデータごとに2行で出力し、データ間に空行を挿入する。
データの第1フィールドは、1行目が操作前の値、2行目が操作後の値である。
第2フィールドと第3フィールドは1行目と2行目で同じものとなっており、
第2フィールドは、この操作が何回目の操作かを示し、
第3フィールドは、この操作が半分にする操作なら0、3倍して1を加える操作なら1になる。

gnuplot> set xlabel 'initial number'
gnuplot> set ylabel 'number of iterations'
gnuplot> unset key
gnuplot> plot [][0:10] 'p6.dat' w l


プロットは下の線から6->3, 3->10, ..., 2->1の各操作を表している。
これでは操作前後の値がどちらか分かりにくいため、
これまでのグラフと同様に半分にする操作は青、3倍して1を加える操作は赤で表示する。

gnuplot> plot [][0:10] 'p6.dat' using 1:($3 == 1 ? $2 : (1/0)) w l, '' using 1:($3 != 1 ? $2 : (1/0)) w l lt 3


第3フィールドが1のものだけと0のものだけとをそれぞれラインタイプを変えてプロットしている。
空のデータファイル名は最新のデータファイル、つまりこの場合ならp6.datを再使用することを意味する。

さらに今までの散布図に合わせるために、
縦軸は各操作前の値が何回で1に至るかの回数の表示に、
横軸は対数表示に、それぞれ変更する。

gnuplot> set xlabel 'initial number'
gnuplot> set ylabel 'number of iterations'
gnuplot> unset key
gnuplot> set logscale x
gnuplot> m = 8
gnuplot> plot [][0:10] 'p6.dat' using 1:($3 == 1 ? (m - $2 + 1) : (1/0)) w l, '' using 1:($3 != 1 ? (m - $2 + 1) : (1/0)) w l lt 3


ユーザ定義変数mは最大操作回数、つまり最後の操作データの操作回数を与える。
プロットは上の線から6->3, 3->10, ..., 2->1の各操作を表している。
半分にする操作に対応する線は当然同じ長さとなり、
2回ある3倍して1を加える操作もほぼ同じ長さに見える。