经典的 Clenshaw 递归（请参阅此处的算法 3.1 ）不如准确。因此，Reich 提出了对其进行修改，并在算法 3.2以及Oliver中进行了讨论。$x\to \pm 1$

在尝试推导两种算法之间的过渡点时，我为这两种算法生成了 ULP 图，并且惊讶地发现我在任何地方都看不到任何差异。

首先，我尝试了 64 个 Unif(0,1) 系数，带有扰动横坐标和扰动系数：

接下来，64 个 Unif(0,1) 精确系数，扰动横坐标：

我还尝试了 Unif(-1,1) 系数，扰动和非扰动，以及更能代表平滑函数的衰减系数。但每次 ULP 图的外观都没有差异。

哪些数据能真正证明 Reinch 修正对 Clenshaw 递归的优越性？

我的代码，供参考：

template<class Real>
inline Real chebyshev_clenshaw_recurrence(const Real* const c, size_t length, const Real& x)
{
    if (length < 2)
    {
        if (length == 0)
        {
            return 0;
        }
        return c[0]/2;
    }
    Real b2 = 0;
    Real b1 = c[length -1];
    for(size_t j = length - 2; j >= 1; --j)
    {
        Real tmp = 2*x*b1 - b2 + c[j];
        b2 = b1;
        b1 = tmp;
    }
    return x*b1 - b2 + c[0]/2;
}

template<class Real>
inline Real modified_chebyshev_clenshaw_recurrence(const Real* const c, size_t length, const Real& x)
{
    using std::abs;
    if (length < 2)
    {
        if (length == 0)
        {
            return 0;
        }
        return c[0]/2;
    }
    Real cutoff = 0.0;
    if (abs(x) <= cutoff)
    {
        Real b2 = 0;
        Real b1 = c[length -1];
        for(size_t j = length - 2; j >= 1; --j)
        {
            Real tmp = 2*x*b1 - b2 + c[j];
            b2 = b1;
            b1 = tmp;
        }
        return x*b1 - b2 + c[0]/2;
    }
    else if (x > cutoff)
    {
        Real b = c[length -1];
        Real d = b;
        Real b2 = b;
        for (size_t r = length - 2; r >= 1; --r)
        {
            d = 2*(x-1)*b + d + c[r];
            b2 = b;
            b = d + b;
        }
        return x*b - b2 + c[0]/2;
    }
    else
    {
        Real b = c[length -1];
        Real d = b;
        Real b2 = b;
        for (size_t r = length - 2; r >= 1; --r)
        {
            d = 2*(x+1)*b - d + c[r];
            b2 = b;
            b = d - b;
        }
        return x*b - b2 + c[0]/2;
    }
}

此外，我通过在 long double 中生成系数和横坐标并将它们四舍五入为浮点数来扰动系数和横坐标，并计算了相对于 long double 精度的“精确”数字。