使用泰勒级数逼近cos

html5 • • 问答

我正在使用泰勒级数来计算一个数的 cos，对于小数，该函数返回准确的结果，例如cos(5)给出0.28366218546322663. 但是对于较大的数字，它会返回不准确的结果，例如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)put 它仍然返回一个用科学记数法书写的数字，1.2194074101485173e+225并带有 e+ 部分。math.cos(1000)给0.5623790762907029，我怎样才能四舍五入我的数字，使它们与 math.cos 方法相同？

回答

麦克劳林级数使用欧拉的思想通过适当的多项式来近似函数的值。多项式显然与函数不同，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|。由于以下术语切换符号，因此它也是结果整体精度的上限。

以上是使用泰勒级数逼近cos的全部内容。
THE END
分享
二维码
< <上一篇
下一篇>>