Diagonal.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) 2007-2009 Benoit Jacob <jacob.benoit.1@gmail.com>
00005 // Copyright (C) 2009-2010 Gael Guennebaud <gael.guennebaud@inria.fr>
00006 //
00007 // Eigen is free software; you can redistribute it and/or
00008 // modify it under the terms of the GNU Lesser General Public
00009 // License as published by the Free Software Foundation; either
00010 // version 3 of the License, or (at your option) any later version.
00011 //
00012 // Alternatively, you can redistribute it and/or
00013 // modify it under the terms of the GNU General Public License as
00014 // published by the Free Software Foundation; either version 2 of
00015 // the License, or (at your option) any later version.
00016 //
00017 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
00018 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
00019 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
00020 // GNU General Public License for more details.
00021 //
00022 // You should have received a copy of the GNU Lesser General Public
00023 // License and a copy of the GNU General Public License along with
00024 // Eigen. If not, see <http://www.gnu.org/licenses/>.
00025 
00026 #ifndef EIGEN_DIAGONAL_H
00027 #define EIGEN_DIAGONAL_H
00028 
00029 namespace Eigen { 
00030 
00050 namespace internal {
00051 template<typename MatrixType, int DiagIndex>
00052 struct traits<Diagonal<MatrixType,DiagIndex> >
00053  : traits<MatrixType>
00054 {
00055   typedef typename nested<MatrixType>::type MatrixTypeNested;
00056   typedef typename remove_reference<MatrixTypeNested>::type _MatrixTypeNested;
00057   typedef typename MatrixType::StorageKind StorageKind;
00058   enum {
00059     AbsDiagIndex = DiagIndex<0 ? -DiagIndex : DiagIndex, // only used if DiagIndex != Dynamic
00060     // FIXME these computations are broken in the case where the matrix is rectangular and DiagIndex!=0
00061     RowsAtCompileTime = (int(DiagIndex) == Dynamic || int(MatrixType::SizeAtCompileTime) == Dynamic) ? Dynamic
00062                       : (EIGEN_SIZE_MIN_PREFER_DYNAMIC(MatrixType::RowsAtCompileTime,
00063                                         MatrixType::ColsAtCompileTime) - AbsDiagIndex),
00064     ColsAtCompileTime = 1,
00065     MaxRowsAtCompileTime = int(MatrixType::MaxSizeAtCompileTime) == Dynamic ? Dynamic
00066                          : DiagIndex == Dynamic ? EIGEN_SIZE_MIN_PREFER_FIXED(MatrixType::MaxRowsAtCompileTime,
00067                                                                     MatrixType::MaxColsAtCompileTime)
00068                          : (EIGEN_SIZE_MIN_PREFER_FIXED(MatrixType::MaxRowsAtCompileTime, MatrixType::MaxColsAtCompileTime) - AbsDiagIndex),
00069     MaxColsAtCompileTime = 1,
00070     MaskLvalueBit = is_lvalue<MatrixType>::value ? LvalueBit : 0,
00071     Flags = (unsigned int)_MatrixTypeNested::Flags & (HereditaryBits | LinearAccessBit | MaskLvalueBit | DirectAccessBit) & ~RowMajorBit,
00072     CoeffReadCost = _MatrixTypeNested::CoeffReadCost,
00073     MatrixTypeOuterStride = outer_stride_at_compile_time<MatrixType>::ret,
00074     InnerStrideAtCompileTime = MatrixTypeOuterStride == Dynamic ? Dynamic : MatrixTypeOuterStride+1,
00075     OuterStrideAtCompileTime = 0
00076   };
00077 };
00078 }
00079 
00080 template<typename MatrixType, int DiagIndex> class Diagonal
00081    : public internal::dense_xpr_base< Diagonal<MatrixType,DiagIndex> >::type
00082 {
00083   public:
00084 
00085     typedef typename internal::dense_xpr_base<Diagonal>::type Base;
00086     EIGEN_DENSE_PUBLIC_INTERFACE(Diagonal)
00087 
00088     inline Diagonal(MatrixType& matrix, Index index = DiagIndex) : m_matrix(matrix), m_index(index) {}
00089 
00090     EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Diagonal)
00091 
00092     inline Index rows() const
00093     { return m_index.value()<0 ? (std::min)(m_matrix.cols(),m_matrix.rows()+m_index.value()) : (std::min)(m_matrix.rows(),m_matrix.cols()-m_index.value()); }
00094 
00095     inline Index cols() const { return 1; }
00096 
00097     inline Index innerStride() const
00098     {
00099       return m_matrix.outerStride() + 1;
00100     }
00101 
00102     inline Index outerStride() const
00103     {
00104       return 0;
00105     }
00106 
00107     typedef typename internal::conditional<
00108                        internal::is_lvalue<MatrixType>::value,
00109                        Scalar,
00110                        const Scalar
00111                      >::type ScalarWithConstIfNotLvalue;
00112 
00113     inline ScalarWithConstIfNotLvalue* data() { return &(m_matrix.const_cast_derived().coeffRef(rowOffset(), colOffset())); }
00114     inline const Scalar* data() const { return &(m_matrix.const_cast_derived().coeffRef(rowOffset(), colOffset())); }
00115 
00116     inline Scalar& coeffRef(Index row, Index)
00117     {
00118       EIGEN_STATIC_ASSERT_LVALUE(MatrixType)
00119       return m_matrix.const_cast_derived().coeffRef(row+rowOffset(), row+colOffset());
00120     }
00121 
00122     inline const Scalar& coeffRef(Index row, Index) const
00123     {
00124       return m_matrix.const_cast_derived().coeffRef(row+rowOffset(), row+colOffset());
00125     }
00126 
00127     inline CoeffReturnType coeff(Index row, Index) const
00128     {
00129       return m_matrix.coeff(row+rowOffset(), row+colOffset());
00130     }
00131 
00132     inline Scalar& coeffRef(Index index)
00133     {
00134       EIGEN_STATIC_ASSERT_LVALUE(MatrixType)
00135       return m_matrix.const_cast_derived().coeffRef(index+rowOffset(), index+colOffset());
00136     }
00137 
00138     inline const Scalar& coeffRef(Index index) const
00139     {
00140       return m_matrix.const_cast_derived().coeffRef(index+rowOffset(), index+colOffset());
00141     }
00142 
00143     inline CoeffReturnType coeff(Index index) const
00144     {
00145       return m_matrix.coeff(index+rowOffset(), index+colOffset());
00146     }
00147 
00148     const typename internal::remove_all<typename MatrixType::Nested>::type& 
00149     nestedExpression() const 
00150     {
00151       return m_matrix;
00152     }
00153 
00154     int index() const
00155     {
00156       return m_index.value();
00157     }
00158 
00159   protected:
00160     typename MatrixType::Nested m_matrix;
00161     const internal::variable_if_dynamic<Index, DiagIndex> m_index;
00162 
00163   private:
00164     // some compilers may fail to optimize std::max etc in case of compile-time constants...
00165     EIGEN_STRONG_INLINE Index absDiagIndex() const { return m_index.value()>0 ? m_index.value() : -m_index.value(); }
00166     EIGEN_STRONG_INLINE Index rowOffset() const { return m_index.value()>0 ? 0 : -m_index.value(); }
00167     EIGEN_STRONG_INLINE Index colOffset() const { return m_index.value()>0 ? m_index.value() : 0; }
00168     // triger a compile time error is someone try to call packet
00169     template<int LoadMode> typename MatrixType::PacketReturnType packet(Index) const;
00170     template<int LoadMode> typename MatrixType::PacketReturnType packet(Index,Index) const;
00171 };
00172 
00181 template<typename Derived>
00182 inline typename MatrixBase<Derived>::DiagonalReturnType
00183 MatrixBase<Derived>::diagonal()
00184 {
00185   return derived();
00186 }
00187 
00189 template<typename Derived>
00190 inline const typename MatrixBase<Derived>::ConstDiagonalReturnType
00191 MatrixBase<Derived>::diagonal() const
00192 {
00193   return ConstDiagonalReturnType(derived());
00194 }
00195 
00207 template<typename Derived>
00208 inline typename MatrixBase<Derived>::template DiagonalIndexReturnType<Dynamic>::Type
00209 MatrixBase<Derived>::diagonal(Index index)
00210 {
00211   return typename DiagonalIndexReturnType<Dynamic>::Type(derived(), index);
00212 }
00213 
00215 template<typename Derived>
00216 inline typename MatrixBase<Derived>::template ConstDiagonalIndexReturnType<Dynamic>::Type
00217 MatrixBase<Derived>::diagonal(Index index) const
00218 {
00219   return typename ConstDiagonalIndexReturnType<Dynamic>::Type(derived(), index);
00220 }
00221 
00233 template<typename Derived>
00234 template<int Index>
00235 inline typename MatrixBase<Derived>::template DiagonalIndexReturnType<Index>::Type
00236 MatrixBase<Derived>::diagonal()
00237 {
00238   return derived();
00239 }
00240 
00242 template<typename Derived>
00243 template<int Index>
00244 inline typename MatrixBase<Derived>::template ConstDiagonalIndexReturnType<Index>::Type
00245 MatrixBase<Derived>::diagonal() const
00246 {
00247   return derived();
00248 }
00249 
00250 } // end namespace Eigen
00251 
00252 #endif // EIGEN_DIAGONAL_H