Tapkee
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 #ifdef TAPKEE_WITH_VIENNACL
15  #define VIENNACL_HAVE_EIGEN
16  #include <viennacl/matrix.hpp>
17  #include <viennacl/vector.hpp>
18  #include <viennacl/linalg/prod.hpp>
19 #endif
20 
21 namespace tapkee
22 {
23 namespace tapkee_internal
24 {
25 
33 {
35  {
36  solver.compute(matrix);
37  }
40  inline DenseMatrix operator()(const DenseMatrix& operatee)
41  {
42  return solver.solve(operatee);
43  }
45  static const char* ARPACK_CODE;
46  static const bool largest;
47 };
49 const bool SparseInverseMatrixOperation::largest = false;
50 
58 {
60  {
61  solver.compute(matrix);
62  }
65  inline DenseMatrix operator()(const DenseMatrix& operatee)
66  {
67  return solver.solve(operatee);
68  }
70  static const char* ARPACK_CODE;
71  static const bool largest;
72 };
74 const bool DenseInverseMatrixOperation::largest = false;
75 
83 {
84  DenseMatrixOperation(const DenseMatrix& matrix) : _matrix(matrix)
85  {
86  }
92  inline DenseMatrix operator()(const DenseMatrix& rhs)
93  {
94  return _matrix.selfadjointView<Eigen::Upper>()*rhs;
95  }
97  static const char* ARPACK_CODE;
98  static const bool largest;
99 };
100 const char* DenseMatrixOperation::ARPACK_CODE = "LA";
101 const bool DenseMatrixOperation::largest = true;
102 
111 {
113  {
114  }
120  inline DenseMatrix operator()(const DenseMatrix& rhs)
121  {
122  return _matrix.selfadjointView<Eigen::Upper>()*(_matrix.selfadjointView<Eigen::Upper>()*rhs);
123  }
125  static const char* ARPACK_CODE;
126  static const bool largest;
127 };
130 
139 {
141  {
142  }
148  inline DenseMatrix operator()(const DenseMatrix& rhs)
149  {
150  return _matrix*(_matrix.transpose()*rhs);
151  }
153  static const char* ARPACK_CODE;
154  static const bool largest;
155 };
158 
159 #ifdef TAPKEE_WITH_VIENNACL
160 struct GPUDenseImplicitSquareMatrixOperation
161 {
162  GPUDenseImplicitSquareMatrixOperation(const DenseMatrix& matrix)
163  {
164  mat = viennacl::matrix<ScalarType>(matrix.cols(),matrix.rows());
165  vec = viennacl::vector<ScalarType>(matrix.cols());
166  res = viennacl::vector<ScalarType>(matrix.cols());
167  viennacl::copy(matrix,mat);
168  }
174  inline DenseVector operator()(const DenseVector& rhs)
175  {
176  viennacl::copy(rhs,vec);
177  res = viennacl::linalg::prod(mat, vec);
178  vec = res;
179  res = viennacl::linalg::prod(mat, vec);
180  DenseVector result(rhs);
181  viennacl::copy(res,result);
182  return result;
183  }
184  viennacl::matrix<ScalarType> mat;
185  viennacl::vector<ScalarType> vec;
186  viennacl::vector<ScalarType> res;
187  static const char* ARPACK_CODE;
188  static const bool largest;
189 };
190 const char* GPUDenseImplicitSquareMatrixOperation::ARPACK_CODE = "LA";
191 const bool GPUDenseImplicitSquareMatrixOperation::largest = true;
192 
193 struct GPUDenseMatrixOperation
194 {
195  GPUDenseMatrixOperation(const DenseMatrix& matrix)
196  {
197  mat = viennacl::matrix<ScalarType>(matrix.cols(),matrix.rows());
198  vec = viennacl::vector<ScalarType>(matrix.cols());
199  res = viennacl::vector<ScalarType>(matrix.cols());
200  viennacl::copy(matrix,mat);
201  }
207  inline DenseVector operator()(const DenseVector& rhs)
208  {
209  viennacl::copy(rhs,vec);
210  res = viennacl::linalg::prod(mat, vec);
211  DenseVector result(rhs);
212  viennacl::copy(res,result);
213  return result;
214  }
215  viennacl::matrix<ScalarType> mat;
216  viennacl::vector<ScalarType> vec;
217  viennacl::vector<ScalarType> res;
218  static const char* ARPACK_CODE;
219  static const bool largest;
220 };
221 const char* GPUDenseMatrixOperation::ARPACK_CODE = "LA";
222 const bool GPUDenseMatrixOperation::largest = true;
223 #endif
224 
225 }
226 }
227 
228 #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)
Eigen::Matrix< tapkee::ScalarType, Eigen::Dynamic, 1 > DenseVector
dense vector type (non-overridable)
Definition: types.hpp:21
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.