[136] | 1 | // This file is part of Eigen, a lightweight C++ template library
|
---|
| 2 | // for linear algebra.
|
---|
| 3 | //
|
---|
| 4 | // Copyright (C) 2011 Gael Guennebaud <gael.guennebaud@inria.fr>
|
---|
| 5 | // Copyright (C) 2012 Desire NUENTSA WAKAM <desire.nuentsa_wakam@inria.fr>
|
---|
| 6 | //
|
---|
| 7 | // This Source Code Form is subject to the terms of the Mozilla
|
---|
| 8 | // Public License v. 2.0. If a copy of the MPL was not distributed
|
---|
| 9 | // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
|
---|
| 10 |
|
---|
| 11 | #ifndef EIGEN_SPARSE_MARKET_IO_H
|
---|
| 12 | #define EIGEN_SPARSE_MARKET_IO_H
|
---|
| 13 |
|
---|
| 14 | #include <iostream>
|
---|
| 15 |
|
---|
| 16 | namespace Eigen {
|
---|
| 17 |
|
---|
| 18 | namespace internal
|
---|
| 19 | {
|
---|
| 20 | template <typename Scalar>
|
---|
| 21 | inline bool GetMarketLine (std::stringstream& line, int& M, int& N, int& i, int& j, Scalar& value)
|
---|
| 22 | {
|
---|
| 23 | line >> i >> j >> value;
|
---|
| 24 | i--;
|
---|
| 25 | j--;
|
---|
| 26 | if(i>=0 && j>=0 && i<M && j<N)
|
---|
| 27 | {
|
---|
| 28 | return true;
|
---|
| 29 | }
|
---|
| 30 | else
|
---|
| 31 | return false;
|
---|
| 32 | }
|
---|
| 33 | template <typename Scalar>
|
---|
| 34 | inline bool GetMarketLine (std::stringstream& line, int& M, int& N, int& i, int& j, std::complex<Scalar>& value)
|
---|
| 35 | {
|
---|
| 36 | Scalar valR, valI;
|
---|
| 37 | line >> i >> j >> valR >> valI;
|
---|
| 38 | i--;
|
---|
| 39 | j--;
|
---|
| 40 | if(i>=0 && j>=0 && i<M && j<N)
|
---|
| 41 | {
|
---|
| 42 | value = std::complex<Scalar>(valR, valI);
|
---|
| 43 | return true;
|
---|
| 44 | }
|
---|
| 45 | else
|
---|
| 46 | return false;
|
---|
| 47 | }
|
---|
| 48 |
|
---|
| 49 | template <typename RealScalar>
|
---|
| 50 | inline void GetVectorElt (const std::string& line, RealScalar& val)
|
---|
| 51 | {
|
---|
| 52 | std::istringstream newline(line);
|
---|
| 53 | newline >> val;
|
---|
| 54 | }
|
---|
| 55 |
|
---|
| 56 | template <typename RealScalar>
|
---|
| 57 | inline void GetVectorElt (const std::string& line, std::complex<RealScalar>& val)
|
---|
| 58 | {
|
---|
| 59 | RealScalar valR, valI;
|
---|
| 60 | std::istringstream newline(line);
|
---|
| 61 | newline >> valR >> valI;
|
---|
| 62 | val = std::complex<RealScalar>(valR, valI);
|
---|
| 63 | }
|
---|
| 64 |
|
---|
| 65 | template<typename Scalar>
|
---|
| 66 | inline void putMarketHeader(std::string& header,int sym)
|
---|
| 67 | {
|
---|
| 68 | header= "%%MatrixMarket matrix coordinate ";
|
---|
| 69 | if(internal::is_same<Scalar, std::complex<float> >::value || internal::is_same<Scalar, std::complex<double> >::value)
|
---|
| 70 | {
|
---|
| 71 | header += " complex";
|
---|
| 72 | if(sym == Symmetric) header += " symmetric";
|
---|
| 73 | else if (sym == SelfAdjoint) header += " Hermitian";
|
---|
| 74 | else header += " general";
|
---|
| 75 | }
|
---|
| 76 | else
|
---|
| 77 | {
|
---|
| 78 | header += " real";
|
---|
| 79 | if(sym == Symmetric) header += " symmetric";
|
---|
| 80 | else header += " general";
|
---|
| 81 | }
|
---|
| 82 | }
|
---|
| 83 |
|
---|
| 84 | template<typename Scalar>
|
---|
| 85 | inline void PutMatrixElt(Scalar value, int row, int col, std::ofstream& out)
|
---|
| 86 | {
|
---|
| 87 | out << row << " "<< col << " " << value << "\n";
|
---|
| 88 | }
|
---|
| 89 | template<typename Scalar>
|
---|
| 90 | inline void PutMatrixElt(std::complex<Scalar> value, int row, int col, std::ofstream& out)
|
---|
| 91 | {
|
---|
| 92 | out << row << " " << col << " " << value.real() << " " << value.imag() << "\n";
|
---|
| 93 | }
|
---|
| 94 |
|
---|
| 95 |
|
---|
| 96 | template<typename Scalar>
|
---|
| 97 | inline void putVectorElt(Scalar value, std::ofstream& out)
|
---|
| 98 | {
|
---|
| 99 | out << value << "\n";
|
---|
| 100 | }
|
---|
| 101 | template<typename Scalar>
|
---|
| 102 | inline void putVectorElt(std::complex<Scalar> value, std::ofstream& out)
|
---|
| 103 | {
|
---|
| 104 | out << value.real << " " << value.imag()<< "\n";
|
---|
| 105 | }
|
---|
| 106 |
|
---|
| 107 | } // end namepsace internal
|
---|
| 108 |
|
---|
| 109 | inline bool getMarketHeader(const std::string& filename, int& sym, bool& iscomplex, bool& isvector)
|
---|
| 110 | {
|
---|
| 111 | sym = 0;
|
---|
| 112 | isvector = false;
|
---|
| 113 | std::ifstream in(filename.c_str(),std::ios::in);
|
---|
| 114 | if(!in)
|
---|
| 115 | return false;
|
---|
| 116 |
|
---|
| 117 | std::string line;
|
---|
| 118 | // The matrix header is always the first line in the file
|
---|
| 119 | std::getline(in, line); eigen_assert(in.good());
|
---|
| 120 |
|
---|
| 121 | std::stringstream fmtline(line);
|
---|
| 122 | std::string substr[5];
|
---|
| 123 | fmtline>> substr[0] >> substr[1] >> substr[2] >> substr[3] >> substr[4];
|
---|
| 124 | if(substr[2].compare("array") == 0) isvector = true;
|
---|
| 125 | if(substr[3].compare("complex") == 0) iscomplex = true;
|
---|
| 126 | if(substr[4].compare("symmetric") == 0) sym = Symmetric;
|
---|
| 127 | else if (substr[4].compare("Hermitian") == 0) sym = SelfAdjoint;
|
---|
| 128 |
|
---|
| 129 | return true;
|
---|
| 130 | }
|
---|
| 131 |
|
---|
| 132 | template<typename SparseMatrixType>
|
---|
| 133 | bool loadMarket(SparseMatrixType& mat, const std::string& filename)
|
---|
| 134 | {
|
---|
| 135 | typedef typename SparseMatrixType::Scalar Scalar;
|
---|
| 136 | std::ifstream input(filename.c_str(),std::ios::in);
|
---|
| 137 | if(!input)
|
---|
| 138 | return false;
|
---|
| 139 |
|
---|
| 140 | const int maxBuffersize = 2048;
|
---|
| 141 | char buffer[maxBuffersize];
|
---|
| 142 |
|
---|
| 143 | bool readsizes = false;
|
---|
| 144 |
|
---|
| 145 | typedef Triplet<Scalar,int> T;
|
---|
| 146 | std::vector<T> elements;
|
---|
| 147 |
|
---|
| 148 | int M(-1), N(-1), NNZ(-1);
|
---|
| 149 | int count = 0;
|
---|
| 150 | while(input.getline(buffer, maxBuffersize))
|
---|
| 151 | {
|
---|
| 152 | // skip comments
|
---|
| 153 | //NOTE An appropriate test should be done on the header to get the symmetry
|
---|
| 154 | if(buffer[0]=='%')
|
---|
| 155 | continue;
|
---|
| 156 |
|
---|
| 157 | std::stringstream line(buffer);
|
---|
| 158 |
|
---|
| 159 | if(!readsizes)
|
---|
| 160 | {
|
---|
| 161 | line >> M >> N >> NNZ;
|
---|
| 162 | if(M > 0 && N > 0 && NNZ > 0)
|
---|
| 163 | {
|
---|
| 164 | readsizes = true;
|
---|
| 165 | std::cout << "sizes: " << M << "," << N << "," << NNZ << "\n";
|
---|
| 166 | mat.resize(M,N);
|
---|
| 167 | mat.reserve(NNZ);
|
---|
| 168 | }
|
---|
| 169 | }
|
---|
| 170 | else
|
---|
| 171 | {
|
---|
| 172 | int i(-1), j(-1);
|
---|
| 173 | Scalar value;
|
---|
| 174 | if( internal::GetMarketLine(line, M, N, i, j, value) )
|
---|
| 175 | {
|
---|
| 176 | ++ count;
|
---|
| 177 | elements.push_back(T(i,j,value));
|
---|
| 178 | }
|
---|
| 179 | else
|
---|
| 180 | std::cerr << "Invalid read: " << i << "," << j << "\n";
|
---|
| 181 | }
|
---|
| 182 | }
|
---|
| 183 | mat.setFromTriplets(elements.begin(), elements.end());
|
---|
| 184 | if(count!=NNZ)
|
---|
| 185 | std::cerr << count << "!=" << NNZ << "\n";
|
---|
| 186 |
|
---|
| 187 | input.close();
|
---|
| 188 | return true;
|
---|
| 189 | }
|
---|
| 190 |
|
---|
| 191 | template<typename VectorType>
|
---|
| 192 | bool loadMarketVector(VectorType& vec, const std::string& filename)
|
---|
| 193 | {
|
---|
| 194 | typedef typename VectorType::Scalar Scalar;
|
---|
| 195 | std::ifstream in(filename.c_str(), std::ios::in);
|
---|
| 196 | if(!in)
|
---|
| 197 | return false;
|
---|
| 198 |
|
---|
| 199 | std::string line;
|
---|
| 200 | int n(0), col(0);
|
---|
| 201 | do
|
---|
| 202 | { // Skip comments
|
---|
| 203 | std::getline(in, line); eigen_assert(in.good());
|
---|
| 204 | } while (line[0] == '%');
|
---|
| 205 | std::istringstream newline(line);
|
---|
| 206 | newline >> n >> col;
|
---|
| 207 | eigen_assert(n>0 && col>0);
|
---|
| 208 | vec.resize(n);
|
---|
| 209 | int i = 0;
|
---|
| 210 | Scalar value;
|
---|
| 211 | while ( std::getline(in, line) && (i < n) ){
|
---|
| 212 | internal::GetVectorElt(line, value);
|
---|
| 213 | vec(i++) = value;
|
---|
| 214 | }
|
---|
| 215 | in.close();
|
---|
| 216 | if (i!=n){
|
---|
| 217 | std::cerr<< "Unable to read all elements from file " << filename << "\n";
|
---|
| 218 | return false;
|
---|
| 219 | }
|
---|
| 220 | return true;
|
---|
| 221 | }
|
---|
| 222 |
|
---|
| 223 | template<typename SparseMatrixType>
|
---|
| 224 | bool saveMarket(const SparseMatrixType& mat, const std::string& filename, int sym = 0)
|
---|
| 225 | {
|
---|
| 226 | typedef typename SparseMatrixType::Scalar Scalar;
|
---|
| 227 | std::ofstream out(filename.c_str(),std::ios::out);
|
---|
| 228 | if(!out)
|
---|
| 229 | return false;
|
---|
| 230 |
|
---|
| 231 | out.flags(std::ios_base::scientific);
|
---|
| 232 | out.precision(64);
|
---|
| 233 | std::string header;
|
---|
| 234 | internal::putMarketHeader<Scalar>(header, sym);
|
---|
| 235 | out << header << std::endl;
|
---|
| 236 | out << mat.rows() << " " << mat.cols() << " " << mat.nonZeros() << "\n";
|
---|
| 237 | int count = 0;
|
---|
| 238 | for(int j=0; j<mat.outerSize(); ++j)
|
---|
| 239 | for(typename SparseMatrixType::InnerIterator it(mat,j); it; ++it)
|
---|
| 240 | {
|
---|
| 241 | ++ count;
|
---|
| 242 | internal::PutMatrixElt(it.value(), it.row()+1, it.col()+1, out);
|
---|
| 243 | // out << it.row()+1 << " " << it.col()+1 << " " << it.value() << "\n";
|
---|
| 244 | }
|
---|
| 245 | out.close();
|
---|
| 246 | return true;
|
---|
| 247 | }
|
---|
| 248 |
|
---|
| 249 | template<typename VectorType>
|
---|
| 250 | bool saveMarketVector (const VectorType& vec, const std::string& filename)
|
---|
| 251 | {
|
---|
| 252 | typedef typename VectorType::Scalar Scalar;
|
---|
| 253 | std::ofstream out(filename.c_str(),std::ios::out);
|
---|
| 254 | if(!out)
|
---|
| 255 | return false;
|
---|
| 256 |
|
---|
| 257 | out.flags(std::ios_base::scientific);
|
---|
| 258 | out.precision(64);
|
---|
| 259 | if(internal::is_same<Scalar, std::complex<float> >::value || internal::is_same<Scalar, std::complex<double> >::value)
|
---|
| 260 | out << "%%MatrixMarket matrix array complex general\n";
|
---|
| 261 | else
|
---|
| 262 | out << "%%MatrixMarket matrix array real general\n";
|
---|
| 263 | out << vec.size() << " "<< 1 << "\n";
|
---|
| 264 | for (int i=0; i < vec.size(); i++){
|
---|
| 265 | internal::putVectorElt(vec(i), out);
|
---|
| 266 | }
|
---|
| 267 | out.close();
|
---|
| 268 | return true;
|
---|
| 269 | }
|
---|
| 270 |
|
---|
| 271 | } // end namespace Eigen
|
---|
| 272 |
|
---|
| 273 | #endif // EIGEN_SPARSE_MARKET_IO_H
|
---|