diff options
author | Gael Guennebaud <g.gael@free.fr> | 2010-01-05 13:07:32 +0100 |
---|---|---|
committer | Gael Guennebaud <g.gael@free.fr> | 2010-01-05 13:07:32 +0100 |
commit | 9d9e00b6080153ddaa26ccfce922d7814811a1ae (patch) | |
tree | a18d77d660e3734a21daec2637c2066afab9021d /Eigen/src/QR | |
parent | 90d2ae7fec1000c244472c94af24126c5f2ca2a2 (diff) | |
parent | 51b8f014f30d0f64fcd4f6dff4b1afa64f8ace48 (diff) |
merge and add start/end to Eigen2Support
Diffstat (limited to 'Eigen/src/QR')
-rw-r--r-- | Eigen/src/QR/ColPivHouseholderQR.h | 14 | ||||
-rw-r--r-- | Eigen/src/QR/FullPivHouseholderQR.h | 10 | ||||
-rw-r--r-- | Eigen/src/QR/HouseholderQR.h | 6 |
3 files changed, 15 insertions, 15 deletions
diff --git a/Eigen/src/QR/ColPivHouseholderQR.h b/Eigen/src/QR/ColPivHouseholderQR.h index b4c1a5fcc..8b705de3c 100644 --- a/Eigen/src/QR/ColPivHouseholderQR.h +++ b/Eigen/src/QR/ColPivHouseholderQR.h @@ -359,14 +359,14 @@ ColPivHouseholderQR<MatrixType>& ColPivHouseholderQR<MatrixType>::compute(const { // first, we look up in our table colSqNorms which column has the biggest squared norm int biggest_col_index; - RealScalar biggest_col_sq_norm = colSqNorms.end(cols-k).maxCoeff(&biggest_col_index); + RealScalar biggest_col_sq_norm = colSqNorms.tail(cols-k).maxCoeff(&biggest_col_index); biggest_col_index += k; // since our table colSqNorms accumulates imprecision at every step, we must now recompute // the actual squared norm of the selected column. // Note that not doing so does result in solve() sometimes returning inf/nan values // when running the unit test with 1000 repetitions. - biggest_col_sq_norm = m_qr.col(biggest_col_index).end(rows-k).squaredNorm(); + biggest_col_sq_norm = m_qr.col(biggest_col_index).tail(rows-k).squaredNorm(); // we store that back into our table: it can't hurt to correct our table. colSqNorms.coeffRef(biggest_col_index) = biggest_col_sq_norm; @@ -379,7 +379,7 @@ ColPivHouseholderQR<MatrixType>& ColPivHouseholderQR<MatrixType>::compute(const if(biggest_col_sq_norm < threshold_helper * (rows-k)) { m_nonzero_pivots = k; - m_hCoeffs.end(size-k).setZero(); + m_hCoeffs.tail(size-k).setZero(); m_qr.corner(BottomRight,rows-k,cols-k) .template triangularView<StrictlyLowerTriangular>() .setZero(); @@ -396,7 +396,7 @@ ColPivHouseholderQR<MatrixType>& ColPivHouseholderQR<MatrixType>::compute(const // generate the householder vector, store it below the diagonal RealScalar beta; - m_qr.col(k).end(rows-k).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta); + m_qr.col(k).tail(rows-k).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta); // apply the householder transformation to the diagonal coefficient m_qr.coeffRef(k,k) = beta; @@ -406,10 +406,10 @@ ColPivHouseholderQR<MatrixType>& ColPivHouseholderQR<MatrixType>::compute(const // apply the householder transformation m_qr.corner(BottomRight, rows-k, cols-k-1) - .applyHouseholderOnTheLeft(m_qr.col(k).end(rows-k-1), m_hCoeffs.coeffRef(k), &temp.coeffRef(k+1)); + .applyHouseholderOnTheLeft(m_qr.col(k).tail(rows-k-1), m_hCoeffs.coeffRef(k), &temp.coeffRef(k+1)); // update our table of squared norms of the columns - colSqNorms.end(cols-k-1) -= m_qr.row(k).end(cols-k-1).cwiseAbs2(); + colSqNorms.tail(cols-k-1) -= m_qr.row(k).tail(cols-k-1).cwiseAbs2(); } m_cols_permutation.setIdentity(cols); @@ -427,7 +427,7 @@ struct ei_solve_retval<ColPivHouseholderQR<_MatrixType>, Rhs> : ei_solve_retval_base<ColPivHouseholderQR<_MatrixType>, Rhs> { EIGEN_MAKE_SOLVE_HELPERS(ColPivHouseholderQR<_MatrixType>,Rhs) - + template<typename Dest> void evalTo(Dest& dst) const { const int rows = dec().rows(), cols = dec().cols(), diff --git a/Eigen/src/QR/FullPivHouseholderQR.h b/Eigen/src/QR/FullPivHouseholderQR.h index 51609ca1a..598f44dc3 100644 --- a/Eigen/src/QR/FullPivHouseholderQR.h +++ b/Eigen/src/QR/FullPivHouseholderQR.h @@ -306,7 +306,7 @@ FullPivHouseholderQR<MatrixType>& FullPivHouseholderQR<MatrixType>::compute(cons m_rows_transpositions.coeffRef(k) = row_of_biggest_in_corner; cols_transpositions.coeffRef(k) = col_of_biggest_in_corner; if(k != row_of_biggest_in_corner) { - m_qr.row(k).end(cols-k).swap(m_qr.row(row_of_biggest_in_corner).end(cols-k)); + m_qr.row(k).tail(cols-k).swap(m_qr.row(row_of_biggest_in_corner).tail(cols-k)); ++number_of_transpositions; } if(k != col_of_biggest_in_corner) { @@ -315,11 +315,11 @@ FullPivHouseholderQR<MatrixType>& FullPivHouseholderQR<MatrixType>::compute(cons } RealScalar beta; - m_qr.col(k).end(rows-k).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta); + 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) - .applyHouseholderOnTheLeft(m_qr.col(k).end(rows-k-1), m_hCoeffs.coeffRef(k), &temp.coeffRef(k+1)); + .applyHouseholderOnTheLeft(m_qr.col(k).tail(rows-k-1), m_hCoeffs.coeffRef(k), &temp.coeffRef(k+1)); } m_cols_permutation.setIdentity(cols); @@ -360,7 +360,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()) - .applyHouseholderOnTheLeft(dec().matrixQR().col(k).end(remainingSize-1), + .applyHouseholderOnTheLeft(dec().matrixQR().col(k).tail(remainingSize-1), dec().hCoeffs().coeff(k), &temp.coeffRef(0)); } @@ -400,7 +400,7 @@ typename FullPivHouseholderQR<MatrixType>::MatrixQType FullPivHouseholderQR<Matr for (int k = size-1; k >= 0; k--) { res.block(k, k, rows-k, rows-k) - .applyHouseholderOnTheLeft(m_qr.col(k).end(rows-k-1), ei_conj(m_hCoeffs.coeff(k)), &temp.coeffRef(k)); + .applyHouseholderOnTheLeft(m_qr.col(k).tail(rows-k-1), ei_conj(m_hCoeffs.coeff(k)), &temp.coeffRef(k)); res.row(k).swap(res.row(m_rows_transpositions.coeff(k))); } return res; diff --git a/Eigen/src/QR/HouseholderQR.h b/Eigen/src/QR/HouseholderQR.h index 895ae046a..24aa96e04 100644 --- a/Eigen/src/QR/HouseholderQR.h +++ b/Eigen/src/QR/HouseholderQR.h @@ -197,12 +197,12 @@ HouseholderQR<MatrixType>& HouseholderQR<MatrixType>::compute(const MatrixType& int remainingCols = cols - k - 1; RealScalar beta; - m_qr.col(k).end(remainingRows).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta); + m_qr.col(k).tail(remainingRows).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta); m_qr.coeffRef(k,k) = beta; // apply H to remaining part of m_qr from the left m_qr.corner(BottomRight, remainingRows, remainingCols) - .applyHouseholderOnTheLeft(m_qr.col(k).end(remainingRows-1), m_hCoeffs.coeffRef(k), &temp.coeffRef(k+1)); + .applyHouseholderOnTheLeft(m_qr.col(k).tail(remainingRows-1), m_hCoeffs.coeffRef(k), &temp.coeffRef(k+1)); } m_isInitialized = true; return *this; @@ -226,7 +226,7 @@ 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().hCoeffs().start(rank)).transpose() + dec().hCoeffs().head(rank)).transpose() ); dec().matrixQR() |