11#ifndef EIGEN_JACOBISVD_H
12#define EIGEN_JACOBISVD_H
14#include "./InternalHeaderCheck.h"
21template<
typename MatrixType,
int QRPreconditioner,
22 bool IsComplex = NumTraits<typename MatrixType::Scalar>::IsComplex>
23struct svd_precondition_2x2_block_to_be_real {};
32enum { PreconditionIfMoreColsThanRows, PreconditionIfMoreRowsThanCols };
34template<
typename MatrixType,
int QRPreconditioner,
int Case>
35struct qr_preconditioner_should_do_anything
37 enum { a = MatrixType::RowsAtCompileTime !=
Dynamic &&
38 MatrixType::ColsAtCompileTime !=
Dynamic &&
39 MatrixType::ColsAtCompileTime <= MatrixType::RowsAtCompileTime,
40 b = MatrixType::RowsAtCompileTime !=
Dynamic &&
41 MatrixType::ColsAtCompileTime !=
Dynamic &&
42 MatrixType::RowsAtCompileTime <= MatrixType::ColsAtCompileTime,
44 (Case == PreconditionIfMoreColsThanRows &&
bool(a)) ||
45 (Case == PreconditionIfMoreRowsThanCols &&
bool(b)) )
49template<
typename MatrixType,
int QRPreconditioner,
int Case,
50 bool DoAnything = qr_preconditioner_should_do_anything<MatrixType, QRPreconditioner, Case>::ret
51>
struct qr_preconditioner_impl {};
53template<
typename MatrixType,
int QRPreconditioner,
int Case>
54class qr_preconditioner_impl<MatrixType, QRPreconditioner, Case, false>
57 void allocate(
const JacobiSVD<MatrixType, QRPreconditioner>&) {}
58 bool run(JacobiSVD<MatrixType, QRPreconditioner>&,
const MatrixType&)
66template<
typename MatrixType>
70 typedef typename MatrixType::Scalar Scalar;
73 RowsAtCompileTime = MatrixType::RowsAtCompileTime,
74 MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime
76 typedef Matrix<Scalar, 1, RowsAtCompileTime, RowMajor, 1, MaxRowsAtCompileTime> WorkspaceType;
78 void allocate(
const JacobiSVD<MatrixType, FullPivHouseholderQRPreconditioner>& svd)
80 if (svd.rows() != m_qr.rows() || svd.cols() != m_qr.cols())
83 ::new (&m_qr) QRType(svd.rows(), svd.cols());
85 if (svd.m_computeFullU) m_workspace.resize(svd.rows());
88 bool run(JacobiSVD<MatrixType, FullPivHouseholderQRPreconditioner>& svd,
const MatrixType& matrix)
90 if(matrix.rows() > matrix.cols())
93 svd.m_workMatrix = m_qr.matrixQR().block(0,0,matrix.cols(),matrix.cols()).template triangularView<Upper>();
94 if(svd.m_computeFullU) m_qr.matrixQ().evalTo(svd.m_matrixU, m_workspace);
95 if(svd.computeV()) svd.m_matrixV = m_qr.colsPermutation();
101 typedef FullPivHouseholderQR<MatrixType> QRType;
103 WorkspaceType m_workspace;
106template<
typename MatrixType>
110 typedef typename MatrixType::Scalar Scalar;
113 RowsAtCompileTime = MatrixType::RowsAtCompileTime,
114 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
115 MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
116 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
117 Options = MatrixType::Options
120 typedef typename internal::make_proper_matrix_type<
121 Scalar, ColsAtCompileTime, RowsAtCompileTime, Options, MaxColsAtCompileTime, MaxRowsAtCompileTime
122 >::type TransposeTypeWithSameStorageOrder;
124 void allocate(
const JacobiSVD<MatrixType, FullPivHouseholderQRPreconditioner>& svd)
126 if (svd.cols() != m_qr.rows() || svd.rows() != m_qr.cols())
129 ::new (&m_qr) QRType(svd.cols(), svd.rows());
131 m_adjoint.resize(svd.cols(), svd.rows());
132 if (svd.m_computeFullV) m_workspace.resize(svd.cols());
135 bool run(JacobiSVD<MatrixType, FullPivHouseholderQRPreconditioner>& svd,
const MatrixType& matrix)
137 if(matrix.cols() > matrix.rows())
139 m_adjoint = matrix.adjoint();
140 m_qr.compute(m_adjoint);
141 svd.m_workMatrix = m_qr.matrixQR().block(0,0,matrix.rows(),matrix.rows()).template triangularView<Upper>().adjoint();
142 if(svd.m_computeFullV) m_qr.matrixQ().evalTo(svd.m_matrixV, m_workspace);
143 if(svd.computeU()) svd.m_matrixU = m_qr.colsPermutation();
149 typedef FullPivHouseholderQR<TransposeTypeWithSameStorageOrder> QRType;
151 TransposeTypeWithSameStorageOrder m_adjoint;
152 typename internal::plain_row_type<MatrixType>::type m_workspace;
157template<
typename MatrixType>
161 void allocate(
const JacobiSVD<MatrixType, ColPivHouseholderQRPreconditioner>& svd)
163 if (svd.rows() != m_qr.rows() || svd.cols() != m_qr.cols())
166 ::new (&m_qr) QRType(svd.rows(), svd.cols());
168 if (svd.m_computeFullU) m_workspace.resize(svd.rows());
169 else if (svd.m_computeThinU) m_workspace.resize(svd.cols());
172 bool run(JacobiSVD<MatrixType, ColPivHouseholderQRPreconditioner>& svd,
const MatrixType& matrix)
174 if(matrix.rows() > matrix.cols())
176 m_qr.compute(matrix);
177 svd.m_workMatrix = m_qr.matrixQR().block(0,0,matrix.cols(),matrix.cols()).template triangularView<Upper>();
178 if(svd.m_computeFullU) m_qr.householderQ().evalTo(svd.m_matrixU, m_workspace);
179 else if(svd.m_computeThinU)
181 svd.m_matrixU.setIdentity(matrix.rows(), matrix.cols());
182 m_qr.householderQ().applyThisOnTheLeft(svd.m_matrixU, m_workspace);
184 if(svd.computeV()) svd.m_matrixV = m_qr.colsPermutation();
191 typedef ColPivHouseholderQR<MatrixType> QRType;
193 typename internal::plain_col_type<MatrixType>::type m_workspace;
196template<
typename MatrixType>
200 typedef typename MatrixType::Scalar Scalar;
203 RowsAtCompileTime = MatrixType::RowsAtCompileTime,
204 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
205 MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
206 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
207 Options = MatrixType::Options
210 typedef typename internal::make_proper_matrix_type<
211 Scalar, ColsAtCompileTime, RowsAtCompileTime, Options, MaxColsAtCompileTime, MaxRowsAtCompileTime
212 >::type TransposeTypeWithSameStorageOrder;
214 void allocate(
const JacobiSVD<MatrixType, ColPivHouseholderQRPreconditioner>& svd)
216 if (svd.cols() != m_qr.rows() || svd.rows() != m_qr.cols())
219 ::new (&m_qr) QRType(svd.cols(), svd.rows());
221 if (svd.m_computeFullV) m_workspace.resize(svd.cols());
222 else if (svd.m_computeThinV) m_workspace.resize(svd.rows());
223 m_adjoint.resize(svd.cols(), svd.rows());
226 bool run(JacobiSVD<MatrixType, ColPivHouseholderQRPreconditioner>& svd,
const MatrixType& matrix)
228 if(matrix.cols() > matrix.rows())
230 m_adjoint = matrix.adjoint();
231 m_qr.compute(m_adjoint);
233 svd.m_workMatrix = m_qr.matrixQR().block(0,0,matrix.rows(),matrix.rows()).template triangularView<Upper>().adjoint();
234 if(svd.m_computeFullV) m_qr.householderQ().evalTo(svd.m_matrixV, m_workspace);
235 else if(svd.m_computeThinV)
237 svd.m_matrixV.setIdentity(matrix.cols(), matrix.rows());
238 m_qr.householderQ().applyThisOnTheLeft(svd.m_matrixV, m_workspace);
240 if(svd.computeU()) svd.m_matrixU = m_qr.colsPermutation();
247 typedef ColPivHouseholderQR<TransposeTypeWithSameStorageOrder> QRType;
249 TransposeTypeWithSameStorageOrder m_adjoint;
250 typename internal::plain_row_type<MatrixType>::type m_workspace;
255template<
typename MatrixType>
259 void allocate(
const JacobiSVD<MatrixType, HouseholderQRPreconditioner>& svd)
261 if (svd.rows() != m_qr.rows() || svd.cols() != m_qr.cols())
264 ::new (&m_qr) QRType(svd.rows(), svd.cols());
266 if (svd.m_computeFullU) m_workspace.resize(svd.rows());
267 else if (svd.m_computeThinU) m_workspace.resize(svd.cols());
270 bool run(JacobiSVD<MatrixType, HouseholderQRPreconditioner>& svd,
const MatrixType& matrix)
272 if(matrix.rows() > matrix.cols())
274 m_qr.compute(matrix);
275 svd.m_workMatrix = m_qr.matrixQR().block(0,0,matrix.cols(),matrix.cols()).template triangularView<Upper>();
276 if(svd.m_computeFullU) m_qr.householderQ().evalTo(svd.m_matrixU, m_workspace);
277 else if(svd.m_computeThinU)
279 svd.m_matrixU.setIdentity(matrix.rows(), matrix.cols());
280 m_qr.householderQ().applyThisOnTheLeft(svd.m_matrixU, m_workspace);
282 if(svd.computeV()) svd.m_matrixV.setIdentity(matrix.cols(), matrix.cols());
288 typedef HouseholderQR<MatrixType> QRType;
290 typename internal::plain_col_type<MatrixType>::type m_workspace;
293template<
typename MatrixType>
297 typedef typename MatrixType::Scalar Scalar;
300 RowsAtCompileTime = MatrixType::RowsAtCompileTime,
301 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
302 MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
303 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
304 Options = MatrixType::Options
307 typedef typename internal::make_proper_matrix_type<
308 Scalar, ColsAtCompileTime, RowsAtCompileTime, Options, MaxColsAtCompileTime, MaxRowsAtCompileTime
309 >::type TransposeTypeWithSameStorageOrder;
311 void allocate(
const JacobiSVD<MatrixType, HouseholderQRPreconditioner>& svd)
313 if (svd.cols() != m_qr.rows() || svd.rows() != m_qr.cols())
316 ::new (&m_qr) QRType(svd.cols(), svd.rows());
318 if (svd.m_computeFullV) m_workspace.resize(svd.cols());
319 else if (svd.m_computeThinV) m_workspace.resize(svd.rows());
320 m_adjoint.resize(svd.cols(), svd.rows());
323 bool run(JacobiSVD<MatrixType, HouseholderQRPreconditioner>& svd,
const MatrixType& matrix)
325 if(matrix.cols() > matrix.rows())
327 m_adjoint = matrix.adjoint();
328 m_qr.compute(m_adjoint);
330 svd.m_workMatrix = m_qr.matrixQR().block(0,0,matrix.rows(),matrix.rows()).template triangularView<Upper>().adjoint();
331 if(svd.m_computeFullV) m_qr.householderQ().evalTo(svd.m_matrixV, m_workspace);
332 else if(svd.m_computeThinV)
334 svd.m_matrixV.setIdentity(matrix.cols(), matrix.rows());
335 m_qr.householderQ().applyThisOnTheLeft(svd.m_matrixV, m_workspace);
337 if(svd.computeU()) svd.m_matrixU.setIdentity(matrix.rows(), matrix.rows());
344 typedef HouseholderQR<TransposeTypeWithSameStorageOrder> QRType;
346 TransposeTypeWithSameStorageOrder m_adjoint;
347 typename internal::plain_row_type<MatrixType>::type m_workspace;
355template<
typename MatrixType,
int QRPreconditioner>
356struct svd_precondition_2x2_block_to_be_real<MatrixType, QRPreconditioner, false>
358 typedef JacobiSVD<MatrixType, QRPreconditioner> SVD;
359 typedef typename MatrixType::RealScalar RealScalar;
360 static bool run(
typename SVD::WorkMatrixType&, SVD&,
Index,
Index, RealScalar&) {
return true; }
363template<
typename MatrixType,
int QRPreconditioner>
364struct svd_precondition_2x2_block_to_be_real<MatrixType, QRPreconditioner, true>
366 typedef JacobiSVD<MatrixType, QRPreconditioner> SVD;
367 typedef typename MatrixType::Scalar Scalar;
368 typedef typename MatrixType::RealScalar RealScalar;
369 static bool run(
typename SVD::WorkMatrixType& work_matrix, SVD& svd,
Index p,
Index q, RealScalar& maxDiagEntry)
375 RealScalar n =
sqrt(numext::abs2(work_matrix.coeff(p,p)) + numext::abs2(work_matrix.coeff(q,p)));
377 const RealScalar considerAsZero = (std::numeric_limits<RealScalar>::min)();
383 work_matrix.coeffRef(p,p) = work_matrix.coeffRef(q,p) =
Scalar(0);
385 if(
abs(numext::imag(work_matrix.coeff(p,q)))>considerAsZero)
388 z =
abs(work_matrix.coeff(p,q)) / work_matrix.coeff(p,q);
389 work_matrix.row(p) *= z;
390 if(svd.computeU()) svd.m_matrixU.col(p) *=
conj(z);
392 if(
abs(numext::imag(work_matrix.coeff(q,q)))>considerAsZero)
394 z =
abs(work_matrix.coeff(q,q)) / work_matrix.coeff(q,q);
395 work_matrix.row(q) *= z;
396 if(svd.computeU()) svd.m_matrixU.col(q) *=
conj(z);
402 rot.c() =
conj(work_matrix.coeff(p,p)) / n;
403 rot.s() = work_matrix.coeff(q,p) / n;
404 work_matrix.applyOnTheLeft(p,q,rot);
405 if(svd.computeU()) svd.m_matrixU.applyOnTheRight(p,q,rot.
adjoint());
406 if(
abs(numext::imag(work_matrix.coeff(p,q)))>considerAsZero)
408 z =
abs(work_matrix.coeff(p,q)) / work_matrix.coeff(p,q);
409 work_matrix.col(q) *= z;
410 if(svd.computeV()) svd.m_matrixV.col(q) *= z;
412 if(
abs(numext::imag(work_matrix.coeff(q,q)))>considerAsZero)
414 z =
abs(work_matrix.coeff(q,q)) / work_matrix.coeff(q,q);
415 work_matrix.row(q) *= z;
416 if(svd.computeU()) svd.m_matrixU.col(q) *=
conj(z);
421 maxDiagEntry = numext::maxi<RealScalar>(maxDiagEntry,numext::maxi<RealScalar>(
abs(work_matrix.coeff(p,p)),
abs(work_matrix.coeff(q,q))));
423 RealScalar threshold = numext::maxi<RealScalar>(considerAsZero, precision * maxDiagEntry);
424 return abs(work_matrix.coeff(p,q))>threshold ||
abs(work_matrix.coeff(q,p)) > threshold;
428template<
typename MatrixType_,
int QRPreconditioner>
429struct traits<JacobiSVD<MatrixType_,QRPreconditioner> >
430 : traits<MatrixType_>
432 typedef MatrixType_ MatrixType;
490template<
typename MatrixType_,
int QRPreconditioner>
class JacobiSVD
491 :
public SVDBase<JacobiSVD<MatrixType_,QRPreconditioner> >
496 typedef MatrixType_ MatrixType;
497 typedef typename MatrixType::Scalar Scalar;
498 typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar;
500 RowsAtCompileTime = MatrixType::RowsAtCompileTime,
501 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
502 DiagSizeAtCompileTime = internal::min_size_prefer_dynamic(RowsAtCompileTime,ColsAtCompileTime),
503 MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
504 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
505 MaxDiagSizeAtCompileTime = internal::min_size_prefer_fixed(MaxRowsAtCompileTime,MaxColsAtCompileTime),
506 MatrixOptions = MatrixType::Options
509 typedef typename Base::MatrixUType MatrixUType;
510 typedef typename Base::MatrixVType MatrixVType;
511 typedef typename Base::SingularValuesType SingularValuesType;
513 typedef typename internal::plain_row_type<MatrixType>::type RowType;
514 typedef typename internal::plain_col_type<MatrixType>::type ColType;
515 typedef Matrix<Scalar, DiagSizeAtCompileTime, DiagSizeAtCompileTime,
516 MatrixOptions, MaxDiagSizeAtCompileTime, MaxDiagSizeAtCompileTime>
536 allocate(
rows,
cols, computationOptions);
549 explicit JacobiSVD(
const MatrixType& matrix,
unsigned int computationOptions = 0)
551 compute(matrix, computationOptions);
564 JacobiSVD&
compute(
const MatrixType& matrix,
unsigned int computationOptions);
574 return compute(matrix, m_computationOptions);
587 using Base::m_matrixU;
588 using Base::m_matrixV;
589 using Base::m_singularValues;
591 using Base::m_isInitialized;
592 using Base::m_isAllocated;
593 using Base::m_usePrescribedThreshold;
594 using Base::m_computeFullU;
595 using Base::m_computeThinU;
596 using Base::m_computeFullV;
597 using Base::m_computeThinV;
598 using Base::m_computationOptions;
599 using Base::m_nonzeroSingularValues;
602 using Base::m_diagSize;
603 using Base::m_prescribedThreshold;
604 WorkMatrixType m_workMatrix;
606 template<
typename MatrixType__,
int QRPreconditioner_,
bool IsComplex_>
607 friend struct internal::svd_precondition_2x2_block_to_be_real;
608 template<
typename MatrixType__,
int QRPreconditioner_,
int Case_,
bool DoAnything_>
609 friend struct internal::qr_preconditioner_impl;
611 internal::qr_preconditioner_impl<MatrixType, QRPreconditioner, internal::PreconditionIfMoreColsThanRows> m_qr_precond_morecols;
612 internal::qr_preconditioner_impl<MatrixType, QRPreconditioner, internal::PreconditionIfMoreRowsThanCols> m_qr_precond_morerows;
613 MatrixType m_scaledMatrix;
616template<
typename MatrixType,
int QRPreconditioner>
617void JacobiSVD<MatrixType, QRPreconditioner>::allocate(
Eigen::Index rows,
Eigen::Index cols,
unsigned int computationOptions)
619 eigen_assert(rows >= 0 && cols >= 0);
624 computationOptions == m_computationOptions)
632 m_isInitialized =
false;
633 m_isAllocated =
true;
634 m_computationOptions = computationOptions;
635 m_computeFullU = (computationOptions &
ComputeFullU) != 0;
636 m_computeThinU = (computationOptions &
ComputeThinU) != 0;
637 m_computeFullV = (computationOptions &
ComputeFullV) != 0;
638 m_computeThinV = (computationOptions &
ComputeThinV) != 0;
639 eigen_assert(!(m_computeFullU && m_computeThinU) &&
"JacobiSVD: you can't ask for both full and thin U");
640 eigen_assert(!(m_computeFullV && m_computeThinV) &&
"JacobiSVD: you can't ask for both full and thin V");
641 eigen_assert(internal::check_implication(m_computeThinU || m_computeThinV, MatrixType::ColsAtCompileTime==
Dynamic) &&
642 "JacobiSVD: thin U and V are only available when your matrix has a dynamic number of columns.");
645 eigen_assert(!(m_computeThinU || m_computeThinV) &&
646 "JacobiSVD: can't compute thin U or thin V with the FullPivHouseholderQR preconditioner. "
647 "Use the ColPivHouseholderQR preconditioner instead.");
649 m_diagSize = (std::min)(m_rows, m_cols);
650 m_singularValues.resize(m_diagSize);
652 m_matrixU.resize(m_rows, m_computeFullU ? m_rows
653 : m_computeThinU ? m_diagSize
656 m_matrixV.resize(m_cols, m_computeFullV ? m_cols
657 : m_computeThinV ? m_diagSize
659 m_workMatrix.resize(m_diagSize, m_diagSize);
661 if(m_cols>m_rows) m_qr_precond_morecols.allocate(*
this);
662 if(m_rows>m_cols) m_qr_precond_morerows.allocate(*
this);
663 if(m_rows!=m_cols) m_scaledMatrix.resize(rows,cols);
666template<
typename MatrixType,
int QRPreconditioner>
667JacobiSVD<MatrixType, QRPreconditioner>&
671 allocate(matrix.rows(), matrix.cols(), computationOptions);
678 const RealScalar considerAsZero = (std::numeric_limits<RealScalar>::min)();
681 RealScalar scale = matrix.cwiseAbs().template maxCoeff<PropagateNaN>();
682 if (!(numext::isfinite)(scale)) {
683 m_isInitialized =
true;
687 if(scale==RealScalar(0)) scale = RealScalar(1);
693 m_scaledMatrix = matrix / scale;
694 m_qr_precond_morecols.run(*
this, m_scaledMatrix);
695 m_qr_precond_morerows.run(*
this, m_scaledMatrix);
699 m_workMatrix = matrix.block(0,0,m_diagSize,m_diagSize) / scale;
700 if(m_computeFullU) m_matrixU.setIdentity(m_rows,m_rows);
701 if(m_computeThinU) m_matrixU.setIdentity(m_rows,m_diagSize);
702 if(m_computeFullV) m_matrixV.setIdentity(m_cols,m_cols);
703 if(m_computeThinV) m_matrixV.setIdentity(m_cols, m_diagSize);
707 RealScalar maxDiagEntry = m_workMatrix.cwiseAbs().diagonal().maxCoeff();
709 bool finished =
false;
716 for(
Index p = 1; p < m_diagSize; ++p)
718 for(
Index q = 0; q < p; ++q)
723 RealScalar threshold = numext::maxi<RealScalar>(considerAsZero, precision * maxDiagEntry);
724 if(
abs(m_workMatrix.coeff(p,q))>threshold ||
abs(m_workMatrix.coeff(q,p)) > threshold)
729 if(internal::svd_precondition_2x2_block_to_be_real<MatrixType, QRPreconditioner>::run(m_workMatrix, *
this, p, q, maxDiagEntry))
732 internal::real_2x2_jacobi_svd(m_workMatrix, p, q, &j_left, &j_right);
735 m_workMatrix.applyOnTheLeft(p,q,j_left);
736 if(computeU()) m_matrixU.applyOnTheRight(p,q,j_left.
transpose());
738 m_workMatrix.applyOnTheRight(p,q,j_right);
739 if(computeV()) m_matrixV.applyOnTheRight(p,q,j_right);
742 maxDiagEntry = numext::maxi<RealScalar>(maxDiagEntry,numext::maxi<RealScalar>(
abs(m_workMatrix.coeff(p,p)),
abs(m_workMatrix.coeff(q,q))));
751 for(
Index i = 0; i < m_diagSize; ++i)
758 RealScalar a =
abs(m_workMatrix.coeff(i,i));
759 m_singularValues.coeffRef(i) =
abs(a);
760 if(computeU()) m_matrixU.col(i) *= m_workMatrix.coeff(i,i)/a;
765 RealScalar a = numext::real(m_workMatrix.coeff(i,i));
766 m_singularValues.coeffRef(i) =
abs(a);
767 if(computeU() && (a<RealScalar(0))) m_matrixU.col(i) = -m_matrixU.col(i);
771 m_singularValues *= scale;
775 m_nonzeroSingularValues = m_diagSize;
776 for(
Index i = 0; i < m_diagSize; i++)
779 RealScalar maxRemainingSingularValue = m_singularValues.tail(m_diagSize-i).maxCoeff(&pos);
780 if(maxRemainingSingularValue == RealScalar(0))
782 m_nonzeroSingularValues = i;
788 std::swap(m_singularValues.coeffRef(i), m_singularValues.coeffRef(pos));
789 if(computeU()) m_matrixU.col(pos).swap(m_matrixU.col(i));
790 if(computeV()) m_matrixV.col(pos).swap(m_matrixV.col(i));
794 m_isInitialized =
true;
805template<
typename Derived>
internal::traits< Homogeneous< MatrixType, Direction_ > >::Scalar Scalar
Definition: DenseBase.h:61
Rotation given by a cosine-sine pair.
Definition: Jacobi.h:37
JacobiRotation adjoint() const
Definition: Jacobi.h:69
JacobiRotation transpose() const
Definition: Jacobi.h:65
Two-sided Jacobi SVD decomposition of a rectangular matrix.
Definition: JacobiSVD.h:492
EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
Definition: EigenBase.h:65
JacobiSVD & compute(const MatrixType &matrix, unsigned int computationOptions)
Method performing the decomposition of given matrix using custom options.
Definition: JacobiSVD.h:668
JacobiSVD(Index rows, Index cols, unsigned int computationOptions=0)
Default Constructor with memory preallocation.
Definition: JacobiSVD.h:534
EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Definition: EigenBase.h:62
JacobiSVD & compute(const MatrixType &matrix)
Method performing the decomposition of given matrix using current options.
Definition: JacobiSVD.h:572
JacobiSVD()
Default Constructor.
Definition: JacobiSVD.h:524
JacobiSVD(const MatrixType &matrix, unsigned int computationOptions=0)
Constructor performing the decomposition of given matrix.
Definition: JacobiSVD.h:549
JacobiSVD< PlainObject > jacobiSvd(unsigned int computationOptions=0) const
Definition: JacobiSVD.h:807
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:182
Base class of SVD algorithms.
Definition: SVDBase.h:66
Index rank() const
Definition: SVDBase.h:150
bool computeV() const
Definition: SVDBase.h:212
bool computeU() const
Definition: SVDBase.h:210
@ NoQRPreconditioner
Definition: Constants.h:427
@ HouseholderQRPreconditioner
Definition: Constants.h:429
@ ColPivHouseholderQRPreconditioner
Definition: Constants.h:431
@ FullPivHouseholderQRPreconditioner
Definition: Constants.h:433
@ InvalidInput
Definition: Constants.h:451
@ Success
Definition: Constants.h:444
@ ComputeFullV
Definition: Constants.h:399
@ ComputeThinV
Definition: Constants.h:401
@ ComputeFullU
Definition: Constants.h:395
@ ComputeThinU
Definition: Constants.h:397
Namespace containing all symbols from the Eigen library.
Definition: B01_Experimental.dox:1
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_sqrt_op< typename Derived::Scalar >, const Derived > sqrt(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_conjugate_op< typename Derived::Scalar >, const Derived > conj(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_abs_op< typename Derived::Scalar >, const Derived > abs(const Eigen::ArrayBase< Derived > &x)
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:59
const int Dynamic
Definition: Constants.h:24
Definition: EigenBase.h:32
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:41
Holds information about the various numeric (i.e. scalar) types allowed by Eigen.
Definition: NumTraits.h:235