c++ - exp() 函数数值计算中的系数

Question

我试图了解 http://gruntthepeon.free.fr/ssemath/sse_mathfun.h 中 exp_ps() 的实现或 http://software-lisc.fbk.eu/avx_mathfun/avx_mathfun.h 中的 exp256_ps() .
我几乎了解计算中的所有内容，除了常数 cephes_exp_C2 是如何确定的。看来它提高了计算的准确性。如果将其从计算中删除，则生成的函数会明显更快，但精确度稍差(对于 +/- 10 左右的值，相对误差仍低于 1%)。我在其他数值库中找到了这样的系数，但没有更详细的解释。

Answer 1

经过一番搜索后 Cephes来源，我认为这是Pommier翻译中的错误。这不是我第一次看到 Pommier 代码中的错误。我建议在 Gromacs 中使用数学库.

<小时/>

来自 Cephe 的 exp.c，

static double C1 = 6.93145751953125E-1;
static double C2 = 1.42860682030941723212E-6;
....
px = floor( LOG2E * x + 0.5 );
n = px;
x -= px * C1;
x -= px * C2;

来自波米尔，

_PS_CONST(cephes_exp_C1, 0.693359375);
_PS_CONST(cephes_exp_C2, -2.12194440e-4);  <-- Wrong value
....

//
//  fx = LOG2E * x + 0.5
//
fx = _mm_mul_ps(x, *(v4sf*)_ps_cephes_LOG2EF);
fx = _mm_add_ps(fx, *(v4sf*)_ps_0p5);

//
//  fx = floor(fx)
//
emm0 = _mm_cvttps_epi32(fx);
tmp  = _mm_cvtepi32_ps(emm0);
v4sf mask = _mm_cmpgt_ps(tmp, fx);    
mask = _mm_and_ps(mask, one);
fx = _mm_sub_ps(tmp, mask);

//
//  x -= fx * C1;
//  x -= fx * C2;  (Using z allows for better ILP in this step)
//
tmp = _mm_mul_ps(fx, *(v4sf*)_ps_cephes_exp_C1);
v4sf z = _mm_mul_ps(fx, *(v4sf*)_ps_cephes_exp_C2);
x = _mm_sub_ps(x, tmp);
x = _mm_sub_ps(x, z);