diff --git a/unsupported/Eigen/src/SparseExtra/MarketIO.h b/unsupported/Eigen/src/SparseExtra/MarketIO.h index 100e617b2..cdc14f86e 100644 --- a/unsupported/Eigen/src/SparseExtra/MarketIO.h +++ b/unsupported/Eigen/src/SparseExtra/MarketIO.h @@ -163,7 +163,7 @@ bool loadMarket(SparseMatrixType& mat, const std::string& filename) if(M > 0 && N > 0 && NNZ > 0) { readsizes = true; - std::cout << "sizes: " << M << "," << N << "," << NNZ << "\n"; + //std::cout << "sizes: " << M << "," << N << "," << NNZ << "\n"; mat.resize(M,N); mat.reserve(NNZ); } diff --git a/unsupported/Eigen/src/SparseExtra/MatrixMarketIterator.h b/unsupported/Eigen/src/SparseExtra/MatrixMarketIterator.h index bf13cf21f..4a40b26be 100644 --- a/unsupported/Eigen/src/SparseExtra/MatrixMarketIterator.h +++ b/unsupported/Eigen/src/SparseExtra/MatrixMarketIterator.h @@ -41,6 +41,7 @@ enum { template class MatrixMarketIterator { + typedef typename NumTraits::Real RealScalar; public: typedef Matrix VectorType; typedef SparseMatrix MatrixType; @@ -81,16 +82,30 @@ class MatrixMarketIterator std::string matrix_file = m_folder + "/" + m_matname + ".mtx"; if ( !loadMarket(m_mat, matrix_file)) { + std::cerr << "Warning loadMarket failed when loading \"" << matrix_file << "\"\n"; m_matIsLoaded = false; return m_mat; } m_matIsLoaded = true; - + if (m_sym != NonSymmetric) - { // Store the upper part of the matrix. It is needed by the solvers dealing with nonsymmetric matrices ?? - MatrixType B; - B = m_mat; - m_mat = B.template selfadjointView(); + { + // Check whether we need to restore a full matrix: + RealScalar diag_norm = m_mat.diagonal().norm(); + RealScalar lower_norm = m_mat.template triangularView().norm(); + RealScalar upper_norm = m_mat.template triangularView().norm(); + if(lower_norm>diag_norm && upper_norm==diag_norm) + { + // only the lower part is stored + MatrixType tmp(m_mat); + m_mat = tmp.template selfadjointView(); + } + else if(upper_norm>diag_norm && lower_norm==diag_norm) + { + // only the upper part is stored + MatrixType tmp(m_mat); + m_mat = tmp.template selfadjointView(); + } } return m_mat; } @@ -143,6 +158,8 @@ class MatrixMarketIterator m_refX.resize(m_mat.cols()); m_hasrefX = loadMarketVector(m_refX, lhs_file); } + else + m_refX.resize(0); return m_refX; }