快速弧算法算法？

我有自己的，非常快的cos函数：

float sine(float x) { const float B = 4/pi; const float C = -4/(pi*pi); float y = B * x + C * x * abs(x); // const float Q = 0.775; const float P = 0.225; y = P * (y * abs(y) - y) + y; // Q * y + P * y * abs(y) return y; } float cosine(float x) { return sine(x + (pi / 2)); } 

但是现在当我描述时，我发现acos（）正在杀死处理器。 我不需要高度精确。 什么是计算acos的快速方法（x）谢谢。

Answer 1

一个简单的三次近似，x∈{-1，-½，0，½，1}的拉格朗日多项式是：

 double acos(x) { return (-0.69813170079773212 * x * x - 0.87266462599716477) * x + 1.5707963267948966; } 

它的最大误差约为0.18弧度。

Answer 2

有多余的记忆？ 查找表（如果需要，带插值）将是最快的。

Answer 3

nVidia有一些很好的资源 ，展示了如何近似非常昂贵的数学函数，例如： acos asin atan2等等……

当执行速度比精度更重要（在合理范围内）时，这些算法会产生良好的结果。 这是他们的acosfunction：

 // Absolute error <= 6.7e-5 float acos(float x) { float negate = float(x < 0); x = abs(x); float ret = -0.0187293; ret = ret * x; ret = ret + 0.0742610; ret = ret * x; ret = ret - 0.2121144; ret = ret * x; ret = ret + 1.5707288; ret = ret * sqrt(1.0-x); ret = ret - 2 * negate * ret; return negate * 3.14159265358979 + ret; } 

以下是计算acos（0.5）时的结果：

 nVidia: result: 1.0471513828611643 math.h: result: 1.0471975511965976 

那非常接近！ 根据您所需的精确度，这可能是一个很好的选择。

Answer 4

我有我自己的。 它非常准确，速度很快。 它基于我围绕四次收敛建立的定理。 它真的很有趣，你可以看到方程式以及它使我的自然对数近似收敛的速度： https ： //www.desmos.com/calculator/yb04qt8jx4

这是我的arccos代码：

 function acos(x) local a=1.43+0.59*xa=(a+(2+2*x)/a)/2 local b=1.65-1.41*xb=(b+(2-2*x)/b)/2 local c=0.88-0.77*xc=(c+(2-a)/c)/2 return (8*(c+(2-a)/c)-(b+(2-2*x)/b))/6 end 

很多只是平方根近似。 它的效果也非常好，除非你太接近于取0的平方根。它的平均误差（不包括x = 0.99到1）为0.0003。 但问题是，在0.99时，它开始变得很糟糕，而在x = 1时，精度的差异变为0.05。 当然，这可以通过在平方根（lol nope）上进行更多迭代来解决，或者只是一点点，如果x> 0.99然后使用一组不同的平方根线性化，但这会使代码变得冗长而丑陋。

如果你不那么关心准确性，那么你可以每平方根进行一次迭代，这应该仍然可以保持在0.0162范围内的某个位置，或者精度如下：

 function acos(x) local a=1.43+0.59*xa=(a+(2+2*x)/a)/2 local b=1.65-1.41*xb=(b+(2-2*x)/b)/2 local c=0.88-0.77*xc=(c+(2-a)/c)/2 return 8/3*cb/3 end 

如果您对它没问题，可以使用预先存在的平方根代码。 它会摆脱x = 1时有点疯狂的等式：

 function acos(x) local a = math.sqrt(2+2*x) local b = math.sqrt(2-2*x) local c = math.sqrt(2-a) return 8/3*db/3 end 

坦率地说，如果你真的很紧张，请记住你可以将arccos线性化为3.14159-1.57079x并且只做：

 function acos(x) return 1.57079-1.57079*x end 

无论如何，如果你想查看我的arccos近似方程式列表，你可以访问https://www.desmos.com/calculator/tcaty2sv8l我知道我的近似值对于某些事情不是最好的，但是如果你’如果我的近似值很有用，请使用它们，但请尽量给予我信任。

Answer 5

您可以采用的另一种方法是使用复数。 来自de Moivre的公式 ，

ⅈx= cos（π/ 2 * x）+ⅈ* sin（π/ 2 * x）

设θ=π/ 2 * x。 那么x =2θ/π，那么

• sin（θ）=ℑ（ⅈ^ ^2θ/π ）
• cos（θ）=ℜ（ⅈ^ ^2θ/π ）

你如何计算没有罪和余的powers的权力？ 从预先计算的表开始，功率为2：

• ⅈ4= 1
• ⅈ2= -1
• ⅈ1=ⅈ
• ⅈ1 ^{/ 2} = 0.7071067811865476 + 0.7071067811865475 *ⅈ
• ⅈ1 ^{/ 4} = 0.9238795325112867 + 0.3826834323650898 *ⅈ
• ⅈ1 ^{/ 8} = 0.9807852804032304 + 0.19509032201612825 *ⅈ
• ⅈ1/ ¹⁶ = 0.9951847266721969 + 0.0980171403295606 *ⅈ
• ⅈ1 ^/ 32 = 0.9987954562051724 + 0.049067674327418015 *ⅈ
• ⅈ1 ^{/ 64} = 0.9996988186962042 + 0.024541228522912288 *ⅈ
• ⅈ1/ 128 = 0.9999247018391445 + 0.012271538285719925 *ⅈ
• ⅈ1/ 256 = 0.9999811752826011 + 0.006135884649154475 *ⅈ

要计算ⅈx的任意值，请将指数近似为二进制分数，然后将表中的相应值相乘。

例如，要找到72°=0.8π/ 2的sin和cos：

ⅈ0.8≈ⅈ205/ ²⁵⁶ =ⅈ0b11001101=ⅈ1 ^{/ 2} *ⅈ1 ^{/ 4} *ⅈ1 ^/ 32 *ⅈ1 ^/ 64 *ⅈ1 ^/ 256
= 0.3078496400415349 + 0.9514350209690084 *ⅈ

• sin（72°）≈0.9514350209690084（“精确”值为0.9510565162951535）
• cos（72°）≈0.3078496400415349（“精确”值为0.30901699437494745）。

要查找asin和acos，可以将此表与Bisection方法一起使用：

例如，要找到asin（0.6）（3-4-5三角形中的最小角度）：

• ⅈ0= 1 + 0 *ⅈ。 罪太小，所以将x增加1/2。
• ⅈ1 ^{/ 2} = 0.7071067811865476 + 0.7071067811865475 *ⅈ。 罪太大，所以将x减少1/4。
• ⅈ1 ^{/ 4} = 0.9238795325112867 + 0.3826834323650898 *ⅈ。 罪太小，所以将x增加1/8。
• ⅈ3 ^{/ 8} = 0.8314696123025452 + 0.5555702330196022 *ⅈ。 罪仍然太小，所以将x增加1/16。
• ⅈ7/ ¹⁶ = 0.773010453362737 + 0.6343932841636455 *ⅈ。 罪过大，所以将x减少1/32。
• ⅈ13/ ³² = 0.8032075314806449 + 0.5956993044924334 *ⅈ。

每次增加x，乘以相应的power功率。 每次减小x，除以相应的power功率。

如果我们在这里停止，我们得到acos（0.6）≈13/ 32 *π/ 2 = 0.6381360077604268（“精确”值是0.6435011087932844。）

当然，准确性取决于迭代次数。 对于快速和脏的近似，使用10次迭代。 对于“高精度”，使用50-60次迭代。

Answer 6

您可以使用dan04建议的多项式来近似反余弦，但多项式是-1和1附近的非常差的近似，其中反余弦的导数变为无穷大。 当你增加多项式的次数时，你会迅速达到递减收益，并且仍然很难在端点周围得到一个很好的近似值。 在这种情况下， 有理函数 （两个多项式的商）可以给出更好的近似。

 acos(x) ≈ π/2 + (ax + bx³) / (1 + cx² + dx⁴) 

哪里

 a = -0.939115566365855 b = 0.9217841528914573 c = -1.2845906244690837 d = 0.295624144969963174 

在区间（-1,1）上的最大绝对误差为0.017弧度（0.96度）。 这是一个图 （黑色的反余弦，红色的三次多项式近似，上面的蓝色函数）用于比较：

选择上述系数以最小化整个域上的最大绝对误差。 如果您愿意在端点处允许更大的错误，则可以使间隔（-0.98,0.98）上的错误更小。 度数为5的分子和度数为2的分母大约与上述函数一样快，但稍微不准确。 以性能为代价，您可以通过使用更高次多项式来提高精度。

关于性能的注意事项：计算两个多项式仍然非常便宜，您可以使用融合乘法 – 加法指令。 划分并不是那么糟糕，因为你可以使用硬件倒数近似和乘法。 与acos近似中的误差相比，倒数近似中的误差可以忽略不计。 在2.6 GHz Skylake i7上，这种近似可以使用AVX每6个周期执行大约8个反余弦。 （即吞吐量，延迟超过6个周期。）

Answer 7

这是一个很棒的网站，有很多选项： https ： //www.ecse.rpi.edu/Homepages/wrf/Research/Short_Notes/arcsin/onlyelem.html

我个人用以下代码去了Chebyshev-Pade商的近似值：

 double arccos(double x) { const double pi = 3.141592653; return pi / 2 - (.5689111419 - .2644381021*x - .4212611542*(2*x - 1)*(2*x - 1) + .1475622352*(2*x - 1)*(2*x - 1)*(2*x - 1)) / (2.006022274 - 2.343685222*x + .3316406750*(2*x - 1)*(2*x - 1) + .02607135626*(2*x - 1)*(2*x - 1)*(2*x - 1)); } 

Answer 8

精确到约0.5度的快速反余弦实现可以基于对于[0,1]中的x，acos（x）≈√（2 *（1-x））的观察 。 额外的比例因子提高了接近零的准确度。 可以通过简单的二进制搜索找到最佳因子。 负参数根据acos（-x）=π – acos（x）处理。

 #include  #include  #include  #include  #include  // Approximate acos(a) with relative error < 5.15e-3 // This uses an idea from Robert Harley's posting in comp.arch.arithmetic on 1996/07/12 // https://groups.google.com/forum/#!original/comp.arch.arithmetic/wqCPkCCXqWs/T9qCkHtGE2YJ float fast_acos (float a) { const float PI = 3.14159265f; const float C = 0.10501094f; float r, s, t, u; t = (a < 0) ? (-a) : a; // handle negative arguments u = 1.0f - t; s = sqrtf (u + u); r = C * u * s + s; // or fmaf (C * u, s, s) if FMA support in hardware if (a < 0) r = PI - r; // handle negative arguments return r; } float uint_as_float (uint32_t a) { float r; memcpy (&r, &a, sizeof(r)); return r; } int main (void) { double maxrelerr = 0.0; uint32_t a = 0; do { float x = uint_as_float (a); float r = fast_acos (x); double xx = (double)x; double res = (double)r; double ref = acos (xx); double relerr = (res - ref) / ref; if (fabs (relerr) > maxrelerr) { maxrelerr = fabs (relerr); printf ("xx=% 15.8e res=% 15.8e ref=% 15.8e rel.err=% 15.8e\n", xx, res, ref, relerr); } a++; } while (a); printf ("maximum relative error = %15.8e\n", maxrelerr); return EXIT_SUCCESS; } 

上述测试脚手架的输出应该类似于：

 xx= 0.00000000e+000 res= 1.56272149e+000 ref= 1.57079633e+000 rel.err=-5.14060021e-003 xx= 2.98023259e-008 res= 1.56272137e+000 ref= 1.57079630e+000 rel.err=-5.14065723e-003 xx= 8.94069672e-008 res= 1.56272125e+000 ref= 1.57079624e+000 rel.err=-5.14069537e-003 xx=-2.98023259e-008 res= 1.57887137e+000 ref= 1.57079636e+000 rel.err= 5.14071269e-003 xx=-8.94069672e-008 res= 1.57887149e+000 ref= 1.57079642e+000 rel.err= 5.14075044e-003 maximum relative error = 5.14075044e-003