diff --git a/Eigen/src/Sparse/SparseMatrix.h b/Eigen/src/Sparse/SparseMatrix.h index 2504d9fd8..b34ed8fd2 100644 --- a/Eigen/src/Sparse/SparseMatrix.h +++ b/Eigen/src/Sparse/SparseMatrix.h @@ -84,16 +84,18 @@ class SparseMatrix : public SparseMatrixBase > const int outer = RowMajor ? row : col; const int inner = RowMajor ? col : row; - int id = m_outerIndex[outer]; + int start = m_outerIndex[outer]; int end = m_outerIndex[outer+1]; - // optimization: let's first check if it is the last coefficient - // (very common in high level algorithms) - if (end>0 && inner==m_data.index(end-1)) - return m_data.value(end-1); - else if (id==end) + if (start==end) return Scalar(0); - const int* r = std::lower_bound(&m_data.index(id),&m_data.index(end),inner); - return (*r==inner) ? m_data.value(r-&m_data.index(0)) : Scalar(0); + else if (end>0 && inner==m_data.index(end-1)) + return m_data.value(end-1); + // ^^ optimization: let's first check if it is the last coefficient + // (very common in high level algorithms) + + const int* r = std::lower_bound(&m_data.index(start),&m_data.index(end),inner); + const int id = r-&m_data.index(0); + return ((*r==inner) && (id > const int outer = RowMajor ? row : col; const int inner = RowMajor ? col : row; - int id = m_outerIndex[outer]; + int start = m_outerIndex[outer]; int end = m_outerIndex[outer+1]; - ei_assert(end>=id && "you probably called coeffRef on a non finalized matrix"); - ei_assert(end>id && "coeffRef cannot be called on a zero coefficient"); - int* r = std::lower_bound(&m_data.index(id),&m_data.index(end),inner); - ei_assert(*r==inner && "coeffRef cannot be called on a zero coefficient"); - return m_data.value(r-&m_data.index(0)); + ei_assert(end>=start && "you probably called coeffRef on a non finalized matrix"); + ei_assert(end>start && "coeffRef cannot be called on a zero coefficient"); + int* r = std::lower_bound(&m_data.index(start),&m_data.index(end),inner); + const int id = r-&m_data.index(0); + ei_assert((*r==inner) && (id { for ( ; col +template void sparse() +{ + int rows = 8, cols = 8; + double density = std::max(8./(rows*cols), 0.01); + typedef Matrix DenseMatrix; + Scalar eps = 1e-6; + + SparseMatrix m(rows, cols); + DenseMatrix refMat = DenseMatrix::Zero(rows, cols); + + std::vector zeroCoords; + std::vector nonzeroCoords; + m.startFill(rows*cols*density); + for(int j=0; j(0,1) < density) ? ei_random() : 0; + if (v!=0) + { + m.fill(i,j) = v; + nonzeroCoords.push_back(Vector2i(i,j)); + } + else + { + zeroCoords.push_back(Vector2i(i,j)); + } + refMat(i,j) = v; + } + } + m.endFill(); + + VERIFY(zeroCoords.size()>0 && "re-run the test"); + VERIFY(nonzeroCoords.size()>0 && "re-run the test"); + + // test coeff and coeffRef + for (int i=0; i, FullyCoherentAccessPattern> w(m); +// for (int i=0; icoeffRef(nonzeroCoords[i].x(),nonzeroCoords[i].y()) = refMat.coeff(nonzeroCoords[i].x(),nonzeroCoords[i].y()); +// } +// } +// VERIFY_IS_APPROX(m, refMat); + + // random setter + { + m.setZero(); + VERIFY_IS_NOT_APPROX(m, refMat); + SparseSetter, RandomAccessPattern> w(m); + std::vector remaining = nonzeroCoords; + while(!remaining.empty()) + { + int i = ei_random(0,remaining.size()-1); + w->coeffRef(remaining[i].x(),remaining[i].y()) = refMat.coeff(remaining[i].x(),remaining[i].y()); + remaining[i] = remaining.back(); + remaining.pop_back(); + } + } + VERIFY_IS_APPROX(m, refMat); +} + void test_sparse() { - int rows = 4, cols = 4; - SparseMatrix m(rows, cols); - - m.startFill(rows); - m.fill(0, 2) = 2; - m.fill(1, 2) = 1; - m.fill(0, 3) = 5; - m.endFill(); - - m.coeffRef(0, 2) = 3; - VERIFY_RAISES_ASSERT( m.coeffRef(0, 0) = 5 ); - VERIFY_IS_MUCH_SMALLER_THAN( m.coeff(0, 0), 0.000001 ); - VERIFY_IS_MUCH_SMALLER_THAN( m.coeff(0, 1), 0.000001 ); - VERIFY_IS_MUCH_SMALLER_THAN( m.coeff(2, 1), 0.000001 ); - VERIFY_IS_APPROX( m.coeff(0, 2), 3.0 ); - VERIFY_IS_APPROX( m.coeff(1, 2), 1.0 ); - VERIFY_IS_APPROX( m.coeff(0, 3), 5.0 ); - - Matrix4d dm; - double r; - m.startFill(rows*cols); - for(int i=0; i(); }