泰勒级数抵消误差累积：如何优化代码适配更小阈值？

我在研究泰勒级数项数与开始发散的xᵢ之间的关系：以10⁻⁴作为误差阈值时，拟合效果良好；但使用10⁻⁸这类更高精度阈值时，较大的n会出现异常表现。推测是大n时大数相减的抵消误差（"大数"-"超大数"≈1）导致Python精度不足，求代码优化方案以适配更小的阈值。

原代码：

import math
import matplotlib.pyplot as plt
import numpy as np

m=50

def taylor_cos(x,n):
    sum=0
    for i in range(0, n):
        sum=sum+ (((-1)**i) * (x**(2*i)) / math.factorial(2*i))
    return sum 

X=[]
Y=[]
arr=[]

for N in range(1,m):
    for x_i in np.linspace(0,N*5,500):
        if abs(math.cos(x_i)-taylor_cos(x_i,N))>10**-4:
            arr.append((N,x_i))
            X.append(N)
            Y.append(x_i)
            break

coeffs = np.polyfit(X, Y, 1)
fit = np.poly1d(coeffs)
print(coeffs)

plt.figure(figsize=(10,10))
plt.grid(True,alpha=0.3)
plt.axhline(0,color='gray',linewidth=0.5)

plt.scatter(X,Y,color='black',s=10)

plt.plot(X, fit(X), color='red', linewidth=1, linestyle='--')

plt.show()

1. 改用递推式计算泰勒项，根除大数运算 #

原代码直接计算x^(2i)和math.factorial(2i)会产生极大中间值，触发浮点数精度损失。递推式可基于前一项迭代计算，完全避免超大数：

• 第i项 = 第i-1项 × (-x²) / [(2i-1)×2i]
• 初始项（i=0）为1

2. 可选：用相对误差替代绝对误差 #

对于大x，cos(x)的绝对值可能极小（比如x接近π/2的奇数倍），固定绝对误差阈值10⁻⁸会不合理。可改用相对误差判断：abs((cos(x)-taylor_cos(x,n))/cos(x))>1e-8，注意避开cos(x)接近0的点（加小偏移量或直接跳过）。

3. 扩展浮点数精度（按需使用） #

Python默认双精度浮点数不够时，可改用decimal模块的高精度浮点数，或numpy的longdouble类型，代价是计算速度下降。

import math
import matplotlib.pyplot as plt
import numpy as np

m = 50

def taylor_cos(x, n):
    if n == 0:
        return 0.0
    current_term = 1.0  # i=0的初始项
    sum_val = current_term
    x_squared = -x * x  # 提前计算-x²，减少重复运算
    for i in range(1, n):
        # 递推计算当前项，避免超大数
        denominator = (2*i - 1) * (2*i)
        current_term *= x_squared / denominator
        sum_val += current_term
    return sum_val

X = []
Y = []
error_threshold = 1e-8  # 更小的精度阈值

for N in range(1, m):
    x_start, x_end = 0, N * 5
    x_list = np.linspace(x_start, x_end, 1000)  # 增加采样点，提升搜索精度
    for x_i in x_list:
        cos_true = math.cos(x_i)
        cos_taylor = taylor_cos(x_i, N)
        # 避开cos(x)接近0的点，避免相对误差分母为0
        if abs(cos_true) < 1e-10:
            continue
        # 用绝对误差判断，也可替换为相对误差
        error = abs(cos_true - cos_taylor)
        # error = abs((cos_true - cos_taylor) / cos_true)  # 相对误差
        if error > error_threshold:
            X.append(N)
            Y.append(x_i)
            break

coeffs = np.polyfit(X, Y, 1)
fit = np.poly1d(coeffs)
print(f"拟合系数：{coeffs}")

plt.figure(figsize=(10,10))
plt.grid(True, alpha=0.3)
plt.axhline(0, color='gray', linewidth=0.5)
plt.scatter(X, Y, color='black', s=10)
plt.plot(X, fit(X), color='red', linewidth=1, linestyle='--')
plt.xlabel("泰勒级数项数N")
plt.ylabel("开始发散的xᵢ")
plt.title("泰勒级数项数与发散点的关系（误差阈值1e-8）")
plt.show()

如果上述优化仍不够，可使用decimal模块提升精度，示例核心代码：

from decimal import Decimal, getcontext

getcontext().prec = 50  # 设置50位高精度

def taylor_cos_decimal(x, n):
    x = Decimal(x)
    if n == 0:
        return Decimal(0)
    current_term = Decimal(1)
    sum_val = current_term
    x_squared = -x * x
    for i in range(1, n):
        denominator = Decimal((2*i -1)*(2*i))
        current_term *= x_squared / denominator
        sum_val += current_term
    return sum_val

内容的提问来源于stack exchange，提问作者Artemii Karapetian