首页 > 解决方案 > 如何在 Python 中获得高精度的数学计算

问题描述

我什至不确定这是否是我想要实现的最佳措辞,但我正在尝试评估 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、阶乘等)中都需要非常精确,这样我就不用自己解决问题了?

标签: python

解决方案


此代码适用于 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

推荐阅读