SHOGUN
v3.2.0
|
00001 /* 00002 * This program is free software; you can redistribute it and/or modify 00003 * it under the terms of the GNU General Public License as published by 00004 * the Free Software Foundation; either version 3 of the License, or 00005 * (at your option) any later version. 00006 * 00007 * Written (W) 2013 Soumyajit De 00008 */ 00009 00010 #include <shogun/lib/config.h> 00011 00012 #ifdef HAVE_EIGEN3 00013 #include <shogun/lib/SGVector.h> 00014 #include <shogun/lib/SGMatrix.h> 00015 #include <shogun/base/Parameter.h> 00016 #include <shogun/mathematics/eigen3.h> 00017 #include <shogun/mathematics/linalg/linop/DenseMatrixOperator.h> 00018 00019 using namespace Eigen; 00020 00021 namespace shogun 00022 { 00023 00024 template<class T> 00025 CDenseMatrixOperator<T>::CDenseMatrixOperator() 00026 : CMatrixOperator<T>() 00027 { 00028 init(); 00029 00030 SG_SGCDEBUG("%s created (%p)\n", this->get_name(), this); 00031 } 00032 00033 template<class T> 00034 CDenseMatrixOperator<T>::CDenseMatrixOperator(SGMatrix<T> op) 00035 : CMatrixOperator<T>(op.num_cols), 00036 m_operator(op) 00037 { 00038 init(); 00039 00040 SG_SGCDEBUG("%s created (%p)\n", this->get_name(), this); 00041 } 00042 00043 template<class T> 00044 CDenseMatrixOperator<T>::CDenseMatrixOperator( 00045 const CDenseMatrixOperator<T>& orig) 00046 : CMatrixOperator<T>(orig.get_dimension()) 00047 { 00048 init(); 00049 00050 m_operator=SGMatrix<T>(orig.m_operator.num_rows, orig.m_operator.num_cols); 00051 for (index_t i=0; i<m_operator.num_cols; ++i) 00052 { 00053 for (index_t j=0; j<m_operator.num_rows; ++j) 00054 m_operator(j,i)=orig.m_operator(j,i); 00055 } 00056 00057 SG_SGCDEBUG("%s deep copy created (%p)\n", this->get_name(), this); 00058 } 00059 00060 template<class T> 00061 void CDenseMatrixOperator<T>::init() 00062 { 00063 CSGObject::set_generic<T>(); 00064 00065 this->m_parameters->add(&m_operator, "dense_matrix", 00066 "The dense matrix of the linear operator"); 00067 } 00068 00069 template<class T> 00070 CDenseMatrixOperator<T>::~CDenseMatrixOperator() 00071 { 00072 SG_SGCDEBUG("%s destroyed (%p)\n", this->get_name(), this); 00073 } 00074 00075 template<class T> 00076 SGMatrix<T> CDenseMatrixOperator<T>::get_matrix_operator() const 00077 { 00078 return m_operator; 00079 } 00080 00081 template<class T> 00082 SGVector<T> CDenseMatrixOperator<T>::get_diagonal() const 00083 { 00084 REQUIRE(m_operator.matrix, "Operator not initialized!\n"); 00085 00086 typedef Matrix<T, Dynamic, 1> VectorXt; 00087 typedef Matrix<T, Dynamic, Dynamic> MatrixXt; 00088 00089 Map<MatrixXt> _op(m_operator.matrix, m_operator.num_rows, 00090 m_operator.num_cols); 00091 00092 SGVector<T> diag(static_cast<int32_t>(_op.diagonalSize())); 00093 Map<VectorXt> _diag(diag.vector, diag.vlen); 00094 _diag=_op.diagonal(); 00095 00096 return diag; 00097 } 00098 00099 template<class T> 00100 void CDenseMatrixOperator<T>::set_diagonal(SGVector<T> diag) 00101 { 00102 REQUIRE(m_operator.matrix, "Operator not initialized!\n"); 00103 REQUIRE(diag.vector, "Diagonal not initialized!\n"); 00104 00105 typedef Matrix<T, Dynamic, 1> VectorXt; 00106 typedef Matrix<T, Dynamic, Dynamic> MatrixXt; 00107 00108 Map<MatrixXt> _op(m_operator.matrix, m_operator.num_rows, 00109 m_operator.num_cols); 00110 00111 REQUIRE(static_cast<int32_t>(_op.diagonalSize())==diag.vlen, 00112 "Dimension mismatch!\n"); 00113 00114 Map<VectorXt> _diag(diag.vector, diag.vlen); 00115 _op.diagonal()=_diag; 00116 } 00117 00118 template<class T> 00119 SGVector<T> CDenseMatrixOperator<T>::apply(SGVector<T> b) const 00120 { 00121 REQUIRE(m_operator.matrix, "Operator not initialized!\n"); 00122 REQUIRE(this->get_dimension()==b.vlen, 00123 "Number of rows of vector must be equal to the " 00124 "number of cols of the operator!\n"); 00125 00126 typedef Matrix<T, Dynamic, 1> VectorXt; 00127 typedef Matrix<T, Dynamic, Dynamic> MatrixXt; 00128 00129 Map<VectorXt> _b(b.vector, b.vlen); 00130 Map<MatrixXt> _op(m_operator.matrix, m_operator.num_rows, 00131 m_operator.num_cols); 00132 00133 SGVector<T> result(m_operator.num_rows); 00134 Map<VectorXt> _result(result.vector, result.vlen); 00135 _result=_op*_b; 00136 00137 return result; 00138 } 00139 00140 #define UNDEFINED(type) \ 00141 template<> \ 00142 SGVector<type> CDenseMatrixOperator<type>::apply(SGVector<type> b) const \ 00143 { \ 00144 SG_SERROR("Not supported for %s\n", #type);\ 00145 return b; \ 00146 } 00147 00148 UNDEFINED(bool) 00149 UNDEFINED(char) 00150 UNDEFINED(int8_t) 00151 UNDEFINED(uint8_t) 00152 #undef UNDEFINED 00153 00154 template class CDenseMatrixOperator<bool>; 00155 template class CDenseMatrixOperator<char>; 00156 template class CDenseMatrixOperator<int8_t>; 00157 template class CDenseMatrixOperator<uint8_t>; 00158 template class CDenseMatrixOperator<int16_t>; 00159 template class CDenseMatrixOperator<uint16_t>; 00160 template class CDenseMatrixOperator<int32_t>; 00161 template class CDenseMatrixOperator<uint32_t>; 00162 template class CDenseMatrixOperator<int64_t>; 00163 template class CDenseMatrixOperator<uint64_t>; 00164 template class CDenseMatrixOperator<float32_t>; 00165 template class CDenseMatrixOperator<float64_t>; 00166 template class CDenseMatrixOperator<floatmax_t>; 00167 template class CDenseMatrixOperator<complex128_t>; 00168 } 00169 #endif // HAVE_EIGEN3