理论上阿克曼函数可以优化吗?

Ξέν*_*νος 42 python algorithm python-3.x ackermann

我想知道是否可以有一个阿克曼函数的版本,其时间复杂度比标准变体更好。

\n

这不是作业,我只是好奇。我知道阿克曼函数除了作为性能基准之外没有任何实际用途,因为它具有深度递归。我知道数字增长得非常快,而且我对计算它不感兴趣。

\n

尽管我使用Python 3并且整数不会溢出,但我的时间确实有限,但我根据维基百科上的定义自己实现了它的一个版本,并计算了极小的值的输出,只是为了使确保输出正确。

\n

在此输入图像描述

\n
def A(m, n):\n    if not m:\n        return n + 1\n    return A(m - 1, A(m, n - 1)) if n else A(m - 1, 1)\n
Run Code Online (Sandbox Code Playgroud)\n

上面的代码是直接翻译图片,速度极慢,不知道如何优化,难道就不能优化了吗?

\n

我能想到的一件事是记住它,但是递归是向后运行的,每次递归调用函数时,参数都是以前没有遇到过的,每次连续的函数调用参数都会减少而不是增加,因此函数的每个返回值都需要要计算,当您第一次使用不同的参数调用函数时,记忆并没有帮助。

\n

仅当您再次使用相同的参数调用它时,记忆化才有帮助,它不会计算结果,而是检索缓存的结果,但如果您使用 (m, n) >= (4, 2 )无论如何它都会使解释器崩溃。

\n

我还根据这个答案实现了另一个版本:

\n
def ack(x, y):\n    for i in range(x, 0, -1):\n        y = ack(i, y - 1) if y else 1\n    return y + 1\n
Run Code Online (Sandbox Code Playgroud)\n

但实际上它更慢:

\n
In [2]: %timeit A(3, 4)\n1.3 ms \xc2\xb1 9.75 \xc2\xb5s per loop (mean \xc2\xb1 std. dev. of 7 runs, 1,000 loops each)\n\nIn [3]: %timeit ack(3, 4)\n2 ms \xc2\xb1 59.9 \xc2\xb5s per loop (mean \xc2\xb1 std. dev. of 7 runs, 1,000 loops each)\n
Run Code Online (Sandbox Code Playgroud)\n

理论上阿克曼函数可以优化吗?如果不是,是否可以肯定证明它的时间复杂度不能降低?

\n
\n

我刚刚测试过A(3, 9)A(4, 1)将使解释器崩溃,以及这两个函数的性能A(3, 8)

\n
In [2]: %timeit A(3, 8)\n432 ms \xc2\xb1 4.63 ms per loop (mean \xc2\xb1 std. dev. of 7 runs, 1 loop each)\n\nIn [3]: %timeit ack(3, 8)\n588 ms \xc2\xb1 10.4 ms per loop (mean \xc2\xb1 std. dev. of 7 runs, 1 loop each)\n
Run Code Online (Sandbox Code Playgroud)\n
\n

我又做了一些实验:

\n
from collections import Counter\nfrom functools import cache\n\nc = Counter()\ndef A1(m, n):\n    c[(m, n)] += 1\n    if not m:\n        return n + 1\n    return A(m - 1, A(m, n - 1)) if n else A(m - 1, 1)\n\ndef test(m, n):\n    c.clear()\n    A1(m, n)\n    return c\n
Run Code Online (Sandbox Code Playgroud)\n

这些论点确实重复了。

\n

但令人惊讶的是,缓存根本没有帮助:

\n
In [9]: %timeit Ackermann = cache(A); Ackermann(3, 4)\n1.3 ms \xc2\xb1 10.1 \xc2\xb5s per loop (mean \xc2\xb1 std. dev. of 7 runs, 1,000 loops each)\n
Run Code Online (Sandbox Code Playgroud)\n

缓存仅在再次使用相同参数调用函数时才有帮助,如下所示:

\n
In [14]: %timeit Ackermann(3, 2)\n101 ns \xc2\xb1 0.47 ns per loop (mean \xc2\xb1 std. dev. of 7 runs, 10,000,000 loops each)\n
Run Code Online (Sandbox Code Playgroud)\n

我已经用不同的参数对其进行了多次测试,并且它总是给出相同的效率提升(实际上没有)。

\n

Ste*_*ann 27

解决方案

\n

我最近根据 templatetypedef 提到的同一篇论文写了一堆解决方案。许多人使用生成器,每个 m 值一个,生成 n=0、n=1、n=2 等的值。这可能是我最喜欢的:

\n
def A_Stefan_generator_stack3(m, n):\n    def a(m):\n        if not m:\n            yield from count(1)\n        x = 1\n        for i, ai in enumerate(a(m-1)):\n            if i == x:\n                x = ai\n                yield x\n    return next(islice(a(m), n, None))\n
Run Code Online (Sandbox Code Playgroud)\n

解释

\n

考虑发电机a(m)。它产生 A(m,0)、A(m,1)、A(m,2) 等。 A(m,n) 的定义使用 A(m-1, A(m, n-1)) 。因此a(m),在其索引 n 处产生 A(m,n),计算如下:

\n
    \n
  • a(m)A(m,n-1) 由生成器本身在索引 n-1 处生成 。这正是该生成器生成的先前值 (x)。
  • \n
  • a(m-1)A(m-1, A(m, n-1)) = A(m-1, x) 由索引 x 处的生成器生成。因此a(m)生成器迭代a(m-1)生成器并获取索引 i == x 处的值。
  • \n
\n

基准

\n

以下是计算 m\xe2\x89\xa43 和 n\xe2\x89\xa417 的所有 A(m,n) 的时间,还包括 templatetypedef 的解决方案:

\n
 1325 ms  A_Stefan_row_class\n 1228 ms  A_Stefan_row_lists\n  544 ms  A_Stefan_generators\n 1363 ms  A_Stefan_paper\n  459 ms  A_Stefan_generators_2\n  866 ms  A_Stefan_m_recursion\n  704 ms  A_Stefan_function_stack\n  468 ms  A_Stefan_generator_stack\n  945 ms  A_Stefan_generator_stack2\n  582 ms  A_Stefan_generator_stack3\n  467 ms  A_Stefan_generator_stack4\n 1652 ms  A_templatetypedef\n
Run Code Online (Sandbox Code Playgroud)\n

注意:使用数学见解/公式甚至更快(更快)的解决方案是可能的,请参阅的评论pts 的答案。我故意没有这样做,因为我对编码技术感兴趣,以避免深度递归并避免重新计算。我的印象是,这也是问题/OP想要的,他们现在确认了这一点(在已删除的答案下,如果您有足够的声誉,则可见)。

\n

代码

\n
def A_Stefan_row_class(m, n):\n    class A0:\n        def __getitem__(self, n):\n            return n + 1\n    class A:\n        def __init__(self, a):\n            self.a = a\n            self.n = 0\n            self.value = a[1]\n        def __getitem__(self, n):\n            while self.n < n:\n                self.value = self.a[self.value]\n                self.n += 1\n            return self.value\n    a = A0()\n    for _ in range(m):\n        a = A(a)\n    return a[n]\n\n\nfrom collections import defaultdict\n\ndef A_Stefan_row_lists(m, n):\n    memo = defaultdict(list)\n    def a(m, n):\n        if not m:\n            return n + 1\n        if m not in memo:\n            memo[m] = [a(m-1, 1)]\n        Am = memo[m]\n        while len(Am) <= n:\n            Am.append(a(m-1, Am[-1]))\n        return Am[n]\n    return a(m, n)\n\n\nfrom itertools import count\n\ndef A_Stefan_generators(m, n):\n    a = count(1)\n    def up(a, x=1):\n        for i, ai in enumerate(a):\n            if i == x:\n                x = ai\n                yield x\n    for _ in range(m):\n        a = up(a)\n    return next(up(a, n))\n\n\ndef A_Stefan_paper(m, n):\n    next = [0] * (m + 1)\n    goal = [1] * m + [-1]\n    while True:\n        value = next[0] + 1\n        transferring = True\n        i = 0\n        while transferring:\n            if next[i] == goal[i]:\n                goal[i] = value\n            else:\n                transferring = False\n            next[i] += 1\n            i += 1\n        if next[m] == n + 1:\n            return value\n\n\ndef A_Stefan_generators_2(m, n):\n    def a0():\n        n = yield\n        while True:\n            n = yield n + 1\n    def up(a):\n        next(a)\n        a = a.send\n        i, x = -1, 1\n        n = yield\n        while True:\n            while i < n:\n                x = a(x)\n                i += 1\n            n = yield x\n    a = a0()\n    for _ in range(m):\n        a = up(a)\n    next(a)\n    return a.send(n)\n\n\ndef A_Stefan_m_recursion(m, n):\n    ix = [None] + [(-1, 1)] * m\n    def a(m, n):\n        if not m:\n            return n + 1\n        i, x = ix[m]\n        while i < n:\n            x = a(m-1, x)\n            i += 1\n        ix[m] = i, x\n        return x\n    return a(m, n)\n\n\ndef A_Stefan_function_stack(m, n):\n    def a(n):\n        return n + 1\n    for _ in range(m):\n        def a(n, a=a, ix=[-1, 1]):\n            i, x = ix\n            while i < n:\n                x = a(x)\n                i += 1\n            ix[:] = i, x\n            return x\n    return a(n)\n\n\nfrom itertools import count, islice\n\ndef A_Stefan_generator_stack(m, n):\n    a = count(1)\n    for _ in range(m):\n        a = (\n            x\n            for a, x in [(a, 1)]\n            for i, ai in enumerate(a)\n            if i == x\n            for x in [ai]\n        )\n    return next(islice(a, n, None))\n\n\nfrom itertools import count, islice\n\ndef A_Stefan_generator_stack2(m, n):\n    a = count(1)\n    def up(a):\n        i, x = 0, 1\n        while True:\n            i, x = x+1, next(islice(a, x-i, None))\n            yield x\n    for _ in range(m):\n        a = up(a)\n    return next(islice(a, n, None))\n\n\ndef A_Stefan_generator_stack3(m, n):\n    def a(m):\n        if not m:\n            yield from count(1)\n        x = 1\n        for i, ai in enumerate(a(m-1)):\n            if i == x:\n                x = ai\n                yield x\n    return next(islice(a(m), n, None))\n\n\ndef A_Stefan_generator_stack4(m, n):\n    def a(m):\n        if not m:\n            return count(1)\n        return (\n            x\n            for x in [1]\n            for i, ai in enumerate(a(m-1))\n            if i == x\n            for x in [ai]\n        )\n    return next(islice(a(m), n, None))\n\n\ndef A_templatetypedef(i, n):\n    positions = [-1] * (i + 1)\n    values = [0] + [1] * i\n    \n    while positions[i] != n:       \n        values[0]    += 1\n        positions[0] += 1\n            \n        j = 1\n        while j <= i and positions[j - 1] == values[j]:\n            values[j] = values[j - 1]\n            positions[j] += 1\n            j += 1\n\n    return values[i]\n\n\nfuncs = [\n    A_Stefan_row_class,\n    A_Stefan_row_lists,\n    A_Stefan_generators,\n    A_Stefan_paper,\n    A_Stefan_generators_2,\n    A_Stefan_m_recursion,\n    A_Stefan_function_stack,\n    A_Stefan_generator_stack,\n    A_Stefan_generator_stack2,\n    A_Stefan_generator_stack3,\n    A_Stefan_generator_stack4,\n    A_templatetypedef,\n]\n\nN = 18\nargs = (\n    [(0, n) for n in range(N)] +\n    [(1, n) for n in range(N)] +\n    [(2, n) for n in range(N)] +\n    [(3, n) for n in range(N)]\n)\n\nfrom time import time\n\ndef print(*args, print=print, file=open(\'out.txt\', \'w\')):\n    print(*args)\n    print(*args, file=file, flush=True)\n    \nexpect = none = object()\nfor _ in range(3):\n  for f in funcs:\n    t = time()\n    result = [f(m, n) for m, n in args]\n    # print(f\'{(time()-t) * 1e3 :5.1f} ms \', f.__name__)\n    print(f\'{(time()-t) * 1e3 :5.0f} ms \', f.__name__)\n    if expect is none:\n        expect = result\n    elif result != expect:\n        raise Exception(f\'{f.__name__} failed\')\n    del result\n  print()\n
Run Code Online (Sandbox Code Playgroud)\n

  • 关于如何改进这一点的思考:我们知道对于所有 n,A(1, n) = 2n + 3。也许每当请求 A(1, n) 时,您都可以通过直接返回 2n+3 来消除最后两级递归? (2认同)

Ano*_*noE 19

您在文字中问了这两个问题:

理论上阿克曼函数可以优化吗?

当然,您可以天真地实现它,然后以技术方式对其进行优化(例如记忆或存储中间值等)。但我认为这不是你问的实际问题。

如果不是,是否可以肯定证明它的时间复杂度不能降低?

优化与时间复杂度无关。“O”符号抽象了常数乘法因子。您可以有两种算法,其中一种算法计算 Ackermann 函数的时间或空间是另一种算法的百万分之一,但它们仍然具有相同的 O 复杂度。

引用维基百科的话,

阿克曼函数以威廉·阿克曼 (Wilhelm Ackermann) 的名字命名,是非原始递归的全可计算函数中最简单和最早发现的示例之一

该函数不是原始递归函数,而且,它

比任何原始递归函数增长得更快。

“原始递归”意味着您可以使用预先已知边界的循环来实现算法;即您不需要可能无限重复的while循环。当然,这是一个有点抽象的概念,再次引用维基百科

原始递归函数的重要性在于,数论(更普遍的是数学)中研究的大多数可计算函数都是原始递归函数。

是的,已经证明阿克曼不是原始递归。发现事实并非如此,可能不会让你赚到任何钱,但肯定会让你在理论计算机科学界享有盛誉。

你无法优化这种复杂性——考虑你的程序只是编写证明阿克曼确实是原始递归的不同格式;你只需要把它转换回数学/逻辑术语,瞧。事实上,这种情况已经很多年没有发生了,再加上链接中的“积极证据”的存在,告诉你,它的行为实际上正如广告中所宣传的那样。

注意:最后,所有这一切都必须从阿克曼函数可能被设计为具有此属性的角度来看待- 正如维基百科页面提到的,在它被发现或创建之前,有些人认为所有可计算函数都是原始的递归的。虽然我不知道是什么促使阿克曼先生在 1920 年代这样做,但阿克曼函数现在确实是理论计算机科学中复杂性理论的基石;一个非常有趣的故事。

  • 朴素的递归斐波那契在“O(Fib(n))”中运行。一个简单的迭代斐波那契在“O(n)”中运行,与记忆递归斐波那契相同。通过算法优化或仅记忆重复计算的子问题,可以在复杂性方面取得巨大的改进。当然,事实证明,您可以为 Ackermann 做的最好的事情并不是设计得那么快,但它可能是与最简单的实现不同的 O(f(n)) 复杂度类别。 (4认同)
  • 原始递归属性意味着如果您需要循环,您可以提前知道需要多少次迭代。一般来说,这是不可能的,退出循环的条件可能只能在执行循环时计算。 (2认同)

小智 15

v0 几乎就是你的代码:

\n
def ack(m, n):\n  if m == 0: return n + 1\n  return ack(m - 1, 1) if n == 0 else ack(m - 1, ack(m, n - 1))\n
Run Code Online (Sandbox Code Playgroud)\n

这将在 2.49 秒内计算出 ack(3, 9)。ack() 被调用 11164370 次。当然,我们可以缓存已经计算出的值,从而节省对函数的大量调用。

\n

v1 使用字典作为缓存,并且仅在结果不在缓存中时才计算:

\n
c = {}\n\ndef ack(m, n):\n  global c\n\n  if "{}-{}".format(m, n) in c: return c["{}-{}".format(m, n)]\n  else:\n    if m == 0: ret = n + 1\n    else: ret = ack(m - 1, 1) if n == 0 else ack(m - 1, ack(m, n - 1))\n\n    c["{}-{}".format(m, n)] = ret\n    return ret\n
Run Code Online (Sandbox Code Playgroud)\n

这会在 0.03 秒内计算 ack(3, 9),因此 ack(3, 9) 不再适合测量执行时间。这次ack()被调用了12294次,节省的空间是巨大的。但我们可以做得更好。从现在开始,我们使用 ack(3, 13),当前运行时间为 0.23 秒。

\n

v2 只关心 m > 0 的缓存,因为 m == 0 的情况是微不足道的。这样空间复杂度就有所降低:

\n
c = {}\n\ndef ack(m, n):\n  global c\n\n  if m == 0: return n + 1\n  else:\n    if "{}-{}".format(m, n) in c: return c["{}-{}".format(m, n)]\n    else: ret = ack(m - 1, 1) if n == 0 else ack(m - 1, ack(m, n - 1))\n\n    c["{}-{}".format(m, n)] = ret\n    return ret\n
Run Code Online (Sandbox Code Playgroud)\n

这会在 0.18 秒内完成 ack(3, 13)。但我们可以做得更好。

\n

v3 通过每次迭代仅计算一次缓存中的键来节省一些时间:

\n
c = {}\n\ndef ack(m, n):\n  global c\n\n  if m == 0: return n + 1\n  else:\n    key = "{}-{}".format(m, n)\n    if key in c: return c[key]\n    else: ret = ack(m - 1, 1) if n == 0 else ack(m - 1, ack(m, n - 1))\n\n    c[key] = ret\n    return ret\n
Run Code Online (Sandbox Code Playgroud)\n

这次运行了0.14s。午夜时分我当然不能做更多的事情,但我会再考虑一下。Ez j\xc3\xb3 mulats\xc3\xa1g, f\xc3\xa9rfi munka volt - 对于那些知道这意味着什么的人。

\n

  • 您可以使用元组作为字典键。 (5认同)

tem*_*def 12

Here\xe2\x80\x99是我对Grossman-Zeitman 算法的 Python 实现的尝试,该算法在 O(A(i, n)) 时间内迭代计算 A(i, n) 。有关该算法如何工作的描述,请检查链接的问题。

\n
def ackermann(i, n):\n    positions = [-1] * (i + 1)\n    values = [0] + [1] * i\n    \n    while positions[i] != n:       \n        values[0]    += 1\n        positions[0] += 1\n            \n        j = 1\n        while j <= i and positions[j - 1] == values[j]:\n            values[j] = values[j - 1]\n            positions[j] += 1\n            j += 1\n\n    return values[i]\n
Run Code Online (Sandbox Code Playgroud)\n

鉴于内部循环非常紧密并且 \xe2\x80\x99s 没有递归,我怀疑这可能会优于您最初发布的基本递归解决方案。然而,我还没有\xe2\x80\x99t 测试过这个实现的正确性或时间;it\xe2\x80\x99s 基于我在链接问题中编写的 Java 代码。

\n


Kel*_*ndy 11

\n

但令人惊讶的是,缓存根本没有帮助:

\n
In [9]: %timeit Ackermann = cache(A); Ackermann(3, 4)\n1.3 ms \xc2\xb1 10.1 \xc2\xb5s per loop (mean \xc2\xb1 std. dev. of 7 runs, 1,000 loops each)\n
Run Code Online (Sandbox Code Playgroud)\n
\n

那是因为你做错了。仅Ackermann记住您的,而不记住您的A。并且递归调用全部转到A.

\n

m,n = 3,8 的时间,包括正确记忆的版本B

\n
440.30 ms  A(m, n)\n431.11 ms  Ackermann = cache(A); Ackermann(m, n)\n  1.74 ms  B.cache_clear(); B(m, n)\n
Run Code Online (Sandbox Code Playgroud)\n

关于B

\n
    \n
  • 之后打印B.cache_info())以查看缓存的效果:CacheInfo(hits=1029, misses=5119, maxsize=None, currsize=5119). 因此,B有 5,119 个“真实”调用(它必须工作)和 1,029 个从缓存应答的调用。如果没有记忆,它会被调用 2,785,999 次。
  • \n
  • 对于 m,n = 3,12,大约需要 32 ms。
  • \n
  • 对于 m,n = 3,13 ,它会由于深度递归而崩溃。
  • \n
\n

代码:

\n
from timeit import repeat\nimport sys\n\nsys.setrecursionlimit(999999)\n\nsetup = '''\nfrom functools import cache\n\ndef A(m, n):\n    if not m:\n        return n + 1\n    return A(m - 1, A(m, n - 1)) if n else A(m - 1, 1)\n\n@cache\ndef B(m, n):\n    if not m:\n        return n + 1\n    return B(m - 1, B(m, n - 1)) if n else B(m - 1, 1)\n\nm, n = 3, 8\n'''\n\ncodes = [\n    'A(m, n)',\n    'Ackermann = cache(A); Ackermann(m, n)',\n    'B.cache_clear(); B(m, n)',\n]\n\nfor code in codes:\n    t = min(repeat(code, setup, number=1))\n    print(f'{t*1e3:6.2f} ms ', code)\n
Run Code Online (Sandbox Code Playgroud)\n


pts*_*pts 9

TL;DR Ackermann 函数增长如此之快,以至于对于 (m >= 4 且 n >= 3),结果无法适合 RAM。对于较小的 m 或 n 值,直接(无需递归)快速计算值非常容易。

我知道这个答案无助于优化递归函数调用,尽管如此,它为在真实的当代计算机上计算阿克曼函数的值提供了一种快速解决方案(通过分析),并且它提供了问题第一段的直接答案。

我们希望将结果存储在计算机上的无符号二进制大整数变量中。为了存储值 2 ** b,我们需要 (b + 1) 位数据和 (c1 + c2 * ceil(log2(b))) 位标头来存储长度 b 本身。(c1 是一个非负整数,c2 是一个正整数。)因此我们需要超过 b 位的 RAM 来存储 2 ** b。

具有大量 RAM 的计算机:

  • 到 2023 年,消费级 PC 的 RAM 高达 128 GiB,并且商用服务器具有 2 TiB RAM ( https://ramforyou.com/is-it-possible-for-a-computer-have-1tb-内存)。
  • 2020 年,构建了具有 6 TiB RAM 的单机架服务器(https://www.youtube.com/watch?v=uHAfTty9UWY)。
  • 2017 年,构建了具有 160 TiB RAM 的大型服务器 ( https://www.forbes.com/sites/aarontilley/2017/05/16/hpe-160-terabytes-memory/ )。
  • 所以我们可以说,在 2023 年建造一台具有 1 PiB == 1024 TiB RAM 的计算机是不可行的,并且 1 EiB == 1024 PiB == 1048576 TiB == (2 ** 63) 位在 2023 年是完全不可能的,不久的将来。
  • 目前对宇宙中原子数量的估计为 10 ** 82 == 10 * 10 ** 81 < 16 * 2 ** 270 < 2 ** 274。
  • 让我们想象一下最大的计算机。即使有更多的原子,比如说 2 ** 300 个,我们也可以在一台计算机中使用所有这些原子,并且我们可以在单个原子中存储 1 EiB 的数据,因此我们有 2 ** 363 位的内存;它仍然太小,无法存储Big值 2 ** (2 ** 363)。

让我们看看 Knuth 的向上箭头表示法 ( https://en.wikipedia.org/wiki/Knuth%27s_up-arrow_notation ) 中的哪个值小于Big

  • 2 ^ b == 2 ** b:适用于 1 <= b <= (2 ** 363 - 1)。

    2 ^ (2 ** 363) == 2 ** (2 ** 363) == 大 >= 大。

  • 2 ^^ b:适用于 1 <= b <= 5。

    2 ^^ 1 == 2;

    2 ^^ 2 == 2 ** 2 == 4;

    2 ^^ 3 == 2 ** (2 ** 2) == 16;

    2 ^^ 4 == 2 ** (2 ** (2 ** 2)) == 65536;

    2 ^^ 5 == 2 ** (2 ** (2 ** (2 ** 2))) == 2 ** 65536 == 2 ** (2 ** 16) < 大;

    2 ^^ 6 == 2 ** (2 ** (2 ** (2 ** (2 ** 2)))) == 2 ** (2 ** 65536) >= 大。

  • 2 ^^^ b:适用于 1 <= b <= 3。

    2 ^^^ 1 == 2;

    2 ^^^ 2 == 2 ^^ 2 == 4;

    2 ^^^ 3 == 2 ^^ (2 ^^ 2) == 65536;

    2 ^^^ 4 == 2 ^^ (2 ^^ (2 ^^ 2)) == 2 ^^ 65536 >= 2 ^^ 6 >= 大。

  • 2 ^^^^ b:适用于 1 <= b <= 2。

    2 ^^^^ 1 == 2;

    2 ^^^^ 2 == 2 ^^^ 2 == 4;

    2 ^^^^ 3 == 2 ^^^ (2 ^^^ 2) == 2 ^^^ 4 == 2 ^^ 65536 >= 2 ^^ 6 >= 大。

  • 2 ^^^^^ b 甚至更多箭头:它适用于 1 <= b <= 2。该值至少为 2 ^^^^ b,请参阅此处。

因此,以下是如何在 Python 中计算可行值:

def up_arrow(a, b):
  if b <= 2:
    if b < 0:
      raise ValueError
    return (1, 2, 4)[b]
  elif a == 1:
    if b >> 363:
      raise ValueError
    return 1 << b  # This may run out of memory for large b.
  elif a == 2:
    if b > 5:
      raise ValueError
    if b == 5:
      return 1 << 65536
    return (16, 65536)[b - 3]
  elif a == 3:
    if b > 3:
      raise ValueError
    return 65536
  else:
    raise ValueError
Run Code Online (Sandbox Code Playgroud)

鉴于 m >= 3,ack(m, n) == up_arrow(m - 2, n + 3) - 3 (另请参阅https://en.wikipedia.org/wiki/Ackermann_function#Table_of_values),我们可以计算阿克曼函数的可行值:

def ack(m, n):
  if n < 0:
    raise ValueError
  if m in (0, 1):
    return n + (m + 1)  # This may run out of memory for large n.
  elif m == 2:
    return (n << 1) + 3  # This may run out of memory for large n.
  elif m == 3:
    if n >> 360:
      raise ValueError
    return (1 << (n + 3)) - 3  # This may run out of memory for large n.
  elif m == 4:
    if n > 2:
      raise ValueError
    if n == 2:
      return (1 << 65536) - 3
    return (13, 65533)[n]
  elif m == 5:
    if n > 0:
      raise ValueError
    return 65533
  else:
    raise ValueError

print([ack(m, 0) for m in range(6)])
print([ack(m, 1) for m in range(5)])
print([ack(m, 2) for m in range(5)])
print([ack(m, 3) for m in range(4)])
print([ack(m, 4) for m in range(4)])
Run Code Online (Sandbox Code Playgroud)


Huy*_* Le 7

问题有标签Python,但是用来ctypes调用C++代码也是Python的特性。它包含 3 个版本(我测量了 的时间成本result = [ackermann(m,n) for m<=3, n <= 17],就像最上面的答案一样)

  • 使用 int64_t 和带有 memoization => 32ms 的递归。它仅适用于小输入(就像顶部答案中的所有算法一样)。它专注于优化递归函数。
  • 使用bignum和与上面相同的算法=> 32ms由于我使用的bignum lib是使用64位基本类型实现的,如果结果适合小数,那么与int64_t直接使用没有什么不同。我添加这个是因为评论提到比较 C++int64_t与 Python 无限长度整数是不公平的。
  • 使用 bignum + 数学公式 => ~~0.14 ms这也适用于更大的输入。

摘要:如果使用数学,则比顶级 Python 答案快 3000 倍,如果不使用数学,则快 13-15 倍。我使用这个单头文件 bignum 库。Stackoverflow 有答案长度限制,因此我无法在此处复制文件。

// ackermann.cpp
#include <iostream>
#include <vector>
#include <chrono>
#include <string>
#include "num.hpp"
using namespace std;

class MyTimer {
    std::chrono::time_point<std::chrono::system_clock> start;

public:
    void startCounter() {
        start = std::chrono::system_clock::now();
    }

    int64_t getCounterNs() {
        return std::chrono::duration_cast<std::chrono::nanoseconds>(std::chrono::system_clock::now() - start).count();
    }

    int64_t getCounterMs() {
        return std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::system_clock::now() - start).count();
    }

    double getCounterMsPrecise() {
        return std::chrono::duration_cast<std::chrono::nanoseconds>(std::chrono::system_clock::now() - start).count()
                / 1000000.0;
    }
};

extern "C" {
  int64_t ackermann(int64_t m, int64_t n) {
    static std::vector<int64_t> cache[4];
    // special signal to clear cache
    if (m < 0 && n < 0) {      
      for (int i = 0; i < 4; i++) {
        cache[i].resize(0);
        cache[i].shrink_to_fit();      
      }
      return -1;
    }
    
    if (n >= cache[m].size()) {
      int cur = cache[m].size();
      cache[m].resize(n + 1);
      for (int i = cur; i < n; i++) cache[m][i] = ackermann(m, i);
    }

    if (cache[m][n]) return cache[m][n];
    if (m == 0) return cache[m][n] = n + 1;
    
    // These 3 lines are kinda cheating, since it uses math formula for special case
    // So I commented them out because the question is about optimizing recursion.
    // if (m == 1) return cache[m][n] = n + 2;
    // if (m == 2) return cache[m][n] = 2 * n + 3;
    // if (m == 3) return cache[m][n] = (1LL << (n + 3)) - 3;

    if (n == 0) return cache[m][n] = ackermann(m - 1, 1);

    return cache[m][n] = ackermann(m - 1, ackermann(m, n - 1));        
  }

  Num ackermann_bignum_smallres(int64_t m, int64_t n) {
    static std::vector<Num> cache[4];
    // special signal to clear cache
    if (m < 0 && n < 0) {      
      for (int i = 0; i < 4; i++) {
        cache[i].resize(0);
        cache[i].shrink_to_fit();      
      }
      return -1;
    }
    
    if (n >= cache[m].size()) {
      int cur = cache[m].size();
      cache[m].resize(n + 1);
      for (int i = cur; i < n; i++) cache[m][i] = ackermann(m, i);
    }

    if (cache[m][n] > 0) return cache[m][n];
    if (m == 0) return cache[m][n] = n + 1;

    if (n == 0) return cache[m][n] = ackermann(m - 1, 1);

    return cache[m][n] = ackermann(m - 1, ackermann(m, n - 1));
  }

  //-----
  Num bignum_pow(const Num& x, const Num& n) {
    if (n == 0) return 1;
    Num mid = bignum_pow(x, n / 2);
    if (n % 2 == 0) return mid * mid;
    else return mid * mid * x;
  }

  Num ackermann_bignum(Num m, Num n) {
    if (m <= 1) return n + (m + 1);
    else if (m == 2) return Num(2) * n + 3;
    else if (m == 3) return bignum_pow(2, n + 3) - 3;
    else {
      cout << "Don't put m >= 4\n";
      exit(1);
    } 
  }
}

Num dummy = 0;
int main(int argc, char* argv[])
{
  int test_type = 0;
  if (argc > 1) {
    try {
      test_type = std::stoi(string(argv[1]));
    } catch (...) {
      test_type = 0;
    }
  }
  int lim = (test_type == 0) ? 63 : 17;

  MyTimer timer;
  timer.startCounter();
    
  for (int m = 0; m <= 3; m++)
  for (int n = 0; n <= lim; n++) {
    if (test_type == 0) {
      dummy = ackermann_bignum(m, n);      
    } else if (test_type == 1) {
      dummy = ackermann_bignum_smallres(m, n);
    } else {
      dummy = ackermann(m, n);      
    }
    cout << "ackermann(" << m << ", " << n << ") = " << dummy << "\n";    
  }

  cout << "ackermann cost = " << timer.getCounterMsPrecise() << "\n";
}

Run Code Online (Sandbox Code Playgroud)

使用以下命令编译上面的内容g++ ackermann.cpp -shared -fPIC -O3 -std=c++17 -o ackermann.so

要单独运行,请使用

g++ -o main_cpp ackermann.cpp -O3 -std=c++17
./main_cpp
Run Code Online (Sandbox Code Playgroud)

然后在同一文件夹中运行此脚本(修改自@Stefan Pochmann 答案)

def A_Stefan_row_class(m, n):
    class A0:
        def __getitem__(self, n):
            return n + 1
    class A:
        def __init__(self, a):
            self.a = a
            self.n = 0
            self.value = a[1]
        def __getitem__(self, n):
            while self.n < n:
                self.value = self.a[self.value]
                self.n += 1
            return self.value
    a = A0()
    for _ in range(m):
        a = A(a)
    return a[n]


from collections import defaultdict

def A_Stefan_row_lists(m, n):
    memo = defaultdict(list)
    def a(m, n):
        if not m:
            return n + 1
        if m not in memo:
            memo[m] = [a(m-1, 1)]
        Am = memo[m]
        while len(Am) <= n:
            Am.append(a(m-1, Am[-1]))
        return Am[n]
    return a(m, n)


from itertools import count

def A_Stefan_generators(m, n):
    a = count(1)
    def up(a, x=1):
        for i, ai in enumerate(a):
            if i == x:
                x = ai
                yield x
    for _ in range(m):
        a = up(a)
    return next(up(a, n))


def A_Stefan_paper(m, n):
    next = [0] * (m + 1)
    goal = [1] * m + [-1]
    while True:
        value = next[0] + 1
        transferring = True
        i = 0
        while transferring:
            if next[i] == goal[i]:
                goal[i] = value
            else:
                transferring = False
            next[i] += 1
            i += 1
        if next[m] == n + 1:
            return value


def A_Stefan_generators_2(m, n):
    def a0():
        n = yield
        while True:
            n = yield n + 1
    def up(a):
        next(a)
        a = a.send
        i, x = -1, 1
        n = yield
        while True:
            while i < n:
                x = a(x)
                i += 1
            n = yield x
    a = a0()
    for _ in range(m):
        a = up(a)
    next(a)
    return a.send(n)


def A_Stefan_m_recursion(m, n):
    ix = [None] + [(-1, 1)] * m
    def a(m, n):
        if not m:
            return n + 1
        i, x = ix[m]
        while i < n:
            x = a(m-1, x)
            i += 1
        ix[m] = i, x
        return x
    return a(m, n)


def A_Stefan_function_stack(m, n):
    def a(n):
        return n + 1
    for _ in range(m):
        def a(n, a=a, ix=[-1, 1]):
            i, x = ix
            while i < n:
                x = a(x)
                i += 1
            ix[:] = i, x
            return x
    return a(n)


from itertools import count, islice

def A_Stefan_generator_stack(m, n):
    a = count(1)
    for _ in range(m):
        a = (
            x
            for a, x in [(a, 1)]
            for i, ai in enumerate(a)
            if i == x
            for x in [ai]
        )
    return next(islice(a, n, None))


from itertools import count, islice

def A_Stefan_generator_stack2(m, n):
    a = count(1)
    def up(a):
        i, x = 0, 1
        while True:
            i, x = x+1, next(islice(a, x-i, None))
            yield x
    for _ in range(m):
        a = up(a)
    return next(islice(a, n, None))


def A_Stefan_generator_stack3(m, n):
    def a(m):
        if not m:
            yield from count(1)
        x = 1
        for i, ai in enumerate(a(m-1)):
            if i == x:
                x = ai
                yield x
    return next(islice(a(m), n, None))


def A_Stefan_generator_stack4(m, n):
    def a(m):
        if not m:
            return count(1)
        return (
            x
            for x in [1]
            for i, ai in enumerate(a(m-1))
            if i == x
            for x in [ai]
        )
    return next(islice(a(m), n, None))


def A_templatetypedef(i, n):
    positions = [-1] * (i + 1)
    values = [0] + [1] * i
    
    while positions[i] != n:       
        values[0]    += 1
        positions[0] += 1
            
        j = 1
        while j <= i and positions[j - 1] == values[j]:
            values[j] = values[j - 1]
            positions[j] += 1
            j += 1

    return values[i]

import ctypes
mylib = ctypes.CDLL('./ackermann.so')
mylib.ackermann.argtypes = [ctypes.c_int64, ctypes.c_int64]
mylib.ackermann.restype = ctypes.c_int64

def c_ackermann(m, n):
    return mylib.ackermann(m,n)

funcs = [
    c_ackermann,
    A_Stefan_row_class,
    A_Stefan_row_lists,
    A_Stefan_generators,
    A_Stefan_paper,
    A_Stefan_generators_2,
    A_Stefan_m_recursion,
    A_Stefan_function_stack,
    A_Stefan_generator_stack,
    A_Stefan_generator_stack2,
    A_Stefan_generator_stack3,
    A_Stefan_generator_stack4,
    A_templatetypedef
]

N = 18
args = (
    [(0, n) for n in range(N)] +
    [(1, n) for n in range(N)] +
    [(2, n) for n in range(N)] +
    [(3, n) for n in range(N)]
)

from time import time

def print(*args, print=print, file=open('out.txt', 'w')):
    print(*args)
    print(*args, file=file, flush=True)
    
expect = none = object()
for _ in range(3):
  for f in funcs:
    t = time()
    result = [f(m, n) for m, n in args]
    # print(f'{(time()-t) * 1e3 :5.1f} ms ', f.__name__)
    print(f'{(time()-t) * 1e3 :5.0f} ms ', f.__name__)
    if expect is none:
        expect = result
    elif result != expect:
        raise Exception(f'{f.__name__} failed')
    del result
  print()

  c_ackermann(-1, -1)

Run Code Online (Sandbox Code Playgroud)

输出:

   32 ms  c_ackermann
 1897 ms  A_Stefan_row_class
 1427 ms  A_Stefan_row_lists
  437 ms  A_Stefan_generators
 1366 ms  A_Stefan_paper
  479 ms  A_Stefan_generators_2
  801 ms  A_Stefan_m_recursion
  725 ms  A_Stefan_function_stack
  716 ms  A_Stefan_generator_stack
 1113 ms  A_Stefan_generator_stack2
  551 ms  A_Stefan_generator_stack3
  682 ms  A_Stefan_generator_stack4
 1622 ms  A_templatetypedef
Run Code Online (Sandbox Code Playgroud)

  • 此解决方案对于“mylib.ackermann(3, 63)”(以及许多其他解决方案)来说是不正确的。它应该返回 9223372036854775808,但它返回一些特定于体系结构的错误值,例如 amd64 上的 1。这是因为该值不适合 int64_t。 (2认同)