Remove Skyline.
diff --git a/unsupported/Eigen/CMakeLists.txt b/unsupported/Eigen/CMakeLists.txt index 3283f40..1517ba9 100644 --- a/unsupported/Eigen/CMakeLists.txt +++ b/unsupported/Eigen/CMakeLists.txt
@@ -16,7 +16,6 @@ NumericalDiff OpenGLSupport Polynomials - Skyline SparseExtra SpecialFunctions Splines
diff --git a/unsupported/Eigen/Skyline b/unsupported/Eigen/Skyline deleted file mode 100644 index ebd7c38..0000000 --- a/unsupported/Eigen/Skyline +++ /dev/null
@@ -1,40 +0,0 @@ -// This file is part of Eigen, a lightweight C++ template library -// for linear algebra. -// -// -// 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_SKYLINE_MODULE_H -#define EIGEN_SKYLINE_MODULE_H - -#include "../../Eigen/Core" - -#include "../../Eigen/src/Core/util/DisableStupidWarnings.h" - -#include <map> -#include <cstdlib> -#include <cstring> -#include <algorithm> - -/** - * \defgroup Skyline_Module Skyline module - * - * - * - * - */ - -// IWYU pragma: begin_exports -#include "src/Skyline/SkylineUtil.h" -#include "src/Skyline/SkylineMatrixBase.h" -#include "src/Skyline/SkylineStorage.h" -#include "src/Skyline/SkylineMatrix.h" -#include "src/Skyline/SkylineInplaceLU.h" -#include "src/Skyline/SkylineProduct.h" -// IWYU pragma: end_exports - -#include "../../Eigen/src/Core/util/ReenableStupidWarnings.h" - -#endif // EIGEN_SKYLINE_MODULE_H
diff --git a/unsupported/Eigen/src/Skyline/InternalHeaderCheck.h b/unsupported/Eigen/src/Skyline/InternalHeaderCheck.h deleted file mode 100644 index bddf0b0..0000000 --- a/unsupported/Eigen/src/Skyline/InternalHeaderCheck.h +++ /dev/null
@@ -1,3 +0,0 @@ -#ifndef EIGEN_SKYLINE_MODULE_H -#error "Please include unsupported/Eigen/Skyline instead of including headers inside the src directory directly." -#endif
diff --git a/unsupported/Eigen/src/Skyline/SkylineInplaceLU.h b/unsupported/Eigen/src/Skyline/SkylineInplaceLU.h deleted file mode 100644 index 50ac882..0000000 --- a/unsupported/Eigen/src/Skyline/SkylineInplaceLU.h +++ /dev/null
@@ -1,326 +0,0 @@ -// This file is part of Eigen, a lightweight C++ template library -// for linear algebra. -// -// Copyright (C) 2008 Guillaume Saupin <guillaume.saupin@cea.fr> -// -// 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_SKYLINEINPLACELU_H -#define EIGEN_SKYLINEINPLACELU_H - -// IWYU pragma: private -#include "./InternalHeaderCheck.h" - -namespace Eigen { - -/** \ingroup Skyline_Module - * - * \class SkylineInplaceLU - * - * \brief Inplace LU decomposition of a skyline matrix and associated features - * - * \param MatrixType the type of the matrix of which we are computing the LU factorization - * - */ -template <typename MatrixType> -class SkylineInplaceLU { - protected: - typedef typename MatrixType::Scalar Scalar; - typedef typename MatrixType::Index Index; - - typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar; - - public: - /** Creates a LU object and compute the respective factorization of \a matrix using - * flags \a flags. */ - SkylineInplaceLU(MatrixType& matrix, int flags = 0) - : /*m_matrix(matrix.rows(), matrix.cols()),*/ m_flags(flags), m_status(0), m_lu(matrix) { - m_precision = RealScalar(0.1) * Eigen::dummy_precision<RealScalar>(); - m_lu.IsRowMajor ? computeRowMajor() : compute(); - } - - /** Sets the relative threshold value used to prune zero coefficients during the decomposition. - * - * Setting a value greater than zero speeds up computation, and yields to an incomplete - * factorization with fewer non zero coefficients. Such approximate factors are especially - * useful to initialize an iterative solver. - * - * Note that the exact meaning of this parameter might depends on the actual - * backend. Moreover, not all backends support this feature. - * - * \sa precision() */ - void setPrecision(RealScalar v) { m_precision = v; } - - /** \returns the current precision. - * - * \sa setPrecision() */ - RealScalar precision() const { return m_precision; } - - /** Sets the flags. Possible values are: - * - CompleteFactorization - * - IncompleteFactorization - * - MemoryEfficient - * - one of the ordering methods - * - etc... - * - * \sa flags() */ - void setFlags(int f) { m_flags = f; } - - /** \returns the current flags */ - int flags() const { return m_flags; } - - void setOrderingMethod(int m) { m_flags = m; } - - int orderingMethod() const { return m_flags; } - - /** Computes/re-computes the LU factorization */ - void compute(); - void computeRowMajor(); - - /** \returns the lower triangular matrix L */ - // inline const MatrixType& matrixL() const { return m_matrixL; } - - /** \returns the upper triangular matrix U */ - // inline const MatrixType& matrixU() const { return m_matrixU; } - - template <typename BDerived, typename XDerived> - bool solve(const MatrixBase<BDerived>& b, MatrixBase<XDerived>* x, const int transposed = 0) const; - - /** \returns true if the factorization succeeded */ - inline bool succeeded(void) const { return m_succeeded; } - - protected: - RealScalar m_precision; - int m_flags; - mutable int m_status; - bool m_succeeded; - MatrixType& m_lu; -}; - -/** Computes / recomputes the in place LU decomposition of the SkylineInplaceLU. - * using the default algorithm. - */ -template <typename MatrixType> -// template<typename Scalar_> -void SkylineInplaceLU<MatrixType>::compute() { - const size_t rows = m_lu.rows(); - const size_t cols = m_lu.cols(); - - eigen_assert(rows == cols && "We do not (yet) support rectangular LU."); - eigen_assert(!m_lu.IsRowMajor && "LU decomposition does not work with rowMajor Storage"); - - for (Index row = 0; row < rows; row++) { - const double pivot = m_lu.coeffDiag(row); - - // Lower matrix Columns update - const Index& col = row; - for (typename MatrixType::InnerLowerIterator lIt(m_lu, col); lIt; ++lIt) { - lIt.valueRef() /= pivot; - } - - // Upper matrix update -> contiguous memory access - typename MatrixType::InnerLowerIterator lIt(m_lu, col); - for (Index rrow = row + 1; rrow < m_lu.rows(); rrow++) { - typename MatrixType::InnerUpperIterator uItPivot(m_lu, row); - typename MatrixType::InnerUpperIterator uIt(m_lu, rrow); - const double coef = lIt.value(); - - uItPivot += (rrow - row - 1); - - // update upper part -> contiguous memory access - for (++uItPivot; uIt && uItPivot;) { - uIt.valueRef() -= uItPivot.value() * coef; - - ++uIt; - ++uItPivot; - } - ++lIt; - } - - // Upper matrix update -> non contiguous memory access - typename MatrixType::InnerLowerIterator lIt3(m_lu, col); - for (Index rrow = row + 1; rrow < m_lu.rows(); rrow++) { - typename MatrixType::InnerUpperIterator uItPivot(m_lu, row); - const double coef = lIt3.value(); - - // update lower part -> non contiguous memory access - for (Index i = 0; i < rrow - row - 1; i++) { - m_lu.coeffRefLower(rrow, row + i + 1) -= uItPivot.value() * coef; - ++uItPivot; - } - ++lIt3; - } - // update diag -> contiguous - typename MatrixType::InnerLowerIterator lIt2(m_lu, col); - for (Index rrow = row + 1; rrow < m_lu.rows(); rrow++) { - typename MatrixType::InnerUpperIterator uItPivot(m_lu, row); - typename MatrixType::InnerUpperIterator uIt(m_lu, rrow); - const double coef = lIt2.value(); - - uItPivot += (rrow - row - 1); - m_lu.coeffRefDiag(rrow) -= uItPivot.value() * coef; - ++lIt2; - } - } -} - -template <typename MatrixType> -void SkylineInplaceLU<MatrixType>::computeRowMajor() { - const size_t rows = m_lu.rows(); - const size_t cols = m_lu.cols(); - - eigen_assert(rows == cols && "We do not (yet) support rectangular LU."); - eigen_assert(m_lu.IsRowMajor && "You're trying to apply rowMajor decomposition on a ColMajor matrix !"); - - for (Index row = 0; row < rows; row++) { - typename MatrixType::InnerLowerIterator llIt(m_lu, row); - - for (Index col = llIt.col(); col < row; col++) { - if (m_lu.coeffExistLower(row, col)) { - const double diag = m_lu.coeffDiag(col); - - typename MatrixType::InnerLowerIterator lIt(m_lu, row); - typename MatrixType::InnerUpperIterator uIt(m_lu, col); - - const Index offset = lIt.col() - uIt.row(); - - Index stop = offset > 0 ? col - lIt.col() : col - uIt.row(); - - // #define VECTORIZE -#ifdef VECTORIZE - Map<VectorXd> rowVal(lIt.valuePtr() + (offset > 0 ? 0 : -offset), stop); - Map<VectorXd> colVal(uIt.valuePtr() + (offset > 0 ? offset : 0), stop); - - Scalar newCoeff = m_lu.coeffLower(row, col) - rowVal.dot(colVal); -#else - if (offset > 0) // Skip zero value of lIt - uIt += offset; - else // Skip zero values of uIt - lIt += -offset; - Scalar newCoeff = m_lu.coeffLower(row, col); - - for (Index k = 0; k < stop; ++k) { - const Scalar tmp = newCoeff; - newCoeff = tmp - lIt.value() * uIt.value(); - ++lIt; - ++uIt; - } -#endif - - m_lu.coeffRefLower(row, col) = newCoeff / diag; - } - } - - // Upper matrix update - const Index col = row; - typename MatrixType::InnerUpperIterator uuIt(m_lu, col); - for (Index rrow = uuIt.row(); rrow < col; rrow++) { - typename MatrixType::InnerLowerIterator lIt(m_lu, rrow); - typename MatrixType::InnerUpperIterator uIt(m_lu, col); - const Index offset = lIt.col() - uIt.row(); - - Index stop = offset > 0 ? rrow - lIt.col() : rrow - uIt.row(); - -#ifdef VECTORIZE - Map<VectorXd> rowVal(lIt.valuePtr() + (offset > 0 ? 0 : -offset), stop); - Map<VectorXd> colVal(uIt.valuePtr() + (offset > 0 ? offset : 0), stop); - - Scalar newCoeff = m_lu.coeffUpper(rrow, col) - rowVal.dot(colVal); -#else - if (offset > 0) // Skip zero value of lIt - uIt += offset; - else // Skip zero values of uIt - lIt += -offset; - Scalar newCoeff = m_lu.coeffUpper(rrow, col); - for (Index k = 0; k < stop; ++k) { - const Scalar tmp = newCoeff; - newCoeff = tmp - lIt.value() * uIt.value(); - - ++lIt; - ++uIt; - } -#endif - m_lu.coeffRefUpper(rrow, col) = newCoeff; - } - - // Diag matrix update - typename MatrixType::InnerLowerIterator lIt(m_lu, row); - typename MatrixType::InnerUpperIterator uIt(m_lu, row); - - const Index offset = lIt.col() - uIt.row(); - - Index stop = offset > 0 ? lIt.size() : uIt.size(); -#ifdef VECTORIZE - Map<VectorXd> rowVal(lIt.valuePtr() + (offset > 0 ? 0 : -offset), stop); - Map<VectorXd> colVal(uIt.valuePtr() + (offset > 0 ? offset : 0), stop); - Scalar newCoeff = m_lu.coeffDiag(row) - rowVal.dot(colVal); -#else - if (offset > 0) // Skip zero value of lIt - uIt += offset; - else // Skip zero values of uIt - lIt += -offset; - Scalar newCoeff = m_lu.coeffDiag(row); - for (Index k = 0; k < stop; ++k) { - const Scalar tmp = newCoeff; - newCoeff = tmp - lIt.value() * uIt.value(); - ++lIt; - ++uIt; - } -#endif - m_lu.coeffRefDiag(row) = newCoeff; - } -} - -/** Computes *x = U^-1 L^-1 b - * - * If \a transpose is set to SvTranspose or SvAdjoint, the solution - * of the transposed/adjoint system is computed instead. - * - * Not all backends implement the solution of the transposed or - * adjoint system. - */ -template <typename MatrixType> -template <typename BDerived, typename XDerived> -bool SkylineInplaceLU<MatrixType>::solve(const MatrixBase<BDerived>& b, MatrixBase<XDerived>* x, - const int transposed) const { - const size_t rows = m_lu.rows(); - const size_t cols = m_lu.cols(); - - for (Index row = 0; row < rows; row++) { - x->coeffRef(row) = b.coeff(row); - Scalar newVal = x->coeff(row); - typename MatrixType::InnerLowerIterator lIt(m_lu, row); - - Index col = lIt.col(); - while (lIt.col() < row) { - newVal -= x->coeff(col++) * lIt.value(); - ++lIt; - } - - x->coeffRef(row) = newVal; - } - - for (Index col = rows - 1; col > 0; col--) { - x->coeffRef(col) = x->coeff(col) / m_lu.coeffDiag(col); - - const Scalar x_col = x->coeff(col); - - typename MatrixType::InnerUpperIterator uIt(m_lu, col); - uIt += uIt.size() - 1; - - while (uIt) { - x->coeffRef(uIt.row()) -= x_col * uIt.value(); - // TODO : introduce --operator - uIt += -1; - } - } - x->coeffRef(0) = x->coeff(0) / m_lu.coeffDiag(0); - - return true; -} - -} // end namespace Eigen - -#endif // EIGEN_SKYLINEINPLACELU_H
diff --git a/unsupported/Eigen/src/Skyline/SkylineMatrix.h b/unsupported/Eigen/src/Skyline/SkylineMatrix.h deleted file mode 100644 index 6615cc7..0000000 --- a/unsupported/Eigen/src/Skyline/SkylineMatrix.h +++ /dev/null
@@ -1,763 +0,0 @@ -// This file is part of Eigen, a lightweight C++ template library -// for linear algebra. -// -// Copyright (C) 2008-2009 Guillaume Saupin <guillaume.saupin@cea.fr> -// -// 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_SKYLINEMATRIX_H -#define EIGEN_SKYLINEMATRIX_H - -#include "SkylineStorage.h" -#include "SkylineMatrixBase.h" - -// IWYU pragma: private -#include "./InternalHeaderCheck.h" - -namespace Eigen { - -/** \ingroup Skyline_Module - * - * \class SkylineMatrix - * - * \brief The main skyline matrix class - * - * This class implements a skyline matrix using the very uncommon storage - * scheme. - * - * \param Scalar_ the scalar type, i.e. the type of the coefficients - * \param Options_ Union of bit flags controlling the storage scheme. Currently the only possibility - * is RowMajor. The default is 0 which means column-major. - * - * - */ -namespace internal { -template <typename Scalar_, int Options_> -struct traits<SkylineMatrix<Scalar_, Options_> > { - typedef Scalar_ Scalar; - typedef Sparse StorageKind; - - enum { - RowsAtCompileTime = Dynamic, - ColsAtCompileTime = Dynamic, - MaxRowsAtCompileTime = Dynamic, - MaxColsAtCompileTime = Dynamic, - Flags = SkylineBit | Options_, - CoeffReadCost = NumTraits<Scalar>::ReadCost, - }; -}; -} // namespace internal - -template <typename Scalar_, int Options_> -class SkylineMatrix : public SkylineMatrixBase<SkylineMatrix<Scalar_, Options_> > { - public: - EIGEN_SKYLINE_GENERIC_PUBLIC_INTERFACE(SkylineMatrix) - EIGEN_SKYLINE_INHERIT_ASSIGNMENT_OPERATOR(SkylineMatrix, +=) - EIGEN_SKYLINE_INHERIT_ASSIGNMENT_OPERATOR(SkylineMatrix, -=) - - using Base::IsRowMajor; - - protected: - typedef SkylineMatrix<Scalar, (Flags & ~RowMajorBit) | (IsRowMajor ? RowMajorBit : 0)> TransposedSkylineMatrix; - - Index m_outerSize; - Index m_innerSize; - - public: - Index* m_colStartIndex; - Index* m_rowStartIndex; - SkylineStorage<Scalar> m_data; - - public: - inline Index rows() const { return IsRowMajor ? m_outerSize : m_innerSize; } - - inline Index cols() const { return IsRowMajor ? m_innerSize : m_outerSize; } - - inline Index innerSize() const { return m_innerSize; } - - inline Index outerSize() const { return m_outerSize; } - - inline Index upperNonZeros() const { return m_data.upperSize(); } - - inline Index lowerNonZeros() const { return m_data.lowerSize(); } - - inline Index upperNonZeros(Index j) const { return m_colStartIndex[j + 1] - m_colStartIndex[j]; } - - inline Index lowerNonZeros(Index j) const { return m_rowStartIndex[j + 1] - m_rowStartIndex[j]; } - - inline const Scalar* _diagPtr() const { return &m_data.diag(0); } - - inline Scalar* _diagPtr() { return &m_data.diag(0); } - - inline const Scalar* _upperPtr() const { return &m_data.upper(0); } - - inline Scalar* _upperPtr() { return &m_data.upper(0); } - - inline const Scalar* _lowerPtr() const { return &m_data.lower(0); } - - inline Scalar* _lowerPtr() { return &m_data.lower(0); } - - inline const Index* _upperProfilePtr() const { return &m_data.upperProfile(0); } - - inline Index* _upperProfilePtr() { return &m_data.upperProfile(0); } - - inline const Index* _lowerProfilePtr() const { return &m_data.lowerProfile(0); } - - inline Index* _lowerProfilePtr() { return &m_data.lowerProfile(0); } - - inline Scalar coeff(Index row, Index col) const { - const Index outer = IsRowMajor ? row : col; - const Index inner = IsRowMajor ? col : row; - - eigen_assert(outer < outerSize()); - eigen_assert(inner < innerSize()); - - if (outer == inner) return this->m_data.diag(outer); - - if (IsRowMajor) { - if (inner > outer) // upper matrix - { - const Index minOuterIndex = inner - m_data.upperProfile(inner); - if (outer >= minOuterIndex) - return this->m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); - else - return Scalar(0); - } - if (inner < outer) // lower matrix - { - const Index minInnerIndex = outer - m_data.lowerProfile(outer); - if (inner >= minInnerIndex) - return this->m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); - else - return Scalar(0); - } - return m_data.upper(m_colStartIndex[inner] + outer - inner); - } else { - if (outer > inner) // upper matrix - { - const Index maxOuterIndex = inner + m_data.upperProfile(inner); - if (outer <= maxOuterIndex) - return this->m_data.upper(m_colStartIndex[inner] + (outer - inner)); - else - return Scalar(0); - } - if (outer < inner) // lower matrix - { - const Index maxInnerIndex = outer + m_data.lowerProfile(outer); - - if (inner <= maxInnerIndex) - return this->m_data.lower(m_rowStartIndex[outer] + (inner - outer)); - else - return Scalar(0); - } - } - } - - inline Scalar& coeffRef(Index row, Index col) { - const Index outer = IsRowMajor ? row : col; - const Index inner = IsRowMajor ? col : row; - - eigen_assert(outer < outerSize()); - eigen_assert(inner < innerSize()); - - if (outer == inner) return this->m_data.diag(outer); - - if (IsRowMajor) { - if (col > row) // upper matrix - { - const Index minOuterIndex = inner - m_data.upperProfile(inner); - eigen_assert(outer >= minOuterIndex && "You tried to access a coeff that does not exist in the storage"); - return this->m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); - } - if (col < row) // lower matrix - { - const Index minInnerIndex = outer - m_data.lowerProfile(outer); - eigen_assert(inner >= minInnerIndex && "You tried to access a coeff that does not exist in the storage"); - return this->m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); - } - } else { - if (outer > inner) // upper matrix - { - const Index maxOuterIndex = inner + m_data.upperProfile(inner); - eigen_assert(outer <= maxOuterIndex && "You tried to access a coeff that does not exist in the storage"); - return this->m_data.upper(m_colStartIndex[inner] + (outer - inner)); - } - if (outer < inner) // lower matrix - { - const Index maxInnerIndex = outer + m_data.lowerProfile(outer); - eigen_assert(inner <= maxInnerIndex && "You tried to access a coeff that does not exist in the storage"); - return this->m_data.lower(m_rowStartIndex[outer] + (inner - outer)); - } - } - } - - inline Scalar coeffDiag(Index idx) const { - eigen_assert(idx < outerSize()); - eigen_assert(idx < innerSize()); - return this->m_data.diag(idx); - } - - inline Scalar coeffLower(Index row, Index col) const { - const Index outer = IsRowMajor ? row : col; - const Index inner = IsRowMajor ? col : row; - - eigen_assert(outer < outerSize()); - eigen_assert(inner < innerSize()); - eigen_assert(inner != outer); - - if (IsRowMajor) { - const Index minInnerIndex = outer - m_data.lowerProfile(outer); - if (inner >= minInnerIndex) - return this->m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); - else - return Scalar(0); - - } else { - const Index maxInnerIndex = outer + m_data.lowerProfile(outer); - if (inner <= maxInnerIndex) - return this->m_data.lower(m_rowStartIndex[outer] + (inner - outer)); - else - return Scalar(0); - } - } - - inline Scalar coeffUpper(Index row, Index col) const { - const Index outer = IsRowMajor ? row : col; - const Index inner = IsRowMajor ? col : row; - - eigen_assert(outer < outerSize()); - eigen_assert(inner < innerSize()); - eigen_assert(inner != outer); - - if (IsRowMajor) { - const Index minOuterIndex = inner - m_data.upperProfile(inner); - if (outer >= minOuterIndex) - return this->m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); - else - return Scalar(0); - } else { - const Index maxOuterIndex = inner + m_data.upperProfile(inner); - if (outer <= maxOuterIndex) - return this->m_data.upper(m_colStartIndex[inner] + (outer - inner)); - else - return Scalar(0); - } - } - - inline Scalar& coeffRefDiag(Index idx) { - eigen_assert(idx < outerSize()); - eigen_assert(idx < innerSize()); - return this->m_data.diag(idx); - } - - inline Scalar& coeffRefLower(Index row, Index col) { - const Index outer = IsRowMajor ? row : col; - const Index inner = IsRowMajor ? col : row; - - eigen_assert(outer < outerSize()); - eigen_assert(inner < innerSize()); - eigen_assert(inner != outer); - - if (IsRowMajor) { - const Index minInnerIndex = outer - m_data.lowerProfile(outer); - eigen_assert(inner >= minInnerIndex && "You tried to access a coeff that does not exist in the storage"); - return this->m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); - } else { - const Index maxInnerIndex = outer + m_data.lowerProfile(outer); - eigen_assert(inner <= maxInnerIndex && "You tried to access a coeff that does not exist in the storage"); - return this->m_data.lower(m_rowStartIndex[outer] + (inner - outer)); - } - } - - inline bool coeffExistLower(Index row, Index col) { - const Index outer = IsRowMajor ? row : col; - const Index inner = IsRowMajor ? col : row; - - eigen_assert(outer < outerSize()); - eigen_assert(inner < innerSize()); - eigen_assert(inner != outer); - - if (IsRowMajor) { - const Index minInnerIndex = outer - m_data.lowerProfile(outer); - return inner >= minInnerIndex; - } else { - const Index maxInnerIndex = outer + m_data.lowerProfile(outer); - return inner <= maxInnerIndex; - } - } - - inline Scalar& coeffRefUpper(Index row, Index col) { - const Index outer = IsRowMajor ? row : col; - const Index inner = IsRowMajor ? col : row; - - eigen_assert(outer < outerSize()); - eigen_assert(inner < innerSize()); - eigen_assert(inner != outer); - - if (IsRowMajor) { - const Index minOuterIndex = inner - m_data.upperProfile(inner); - eigen_assert(outer >= minOuterIndex && "You tried to access a coeff that does not exist in the storage"); - return this->m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); - } else { - const Index maxOuterIndex = inner + m_data.upperProfile(inner); - eigen_assert(outer <= maxOuterIndex && "You tried to access a coeff that does not exist in the storage"); - return this->m_data.upper(m_colStartIndex[inner] + (outer - inner)); - } - } - - inline bool coeffExistUpper(Index row, Index col) { - const Index outer = IsRowMajor ? row : col; - const Index inner = IsRowMajor ? col : row; - - eigen_assert(outer < outerSize()); - eigen_assert(inner < innerSize()); - eigen_assert(inner != outer); - - if (IsRowMajor) { - const Index minOuterIndex = inner - m_data.upperProfile(inner); - return outer >= minOuterIndex; - } else { - const Index maxOuterIndex = inner + m_data.upperProfile(inner); - return outer <= maxOuterIndex; - } - } - - protected: - public: - class InnerUpperIterator; - class InnerLowerIterator; - - class OuterUpperIterator; - class OuterLowerIterator; - - /** Removes all non zeros */ - inline void setZero() { - m_data.clear(); - std::fill_n(m_colStartIndex, m_outerSize + 1, Index(0)); - std::fill_n(m_rowStartIndex, m_outerSize + 1, Index(0)); - } - - /** \returns the number of non zero coefficients */ - inline Index nonZeros() const { return m_data.diagSize() + m_data.upperSize() + m_data.lowerSize(); } - - /** Preallocates \a reserveSize non zeros */ - inline void reserve(Index reserveSize, Index reserveUpperSize, Index reserveLowerSize) { - m_data.reserve(reserveSize, reserveUpperSize, reserveLowerSize); - } - - /** \returns a reference to a novel non zero coefficient with coordinates \a row x \a col. - - * - * \warning This function can be extremely slow if the non zero coefficients - * are not inserted in a coherent order. - * - * After an insertion session, you should call the finalize() function. - */ - EIGEN_DONT_INLINE Scalar& insert(Index row, Index col) { - const Index outer = IsRowMajor ? row : col; - const Index inner = IsRowMajor ? col : row; - - eigen_assert(outer < outerSize()); - eigen_assert(inner < innerSize()); - - if (outer == inner) return m_data.diag(col); - - if (IsRowMajor) { - if (outer < inner) // upper matrix - { - Index minOuterIndex = 0; - minOuterIndex = inner - m_data.upperProfile(inner); - - if (outer < minOuterIndex) // The value does not yet exist - { - const Index previousProfile = m_data.upperProfile(inner); - - m_data.upperProfile(inner) = inner - outer; - - const Index bandIncrement = m_data.upperProfile(inner) - previousProfile; - // shift data stored after this new one - const Index stop = m_colStartIndex[cols()]; - const Index start = m_colStartIndex[inner]; - - for (Index innerIdx = stop; innerIdx >= start; innerIdx--) { - m_data.upper(innerIdx + bandIncrement) = m_data.upper(innerIdx); - } - - for (Index innerIdx = cols(); innerIdx > inner; innerIdx--) { - m_colStartIndex[innerIdx] += bandIncrement; - } - - // zeros new data - std::fill_n(this->_upperPtr() + start, bandIncrement - 1, Scalar(0)); - - return m_data.upper(m_colStartIndex[inner]); - } else { - return m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); - } - } - - if (outer > inner) // lower matrix - { - const Index minInnerIndex = outer - m_data.lowerProfile(outer); - if (inner < minInnerIndex) // The value does not yet exist - { - const Index previousProfile = m_data.lowerProfile(outer); - m_data.lowerProfile(outer) = outer - inner; - - const Index bandIncrement = m_data.lowerProfile(outer) - previousProfile; - // shift data stored after this new one - const Index stop = m_rowStartIndex[rows()]; - const Index start = m_rowStartIndex[outer]; - - for (Index innerIdx = stop; innerIdx >= start; innerIdx--) { - m_data.lower(innerIdx + bandIncrement) = m_data.lower(innerIdx); - } - - for (Index innerIdx = rows(); innerIdx > outer; innerIdx--) { - m_rowStartIndex[innerIdx] += bandIncrement; - } - - // zeros new data - std::fill_n(this->_lowerPtr() + start, bandIncrement - 1, Scalar(0)); - return m_data.lower(m_rowStartIndex[outer]); - } else { - return m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); - } - } - } else { - if (outer > inner) // upper matrix - { - const Index maxOuterIndex = inner + m_data.upperProfile(inner); - if (outer > maxOuterIndex) // The value does not yet exist - { - const Index previousProfile = m_data.upperProfile(inner); - m_data.upperProfile(inner) = outer - inner; - - const Index bandIncrement = m_data.upperProfile(inner) - previousProfile; - // shift data stored after this new one - const Index stop = m_rowStartIndex[rows()]; - const Index start = m_rowStartIndex[inner + 1]; - - for (Index innerIdx = stop; innerIdx >= start; innerIdx--) { - m_data.upper(innerIdx + bandIncrement) = m_data.upper(innerIdx); - } - - for (Index innerIdx = inner + 1; innerIdx < outerSize() + 1; innerIdx++) { - m_rowStartIndex[innerIdx] += bandIncrement; - } - std::fill_n(this->_upperPtr() + m_rowStartIndex[inner] + previousProfile + 1, bandIncrement - 1, Scalar(0)); - return m_data.upper(m_rowStartIndex[inner] + m_data.upperProfile(inner)); - } else { - return m_data.upper(m_rowStartIndex[inner] + (outer - inner)); - } - } - - if (outer < inner) // lower matrix - { - const Index maxInnerIndex = outer + m_data.lowerProfile(outer); - if (inner > maxInnerIndex) // The value does not yet exist - { - const Index previousProfile = m_data.lowerProfile(outer); - m_data.lowerProfile(outer) = inner - outer; - - const Index bandIncrement = m_data.lowerProfile(outer) - previousProfile; - // shift data stored after this new one - const Index stop = m_colStartIndex[cols()]; - const Index start = m_colStartIndex[outer + 1]; - - for (Index innerIdx = stop; innerIdx >= start; innerIdx--) { - m_data.lower(innerIdx + bandIncrement) = m_data.lower(innerIdx); - } - - for (Index innerIdx = outer + 1; innerIdx < outerSize() + 1; innerIdx++) { - m_colStartIndex[innerIdx] += bandIncrement; - } - std::fill_n(this->_lowerPtr() + m_colStartIndex[outer] + previousProfile + 1, bandIncrement - 1, Scalar(0)); - return m_data.lower(m_colStartIndex[outer] + m_data.lowerProfile(outer)); - } else { - return m_data.lower(m_colStartIndex[outer] + (inner - outer)); - } - } - } - } - - /** Must be called after inserting a set of non zero entries. - */ - inline void finalize() { - if (IsRowMajor) { - if (rows() > cols()) - m_data.resize(cols(), cols(), rows(), m_colStartIndex[cols()] + 1, m_rowStartIndex[rows()] + 1); - else - m_data.resize(rows(), cols(), rows(), m_colStartIndex[cols()] + 1, m_rowStartIndex[rows()] + 1); - - // eigen_assert(rows() == cols() && "memory reorganisatrion only works with suare matrix"); - // - // Scalar* newArray = new Scalar[m_colStartIndex[cols()] + 1 + m_rowStartIndex[rows()] + 1]; - // Index dataIdx = 0; - // for (Index row = 0; row < rows(); row++) { - // - // const Index nbLowerElts = m_rowStartIndex[row + 1] - m_rowStartIndex[row]; - // // std::cout << "nbLowerElts" << nbLowerElts << std::endl; - // memcpy(newArray + dataIdx, m_data.m_lower + m_rowStartIndex[row], nbLowerElts * sizeof - // (Scalar)); m_rowStartIndex[row] = dataIdx; dataIdx += nbLowerElts; - // - // const Index nbUpperElts = m_colStartIndex[row + 1] - m_colStartIndex[row]; - // memcpy(newArray + dataIdx, m_data.m_upper + m_colStartIndex[row], nbUpperElts * sizeof - // (Scalar)); m_colStartIndex[row] = dataIdx; dataIdx += nbUpperElts; - // - // - // } - // //todo : don't access m_data profile directly : add an accessor from SkylineMatrix - // m_rowStartIndex[rows()] = m_rowStartIndex[rows()-1] + m_data.lowerProfile(rows()-1); - // m_colStartIndex[cols()] = m_colStartIndex[cols()-1] + m_data.upperProfile(cols()-1); - // - // delete[] m_data.m_lower; - // delete[] m_data.m_upper; - // - // m_data.m_lower = newArray; - // m_data.m_upper = newArray; - } else { - if (rows() > cols()) - m_data.resize(cols(), rows(), cols(), m_rowStartIndex[cols()] + 1, m_colStartIndex[cols()] + 1); - else - m_data.resize(rows(), rows(), cols(), m_rowStartIndex[rows()] + 1, m_colStartIndex[rows()] + 1); - } - } - - inline void squeeze() { - finalize(); - m_data.squeeze(); - } - - void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar>()) { - // TODO - } - - /** Resizes the matrix to a \a rows x \a cols matrix and initializes it to zero - * \sa resizeNonZeros(Index), reserve(), setZero() - */ - void resize(size_t rows, size_t cols) { - const Index diagSize = rows > cols ? cols : rows; - m_innerSize = IsRowMajor ? cols : rows; - - eigen_assert(rows == cols && "Skyline matrix must be square matrix"); - - if (diagSize % 2) { // diagSize is odd - const Index k = (diagSize - 1) / 2; - - m_data.resize(diagSize, IsRowMajor ? cols : rows, IsRowMajor ? rows : cols, 2 * k * k + k + 1, 2 * k * k + k + 1); - - } else // diagSize is even - { - const Index k = diagSize / 2; - m_data.resize(diagSize, IsRowMajor ? cols : rows, IsRowMajor ? rows : cols, 2 * k * k - k + 1, 2 * k * k - k + 1); - } - - if (m_colStartIndex && m_rowStartIndex) { - delete[] m_colStartIndex; - delete[] m_rowStartIndex; - } - m_colStartIndex = new Index[cols + 1]; - m_rowStartIndex = new Index[rows + 1]; - m_outerSize = diagSize; - - m_data.reset(); - m_data.clear(); - - m_outerSize = diagSize; - std::fill_n(m_colStartIndex, cols + 1, Index(0)); - std::fill_n(m_rowStartIndex, rows + 1, Index(0)); - } - - void resizeNonZeros(Index size) { m_data.resize(size); } - - inline SkylineMatrix() : m_outerSize(-1), m_innerSize(0), m_colStartIndex(0), m_rowStartIndex(0) { resize(0, 0); } - - inline SkylineMatrix(size_t rows, size_t cols) - : m_outerSize(0), m_innerSize(0), m_colStartIndex(0), m_rowStartIndex(0) { - resize(rows, cols); - } - - template <typename OtherDerived> - inline SkylineMatrix(const SkylineMatrixBase<OtherDerived>& other) - : m_outerSize(0), m_innerSize(0), m_colStartIndex(0), m_rowStartIndex(0) { - *this = other.derived(); - } - - inline SkylineMatrix(const SkylineMatrix& other) - : Base(), m_outerSize(0), m_innerSize(0), m_colStartIndex(0), m_rowStartIndex(0) { - *this = other.derived(); - } - - inline void swap(SkylineMatrix& other) { - // EIGEN_DBG_SKYLINE(std::cout << "SkylineMatrix:: swap\n"); - std::swap(m_colStartIndex, other.m_colStartIndex); - std::swap(m_rowStartIndex, other.m_rowStartIndex); - std::swap(m_innerSize, other.m_innerSize); - std::swap(m_outerSize, other.m_outerSize); - m_data.swap(other.m_data); - } - - inline SkylineMatrix& operator=(const SkylineMatrix& other) { - std::cout << "SkylineMatrix& operator=(const SkylineMatrix& other)\n"; - if (other.isRValue()) { - swap(other.const_cast_derived()); - } else { - resize(other.rows(), other.cols()); - memcpy(m_colStartIndex, other.m_colStartIndex, (m_outerSize + 1) * sizeof(Index)); - memcpy(m_rowStartIndex, other.m_rowStartIndex, (m_outerSize + 1) * sizeof(Index)); - m_data = other.m_data; - } - return *this; - } - - template <typename OtherDerived> - inline SkylineMatrix& operator=(const SkylineMatrixBase<OtherDerived>& other) { - const bool needToTranspose = (Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit); - if (needToTranspose) { - // TODO - // return *this; - } else { - // there is no special optimization - return SkylineMatrixBase<SkylineMatrix>::operator=(other.derived()); - } - } - - friend std::ostream& operator<<(std::ostream& s, const SkylineMatrix& m) { - EIGEN_DBG_SKYLINE( - std::cout << "upper elements : " << std::endl; - for (Index i = 0; i < m.m_data.upperSize(); i++) std::cout << m.m_data.upper(i) << "\t"; std::cout << std::endl; - std::cout << "upper profile : " << std::endl; - for (Index i = 0; i < m.m_data.upperProfileSize(); i++) std::cout << m.m_data.upperProfile(i) << "\t"; - std::cout << std::endl; std::cout << "lower startIdx : " << std::endl; - for (Index i = 0; i < m.m_data.upperProfileSize(); i++) std::cout - << (IsRowMajor ? m.m_colStartIndex[i] : m.m_rowStartIndex[i]) << "\t"; - std::cout << std::endl; - - std::cout << "lower elements : " << std::endl; - for (Index i = 0; i < m.m_data.lowerSize(); i++) std::cout << m.m_data.lower(i) << "\t"; std::cout << std::endl; - std::cout << "lower profile : " << std::endl; - for (Index i = 0; i < m.m_data.lowerProfileSize(); i++) std::cout << m.m_data.lowerProfile(i) << "\t"; - std::cout << std::endl; std::cout << "lower startIdx : " << std::endl; - for (Index i = 0; i < m.m_data.lowerProfileSize(); i++) std::cout - << (IsRowMajor ? m.m_rowStartIndex[i] : m.m_colStartIndex[i]) << "\t"; - std::cout << std::endl;); - for (Index rowIdx = 0; rowIdx < m.rows(); rowIdx++) { - for (Index colIdx = 0; colIdx < m.cols(); colIdx++) { - s << m.coeff(rowIdx, colIdx) << "\t"; - } - s << std::endl; - } - return s; - } - - /** Destructor */ - inline ~SkylineMatrix() { - delete[] m_colStartIndex; - delete[] m_rowStartIndex; - } - - /** Overloaded for performance */ - Scalar sum() const; -}; - -template <typename Scalar, int Options_> -class SkylineMatrix<Scalar, Options_>::InnerUpperIterator { - public: - InnerUpperIterator(const SkylineMatrix& mat, Index outer) - : m_matrix(mat), - m_outer(outer), - m_id(Options_ == RowMajor ? mat.m_colStartIndex[outer] : mat.m_rowStartIndex[outer] + 1), - m_start(m_id), - m_end(Options_ == RowMajor ? mat.m_colStartIndex[outer + 1] : mat.m_rowStartIndex[outer + 1] + 1) {} - - inline InnerUpperIterator& operator++() { - m_id++; - return *this; - } - - inline InnerUpperIterator& operator+=(Index shift) { - m_id += shift; - return *this; - } - - inline Scalar value() const { return m_matrix.m_data.upper(m_id); } - - inline Scalar* valuePtr() { return const_cast<Scalar*>(&(m_matrix.m_data.upper(m_id))); } - - inline Scalar& valueRef() { return const_cast<Scalar&>(m_matrix.m_data.upper(m_id)); } - - inline Index index() const { - return IsRowMajor ? m_outer - m_matrix.m_data.upperProfile(m_outer) + (m_id - m_start) - : m_outer + (m_id - m_start) + 1; - } - - inline Index row() const { return IsRowMajor ? index() : m_outer; } - - inline Index col() const { return IsRowMajor ? m_outer : index(); } - - inline size_t size() const { return m_matrix.m_data.upperProfile(m_outer); } - - inline operator bool() const { return (m_id < m_end) && (m_id >= m_start); } - - protected: - const SkylineMatrix& m_matrix; - const Index m_outer; - Index m_id; - const Index m_start; - const Index m_end; -}; - -template <typename Scalar, int Options_> -class SkylineMatrix<Scalar, Options_>::InnerLowerIterator { - public: - InnerLowerIterator(const SkylineMatrix& mat, Index outer) - : m_matrix(mat), - m_outer(outer), - m_id(Options_ == RowMajor ? mat.m_rowStartIndex[outer] : mat.m_colStartIndex[outer] + 1), - m_start(m_id), - m_end(Options_ == RowMajor ? mat.m_rowStartIndex[outer + 1] : mat.m_colStartIndex[outer + 1] + 1) {} - - inline InnerLowerIterator& operator++() { - m_id++; - return *this; - } - - inline InnerLowerIterator& operator+=(Index shift) { - m_id += shift; - return *this; - } - - inline Scalar value() const { return m_matrix.m_data.lower(m_id); } - - inline Scalar* valuePtr() { return const_cast<Scalar*>(&(m_matrix.m_data.lower(m_id))); } - - inline Scalar& valueRef() { return const_cast<Scalar&>(m_matrix.m_data.lower(m_id)); } - - inline Index index() const { - return IsRowMajor ? m_outer - m_matrix.m_data.lowerProfile(m_outer) + (m_id - m_start) - : m_outer + (m_id - m_start) + 1; - ; - } - - inline Index row() const { return IsRowMajor ? m_outer : index(); } - - inline Index col() const { return IsRowMajor ? index() : m_outer; } - - inline size_t size() const { return m_matrix.m_data.lowerProfile(m_outer); } - - inline operator bool() const { return (m_id < m_end) && (m_id >= m_start); } - - protected: - const SkylineMatrix& m_matrix; - const Index m_outer; - Index m_id; - const Index m_start; - const Index m_end; -}; - -} // end namespace Eigen - -#endif // EIGEN_SKYLINEMATRIX_H
diff --git a/unsupported/Eigen/src/Skyline/SkylineMatrixBase.h b/unsupported/Eigen/src/Skyline/SkylineMatrixBase.h deleted file mode 100644 index 1927458..0000000 --- a/unsupported/Eigen/src/Skyline/SkylineMatrixBase.h +++ /dev/null
@@ -1,188 +0,0 @@ -// This file is part of Eigen, a lightweight C++ template library -// for linear algebra. -// -// Copyright (C) 2008-2009 Guillaume Saupin <guillaume.saupin@cea.fr> -// -// 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_SKYLINEMATRIXBASE_H -#define EIGEN_SKYLINEMATRIXBASE_H - -#include "SkylineUtil.h" - -// IWYU pragma: private -#include "./InternalHeaderCheck.h" - -namespace Eigen { - -/** \ingroup Skyline_Module - * - * \class SkylineMatrixBase - * - * \brief Base class of any skyline matrices or skyline expressions - * - * \param Derived - * - */ -template <typename Derived> -class SkylineMatrixBase : public EigenBase<Derived> { - public: - typedef typename internal::traits<Derived>::Scalar Scalar; - typedef typename internal::traits<Derived>::StorageKind StorageKind; - typedef typename internal::index<StorageKind>::type Index; - - enum { - RowsAtCompileTime = internal::traits<Derived>::RowsAtCompileTime, - /**< The number of rows at compile-time. This is just a copy of the value provided - * by the \a Derived type. If a value is not known at compile-time, - * it is set to the \a Dynamic constant. - * \sa MatrixBase::rows(), MatrixBase::cols(), ColsAtCompileTime, SizeAtCompileTime */ - - ColsAtCompileTime = internal::traits<Derived>::ColsAtCompileTime, - /**< The number of columns at compile-time. This is just a copy of the value provided - * by the \a Derived type. If a value is not known at compile-time, - * it is set to the \a Dynamic constant. - * \sa MatrixBase::rows(), MatrixBase::cols(), RowsAtCompileTime, SizeAtCompileTime */ - - SizeAtCompileTime = (internal::size_of_xpr_at_compile_time<Derived>::ret), - /**< This is equal to the number of coefficients, i.e. the number of - * rows times the number of columns, or to \a Dynamic if this is not - * known at compile-time. \sa RowsAtCompileTime, ColsAtCompileTime */ - - MaxRowsAtCompileTime = RowsAtCompileTime, - MaxColsAtCompileTime = ColsAtCompileTime, - - MaxSizeAtCompileTime = (internal::size_at_compile_time(MaxRowsAtCompileTime, MaxColsAtCompileTime)), - - IsVectorAtCompileTime = RowsAtCompileTime == 1 || ColsAtCompileTime == 1, - /**< This is set to true if either the number of rows or the number of - * columns is known at compile-time to be equal to 1. Indeed, in that case, - * we are dealing with a column-vector (if there is only one column) or with - * a row-vector (if there is only one row). */ - - Flags = internal::traits<Derived>::Flags, - /**< This stores expression \ref flags flags which may or may not be inherited by new expressions - * constructed from this one. See the \ref flags "list of flags". - */ - - CoeffReadCost = internal::traits<Derived>::CoeffReadCost, - /**< This is a rough measure of how expensive it is to read one coefficient from - * this expression. - */ - - IsRowMajor = Flags & RowMajorBit ? 1 : 0 - }; - -#ifndef EIGEN_PARSED_BY_DOXYGEN - /** This is the "real scalar" type; if the \a Scalar type is already real numbers - * (e.g. int, float or double) then \a RealScalar is just the same as \a Scalar. If - * \a Scalar is \a std::complex<T> then RealScalar is \a T. - * - * \sa class NumTraits - */ - typedef typename NumTraits<Scalar>::Real RealScalar; - - /** type of the equivalent square matrix */ - typedef Matrix<Scalar, internal::max_size_prefer_dynamic(RowsAtCompileTime, ColsAtCompileTime), - internal::max_size_prefer_dynamic(RowsAtCompileTime, ColsAtCompileTime)> - SquareMatrixType; - - inline const Derived& derived() const { return *static_cast<const Derived*>(this); } - - inline Derived& derived() { return *static_cast<Derived*>(this); } - - inline Derived& const_cast_derived() const { return *static_cast<Derived*>(const_cast<SkylineMatrixBase*>(this)); } -#endif // not EIGEN_PARSED_BY_DOXYGEN - - /** \returns the number of rows. \sa cols(), RowsAtCompileTime */ - inline EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT { return derived().rows(); } - - /** \returns the number of columns. \sa rows(), ColsAtCompileTime*/ - inline EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT { return derived().cols(); } - - /** \returns the number of coefficients, which is \a rows()*cols(). - * \sa rows(), cols(), SizeAtCompileTime. */ - inline EIGEN_CONSTEXPR Index size() const EIGEN_NOEXCEPT { return rows() * cols(); } - - /** \returns the number of nonzero coefficients which is in practice the number - * of stored coefficients. */ - inline Index nonZeros() const { return derived().nonZeros(); } - - /** \returns the size of the storage major dimension, - * i.e., the number of columns for a columns major matrix, and the number of rows otherwise */ - Index outerSize() const { return (int(Flags) & RowMajorBit) ? this->rows() : this->cols(); } - - /** \returns the size of the inner dimension according to the storage order, - * i.e., the number of rows for a columns major matrix, and the number of cols otherwise */ - Index innerSize() const { return (int(Flags) & RowMajorBit) ? this->cols() : this->rows(); } - - bool isRValue() const { return m_isRValue; } - - Derived& markAsRValue() { - m_isRValue = true; - return derived(); - } - - SkylineMatrixBase() : m_isRValue(false) { /* TODO check flags */ - } - - inline Derived& operator=(const Derived& other) { - this->operator= <Derived>(other); - return derived(); - } - - template <typename OtherDerived> - inline void assignGeneric(const OtherDerived& other) { - derived().resize(other.rows(), other.cols()); - for (Index row = 0; row < rows(); row++) - for (Index col = 0; col < cols(); col++) { - if (other.coeff(row, col) != Scalar(0)) derived().insert(row, col) = other.coeff(row, col); - } - derived().finalize(); - } - - template <typename OtherDerived> - inline Derived& operator=(const SkylineMatrixBase<OtherDerived>& other) { - // TODO - } - - template <typename Lhs, typename Rhs> - inline Derived& operator=(const SkylineProduct<Lhs, Rhs, SkylineTimeSkylineProduct>& product); - - friend std::ostream& operator<<(std::ostream& s, const SkylineMatrixBase& m) { - s << m.derived(); - return s; - } - - template <typename OtherDerived> - const typename SkylineProductReturnType<Derived, OtherDerived>::Type operator*( - const MatrixBase<OtherDerived>& other) const; - - /** \internal use operator= */ - template <typename DenseDerived> - void evalTo(MatrixBase<DenseDerived>& dst) const { - dst.setZero(); - for (Index i = 0; i < rows(); i++) - for (Index j = 0; j < rows(); j++) dst(i, j) = derived().coeff(i, j); - } - - Matrix<Scalar, RowsAtCompileTime, ColsAtCompileTime> toDense() const { return derived(); } - - /** \returns the matrix or vector obtained by evaluating this expression. - * - * Notice that in the case of a plain matrix or vector (not an expression) this function just returns - * a const reference, in order to avoid a useless copy. - */ - EIGEN_STRONG_INLINE const typename internal::eval<Derived, IsSkyline>::type eval() const { - return typename internal::eval<Derived>::type(derived()); - } - - protected: - bool m_isRValue; -}; - -} // end namespace Eigen - -#endif // EIGEN_SKYLINEMATRIXBASE_H
diff --git a/unsupported/Eigen/src/Skyline/SkylineProduct.h b/unsupported/Eigen/src/Skyline/SkylineProduct.h deleted file mode 100644 index 50be597..0000000 --- a/unsupported/Eigen/src/Skyline/SkylineProduct.h +++ /dev/null
@@ -1,262 +0,0 @@ -// This file is part of Eigen, a lightweight C++ template library -// for linear algebra. -// -// Copyright (C) 2008-2009 Guillaume Saupin <guillaume.saupin@cea.fr> -// -// 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_SKYLINEPRODUCT_H -#define EIGEN_SKYLINEPRODUCT_H - -// IWYU pragma: private -#include "./InternalHeaderCheck.h" - -namespace Eigen { - -template <typename Lhs, typename Rhs, int ProductMode> -struct SkylineProductReturnType { - typedef const typename internal::nested_eval<Lhs, Rhs::RowsAtCompileTime>::type LhsNested; - typedef const typename internal::nested_eval<Rhs, Lhs::RowsAtCompileTime>::type RhsNested; - - typedef SkylineProduct<LhsNested, RhsNested, ProductMode> Type; -}; - -template <typename LhsNested, typename RhsNested, int ProductMode> -struct internal::traits<SkylineProduct<LhsNested, RhsNested, ProductMode>> { - // clean the nested types: - typedef internal::remove_all_t<LhsNested> LhsNested_; - typedef internal::remove_all_t<RhsNested> RhsNested_; - typedef typename LhsNested_::Scalar Scalar; - - enum { - LhsCoeffReadCost = LhsNested_::CoeffReadCost, - RhsCoeffReadCost = RhsNested_::CoeffReadCost, - LhsFlags = LhsNested_::Flags, - RhsFlags = RhsNested_::Flags, - - RowsAtCompileTime = LhsNested_::RowsAtCompileTime, - ColsAtCompileTime = RhsNested_::ColsAtCompileTime, - InnerSize = internal::min_size_prefer_fixed(LhsNested_::ColsAtCompileTime, RhsNested_::RowsAtCompileTime), - - MaxRowsAtCompileTime = LhsNested_::MaxRowsAtCompileTime, - MaxColsAtCompileTime = RhsNested_::MaxColsAtCompileTime, - - EvalToRowMajor = (RhsFlags & LhsFlags & RowMajorBit), - ResultIsSkyline = ProductMode == SkylineTimeSkylineProduct, - - RemovedBits = ~((EvalToRowMajor ? 0 : RowMajorBit) | (ResultIsSkyline ? 0 : SkylineBit)), - - Flags = (int(LhsFlags | RhsFlags) & HereditaryBits & RemovedBits) | EvalBeforeAssigningBit | EvalBeforeNestingBit, - - CoeffReadCost = HugeCost - }; - - typedef std::conditional_t<ResultIsSkyline, SkylineMatrixBase<SkylineProduct<LhsNested, RhsNested, ProductMode>>, - MatrixBase<SkylineProduct<LhsNested, RhsNested, ProductMode>>> - Base; -}; - -namespace internal { -template <typename LhsNested, typename RhsNested, int ProductMode> -class SkylineProduct : no_assignment_operator, public traits<SkylineProduct<LhsNested, RhsNested, ProductMode>>::Base { - public: - EIGEN_GENERIC_PUBLIC_INTERFACE(SkylineProduct) - - private: - typedef typename traits<SkylineProduct>::LhsNested_ LhsNested_; - typedef typename traits<SkylineProduct>::RhsNested_ RhsNested_; - - public: - template <typename Lhs, typename Rhs> - EIGEN_STRONG_INLINE SkylineProduct(const Lhs& lhs, const Rhs& rhs) : m_lhs(lhs), m_rhs(rhs) { - eigen_assert(lhs.cols() == rhs.rows()); - - enum { - ProductIsValid = LhsNested_::ColsAtCompileTime == Dynamic || RhsNested_::RowsAtCompileTime == Dynamic || - int(LhsNested_::ColsAtCompileTime) == int(RhsNested_::RowsAtCompileTime), - AreVectors = LhsNested_::IsVectorAtCompileTime && RhsNested_::IsVectorAtCompileTime, - SameSizes = EIGEN_PREDICATE_SAME_MATRIX_SIZE(LhsNested_, RhsNested_) - }; - // note to the lost user: - // * for a dot product use: v1.dot(v2) - // * for a coeff-wise product use: v1.cwise()*v2 - EIGEN_STATIC_ASSERT( - ProductIsValid || !(AreVectors && SameSizes), - INVALID_VECTOR_VECTOR_PRODUCT__IF_YOU_WANTED_A_DOT_OR_COEFF_WISE_PRODUCT_YOU_MUST_USE_THE_EXPLICIT_FUNCTIONS) - EIGEN_STATIC_ASSERT(ProductIsValid || !(SameSizes && !AreVectors), - INVALID_MATRIX_PRODUCT__IF_YOU_WANTED_A_COEFF_WISE_PRODUCT_YOU_MUST_USE_THE_EXPLICIT_FUNCTION) - EIGEN_STATIC_ASSERT(ProductIsValid || SameSizes, INVALID_MATRIX_PRODUCT) - } - - EIGEN_STRONG_INLINE Index rows() const { return m_lhs.rows(); } - - EIGEN_STRONG_INLINE Index cols() const { return m_rhs.cols(); } - - EIGEN_STRONG_INLINE const LhsNested_& lhs() const { return m_lhs; } - - EIGEN_STRONG_INLINE const RhsNested_& rhs() const { return m_rhs; } - - protected: - LhsNested m_lhs; - RhsNested m_rhs; -}; - -// dense = skyline * dense -// Note that here we force no inlining and separate the setZero() because GCC messes up otherwise - -template <typename Lhs, typename Rhs, typename Dest> -EIGEN_DONT_INLINE void skyline_row_major_time_dense_product(const Lhs& lhs, const Rhs& rhs, Dest& dst) { - typedef remove_all_t<Lhs> Lhs_; - typedef remove_all_t<Rhs> Rhs_; - typedef typename traits<Lhs>::Scalar Scalar; - - enum { - LhsIsRowMajor = (Lhs_::Flags & RowMajorBit) == RowMajorBit, - LhsIsSelfAdjoint = (Lhs_::Flags & SelfAdjointBit) == SelfAdjointBit, - ProcessFirstHalf = LhsIsSelfAdjoint && (((Lhs_::Flags & (UpperTriangularBit | LowerTriangularBit)) == 0) || - ((Lhs_::Flags & UpperTriangularBit) && !LhsIsRowMajor) || - ((Lhs_::Flags & LowerTriangularBit) && LhsIsRowMajor)), - ProcessSecondHalf = LhsIsSelfAdjoint && (!ProcessFirstHalf) - }; - - // Use matrix diagonal part <- Improvement : use inner iterator on dense matrix. - for (Index col = 0; col < rhs.cols(); col++) { - for (Index row = 0; row < lhs.rows(); row++) { - dst(row, col) = lhs.coeffDiag(row) * rhs(row, col); - } - } - // Use matrix lower triangular part - for (Index row = 0; row < lhs.rows(); row++) { - typename Lhs_::InnerLowerIterator lIt(lhs, row); - const Index stop = lIt.col() + lIt.size(); - for (Index col = 0; col < rhs.cols(); col++) { - Index k = lIt.col(); - Scalar tmp = 0; - while (k < stop) { - tmp += lIt.value() * rhs(k++, col); - ++lIt; - } - dst(row, col) += tmp; - lIt += -lIt.size(); - } - } - - // Use matrix upper triangular part - for (Index lhscol = 0; lhscol < lhs.cols(); lhscol++) { - typename Lhs_::InnerUpperIterator uIt(lhs, lhscol); - const Index stop = uIt.size() + uIt.row(); - for (Index rhscol = 0; rhscol < rhs.cols(); rhscol++) { - const Scalar rhsCoeff = rhs.coeff(lhscol, rhscol); - Index k = uIt.row(); - while (k < stop) { - dst(k++, rhscol) += uIt.value() * rhsCoeff; - ++uIt; - } - uIt += -uIt.size(); - } - } -} - -template <typename Lhs, typename Rhs, typename Dest> -EIGEN_DONT_INLINE void skyline_col_major_time_dense_product(const Lhs& lhs, const Rhs& rhs, Dest& dst) { - typedef remove_all_t<Lhs> Lhs_; - typedef remove_all_t<Rhs> Rhs_; - typedef typename traits<Lhs>::Scalar Scalar; - - enum { - LhsIsRowMajor = (Lhs_::Flags & RowMajorBit) == RowMajorBit, - LhsIsSelfAdjoint = (Lhs_::Flags & SelfAdjointBit) == SelfAdjointBit, - ProcessFirstHalf = LhsIsSelfAdjoint && (((Lhs_::Flags & (UpperTriangularBit | LowerTriangularBit)) == 0) || - ((Lhs_::Flags & UpperTriangularBit) && !LhsIsRowMajor) || - ((Lhs_::Flags & LowerTriangularBit) && LhsIsRowMajor)), - ProcessSecondHalf = LhsIsSelfAdjoint && (!ProcessFirstHalf) - }; - - // Use matrix diagonal part <- Improvement : use inner iterator on dense matrix. - for (Index col = 0; col < rhs.cols(); col++) { - for (Index row = 0; row < lhs.rows(); row++) { - dst(row, col) = lhs.coeffDiag(row) * rhs(row, col); - } - } - - // Use matrix upper triangular part - for (Index row = 0; row < lhs.rows(); row++) { - typename Lhs_::InnerUpperIterator uIt(lhs, row); - const Index stop = uIt.col() + uIt.size(); - for (Index col = 0; col < rhs.cols(); col++) { - Index k = uIt.col(); - Scalar tmp = 0; - while (k < stop) { - tmp += uIt.value() * rhs(k++, col); - ++uIt; - } - - dst(row, col) += tmp; - uIt += -uIt.size(); - } - } - - // Use matrix lower triangular part - for (Index lhscol = 0; lhscol < lhs.cols(); lhscol++) { - typename Lhs_::InnerLowerIterator lIt(lhs, lhscol); - const Index stop = lIt.size() + lIt.row(); - for (Index rhscol = 0; rhscol < rhs.cols(); rhscol++) { - const Scalar rhsCoeff = rhs.coeff(lhscol, rhscol); - Index k = lIt.row(); - while (k < stop) { - dst(k++, rhscol) += lIt.value() * rhsCoeff; - ++lIt; - } - lIt += -lIt.size(); - } - } -} - -template <typename Lhs, typename Rhs, typename ResultType, int LhsStorageOrder = traits<Lhs>::Flags & RowMajorBit> -struct skyline_product_selector; - -template <typename Lhs, typename Rhs, typename ResultType> -struct skyline_product_selector<Lhs, Rhs, ResultType, RowMajor> { - typedef typename traits<remove_all_t<Lhs>>::Scalar Scalar; - - static void run(const Lhs& lhs, const Rhs& rhs, ResultType& res) { - skyline_row_major_time_dense_product<Lhs, Rhs, ResultType>(lhs, rhs, res); - } -}; - -template <typename Lhs, typename Rhs, typename ResultType> -struct skyline_product_selector<Lhs, Rhs, ResultType, ColMajor> { - typedef typename traits<remove_all_t<Lhs>>::Scalar Scalar; - - static void run(const Lhs& lhs, const Rhs& rhs, ResultType& res) { - skyline_col_major_time_dense_product<Lhs, Rhs, ResultType>(lhs, rhs, res); - } -}; - -} // end namespace internal - -// template<typename Derived> -// template<typename Lhs, typename Rhs > -// Derived & MatrixBase<Derived>::lazyAssign(const SkylineProduct<Lhs, Rhs, SkylineTimeDenseProduct>& product) { -// typedef internal::remove_all_t<Lhs> Lhs_; -// internal::skyline_product_selector<internal::remove_all_t<Lhs>, -// internal::remove_all_t<Rhs>, -// Derived>::run(product.lhs(), product.rhs(), derived()); -// -// return derived(); -// } - -// skyline * dense - -template <typename Derived> -template <typename OtherDerived> -EIGEN_STRONG_INLINE const typename SkylineProductReturnType<Derived, OtherDerived>::Type -SkylineMatrixBase<Derived>::operator*(const MatrixBase<OtherDerived>& other) const { - return typename SkylineProductReturnType<Derived, OtherDerived>::Type(derived(), other.derived()); -} - -} // end namespace Eigen - -#endif // EIGEN_SKYLINEPRODUCT_H
diff --git a/unsupported/Eigen/src/Skyline/SkylineStorage.h b/unsupported/Eigen/src/Skyline/SkylineStorage.h deleted file mode 100644 index 893fb4f..0000000 --- a/unsupported/Eigen/src/Skyline/SkylineStorage.h +++ /dev/null
@@ -1,224 +0,0 @@ -// This file is part of Eigen, a lightweight C++ template library -// for linear algebra. -// -// Copyright (C) 2008-2009 Guillaume Saupin <guillaume.saupin@cea.fr> -// -// 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_SKYLINE_STORAGE_H -#define EIGEN_SKYLINE_STORAGE_H - -// IWYU pragma: private -#include "./InternalHeaderCheck.h" - -namespace Eigen { - -/** Stores a skyline set of values in three structures : - * The diagonal elements - * The upper elements - * The lower elements - * - */ -template <typename Scalar> -class SkylineStorage { - typedef typename NumTraits<Scalar>::Real RealScalar; - typedef SparseIndex Index; - - public: - SkylineStorage() - : m_diag(0), - m_lower(0), - m_upper(0), - m_lowerProfile(0), - m_upperProfile(0), - m_diagSize(0), - m_upperSize(0), - m_lowerSize(0), - m_upperProfileSize(0), - m_lowerProfileSize(0), - m_allocatedSize(0) {} - - SkylineStorage(const SkylineStorage& other) - : m_diag(0), - m_lower(0), - m_upper(0), - m_lowerProfile(0), - m_upperProfile(0), - m_diagSize(0), - m_upperSize(0), - m_lowerSize(0), - m_upperProfileSize(0), - m_lowerProfileSize(0), - m_allocatedSize(0) { - *this = other; - } - - SkylineStorage& operator=(const SkylineStorage& other) { - resize(other.diagSize(), other.m_upperProfileSize, other.m_lowerProfileSize, other.upperSize(), other.lowerSize()); - memcpy(m_diag, other.m_diag, m_diagSize * sizeof(Scalar)); - memcpy(m_upper, other.m_upper, other.upperSize() * sizeof(Scalar)); - memcpy(m_lower, other.m_lower, other.lowerSize() * sizeof(Scalar)); - memcpy(m_upperProfile, other.m_upperProfile, m_upperProfileSize * sizeof(Index)); - memcpy(m_lowerProfile, other.m_lowerProfile, m_lowerProfileSize * sizeof(Index)); - return *this; - } - - void swap(SkylineStorage& other) { - std::swap(m_diag, other.m_diag); - std::swap(m_upper, other.m_upper); - std::swap(m_lower, other.m_lower); - std::swap(m_upperProfile, other.m_upperProfile); - std::swap(m_lowerProfile, other.m_lowerProfile); - std::swap(m_diagSize, other.m_diagSize); - std::swap(m_upperSize, other.m_upperSize); - std::swap(m_lowerSize, other.m_lowerSize); - std::swap(m_allocatedSize, other.m_allocatedSize); - } - - ~SkylineStorage() { - delete[] m_diag; - delete[] m_upper; - if (m_upper != m_lower) delete[] m_lower; - delete[] m_upperProfile; - delete[] m_lowerProfile; - } - - void reserve(Index size, Index upperProfileSize, Index lowerProfileSize, Index upperSize, Index lowerSize) { - Index newAllocatedSize = size + upperSize + lowerSize; - if (newAllocatedSize > m_allocatedSize) reallocate(size, upperProfileSize, lowerProfileSize, upperSize, lowerSize); - } - - void squeeze() { - if (m_allocatedSize > m_diagSize + m_upperSize + m_lowerSize) - reallocate(m_diagSize, m_upperProfileSize, m_lowerProfileSize, m_upperSize, m_lowerSize); - } - - void resize(Index diagSize, Index upperProfileSize, Index lowerProfileSize, Index upperSize, Index lowerSize, - float reserveSizeFactor = 0) { - if (m_allocatedSize < diagSize + upperSize + lowerSize) - reallocate(diagSize, upperProfileSize, lowerProfileSize, upperSize + Index(reserveSizeFactor * upperSize), - lowerSize + Index(reserveSizeFactor * lowerSize)); - m_diagSize = diagSize; - m_upperSize = upperSize; - m_lowerSize = lowerSize; - m_upperProfileSize = upperProfileSize; - m_lowerProfileSize = lowerProfileSize; - } - - inline Index diagSize() const { return m_diagSize; } - - inline Index upperSize() const { return m_upperSize; } - - inline Index lowerSize() const { return m_lowerSize; } - - inline Index upperProfileSize() const { return m_upperProfileSize; } - - inline Index lowerProfileSize() const { return m_lowerProfileSize; } - - inline Index allocatedSize() const { return m_allocatedSize; } - - inline void clear() { m_diagSize = 0; } - - inline Scalar& diag(Index i) { return m_diag[i]; } - - inline const Scalar& diag(Index i) const { return m_diag[i]; } - - inline Scalar& upper(Index i) { return m_upper[i]; } - - inline const Scalar& upper(Index i) const { return m_upper[i]; } - - inline Scalar& lower(Index i) { return m_lower[i]; } - - inline const Scalar& lower(Index i) const { return m_lower[i]; } - - inline Index& upperProfile(Index i) { return m_upperProfile[i]; } - - inline const Index& upperProfile(Index i) const { return m_upperProfile[i]; } - - inline Index& lowerProfile(Index i) { return m_lowerProfile[i]; } - - inline const Index& lowerProfile(Index i) const { return m_lowerProfile[i]; } - - static SkylineStorage Map(Index* upperProfile, Index* lowerProfile, Scalar* diag, Scalar* upper, Scalar* lower, - Index size, Index upperSize, Index lowerSize) { - SkylineStorage res; - res.m_upperProfile = upperProfile; - res.m_lowerProfile = lowerProfile; - res.m_diag = diag; - res.m_upper = upper; - res.m_lower = lower; - res.m_allocatedSize = res.m_diagSize = size; - res.m_upperSize = upperSize; - res.m_lowerSize = lowerSize; - return res; - } - - inline void reset() { - std::fill_n(m_diag, m_diagSize, Scalar(0)); - std::fill_n(m_upper, m_upperSize, Scalar(0)); - std::fill_n(m_lower, m_lowerSize, Scalar(0)); - std::fill_n(m_upperProfile, m_diagSize, Index(0)); - std::fill_n(m_lowerProfile, m_diagSize, Index(0)); - } - - void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar>()) { - // TODO - } - - protected: - inline void reallocate(Index diagSize, Index upperProfileSize, Index lowerProfileSize, Index upperSize, - Index lowerSize) { - Scalar* diag = new Scalar[diagSize]; - Scalar* upper = new Scalar[upperSize]; - Scalar* lower = new Scalar[lowerSize]; - Index* upperProfile = new Index[upperProfileSize]; - Index* lowerProfile = new Index[lowerProfileSize]; - - Index copyDiagSize = (std::min)(diagSize, m_diagSize); - Index copyUpperSize = (std::min)(upperSize, m_upperSize); - Index copyLowerSize = (std::min)(lowerSize, m_lowerSize); - Index copyUpperProfileSize = (std::min)(upperProfileSize, m_upperProfileSize); - Index copyLowerProfileSize = (std::min)(lowerProfileSize, m_lowerProfileSize); - - // copy - memcpy(diag, m_diag, copyDiagSize * sizeof(Scalar)); - memcpy(upper, m_upper, copyUpperSize * sizeof(Scalar)); - memcpy(lower, m_lower, copyLowerSize * sizeof(Scalar)); - memcpy(upperProfile, m_upperProfile, copyUpperProfileSize * sizeof(Index)); - memcpy(lowerProfile, m_lowerProfile, copyLowerProfileSize * sizeof(Index)); - - // delete old stuff - delete[] m_diag; - delete[] m_upper; - delete[] m_lower; - delete[] m_upperProfile; - delete[] m_lowerProfile; - m_diag = diag; - m_upper = upper; - m_lower = lower; - m_upperProfile = upperProfile; - m_lowerProfile = lowerProfile; - m_allocatedSize = diagSize + upperSize + lowerSize; - m_upperSize = upperSize; - m_lowerSize = lowerSize; - } - - public: - Scalar* m_diag; - Scalar* m_upper; - Scalar* m_lower; - Index* m_upperProfile; - Index* m_lowerProfile; - Index m_diagSize; - Index m_upperSize; - Index m_lowerSize; - Index m_upperProfileSize; - Index m_lowerProfileSize; - Index m_allocatedSize; -}; - -} // end namespace Eigen - -#endif // EIGEN_SKYLINE_STORAGE_H
diff --git a/unsupported/Eigen/src/Skyline/SkylineUtil.h b/unsupported/Eigen/src/Skyline/SkylineUtil.h deleted file mode 100644 index e9c828f..0000000 --- a/unsupported/Eigen/src/Skyline/SkylineUtil.h +++ /dev/null
@@ -1,94 +0,0 @@ -// This file is part of Eigen, a lightweight C++ template library -// for linear algebra. -// -// Copyright (C) 2009 Guillaume Saupin <guillaume.saupin@cea.fr> -// -// 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_SKYLINEUTIL_H -#define EIGEN_SKYLINEUTIL_H - -// IWYU pragma: private -#include "./InternalHeaderCheck.h" - -namespace Eigen { - -#ifdef NDEBUG -#define EIGEN_DBG_SKYLINE(X) -#else -#define EIGEN_DBG_SKYLINE(X) X -#endif - -const unsigned int SkylineBit = 0x1200; -template <typename Lhs, typename Rhs, int ProductMode> -class SkylineProduct; -enum AdditionalProductEvaluationMode { SkylineTimeDenseProduct, SkylineTimeSkylineProduct, DenseTimeSkylineProduct }; -enum { IsSkyline = SkylineBit }; - -#define EIGEN_SKYLINE_INHERIT_ASSIGNMENT_OPERATOR(Derived, Op) \ - template <typename OtherDerived> \ - EIGEN_STRONG_INLINE Derived& operator Op(const Eigen::SkylineMatrixBase<OtherDerived>& other) { \ - return Base::operator Op(other.derived()); \ - } \ - EIGEN_STRONG_INLINE Derived& operator Op(const Derived & other) { return Base::operator Op(other); } - -#define EIGEN_SKYLINE_INHERIT_SCALAR_ASSIGNMENT_OPERATOR(Derived, Op) \ - template <typename Other> \ - EIGEN_STRONG_INLINE Derived& operator Op(const Other & scalar) { \ - return Base::operator Op(scalar); \ - } - -#define EIGEN_SKYLINE_INHERIT_ASSIGNMENT_OPERATORS(Derived) \ - EIGEN_SKYLINE_INHERIT_ASSIGNMENT_OPERATOR(Derived, =) \ - EIGEN_SKYLINE_INHERIT_ASSIGNMENT_OPERATOR(Derived, +=) \ - EIGEN_SKYLINE_INHERIT_ASSIGNMENT_OPERATOR(Derived, -=) \ - EIGEN_SKYLINE_INHERIT_SCALAR_ASSIGNMENT_OPERATOR(Derived, *=) \ - EIGEN_SKYLINE_INHERIT_SCALAR_ASSIGNMENT_OPERATOR(Derived, /=) - -#define EIGEN_SKYLINE_GENERIC_PUBLIC_INTERFACE_(Derived, BaseClass) \ - typedef BaseClass Base; \ - typedef typename Eigen::internal::traits<Derived>::Scalar Scalar; \ - typedef typename Eigen::NumTraits<Scalar>::Real RealScalar; \ - typedef typename Eigen::internal::traits<Derived>::StorageKind StorageKind; \ - typedef typename Eigen::internal::index<StorageKind>::type Index; \ - enum { \ - Flags = Eigen::internal::traits<Derived>::Flags, \ - }; - -#define EIGEN_SKYLINE_GENERIC_PUBLIC_INTERFACE(Derived) \ - EIGEN_SKYLINE_GENERIC_PUBLIC_INTERFACE_(Derived, Eigen::SkylineMatrixBase<Derived>) - -template <typename Derived> -class SkylineMatrixBase; -template <typename Scalar_, int Flags_ = 0> -class SkylineMatrix; -template <typename Scalar_, int Flags_ = 0> -class DynamicSkylineMatrix; -template <typename Scalar_, int Flags_ = 0> -class SkylineVector; -template <typename Scalar_, int Flags_ = 0> -class MappedSkylineMatrix; - -namespace internal { - -template <typename Lhs, typename Rhs> -struct skyline_product_mode; -template <typename Lhs, typename Rhs, int ProductMode = skyline_product_mode<Lhs, Rhs>::value> -struct SkylineProductReturnType; - -template <typename T> -class eval<T, IsSkyline> { - typedef typename traits<T>::Scalar Scalar_; - enum { Flags_ = traits<T>::Flags }; - - public: - typedef SkylineMatrix<Scalar_, Flags_> type; -}; - -} // end namespace internal - -} // end namespace Eigen - -#endif // EIGEN_SKYLINEUTIL_H