遗忘河

一个人的烦恼是因为记性太好.

首页 新随笔 联系 管理
  11 Posts :: 4 Stories :: 13 Comments :: 0 Trackbacks
inline bool LMatrix44::GetInverse(LMatrix44& matInvDest) const
{
    
// 事实上,我们绝大多数的逆操作都能成功. 
    
// 注意到我们GetDeterm() 操作中计算了一部分余子式,如果逆操作成功,那么GetDeterm就浪费了.
    
// 因此我们选用一种对以可逆为前提的优化.
    assert(&matInvDest != this);
    
float tmp[12];
    tmp[
0= m33 * m44; tmp[1= m43 * m34; tmp[2= m23 * m44; tmp[3= m24 * m43;
    tmp[
4= m23 * m34; tmp[5= m33 * m24; tmp[6= m13 * m44; tmp[7= m43 * m14;
    tmp[
8= m13 * m34; tmp[9= m33 * m14; tmp[10]= m13 * m24; tmp[11]= m23 * m14;    // # 12*

    
// 计算余子式 (倒置伴随矩阵)
    matInvDest.m11 = m22*(tmp[0- tmp[1]) + m32*(tmp[3- tmp[2]) + m42*(tmp[4- tmp[5]);  // 3* 3- 2+
    matInvDest.m12 = m12*(tmp[1- tmp[0]) + m32*(tmp[6- tmp[7]) + m42*(tmp[9- tmp[8]);
    matInvDest.m13 
= m12*(tmp[2- tmp[3]) + m22*(tmp[7- tmp[6]) + m42*(tmp[10- tmp[11]);
    matInvDest.m14 
= m12*(tmp[5- tmp[4]) + m22*(tmp[8- tmp[9]) + m32*(tmp[11- tmp[10]);
    matInvDest.m21 
= m21*(tmp[1- tmp[0]) + m31*(tmp[2- tmp[3]) + m41*(tmp[5- tmp[4]);
    matInvDest.m22 
= m11*(tmp[0- tmp[1]) + m31*(tmp[7- tmp[6]) + m41*(tmp[8- tmp[9]);
    matInvDest.m23 
= m11*(tmp[3- tmp[2]) + m21*(tmp[6- tmp[7]) + m41*(tmp[11- tmp[10]);
    matInvDest.m24 
= m11*(tmp[4- tmp[5]) + m21*(tmp[9- tmp[8]) + m31*(tmp[10- tmp[11]);
                                                                                        
//36* 24- 16+
    tmp[0= m31 * m42; tmp[1= m41 * m32; tmp[2= m21 * m42; tmp[3= m41 * m22;
    tmp[
4= m21 * m32; tmp[5= m31 * m22; tmp[6= m11 * m42; tmp[7= m41 * m12;
    tmp[
8= m11 * m32; tmp[9= m31 * m12; tmp[10]= m11 * m22; tmp[11]= m21 * m12;
    
    matInvDest.m31 
= m24*(tmp[0- tmp[1]) + m34*(tmp[3- tmp[2]) + m44*(tmp[4- tmp[5]);
    matInvDest.m32 
= m14*(tmp[1- tmp[0]) + m34*(tmp[6- tmp[7]) + m44*(tmp[9- tmp[8]);
    matInvDest.m33 
= m14*(tmp[2- tmp[3]) + m23*(tmp[7- tmp[6]) + m44*(tmp[10- tmp[11]);
    matInvDest.m34 
= m14*(tmp[5- tmp[4]) + m24*(tmp[8- tmp[9]) + m34*(tmp[11- tmp[10]);
    matInvDest.m41 
= m23*(tmp[1- tmp[0]) + m33*(tmp[2- tmp[3]) + m43*(tmp[5- tmp[4]);
    matInvDest.m42 
= m33*(tmp[7- tmp[6]) + m13*(tmp[0- tmp[1]) + m43*(tmp[8- tmp[9]);
    matInvDest.m43 
= m23*(tmp[6- tmp[7]) + m13*(tmp[3- tmp[2]) + m43*(tmp[11- tmp[10]);
    matInvDest.m44 
= m13*(tmp[4- tmp[5]) + m33*(tmp[10- tmp[11]) + m23*(tmp[9- tmp[8]);
                                                                                        
// 72* 48- 32+
    float fDet=(m11*matInvDest.m11 + m21*matInvDest.m12 + m31*matInvDest.m13 + m41*matInvDest.m14);        // 4* 3+
    if (LAbs(fDet) < 0.000001f)
    
{
        matInvDest.SetIdentity();
        
return false;
    }


    fDet 
= 1.0f/fDet;
    matInvDest.m11
*= fDet;matInvDest.m12*= fDet;matInvDest.m13*= fDet;matInvDest.m14*= fDet;
    matInvDest.m21
*= fDet;matInvDest.m22*= fDet;matInvDest.m23*= fDet;matInvDest.m24*= fDet;
    matInvDest.m31
*= fDet;matInvDest.m32*= fDet;matInvDest.m33*= fDet;matInvDest.m34*= fDet;
    matInvDest.m41
*= fDet;matInvDest.m42*= fDet;matInvDest.m43*= fDet;matInvDest.m44*= fDet;//16*
    return true;
}

以上是优化实现,以12个浮点的零时变量消耗,通过计算倒置伴随矩阵,最终获得逆矩阵. 在实现中,以矩阵可逆为前提,避免了大部分的cofactor 计算.
比常规通过求模,然后求逆要少60次乘法.

缺点: 没有形成很好的16字节排列,依旧不便于SSE/SSE2指令优化.
posted on 2007-08-15 16:20 大悟终空 阅读(1437) 评论(0)  编辑 收藏 引用
只有注册用户登录后才能发表评论。