Thư viện tăng cường, làm thế nào để có được yếu tố quyết định từ lu_factorize()?

Question

Tôi đang cố gắng tính toán một yếu tố quyết định bằng cách sử dụng thư viện boost C++. Tôi tìm thấy mã cho hàm InvertMatrix() mà tôi đã sao chép bên dưới. Mỗi khi tôi tính toán nghịch đảo này, tôi cũng muốn định thức. Tôi có một ý tưởng tốt như thế nào để tính toán, bằng cách nhân xuống đường chéo của ma trận U từ sự phân hủy LU. Có một vấn đề, tôi có thể tính toán định thức đúng cách, ngoại trừ dấu hiệu. Tùy thuộc vào việc xoay vòng tôi nhận được dấu hiệu không chính xác một nửa thời gian. Có ai có một gợi ý về làm thế nào để có được những dấu hiệu đúng mỗi lần? Cảm ơn trước.

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; 

Đây là nơi tôi chèn ảnh tốt nhất của mình khi tính toán yếu tố quyết định.

T determinant = 1; 

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

Kết thúc phần mã của tô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

Các hoán vị ma trận pm chứa các thông tin mà bạn sẽ cần phải xác định sự thay đổi dấu: bạn sẽ muốn nhân tố quyết định của bạn bằng cách định thức của ma trận hoán vị.

Perusing tệp nguồn lu.hpp chúng tôi tìm thấy hàm có tên swap_rows cho biết cách áp dụng ma trận hoán vị cho ma trận. Nó dễ dàng sửa đổi để mang yếu tố quyết định của ma trận hoán vị (dấu hiệu của sự hoán vị), cho rằng mỗi hoán đổi thực tế đóng góp một yếu tố của -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; 
} 

Một lựa chọn khác là sử dụng các phương pháp lu_factorize và lu_substitute mà không thực hiện bất kỳ sự xoay vòng nào (tham khảo ý kiến ​​nguồn, nhưng về cơ bản, hãy thả pm trong các cuộc gọi đến lu_factorize và lu_substitute). Sự thay đổi đó sẽ làm cho công việc tính toán quyết định của bạn như hiện tại. Hãy cẩn thận, tuy nhiên: loại bỏ xoay vòng sẽ làm cho thuật toán ít ổn định về số lượng.

Answer 2

Theo http://qiangsong.wordpress.com/2011/07/16/lu-factorisation-in-ublas/:

Chỉ cần thay đổi determinant *= A(i,i) để determinant *= (pm(i) == i ? 1 : -1) * A(i,i). Tôi đã thử cách này và nó hoạt động.

Tôi biết, nó thực sự rất giống với câu trả lời của Managu và ý tưởng là như nhau, nhưng tôi tin rằng nó đơn giản hơn (và "2 trong 1" nếu được sử dụng trong chức năng InvertMatrix).