Collatz問題

Haskellのプログラムを参考にCollatz問題のプログラムを書いてみました。
Fortran2003のREALLOCATE機能を使うことで、配列生成子を使うことで、配列をリストの代用とし、再帰を使って自動で伸ばしてゆきます。

Collatz問題とは何かは、詳しくはここを見てください。

簡単に言えば、ある正整数nが与えられたときに、n偶数ならn/2,n奇数なら3*n+1としてn=1になるまで繰り返してゆき、これが停止するかどうかを問います。かなり大きな数まで必ず1になって停止することが知られています。

■実行結果

１０万までの自然数に対して停止するかを実行しました。最長の数列が更新される毎に、その数列を出力しました。ここに示したのは、１０万までで最長のもの。77031で初期値を含め351ステップで停止しています。

■ソースコード

オプション　/assume:realloc_lhs が必要です。
テンポラリ配列 m(:) に代入しているのは、FUNCTION 返り値の ALLOCATABLE 配列が代入操作をしないと解放されないためです。これはバグではないかという気がします。
（しかし、おぼろげな記憶では Fortran95/2003 Explained には、FUNCTION 返り値の ALLOCATABLE 配列は、代入操作で解放されるとはありましたが、その他どういう場合に解放されるのかわからないので、もう一度 Fortran95/2003 Explained なり規格書を読まないとなんともいえません。）

MODULE m_Collatz
   IMPLICIT NONE
 CONTAINS
   RECURSIVE FUNCTION iCollatz(n) RESULT(ires)
    INTEGER, INTENT(IN) :: n
    INTEGER, ALLOCATABLE :: ires(:)
    IF (n == 1) THEN
     ires = [1]
    ELSE IF ( MOD(n, 2) == 0 ) THEN ! even
     ires = [n, iCollatz(n / 2)]
    ELSE                             ! odd  
     ires = [n, iCollatz(3 * n + 1)]     
    END IF
    RETURN   
   END FUNCTION iCollatz
END MODULE m_Collatz
!===========================
PROGRAM test
   USE m_Collatz
   IMPLICIT NONE
   INTEGER :: i, k
   INTEGER, ALLOCATABLE :: m(:)
   k = 0
   DO i = 1, 10**5
    m = iCollatz(i)
    IF ( k < SIZE(m) ) THEN 
     PRINT '(2i10)', i, SIZE(m)
     PRINT '(8i9)', m
     PRINT *
     k = max(k, SIZE(m))
    END IF
   END DO 
   STOP
END PROGRAM test