aboutsummaryrefslogtreecommitdiffhomepage
path: root/Eigen
diff options
context:
space:
mode:
authorGravatar Benoit Jacob <jacob.benoit.1@gmail.com>2010-04-22 14:11:18 -0400
committerGravatar Benoit Jacob <jacob.benoit.1@gmail.com>2010-04-22 14:11:18 -0400
commit9962c59b56960569c8df332144190e62c1eb3b01 (patch)
treea3efa574460c6a08f4ed17a3896b497d5bfc374f /Eigen
parent28dde19e40a3d758faa94f0fc228857f7b3192ea (diff)
* implement the corner() API change: new methods topLeftCorner() etc
* get rid of BlockReturnType: it was not needed, and code was not always using it consistently anyway * add topRows(), leftCols(), bottomRows(), rightCols() * add corners unit-test covering all of that * adapt docs, expand "porting from eigen 2 to 3" * adapt Eigen2Support
Diffstat (limited to 'Eigen')
-rw-r--r--Eigen/src/Core/Block.h604
-rw-r--r--Eigen/src/Core/DenseBase.h74
-rw-r--r--Eigen/src/Core/DenseStorageBase.h4
-rw-r--r--Eigen/src/Core/util/Constants.h4
-rw-r--r--Eigen/src/Core/util/XprHelper.h7
-rw-r--r--Eigen/src/Eigenvalues/HessenbergDecomposition.h8
-rw-r--r--Eigen/src/Eigenvalues/RealSchur.h2
-rw-r--r--Eigen/src/Eigenvalues/Tridiagonalization.h12
-rw-r--r--Eigen/src/Geometry/Transform.h13
-rw-r--r--Eigen/src/Geometry/Umeyama.h2
-rw-r--r--Eigen/src/Householder/HouseholderSequence.h8
-rw-r--r--Eigen/src/LU/FullPivLU.h24
-rw-r--r--Eigen/src/LU/PartialPivLU.h2
-rw-r--r--Eigen/src/QR/ColPivHouseholderQR.h14
-rw-r--r--Eigen/src/QR/FullPivHouseholderQR.h14
-rw-r--r--Eigen/src/QR/HouseholderQR.h12
-rw-r--r--Eigen/src/SVD/UpperBidiagonalization.h4
-rw-r--r--Eigen/src/Sparse/SparseMatrixBase.h8
18 files changed, 647 insertions, 169 deletions
diff --git a/Eigen/src/Core/Block.h b/Eigen/src/Core/Block.h
index d0e6d4ee8..32f4b0ab2 100644
--- a/Eigen/src/Core/Block.h
+++ b/Eigen/src/Core/Block.h
@@ -328,128 +328,456 @@ class Block<XprType,BlockRows,BlockCols,HasDirectAccess>
* \sa class Block, block(int,int)
*/
template<typename Derived>
-inline typename BlockReturnType<Derived>::Type DenseBase<Derived>
+inline Block<Derived> DenseBase<Derived>
::block(int startRow, int startCol, int blockRows, int blockCols)
{
- return typename BlockReturnType<Derived>::Type(derived(), startRow, startCol, blockRows, blockCols);
+ return Block<Derived>(derived(), startRow, startCol, blockRows, blockCols);
}
/** This is the const version of block(int,int,int,int). */
template<typename Derived>
-inline const typename BlockReturnType<Derived>::Type DenseBase<Derived>
+inline const Block<Derived> DenseBase<Derived>
::block(int startRow, int startCol, int blockRows, int blockCols) const
{
- return typename BlockReturnType<Derived>::Type(derived(), startRow, startCol, blockRows, blockCols);
+ return Block<Derived>(derived(), startRow, startCol, blockRows, blockCols);
}
-/** \returns a dynamic-size expression of a corner of *this.
+
+
+
+/** \returns a dynamic-size expression of a top-right corner of *this.
*
- * \param type the type of corner. Can be \a Eigen::TopLeft, \a Eigen::TopRight,
- * \a Eigen::BottomLeft, \a Eigen::BottomRight.
* \param cRows the number of rows in the corner
* \param cCols the number of columns in the corner
*
- * Example: \include MatrixBase_corner_enum_int_int.cpp
- * Output: \verbinclude MatrixBase_corner_enum_int_int.out
+ * Example: \include MatrixBase_topRightCorner_int_int.cpp
+ * Output: \verbinclude MatrixBase_topRightCorner_int_int.out
*
- * \note Even though the returned expression has dynamic size, in the case
- * when it is applied to a fixed-size matrix, it inherits a fixed maximal size,
- * which means that evaluating it does not cause a dynamic memory allocation.
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+inline Block<Derived> DenseBase<Derived>
+ ::topRightCorner(int cRows, int cCols)
+{
+ return Block<Derived>(derived(), 0, cols() - cCols, cRows, cCols);
+}
+
+/** This is the const version of topRightCorner(int, int).*/
+template<typename Derived>
+inline const Block<Derived>
+DenseBase<Derived>::topRightCorner(int cRows, int cCols) const
+{
+ return Block<Derived>(derived(), 0, cols() - cCols, cRows, cCols);
+}
+
+/** \returns an expression of a fixed-size top-right corner of *this.
+ *
+ * The template parameters CRows and CCols are the number of rows and columns in the corner.
+ *
+ * Example: \include MatrixBase_template_int_int_topRightCorner.cpp
+ * Output: \verbinclude MatrixBase_template_int_int_topRightCorner.out
*
* \sa class Block, block(int,int,int,int)
*/
template<typename Derived>
-inline typename BlockReturnType<Derived>::Type DenseBase<Derived>
- ::corner(CornerType type, int cRows, int cCols)
+template<int CRows, int CCols>
+inline Block<Derived, CRows, CCols>
+DenseBase<Derived>::topRightCorner()
{
- switch(type)
- {
- default:
- ei_assert(false && "Bad corner type.");
- case TopLeft:
- return typename BlockReturnType<Derived>::Type(derived(), 0, 0, cRows, cCols);
- case TopRight:
- return typename BlockReturnType<Derived>::Type(derived(), 0, cols() - cCols, cRows, cCols);
- case BottomLeft:
- return typename BlockReturnType<Derived>::Type(derived(), rows() - cRows, 0, cRows, cCols);
- case BottomRight:
- return typename BlockReturnType<Derived>::Type(derived(), rows() - cRows, cols() - cCols, cRows, cCols);
- }
+ return Block<Derived, CRows, CCols>(derived(), 0, cols() - CCols);
}
-/** This is the const version of corner(CornerType, int, int).*/
+/** This is the const version of topRightCorner<int, int>().*/
template<typename Derived>
-inline const typename BlockReturnType<Derived>::Type
-DenseBase<Derived>::corner(CornerType type, int cRows, int cCols) const
+template<int CRows, int CCols>
+inline const Block<Derived, CRows, CCols>
+DenseBase<Derived>::topRightCorner() const
{
- switch(type)
- {
- default:
- ei_assert(false && "Bad corner type.");
- case TopLeft:
- return typename BlockReturnType<Derived>::Type(derived(), 0, 0, cRows, cCols);
- case TopRight:
- return typename BlockReturnType<Derived>::Type(derived(), 0, cols() - cCols, cRows, cCols);
- case BottomLeft:
- return typename BlockReturnType<Derived>::Type(derived(), rows() - cRows, 0, cRows, cCols);
- case BottomRight:
- return typename BlockReturnType<Derived>::Type(derived(), rows() - cRows, cols() - cCols, cRows, cCols);
- }
+ return Block<Derived, CRows, CCols>(derived(), 0, cols() - CCols);
}
-/** \returns a fixed-size expression of a corner of *this.
+
+
+
+/** \returns a dynamic-size expression of a top-left corner of *this.
*
- * \param type the type of corner. Can be \a Eigen::TopLeft, \a Eigen::TopRight,
- * \a Eigen::BottomLeft, \a Eigen::BottomRight.
+ * \param cRows the number of rows in the corner
+ * \param cCols the number of columns in the corner
*
- * The template parameters CRows and CCols arethe number of rows and columns in the corner.
+ * Example: \include MatrixBase_topLeftCorner_int_int.cpp
+ * Output: \verbinclude MatrixBase_topLeftCorner_int_int.out
*
- * Example: \include MatrixBase_template_int_int_corner_enum.cpp
- * Output: \verbinclude MatrixBase_template_int_int_corner_enum.out
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+inline Block<Derived> DenseBase<Derived>
+ ::topLeftCorner(int cRows, int cCols)
+{
+ return Block<Derived>(derived(), 0, 0, cRows, cCols);
+}
+
+/** This is the const version of topLeftCorner(int, int).*/
+template<typename Derived>
+inline const Block<Derived>
+DenseBase<Derived>::topLeftCorner(int cRows, int cCols) const
+{
+ return Block<Derived>(derived(), 0, 0, cRows, cCols);
+}
+
+/** \returns an expression of a fixed-size top-left corner of *this.
+ *
+ * The template parameters CRows and CCols are the number of rows and columns in the corner.
+ *
+ * Example: \include MatrixBase_template_int_int_topLeftCorner.cpp
+ * Output: \verbinclude MatrixBase_template_int_int_topLeftCorner.out
*
* \sa class Block, block(int,int,int,int)
*/
template<typename Derived>
template<int CRows, int CCols>
-inline typename BlockReturnType<Derived, CRows, CCols>::Type
-DenseBase<Derived>::corner(CornerType type)
+inline Block<Derived, CRows, CCols>
+DenseBase<Derived>::topLeftCorner()
{
- switch(type)
- {
- default:
- ei_assert(false && "Bad corner type.");
- case TopLeft:
- return Block<Derived, CRows, CCols>(derived(), 0, 0);
- case TopRight:
- return Block<Derived, CRows, CCols>(derived(), 0, cols() - CCols);
- case BottomLeft:
- return Block<Derived, CRows, CCols>(derived(), rows() - CRows, 0);
- case BottomRight:
- return Block<Derived, CRows, CCols>(derived(), rows() - CRows, cols() - CCols);
- }
+ return Block<Derived, CRows, CCols>(derived(), 0, 0);
}
-/** This is the const version of corner<int, int>(CornerType).*/
+/** This is the const version of topLeftCorner<int, int>().*/
template<typename Derived>
template<int CRows, int CCols>
-inline const typename BlockReturnType<Derived, CRows, CCols>::Type
-DenseBase<Derived>::corner(CornerType type) const
+inline const Block<Derived, CRows, CCols>
+DenseBase<Derived>::topLeftCorner() const
{
- switch(type)
- {
- default:
- ei_assert(false && "Bad corner type.");
- case TopLeft:
- return Block<Derived, CRows, CCols>(derived(), 0, 0);
- case TopRight:
- return Block<Derived, CRows, CCols>(derived(), 0, cols() - CCols);
- case BottomLeft:
- return Block<Derived, CRows, CCols>(derived(), rows() - CRows, 0);
- case BottomRight:
- return Block<Derived, CRows, CCols>(derived(), rows() - CRows, cols() - CCols);
- }
+ return Block<Derived, CRows, CCols>(derived(), 0, 0);
+}
+
+
+
+
+
+
+/** \returns a dynamic-size expression of a bottom-right corner of *this.
+ *
+ * \param cRows the number of rows in the corner
+ * \param cCols the number of columns in the corner
+ *
+ * Example: \include MatrixBase_bottomRightCorner_int_int.cpp
+ * Output: \verbinclude MatrixBase_bottomRightCorner_int_int.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+inline Block<Derived> DenseBase<Derived>
+ ::bottomRightCorner(int cRows, int cCols)
+{
+ return Block<Derived>(derived(), rows() - cRows, cols() - cCols, cRows, cCols);
+}
+
+/** This is the const version of bottomRightCorner(int, int).*/
+template<typename Derived>
+inline const Block<Derived>
+DenseBase<Derived>::bottomRightCorner(int cRows, int cCols) const
+{
+ return Block<Derived>(derived(), rows() - cRows, cols() - cCols, cRows, cCols);
+}
+
+/** \returns an expression of a fixed-size bottom-right corner of *this.
+ *
+ * The template parameters CRows and CCols are the number of rows and columns in the corner.
+ *
+ * Example: \include MatrixBase_template_int_int_bottomRightCorner.cpp
+ * Output: \verbinclude MatrixBase_template_int_int_bottomRightCorner.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+template<int CRows, int CCols>
+inline Block<Derived, CRows, CCols>
+DenseBase<Derived>::bottomRightCorner()
+{
+ return Block<Derived, CRows, CCols>(derived(), rows() - CRows, cols() - CCols);
+}
+
+/** This is the const version of bottomRightCorner<int, int>().*/
+template<typename Derived>
+template<int CRows, int CCols>
+inline const Block<Derived, CRows, CCols>
+DenseBase<Derived>::bottomRightCorner() const
+{
+ return Block<Derived, CRows, CCols>(derived(), rows() - CRows, cols() - CCols);
+}
+
+
+
+
+/** \returns a dynamic-size expression of a bottom-left corner of *this.
+ *
+ * \param cRows the number of rows in the corner
+ * \param cCols the number of columns in the corner
+ *
+ * Example: \include MatrixBase_bottomLeftCorner_int_int.cpp
+ * Output: \verbinclude MatrixBase_bottomLeftCorner_int_int.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+inline Block<Derived> DenseBase<Derived>
+ ::bottomLeftCorner(int cRows, int cCols)
+{
+ return Block<Derived>(derived(), rows() - cRows, 0, cRows, cCols);
+}
+
+/** This is the const version of bottomLeftCorner(int, int).*/
+template<typename Derived>
+inline const Block<Derived>
+DenseBase<Derived>::bottomLeftCorner(int cRows, int cCols) const
+{
+ return Block<Derived>(derived(), rows() - cRows, 0, cRows, cCols);
+}
+
+/** \returns an expression of a fixed-size bottom-left corner of *this.
+ *
+ * The template parameters CRows and CCols are the number of rows and columns in the corner.
+ *
+ * Example: \include MatrixBase_template_int_int_bottomLeftCorner.cpp
+ * Output: \verbinclude MatrixBase_template_int_int_bottomLeftCorner.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+template<int CRows, int CCols>
+inline Block<Derived, CRows, CCols>
+DenseBase<Derived>::bottomLeftCorner()
+{
+ return Block<Derived, CRows, CCols>(derived(), rows() - CRows, 0);
+}
+
+/** This is the const version of bottomLeftCorner<int, int>().*/
+template<typename Derived>
+template<int CRows, int CCols>
+inline const Block<Derived, CRows, CCols>
+DenseBase<Derived>::bottomLeftCorner() const
+{
+ return Block<Derived, CRows, CCols>(derived(), rows() - CRows, 0);
+}
+
+
+
+/** \returns a block consisting of the top rows of *this.
+ *
+ * \param n the number of rows in the block
+ *
+ * Example: \include MatrixBase_topRows_int.cpp
+ * Output: \verbinclude MatrixBase_topRows_int.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+inline typename DenseBase<Derived>::RowsBlockXpr DenseBase<Derived>
+ ::topRows(int n)
+{
+ return RowsBlockXpr(derived(), 0, 0, n, cols());
+}
+
+/** This is the const version of topRows(int).*/
+template<typename Derived>
+inline const typename DenseBase<Derived>::RowsBlockXpr
+DenseBase<Derived>::topRows(int n) const
+{
+ return RowsBlockXpr(derived(), 0, 0, n, cols());
+}
+
+/** \returns a block consisting of the top rows of *this.
+ *
+ * \param N the number of rows in the block
+ *
+ * Example: \include MatrixBase_template_int_topRows.cpp
+ * Output: \verbinclude MatrixBase_template_int_topRows.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+template<int N>
+inline typename DenseBase<Derived>::template NRowsBlockXpr<N>::Type
+DenseBase<Derived>::topRows()
+{
+ return typename DenseBase<Derived>::template NRowsBlockXpr<N>::Type(derived(), 0, 0, N, cols());
+}
+
+/** This is the const version of topRows<int>().*/
+template<typename Derived>
+template<int N>
+inline const typename DenseBase<Derived>::template NRowsBlockXpr<N>::Type
+DenseBase<Derived>::topRows() const
+{
+ return typename DenseBase<Derived>::template NRowsBlockXpr<N>::Type(derived(), 0, 0, N, cols());
}
+
+
+
+
+/** \returns a block consisting of the bottom rows of *this.
+ *
+ * \param n the number of rows in the block
+ *
+ * Example: \include MatrixBase_bottomRows_int.cpp
+ * Output: \verbinclude MatrixBase_bottomRows_int.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+inline typename DenseBase<Derived>::RowsBlockXpr DenseBase<Derived>
+ ::bottomRows(int n)
+{
+ return RowsBlockXpr(derived(), rows() - n, 0, n, cols());
+}
+
+/** This is the const version of bottomRows(int).*/
+template<typename Derived>
+inline const typename DenseBase<Derived>::RowsBlockXpr
+DenseBase<Derived>::bottomRows(int n) const
+{
+ return RowsBlockXpr(derived(), rows() - n, 0, n, cols());
+}
+
+/** \returns a block consisting of the bottom rows of *this.
+ *
+ * \param N the number of rows in the block
+ *
+ * Example: \include MatrixBase_template_int_bottomRows.cpp
+ * Output: \verbinclude MatrixBase_template_int_bottomRows.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+template<int N>
+inline typename DenseBase<Derived>::template NRowsBlockXpr<N>::Type
+DenseBase<Derived>::bottomRows()
+{
+ return typename NRowsBlockXpr<N>::Type(derived(), rows() - N, 0, N, cols());
+}
+
+/** This is the const version of bottomRows<int>().*/
+template<typename Derived>
+template<int N>
+inline const typename DenseBase<Derived>::template NRowsBlockXpr<N>::Type
+DenseBase<Derived>::bottomRows() const
+{
+ return typename NRowsBlockXpr<N>::Type(derived(), rows() - N, 0, N, cols());
+}
+
+
+
+
+
+/** \returns a block consisting of the top columns of *this.
+ *
+ * \param n the number of columns in the block
+ *
+ * Example: \include MatrixBase_leftCols_int.cpp
+ * Output: \verbinclude MatrixBase_leftCols_int.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+inline typename DenseBase<Derived>::ColsBlockXpr DenseBase<Derived>
+ ::leftCols(int n)
+{
+ return ColsBlockXpr(derived(), 0, 0, rows(), n);
+}
+
+/** This is the const version of leftCols(int).*/
+template<typename Derived>
+inline const typename DenseBase<Derived>::ColsBlockXpr
+DenseBase<Derived>::leftCols(int n) const
+{
+ return ColsBlockXpr(derived(), 0, 0, rows(), n);
+}
+
+/** \returns a block consisting of the top columns of *this.
+ *
+ * \param N the number of columns in the block
+ *
+ * Example: \include MatrixBase_template_int_leftCols.cpp
+ * Output: \verbinclude MatrixBase_template_int_leftCols.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+template<int N>
+inline typename DenseBase<Derived>::template NColsBlockXpr<N>::Type
+DenseBase<Derived>::leftCols()
+{
+ return typename NColsBlockXpr<N>::Type(derived(), 0, 0, rows(), N);
+}
+
+/** This is the const version of leftCols<int>().*/
+template<typename Derived>
+template<int N>
+inline const typename DenseBase<Derived>::template NColsBlockXpr<N>::Type
+DenseBase<Derived>::leftCols() const
+{
+ return typename NColsBlockXpr<N>::Type(derived(), 0, 0, rows(), N);
+}
+
+
+
+
+
+/** \returns a block consisting of the top columns of *this.
+ *
+ * \param n the number of columns in the block
+ *
+ * Example: \include MatrixBase_rightCols_int.cpp
+ * Output: \verbinclude MatrixBase_rightCols_int.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+inline typename DenseBase<Derived>::ColsBlockXpr DenseBase<Derived>
+ ::rightCols(int n)
+{
+ return ColsBlockXpr(derived(), 0, cols() - n, rows(), n);
+}
+
+/** This is the const version of rightCols(int).*/
+template<typename Derived>
+inline const typename DenseBase<Derived>::ColsBlockXpr
+DenseBase<Derived>::rightCols(int n) const
+{
+ return ColsBlockXpr(derived(), 0, cols() - n, rows(), n);
+}
+
+/** \returns a block consisting of the top columns of *this.
+ *
+ * \param N the number of columns in the block
+ *
+ * Example: \include MatrixBase_template_int_rightCols.cpp
+ * Output: \verbinclude MatrixBase_template_int_rightCols.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+template<int N>
+inline typename DenseBase<Derived>::template NColsBlockXpr<N>::Type
+DenseBase<Derived>::rightCols()
+{
+ return typename DenseBase<Derived>::template NColsBlockXpr<N>::Type(derived(), 0, cols() - N, rows(), N);
+}
+
+/** This is the const version of rightCols<int>().*/
+template<typename Derived>
+template<int N>
+inline const typename DenseBase<Derived>::template NColsBlockXpr<N>::Type
+DenseBase<Derived>::rightCols() const
+{
+ return typename DenseBase<Derived>::template NColsBlockXpr<N>::Type(derived(), 0, cols() - N, rows(), N);
+}
+
+
+
+
+
/** \returns a fixed-size expression of a block in *this.
*
* The template parameters \a BlockRows and \a BlockCols are the number of
@@ -468,7 +796,7 @@ DenseBase<Derived>::corner(CornerType type) const
*/
template<typename Derived>
template<int BlockRows, int BlockCols>
-inline typename BlockReturnType<Derived, BlockRows, BlockCols>::Type
+inline Block<Derived, BlockRows, BlockCols>
DenseBase<Derived>::block(int startRow, int startCol)
{
return Block<Derived, BlockRows, BlockCols>(derived(), startRow, startCol);
@@ -477,7 +805,7 @@ DenseBase<Derived>::block(int startRow, int startCol)
/** This is the const version of block<>(int, int). */
template<typename Derived>
template<int BlockRows, int BlockCols>
-inline const typename BlockReturnType<Derived, BlockRows, BlockCols>::Type
+inline const Block<Derived, BlockRows, BlockCols>
DenseBase<Derived>::block(int startRow, int startCol) const
{
return Block<Derived, BlockRows, BlockCols>(derived(), startRow, startCol);
@@ -525,4 +853,116 @@ DenseBase<Derived>::row(int i) const
return RowXpr(derived(), i);
}
+#ifdef EIGEN2_SUPPORT
+
+/** \returns a dynamic-size expression of a corner of *this.
+ *
+ * \param type the type of corner. Can be \a Eigen::TopLeft, \a Eigen::TopRight,
+ * \a Eigen::BottomLeft, \a Eigen::BottomRight.
+ * \param cRows the number of rows in the corner
+ * \param cCols the number of columns in the corner
+ *
+ * Example: \include MatrixBase_corner_enum_int_int.cpp
+ * Output: \verbinclude MatrixBase_corner_enum_int_int.out
+ *
+ * \note Even though the returned expression has dynamic size, in the case
+ * when it is applied to a fixed-size matrix, it inherits a fixed maximal size,
+ * which means that evaluating it does not cause a dynamic memory allocation.
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+inline Block<Derived> DenseBase<Derived>
+ ::corner(CornerType type, int cRows, int cCols)
+{
+ switch(type)
+ {
+ default:
+ ei_assert(false && "Bad corner type.");
+ case TopLeft:
+ return Block<Derived>(derived(), 0, 0, cRows, cCols);
+ case TopRight:
+ return Block<Derived>(derived(), 0, cols() - cCols, cRows, cCols);
+ case BottomLeft:
+ return Block<Derived>(derived(), rows() - cRows, 0, cRows, cCols);
+ case BottomRight:
+ return Block<Derived>(derived(), rows() - cRows, cols() - cCols, cRows, cCols);
+ }
+}
+
+/** This is the const version of corner(CornerType, int, int).*/
+template<typename Derived>
+inline const Block<Derived>
+DenseBase<Derived>::corner(CornerType type, int cRows, int cCols) const
+{
+ switch(type)
+ {
+ default:
+ ei_assert(false && "Bad corner type.");
+ case TopLeft:
+ return Block<Derived>(derived(), 0, 0, cRows, cCols);
+ case TopRight:
+ return Block<Derived>(derived(), 0, cols() - cCols, cRows, cCols);
+ case BottomLeft:
+ return Block<Derived>(derived(), rows() - cRows, 0, cRows, cCols);
+ case BottomRight:
+ return Block<Derived>(derived(), rows() - cRows, cols() - cCols, cRows, cCols);
+ }
+}
+
+/** \returns a fixed-size expression of a corner of *this.
+ *
+ * \param type the type of corner. Can be \a Eigen::TopLeft, \a Eigen::TopRight,
+ * \a Eigen::BottomLeft, \a Eigen::BottomRight.
+ *
+ * The template parameters CRows and CCols arethe number of rows and columns in the corner.
+ *
+ * Example: \include MatrixBase_template_int_int_corner_enum.cpp
+ * Output: \verbinclude MatrixBase_template_int_int_corner_enum.out
+ *
+ * \sa class Block, block(int,int,int,int)
+ */
+template<typename Derived>
+template<int CRows, int CCols>
+inline Block<Derived, CRows, CCols>
+DenseBase<Derived>::corner(CornerType type)
+{
+ switch(type)
+ {
+ default:
+ ei_assert(false && "Bad corner type.");
+ case TopLeft:
+ return Block<Derived, CRows, CCols>(derived(), 0, 0);
+ case TopRight:
+ return Block<Derived, CRows, CCols>(derived(), 0, cols() - CCols);
+ case BottomLeft:
+ return Block<Derived, CRows, CCols>(derived(), rows() - CRows, 0);
+ case BottomRight:
+ return Block<Derived, CRows, CCols>(derived(), rows() - CRows, cols() - CCols);
+ }
+}
+
+/** This is the const version of corner<int, int>(CornerType).*/
+template<typename Derived>
+template<int CRows, int CCols>
+inline const Block<Derived, CRows, CCols>
+DenseBase<Derived>::corner(CornerType type) const
+{
+ switch(type)
+ {
+ default:
+ ei_assert(false && "Bad corner type.");
+ case TopLeft:
+ return Block<Derived, CRows, CCols>(derived(), 0, 0);
+ case TopRight:
+ return Block<Derived, CRows, CCols>(derived(), 0, cols() - CCols);
+ case BottomLeft:
+ return Block<Derived, CRows, CCols>(derived(), rows() - CRows, 0);
+ case BottomRight:
+ return Block<Derived, CRows, CCols>(derived(), rows() - CRows, cols() - CCols);
+ }
+}
+
+#endif // EIGEN2_SUPPORT
+
#endif // EIGEN_BLOCK_H
diff --git a/Eigen/src/Core/DenseBase.h b/Eigen/src/Core/DenseBase.h
index 566b4b410..da464067b 100644
--- a/Eigen/src/Core/DenseBase.h
+++ b/Eigen/src/Core/DenseBase.h
@@ -220,8 +220,18 @@ template<typename Derived> class DenseBase
typedef Matrix<typename NumTraits<typename ei_traits<Derived>::Scalar>::Real, ei_traits<Derived>::ColsAtCompileTime, 1> EigenvaluesReturnType;
/** \internal expression type of a column */
typedef Block<Derived, ei_traits<Derived>::RowsAtCompileTime, 1> ColXpr;
- /** \internal expression type of a column */
+ /** \internal expression type of a row */
typedef Block<Derived, 1, ei_traits<Derived>::ColsAtCompileTime> RowXpr;
+ /** \internal expression type of a block of whole columns */
+ typedef Block<Derived, ei_traits<Derived>::RowsAtCompileTime, Dynamic> ColsBlockXpr;
+ /** \internal expression type of a block of whole rows */
+ typedef Block<Derived, Dynamic, ei_traits<Derived>::ColsAtCompileTime> RowsBlockXpr;
+ /** \internal expression type of a block of whole columns */
+ template<int N> struct NColsBlockXpr { typedef Block<Derived, ei_traits<Derived>::RowsAtCompileTime, N> Type; };
+ /** \internal expression type of a block of whole rows */
+ template<int N> struct NRowsBlockXpr { typedef Block<Derived, N, ei_traits<Derived>::ColsAtCompileTime> Type; };
+
+
#endif // not EIGEN_PARSED_BY_DOXYGEN
const CoeffReturnType x() const;
@@ -334,9 +344,8 @@ template<typename Derived> class DenseBase
ColXpr col(int i);
const ColXpr col(int i) const;
- typename BlockReturnType<Derived>::Type block(int startRow, int startCol, int blockRows, int blockCols);
- const typename BlockReturnType<Derived>::Type
- block(int startRow, int startCol, int blockRows, int blockCols) const;
+ Block<Derived> block(int startRow, int startCol, int blockRows, int blockCols);
+ const Block<Derived> block(int startRow, int startCol, int blockRows, int blockCols) const;
VectorBlock<Derived> segment(int start, int size);
const VectorBlock<Derived> segment(int start, int size) const;
@@ -347,18 +356,46 @@ template<typename Derived> class DenseBase
VectorBlock<Derived> tail(int size);
const VectorBlock<Derived> tail(int size) const;
- typename BlockReturnType<Derived>::Type corner(CornerType type, int cRows, int cCols);
- const typename BlockReturnType<Derived>::Type corner(CornerType type, int cRows, int cCols) const;
+ Block<Derived> topLeftCorner(int cRows, int cCols);
+ const Block<Derived> topLeftCorner(int cRows, int cCols) const;
+ Block<Derived> topRightCorner(int cRows, int cCols);
+ const Block<Derived> topRightCorner(int cRows, int cCols) const;
+ Block<Derived> bottomLeftCorner(int cRows, int cCols);
+ const Block<Derived> bottomLeftCorner(int cRows, int cCols) const;
+ Block<Derived> bottomRightCorner(int cRows, int cCols);
+ const Block<Derived> bottomRightCorner(int cRows, int cCols) const;
+
+ RowsBlockXpr topRows(int n);
+ const RowsBlockXpr topRows(int n) const;
+ RowsBlockXpr bottomRows(int n);
+ const RowsBlockXpr bottomRows(int n) const;
+ ColsBlockXpr leftCols(int n);
+ const ColsBlockXpr leftCols(int n) const;
+ ColsBlockXpr rightCols(int n);
+ const ColsBlockXpr rightCols(int n) const;
+
+ template<int CRows, int CCols> Block<Derived, CRows, CCols> topLeftCorner();
+ template<int CRows, int CCols> const Block<Derived, CRows, CCols> topLeftCorner() const;
+ template<int CRows, int CCols> Block<Derived, CRows, CCols> topRightCorner();
+ template<int CRows, int CCols> const Block<Derived, CRows, CCols> topRightCorner() const;
+ template<int CRows, int CCols> Block<Derived, CRows, CCols> bottomLeftCorner();
+ template<int CRows, int CCols> const Block<Derived, CRows, CCols> bottomLeftCorner() const;
+ template<int CRows, int CCols> Block<Derived, CRows, CCols> bottomRightCorner();
+ template<int CRows, int CCols> const Block<Derived, CRows, CCols> bottomRightCorner() const;
+
+ template<int NRows> typename NRowsBlockXpr<NRows>::Type topRows();
+ template<int NRows> const typename NRowsBlockXpr<NRows>::Type topRows() const;
+ template<int NRows> typename NRowsBlockXpr<NRows>::Type bottomRows();
+ template<int NRows> const typename NRowsBlockXpr<NRows>::Type bottomRows() const;
+ template<int NCols> typename NColsBlockXpr<NCols>::Type leftCols();
+ template<int NCols> const typename NColsBlockXpr<NCols>::Type leftCols() const;
+ template<int NCols> typename NColsBlockXpr<NCols>::Type rightCols();
+ template<int NCols> const typename NColsBlockXpr<NCols>::Type rightCols() const;
template<int BlockRows, int BlockCols>
- typename BlockReturnType<Derived, BlockRows, BlockCols>::Type block(int startRow, int startCol);
+ Block<Derived, BlockRows, BlockCols> block(int startRow, int startCol);
template<int BlockRows, int BlockCols>
- const typename BlockReturnType<Derived, BlockRows, BlockCols>::Type block(int startRow, int startCol) const;
-
- template<int CRows, int CCols>
- typename BlockReturnType<Derived, CRows, CCols>::Type corner(CornerType type);
- template<int CRows, int CCols>
- const typename BlockReturnType<Derived, CRows, CCols>::Type corner(CornerType type) const;
+ const Block<Derived, BlockRows, BlockCols> block(int startRow, int startCol) const;
template<int Size> VectorBlock<Derived,Size> head(void);
template<int Size> const VectorBlock<Derived,Size> head() const;
@@ -523,6 +560,17 @@ template<typename Derived> class DenseBase
const Eigen::Reverse<Derived, BothDirections> reverse() const;
void reverseInPlace();
+#ifdef EIGEN2_SUPPORT
+
+ Block<Derived> corner(CornerType type, int cRows, int cCols);
+ const Block<Derived> corner(CornerType type, int cRows, int cCols) const;
+ template<int CRows, int CCols>
+ Block<Derived, CRows, CCols> corner(CornerType type);
+ template<int CRows, int CCols>
+ const Block<Derived, CRows, CCols> corner(CornerType type) const;
+
+#endif // EIGEN2_SUPPORT
+
#ifdef EIGEN_DENSEBASE_PLUGIN
#include EIGEN_DENSEBASE_PLUGIN
#endif
diff --git a/Eigen/src/Core/DenseStorageBase.h b/Eigen/src/Core/DenseStorageBase.h
index 687c76cf8..220cd86c7 100644
--- a/Eigen/src/Core/DenseStorageBase.h
+++ b/Eigen/src/Core/DenseStorageBase.h
@@ -566,9 +566,9 @@ struct ei_conservative_resize_like_impl
const int new_cols = other.cols() - _this.cols();
_this.derived().m_storage.conservativeResize(other.size(),other.rows(),other.cols());
if (new_rows>0)
- _this.corner(BottomRight, new_rows, other.cols()) = other.corner(BottomRight, new_rows, other.cols());
+ _this.bottomRightCorner(new_rows, other.cols()) = other.bottomRows(new_rows);
else if (new_cols>0)
- _this.corner(BottomRight, other.rows(), new_cols) = other.corner(BottomRight, other.rows(), new_cols);
+ _this.bottomRightCorner(other.rows(), new_cols) = other.rightCols(new_cols);
}
else
{
diff --git a/Eigen/src/Core/util/Constants.h b/Eigen/src/Core/util/Constants.h
index 51c3dd287..023ef3852 100644
--- a/Eigen/src/Core/util/Constants.h
+++ b/Eigen/src/Core/util/Constants.h
@@ -165,7 +165,11 @@ enum {
enum { Unaligned=0, Aligned=1 };
enum { ConditionalJumpCost = 5 };
+
+// FIXME after the corner() API change, this was not needed anymore, except by AlignedBox
+// TODO: find out what to do with that. Adapt the AlignedBox API ?
enum CornerType { TopLeft, TopRight, BottomLeft, BottomRight };
+
enum DirectionType { Vertical, Horizontal, BothDirections };
enum ProductEvaluationMode { NormalProduct, CacheFriendlyProduct };
diff --git a/Eigen/src/Core/util/XprHelper.h b/Eigen/src/Core/util/XprHelper.h
index 4eb132cf7..f895dad1c 100644
--- a/Eigen/src/Core/util/XprHelper.h
+++ b/Eigen/src/Core/util/XprHelper.h
@@ -329,13 +329,6 @@ struct ei_special_scalar_op_base<Derived,Scalar,OtherScalar,true> : public Eige
{ return static_cast<const ei_special_scalar_op_base&>(matrix).operator*(scalar); }
};
-/** \internal Gives the type of a sub-matrix or sub-vector of a matrix of type \a ExpressionType and size \a Size
- * TODO: could be a good idea to define a big ReturnType struct ??
- */
-template<typename ExpressionType, int RowsOrSize=Dynamic, int Cols=Dynamic> struct BlockReturnType {
- typedef Block<ExpressionType, RowsOrSize, Cols> Type;
-};
-
template<typename ExpressionType> struct HNormalizedReturnType {
enum {
diff --git a/Eigen/src/Eigenvalues/HessenbergDecomposition.h b/Eigen/src/Eigenvalues/HessenbergDecomposition.h
index 451604252..a1ac31981 100644
--- a/Eigen/src/Eigenvalues/HessenbergDecomposition.h
+++ b/Eigen/src/Eigenvalues/HessenbergDecomposition.h
@@ -272,11 +272,11 @@ void HessenbergDecomposition<MatrixType>::_compute(MatrixType& matA, CoeffVector
// i.e., compute A = H A H'
// A = H A
- matA.corner(BottomRight, remainingSize, remainingSize)
+ matA.bottomRightCorner(remainingSize, remainingSize)
.applyHouseholderOnTheLeft(matA.col(i).tail(remainingSize-1), h, &temp.coeffRef(0));
// A = A H'
- matA.corner(BottomRight, n, remainingSize)
+ matA.rightCols(remainingSize)
.applyHouseholderOnTheRight(matA.col(i).tail(remainingSize-1).conjugate(), ei_conj(h), &temp.coeffRef(0));
}
}
@@ -290,7 +290,7 @@ HessenbergDecomposition<MatrixType>::matrixQ() const
VectorType temp(n);
for (int i = n-2; i>=0; i--)
{
- matQ.corner(BottomRight,n-i-1,n-i-1)
+ matQ.bottomRightCorner(n-i-1,n-i-1)
.applyHouseholderOnTheLeft(m_matrix.col(i).tail(n-i-2), ei_conj(m_hCoeffs.coeff(i)), &temp.coeffRef(0,0));
}
return matQ;
@@ -307,7 +307,7 @@ HessenbergDecomposition<MatrixType>::matrixH() const
int n = m_matrix.rows();
MatrixType matH = m_matrix;
if (n>2)
- matH.corner(BottomLeft,n-2, n-2).template triangularView<Lower>().setZero();
+ matH.bottomLeftCorner(n-2, n-2).template triangularView<Lower>().setZero();
return matH;
}
diff --git a/Eigen/src/Eigenvalues/RealSchur.h b/Eigen/src/Eigenvalues/RealSchur.h
index dd4a1ab7e..c2835f1fb 100644
--- a/Eigen/src/Eigenvalues/RealSchur.h
+++ b/Eigen/src/Eigenvalues/RealSchur.h
@@ -249,7 +249,7 @@ inline typename MatrixType::Scalar RealSchur<MatrixType>::computeNormOfT()
const int size = m_matU.cols();
// FIXME to be efficient the following would requires a triangular reduxion code
// Scalar norm = m_matT.upper().cwiseAbs().sum()
- // + m_matT.corner(BottomLeft,size-1,size-1).diagonal().cwiseAbs().sum();
+ // + m_matT.bottomLeftCorner(size-1,size-1).diagonal().cwiseAbs().sum();
Scalar norm = 0.0;
for (int j = 0; j < size; ++j)
norm += m_matT.row(j).segment(std::max(j-1,0), size-std::max(j-1,0)).cwiseAbs().sum();
diff --git a/Eigen/src/Eigenvalues/Tridiagonalization.h b/Eigen/src/Eigenvalues/Tridiagonalization.h
index c320242ef..a41a5f670 100644
--- a/Eigen/src/Eigenvalues/Tridiagonalization.h
+++ b/Eigen/src/Eigenvalues/Tridiagonalization.h
@@ -171,11 +171,11 @@ Tridiagonalization<MatrixType>::matrixT(void) const
// and fill it ? (to avoid temporaries)
int n = m_matrix.rows();
MatrixType matT = m_matrix;
- matT.corner(TopRight,n-1, n-1).diagonal() = subDiagonal().template cast<Scalar>().conjugate();
+ matT.topRightCorner(n-1, n-1).diagonal() = subDiagonal().template cast<Scalar>().conjugate();
if (n>2)
{
- matT.corner(TopRight,n-2, n-2).template triangularView<Upper>().setZero();
- matT.corner(BottomLeft,n-2, n-2).template triangularView<Lower>().setZero();
+ matT.topRightCorner(n-2, n-2).template triangularView<Upper>().setZero();
+ matT.bottomLeftCorner(n-2, n-2).template triangularView<Lower>().setZero();
}
return matT;
}
@@ -210,12 +210,12 @@ void Tridiagonalization<MatrixType>::_compute(MatrixType& matA, CoeffVectorType&
// i.e., A = H A H' where H = I - h v v' and v = matA.col(i).tail(n-i-1)
matA.col(i).coeffRef(i+1) = 1;
- hCoeffs.tail(n-i-1).noalias() = (matA.corner(BottomRight,remainingSize,remainingSize).template selfadjointView<Lower>()
+ hCoeffs.tail(n-i-1).noalias() = (matA.bottomRightCorner(remainingSize,remainingSize).template selfadjointView<Lower>()
* (ei_conj(h) * matA.col(i).tail(remainingSize)));
hCoeffs.tail(n-i-1) += (ei_conj(h)*Scalar(-0.5)*(hCoeffs.tail(remainingSize).dot(matA.col(i).tail(remainingSize)))) * matA.col(i).tail(n-i-1);
- matA.corner(BottomRight, remainingSize, remainingSize).template selfadjointView<Lower>()
+ matA.bottomRightCorner(remainingSize, remainingSize).template selfadjointView<Lower>()
.rankUpdate(matA.col(i).tail(remainingSize), hCoeffs.tail(remainingSize), -1);
matA.col(i).coeffRef(i+1) = beta;
@@ -243,7 +243,7 @@ void Tridiagonalization<MatrixType>::matrixQInPlace(MatrixBase<QDerived>* q) con
RowVectorType aux(n);
for (int i = n-2; i>=0; i--)
{
- matQ.corner(BottomRight,n-i-1,n-i-1)
+ matQ.bottomRightCorner(n-i-1,n-i-1)
.applyHouseholderOnTheLeft(m_matrix.col(i).tail(n-i-2), ei_conj(m_hCoeffs.coeff(i)), &aux.coeffRef(0,0));
}
}
diff --git a/Eigen/src/Geometry/Transform.h b/Eigen/src/Geometry/Transform.h
index 1240a95bc..1f4db8098 100644
--- a/Eigen/src/Geometry/Transform.h
+++ b/Eigen/src/Geometry/Transform.h
@@ -653,7 +653,7 @@ Transform<Scalar,Dim,Mode>::prescale(const MatrixBase<OtherDerived> &other)
template<typename Scalar, int Dim, int Mode>
inline Transform<Scalar,Dim,Mode>& Transform<Scalar,Dim,Mode>::prescale(Scalar s)
{
- m_matrix.template corner<Dim,HDim>(TopLeft) *= s;
+ m_matrix.template topRows<Dim>() *= s;
return *this;
}
@@ -940,18 +940,19 @@ Transform<Scalar,Dim,Mode>::inverse(TransformTraits hint) const
{
if (hint == Isometry)
{
- res.matrix().template corner<Dim,Dim>(TopLeft) = linear().transpose();
+ res.matrix().template topLeftCorner<Dim,Dim>() = linear().transpose();
}
else if(hint&Affine)
{
- res.matrix().template corner<Dim,Dim>(TopLeft) = linear().inverse();
+ res.matrix().template topLeftCorner<Dim,Dim>() = linear().inverse();
}
else
{
ei_assert(false && "Invalid transform traits in Transform::Inverse");
}
// translation and remaining parts
- res.matrix().template corner<Dim,1>(TopRight) = - res.matrix().template corner<Dim,Dim>(TopLeft) * translation();
+ res.matrix().template topRightCorner<Dim,1>()
+ = - res.matrix().template topLeftCorner<Dim,Dim>() * translation();
if(int(Mode)!=int(AffineCompact))
{
res.matrix().template block<1,Dim>(Dim,0).setZero();
@@ -1247,8 +1248,8 @@ struct ei_transform_left_product_impl<Other,Mode,Dim,HDim, Dim,Dim>
TransformType res;
if(Mode!=AffineCompact)
res.matrix().row(Dim) = tr.matrix().row(Dim);
- res.matrix().template corner<Dim,HDim>(TopLeft).noalias()
- = other * tr.matrix().template corner<Dim,HDim>(TopLeft);
+ res.matrix().template topRows<Dim>().noalias()
+ = other * tr.matrix().template topRows<Dim>();
return res;
}
};
diff --git a/Eigen/src/Geometry/Umeyama.h b/Eigen/src/Geometry/Umeyama.h
index c5d99d533..262d27aa3 100644
--- a/Eigen/src/Geometry/Umeyama.h
+++ b/Eigen/src/Geometry/Umeyama.h
@@ -171,7 +171,7 @@ umeyama(const MatrixBase<Derived>& src, const MatrixBase<OtherDerived>& dst, boo
// Note that we first assign dst_mean to the destination so that there no need
// for a temporary.
Rt.col(m).head(m) = dst_mean;
- Rt.col(m).head(m).noalias() -= c*Rt.corner(TopLeft,m,m)*src_mean;
+ Rt.col(m).head(m).noalias() -= c*Rt.topLeftCorner(m,m)*src_mean;
if (with_scaling) Rt.block(0,0,m,m) *= c;
diff --git a/Eigen/src/Householder/HouseholderSequence.h b/Eigen/src/Householder/HouseholderSequence.h
index 05e883168..4bac484eb 100644
--- a/Eigen/src/Householder/HouseholderSequence.h
+++ b/Eigen/src/Householder/HouseholderSequence.h
@@ -161,10 +161,10 @@ template<typename VectorsType, typename CoeffsType, int Side> class HouseholderS
{
int cornerSize = rows() - k - m_shift;
if(m_trans)
- dst.corner(BottomRight, cornerSize, cornerSize)
+ dst.bottomRightCorner(cornerSize, cornerSize)
.applyHouseholderOnTheRight(essentialVector(k), m_coeffs.coeff(k), &temp.coeffRef(0));
else
- dst.corner(BottomRight, cornerSize, cornerSize)
+ dst.bottomRightCorner(cornerSize, cornerSize)
.applyHouseholderOnTheLeft(essentialVector(k), m_coeffs.coeff(k), &temp.coeffRef(0));
}
}
@@ -176,7 +176,7 @@ template<typename VectorsType, typename CoeffsType, int Side> class HouseholderS
for(int k = 0; k < m_actualVectors; ++k)
{
int actual_k = m_trans ? m_actualVectors-k-1 : k;
- dst.corner(BottomRight, dst.rows(), rows()-m_shift-actual_k)
+ dst.rightCols(rows()-m_shift-actual_k)
.applyHouseholderOnTheRight(essentialVector(actual_k), m_coeffs.coeff(actual_k), &temp.coeffRef(0));
}
}
@@ -188,7 +188,7 @@ template<typename VectorsType, typename CoeffsType, int Side> class HouseholderS
for(int k = 0; k < m_actualVectors; ++k)
{
int actual_k = m_trans ? k : m_actualVectors-k-1;
- dst.corner(BottomRight, rows()-m_shift-actual_k, dst.cols())
+ dst.bottomRows(rows()-m_shift-actual_k)
.applyHouseholderOnTheLeft(essentialVector(actual_k), m_coeffs.coeff(actual_k), &temp.coeffRef(0));
}
}
diff --git a/Eigen/src/LU/FullPivLU.h b/Eigen/src/LU/FullPivLU.h
index 241fb9d16..1b0e67b3a 100644
--- a/Eigen/src/LU/FullPivLU.h
+++ b/Eigen/src/LU/FullPivLU.h
@@ -450,7 +450,7 @@ FullPivLU<MatrixType>& FullPivLU<MatrixType>::compute(const MatrixType& matrix)
// biggest coefficient in the remaining bottom-right corner (starting at row k, col k)
int row_of_biggest_in_corner, col_of_biggest_in_corner;
RealScalar biggest_in_corner;
- biggest_in_corner = m_lu.corner(Eigen::BottomRight, rows-k, cols-k)
+ biggest_in_corner = m_lu.bottomRightCorner(rows-k, cols-k)
.cwiseAbs()
.maxCoeff(&row_of_biggest_in_corner, &col_of_biggest_in_corner);
row_of_biggest_in_corner += k; // correct the values! since they were computed in the corner,
@@ -535,9 +535,9 @@ MatrixType FullPivLU<MatrixType>::reconstructedMatrix() const
// LU
MatrixType res(m_lu.rows(),m_lu.cols());
// FIXME the .toDenseMatrix() should not be needed...
- res = m_lu.corner(TopLeft,m_lu.rows(),smalldim)
+ res = m_lu.leftCols(smalldim)
.template triangularView<UnitLower>().toDenseMatrix()
- * m_lu.corner(TopLeft,smalldim,m_lu.cols())
+ * m_lu.topRows(smalldim)
.template triangularView<Upper>().toDenseMatrix();
// P^{-1}(LU)
@@ -618,9 +618,9 @@ struct ei_kernel_retval<FullPivLU<_MatrixType> >
// ok, we have our trapezoid matrix, we can apply the triangular solver.
// notice that the math behind this suggests that we should apply this to the
// negative of the RHS, but for performance we just put the negative sign elsewhere, see below.
- m.corner(TopLeft, rank(), rank())
+ m.topLeftCorner(rank(), rank())
.template triangularView<Upper>().solveInPlace(
- m.corner(TopRight, rank(), dimker)
+ m.topRightCorner(rank(), dimker)
);
// now we must undo the column permutation that we had applied!
@@ -707,21 +707,21 @@ struct ei_solve_retval<FullPivLU<_MatrixType>, Rhs>
// Step 2
dec().matrixLU()
- .corner(Eigen::TopLeft,smalldim,smalldim)
+ .topLeftCorner(smalldim,smalldim)
.template triangularView<UnitLower>()
- .solveInPlace(c.corner(Eigen::TopLeft, smalldim, c.cols()));
+ .solveInPlace(c.topRows(smalldim));
if(rows>cols)
{
- c.corner(Eigen::BottomLeft, rows-cols, c.cols())
- -= dec().matrixLU().corner(Eigen::BottomLeft, rows-cols, cols)
- * c.corner(Eigen::TopLeft, cols, c.cols());
+ c.bottomRows(rows-cols)
+ -= dec().matrixLU().bottomRows(rows-cols)
+ * c.topRows(cols);
}
// Step 3
dec().matrixLU()
- .corner(TopLeft, nonzero_pivots, nonzero_pivots)
+ .topLeftCorner(nonzero_pivots, nonzero_pivots)
.template triangularView<Upper>()
- .solveInPlace(c.corner(TopLeft, nonzero_pivots, c.cols()));
+ .solveInPlace(c.topRows(nonzero_pivots));
// Step 4
for(int i = 0; i < nonzero_pivots; ++i)
diff --git a/Eigen/src/LU/PartialPivLU.h b/Eigen/src/LU/PartialPivLU.h
index 42c74ee53..695b7d75c 100644
--- a/Eigen/src/LU/PartialPivLU.h
+++ b/Eigen/src/LU/PartialPivLU.h
@@ -283,7 +283,7 @@ struct ei_partial_lu_impl
int rrows = rows-k-1;
int rsize = size-k-1;
lu.col(k).tail(rrows) /= lu.coeff(k,k);
- lu.corner(BottomRight,rrows,rsize).noalias() -= lu.col(k).tail(rrows) * lu.row(k).tail(rsize);
+ lu.bottomRightCorner(rrows,rsize).noalias() -= lu.col(k).tail(rrows) * lu.row(k).tail(rsize);
}
}
return true;
diff --git a/Eigen/src/QR/ColPivHouseholderQR.h b/Eigen/src/QR/ColPivHouseholderQR.h
index 5893125cd..8c9c8840b 100644
--- a/Eigen/src/QR/ColPivHouseholderQR.h
+++ b/Eigen/src/QR/ColPivHouseholderQR.h
@@ -411,7 +411,7 @@ ColPivHouseholderQR<MatrixType>& ColPivHouseholderQR<MatrixType>::compute(const
{
m_nonzero_pivots = k;
m_hCoeffs.tail(size-k).setZero();
- m_qr.corner(BottomRight,rows-k,cols-k)
+ m_qr.bottomRightCorner(rows-k,cols-k)
.template triangularView<StrictlyLower>()
.setZero();
break;
@@ -436,7 +436,7 @@ ColPivHouseholderQR<MatrixType>& ColPivHouseholderQR<MatrixType>::compute(const
if(ei_abs(beta) > m_maxpivot) m_maxpivot = ei_abs(beta);
// apply the householder transformation
- m_qr.corner(BottomRight, rows-k, cols-k-1)
+ m_qr.bottomRightCorner(rows-k, cols-k-1)
.applyHouseholderOnTheLeft(m_qr.col(k).tail(rows-k-1), m_hCoeffs.coeffRef(k), &m_temp.coeffRef(k+1));
// update our table of squared norms of the columns
@@ -483,17 +483,17 @@ struct ei_solve_retval<ColPivHouseholderQR<_MatrixType>, Rhs>
));
dec().matrixQR()
- .corner(TopLeft, nonzero_pivots, nonzero_pivots)
+ .topLeftCorner(nonzero_pivots, nonzero_pivots)
.template triangularView<Upper>()
- .solveInPlace(c.corner(TopLeft, nonzero_pivots, c.cols()));
+ .solveInPlace(c.topRows(nonzero_pivots));
typename Rhs::PlainObject d(c);
- d.corner(TopLeft, nonzero_pivots, c.cols())
+ d.topRows(nonzero_pivots)
= dec().matrixQR()
- .corner(TopLeft, nonzero_pivots, nonzero_pivots)
+ .topLeftCorner(nonzero_pivots, nonzero_pivots)
.template triangularView<Upper>()
- * c.corner(TopLeft, nonzero_pivots, c.cols());
+ * c.topRows(nonzero_pivots);
for(int i = 0; i < nonzero_pivots; ++i) dst.row(dec().colsPermutation().indices().coeff(i)) = c.row(i);
for(int i = nonzero_pivots; i < cols; ++i) dst.row(dec().colsPermutation().indices().coeff(i)).setZero();
diff --git a/Eigen/src/QR/FullPivHouseholderQR.h b/Eigen/src/QR/FullPivHouseholderQR.h
index 8df2ed1e0..0195e1330 100644
--- a/Eigen/src/QR/FullPivHouseholderQR.h
+++ b/Eigen/src/QR/FullPivHouseholderQR.h
@@ -314,7 +314,7 @@ FullPivHouseholderQR<MatrixType>& FullPivHouseholderQR<MatrixType>::compute(cons
int row_of_biggest_in_corner, col_of_biggest_in_corner;
RealScalar biggest_in_corner;
- biggest_in_corner = m_qr.corner(Eigen::BottomRight, rows-k, cols-k)
+ biggest_in_corner = m_qr.bottomRightCorner(rows-k, cols-k)
.cwiseAbs()
.maxCoeff(&row_of_biggest_in_corner, &col_of_biggest_in_corner);
row_of_biggest_in_corner += k;
@@ -349,7 +349,7 @@ FullPivHouseholderQR<MatrixType>& FullPivHouseholderQR<MatrixType>::compute(cons
m_qr.col(k).tail(rows-k).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta);
m_qr.coeffRef(k,k) = beta;
- m_qr.corner(BottomRight, rows-k, cols-k-1)
+ m_qr.bottomRightCorner(rows-k, cols-k-1)
.applyHouseholderOnTheLeft(m_qr.col(k).tail(rows-k-1), m_hCoeffs.coeffRef(k), &m_temp.coeffRef(k+1));
}
@@ -389,7 +389,7 @@ struct ei_solve_retval<FullPivHouseholderQR<_MatrixType>, Rhs>
{
int remainingSize = rows-k;
c.row(k).swap(c.row(dec().rowsTranspositions().coeff(k)));
- c.corner(BottomRight, remainingSize, rhs().cols())
+ c.bottomRightCorner(remainingSize, rhs().cols())
.applyHouseholderOnTheLeft(dec().matrixQR().col(k).tail(remainingSize-1),
dec().hCoeffs().coeff(k), &temp.coeffRef(0));
}
@@ -397,17 +397,17 @@ struct ei_solve_retval<FullPivHouseholderQR<_MatrixType>, Rhs>
if(!dec().isSurjective())
{
// is c is in the image of R ?
- RealScalar biggest_in_upper_part_of_c = c.corner(TopLeft, dec().rank(), c.cols()).cwiseAbs().maxCoeff();
- RealScalar biggest_in_lower_part_of_c = c.corner(BottomLeft, rows-dec().rank(), c.cols()).cwiseAbs().maxCoeff();
+ RealScalar biggest_in_upper_part_of_c = c.topRows( dec().rank() ).cwiseAbs().maxCoeff();
+ RealScalar biggest_in_lower_part_of_c = c.bottomRows(rows-dec().rank()).cwiseAbs().maxCoeff();
// FIXME brain dead
const RealScalar m_precision = NumTraits<Scalar>::epsilon() * std::min(rows,cols);
if(!ei_isMuchSmallerThan(biggest_in_lower_part_of_c, biggest_in_upper_part_of_c, m_precision))
return;
}
dec().matrixQR()
- .corner(TopLeft, dec().rank(), dec().rank())
+ .topLeftCorner(dec().rank(), dec().rank())
.template triangularView<Upper>()
- .solveInPlace(c.corner(TopLeft, dec().rank(), c.cols()));
+ .solveInPlace(c.topRows(dec().rank()));
for(int i = 0; i < dec().rank(); ++i) dst.row(dec().colsPermutation().indices().coeff(i)) = c.row(i);
for(int i = dec().rank(); i < cols; ++i) dst.row(dec().colsPermutation().indices().coeff(i)).setZero();
diff --git a/Eigen/src/QR/HouseholderQR.h b/Eigen/src/QR/HouseholderQR.h
index 10b9fb93f..6a2883939 100644
--- a/Eigen/src/QR/HouseholderQR.h
+++ b/Eigen/src/QR/HouseholderQR.h
@@ -216,7 +216,7 @@ HouseholderQR<MatrixType>& HouseholderQR<MatrixType>::compute(const MatrixType&
m_qr.coeffRef(k,k) = beta;
// apply H to remaining part of m_qr from the left
- m_qr.corner(BottomRight, remainingRows, remainingCols)
+ m_qr.bottomRightCorner(remainingRows, remainingCols)
.applyHouseholderOnTheLeft(m_qr.col(k).tail(remainingRows-1), m_hCoeffs.coeffRef(k), &m_temp.coeffRef(k+1));
}
m_isInitialized = true;
@@ -239,17 +239,17 @@ struct ei_solve_retval<HouseholderQR<_MatrixType>, Rhs>
// Note that the matrix Q = H_0^* H_1^*... so its inverse is Q^* = (H_0 H_1 ...)^T
c.applyOnTheLeft(householderSequence(
- dec().matrixQR().corner(TopLeft,rows,rank),
+ dec().matrixQR().leftCols(rank),
dec().hCoeffs().head(rank)).transpose()
);
dec().matrixQR()
- .corner(TopLeft, rank, rank)
+ .topLeftCorner(rank, rank)
.template triangularView<Upper>()
- .solveInPlace(c.corner(TopLeft, rank, c.cols()));
+ .solveInPlace(c.topRows(rank));
- dst.corner(TopLeft, rank, c.cols()) = c.corner(TopLeft, rank, c.cols());
- dst.corner(BottomLeft, cols-rank, c.cols()).setZero();
+ dst.topRows(rank) = c.topRows(rank);
+ dst.bottomRows(cols-rank).setZero();
}
};
diff --git a/Eigen/src/SVD/UpperBidiagonalization.h b/Eigen/src/SVD/UpperBidiagonalization.h
index 0a63b4265..53e04076a 100644
--- a/Eigen/src/SVD/UpperBidiagonalization.h
+++ b/Eigen/src/SVD/UpperBidiagonalization.h
@@ -114,7 +114,7 @@ UpperBidiagonalization<_MatrixType>& UpperBidiagonalization<_MatrixType>::comput
.makeHouseholderInPlace(m_householder.coeffRef(k,k),
m_bidiagonal.template diagonal<0>().coeffRef(k));
// apply householder transform to remaining part of m_householder on the left
- m_householder.corner(BottomRight, remainingRows, remainingCols)
+ m_householder.bottomRightCorner(remainingRows, remainingCols)
.applyHouseholderOnTheLeft(m_householder.col(k).tail(remainingRows-1),
m_householder.coeff(k,k),
temp.data());
@@ -126,7 +126,7 @@ UpperBidiagonalization<_MatrixType>& UpperBidiagonalization<_MatrixType>::comput
.makeHouseholderInPlace(m_householder.coeffRef(k,k+1),
m_bidiagonal.template diagonal<1>().coeffRef(k));
// apply householder transform to remaining part of m_householder on the left
- m_householder.corner(BottomRight, remainingRows-1, remainingCols)
+ m_householder.bottomRightCorner(remainingRows-1, remainingCols)
.applyHouseholderOnTheRight(m_householder.row(k).tail(remainingCols-1).transpose(),
m_householder.coeff(k,k+1),
temp.data());
diff --git a/Eigen/src/Sparse/SparseMatrixBase.h b/Eigen/src/Sparse/SparseMatrixBase.h
index d269ce604..65fa19a79 100644
--- a/Eigen/src/Sparse/SparseMatrixBase.h
+++ b/Eigen/src/Sparse/SparseMatrixBase.h
@@ -433,19 +433,11 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived>
// typename BlockReturnType<Derived,Dynamic>::SubVectorType end(int size);
// const typename BlockReturnType<Derived,Dynamic>::SubVectorType end(int size) const;
//
-// typename BlockReturnType<Derived>::Type corner(CornerType type, int cRows, int cCols);
-// const typename BlockReturnType<Derived>::Type corner(CornerType type, int cRows, int cCols) const;
-//
// template<int BlockRows, int BlockCols>
// typename BlockReturnType<Derived, BlockRows, BlockCols>::Type block(int startRow, int startCol);
// template<int BlockRows, int BlockCols>
// const typename BlockReturnType<Derived, BlockRows, BlockCols>::Type block(int startRow, int startCol) const;
-// template<int CRows, int CCols>
-// typename BlockReturnType<Derived, CRows, CCols>::Type corner(CornerType type);
-// template<int CRows, int CCols>
-// const typename BlockReturnType<Derived, CRows, CCols>::Type corner(CornerType type) const;
-
// template<int Size> typename BlockReturnType<Derived,Size>::SubVectorType start(void);
// template<int Size> const typename BlockReturnType<Derived,Size>::SubVectorType start() const;