diff options
Diffstat (limited to 'unsupported/Eigen/src/SparseExtra/MarketIO.h')
-rw-r--r-- | unsupported/Eigen/src/SparseExtra/MarketIO.h | 273 |
1 files changed, 273 insertions, 0 deletions
diff --git a/unsupported/Eigen/src/SparseExtra/MarketIO.h b/unsupported/Eigen/src/SparseExtra/MarketIO.h new file mode 100644 index 000000000..de958de9f --- /dev/null +++ b/unsupported/Eigen/src/SparseExtra/MarketIO.h @@ -0,0 +1,273 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. +// +// Copyright (C) 2011 Gael Guennebaud <gael.guennebaud@inria.fr> +// Copyright (C) 2012 Desire NUENTSA WAKAM <desire.nuentsa_wakam@inria.fr> +// +// This Source Code Form is subject to the terms of the Mozilla +// Public License v. 2.0. If a copy of the MPL was not distributed +// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. + +#ifndef EIGEN_SPARSE_MARKET_IO_H +#define EIGEN_SPARSE_MARKET_IO_H + +#include <iostream> + +namespace Eigen { + +namespace internal +{ + template <typename Scalar> + inline bool GetMarketLine (std::stringstream& line, int& M, int& N, int& i, int& j, Scalar& value) + { + line >> i >> j >> value; + i--; + j--; + if(i>=0 && j>=0 && i<M && j<N) + { + return true; + } + else + return false; + } + template <typename Scalar> + inline bool GetMarketLine (std::stringstream& line, int& M, int& N, int& i, int& j, std::complex<Scalar>& value) + { + Scalar valR, valI; + line >> i >> j >> valR >> valI; + i--; + j--; + if(i>=0 && j>=0 && i<M && j<N) + { + value = std::complex<Scalar>(valR, valI); + return true; + } + else + return false; + } + + template <typename RealScalar> + inline void GetVectorElt (const std::string& line, RealScalar& val) + { + std::istringstream newline(line); + newline >> val; + } + + template <typename RealScalar> + inline void GetVectorElt (const std::string& line, std::complex<RealScalar>& val) + { + RealScalar valR, valI; + std::istringstream newline(line); + newline >> valR >> valI; + val = std::complex<RealScalar>(valR, valI); + } + + template<typename Scalar> + inline void putMarketHeader(std::string& header,int sym) + { + header= "%%MatrixMarket matrix coordinate "; + if(internal::is_same<Scalar, std::complex<float> >::value || internal::is_same<Scalar, std::complex<double> >::value) + { + header += " complex"; + if(sym == Symmetric) header += " symmetric"; + else if (sym == SelfAdjoint) header += " Hermitian"; + else header += " general"; + } + else + { + header += " real"; + if(sym == Symmetric) header += " symmetric"; + else header += " general"; + } + } + + template<typename Scalar> + inline void PutMatrixElt(Scalar value, int row, int col, std::ofstream& out) + { + out << row << " "<< col << " " << value << "\n"; + } + template<typename Scalar> + inline void PutMatrixElt(std::complex<Scalar> value, int row, int col, std::ofstream& out) + { + out << row << " " << col << " " << value.real() << " " << value.imag() << "\n"; + } + + + template<typename Scalar> + inline void putVectorElt(Scalar value, std::ofstream& out) + { + out << value << "\n"; + } + template<typename Scalar> + inline void putVectorElt(std::complex<Scalar> value, std::ofstream& out) + { + out << value.real << " " << value.imag()<< "\n"; + } + +} // end namepsace internal + +inline bool getMarketHeader(const std::string& filename, int& sym, bool& iscomplex, bool& isvector) +{ + sym = 0; + isvector = false; + std::ifstream in(filename.c_str(),std::ios::in); + if(!in) + return false; + + std::string line; + // The matrix header is always the first line in the file + std::getline(in, line); assert(in.good()); + + std::stringstream fmtline(line); + std::string substr[5]; + fmtline>> substr[0] >> substr[1] >> substr[2] >> substr[3] >> substr[4]; + if(substr[2].compare("array") == 0) isvector = true; + if(substr[3].compare("complex") == 0) iscomplex = true; + if(substr[4].compare("symmetric") == 0) sym = Symmetric; + else if (substr[4].compare("Hermitian") == 0) sym = SelfAdjoint; + + return true; +} + +template<typename SparseMatrixType> +bool loadMarket(SparseMatrixType& mat, const std::string& filename) +{ + typedef typename SparseMatrixType::Scalar Scalar; + std::ifstream input(filename.c_str(),std::ios::in); + if(!input) + return false; + + const int maxBuffersize = 2048; + char buffer[maxBuffersize]; + + bool readsizes = false; + + typedef Triplet<Scalar,int> T; + std::vector<T> elements; + + int M(-1), N(-1), NNZ(-1); + int count = 0; + while(input.getline(buffer, maxBuffersize)) + { + // skip comments + //NOTE An appropriate test should be done on the header to get the symmetry + if(buffer[0]=='%') + continue; + + std::stringstream line(buffer); + + if(!readsizes) + { + line >> M >> N >> NNZ; + if(M > 0 && N > 0 && NNZ > 0) + { + readsizes = true; + std::cout << "sizes: " << M << "," << N << "," << NNZ << "\n"; + mat.resize(M,N); + mat.reserve(NNZ); + } + } + else + { + int i(-1), j(-1); + Scalar value; + if( internal::GetMarketLine(line, M, N, i, j, value) ) + { + ++ count; + elements.push_back(T(i,j,value)); + } + else + std::cerr << "Invalid read: " << i << "," << j << "\n"; + } + } + mat.setFromTriplets(elements.begin(), elements.end()); + if(count!=NNZ) + std::cerr << count << "!=" << NNZ << "\n"; + + input.close(); + return true; +} + +template<typename VectorType> +bool loadMarketVector(VectorType& vec, const std::string& filename) +{ + typedef typename VectorType::Scalar Scalar; + std::ifstream in(filename.c_str(), std::ios::in); + if(!in) + return false; + + std::string line; + int n(0), col(0); + do + { // Skip comments + std::getline(in, line); assert(in.good()); + } while (line[0] == '%'); + std::istringstream newline(line); + newline >> n >> col; + assert(n>0 && col>0); + vec.resize(n); + int i = 0; + Scalar value; + while ( std::getline(in, line) && (i < n) ){ + internal::GetVectorElt(line, value); + vec(i++) = value; + } + in.close(); + if (i!=n){ + std::cerr<< "Unable to read all elements from file " << filename << "\n"; + return false; + } + return true; +} + +template<typename SparseMatrixType> +bool saveMarket(const SparseMatrixType& mat, const std::string& filename, int sym = 0) +{ + typedef typename SparseMatrixType::Scalar Scalar; + std::ofstream out(filename.c_str(),std::ios::out); + if(!out) + return false; + + out.flags(std::ios_base::scientific); + out.precision(64); + std::string header; + internal::putMarketHeader<Scalar>(header, sym); + out << header << std::endl; + out << mat.rows() << " " << mat.cols() << " " << mat.nonZeros() << "\n"; + int count = 0; + for(int j=0; j<mat.outerSize(); ++j) + for(typename SparseMatrixType::InnerIterator it(mat,j); it; ++it) + { + ++ count; + internal::PutMatrixElt(it.value(), it.row()+1, it.col()+1, out); + // out << it.row()+1 << " " << it.col()+1 << " " << it.value() << "\n"; + } + out.close(); + return true; +} + +template<typename VectorType> +bool saveMarketVector (const VectorType& vec, const std::string& filename) +{ + typedef typename VectorType::Scalar Scalar; + std::ofstream out(filename.c_str(),std::ios::out); + if(!out) + return false; + + out.flags(std::ios_base::scientific); + out.precision(64); + if(internal::is_same<Scalar, std::complex<float> >::value || internal::is_same<Scalar, std::complex<double> >::value) + out << "%%MatrixMarket matrix array complex general\n"; + else + out << "%%MatrixMarket matrix array real general\n"; + out << vec.size() << " "<< 1 << "\n"; + for (int i=0; i < vec.size(); i++){ + internal::putVectorElt(vec(i), out); + } + out.close(); + return true; +} + +} // end namespace Eigen + +#endif // EIGEN_SPARSE_MARKET_IO_H |