SparseTriangularView.h
Go to the documentation of this file.
00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr>
00005 //
00006 // Eigen is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 3 of the License, or (at your option) any later version.
00010 //
00011 // Alternatively, you can redistribute it and/or
00012 // modify it under the terms of the GNU General Public License as
00013 // published by the Free Software Foundation; either version 2 of
00014 // the License, or (at your option) any later version.
00015 //
00016 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
00017 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
00018 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
00019 // GNU General Public License for more details.
00020 //
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License and a copy of the GNU General Public License along with
00023 // Eigen. If not, see <http://www.gnu.org/licenses/>.
00024 
00025 #ifndef EIGEN_SPARSE_TRIANGULARVIEW_H
00026 #define EIGEN_SPARSE_TRIANGULARVIEW_H
00027 
00028 namespace Eigen { 
00029 
00030 namespace internal {
00031   
00032 template<typename MatrixType, int Mode>
00033 struct traits<SparseTriangularView<MatrixType,Mode> >
00034 : public traits<MatrixType>
00035 {};
00036 
00037 } // namespace internal
00038 
00039 template<typename MatrixType, int Mode> class SparseTriangularView
00040   : public SparseMatrixBase<SparseTriangularView<MatrixType,Mode> >
00041 {
00042     enum { SkipFirst = ((Mode&Lower) && !(MatrixType::Flags&RowMajorBit))
00043                     || ((Mode&Upper) &&  (MatrixType::Flags&RowMajorBit)),
00044            SkipLast = !SkipFirst,
00045            HasUnitDiag = (Mode&UnitDiag) ? 1 : 0
00046     };
00047 
00048   public:
00049     
00050     EIGEN_SPARSE_PUBLIC_INTERFACE(SparseTriangularView)
00051 
00052     class InnerIterator;
00053     class ReverseInnerIterator;
00054 
00055     inline Index rows() const { return m_matrix.rows(); }
00056     inline Index cols() const { return m_matrix.cols(); }
00057 
00058     typedef typename MatrixType::Nested MatrixTypeNested;
00059     typedef typename internal::remove_reference<MatrixTypeNested>::type MatrixTypeNestedNonRef;
00060     typedef typename internal::remove_all<MatrixTypeNested>::type MatrixTypeNestedCleaned;
00061 
00062     inline SparseTriangularView(const MatrixType& matrix) : m_matrix(matrix) {}
00063 
00065     inline const MatrixTypeNestedCleaned& nestedExpression() const { return m_matrix; }
00066 
00067     template<typename OtherDerived>
00068     typename internal::plain_matrix_type_column_major<OtherDerived>::type
00069     solve(const MatrixBase<OtherDerived>& other) const;
00070 
00071     template<typename OtherDerived> void solveInPlace(MatrixBase<OtherDerived>& other) const;
00072     template<typename OtherDerived> void solveInPlace(SparseMatrixBase<OtherDerived>& other) const;
00073 
00074   protected:
00075     MatrixTypeNested m_matrix;
00076 };
00077 
00078 template<typename MatrixType, int Mode>
00079 class SparseTriangularView<MatrixType,Mode>::InnerIterator : public MatrixTypeNestedCleaned::InnerIterator
00080 {
00081     typedef typename MatrixTypeNestedCleaned::InnerIterator Base;
00082   public:
00083 
00084     EIGEN_STRONG_INLINE InnerIterator(const SparseTriangularView& view, Index outer)
00085       : Base(view.nestedExpression(), outer), m_returnOne(false)
00086     {
00087       if(SkipFirst)
00088       {
00089         while((*this) && (HasUnitDiag ? this->index()<=outer : this->index()<outer))
00090           Base::operator++();
00091         if(HasUnitDiag)
00092           m_returnOne = true;
00093       }
00094       else if(HasUnitDiag && ((!Base::operator bool()) || Base::index()>=Base::outer()))
00095       {
00096         if((!SkipFirst) && Base::operator bool())
00097           Base::operator++();
00098         m_returnOne = true;
00099       }
00100     }
00101 
00102     EIGEN_STRONG_INLINE InnerIterator& operator++()
00103     {
00104       if(HasUnitDiag && m_returnOne)
00105         m_returnOne = false;
00106       else
00107       {
00108         Base::operator++();
00109         if(HasUnitDiag && (!SkipFirst) && ((!Base::operator bool()) || Base::index()>=Base::outer()))
00110         {
00111           if((!SkipFirst) && Base::operator bool())
00112             Base::operator++();
00113           m_returnOne = true;
00114         }
00115       }
00116       return *this;
00117     }
00118 
00119     inline Index row() const { return Base::row(); }
00120     inline Index col() const { return Base::col(); }
00121     inline Index index() const
00122     {
00123       if(HasUnitDiag && m_returnOne)  return Base::outer();
00124       else                            return Base::index();
00125     }
00126     inline Scalar value() const
00127     {
00128       if(HasUnitDiag && m_returnOne)  return Scalar(1);
00129       else                            return Base::value();
00130     }
00131 
00132     EIGEN_STRONG_INLINE operator bool() const
00133     {
00134       if(HasUnitDiag && m_returnOne)
00135         return true;
00136       return (SkipFirst ? Base::operator bool() : (Base::operator bool() && this->index() <= this->outer()));
00137     }
00138   protected:
00139     bool m_returnOne;
00140 };
00141 
00142 template<typename MatrixType, int Mode>
00143 class SparseTriangularView<MatrixType,Mode>::ReverseInnerIterator : public MatrixTypeNestedCleaned::ReverseInnerIterator
00144 {
00145     typedef typename MatrixTypeNestedCleaned::ReverseInnerIterator Base;
00146   public:
00147 
00148     EIGEN_STRONG_INLINE ReverseInnerIterator(const SparseTriangularView& view, Index outer)
00149       : Base(view.nestedExpression(), outer)
00150     {
00151       eigen_assert((!HasUnitDiag) && "ReverseInnerIterator does not support yet triangular views with a unit diagonal");
00152       if(SkipLast)
00153         while((*this) && this->index()>outer)
00154           --(*this);
00155     }
00156 
00157     EIGEN_STRONG_INLINE InnerIterator& operator--()
00158     { Base::operator--(); return *this; }
00159 
00160     inline Index row() const { return Base::row(); }
00161     inline Index col() const { return Base::col(); }
00162 
00163     EIGEN_STRONG_INLINE operator bool() const
00164     {
00165       return SkipLast ? Base::operator bool() : (Base::operator bool() && this->index() >= this->outer());
00166     }
00167 };
00168 
00169 template<typename Derived>
00170 template<int Mode>
00171 inline const SparseTriangularView<Derived, Mode>
00172 SparseMatrixBase<Derived>::triangularView() const
00173 {
00174   return derived();
00175 }
00176 
00177 } // end namespace Eigen
00178 
00179 #endif // EIGEN_SPARSE_TRIANGULARVIEW_H