SSE优化在数学库中的应用之二

 引自：http://hi.baidu.com/sige_online/blog/item/d8fdfffc8f0033f7fd037fac.html

下面将通过几个简单的运算例子介绍SSE intrinsic的使用。首先，使用SSE需要一个新的头文件

#include <xmmintrin.h>

里面定义了一个新的数据类型，__m128，这是一个128位、４个32位单精度浮点数的结构，如果你正在使用VC.net，你会看到它是一个关键字，被当作一种基本数据类型。要是你不打算使用汇编SSE，那么就没必要深究编译器在幕后到底如何处理__m128类型的数据，你只需要知道里面能存放四个float，而这四个float可以进行并行运算。

 

在定义了__m128后，文件声明一大堆对__m128进行运算的函数，如_mm_add_ps、_mm_sub_ps等等，这就是SSE运算指令的声明。使用SSE优化在这些声明的帮助下变得非常简单，如计算两个向量之和，平时需要每一个元素进行一次加法运算，现在只需要简单地：

__m128 a , b , c;

c = _mm_add_ps( a , b );

这样等价于：

 

float a[4] , b[4] , c[4];

for( int i = 0 ; i < 4 ; ++ i )

    c[i] = a[i] + b[i];

但前者的运算是并行的，在一般情况下效率远比后者要高。况且在描述复杂的运算的时候，如：

 

      a = b * c + d / e;

 

则可以直接写成：

 

__m128 a = _mm_add_ps( _mm_mul_ps( b , c ) , _mm_div_ps( d , e ) );

咋看之下，很多效率至上的人马上就会大叫“昂贵的函数调用啊！Bad smell code!”。其实我正要告诉你，我也是效率至上派的。前面已经说过了，这些看上去貌似函数的调用实际上并非函数，而是所谓intrinsic，它们在编译优化中将被解释为单条或多条SSE指令，而且编译器会自动调节调用顺序以使其最大并行效率。

 

不过除了直接使用这些intrinsic以外，我们还可以把它们封装到类里面，重载运算符，这样就可以把运算写成可读性更强的算术式。如果你不愿意自己动手封装，也可以使用Intel封装好了的F32vec4类，它提供了完备的运算符重载，完全使用SSE，非常方便。

虽然Intel封装好的类已经很完善了，但还有一大堆数学运算需要我们自己动手进行编写，如内积（点积）和外积（叉积）。

首先来看一个比较实用的运算，求倒数。求倒数在很多数学库里都有专门的优化，通常原理都是先求出一个近似值，然后通过Newton-Raphson逼近法求出较精确值，下面的代码摘自NV的fastmath.cpp:

 

#define FP_ONE_BITS 0x3F800000

// r = 1/p

 

#define FP_INV(r,p)                                                   　 \

{                                                                          　　 \

    int _i = 2 * FP_ONE_BITS - *(int *)&(p);                     \

    r = *(float *)&_i;                                                　　 \

    r = r * (2.0f - (p) * r);                                       　　  \

}

而在SSE里也提供了两条求倒数的指令rcpss/rcpps（对应的intrinsic是_mm_rcp_ss与_mm_rcp_ps），不过这两条指令求的并非是精确值，而是近似值，所以我们需要对它的结果进行逼近处理。

 

 

float __rcp<float>( const float& a ) {

    register float r;

    __m128 rcp = _mm_load_ss( &a );

　　rcp = _mm_rcp_ss( rcp );

　　_mm_store_ss( &r, rcp );

　　/* [2 * rcpps(x) - (x * rcpps(x) * rcpps(x))]*/

　　r = 2.0f * r - ( a * r * r );

　　return r;

}

 

 原理一致，只不过我们还可以用_mm_rcp_ps并行求四分量的倒数。如果你还对SSE的威力有所保留，那我建议你设计一个测试单元测试一下使用除法求倒数与使用SSE求倒数，看效率到底是谁更高、高多少。当然，我自己已经测试过很多次了。