行列式计算(模拟递归)

typedef float MatEle;

typedef struct

{

    MatEle* ele;

    int row,col;

}SMatrix;

 

/*//////// n×n矩阵对应行列式的值(非递归算法) ////////////////////////////*/
/*
    使用了I[]模拟堆栈，从1、2、3...n逐渐递增到 n、n-1、...2、1
    使用Vs[i]保存前i个的相乘积，避免了重复的计算，
    比如由1、2、3、4 变为 1、2、4、3 ，只需要计算：
    Vs[3]=Vs[2]*a[3][4]； Vs[4]=Vs[3]*a[4][3]
    而不重新由a[1][1]*a[2][2]*a[3][4]*a[4][3]计算。
*/

MatEle HLSAbs(SMatrix mat)
{
    int    i,j,k,t,N,n;
    MatEle rv;        /* 返回值 */
    int    *I, minus; /* 全排列、负号标识 */
    MatEle *Vs,*V;    /* 数据暂存、 指向全排列的值*/

    n = mat.row;
    N = n + 1;
    I = (int*)    malloc( sizeof(I)*N      ); /* 多出的I[0]作为末尾标识 */
    Vs= (MatEle*) malloc( sizeof(MatEle)*N );
    for(i=0; i<N; i++) I[i] = i-1;
    Vs[0] = 1;

    V = Vs+n;     /* 令V指向全排列相乘的值 */
    i = 1;
    rv = 0;
    minus = 0;
    while(1)
    {
        for( ; i<N; i++)
        {
            Vs[i] = Vs[i-1]*mat.ele[(i-1)*n+I[i]];
        }
        if(minus) rv -= *V;
        else rv += *V;
        i = n;
        while( I[i] < I[--i] ) ;

/* 如果i==0，说明全排列为 (n,n-1,...,2,1)，计算完毕，返回；*/
        if(i<1)
        {
            free(I);
            free(Vs);
            return rv;
        }

        for(k=n; I[k]<I[i]; k--);

        t    = I[i];
        I[i] = I[k];
        I[k] = t;

        for(j=i+1, k=n ; j<k ; j++, k--)
        { /* I[j]和I[k]之间必然是由小到大排列，此举在于使其由大到小排列 */
            t=I[j];
            I[j]=I[k];
            I[k]=t;
        }

  minus ^= ((n-i)>>1^0x01) ;
    }
}