使用泰勒级数近似计算cos函数

Question

使用泰勒级数近似计算cos函数

6

我正在使用泰勒级数计算一个数字的余弦值。对于小的数字，该函数会返回精确的结果，例如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)，但它仍然返回一个科学计数法表示的数字1.2194074101485173e+225，其中包含了e+部分。而math.cos(1000)则返回0.5623790762907029。我应该如何四舍五入使得结果与math.cos方法的结果相同？

- user9321739

1

"round(cos(1000), 8)" 这个函数返回的是一个浮点数，它的大小不是十进制表示法、科学计数法或其他任何表示法。因此，它并不会返回科学计数法的数字 1.2194074101485173e+225。实际上，科学计数法只是一种表示方法。这个结果是由你如何打印这个值决定的。 - Alexander

1

看一下余弦的多项式展开式。1000/6.2... > 100，因此你的多项式无法收敛到离零那么远的答案。你需要对参数取模2pi。 - Mad Physicist

4

你关注科学计数法，但完全忽略了你手头的数字是10的255次方倍。 - Mad Physicist

@Jared。我又更新了我的答案。 - Mad Physicist

2个回答

1

返回结果: 返回的数字只是一个数字；直到你打印它，它才有符号的意义。如果你想控制值的打印方式，比如说你这样打印

print(cos(1000))

然后我们可以使用格式化字符串来控制输出

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

如果您使用的是Python 3.6或更新版本，我们甚至可以将其直接插值到字符串文字中。

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

您可以阅读上述链接以了解有关格式迷你语言的更多详细信息（两个功能之间的语言是相同的）。例如，如果您想要打印特定数量的小数位数，也可以请求。我们可以按如下方式精确地打印三位小数。

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

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

- Silvio Mayolo

函数有问题，你不能神奇地将 1.2*10^255 转换成从负一到正一之间的数字。 - Mad Physicist

糟糕！我以为那是 e-255（因此，是一个非常接近零的数字）。我的错误。 - Silvio Mayolo

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

我会把我的答案放在这里，因为我认为格式信息仍然有帮助作用，但你的答案肯定是对这个问题的正确答案。 - Silvio Mayolo

我会取消我的踩，因为我认为这个信息很有用，而且我现在明白你为什么认为它相关了。 - Mad Physicist

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Mad Physicist · Accepted Answer

一个麦克劳林级数运用欧拉的思想，使用适当的多项式来近似计算函数的值。不过，这些多项式显然会与像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|的限制。由于以下项会切换符号，因此它还是结果总精度的上限。