00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 #ifndef EIGEN_SPARSE_MARKET_IO_H
00027 #define EIGEN_SPARSE_MARKET_IO_H
00028
00029 #include <iostream>
00030
00031 namespace Eigen {
00032
00033 namespace internal
00034 {
00035 template <typename Scalar>
00036 inline bool GetMarketLine (std::stringstream& line, int& M, int& N, int& i, int& j, Scalar& value)
00037 {
00038 line >> i >> j >> value;
00039 i--;
00040 j--;
00041 if(i>=0 && j>=0 && i<M && j<N)
00042 {
00043 return true;
00044 }
00045 else
00046 return false;
00047 }
00048 template <typename Scalar>
00049 inline bool GetMarketLine (std::stringstream& line, int& M, int& N, int& i, int& j, std::complex<Scalar>& value)
00050 {
00051 Scalar valR, valI;
00052 line >> i >> j >> valR >> valI;
00053 i--;
00054 j--;
00055 if(i>=0 && j>=0 && i<M && j<N)
00056 {
00057 value = std::complex<Scalar>(valR, valI);
00058 return true;
00059 }
00060 else
00061 return false;
00062 }
00063
00064 template <typename RealScalar>
00065 inline void GetVectorElt (const std::string& line, RealScalar& val)
00066 {
00067 std::istringstream newline(line);
00068 newline >> val;
00069 }
00070
00071 template <typename RealScalar>
00072 inline void GetVectorElt (const std::string& line, std::complex<RealScalar>& val)
00073 {
00074 RealScalar valR, valI;
00075 std::istringstream newline(line);
00076 newline >> valR >> valI;
00077 val = std::complex<RealScalar>(valR, valI);
00078 }
00079
00080 template<typename Scalar>
00081 inline void putMarketHeader(std::string& header,int sym)
00082 {
00083 header= "%%MatrixMarket matrix coordinate ";
00084 if(internal::is_same<Scalar, std::complex<float> >::value || internal::is_same<Scalar, std::complex<double> >::value)
00085 {
00086 header += " complex";
00087 if(sym == Symmetric) header += " symmetric";
00088 else if (sym == SelfAdjoint) header += " Hermitian";
00089 else header += " general";
00090 }
00091 else
00092 {
00093 header += " real";
00094 if(sym == Symmetric) header += " symmetric";
00095 else header += " general";
00096 }
00097 }
00098
00099 template<typename Scalar>
00100 inline void PutMatrixElt(Scalar value, int row, int col, std::ofstream& out)
00101 {
00102 out << row << " "<< col << " " << value << "\n";
00103 }
00104 template<typename Scalar>
00105 inline void PutMatrixElt(std::complex<Scalar> value, int row, int col, std::ofstream& out)
00106 {
00107 out << row << " " << col << " " << value.real() << " " << value.imag() << "\n";
00108 }
00109
00110
00111 template<typename Scalar>
00112 inline void putVectorElt(Scalar value, std::ofstream& out)
00113 {
00114 out << value << "\n";
00115 }
00116 template<typename Scalar>
00117 inline void putVectorElt(std::complex<Scalar> value, std::ofstream& out)
00118 {
00119 out << value.real << " " << value.imag()<< "\n";
00120 }
00121
00122 }
00123
00124 inline bool getMarketHeader(const std::string& filename, int& sym, bool& iscomplex, bool& isvector)
00125 {
00126 sym = 0;
00127 isvector = false;
00128 std::ifstream in(filename.c_str(),std::ios::in);
00129 if(!in)
00130 return false;
00131
00132 std::string line;
00133
00134 std::getline(in, line); assert(in.good());
00135
00136 std::stringstream fmtline(line);
00137 std::string substr[5];
00138 fmtline>> substr[0] >> substr[1] >> substr[2] >> substr[3] >> substr[4];
00139 if(substr[2].compare("array") == 0) isvector = true;
00140 if(substr[3].compare("complex") == 0) iscomplex = true;
00141 if(substr[4].compare("symmetric") == 0) sym = Symmetric;
00142 else if (substr[4].compare("Hermitian") == 0) sym = SelfAdjoint;
00143
00144 return true;
00145 }
00146
00147 template<typename SparseMatrixType>
00148 bool loadMarket(SparseMatrixType& mat, const std::string& filename)
00149 {
00150 typedef typename SparseMatrixType::Scalar Scalar;
00151 std::ifstream input(filename.c_str(),std::ios::in);
00152 if(!input)
00153 return false;
00154
00155 const int maxBuffersize = 2048;
00156 char buffer[maxBuffersize];
00157
00158 bool readsizes = false;
00159
00160 int M(-1), N(-1), NNZ(-1);
00161 int count = 0;
00162 while(input.getline(buffer, maxBuffersize))
00163 {
00164
00165
00166 if(buffer[0]=='%')
00167 continue;
00168
00169 std::stringstream line(buffer);
00170
00171 if(!readsizes)
00172 {
00173 line >> M >> N >> NNZ;
00174 if(M > 0 && N > 0 && NNZ > 0)
00175 {
00176 readsizes = true;
00177 std::cout << "sizes: " << M << "," << N << "," << NNZ << "\n";
00178 mat.resize(M,N);
00179 mat.reserve(NNZ);
00180 }
00181 }
00182 else
00183 {
00184 int i(-1), j(-1);
00185 Scalar value;
00186 if( internal::GetMarketLine(line, M, N, i, j, value) )
00187 {
00188 ++ count;
00189 mat.insert(i,j) = value;
00190 }
00191 else
00192 std::cerr << "Invalid read: " << i << "," << j << "\n";
00193 }
00194 }
00195 mat.makeCompressed();
00196 if(count!=NNZ)
00197 std::cerr << count << "!=" << NNZ << "\n";
00198
00199 input.close();
00200 return true;
00201 }
00202
00203 template<typename VectorType>
00204 bool loadMarketVector(VectorType& vec, const std::string& filename)
00205 {
00206 typedef typename VectorType::Scalar Scalar;
00207 std::ifstream in(filename.c_str(), std::ios::in);
00208 if(!in)
00209 return false;
00210
00211 std::string line;
00212 int n(0), col(0);
00213 do
00214 {
00215 std::getline(in, line); assert(in.good());
00216 } while (line[0] == '%');
00217 std::istringstream newline(line);
00218 newline >> n >> col;
00219 assert(n>0 && col>0);
00220 vec.resize(n);
00221 int i = 0;
00222 Scalar value;
00223 while ( std::getline(in, line) && (i < n) ){
00224 internal::GetVectorElt(line, value);
00225 vec(i++) = value;
00226 }
00227 in.close();
00228 if (i!=n){
00229 std::cerr<< "Unable to read all elements from file " << filename << "\n";
00230 return false;
00231 }
00232 return true;
00233 }
00234
00235 template<typename SparseMatrixType>
00236 bool saveMarket(const SparseMatrixType& mat, const std::string& filename, int sym = 0)
00237 {
00238 typedef typename SparseMatrixType::Scalar Scalar;
00239 std::ofstream out(filename.c_str(),std::ios::out);
00240 if(!out)
00241 return false;
00242
00243 out.flags(std::ios_base::scientific);
00244 out.precision(64);
00245 std::string header;
00246 internal::putMarketHeader<Scalar>(header, sym);
00247 out << header << std::endl;
00248 out << mat.rows() << " " << mat.cols() << " " << mat.nonZeros() << "\n";
00249 int count = 0;
00250 for(int j=0; j<mat.outerSize(); ++j)
00251 for(typename SparseMatrixType::InnerIterator it(mat,j); it; ++it)
00252 {
00253 ++ count;
00254 internal::PutMatrixElt(it.value(), it.row()+1, it.col()+1, out);
00255
00256 }
00257 out.close();
00258 return true;
00259 }
00260
00261 template<typename VectorType>
00262 bool saveMarketVector (const VectorType& vec, const std::string& filename)
00263 {
00264 typedef typename VectorType::Scalar Scalar;
00265 std::ofstream out(filename.c_str(),std::ios::out);
00266 if(!out)
00267 return false;
00268
00269 out.flags(std::ios_base::scientific);
00270 out.precision(64);
00271 if(internal::is_same<Scalar, std::complex<float> >::value || internal::is_same<Scalar, std::complex<double> >::value)
00272 out << "%%MatrixMarket matrix array complex general\n";
00273 else
00274 out << "%%MatrixMarket matrix array real general\n";
00275 out << vec.size() << " "<< 1 << "\n";
00276 for (int i=0; i < vec.size(); i++){
00277 internal::putVectorElt(vec(i), out);
00278 }
00279 out.close();
00280 return true;
00281 }
00282
00283 }
00284
00285 #endif // EIGEN_SPARSE_MARKET_IO_H