Tapkee
matrix_operations.hpp
Go to the documentation of this file.
00001 /* This software is distributed under BSD 3-clause license (see LICENSE file).
00002  *
00003  * Copyright (c) 2012-2013 Sergey Lisitsyn
00004  *
00005  */
00006 
00007 #ifndef TAPKEE_MATRIX_OPS_H_
00008 #define TAPKEE_MATRIX_OPS_H_
00009 
00010 /* Tapkee includes */
00011 #include <tapkee/defines.hpp>
00012 /* End of Tapkee includes */
00013 
00014 namespace tapkee
00015 {
00016 namespace tapkee_internal
00017 {
00018 
00025 struct SparseInverseMatrixOperation
00026 {
00027     SparseInverseMatrixOperation(const SparseWeightMatrix& matrix) : solver()
00028     {
00029         solver.compute(matrix);
00030     }
00033     inline DenseMatrix operator()(const DenseMatrix& operatee)
00034     {
00035         return solver.solve(operatee);
00036     }
00037     SparseSolver solver;
00038     static const char* ARPACK_CODE;
00039     static const bool largest;
00040 };
00041 const char* SparseInverseMatrixOperation::ARPACK_CODE = "SM";
00042 const bool SparseInverseMatrixOperation::largest = false;
00043 
00050 struct DenseInverseMatrixOperation
00051 {
00052     DenseInverseMatrixOperation(const DenseMatrix& matrix) : solver()
00053     {
00054         solver.compute(matrix);
00055     }
00058     inline DenseMatrix operator()(const DenseMatrix& operatee)
00059     {
00060         return solver.solve(operatee);
00061     }
00062     DenseSolver solver;
00063     static const char* ARPACK_CODE;
00064     static const bool largest;
00065 };
00066 const char* DenseInverseMatrixOperation::ARPACK_CODE = "SM";
00067 const bool DenseInverseMatrixOperation::largest = false;
00068 
00075 struct DenseMatrixOperation
00076 {
00077     DenseMatrixOperation(const DenseMatrix& matrix) : _matrix(matrix)
00078     {
00079     }
00085     inline DenseMatrix operator()(const DenseMatrix& rhs)
00086     {
00087         return _matrix.selfadjointView<Eigen::Upper>()*rhs;
00088     }
00089     const DenseMatrix& _matrix;
00090     static const char* ARPACK_CODE;
00091     static const bool largest;
00092 };
00093 const char* DenseMatrixOperation::ARPACK_CODE = "LM";
00094 const bool DenseMatrixOperation::largest = true;
00095 
00103 struct DenseImplicitSquareSymmetricMatrixOperation
00104 {
00105     DenseImplicitSquareSymmetricMatrixOperation(const DenseMatrix& matrix) : _matrix(matrix)
00106     {
00107     }
00113     inline DenseMatrix operator()(const DenseMatrix& rhs)
00114     {
00115         return _matrix.selfadjointView<Eigen::Upper>()*(_matrix.selfadjointView<Eigen::Upper>()*rhs);
00116     }
00117     const DenseMatrix& _matrix;
00118     static const char* ARPACK_CODE;
00119     static const bool largest;
00120 };
00121 const char* DenseImplicitSquareSymmetricMatrixOperation::ARPACK_CODE = "LM";
00122 const bool DenseImplicitSquareSymmetricMatrixOperation::largest = true;
00123 
00131 struct DenseImplicitSquareMatrixOperation
00132 {
00133     DenseImplicitSquareMatrixOperation(const DenseMatrix& matrix) : _matrix(matrix)
00134     {
00135     }
00141     inline DenseMatrix operator()(const DenseMatrix& rhs)
00142     {
00143         return _matrix*(_matrix.transpose()*rhs);
00144     }
00145     const DenseMatrix& _matrix;
00146     static const char* ARPACK_CODE;
00147     static const bool largest;
00148 };
00149 const char* DenseImplicitSquareMatrixOperation::ARPACK_CODE = "LM";
00150 const bool DenseImplicitSquareMatrixOperation::largest = true;
00151 
00152 #ifdef TAPKEE_GPU
00153 struct GPUDenseImplicitSquareMatrixOperation
00154 {
00155     GPUDenseImplicitSquareMatrixOperation(const DenseMatrix& matrix)
00156     {
00157         timed_context c("Storing matrices");
00158         mat = viennacl::matrix<ScalarType>(matrix.cols(),matrix.rows());
00159         vec = viennacl::matrix<ScalarType>(matrix.cols(),1);
00160         res = viennacl::matrix<ScalarType>(matrix.cols(),1);
00161         viennacl::copy(matrix,mat);
00162     }
00168     inline DenseMatrix operator()(const DenseMatrix& rhs)
00169     {
00170         timed_context c("Computing product");
00171         viennacl::copy(rhs,vec);
00172         res = viennacl::linalg::prod(mat, vec);
00173         vec = res;
00174         res = viennacl::linalg::prod(mat, vec);
00175         DenseMatrix result(rhs);
00176         viennacl::copy(res,result);
00177         return result;
00178     }
00179     viennacl::matrix<ScalarType> mat;
00180     viennacl::matrix<ScalarType> vec;
00181     viennacl::matrix<ScalarType> res;
00182     static const char* ARPACK_CODE;
00183     static bool largest;
00184 };
00185 const char* GPUDenseImplicitSquareMatrixOperation::ARPACK_CODE = "LM";
00186 const bool GPUDenseImplicitSquareMatrixOperation::largest = true;
00187 
00188 struct GPUDenseMatrixOperation
00189 {
00190     GPUDenseMatrixOperation(const DenseMatrix& matrix)
00191     {
00192         mat = viennacl::matrix<ScalarType>(matrix.cols(),matrix.rows());
00193         vec = viennacl::matrix<ScalarType>(matrix.cols(),1);
00194         res = viennacl::matrix<ScalarType>(matrix.cols(),1);
00195         viennacl::copy(matrix,mat);
00196     }
00202     inline DenseMatrix operator()(const DenseMatrix& rhs)
00203     {
00204         viennacl::copy(rhs,vec);
00205         res = viennacl::linalg::prod(mat, vec);
00206         DenseMatrix result(rhs);
00207         viennacl::copy(res,result);
00208         return result;
00209     }
00210     viennacl::matrix<ScalarType> mat;
00211     viennacl::matrix<ScalarType> vec;
00212     viennacl::matrix<ScalarType> res;
00213     static const char* ARPACK_CODE;
00214     static bool largest;
00215 };
00216 const char* GPUDenseMatrixOperation::ARPACK_CODE = "LM";
00217 const bool GPUDenseMatrixOperation::largest = true;
00218 #endif
00219 
00220 }
00221 }
00222 
00223 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines