使用泰勒级数逼近 cos

Posted 2023-03-28

【中文标题】使用泰勒级数逼近 cos

【英文标题】：Approximating cos using the Taylor series

【发布时间】：2021-08-12 08:55:34

【问题描述】：

我使用Taylors series 来计算一个数字的余弦值，对于小数字，该函数会返回准确的结果，例如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 方法相同？

【问题讨论】：

"round(cos(1000), 8) put 它仍然返回一个用科学记数法写的数字1.2194074101485173e+225" 实际上并没有。科学记数法就是……一种记数法。 round 的结果是一个浮点数，其大小不是十进制计数法、科学计数法或任何其他计数法。这是您打印值的方式的结果。

看看余弦的多项式展开。 1000/6.2... > 100，所以你的多项式不能收敛到远离零的答案。您需要将参数取模 2pi

您担心科学记数法，但完全忽略了一个事实，即您的数字是 10^255 的倍数

@Jared。我已经更新了我的答案。

【参考方案1】：

返回的数字只是一个数字；在您打印之前，它没有符号意义。如果您希望控制值的打印方式，假设您是这样打印的

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")

但是，正如 Mad Physicist 指出的那样，您的代码还存在一些数学问题，因此我强烈建议您也阅读他的答案。

【讨论】：

函数不行，不能神奇地将1.2*10^255变成负一到一之间的数。

糟糕！我认为阅读e-255（因此，是一个非常接近于零的数字）。我的错。

你的回答现在更有意义了！

我会在这里留下我的答案，因为我认为格式信息仍然很有帮助，但你的绝对是问题的正确答案。

我将删除我的反对票，因为我同意该信息很有用，并且我理解您为什么认为它现在是相关的。

【参考方案2】：

McLaurin 级数使用 Euler 的思想，使用适当的多项式来近似函数的值。多项式显然与cos(x) 之类的函数不同，因为它们都在某个点趋于无穷大，而cos 则不然。一个 100 阶多项式在零的每一侧最多可以近似函数的 50 个周期。由于 50 * 2pi 不能近似于 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 收敛到零，正如 Euler 所期望的那样。

您可以使用此信息进一步改进您的循环。由于浮点数的精度有限（具体而言，尾数为 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| 上的界限。由于以下术语切换符号，它也是结果整体精度的上限。

【讨论】：

当我将x %= 2.0 * math.pi 添加到函数的开头时，我得到n -= ((x**c) / factorial(c)) OverflowError: int too large to convert to float

@Ironstone1_：您使用的术语太多了。 factorial(172) 已经大于最大的浮点值。

也许泰勒级数不是计算 cos 的最佳方法，我会寻找另一个。

@Ironstone1_。大多数快速实现都使用一堆技巧和称为 Chebyshev 多项式的东西

@Ironstone1_。话虽如此，我已经添加了一个具体的实现，它可能足以满足您的当前需求