Python勉強記(11) アッカーマン関数

はじめに

これまで、少しずつPython言語に必要なパーツを学習してきました。

次に、Pythonプログラムがどのように組み立てられるかを俯瞰しつつ、プログラムが組み立てられるように、少しずつレベルを上げていきます。

ライブラリを駆使して、機械学習・統計処理という段階までには、もう少しかかります。ですが、部品はそろったのですから、これを生かして何かできないか？を考えていきます。

今回は、情報系の学部や仕事をした人ならば一度はみたことがある、有名な計算関数をいじってみたいと思います。

アッカーマン関数とは？

m, n を負でない整数とするとき、以下の式で再帰的（帰納的）に定義されるA(m, n)をアッカーマン関数といいます。

$$
A(0,n) = n + 1 \\
A(m,0) = A(m -1, 1) \\
A(m,n ) = A(m - 1,A(m, n - 1))
$$

A(m,n)の性質

アッカーマン関数 - Wikipedia

m=3までは、Aの一般項は簡単な初等関数で表せます。

$$
A(1,n) = n + 2 = 2 + (n + 3) -3 \\
A(2,n) = 2n + 3  = 2(n+3) - 3 \\
A(3,n) = 2^{n+3}  - 3 \\
A(4,n) = 2↑↑(n+3) - 3
$$

Wikipedia

m=4の「2↑↑(n+3)」は「２の累乗をn+3回繰り返す」という意味です。このような特別な演算子を用いないと表現できませんが、mが超えると爆発的に増加することはわかります。

クヌースの矢印表記

定義に従ったPythonプログラム

Pythonの関数は、再帰的に（自分自身を）呼び出せるため。定義に従って容易に計算するコードを記述できます。以下コード内ではAck(m.n）と書くことにします。
コンソールよりm, nの値をinput関数で自由に入力します。

def Ack(m, n):
    if m == 0:
        return n + 1
    elif n == 0:
        return Ack(m - 1, 1)
    else:
        return Ack(m - 1, Ack(m,  n - 1))
        
m = int(input("m="))
n = int(input("n="))
print(Ack(m, n))

小さな(m,n)では問題なく計算できます（が、アッカーマン関数の性質を先に知ると、m=3, 4あたりで計算が厳しくなってくるのが予想できるでしょう）

m=3として、nの値を0から少しずつ増やしていくと、Ack(3,8) = 2045 あたりから少しもっさりし始め、Ack(3,9) でおそらくエラーが出ます。
当方のPCのJupyter Lab 上でm=3, n=9 とすると

Ack(3,9)は上記プログラムでは計算できない。

「RecursionError: maximum recursion depth exceeded」とは一体どういうことなのでしょうか？

「maximum recurtion depth」について

A(2,2)を手計算していこうとしてみると、関数の定義式より

A(2, 2) = A(1, A(2,1))
           = A(1, A(1, A(2,0)))
　　　= A(1, A(1,2))) 
           = A(1, A(0, A(1,1))))
　　   = …

0

となり（このレベルで手計算では無理、となりますが）、Pythonプログラム上で計算中、右辺の入れ子状態を保存しておかなければいけません。この保存領域をスタックと言います。

スタックの入れ子はどんどん伸びていき（＝再帰の深さ(recursion depth)がどんどん深くなり）、アッカーマン関数では指数的に増加していくために処理系のリミットを超えます。「maximum recursion depth exceed」とはリミットを超えた状態です。

「maximum recursion depth」はPythonの標準ライブラリ関数で取得できます。

import sys
print(sys.getrecursionlimit())

実際に取得してみると、私のjupyter labでは3000と出ました（標準のpythonでは1000とのことです）

再帰をもっと深くできるか

作成したPython コードの先頭に、以下２行を挿入します。

import sys
sys.setrecursionlimit(10000)

めでたくAck(3,9)は計算できました。必要な再帰の深さは3000よりは大きいけれど10000以下だったことがわかります。

【注意】sys.setrecursionlimit(100000)などとやってはいけません。スタックがメモリを食いつぶし、PCがクラッシュする原因となります。
Google Colabで試したところ、30000で「原因不明なエラー」が発生しました。

計算を工夫する。

以下のサイトを参考にしました。
Ackermann 函数と Python で遊んだ
https://gist.github.com/metasta/47414638f6e3e7688124a30b96006e85

「手習」でこういうコード書いてしまう人すげぇ…。

サイト内のアッカーマン関数をリスト操作に置き換えた素晴らしいAckermann_listという関数を拝借・計算してみます。

Ackermann_list（上記サイト内より引用）

難なくA(3,9) = 2^12-3 =4093は計算できました、が
次にA(4,1)=65533 を計算しようとすると…帰ってきません。
プログラムが悪いのではなく、恐ろしい再帰爆発のために、計算量もとんでもないことになっていると予想できます。

そうだ時間測ろう。

上のAckermann_list関数の計算時間を測れるように、少し書き換えます。
標準ライブラリのtime.time()関数で測定時間を測ります。

import time as tm

def ackermann_list(p, q):
  a = [p, q]
  while len(a) > 1:
    m, n = a[-2:]
    if m == 0:
      a[-2:] = [n+1]
    elif n == 0:
      a[-2:] = [m-1, 1]
    else:
      a[-2:] = [m-1, m, n-1]
  return a[0]

m = int(input("m ="))
n = int(input("n ="))
tm0 = tm.time()
ack1 = ackermann_list(m, n)
tm1 = tm.time()
print(ack1, " Time=", tm1-tm0)

まずはJupyter lab上でA(3,9)を計算してみました。

A(3,9) on Jupyter Lab

約3.68秒。結構かかっていますね。
もしかしたらJupyter labのIPython なんかでやってるのが悪いのかも。

次に上記ソースをファイルack_list.py に保存し、コマンドラインでやってみます。(python ack_list.py)

A(3,9) on CPython

約1.97秒。全然比較にならないほど速い。といっても２倍弱ですが。

A(4.1)をコマンドライン版Python（CPython）で流してみます・・・

A(4,1) on CPython

やっと出た。536秒＝8.93分。結構すごい。
これ以上はPCでは計算できない（してはいけない）領域でしょう。