7 #ifndef TAPKEE_MATRIX_OPS_H_
8 #define TAPKEE_MATRIX_OPS_H_
16 namespace tapkee_internal
35 return solver.solve(operatee);
60 return solver.solve(operatee);
87 return _matrix.selfadjointView<Eigen::Upper>()*rhs;
115 return _matrix.selfadjointView<Eigen::Upper>()*(
_matrix.selfadjointView<Eigen::Upper>()*rhs);
153 struct GPUDenseImplicitSquareMatrixOperation
155 GPUDenseImplicitSquareMatrixOperation(
const DenseMatrix& matrix)
158 mat = viennacl::matrix<ScalarType>(matrix.cols(),matrix.rows());
159 vec = viennacl::matrix<ScalarType>(matrix.cols(),1);
160 res = viennacl::matrix<ScalarType>(matrix.cols(),1);
161 viennacl::copy(matrix,mat);
170 timed_context c(
"Computing product");
171 viennacl::copy(rhs,vec);
172 res = viennacl::linalg::prod(mat, vec);
174 res = viennacl::linalg::prod(mat, vec);
176 viennacl::copy(res,result);
179 viennacl::matrix<ScalarType> mat;
180 viennacl::matrix<ScalarType> vec;
181 viennacl::matrix<ScalarType> res;
182 static const char* ARPACK_CODE;
185 const char* GPUDenseImplicitSquareMatrixOperation::ARPACK_CODE =
"LM";
186 const bool GPUDenseImplicitSquareMatrixOperation::largest =
true;
188 struct GPUDenseMatrixOperation
192 mat = viennacl::matrix<ScalarType>(matrix.cols(),matrix.rows());
193 vec = viennacl::matrix<ScalarType>(matrix.cols(),1);
194 res = viennacl::matrix<ScalarType>(matrix.cols(),1);
195 viennacl::copy(matrix,mat);
204 viennacl::copy(rhs,vec);
205 res = viennacl::linalg::prod(mat, vec);
207 viennacl::copy(res,result);
210 viennacl::matrix<ScalarType> mat;
211 viennacl::matrix<ScalarType> vec;
212 viennacl::matrix<ScalarType> res;
213 static const char* ARPACK_CODE;
216 const char* GPUDenseMatrixOperation::ARPACK_CODE =
"LM";
217 const bool GPUDenseMatrixOperation::largest =
true;
DenseMatrixOperation(const DenseMatrix &matrix)
DenseInverseMatrixOperation(const DenseMatrix &matrix)
Matrix-matrix operation used to compute smallest eigenvalues and associated eigenvectors of a dense m...
SparseInverseMatrixOperation(const SparseWeightMatrix &matrix)
Eigen::Matrix< tapkee::ScalarType, Eigen::Dynamic, Eigen::Dynamic > DenseMatrix
dense matrix type (non-overridable)
DenseMatrix operator()(const DenseMatrix &operatee)
static const char * ARPACK_CODE
static const char * ARPACK_CODE
Matrix-matrix operation used to compute largest eigenvalues and associated eigenvectors of X*X^T like...
const DenseMatrix & _matrix
const DenseMatrix & _matrix
static const bool largest
Eigen::SparseMatrix< tapkee::ScalarType > SparseWeightMatrix
sparse weight matrix type (non-overridable)
DenseMatrix operator()(const DenseMatrix &operatee)
static const char * ARPACK_CODE
Matrix-matrix operation used to compute smallest eigenvalues and associated eigenvectors of a sparse ...
const DenseMatrix & _matrix
Eigen::SimplicialLDLT< tapkee::SparseWeightMatrix > SparseSolver
DenseMatrix operator()(const DenseMatrix &rhs)
Computes matrix product of the matrix and provided right-hand side matrix twice.
static const bool largest
DenseImplicitSquareSymmetricMatrixOperation(const DenseMatrix &matrix)
DenseMatrix operator()(const DenseMatrix &rhs)
Computes matrix product of the matrix and provided right-hand side matrix.
static const bool largest
Matrix-matrix operation used to compute largest eigenvalues and associated eigenvectors of X*X^T like...
static const bool largest
Matrix-matrix operation used to compute largest eigenvalues and associated eigenvectors. Essentially computes matrix product with provided right-hand side part.
static const bool largest
Eigen::LDLT< tapkee::DenseMatrix > DenseSolver
dense solver (non-overridable)
DenseImplicitSquareMatrixOperation(const DenseMatrix &matrix)
static const char * ARPACK_CODE
DenseMatrix operator()(const DenseMatrix &rhs)
Computes matrix product of the matrix and provided right-hand side matrix twice.
static const char * ARPACK_CODE