$treeview $search $mathjax
Eigen-unsupported
3.2.5
$projectbrief
|
$projectbrief
|
$searchbox |
00001 00002 // This file is part of Eigen, a lightweight C++ template library 00003 // for linear algebra. 00004 // 00005 // Copyright (C) 2012 Desire NUENTSA WAKAM <desire.nuentsa_wakam@inria.fr> 00006 // 00007 // This Source Code Form is subject to the terms of the Mozilla 00008 // Public License v. 2.0. If a copy of the MPL was not distributed 00009 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00010 00011 #ifndef EIGEN_BROWSE_MATRICES_H 00012 #define EIGEN_BROWSE_MATRICES_H 00013 00014 namespace Eigen { 00015 00016 enum { 00017 SPD = 0x100, 00018 NonSymmetric = 0x0 00019 }; 00020 00041 template <typename Scalar> 00042 class MatrixMarketIterator 00043 { 00044 public: 00045 typedef Matrix<Scalar,Dynamic,1> VectorType; 00046 typedef SparseMatrix<Scalar,ColMajor> MatrixType; 00047 00048 public: 00049 MatrixMarketIterator(const std::string folder):m_sym(0),m_isvalid(false),m_matIsLoaded(false),m_hasRhs(false),m_hasrefX(false),m_folder(folder) 00050 { 00051 m_folder_id = opendir(folder.c_str()); 00052 if (!m_folder_id){ 00053 m_isvalid = false; 00054 std::cerr << "The provided Matrix folder could not be opened \n\n"; 00055 abort(); 00056 } 00057 Getnextvalidmatrix(); 00058 } 00059 00060 ~MatrixMarketIterator() 00061 { 00062 if (m_folder_id) closedir(m_folder_id); 00063 } 00064 00065 inline MatrixMarketIterator& operator++() 00066 { 00067 m_matIsLoaded = false; 00068 m_hasrefX = false; 00069 m_hasRhs = false; 00070 Getnextvalidmatrix(); 00071 return *this; 00072 } 00073 inline operator bool() const { return m_isvalid;} 00074 00076 inline MatrixType& matrix() 00077 { 00078 // Read the matrix 00079 if (m_matIsLoaded) return m_mat; 00080 00081 std::string matrix_file = m_folder + "/" + m_matname + ".mtx"; 00082 if ( !loadMarket(m_mat, matrix_file)) 00083 { 00084 m_matIsLoaded = false; 00085 return m_mat; 00086 } 00087 m_matIsLoaded = true; 00088 00089 if (m_sym != NonSymmetric) 00090 { // Store the upper part of the matrix. It is needed by the solvers dealing with nonsymmetric matrices ?? 00091 MatrixType B; 00092 B = m_mat; 00093 m_mat = B.template selfadjointView<Lower>(); 00094 } 00095 return m_mat; 00096 } 00097 00101 inline VectorType& rhs() 00102 { 00103 // Get the right hand side 00104 if (m_hasRhs) return m_rhs; 00105 00106 std::string rhs_file; 00107 rhs_file = m_folder + "/" + m_matname + "_b.mtx"; // The pattern is matname_b.mtx 00108 m_hasRhs = Fileexists(rhs_file); 00109 if (m_hasRhs) 00110 { 00111 m_rhs.resize(m_mat.cols()); 00112 m_hasRhs = loadMarketVector(m_rhs, rhs_file); 00113 } 00114 if (!m_hasRhs) 00115 { 00116 // Generate a random right hand side 00117 if (!m_matIsLoaded) this->matrix(); 00118 m_refX.resize(m_mat.cols()); 00119 m_refX.setRandom(); 00120 m_rhs = m_mat * m_refX; 00121 m_hasrefX = true; 00122 m_hasRhs = true; 00123 } 00124 return m_rhs; 00125 } 00126 00133 inline VectorType& refX() 00134 { 00135 // Check if a reference solution is provided 00136 if (m_hasrefX) return m_refX; 00137 00138 std::string lhs_file; 00139 lhs_file = m_folder + "/" + m_matname + "_x.mtx"; 00140 m_hasrefX = Fileexists(lhs_file); 00141 if (m_hasrefX) 00142 { 00143 m_refX.resize(m_mat.cols()); 00144 m_hasrefX = loadMarketVector(m_refX, lhs_file); 00145 } 00146 return m_refX; 00147 } 00148 00149 inline std::string& matname() { return m_matname; } 00150 00151 inline int sym() { return m_sym; } 00152 00153 inline bool hasRhs() {return m_hasRhs; } 00154 inline bool hasrefX() {return m_hasrefX; } 00155 00156 protected: 00157 00158 inline bool Fileexists(std::string file) 00159 { 00160 std::ifstream file_id(file.c_str()); 00161 if (!file_id.good() ) 00162 { 00163 return false; 00164 } 00165 else 00166 { 00167 file_id.close(); 00168 return true; 00169 } 00170 } 00171 00172 void Getnextvalidmatrix( ) 00173 { 00174 m_isvalid = false; 00175 // Here, we return with the next valid matrix in the folder 00176 while ( (m_curs_id = readdir(m_folder_id)) != NULL) { 00177 m_isvalid = false; 00178 std::string curfile; 00179 curfile = m_folder + "/" + m_curs_id->d_name; 00180 // Discard if it is a folder 00181 if (m_curs_id->d_type == DT_DIR) continue; //FIXME This may not be available on non BSD systems 00182 // struct stat st_buf; 00183 // stat (curfile.c_str(), &st_buf); 00184 // if (S_ISDIR(st_buf.st_mode)) continue; 00185 00186 // Determine from the header if it is a matrix or a right hand side 00187 bool isvector,iscomplex=false; 00188 if(!getMarketHeader(curfile,m_sym,iscomplex,isvector)) continue; 00189 if(isvector) continue; 00190 if (!iscomplex) 00191 { 00192 if(internal::is_same<Scalar, std::complex<float> >::value || internal::is_same<Scalar, std::complex<double> >::value) 00193 continue; 00194 } 00195 if (iscomplex) 00196 { 00197 if(internal::is_same<Scalar, float>::value || internal::is_same<Scalar, double>::value) 00198 continue; 00199 } 00200 00201 00202 // Get the matrix name 00203 std::string filename = m_curs_id->d_name; 00204 m_matname = filename.substr(0, filename.length()-4); 00205 00206 // Find if the matrix is SPD 00207 size_t found = m_matname.find("SPD"); 00208 if( (found!=std::string::npos) && (m_sym != NonSymmetric) ) 00209 m_sym = SPD; 00210 00211 m_isvalid = true; 00212 break; 00213 } 00214 } 00215 int m_sym; // Symmetry of the matrix 00216 MatrixType m_mat; // Current matrix 00217 VectorType m_rhs; // Current vector 00218 VectorType m_refX; // The reference solution, if exists 00219 std::string m_matname; // Matrix Name 00220 bool m_isvalid; 00221 bool m_matIsLoaded; // Determine if the matrix has already been loaded from the file 00222 bool m_hasRhs; // The right hand side exists 00223 bool m_hasrefX; // A reference solution is provided 00224 std::string m_folder; 00225 DIR * m_folder_id; 00226 struct dirent *m_curs_id; 00227 00228 }; 00229 00230 } // end namespace Eigen 00231 00232 #endif