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 #ifndef SPARSE_MATRIX_OPERATOR_H_ 00011 #define SPARSE_MATRIX_OPERATOR_H_ 00012 00013 #include <shogun/lib/config.h> 00014 #include <shogun/lib/SGSparseVector.h> 00015 #include <shogun/lib/SGSparseMatrix.h> 00016 #include <shogun/mathematics/linalg/linop/MatrixOperator.h> 00017 00018 namespace shogun 00019 { 00020 template<class T> class SGVector; 00021 template<class T> class SGSparseMatrix; 00022 00028 struct SparsityStructure 00029 { 00031 SparsityStructure() : m_num_rows(0), m_ptr(NULL) {} 00032 00040 SparsityStructure(index_t* row_offsets, index_t* column_indices, 00041 index_t num_rows) 00042 : m_num_rows(num_rows), 00043 m_ptr(new int32_t*[num_rows]()) 00044 { 00045 for (index_t i=0; i<m_num_rows; ++i) 00046 { 00047 index_t current_index=row_offsets[i]; 00048 index_t new_index=row_offsets[i+1]; 00049 index_t length_row=(new_index-current_index); 00050 00051 m_ptr[i]=new int32_t[length_row+1](); 00052 m_ptr[i][0]=length_row; 00053 00054 for (index_t j=1; j<=length_row; ++j) 00055 m_ptr[i][j]=column_indices[current_index++]; 00056 } 00057 } 00058 00060 ~SparsityStructure() 00061 { 00062 for (index_t i=0; i<m_num_rows; ++i) 00063 delete [] m_ptr[i]; 00064 delete [] m_ptr; 00065 } 00066 00068 void display_sparsity_structure() 00069 { 00070 for (index_t i=0; i<m_num_rows; ++i) 00071 { 00072 index_t nnzs=m_ptr[i][0]; 00073 SG_SPRINT("Row number %d. Number of Non-zeros %d. Colums ", i, nnzs); 00074 for(index_t j=1; j<=nnzs; ++j) 00075 { 00076 SG_SPRINT("%d", m_ptr[i][j]); 00077 if (j<nnzs) 00078 SG_SPRINT(", "); 00079 } 00080 SG_SPRINT("\n"); 00081 } 00082 } 00083 00085 index_t m_num_rows; 00086 00088 int32_t **m_ptr; 00089 }; 00090 00091 00098 template<class T> class CSparseMatrixOperator : public CMatrixOperator<T> 00099 { 00101 typedef bool supports_complex128_t; 00102 00103 public: 00105 CSparseMatrixOperator(); 00106 00112 explicit CSparseMatrixOperator(SGSparseMatrix<T> op); 00113 00119 CSparseMatrixOperator(const CSparseMatrixOperator<T>& orig); 00120 00122 ~CSparseMatrixOperator(); 00123 00130 virtual SGVector<T> apply(SGVector<T> b) const; 00131 00137 virtual void set_diagonal(SGVector<T> diag); 00138 00144 virtual SGVector<T> get_diagonal() const; 00145 00147 SGSparseMatrix<T> get_matrix_operator() const; 00148 00150 SparsityStructure* get_sparsity_structure(int64_t power=1) const; 00151 00155 template<class Scalar> 00156 inline operator CSparseMatrixOperator<Scalar>*() const 00157 { 00158 REQUIRE(m_operator.sparse_matrix, "Matrix is not initialized!\n"); 00159 typedef SGSparseVector<Scalar> vector; 00160 typedef SGSparseVectorEntry<Scalar> entry; 00161 00162 vector* rows=SG_MALLOC(vector, m_operator.num_vectors); 00163 00164 for (index_t i=0; i<m_operator.num_vectors; ++i) 00165 { 00166 entry* features=SG_MALLOC(entry, m_operator[i].num_feat_entries); 00167 for (index_t j=0; j<m_operator[i].num_feat_entries; ++j) 00168 { 00169 features[j].feat_index=m_operator[i].features[j].feat_index; 00170 features[j].entry=static_cast<Scalar>(m_operator[i].features[j].entry); 00171 } 00172 rows[i].features=features; 00173 rows[i].num_feat_entries=m_operator[i].num_feat_entries; 00174 } 00175 00176 SGSparseMatrix<Scalar> casted_m; 00177 casted_m.sparse_matrix=rows; 00178 casted_m.num_vectors=m_operator.num_vectors; 00179 casted_m.num_features= m_operator.num_features; 00180 00181 return new CSparseMatrixOperator<Scalar>(casted_m); 00182 } 00183 00185 virtual const char* get_name() const 00186 { 00187 return "SparseMatrixOperator"; 00188 } 00189 00190 private: 00192 SGSparseMatrix<T> m_operator; 00193 00195 void init(); 00196 00197 }; 00198 00199 } 00200 #endif // SPARSE_MATRIX_OPERATOR_H_