Matrixmarket extension

unsupported/Eigen/SparseExtra

@@ -31,10 +31,12 @@
 /**
   * \defgroup SparseExtra_Module SparseExtra module
   *
-  * This module contains some experimental features extending the sparse module.
+  * This module contains some experimental features extending the sparse module:
+  * - A RandomSetter which is a wrapper object allowing to set/update a sparse matrix with random access.
+  * - MatrixMarket format(https://math.nist.gov/MatrixMarket/formats.html) readers and writers for sparse and dense matrices.
   *
   * \code
-  * #include <Eigen/SparseExtra>
+  * #include <unsupported/Eigen/SparseExtra>
   * \endcode
   */
 

unsupported/Eigen/src/SparseExtra/MarketIO.h

@@ -47,14 +47,14 @@ namespace internal
   }
 
   template <typename RealScalar>
-  inline void  GetVectorElt (const std::string& line, RealScalar& val)
+  inline void  GetDenseElt (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)
+  inline void GetDenseElt (const std::string& line, std::complex<RealScalar>& val)
   {
     RealScalar valR, valI; 
     std::istringstream newline(line);
@@ -94,23 +94,34 @@ namespace internal
 
 
   template<typename Scalar>
-  inline void putVectorElt(Scalar value, std::ofstream& out)
+  inline void putDenseElt(Scalar value, std::ofstream& out)
   {
     out << value << "\n"; 
   }
   template<typename Scalar>
-  inline void putVectorElt(std::complex<Scalar> value, std::ofstream& out)
+  inline void putDenseElt(std::complex<Scalar> value, std::ofstream& out)
   {
     out << value.real() << " " << value.imag()<< "\n"; 
   }
 
 } // end namespace internal
 
-inline bool getMarketHeader(const std::string& filename, int& sym, bool& iscomplex, bool& isvector)
+
+/**
+ * \ingroup SparseExtra_Module
+ * @brief Reads the header of a matrixmarket file and determines the properties of a matrix
+ * 
+ * @param filename of the file
+ * @param sym if the matrix is hermitian,symmetric or none of the latter (sym=0) 
+ * @param iscomplex if the matrix has complex or real coefficients 
+ * @param isdense if the matrix is dense or sparse
+ * @return true if the file was found
+ */
+inline bool getMarketHeader(const std::string& filename, int& sym, bool& iscomplex, bool& isdense)
 {
   sym = 0; 
   iscomplex = false;
-  isvector = false;
+  isdense = false;
   std::ifstream in(filename.c_str(),std::ios::in);
   if(!in)
     return false;
@@ -122,14 +133,22 @@ inline bool getMarketHeader(const std::string& filename, int& sym, bool& iscompl
   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[2].compare("array") == 0) isdense = 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;
 }
-  
+/**
+ * \ingroup SparseExtra_Module
+ * @brief Loads a sparse matrix from a matrixmarket format file.
+ * 
+ * @tparam SparseMatrixType to read into, symmetries are not supported
+ * @param mat SparseMatrix to read into, current values are overwritten
+ * @param filename to parse matrix from
+ * @return returns true if file exists. Returns false if the parsing did not suceed.
+ */
 template<typename SparseMatrixType>
 bool loadMarket(SparseMatrixType& mat, const std::string& filename)
 {
@@ -184,50 +203,108 @@ bool loadMarket(SparseMatrixType& mat, const std::string& filename)
         elements.push_back(T(i,j,value));
       }
       else
-        std::cerr << "Invalid read: " << i << "," << j << "\n";        
+      {
+        std::cerr << "Invalid read: " << i << "," << j << "\n";   
+        return false;
+      }     
     }
   }
 
   mat.setFromTriplets(elements.begin(), elements.end());
-  if(count!=NNZ)
+  if(count!=NNZ){
     std::cerr << count << "!=" << NNZ << "\n";
-  
+    return false;
+  }
   input.close();
   return true;
 }
 
-template<typename VectorType>
+
-bool loadMarketVector(VectorType& vec, const std::string& filename)
+/**
+ * \ingroup SparseExtra_Module
+ * @brief Loads a dense Matrix or Vector from a matrixmarket file. If a statically sized matrix has to be parsed and the file contains the wrong dimensions it is undefined behaviour.
+ * 
+ * @tparam DenseMatrixType to read into
+ * @param mat DenseMatrix to read into, current values are overwritten, symmetries are not supported
+ * @param filename to parse matrix from
+ * @return true if parsing was successful. Returns false if the parsing did not suceed.
+ */
+template<typename DenseType>
+bool loadMarketDense(DenseType& mat, const std::string& filename)
 {
-   typedef typename VectorType::Scalar Scalar;
+   typedef typename DenseType::Scalar Scalar;
   std::ifstream in(filename.c_str(), std::ios::in);
   if(!in)
     return false;
   
   std::string line; 
-  int n(0), col(0); 
+  Index rows(0), cols(0); 
   do 
   { // Skip comments
     std::getline(in, line); eigen_assert(in.good());
   } while (line[0] == '%');
   std::istringstream newline(line);
-  newline  >> n >> col; 
+  newline  >> rows >> cols; 
-  eigen_assert(n>0 && col>0);
+
-  vec.resize(n);
+  bool sizes_not_positive=(rows<1 || cols<1);
-  int i = 0; 
+  bool wrong_input_rows = (DenseType::MaxRowsAtCompileTime != Dynamic && rows > DenseType::MaxRowsAtCompileTime) ||
+                          (DenseType::RowsAtCompileTime!=Dynamic && rows!=DenseType::RowsAtCompileTime);
+  bool wrong_input_cols = (DenseType::MaxColsAtCompileTime != Dynamic && cols > DenseType::MaxColsAtCompileTime) ||
+                          (DenseType::ColsAtCompileTime!=Dynamic && cols!=DenseType::ColsAtCompileTime);
+
+  if(sizes_not_positive || wrong_input_rows || wrong_input_cols){
+    if(sizes_not_positive){
+      std::cerr<< "non-positive row or column size in file" << filename << "\n";
+    }else{
+      std::cerr<< "Input matrix can not be resized to"<<rows<<" x "<<cols<< "as given in " << filename << "\n";
+    }
+    in.close();
+    return false;
+  }
+
+  mat.resize(rows,cols);
+  Index row = 0;
+  Index col = 0; 
+  Index n=0;
   Scalar value; 
-  while ( std::getline(in, line) && (i < n) ){
+  while ( std::getline(in, line) && (row < rows) && (col < cols)){
-    internal::GetVectorElt(line, value); 
+    internal::GetDenseElt(line, value); 
-    vec(i++) = value; 
+    //matrixmarket format is column major
+    mat(row,col) = value; 
+    row++;
+    if(row==rows){
+      row=0;
+      col++;
+    }
+    n++;
   }
   in.close();
-  if (i!=n){
+  if (n!=mat.size()){
     std::cerr<< "Unable to read all elements from file " << filename << "\n";
     return false;
   }
   return true;
 }
+/**
+ * \ingroup SparseExtra_Module
+ * @brief Same functionality as loadMarketDense, deprecated
+ */
+template<typename VectorType>
+bool loadMarketVector(VectorType& vec, const std::string& filename)
+{
+ return loadMarketDense(vec, filename);
+}
 
+/**
+ * \ingroup SparseExtra_Module
+ * @brief writes a sparse Matrix to a marketmarket format file
+ * 
+ * @tparam SparseMatrixType to write to file
+ * @param mat matrix to write to file
+ * @param filename filename to write to 
+ * @param sym at the moment no symmetry operations are supported
+ * @return true if writing suceeded
+ */
 template<typename SparseMatrixType>
 bool saveMarket(const SparseMatrixType& mat, const std::string& filename, int sym = 0)
 {
@@ -254,11 +331,22 @@ bool saveMarket(const SparseMatrixType& mat, const std::string& filename, int sy
   return true;
 }
 
-template<typename VectorType>
+
-bool saveMarketVector (const VectorType& vec, const std::string& filename)
+/**
+ * \ingroup SparseExtra_Module
+ * @brief writes a dense Matrix or vector to a marketmarket format file
+ * 
+ * @tparam DenseMatrixType to write to file
+ * @param mat matrix to write to file
+ * @param filename filename to write to 
+ * @return true if writing suceeded
+ */
+
+template<typename DenseType>
+bool saveMarketDense (const DenseType& mat, const std::string& filename)
 {
- typedef typename VectorType::Scalar Scalar;
+ typedef typename DenseType::Scalar Scalar;
- typedef typename VectorType::RealScalar RealScalar;
+ typedef typename DenseType::RealScalar RealScalar;
  std::ofstream out(filename.c_str(),std::ios::out);
   if(!out)
     return false;
@@ -269,14 +357,26 @@ bool saveMarketVector (const VectorType& vec, const std::string& filename)
       out << "%%MatrixMarket matrix array complex general\n"; 
   else
     out << "%%MatrixMarket matrix array real general\n"; 
-  out << vec.size() << " "<< 1 << "\n";
+  out << mat.rows() << " "<< mat.cols() << "\n";
-  for (int i=0; i < vec.size(); i++){
+  for (Index i=0; i < mat.cols(); i++){
-    internal::putVectorElt(vec(i), out); 
+    for (Index j=0; j < mat.rows(); j++){
+      internal::putDenseElt(mat(j,i), out); 
+    }
   }
   out.close();
   return true; 
 }
 
+/**
+ * \ingroup SparseExtra_Module
+ * @brief Same functionality as saveMarketDense, deprecated
+ */
+template<typename VectorType>
+bool saveMarketVector (const VectorType& vec, const std::string& filename)
+{
+  return saveMarketDense(vec, filename);
+}
+
 } // end namespace Eigen
 
 #endif // EIGEN_SPARSE_MARKET_IO_H

unsupported/Eigen/src/SparseExtra/RandomSetter.h

@@ -101,7 +101,7 @@ template<typename Scalar> struct GoogleSparseHashMapTraits
 #endif
 
 /** \class RandomSetter
-  *
+  * \ingroup SparseExtra_Module
   * \brief The RandomSetter is a wrapper object allowing to set/update a sparse matrix with random access
   *
   * \tparam SparseMatrixType the type of the sparse matrix we are updating

unsupported/test/sparse_extra.cpp

@@ -173,6 +173,30 @@ void check_marketio_vector()
   VERIFY_IS_EQUAL(v1,v2);
 }
 
+template<typename DenseMatrixType>
+void check_marketio_dense()
+{
+  Index rows=DenseMatrixType::MaxRowsAtCompileTime;
+  if (DenseMatrixType::MaxRowsAtCompileTime==Dynamic){
+    rows=internal::random<Index>(1,100);
+  }else if(DenseMatrixType::RowsAtCompileTime==Dynamic){
+    rows=internal::random<Index>(1,DenseMatrixType::MaxRowsAtCompileTime);
+  }
+
+  Index cols =DenseMatrixType::MaxColsAtCompileTime; 
+  if (DenseMatrixType::MaxColsAtCompileTime==Dynamic){
+    cols=internal::random<Index>(1,100);
+  }else if(DenseMatrixType::ColsAtCompileTime==Dynamic){
+    cols=internal::random<Index>(1,DenseMatrixType::MaxColsAtCompileTime);
+  }
+
+  DenseMatrixType m1, m2;
+  m1= DenseMatrixType::Random(rows,cols);
+  saveMarketDense(m1, "dense_extra.mtx");
+  loadMarketDense(m2, "dense_extra.mtx");
+  VERIFY_IS_EQUAL(m1,m2);
+}
+
 EIGEN_DECLARE_TEST(sparse_extra)
 {
   for(int i = 0; i < g_repeat; i++) {
@@ -181,15 +205,24 @@ EIGEN_DECLARE_TEST(sparse_extra)
     CALL_SUBTEST_2( sparse_extra(SparseMatrix<std::complex<double> >(s, s)) );
     CALL_SUBTEST_1( sparse_extra(SparseMatrix<double>(s, s)) );
 
-    CALL_SUBTEST_4( (check_marketio<SparseMatrix<float,ColMajor,int> >()) );
+    CALL_SUBTEST_3( (check_marketio<SparseMatrix<float,ColMajor,int> >()) );
-    CALL_SUBTEST_4( (check_marketio<SparseMatrix<double,ColMajor,int> >()) );
+    CALL_SUBTEST_3( (check_marketio<SparseMatrix<double,ColMajor,int> >()) );
-    CALL_SUBTEST_4( (check_marketio<SparseMatrix<std::complex<float>,ColMajor,int> >()) );
+    CALL_SUBTEST_3( (check_marketio<SparseMatrix<std::complex<float>,ColMajor,int> >()) );
-    CALL_SUBTEST_4( (check_marketio<SparseMatrix<std::complex<double>,ColMajor,int> >()) );
+    CALL_SUBTEST_3( (check_marketio<SparseMatrix<std::complex<double>,ColMajor,int> >()) );
-    CALL_SUBTEST_4( (check_marketio<SparseMatrix<float,ColMajor,long int> >()) );
+    CALL_SUBTEST_3( (check_marketio<SparseMatrix<float,ColMajor,long int> >()) );
-    CALL_SUBTEST_4( (check_marketio<SparseMatrix<double,ColMajor,long int> >()) );
+    CALL_SUBTEST_3( (check_marketio<SparseMatrix<double,ColMajor,long int> >()) );
-    CALL_SUBTEST_4( (check_marketio<SparseMatrix<std::complex<float>,ColMajor,long int> >()) );
+    CALL_SUBTEST_3( (check_marketio<SparseMatrix<std::complex<float>,ColMajor,long int> >()) );
-    CALL_SUBTEST_4( (check_marketio<SparseMatrix<std::complex<double>,ColMajor,long int> >()) );
+    CALL_SUBTEST_3( (check_marketio<SparseMatrix<std::complex<double>,ColMajor,long int> >()) );
 
+    CALL_SUBTEST_4( (check_marketio_dense<Matrix<float,Dynamic,Dynamic> >()) );
+    CALL_SUBTEST_4( (check_marketio_dense<Matrix<float,Dynamic,Dynamic,RowMajor> >()) );
+    CALL_SUBTEST_4( (check_marketio_dense<Matrix<double,Dynamic,Dynamic> >()) );
+    CALL_SUBTEST_4( (check_marketio_dense<Matrix<std::complex<float>,Dynamic,Dynamic> >()) );
+    CALL_SUBTEST_4( (check_marketio_dense<Matrix<std::complex<double>,Dynamic,Dynamic> >()) );
+    CALL_SUBTEST_4( (check_marketio_dense<Matrix<float,Dynamic,3> >()) );
+    CALL_SUBTEST_4( (check_marketio_dense<Matrix<double,3,Dynamic> >()) );
+    CALL_SUBTEST_4( (check_marketio_dense<Matrix<double,3,4> >()) );
+    CALL_SUBTEST_4( (check_marketio_dense<Matrix<double,Dynamic,Dynamic,ColMajor,5,5> >()) );
 
     CALL_SUBTEST_5( (check_marketio_vector<Matrix<float,1,Dynamic> >()) );
     CALL_SUBTEST_5( (check_marketio_vector<Matrix<double,1,Dynamic> >()) );