Tapkee
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
matrix_operations.hpp
Go to the documentation of this file.
1 /* This software is distributed under BSD 3-clause license (see LICENSE file).
2  *
3  * Copyright (c) 2012-2013 Sergey Lisitsyn
4  *
5  */
6 
7 #ifndef TAPKEE_MATRIX_OPS_H_
8 #define TAPKEE_MATRIX_OPS_H_
9 
10 /* Tapkee includes */
11 #include <tapkee/defines.hpp>
12 /* End of Tapkee includes */
13 
14 namespace tapkee
15 {
16 namespace tapkee_internal
17 {
18 
26 {
28  {
29  solver.compute(matrix);
30  }
33  inline DenseMatrix operator()(const DenseMatrix& operatee)
34  {
35  return solver.solve(operatee);
36  }
38  static const char* ARPACK_CODE;
39  static const bool largest;
40 };
42 const bool SparseInverseMatrixOperation::largest = false;
43 
51 {
53  {
54  solver.compute(matrix);
55  }
58  inline DenseMatrix operator()(const DenseMatrix& operatee)
59  {
60  return solver.solve(operatee);
61  }
63  static const char* ARPACK_CODE;
64  static const bool largest;
65 };
67 const bool DenseInverseMatrixOperation::largest = false;
68 
76 {
77  DenseMatrixOperation(const DenseMatrix& matrix) : _matrix(matrix)
78  {
79  }
85  inline DenseMatrix operator()(const DenseMatrix& rhs)
86  {
87  return _matrix.selfadjointView<Eigen::Upper>()*rhs;
88  }
90  static const char* ARPACK_CODE;
91  static const bool largest;
92 };
93 const char* DenseMatrixOperation::ARPACK_CODE = "LM";
94 const bool DenseMatrixOperation::largest = true;
95 
104 {
106  {
107  }
113  inline DenseMatrix operator()(const DenseMatrix& rhs)
114  {
115  return _matrix.selfadjointView<Eigen::Upper>()*(_matrix.selfadjointView<Eigen::Upper>()*rhs);
116  }
118  static const char* ARPACK_CODE;
119  static const bool largest;
120 };
123 
132 {
134  {
135  }
141  inline DenseMatrix operator()(const DenseMatrix& rhs)
142  {
143  return _matrix*(_matrix.transpose()*rhs);
144  }
146  static const char* ARPACK_CODE;
147  static const bool largest;
148 };
151 
152 #ifdef TAPKEE_GPU
153 struct GPUDenseImplicitSquareMatrixOperation
154 {
155  GPUDenseImplicitSquareMatrixOperation(const DenseMatrix& matrix)
156  {
157  timed_context c("Storing matrices");
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);
162  }
168  inline DenseMatrix operator()(const DenseMatrix& rhs)
169  {
170  timed_context c("Computing product");
171  viennacl::copy(rhs,vec);
172  res = viennacl::linalg::prod(mat, vec);
173  vec = res;
174  res = viennacl::linalg::prod(mat, vec);
175  DenseMatrix result(rhs);
176  viennacl::copy(res,result);
177  return result;
178  }
179  viennacl::matrix<ScalarType> mat;
180  viennacl::matrix<ScalarType> vec;
181  viennacl::matrix<ScalarType> res;
182  static const char* ARPACK_CODE;
183  static bool largest;
184 };
185 const char* GPUDenseImplicitSquareMatrixOperation::ARPACK_CODE = "LM";
186 const bool GPUDenseImplicitSquareMatrixOperation::largest = true;
187 
188 struct GPUDenseMatrixOperation
189 {
190  GPUDenseMatrixOperation(const DenseMatrix& matrix)
191  {
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);
196  }
202  inline DenseMatrix operator()(const DenseMatrix& rhs)
203  {
204  viennacl::copy(rhs,vec);
205  res = viennacl::linalg::prod(mat, vec);
206  DenseMatrix result(rhs);
207  viennacl::copy(res,result);
208  return result;
209  }
210  viennacl::matrix<ScalarType> mat;
211  viennacl::matrix<ScalarType> vec;
212  viennacl::matrix<ScalarType> res;
213  static const char* ARPACK_CODE;
214  static bool largest;
215 };
216 const char* GPUDenseMatrixOperation::ARPACK_CODE = "LM";
217 const bool GPUDenseMatrixOperation::largest = true;
218 #endif
219 
220 }
221 }
222 
223 #endif
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)
Definition: types.hpp:23
DenseMatrix operator()(const DenseMatrix &operatee)
Matrix-matrix operation used to compute largest eigenvalues and associated eigenvectors of X*X^T like...
Eigen::SparseMatrix< tapkee::ScalarType > SparseWeightMatrix
sparse weight matrix type (non-overridable)
Definition: types.hpp:29
DenseMatrix operator()(const DenseMatrix &operatee)
Matrix-matrix operation used to compute smallest eigenvalues and associated eigenvectors of a sparse ...
Eigen::SimplicialLDLT< tapkee::SparseWeightMatrix > SparseSolver
Definition: types.hpp:45
DenseMatrix operator()(const DenseMatrix &rhs)
Computes matrix product of the matrix and provided right-hand side matrix twice.
DenseMatrix operator()(const DenseMatrix &rhs)
Computes matrix product of the matrix and provided right-hand side matrix.
Matrix-matrix operation used to compute largest eigenvalues and associated eigenvectors of X*X^T like...
Matrix-matrix operation used to compute largest eigenvalues and associated eigenvectors. Essentially computes matrix product with provided right-hand side part.
Eigen::LDLT< tapkee::DenseMatrix > DenseSolver
dense solver (non-overridable)
Definition: types.hpp:35
DenseMatrix operator()(const DenseMatrix &rhs)
Computes matrix product of the matrix and provided right-hand side matrix twice.