Boost Library,如何从lu_factorize()获得决定因素？

Question

我试图使用boost c ++库计算一个行列式.我找到了我在下面复制的函数InvertMatrix()的代码.每次我计算这个逆,我也想要行列式.我很清楚如何通过将U矩阵的对角线与LU分解相乘来计算.有一个问题,我能够正确计算行列式,除了符号.根据旋转,我在一半时间内得到的标志不正确.有没有人建议如何每次都正确的标志？提前致谢.

template<class T>
bool InvertMatrix(const ublas::matrix<T>& input, ublas::matrix<T>& inverse)
{
 using namespace boost::numeric::ublas;
 typedef permutation_matrix<std::size_t> pmatrix;
 // create a working copy of the input
 matrix<T> A(input);
 // create a permutation matrix for the LU-factorization
 pmatrix pm(A.size1());

 // perform LU-factorization
 int res = lu_factorize(A,pm);
 if( res != 0 ) return false;

这是我在计算行列式时插入最佳镜头的地方.

 T determinant = 1;

 for(int i = 0; i < A.size1(); i++)
 {
  determinant *= A(i,i);
 }

结束我的部分代码.

 // create identity matrix of "inverse"
 inverse.assign(ublas::identity_matrix<T>(A.size1()));

 // backsubstitute to get the inverse
 lu_substitute(A, pm, inverse);

 return true;
}

Answer 1

置换矩阵pm包含确定符号变化所需的信息：您需要将行列式乘以置换矩阵的行列式。

仔细阅读源文件，lu.hpp我们发现一个名为 的函数swap_rows，它告诉我们如何将置换矩阵应用于矩阵。考虑到每个实际交换贡献的因子为 -1，它很容易修改以产生排列矩阵的行列式（排列的符号）：

template <typename size_type, typename A>
int determinant(const permutation_matrix<size_type,A>& pm)
{
    int pm_sign=1;
    size_type size=pm.size();
    for (size_type i = 0; i < size; ++i)
        if (i != pm(i))
            pm_sign* = -1; // swap_rows would swap a pair of rows here, so we change sign
    return pm_sign;
}

另一种选择是使用不进行任何旋转的lu_factorize和方法（查阅源代码，但基本上删除对和 的调用）。该更改将使您的行列式计算按原样进行。但要小心：删除旋转将使算法在数值上不太稳定。lu_substitutepmlu_factorizelu_substitute