// This file is part of Eigen, a lightweight C++ template library // for linear algebra. // // Copyright (C) 2008 Benoit Jacob <[email protected]> // Copyright (C) 2008-2009 Gael Guennebaud <[email protected]> // // 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 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. #ifndef EIGEN_TRIANGULARMATRIX_H #define EIGEN_TRIANGULARMATRIX_H // IWYU pragma: private #include "./InternalHeaderCheck.h" namespace Eigen { namespace internal { template <int Side, typename TriangularType, typename Rhs> struct triangular_solve_retval; } /** \class TriangularBase * \ingroup Core_Module * * \brief Base class for triangular part in a matrix */ template <typename Derived> class TriangularBase : public EigenBase<Derived> { … }; /** \class TriangularView * \ingroup Core_Module * * \brief Expression of a triangular part in a matrix * * \tparam MatrixType the type of the object in which we are taking the triangular part * \tparam Mode the kind of triangular matrix expression to construct. Can be #Upper, * #Lower, #UnitUpper, #UnitLower, #StrictlyUpper, or #StrictlyLower. * This is in fact a bit field; it must have either #Upper or #Lower, * and additionally it may have #UnitDiag or #ZeroDiag or neither. * * This class represents a triangular part of a matrix, not necessarily square. Strictly speaking, for rectangular * matrices one should speak of "trapezoid" parts. This class is the return type * of MatrixBase::triangularView() and SparseMatrixBase::triangularView(), and most of the time this is the only way it * is used. * * \sa MatrixBase::triangularView() */ namespace internal { traits<TriangularView<MatrixType, Mode_>>; } // namespace internal template <typename MatrixType_, unsigned int Mode_, typename StorageKind> class TriangularViewImpl; template <typename MatrixType_, unsigned int Mode_> class TriangularView : public TriangularViewImpl<MatrixType_, Mode_, typename internal::traits<MatrixType_>::StorageKind> { … }; /** \ingroup Core_Module * * \brief Base class for a triangular part in a \b dense matrix * * This class is an abstract base class of class TriangularView, and objects of type TriangularViewImpl cannot be * instantiated. It extends class TriangularView with additional methods which available for dense expressions only. * * \sa class TriangularView, MatrixBase::triangularView() */ TriangularViewImpl<MatrixType_, Mode_, Dense>; /*************************************************************************** * Implementation of triangular evaluation/assignment ***************************************************************************/ #ifndef EIGEN_PARSED_BY_DOXYGEN // FIXME should we keep that possibility template <typename MatrixType, unsigned int Mode> template <typename OtherDerived> EIGEN_DEVICE_FUNC inline TriangularView<MatrixType, Mode>& TriangularViewImpl<MatrixType, Mode, Dense>::operator=( const MatrixBase<OtherDerived>& other) { … } // FIXME should we keep that possibility template <typename MatrixType, unsigned int Mode> template <typename OtherDerived> EIGEN_DEVICE_FUNC void TriangularViewImpl<MatrixType, Mode, Dense>::lazyAssign(const MatrixBase<OtherDerived>& other) { … } template <typename MatrixType, unsigned int Mode> template <typename OtherDerived> EIGEN_DEVICE_FUNC inline TriangularView<MatrixType, Mode>& TriangularViewImpl<MatrixType, Mode, Dense>::operator=( const TriangularBase<OtherDerived>& other) { … } template <typename MatrixType, unsigned int Mode> template <typename OtherDerived> EIGEN_DEVICE_FUNC void TriangularViewImpl<MatrixType, Mode, Dense>::lazyAssign( const TriangularBase<OtherDerived>& other) { … } #endif /*************************************************************************** * Implementation of TriangularBase methods ***************************************************************************/ /** Assigns a triangular or selfadjoint matrix to a dense matrix. * If the matrix is triangular, the opposite part is set to zero. */ template <typename Derived> template <typename DenseDerived> EIGEN_DEVICE_FUNC void TriangularBase<Derived>::evalTo(MatrixBase<DenseDerived>& other) const { … } /*************************************************************************** * Implementation of TriangularView methods ***************************************************************************/ /*************************************************************************** * Implementation of MatrixBase methods ***************************************************************************/ /** * \returns an expression of a triangular view extracted from the current matrix * * The parameter \a Mode can have the following values: \c #Upper, \c #StrictlyUpper, \c #UnitUpper, * \c #Lower, \c #StrictlyLower, \c #UnitLower. * * Example: \include MatrixBase_triangularView.cpp * Output: \verbinclude MatrixBase_triangularView.out * * \sa class TriangularView */ template <typename Derived> template <unsigned int Mode> EIGEN_DEVICE_FUNC typename MatrixBase<Derived>::template TriangularViewReturnType<Mode>::Type MatrixBase<Derived>::triangularView() { … } /** This is the const version of MatrixBase::triangularView() */ template <typename Derived> template <unsigned int Mode> EIGEN_DEVICE_FUNC typename MatrixBase<Derived>::template ConstTriangularViewReturnType<Mode>::Type MatrixBase<Derived>::triangularView() const { … } /** \returns true if *this is approximately equal to an upper triangular matrix, * within the precision given by \a prec. * * \sa isLowerTriangular() */ template <typename Derived> bool MatrixBase<Derived>::isUpperTriangular(const RealScalar& prec) const { … } /** \returns true if *this is approximately equal to a lower triangular matrix, * within the precision given by \a prec. * * \sa isUpperTriangular() */ template <typename Derived> bool MatrixBase<Derived>::isLowerTriangular(const RealScalar& prec) const { … } /*************************************************************************** **************************************************************************** * Evaluators and Assignment of triangular expressions *************************************************************************** ***************************************************************************/ namespace internal { // TODO currently a triangular expression has the form TriangularView<.,.> // in the future triangular-ness should be defined by the expression traits // such that Transpose<TriangularView<.,.> > is valid. (currently TriangularBase::transpose() is overloaded to make // it work) evaluator_traits<TriangularView<MatrixType, Mode>>; unary_evaluator<TriangularView<MatrixType, Mode>, IndexBased>; // Additional assignment kinds: struct Triangular2Triangular { … }; struct Triangular2Dense { … }; struct Dense2Triangular { … }; template <typename Kernel, unsigned int Mode, int UnrollCount, bool ClearOpposite> struct triangular_assignment_loop; /** \internal Specialization of the dense assignment kernel for triangular matrices. * The main difference is that the triangular, diagonal, and opposite parts are processed through three different * functions. \tparam UpLo must be either Lower or Upper \tparam Mode must be either 0, UnitDiag, ZeroDiag, or * SelfAdjoint */ template <int UpLo, int Mode, int SetOpposite, typename DstEvaluatorTypeT, typename SrcEvaluatorTypeT, typename Functor, int Version = Specialized> class triangular_dense_assignment_kernel : public generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version> { … }; template <int Mode, bool SetOpposite, typename DstXprType, typename SrcXprType, typename Functor> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void call_triangular_assignment_loop(DstXprType& dst, const SrcXprType& src, const Functor& func) { … } template <int Mode, bool SetOpposite, typename DstXprType, typename SrcXprType> EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void call_triangular_assignment_loop(DstXprType& dst, const SrcXprType& src) { … } template <> struct AssignmentKind<TriangularShape, TriangularShape> { … }; template <> struct AssignmentKind<DenseShape, TriangularShape> { … }; template <> struct AssignmentKind<TriangularShape, DenseShape> { … }; Assignment<DstXprType, SrcXprType, Functor, Triangular2Triangular>; Assignment<DstXprType, SrcXprType, Functor, Triangular2Dense>; Assignment<DstXprType, SrcXprType, Functor, Dense2Triangular>; template <typename Kernel, unsigned int Mode, int UnrollCount, bool SetOpposite> struct triangular_assignment_loop { … }; // prevent buggy user code from causing an infinite recursion triangular_assignment_loop<Kernel, Mode, 0, SetOpposite>; // TODO: experiment with a recursive assignment procedure splitting the current // triangular part into one rectangular and two triangular parts. triangular_assignment_loop<Kernel, Mode, Dynamic, SetOpposite>; } // end namespace internal /** Assigns a triangular or selfadjoint matrix to a dense matrix. * If the matrix is triangular, the opposite part is set to zero. */ template <typename Derived> template <typename DenseDerived> EIGEN_DEVICE_FUNC void TriangularBase<Derived>::evalToLazy(MatrixBase<DenseDerived>& other) const { … } namespace internal { // Triangular = Product Assignment<DstXprType, Product<Lhs, Rhs, DefaultProduct>, internal::assign_op<Scalar, typename Product<Lhs, Rhs, DefaultProduct>::Scalar>, Dense2Triangular>; // Triangular += Product Assignment<DstXprType, Product<Lhs, Rhs, DefaultProduct>, internal::add_assign_op<Scalar, typename Product<Lhs, Rhs, DefaultProduct>::Scalar>, Dense2Triangular>; // Triangular -= Product Assignment<DstXprType, Product<Lhs, Rhs, DefaultProduct>, internal::sub_assign_op<Scalar, typename Product<Lhs, Rhs, DefaultProduct>::Scalar>, Dense2Triangular>; } // end namespace internal } // end namespace Eigen #endif // EIGEN_TRIANGULARMATRIX_H