diff options
Diffstat (limited to 'Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h')
-rw-r--r-- | Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h | 42 |
1 files changed, 23 insertions, 19 deletions
diff --git a/Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h b/Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h index 3eae95479..ad6f2183f 100644 --- a/Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h +++ b/Eigen/src/SparseLU/SparseLU_SupernodalMatrix.h @@ -189,13 +189,14 @@ class MappedSuperNodalMatrix<Scalar,Index>::InnerIterator m_idval(mat.colIndexPtr()[outer]), m_startidval(m_idval), m_endidval(mat.colIndexPtr()[outer+1]), - m_idrow(mat.rowIndexPtr()[outer]) + m_idrow(mat.rowIndexPtr()[outer]), + m_endidrow(mat.rowIndexPtr()[outer+1]) {} inline InnerIterator& operator++() { m_idval++; m_idrow++; - return *this; + return *this; } inline Scalar value() const { return m_matrix.valuePtr()[m_idval]; } @@ -209,17 +210,19 @@ class MappedSuperNodalMatrix<Scalar,Index>::InnerIterator inline operator bool() const { - return ( (m_idval < m_endidval) && (m_idval >= m_startidval) ); + return ( (m_idval < m_endidval) && (m_idval >= m_startidval) + && (m_idrow < m_endidrow) ); } protected: const MappedSuperNodalMatrix& m_matrix; // Supernodal lower triangular matrix - const Index m_outer; // Current column - const Index m_supno; // Current SuperNode number - Index m_idval; //Index to browse the values in the current column - const Index m_startidval; // Start of the column value - const Index m_endidval; // End of the column value - Index m_idrow; //Index to browse the row indices + const Index m_outer; // Current column + const Index m_supno; // Current SuperNode number + Index m_idval; // Index to browse the values in the current column + const Index m_startidval; // Start of the column value + const Index m_endidval; // End of the column value + Index m_idrow; // Index to browse the row indices + Index m_endidrow; // End index of row indices of the current column }; /** @@ -232,32 +235,32 @@ void MappedSuperNodalMatrix<Scalar,Index>::solveInPlace( MatrixBase<Dest>&X) con { Index n = X.rows(); Index nrhs = X.cols(); - const Scalar * Lval = valuePtr(); // Nonzero values - Matrix<Scalar,Dynamic,Dynamic> work(n, nrhs); // working vector + const Scalar * Lval = valuePtr(); // Nonzero values + Matrix<Scalar,Dynamic,Dynamic> work(n, nrhs); // working vector work.setZero(); for (Index k = 0; k <= nsuper(); k ++) { - Index fsupc = supToCol()[k]; // First column of the current supernode - Index istart = rowIndexPtr()[fsupc]; // Pointer index to the subscript of the current column + Index fsupc = supToCol()[k]; // First column of the current supernode + Index istart = rowIndexPtr()[fsupc]; // Pointer index to the subscript of the current column Index nsupr = rowIndexPtr()[fsupc+1] - istart; // Number of rows in the current supernode - Index nsupc = supToCol()[k+1] - fsupc; // Number of columns in the current supernode - Index nrow = nsupr - nsupc; // Number of rows in the non-diagonal part of the supernode - Index irow; //Current index row + Index nsupc = supToCol()[k+1] - fsupc; // Number of columns in the current supernode + Index nrow = nsupr - nsupc; // Number of rows in the non-diagonal part of the supernode + Index irow; //Current index row if (nsupc == 1 ) { for (Index j = 0; j < nrhs; j++) { - InnerIterator it(*this, fsupc); + InnerIterator it(*this, fsupc); ++it; // Skip the diagonal element for (; it; ++it) { irow = it.row(); - X(irow, j) -= X(fsupc, j) * it.value(); + X(irow, j) -= X(fsupc, j) * it.value(); } } } - else + else { // The supernode has more than one column Index luptr = colIndexPtr()[fsupc]; @@ -291,4 +294,5 @@ void MappedSuperNodalMatrix<Scalar,Index>::solveInPlace( MatrixBase<Dest>&X) con } // end namespace internal } // end namespace Eigen + #endif // EIGEN_SPARSELU_MATRIX_H |