GitHub-Proxy/eigen
1
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2025-08-12 19:59:05 +08:00
Code Issues Packages Projects Releases Wiki Activity

fix stuff after the PermutationMatrix changes.

Browse Source 
I still have JacobiSVD errors when cols>rows
This commit is contained in:
Benoit Jacob 2009-11-18 12:41:24 -05:00
parent 94c706d04f
commit 8860203e6a
4 changed files with 13 additions and 10 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
Eigen/src/QR/ColPivHouseholderQR.h
View File

@ -314,7 +314,7 @@ ColPivHouseholderQR<MatrixType>& ColPivHouseholderQR<MatrixType>::compute(const
m_cols_permutation.setIdentity(cols); m_cols_permutation.setIdentity(cols);
for(int k = 0; k < size; ++k) for(int k = 0; k < size; ++k)
m_cols_permutation.applyTranspositionOnTheLeft(k, cols_transpositions.coeff(k)); m_cols_permutation.applyTranspositionOnTheRight(k, cols_transpositions.coeff(k));
m_det_pq = (number_of_transpositions%2) ? -1 : 1; m_det_pq = (number_of_transpositions%2) ? -1 : 1;
m_isInitialized = true; m_isInitialized = true;

9
Eigen/src/SVD/JacobiSVD.h
View File

@ -236,9 +236,8 @@ struct ei_svd_precondition_if_more_rows_than_cols<MatrixType, Options, true>
FullPivHouseholderQR<MatrixType> qr(matrix); FullPivHouseholderQR<MatrixType> qr(matrix);
work_matrix = qr.matrixQR().block(0,0,diagSize,diagSize).template triangularView<UpperTriangular>(); work_matrix = qr.matrixQR().block(0,0,diagSize,diagSize).template triangularView<UpperTriangular>();
if(ComputeU) svd.m_matrixU = qr.matrixQ(); if(ComputeU) svd.m_matrixU = qr.matrixQ();
if(ComputeV) if(ComputeV) svd.m_matrixV = qr.colsPermutation();
for(int i = 0; i < cols; i++)
svd.m_matrixV.coeffRef(qr.colsPermutation().coeff(i),i) = Scalar(1);
return true; return true;
} }
else return false; else return false;
@ -281,9 +280,7 @@ struct ei_svd_precondition_if_more_cols_than_rows<MatrixType, Options, true>
FullPivHouseholderQR<TransposeTypeWithSameStorageOrder> qr(matrix.adjoint()); FullPivHouseholderQR<TransposeTypeWithSameStorageOrder> qr(matrix.adjoint());
work_matrix = qr.matrixQR().block(0,0,diagSize,diagSize).template triangularView<UpperTriangular>().adjoint(); work_matrix = qr.matrixQR().block(0,0,diagSize,diagSize).template triangularView<UpperTriangular>().adjoint();
if(ComputeV) svd.m_matrixV = qr.matrixQ(); if(ComputeV) svd.m_matrixV = qr.matrixQ();
if(ComputeU) if(ComputeU) svd.m_matrixU = qr.colsPermutation();
for(int i = 0; i < rows; i++)
svd.m_matrixU.coeffRef(qr.colsPermutation().coeff(i),i) = Scalar(1);
return true; return true;
} }
else return false; else return false;

3
test/jacobisvd.cpp
View File

@ -54,6 +54,9 @@ template<typename MatrixType, unsigned int Options> void svd(const MatrixType& m
MatrixUType u = svd.matrixU(); MatrixUType u = svd.matrixU();
MatrixVType v = svd.matrixV(); MatrixVType v = svd.matrixV();
std::cout << "a\n" << a << std::endl;
std::cout << "b\n" << u * sigma * v.adjoint() << std::endl;
VERIFY_IS_APPROX(a, u * sigma * v.adjoint()); VERIFY_IS_APPROX(a, u * sigma * v.adjoint());
VERIFY_IS_UNITARY(u); VERIFY_IS_UNITARY(u);
VERIFY_IS_UNITARY(v); VERIFY_IS_UNITARY(v);

9
test/qr_colpivoting.cpp
View File

@ -28,7 +28,8 @@
template<typename MatrixType> void qr() template<typename MatrixType> void qr()
{ {
int rows = ei_random<int>(20,200), cols = ei_random<int>(20,200), cols2 = ei_random<int>(20,200); // int rows = ei_random<int>(20,200), cols = ei_random<int>(20,200), cols2 = ei_random<int>(20,200);
int rows=3, cols=3, cols2=3;
int rank = ei_random<int>(1, std::min(rows, cols)-1); int rank = ei_random<int>(1, std::min(rows, cols)-1);
typedef typename MatrixType::Scalar Scalar; typedef typename MatrixType::Scalar Scalar;
@ -55,7 +56,9 @@ template<typename MatrixType> void qr()
MatrixType c = MatrixType::Zero(rows,cols); MatrixType c = MatrixType::Zero(rows,cols);
for(int i = 0; i < cols; ++i) c.col(qr.colsPermutation().coeff(i)) = b.col(i); for(int i = 0; i < cols; ++i) c.col(qr.colsPermutation().indices().coeff(i)) = b.col(i);
std::cout << "m1:\n" << m1 << std::endl;
std::cout << "c:\n" << c << std::endl;
VERIFY_IS_APPROX(m1, c); VERIFY_IS_APPROX(m1, c);
MatrixType m2 = MatrixType::Random(cols,cols2); MatrixType m2 = MatrixType::Random(cols,cols2);
@ -87,7 +90,7 @@ template<typename MatrixType, int Cols2> void qr_fixedsize()
Matrix<Scalar,Rows,Cols> c = MatrixType::Zero(Rows,Cols); Matrix<Scalar,Rows,Cols> c = MatrixType::Zero(Rows,Cols);
for(int i = 0; i < Cols; ++i) c.col(qr.colsPermutation().coeff(i)) = b.col(i); for(int i = 0; i < Cols; ++i) c.col(qr.colsPermutation().indices().coeff(i)) = b.col(i);
VERIFY_IS_APPROX(m1, c); VERIFY_IS_APPROX(m1, c);
Matrix<Scalar,Cols,Cols2> m2 = Matrix<Scalar,Cols,Cols2>::Random(Cols,Cols2); Matrix<Scalar,Cols,Cols2> m2 = Matrix<Scalar,Cols,Cols2>::Random(Cols,Cols2);