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

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

然后我们把注意力放到一条非常特殊的指令shufps(对应intrinsic是_mm_shuffle_ps)上面。这是一条非常有用的指令,它可以把两个操作数的分量以特定的顺序排列并赋予给目标数。比如
__m128 b = _mm_shuffle_ps( a , a , 0 );
  
则 b 的所有分量都是 a 中下标为0的分量。第三个参数控制分量分配,是一个8bit的常量,这个常量的1~8位分别控制了从两个操作数中选择分量的情况,具体怎么控制将在后面讨论SSE汇编中一并说明,而在使用intrinsic的时候,最好使用_MM_SHUFFLE宏,它可以定义分配情况。下面我们来复习一下叉积的求法。
c = a x b
可以写成:
Vector cross(const Vector& a , const Vector& b ) {
    return Vector(
        ( a[1] * b[2] - a[2] * b[1] ) ,
        ( a[2] * b[0] - a[0] * b[2] ) ,
        ( a[0] * b[1] - a[1] * b[0] ) );
}
那么写成SSE intrinsic形式则是:
/* cross */
__m128 _mm_cross_ps( __m128a , __m128 b ) {
    __m128 ea , eb;
    // set to a[1][2][0][3] , b[2][0][1][3]
    ea = _mm_shuffle_ps( a , a , _MM_SHUFFLE( 3 , 0 , 2 , 1 ) );
    eb = _mm_shuffle_ps( b , b , _MM_SHUFFLE( 3 , 1 , 0 , 2 ) );
    // multiply
    __m128 xa = _mm_mul_ps( ea , eb );
    // set to a[2][0][1][3] , b[1][2][0][3]
    a = _mm_shuffle_ps( a , a , _MM_SHUFFLE( 3 , 1 , 0 , 2 ) );
    b = _mm_shuffle_ps( b , b , _MM_SHUFFLE( 3 , 0 , 2 , 1 ) );
    // multiply
    __m128 xb = _mm_mul_ps( a , b );
    // subtract
    return _mm_sub_ps( xa , xb );
}
这就是shuffle强大的地方,它可以直接在寄存器里直接调整分量的顺序。而且配合_mm_movehl_ps,我们可以轻松解决点积的运算。_mm_movehl_ps把操作数高位两个分量赋予目标数的低位两分量,而目标数的高位两分量值不变,相当于:
a[0] = b[2];
a[1] = b[3];
三分量的向量求点积,可以写成:
float dot( const float& a , const float& b ) const {
    return a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
}
则用SSE intrinsic可以写成:
/* x[0] * x[1] + y[0] * y[1] + z[0] * z[1] */
__m128 _mm_dot_ps( __m128 x , __m128 y ) {
    __m128 s , r;
    s = _mm_mul_ps( x , y );
    r = _mm_add_ss( s , _mm_movehl_ps( s , s ) );
    r = _mm_add_ss( r , _mm_shuffle_ps( r , r , 1 ) );
    return r;
}
通过这两个例子,可以留意到向量内元素的垂直相加一般形式,即:
/* x[0] + x[1] + x[2] + x[3] */
__m128 _mm_sum_ps( __m128 x ) {
    __m128 r;
    r = _mm_add_ps( x , _mm_movehl_ps( x , x ) );
    r = _mm_add_ss( r , _mm_shuffle_ps( r , r , 1 ) );
    return r;
}
那么通过扩展,可以得到求向量长度的函数,首先是求分量平方和函数:
/* x[0] * x[0] + y[0] * y[0] + z[0] * z[0] */
__m128 _mm_square_ps( __m128 x ) {
    __m128 s , r;
    s = _mm_mul_ps( x , x );
    r = _mm_add_ss( s , _mm_movehl_ps( s , s ) );
    r = _mm_add_ss( r , _mm_shuffle_ps( r , r , 1 ) );
    return r;
}
然后就可以直接把结果求平方根,可得长度。解决了长度,接下来则是很重要的单位化了。可以说单位化是最重要的一个函数,它经常被调用到,而函数内的陷阱却又最多。求单位化其实并不难,就是分量除以向量长度,可以写成:
void normalize( const Vector& a ) {
    float len = a[0] * a[0] + a[1] * a[1] + a[2] * a[2];
    if( is_zero( len ) )
        return;
    len = 1 / len;
    a[0] *= len;
    a[1] *= len;
    a[2] *= len;
}
原文地址:https://www.cnblogs.com/elvisxu/p/2090830.html