Add support with unit test for off-diagonal sparse matrix views

This commit is contained in:
Desire NUENTSA 2013-06-11 14:42:29 +02:00
parent b3fff170a0
commit 4cd8245c96
2 changed files with 37 additions and 7 deletions

View File

@ -2,6 +2,7 @@
// for linear algebra. // for linear algebra.
// //
// Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr> // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr>
// Copyright (C) 2012 Désiré Nuentsa-Wakam <desire.nuentsa_wakam@inria.fr>
// //
// This Source Code Form is subject to the terms of the Mozilla // 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 // Public License v. 2.0. If a copy of the MPL was not distributed
@ -27,6 +28,7 @@ template<typename MatrixType, int Mode> class SparseTriangularView
enum { SkipFirst = ((Mode&Lower) && !(MatrixType::Flags&RowMajorBit)) enum { SkipFirst = ((Mode&Lower) && !(MatrixType::Flags&RowMajorBit))
|| ((Mode&Upper) && (MatrixType::Flags&RowMajorBit)), || ((Mode&Upper) && (MatrixType::Flags&RowMajorBit)),
SkipLast = !SkipFirst, SkipLast = !SkipFirst,
SkipDiag = (Mode&ZeroDiag) ? 1 : 0,
HasUnitDiag = (Mode&UnitDiag) ? 1 : 0 HasUnitDiag = (Mode&UnitDiag) ? 1 : 0
}; };
@ -71,7 +73,7 @@ class SparseTriangularView<MatrixType,Mode>::InnerIterator : public MatrixTypeNe
{ {
if(SkipFirst) if(SkipFirst)
{ {
while((*this) && (HasUnitDiag ? this->index()<=outer : this->index()<outer)) while((*this) && ((HasUnitDiag||SkipDiag) ? this->index()<=outer : this->index()<outer))
Base::operator++(); Base::operator++();
if(HasUnitDiag) if(HasUnitDiag)
m_returnOne = true; m_returnOne = true;
@ -101,8 +103,8 @@ class SparseTriangularView<MatrixType,Mode>::InnerIterator : public MatrixTypeNe
return *this; return *this;
} }
inline Index row() const { return Base::row(); } inline Index row() const { return (MatrixType::Flags&RowMajorBit ? Base::outer() : this->index()); }
inline Index col() const { return Base::col(); } inline Index col() const { return (MatrixType::Flags&RowMajorBit ? this->index() : Base::outer()); }
inline Index index() const inline Index index() const
{ {
if(HasUnitDiag && m_returnOne) return Base::outer(); if(HasUnitDiag && m_returnOne) return Base::outer();
@ -118,7 +120,12 @@ class SparseTriangularView<MatrixType,Mode>::InnerIterator : public MatrixTypeNe
{ {
if(HasUnitDiag && m_returnOne) if(HasUnitDiag && m_returnOne)
return true; return true;
return (SkipFirst ? Base::operator bool() : (Base::operator bool() && this->index() <= this->outer())); if(SkipFirst) return Base::operator bool();
else
{
if (SkipDiag) return (Base::operator bool() && this->index() < this->outer());
else return (Base::operator bool() && this->index() <= this->outer());
}
} }
protected: protected:
bool m_returnOne; bool m_returnOne;
@ -134,9 +141,10 @@ class SparseTriangularView<MatrixType,Mode>::ReverseInnerIterator : public Matri
: Base(view.nestedExpression(), outer) : Base(view.nestedExpression(), outer)
{ {
eigen_assert((!HasUnitDiag) && "ReverseInnerIterator does not support yet triangular views with a unit diagonal"); eigen_assert((!HasUnitDiag) && "ReverseInnerIterator does not support yet triangular views with a unit diagonal");
if(SkipLast) if(SkipLast) {
while((*this) && this->index()>outer) while((*this) && (SkipDiag ? this->index()>=outer : this->index()>outer))
--(*this); --(*this);
}
} }
EIGEN_STRONG_INLINE ReverseInnerIterator& operator--() EIGEN_STRONG_INLINE ReverseInnerIterator& operator--()
@ -147,7 +155,12 @@ class SparseTriangularView<MatrixType,Mode>::ReverseInnerIterator : public Matri
EIGEN_STRONG_INLINE operator bool() const EIGEN_STRONG_INLINE operator bool() const
{ {
return SkipLast ? Base::operator bool() : (Base::operator bool() && this->index() >= this->outer()); if (SkipLast) return Base::operator bool() ;
else
{
if(SkipDiag) return (Base::operator bool() && this->index() > this->outer());
else return (Base::operator bool() && this->index() >= this->outer());
}
} }
}; };

View File

@ -3,6 +3,7 @@
// //
// Copyright (C) 2008-2011 Gael Guennebaud <gael.guennebaud@inria.fr> // Copyright (C) 2008-2011 Gael Guennebaud <gael.guennebaud@inria.fr>
// Copyright (C) 2008 Daniel Gomez Ferro <dgomezferro@gmail.com> // Copyright (C) 2008 Daniel Gomez Ferro <dgomezferro@gmail.com>
// Copyright (C) 2013 Désiré Nuentsa-Wakam <desire.nuentsa_wakam@inria.fr>
// //
// This Source Code Form is subject to the terms of the Mozilla // 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 // Public License v. 2.0. If a copy of the MPL was not distributed
@ -391,6 +392,14 @@ template<typename SparseMatrixType> void sparse_basic(const SparseMatrixType& re
refMat3 = refMat2.template triangularView<UnitLower>(); refMat3 = refMat2.template triangularView<UnitLower>();
m3 = m2.template triangularView<UnitLower>(); m3 = m2.template triangularView<UnitLower>();
VERIFY_IS_APPROX(m3, refMat3); VERIFY_IS_APPROX(m3, refMat3);
refMat3 = refMat2.template triangularView<StrictlyUpper>();
m3 = m2.template triangularView<StrictlyUpper>();
VERIFY_IS_APPROX(m3, refMat3);
refMat3 = refMat2.template triangularView<StrictlyLower>();
m3 = m2.template triangularView<StrictlyLower>();
VERIFY_IS_APPROX(m3, refMat3);
} }
// test selfadjointView // test selfadjointView
@ -454,6 +463,14 @@ template<typename SparseMatrixType> void sparse_basic(const SparseMatrixType& re
} }
} }
// test Identity matrix
{
DenseMatrix refMat1 = DenseMatrix::Identity(rows, rows);
SparseMatrixType m1(rows, rows);
m1.setIdentity();
VERIFY_IS_APPROX(m1, refMat1);
}
} }
void test_sparse_basic() void test_sparse_basic()