遗忘河

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

首页 新随笔 联系 管理
  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 大悟终空 阅读(1512) 评论(0)  编辑 收藏 引用
只有注册用户登录后才能发表评论。