Eigen::internal::lapacke_helpers::BDCSVD_LAPACKE< MatrixType_, Options > Class Template Reference

#include <BDCSVD_LAPACKE.h>

+ Inheritance diagram for Eigen::internal::lapacke_helpers::BDCSVD_LAPACKE< MatrixType_, Options >:

Public Member Functions

 BDCSVD_LAPACKE (SVD &&svd)
 
void compute_impl_lapacke (const MatrixType &matrix, unsigned int computationOptions)
 
- Public Member Functions inherited from Eigen::BDCSVD< MatrixType_, Options >
 BDCSVD ()
 Default Constructor. More...
 
 BDCSVD (Index rows, Index cols)
 Default Constructor with memory preallocation. More...
 
EIGEN_DEPRECATED BDCSVD (Index rows, Index cols, unsigned int computationOptions)
 Default Constructor with memory preallocation. More...
 
 BDCSVD (const MatrixType &matrix)
 Constructor performing the decomposition of given matrix, using the custom options specified with the Options template parameter. More...
 
EIGEN_DEPRECATED BDCSVD (const MatrixType &matrix, unsigned int computationOptions)
 Constructor performing the decomposition of given matrix using specified options for computing unitaries. More...
 
 ~BDCSVD ()
 
BDCSVDcompute (const MatrixType &matrix)
 Method performing the decomposition of given matrix. Computes Thin/Full unitaries U/V if specified using the Options template parameter or the class constructor. More...
 
EIGEN_DEPRECATED BDCSVDcompute (const MatrixType &matrix, unsigned int computationOptions)
 Method performing the decomposition of given matrix, as specified by the computationOptions parameter. More...
 
void setSwitchSize (int s)
 
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index cols () const EIGEN_NOEXCEPT
 
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index rows () const EIGEN_NOEXCEPT
 
- Public Member Functions inherited from Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >
BDCSVD< MatrixType_, Options_ > & derived ()
 
const BDCSVD< MatrixType_, Options_ > & derived () const
 
const MatrixUTypematrixU () const
 
const MatrixVTypematrixV () const
 
const SingularValuesTypesingularValues () const
 
Index nonzeroSingularValues () const
 
Index rank () const
 
BDCSVD< MatrixType_, Options_ > & setThreshold (const RealScalar &threshold)
 
BDCSVD< MatrixType_, Options_ > & setThreshold (Default_t)
 
RealScalar threshold () const
 
bool computeU () const
 
bool computeV () const
 
Index rows () const
 
Index cols () const
 
Index diagSize () const
 
EIGEN_DEVICE_FUNC ComputationInfo info () const
 Reports whether previous computation was successful. More...
 
void _solve_impl (const RhsType &rhs, DstType &dst) const
 
void _solve_impl_transposed (const RhsType &rhs, DstType &dst) const
 
- Public Member Functions inherited from Eigen::SolverBase< Derived >
 SolverBase ()
 
 ~SolverBase ()
 
template<typename Rhs >
const Solve< Derived, Rhs > solve (const MatrixBase< Rhs > &b) const
 
const ConstTransposeReturnType transpose () const
 
const AdjointReturnType adjoint () const
 
constexpr EIGEN_DEVICE_FUNC Derived & derived ()
 
constexpr EIGEN_DEVICE_FUNC const Derived & derived () const
 
- Public Member Functions inherited from Eigen::EigenBase< Derived >
constexpr EIGEN_DEVICE_FUNC Derived & derived ()
 
constexpr EIGEN_DEVICE_FUNC const Derived & derived () const
 
EIGEN_DEVICE_FUNC Derived & const_cast_derived () const
 
EIGEN_DEVICE_FUNC const Derived & const_derived () const
 
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index rows () const EIGEN_NOEXCEPT
 
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index cols () const EIGEN_NOEXCEPT
 
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index size () const EIGEN_NOEXCEPT
 
template<typename Dest >
EIGEN_DEVICE_FUNC void evalTo (Dest &dst) const
 
template<typename Dest >
EIGEN_DEVICE_FUNC void addTo (Dest &dst) const
 
template<typename Dest >
EIGEN_DEVICE_FUNC void subTo (Dest &dst) const
 
template<typename Dest >
EIGEN_DEVICE_FUNC void applyThisOnTheRight (Dest &dst) const
 
template<typename Dest >
EIGEN_DEVICE_FUNC void applyThisOnTheLeft (Dest &dst) const
 
template<typename Device >
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE DeviceWrapper< Derived, Device > device (Device &device)
 
template<typename Device >
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE DeviceWrapper< const Derived, Device > device (Device &device) const
 

Private Types

typedef BDCSVD< MatrixType_, OptionsSVD
 
typedef SVD::MatrixType MatrixType
 
typedef SVD::Scalar Scalar
 
typedef SVD::RealScalar RealScalar
 

Additional Inherited Members

- Public Types inherited from Eigen::BDCSVD< MatrixType_, Options >
enum  
 
typedef MatrixType_ MatrixType
 
typedef Base::Scalar Scalar
 
typedef Base::RealScalar RealScalar
 
typedef NumTraits< RealScalar >::Literal Literal
 
typedef Base::Index Index
 
typedef Base::MatrixUType MatrixUType
 
typedef Base::MatrixVType MatrixVType
 
typedef Base::SingularValuesType SingularValuesType
 
typedef Matrix< Scalar, Dynamic, Dynamic, ColMajor > MatrixX
 
typedef Matrix< RealScalar, Dynamic, Dynamic, ColMajor > MatrixXr
 
typedef Matrix< RealScalar, Dynamic, 1 > VectorType
 
typedef Array< RealScalar, Dynamic, 1 > ArrayXr
 
typedef Array< Index, 1, Dynamic > ArrayXi
 
typedef Ref< ArrayXrArrayRef
 
typedef Ref< ArrayXiIndicesRef
 
- Public Types inherited from Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >
enum  
 
typedef internal::traits< BDCSVD< MatrixType_, Options_ > >::MatrixType MatrixType
 
typedef MatrixType::Scalar Scalar
 
typedef NumTraits< typename MatrixType::Scalar >::Real RealScalar
 
typedef Eigen::internal::traits< SVDBase >::StorageIndex StorageIndex
 
typedef internal::make_proper_matrix_type< Scalar, RowsAtCompileTime, MatrixUColsAtCompileTime, MatrixOptions, MaxRowsAtCompileTime, MatrixUMaxColsAtCompileTime >::type MatrixUType
 
typedef internal::make_proper_matrix_type< Scalar, ColsAtCompileTime, MatrixVColsAtCompileTime, MatrixOptions, MaxColsAtCompileTime, MatrixVMaxColsAtCompileTime >::type MatrixVType
 
typedef internal::plain_diag_type< MatrixType, RealScalar >::type SingularValuesType
 
- Public Types inherited from Eigen::SolverBase< Derived >
enum  {
  RowsAtCompileTime = internal::traits<Derived>::RowsAtCompileTime , ColsAtCompileTime = internal::traits<Derived>::ColsAtCompileTime , SizeAtCompileTime = (internal::size_of_xpr_at_compile_time<Derived>::ret) , MaxRowsAtCompileTime = internal::traits<Derived>::MaxRowsAtCompileTime ,
  MaxColsAtCompileTime = internal::traits<Derived>::MaxColsAtCompileTime , MaxSizeAtCompileTime , IsVectorAtCompileTime , NumDimensions
}
 
typedef EigenBase< Derived > Base
 
typedef internal::traits< Derived >::Scalar Scalar
 
typedef Scalar CoeffReturnType
 
typedef Transpose< const Derived > ConstTransposeReturnType
 
typedef std::conditional_t< NumTraits< Scalar >::IsComplex, CwiseUnaryOp< internal::scalar_conjugate_op< Scalar >, const ConstTransposeReturnType >, const ConstTransposeReturnTypeAdjointReturnType
 
- Public Types inherited from Eigen::EigenBase< Derived >
typedef Eigen::Index Index
 The interface type of indices. More...
 
typedef internal::traits< Derived >::StorageKind StorageKind
 
- Public Attributes inherited from Eigen::BDCSVD< MatrixType_, Options >
int m_numIters
 
- Static Public Attributes inherited from Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >
static constexpr bool ShouldComputeFullU
 
static constexpr bool ShouldComputeThinU
 
static constexpr bool ShouldComputeFullV
 
static constexpr bool ShouldComputeThinV
 
- Protected Member Functions inherited from Eigen::BDCSVD< MatrixType_, Options >
void allocate (Index rows, Index cols, unsigned int computationOptions)
 
- Protected Member Functions inherited from Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >
void _check_compute_assertions () const
 
void _check_solve_assertion (const Rhs &b) const
 
bool allocate (Index rows, Index cols, unsigned int computationOptions)
 
 SVDBase ()
 Default Constructor. More...
 
- Protected Member Functions inherited from Eigen::SolverBase< Derived >
template<bool Transpose_, typename Rhs >
void _check_solve_assertion (const Rhs &b) const
 
- Protected Attributes inherited from Eigen::BDCSVD< MatrixType_, Options >
MatrixXr m_naiveU
 
MatrixXr m_naiveV
 
MatrixXr m_computed
 
Index m_nRec
 
ArrayXr m_workspace
 
ArrayXi m_workspaceI
 
int m_algoswap
 
bool m_isTranspose
 
bool m_compU
 
bool m_compV
 
bool m_useQrDecomp
 
JacobiSVD< MatrixType, ComputationOptionssmallSvd
 
HouseholderQR< MatrixXqrDecomp
 
internal::UpperBidiagonalization< MatrixXbid
 
MatrixX copyWorkspace
 
MatrixX reducedTriangle
 
- Protected Attributes inherited from Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >
MatrixUType m_matrixU
 
MatrixVType m_matrixV
 
SingularValuesType m_singularValues
 
ComputationInfo m_info
 
bool m_isInitialized
 
bool m_isAllocated
 
bool m_usePrescribedThreshold
 
bool m_computeFullU
 
bool m_computeThinU
 
bool m_computeFullV
 
bool m_computeThinV
 
unsigned int m_computationOptions
 
Index m_nonzeroSingularValues
 
internal::variable_if_dynamic< Index, RowsAtCompileTimem_rows
 
internal::variable_if_dynamic< Index, ColsAtCompileTimem_cols
 
internal::variable_if_dynamic< Index, DiagSizeAtCompileTimem_diagSize
 
RealScalar m_prescribedThreshold
 

Detailed Description

template<typename MatrixType_, int Options>
class Eigen::internal::lapacke_helpers::BDCSVD_LAPACKE< MatrixType_, Options >

Specialization for the data types supported by LAPACKe

Member Typedef Documentation

◆ MatrixType

template<typename MatrixType_ , int Options>
typedef SVD::MatrixType Eigen::internal::lapacke_helpers::BDCSVD_LAPACKE< MatrixType_, Options >::MatrixType
private

◆ RealScalar

template<typename MatrixType_ , int Options>
typedef SVD::RealScalar Eigen::internal::lapacke_helpers::BDCSVD_LAPACKE< MatrixType_, Options >::RealScalar
private

◆ Scalar

template<typename MatrixType_ , int Options>
typedef SVD::Scalar Eigen::internal::lapacke_helpers::BDCSVD_LAPACKE< MatrixType_, Options >::Scalar
private

◆ SVD

template<typename MatrixType_ , int Options>
typedef BDCSVD<MatrixType_, Options> Eigen::internal::lapacke_helpers::BDCSVD_LAPACKE< MatrixType_, Options >::SVD
private

Constructor & Destructor Documentation

◆ BDCSVD_LAPACKE()

template<typename MatrixType_ , int Options>
Eigen::internal::lapacke_helpers::BDCSVD_LAPACKE< MatrixType_, Options >::BDCSVD_LAPACKE ( SVD &&  svd)
inline
59 : SVD(std::move(svd)) {}
cout<< "Here is the matrix m:"<< endl<< m<< endl;JacobiSVD< MatrixXf, ComputeThinU|ComputeThinV > svd(m)
BDCSVD< MatrixType_, Options > SVD
Definition: BDCSVD_LAPACKE.h:52

Member Function Documentation

◆ compute_impl_lapacke()

template<typename MatrixType_ , int Options>
void Eigen::internal::lapacke_helpers::BDCSVD_LAPACKE< MatrixType_, Options >::compute_impl_lapacke ( const MatrixType matrix,
unsigned int  computationOptions 
)
inline
61  {
62  SVD::allocate(matrix.rows(), matrix.cols(), computationOptions);
63 
65 
66  // prepare arguments to ?gesdd
67  const lapack_int matrix_order = lapack_storage_of(matrix);
68  const char jobz = (SVD::m_computeFullU || SVD::m_computeFullV) ? 'A'
70  : 'N';
71  const lapack_int u_cols = (jobz == 'A') ? to_lapack(SVD::rows()) : (jobz == 'S') ? to_lapack(SVD::diagSize()) : 1;
72  const lapack_int vt_rows = (jobz == 'A') ? to_lapack(SVD::cols()) : (jobz == 'S') ? to_lapack(SVD::diagSize()) : 1;
73  lapack_int ldu, ldvt;
74  Scalar *u, *vt, dummy;
75  MatrixType localU;
77  ldu = to_lapack(SVD::m_matrixU.outerStride());
78  u = SVD::m_matrixU.data();
79  } else if (SVD::computeV()) {
80  localU.resize(SVD::rows(), u_cols);
81  ldu = to_lapack(localU.outerStride());
82  u = localU.data();
83  } else {
84  ldu = 1;
85  u = &dummy;
86  }
87  MatrixType localV;
88  if (SVD::computeU() || SVD::computeV()) {
89  localV.resize(vt_rows, SVD::cols());
90  ldvt = to_lapack(localV.outerStride());
91  vt = localV.data();
92  } else {
93  ldvt = 1;
94  vt = &dummy;
95  }
96  MatrixType temp;
97  temp = matrix;
98 
99  // actual call to ?gesdd
100  lapack_int info = gesdd(matrix_order, jobz, to_lapack(SVD::rows()), to_lapack(SVD::cols()), to_lapack(temp.data()),
101  to_lapack(temp.outerStride()), (RealScalar*)SVD::m_singularValues.data(), to_lapack(u), ldu,
102  to_lapack(vt), ldvt);
103 
104  // Check the result of the LAPACK call
105  if (info < 0 || !SVD::m_singularValues.allFinite()) {
106  // this includes info == -4 => NaN entry in A
108  } else if (info > 0) {
110  } else {
113  SVD::m_matrixU = localU.leftCols(SVD::m_matrixU.cols());
114  }
115  if (SVD::computeV()) {
116  SVD::m_matrixV = localV.adjoint().leftCols(SVD::m_matrixV.cols());
117  }
118  }
119  SVD::m_isInitialized = true;
120  }
SCALAR Scalar
Definition: bench_gemm.cpp:45
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:46
MatrixXf MatrixType
Definition: benchmark-blocking-sizes.cpp:52
void allocate(Index rows, Index cols, unsigned int computationOptions)
Definition: BDCSVD.h:268
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
Definition: EigenBase.h:61
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Definition: EigenBase.h:59
constexpr EIGEN_DEVICE_FUNC const Scalar * data() const
Definition: PlainObjectBase.h:273
bool m_computeFullV
Definition: SVDBase.h:337
internal::variable_if_dynamic< Index, DiagSizeAtCompileTime > m_diagSize
Definition: SVDBase.h:342
ComputationInfo m_info
Definition: SVDBase.h:334
bool m_computeThinU
Definition: SVDBase.h:336
EIGEN_DEVICE_FUNC ComputationInfo info() const
Reports whether previous computation was successful.
Definition: SVDBase.h:300
bool computeV() const
Definition: SVDBase.h:275
bool m_isInitialized
Definition: SVDBase.h:335
MatrixVType m_matrixV
Definition: SVDBase.h:332
bool computeU() const
Definition: SVDBase.h:273
Index m_nonzeroSingularValues
Definition: SVDBase.h:339
SingularValuesType m_singularValues
Definition: SVDBase.h:333
bool m_computeThinV
Definition: SVDBase.h:337
bool m_computeFullU
Definition: SVDBase.h:336
Index diagSize() const
Definition: SVDBase.h:279
MatrixUType m_matrixU
Definition: SVDBase.h:331
Eigen::Map< Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic, Eigen::ColMajor >, 0, Eigen::OuterStride<> > matrix(T *data, int rows, int cols, int stride)
Definition: common.h:85
@ InvalidInput
Definition: Constants.h:447
@ Success
Definition: Constants.h:440
@ NoConvergence
Definition: Constants.h:444
#define lapack_int
Definition: lapacke.h:52
EIGEN_ALWAYS_INLINE auto to_lapack(Source value)
Definition: lapacke_helpers.h:61
EIGEN_ALWAYS_INLINE EIGEN_CONSTEXPR lapack_int lapack_storage_of(const EigenBase< Derived > &)
translates storage order of the given Eigen object to the corresponding lapack constant
Definition: lapacke_helpers.h:78

References Eigen::BDCSVD< MatrixType_, Options >::allocate(), Eigen::BDCSVD< MatrixType_, Options >::cols(), Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::computeU(), Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::computeV(), Eigen::PlainObjectBase< Derived >::data(), Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::diagSize(), Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::info(), Eigen::InvalidInput, lapack_int, Eigen::internal::lapacke_helpers::lapack_storage_of(), Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_computeFullU, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_computeFullV, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_computeThinU, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_computeThinV, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_diagSize, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_info, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_isInitialized, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_matrixU, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_matrixV, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_nonzeroSingularValues, Eigen::SVDBase< BDCSVD< MatrixType_, Options_ > >::m_singularValues, matrix(), Eigen::NoConvergence, Eigen::BDCSVD< MatrixType_, Options >::rows(), Eigen::Success, and Eigen::internal::lapacke_helpers::to_lapack().

Referenced by Eigen::internal::lapacke_helpers::BDCSVD_wrapper().


The documentation for this class was generated from the following file: