为什么SSE标量sqrt（x）慢于rsqrt（x）* x？
我一直在分析我们在Intel Core Duo上的一些核心math，并且在研究各种平方根的方法时，我注意到一些奇怪的事情：使用SSE标量运算，取相反的平方根并乘以它得到sqrt，比使用本地的sqrt操作码！
我用一个循环来testing它：
inline float TestSqrtFunction( float in ); void TestFunc() { #define ARRAYSIZE 4096 #define NUMITERS 16386 float flIn[ ARRAYSIZE ]; // filled with random numbers ( 0 .. 2^22 ) float flOut [ ARRAYSIZE ]; // filled with 0 to force fetch into L1 cache cyclecounter.Start(); for ( int i = 0 ; i < NUMITERS ; ++i ) for ( int j = 0 ; j < ARRAYSIZE ; ++j ) { flOut[j] = TestSqrtFunction( flIn[j] ); // unrolling this loop makes no difference -- I tested it. } cyclecounter.Stop(); printf( "%d loops over %d floats took %.3f milliseconds", NUMITERS, ARRAYSIZE, cyclecounter.Milliseconds() ); } 我已经用TestSqrtFunction的几个不同的机构来尝试这个，而且我有一些时间确实在困扰着我。 迄今为止最糟糕的是使用本地sqrt（）函数，并让“智能”编译器“优化”。 在24ns / float时，使用x87 FPU这是一个非常糟糕的情况：
inline float TestSqrtFunction( float in ) { return sqrt(in); }
接下来我尝试的是使用一个内部函数强制编译器使用SSE的标量sqrt操作码：
inline void SSESqrt( float * restrict pOut, float * restrict pIn ) { _mm_store_ss( pOut, _mm_sqrt_ss( _mm_load_ss( pIn ) ) ); // compiles to movss, sqrtss, movss }
这是更好的，在11.9ns /浮动。 我也尝试过Carmack古怪的Newton-Rhapson逼近技术 ，它的运行速度甚至比硬件要高出4.3ns / float，尽pipe误差在1/2（对我来说太过分了）。
当我尝试使用SSE运算倒数平方根，然后使用乘法得到平方根（x * 1 /√x=√x）时，这个工作就很复杂了。 尽pipe这需要两个相关的操作，但它是迄今为止最快的解决scheme，在1.24ns / float和精确到2 -14 ：
inline void SSESqrt_Recip_Times_X( float * restrict pOut, float * restrict pIn ) { __m128 in = _mm_load_ss( pIn ); _mm_store_ss( pOut, _mm_mul_ss( in, _mm_rsqrt_ss( in ) ) ); // compiles to movss, movaps, rsqrtss, mulss, movss }
我的问题基本上是什么给 ？ 为什么SSE内置到硬件的平方根操作码要比从其他两个math运算中综合出来慢呢？
我敢肯定，这真的是操作本身的成本，因为我已经validation：
- 所有数据都适合caching，并且访问是顺序的
- 函数被内联
- 展开循环没有任何区别
- 编译器标志设置为完全优化（和组件是好的，我检查）
（ 编辑 ：stephentyrone正确地指出，对长数字串的操作应该使用向量化的SIMD压缩操作，如rsqrtps – 但这里的数组数据结构仅用于testing目的：我真正要测量的是标量性能不能被vector化的代码。）
sqrtss给出了一个正确的舍入结果。 rsqrtss给出了倒数的近似值 ，精确到11位。
sqrtss正在产生一个更准确的结果，因为什么时候需要精度。 rsqrtss满足的情况， rsqrtss存在，但速度是必需的。 如果您阅读了英特尔的文档，您还可以find一个指令序列（倒数平方根逼近，接着是一个单一的牛顿 – 拉夫逊步骤），可以提供几乎完全的精度（如果我没有记错的话，精度可达23位）快于sqrtss 。
编辑：如果速度至关重要，而且你真的在循环中调用了很多值，那么应该使用这些指令的向量化版本， rsqrtps或sqrtps ，这两个指令每个指令处理四个浮点数。
分工也是如此。 MULSS（a，RCPSS（b））比DIVSS（a，b）快得多。 实际上，即使通过Newton-Rhapson迭代提高精度，速度仍然更快。
英特尔和AMD都在他们的优化手册中推荐了这种技术。 在不需要符合IEEE-754的应用程序中，使用div / sqrt的唯一原因是代码可读性。
而不是提供一个答案，实际上可能是不正确的（我也不会去检查或争论caching和其他的东西，让我们说，他们是相同的）我会尽力指出你可以回答你的问题的来源。
差异可能在于如何计算sqrt和rsqrt。 您可以在这里阅读更多http://www.intel.com/products/processor/manuals/ 。 我build议从阅读你正在使用的处理器函数开始，有一些信息，尤其是关于rsqrt（cpu是使用内部查找表与大量近似，这使得它更简单得到结果）。 看起来，rsqrt比sqrt快得多，另外多一个操作（这不是花钱）可能不会改变这里的情况。
编辑：可能值得一提的几个事实：
1.一旦我为graphics库做了一些微观优化，并使用rsqrt来计算vector的长度。 （而不是sqrt，我已经乘以了rsqrt的平方和，这正是你在testing中所做的），并且它performance更好。
2.使用简单的查找表来计算rsqrt可能更容易，对于rsqrt，当x变为无穷大时，1 / sqrt（x）变为0，所以对于小的x，函数值不会改变（很多），而对于sqrt – 它走向无穷大，所以就是这么简单的情况）。
此外，澄清：我不知道我已经在我已经链接的书中find它，但我敢肯定，我已经读了rsqrt正在使用一些查找表，它应该只使用，当结果不需要确切的，虽然 – 我可能是错误的，因为它是前一段时间:)。
牛顿 – 拉夫逊收敛于f(x)的零点，增量等于-f/f' ，其中f'是导数。
对于x=sqrt(y) ，可以尝试使用f(x) = x^2 - y来求解f(x) = 0 ;
然后，增量为：其中dx = -f/f' = 1/2 (x - y/x) = 1/2 (x^2 - y) / x 。
你可以尝试其他函数（比如f(x) = 1/y - 1/x^2 ），但是它们同样复杂。
现在我们来看1/sqrt(y) 。 你可以尝试f(x) = x^2 - 1/y ，但它会同样复杂：例如， dx = 2xy / (y*x^2 - 1) 。 f(x)一个非显而易见的替代select是： f(x) = y - 1/x^2
那么： dx = -f/f' = (y - 1/x^2) / (2/x^3) = 1/2 * x * (1 - y * x^2)
啊! 这不是一个微不足道的expression，但你只有在它的倍增，没有分裂。 =>更快！
并且：完整的更新步骤new_x = x + dx然后读取：
x *= 3/2 - y/2 * x * x这也很容易。
由于这些指令忽略舍入模式，所以速度更快，并且不处理floatin点exception或非常化数字。 由于这些原因，stream水线，推测和执行其他fp指令失序要容易得多。
