使用泰勒级数逼近 cos

Approximating cos using the Taylor series

我正在使用 Taylors series 计算数字的余弦值,对于较小的数字,函数 returns 的准确结果例如 cos(5) 给出 0.28366218546322663。但是对于更大的数字,它 returns 不准确的结果如 cos(1000) 给出 1.2194074101485173e+225

def factorial(n):
    c = n
    for i in range(n-1, 0, -1):
        c *= i
    return c

def cos(x, i=100):
    c = 2
    n = 0
    for i in range(i):
        if i % 2 == 0:
            n += ((x**c) / factorial(c))
        else:
            n -= ((x**c) / factorial(c))
        c += 2
    return 1 - n

我尝试使用 round(cos(1000), 8) 把它仍然 returns 一个用科学记数法写成的数字 1.2194074101485173e+225 与 e+ 部分。 math.cos(1000) 给出 0.5623790762907029,我如何舍入我的数字,使它们与 math.cos 方法相同?

返回的数字只是一个数字;在您打印它之前,它没有符号感。如果您希望控制值的 打印方式 ,我们假设您的打印方式是这样

print(cos(1000))

然后我们可以用format strings来控制输出

print("{:f}".format(cos(1000)))

如果您使用的是 Python 3.6 或更高版本,我们甚至可以 interpolate 将其直接放入字符串文字中。

print(f"{cos(1000):f}")

您可以阅读以上链接以了解有关格式迷你语言的更多详细信息(两种功能之间的语言相同)。例如,如果你想打印特定数量的小数位,你也可以请求。我们可以打印精确的小数点后三位如下

print("{:.3f}".format(cos(1000)))
print(f"{cos(1000):.3f}")

但是,正如疯狂物理学家指出的那样,您的代码还存在一些更多的数学问题,因此我强烈建议您也阅读他的回答。

麦克劳林级数使用欧拉的思想,使用适当的多项式来逼近函数的值。多项式明显不同于 cos(x) 这样的函数,因为它们在某个点都趋于无穷大,而 cos 则不然。一个 100 阶多项式在零的每一侧最多可以逼近函数的 50 个周期。由于 50 * 2pi << 1000,您的多项式 不能 近似 cos(1000).

为了更接近合理的解决方案,多项式的阶数必须至少为 x / pi。您可以尝试计算 300 阶以上的多项式,但由于浮点数的有限精度和阶乘的巨大性,您很可能 运行 陷入一些主要的数值问题。

相反,使用 cos(x) 的周期性并添加以下内容作为函数的第一行:

x %= 2.0 * math.pi

您还需要限制多项式的阶数,以避免阶乘太大而无法放入浮点数的问题。此外,您可以而且应该通过递增先前的结果来计算阶乘,而不是在每次迭代时都从头开始。这是一个具体的例子:

import math

def cos(x, i=30):
    x %= 2 * math.pi
    c = 2
    n = 0
    f = 2
    for i in range(i):
        if i % 2 == 0:
            n += x**c / f
        else:
            n -= x**c / f
        c += 2
        f *= c * (c - 1)
    return 1 - n
>>> print(cos(5), math.cos(5))
0.28366218546322663 0.28366218546322625

>>> print(cos(1000), math.cos(1000))
0.5623790762906707 0.5623790762907029

>>> print(cos(1000, i=86))
...
OverflowError: int too large to convert to float

您可以通过注意增量产品 x**2 / (c * (c - 1)) 来进一步摆脱数值瓶颈。对于比直接阶乘所能支持的大得多的i,这将保持良好的界限:

import math

def cos(x, i=30):
    x %= 2 * math.pi
    n = 0
    dn = x**2 / 2
    for c in range(2, 2 * i + 2, 2):
        n += dn
        dn *= -x**2 / ((c + 1) * (c + 2))
    return 1 - n
>>> print(cos(5), math.cos(5))
0.28366218546322675 0.28366218546322625
>>> print(cos(1000), math.cos(1000))
0.5623790762906709 0.5623790762907029
>>> print(cos(1000, i=86), math.cos(1000))
0.5623790762906709 0.5623790762907029
>>> print(cos(1000, i=1000), math.cos(1000))
0.5623790762906709 0.5623790762907029

请注意,经过某个点后,无论循环多少次,结果都不会改变。这是因为现在 dn 收敛于零,正如欧拉预期的那样。

您可以使用此信息进一步改进循环。由于浮点数的精度有限(具体来说是尾数中的 53 位),您可以在 |dn / n| < 2**-53:

时停止迭代
import math

def cos(x, conv=2**-53):
    x %= 2 * math.pi
    c = 2
    n = 1.0
    dn = -x**2 / 2.0
    while abs(n / dn) > conv:
        n += dn
        c += 2
        dn *= -x**2 / (c * (c - 1))
    return n
>>> print(cos2(5), math.cos(5))
0.28366218546322675 0.28366218546322625
>>> print(cos(1000), math.cos(1000))
0.5623790762906709 0.5623790762907029
>>> print(cos(1000, 1e-6), math.cos(1000))
0.5623792855306163 0.5623790762907029
>>> print(cos2(1000, 1e-100), math.cos(1000))
0.5623790762906709 0.5623790762907029

参数 conv 不仅仅是 |dn/n| 上的界限。由于以下术语切换符号,因此它也是结果整体精度的上限。