6

我正在尝试使用 boost 进行简单的矩阵求逆运算。但我收到一个错误。基本上我想要找到的是 inversted_matrix = inverse(trans(matrix) * matrix) 但我收到一个错误

Check failed in file boost_1_53_0/boost/numeric/ublas/lu.hpp at line 299: 
detail::expression_type_check (prod (triangular_adaptor<const_matrix_type, 
upper> (m), e), cm2) 
terminate called after throwing an instance of 
'boost::numeric::ublas::internal_logic' 
  what(): internal logic 
Aborted (core dumped)

我的尝试：

#include <boost/numeric/ublas/matrix.hpp> 
#include <boost/numeric/ublas/vector.hpp> 
#include <boost/numeric/ublas/io.hpp> 
#include <boost/numeric/ublas/vector_proxy.hpp> 
#include <boost/numeric/ublas/matrix.hpp> 
#include <boost/numeric/ublas/triangular.hpp> 
#include <boost/numeric/ublas/lu.hpp> 

namespace ublas = boost::numeric::ublas; 
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; 
    // 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; 
}

int main(){ 
    using namespace boost::numeric::ublas; 
    matrix<double> m(4,5); 
    vector<double> v(4); 
    vector<double> thetas; 
    m(0,0) = 1; m(0,1) = 2104; m(0,2) = 5; m(0,3) = 1;m(0,4) = 45; 
    m(1,0) = 1; m(1,1) = 1416; m(1,2) = 3; m(1,3) = 2;m(1,4) = 40; 
    m(2,0) = 1; m(2,1) = 1534; m(2,2) = 3; m(2,3) = 2;m(2,4) = 30; 
    m(3,0) = 1; m(3,1) = 852; m(3,2) = 2; m(3,3) = 1;m(3,4) = 36; 
    std::cout<<m<<std::endl; 
    matrix<double> product = prod(trans(m), m); 
    std::cout<<product<<std::endl; 
    matrix<double> inversion(5,5); 
    bool inverted; 
    inverted = InvertMatrix(product, inversion); 
    std::cout << inversion << std::endl; 
}
4

3 回答 3

7

Boost Ublas 具有运行时检查，以确保数值稳定性等。如果您查看错误的来源，您会看到它试图确保 U*X = B, X = U^-1*B, U*X = B（或类似的东西）在某个 epsilon 内是 coorect 的。如果你在数值上的偏差太大，这可能不会成立。

您可以通过禁用检查-DBOOST_UBLAS_NDEBUG或使用BOOST_UBLAS_TYPE_CHECK_EPSILON, BOOST_UBLAS_TYPE_CHECK_MIN.

于 2013-06-21T18:24:36.977 回答
5

As m has only 4 rows, prod(trans(m), m) cannot have a rank higher than 4, and as the product is a 5x5 matrix, it must be singular (i.e. it has determinant 0) and calculating the inverse of a singular matrix is like division by 0. Add independent rows to m to solve this singularity problem.

于 2013-10-30T10:28:03.730 回答
0

我认为您的矩阵尺寸 4 x 5 导致了错误。就像Maarten Hilferink 提到的那样，您可以尝试使用 5 x 5 之类的方阵。以下是逆矩阵的要求：

  1. 矩阵必须是正方形的（行数和列数相同）。
  2. 矩阵的行列式不能为零（行列式在第 6.4 节中介绍）。这不是实数不为零以具有逆，行列式必须不为零以具有逆。
于 2022-02-26T00:24:38.143 回答