如何在Python中实现高精度数学计算

Question

我甚至不确定这是否是我试图实现的最好的措辞，但我正在尝试用Python计算幂级数，正如我所想象的那样，因为我是Python的新手，而且由于使用浮点数时，幂级数计算很棘手，我的计算正在崩溃。

我的幂级数是简单素数幂计数函数J(x)的近似值。在Mathematica中，我可以计算这个幂级数(我已经计算到了素数151，非常耗时)，但即使在Mathematica中，如果我在公式中使用浮点数(而不是符号，然后转换为十进制)，结果也是错误的。我的计划是创建一个可执行文件，让它在Unix或Windows上的作业远程桌面上运行。

这是我在Python中创建的代码(使用scipy)：

x = 12
M = 9 * x
soma = 0
for i in range(1, M+1):
    termo1 = (-1)**i * x * (2*pi*x)**(2*i) / (2*i+1)
    for j in range(1, i+1):
        termo2 = (-1)**j * log(zeta(2*j))/((2*pi)**(2*j) * factorial(2*i+1-2*j))
        soma = soma + termo1 * termo2

soma = -4 * soma

当我们刚接触一些东西时，弄清楚如何做一些简单的事情的压力可能是很难的(我是一个新手，所以仅仅说使用Decimal对我来说是希腊语)。这是我得到的错误：

File "D:/iTunes/Python/PyCharm/Zeta.py", line 31, in <module>
termo1 = (-1)**i * x * (2*pi*x)**(2*i) / (2*i+1)
OverflowError: (34, 'Result too large')

我该如何解决这个问题？或者更好的是，有没有一个软件包假设我在所有函数(zeta，log，factorial等)中都需要很高的精度。这就省去了我自己弄清楚事情的麻烦？

Answer 1

这段代码适用于Google Colab。请根据需要更新精度。请注意，乘法和除法会使最后几个数字变得不准确。

import math
from math import pi, log, factorial
from scipy.special import zeta
from decimal import Decimal, getcontext

getcontext().prec = 100
x = Decimal(12)
M = int(9 * x)
soma = Decimal(0)
power_of_one = Decimal(-1)
power_of_two_pi_x = (2 * Decimal(pi) * x) * (2 * Decimal(pi) * x)
print(type(power_of_two_pi_x))
for i in range(1, M+1):
    termo1 = power_of_one * x * power_of_two_pi_x / (2*i+1)
    power_of_one *= -1
    power_of_two_pi_x *= (2 * Decimal(pi) * x) * (2 * Decimal(pi) * x)

    power_of_one_j = -1
    power_of_two_pi = Decimal(4 * pi * pi)
    factori = factorial(2 * i - 1)
    for j in range(1, i+1):
        termo2 = power_of_one_j * Decimal(log(zeta(2*j)))/power_of_two_pi * Decimal(factori)
        power_of_one_j *= -1
        power_of_two_pi *= Decimal(4 * pi * pi)
        factori //= 2 * i + 1 - 2 * j
        if i > j:
            factori //= 2 * i - 2* j

        soma += termo1 * termo2

soma = -4 * soma