11 #ifndef EIGEN_SPARSE_MARKET_IO_H
12 #define EIGEN_SPARSE_MARKET_IO_H
17 #include "./InternalHeaderCheck.h"
23 template <
typename Scalar,
typename StorageIndex>
24 inline void GetMarketLine (
const char* line, StorageIndex& i, StorageIndex& j, Scalar& value)
26 std::stringstream sline(line);
27 sline >> i >> j >> value;
30 template<>
inline void GetMarketLine (
const char* line,
int& i,
int& j,
float& value)
31 { std::sscanf(line,
"%d %d %g", &i, &j, &value); }
33 template<>
inline void GetMarketLine (
const char* line,
int& i,
int& j,
double& value)
34 { std::sscanf(line,
"%d %d %lg", &i, &j, &value); }
36 template<>
inline void GetMarketLine (
const char* line,
int& i,
int& j, std::complex<float>& value)
37 { std::sscanf(line,
"%d %d %g %g", &i, &j, &numext::real_ref(value), &numext::imag_ref(value)); }
39 template<>
inline void GetMarketLine (
const char* line,
int& i,
int& j, std::complex<double>& value)
40 { std::sscanf(line,
"%d %d %lg %lg", &i, &j, &numext::real_ref(value), &numext::imag_ref(value)); }
42 template <
typename Scalar,
typename StorageIndex>
43 inline void GetMarketLine (
const char* line, StorageIndex& i, StorageIndex& j, std::complex<Scalar>& value)
45 std::stringstream sline(line);
47 sline >> i >> j >> valR >> valI;
48 value = std::complex<Scalar>(valR,valI);
51 template <
typename RealScalar>
52 inline void GetDenseElt (
const std::string& line, RealScalar& val)
54 std::istringstream newline(line);
58 template <
typename RealScalar>
59 inline void GetDenseElt (
const std::string& line, std::complex<RealScalar>& val)
61 RealScalar valR, valI;
62 std::istringstream newline(line);
63 newline >> valR >> valI;
64 val = std::complex<RealScalar>(valR, valI);
67 template<
typename Scalar>
68 inline void putMarketHeader(std::string& header,
int sym)
70 header=
"%%MatrixMarket matrix coordinate ";
71 if(internal::is_same<Scalar, std::complex<float> >::value || internal::is_same<Scalar, std::complex<double> >::value)
74 if(sym ==
Symmetric) header +=
" symmetric";
75 else if (sym ==
SelfAdjoint) header +=
" Hermitian";
76 else header +=
" general";
81 if(sym ==
Symmetric) header +=
" symmetric";
82 else header +=
" general";
86 template<
typename Scalar,
typename StorageIndex>
87 inline void PutMatrixElt(Scalar value, StorageIndex row, StorageIndex col, std::ofstream& out)
89 out << row <<
" "<< col <<
" " << value <<
"\n";
91 template<
typename Scalar,
typename StorageIndex>
92 inline void PutMatrixElt(std::complex<Scalar> value, StorageIndex row, StorageIndex col, std::ofstream& out)
94 out << row <<
" " << col <<
" " << value.real() <<
" " << value.imag() <<
"\n";
98 template<
typename Scalar>
99 inline void putDenseElt(Scalar value, std::ofstream& out)
101 out << value <<
"\n";
103 template<
typename Scalar>
104 inline void putDenseElt(std::complex<Scalar> value, std::ofstream& out)
106 out << value.real() <<
" " << value.imag()<<
"\n";
122 inline bool getMarketHeader(
const std::string& filename,
int& sym,
bool& iscomplex,
bool& isdense)
127 std::ifstream in(filename.c_str(),std::ios::in);
133 std::getline(in, line); eigen_assert(in.good());
135 std::stringstream fmtline(line);
136 std::string substr[5];
137 fmtline>> substr[0] >> substr[1] >> substr[2] >> substr[3] >> substr[4];
138 if(substr[2].compare(
"array") == 0) isdense =
true;
139 if(substr[3].compare(
"complex") == 0) iscomplex =
true;
140 if(substr[4].compare(
"symmetric") == 0) sym =
Symmetric;
141 else if (substr[4].compare(
"Hermitian") == 0) sym =
SelfAdjoint;
154 template<
typename SparseMatrixType>
155 bool loadMarket(SparseMatrixType& mat,
const std::string& filename)
157 typedef typename SparseMatrixType::Scalar Scalar;
158 typedef typename SparseMatrixType::StorageIndex StorageIndex;
159 std::ifstream input(filename.c_str(),std::ios::in);
164 input.rdbuf()->pubsetbuf(rdbuffer, 4096);
166 const int maxBuffersize = 2048;
167 char buffer[maxBuffersize];
169 bool readsizes =
false;
172 std::vector<T> elements;
174 Index M(-1), N(-1), NNZ(-1);
176 while(input.getline(buffer, maxBuffersize))
185 std::stringstream line(buffer);
186 line >> M >> N >> NNZ;
196 StorageIndex i(-1), j(-1);
198 internal::GetMarketLine(buffer, i, j, value);
202 if(i>=0 && j>=0 && i<M && j<N)
205 elements.push_back(T(i,j,value));
209 std::cerr <<
"Invalid read: " << i <<
"," << j <<
"\n";
215 mat.setFromTriplets(elements.begin(), elements.end());
217 std::cerr << count <<
"!=" << NNZ <<
"\n";
234 template<
typename DenseType>
237 typedef typename DenseType::Scalar Scalar;
238 std::ifstream in(filename.c_str(), std::ios::in);
243 Index rows(0), cols(0);
246 std::getline(in, line); eigen_assert(in.good());
247 }
while (line[0] ==
'%');
248 std::istringstream newline(line);
249 newline >> rows >> cols;
251 bool sizes_not_positive=(rows<1 || cols<1);
252 bool wrong_input_rows = (DenseType::MaxRowsAtCompileTime !=
Dynamic && rows > DenseType::MaxRowsAtCompileTime) ||
253 (DenseType::RowsAtCompileTime!=
Dynamic && rows!=DenseType::RowsAtCompileTime);
254 bool wrong_input_cols = (DenseType::MaxColsAtCompileTime !=
Dynamic && cols > DenseType::MaxColsAtCompileTime) ||
255 (DenseType::ColsAtCompileTime!=
Dynamic && cols!=DenseType::ColsAtCompileTime);
257 if(sizes_not_positive || wrong_input_rows || wrong_input_cols){
258 if(sizes_not_positive){
259 std::cerr<<
"non-positive row or column size in file" << filename <<
"\n";
261 std::cerr<<
"Input matrix can not be resized to"<<rows<<
" x "<<cols<<
"as given in " << filename <<
"\n";
267 mat.resize(rows,cols);
272 while ( std::getline(in, line) && (row < rows) && (col < cols)){
273 internal::GetDenseElt(line, value);
275 mat(row,col) = value;
285 std::cerr<<
"Unable to read all elements from file " << filename <<
"\n";
294 template<
typename VectorType>
310 template<
typename SparseMatrixType>
311 bool saveMarket(
const SparseMatrixType& mat,
const std::string& filename,
int sym = 0)
313 typedef typename SparseMatrixType::Scalar Scalar;
314 typedef typename SparseMatrixType::RealScalar RealScalar;
315 std::ofstream out(filename.c_str(),std::ios::out);
319 out.flags(std::ios_base::scientific);
320 out.precision(std::numeric_limits<RealScalar>::digits10 + 2);
322 internal::putMarketHeader<Scalar>(header, sym);
323 out << header << std::endl;
324 out << mat.rows() <<
" " << mat.cols() <<
" " << mat.nonZeros() <<
"\n";
326 for(
int j=0; j<mat.outerSize(); ++j)
327 for(
typename SparseMatrixType::InnerIterator it(mat,j); it; ++it)
330 internal::PutMatrixElt(it.value(), it.row()+1, it.col()+1, out);
347 template<
typename DenseType>
350 typedef typename DenseType::Scalar Scalar;
351 typedef typename DenseType::RealScalar RealScalar;
352 std::ofstream out(filename.c_str(),std::ios::out);
356 out.flags(std::ios_base::scientific);
357 out.precision(std::numeric_limits<RealScalar>::digits10 + 2);
358 if(internal::is_same<Scalar, std::complex<float> >::value || internal::is_same<Scalar, std::complex<double> >::value)
359 out <<
"%%MatrixMarket matrix array complex general\n";
361 out <<
"%%MatrixMarket matrix array real general\n";
362 out << mat.rows() <<
" "<< mat.cols() <<
"\n";
363 for (
Index i=0; i < mat.cols(); i++){
364 for (
Index j=0; j < mat.rows(); j++){
365 internal::putDenseElt(mat(j,i), out);
376 template<
typename VectorType>
Namespace containing all symbols from the Eigen library.
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index