Eigne Matrix Class Library
Dependents: Eigen_test Odometry_test AttitudeEstimation_usingTicker MPU9250_Quaternion_Binary_Serial ... more
RealSchur.h
00001 // This file is part of Eigen, a lightweight C++ template library 00002 // for linear algebra. 00003 // 00004 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@inria.fr> 00005 // Copyright (C) 2010,2012 Jitse Niesen <jitse@maths.leeds.ac.uk> 00006 // 00007 // This Source Code Form is subject to the terms of the Mozilla 00008 // Public License v. 2.0. If a copy of the MPL was not distributed 00009 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00010 00011 #ifndef EIGEN_REAL_SCHUR_H 00012 #define EIGEN_REAL_SCHUR_H 00013 00014 #include "./HessenbergDecomposition.h" 00015 00016 namespace Eigen { 00017 00018 /** \eigenvalues_module \ingroup Eigenvalues_Module 00019 * 00020 * 00021 * \class RealSchur 00022 * 00023 * \brief Performs a real Schur decomposition of a square matrix 00024 * 00025 * \tparam _MatrixType the type of the matrix of which we are computing the 00026 * real Schur decomposition; this is expected to be an instantiation of the 00027 * Matrix class template. 00028 * 00029 * Given a real square matrix A, this class computes the real Schur 00030 * decomposition: \f$ A = U T U^T \f$ where U is a real orthogonal matrix and 00031 * T is a real quasi-triangular matrix. An orthogonal matrix is a matrix whose 00032 * inverse is equal to its transpose, \f$ U^{-1} = U^T \f$. A quasi-triangular 00033 * matrix is a block-triangular matrix whose diagonal consists of 1-by-1 00034 * blocks and 2-by-2 blocks with complex eigenvalues. The eigenvalues of the 00035 * blocks on the diagonal of T are the same as the eigenvalues of the matrix 00036 * A, and thus the real Schur decomposition is used in EigenSolver to compute 00037 * the eigendecomposition of a matrix. 00038 * 00039 * Call the function compute() to compute the real Schur decomposition of a 00040 * given matrix. Alternatively, you can use the RealSchur(const MatrixType&, bool) 00041 * constructor which computes the real Schur decomposition at construction 00042 * time. Once the decomposition is computed, you can use the matrixU() and 00043 * matrixT() functions to retrieve the matrices U and T in the decomposition. 00044 * 00045 * The documentation of RealSchur(const MatrixType&, bool) contains an example 00046 * of the typical use of this class. 00047 * 00048 * \note The implementation is adapted from 00049 * <a href="http://math.nist.gov/javanumerics/jama/">JAMA</a> (public domain). 00050 * Their code is based on EISPACK. 00051 * 00052 * \sa class ComplexSchur, class EigenSolver, class ComplexEigenSolver 00053 */ 00054 template<typename _MatrixType> class RealSchur 00055 { 00056 public: 00057 typedef _MatrixType MatrixType; 00058 enum { 00059 RowsAtCompileTime = MatrixType::RowsAtCompileTime, 00060 ColsAtCompileTime = MatrixType::ColsAtCompileTime, 00061 Options = MatrixType::Options, 00062 MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime, 00063 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime 00064 }; 00065 typedef typename MatrixType::Scalar Scalar; 00066 typedef std::complex<typename NumTraits<Scalar>::Real> ComplexScalar; 00067 typedef typename MatrixType::Index Index; 00068 00069 typedef Matrix<ComplexScalar, ColsAtCompileTime, 1, Options & ~RowMajor, MaxColsAtCompileTime, 1> EigenvalueType ; 00070 typedef Matrix<Scalar, ColsAtCompileTime, 1, Options & ~RowMajor, MaxColsAtCompileTime, 1> ColumnVectorType ; 00071 00072 /** \brief Default constructor. 00073 * 00074 * \param [in] size Positive integer, size of the matrix whose Schur decomposition will be computed. 00075 * 00076 * The default constructor is useful in cases in which the user intends to 00077 * perform decompositions via compute(). The \p size parameter is only 00078 * used as a hint. It is not an error to give a wrong \p size, but it may 00079 * impair performance. 00080 * 00081 * \sa compute() for an example. 00082 */ 00083 RealSchur(Index size = RowsAtCompileTime==Dynamic ? 1 : RowsAtCompileTime) 00084 : m_matT(size, size), 00085 m_matU(size, size), 00086 m_workspaceVector(size), 00087 m_hess(size), 00088 m_isInitialized(false), 00089 m_matUisUptodate(false), 00090 m_maxIters(-1) 00091 { } 00092 00093 /** \brief Constructor; computes real Schur decomposition of given matrix. 00094 * 00095 * \param[in] matrix Square matrix whose Schur decomposition is to be computed. 00096 * \param[in] computeU If true, both T and U are computed; if false, only T is computed. 00097 * 00098 * This constructor calls compute() to compute the Schur decomposition. 00099 * 00100 * Example: \include RealSchur_RealSchur_MatrixType.cpp 00101 * Output: \verbinclude RealSchur_RealSchur_MatrixType.out 00102 */ 00103 RealSchur(const MatrixType& matrix, bool computeU = true) 00104 : m_matT(matrix.rows(),matrix.cols()), 00105 m_matU(matrix.rows(),matrix.cols()), 00106 m_workspaceVector(matrix.rows()), 00107 m_hess(matrix.rows()), 00108 m_isInitialized(false), 00109 m_matUisUptodate(false), 00110 m_maxIters(-1) 00111 { 00112 compute(matrix, computeU); 00113 } 00114 00115 /** \brief Returns the orthogonal matrix in the Schur decomposition. 00116 * 00117 * \returns A const reference to the matrix U. 00118 * 00119 * \pre Either the constructor RealSchur(const MatrixType&, bool) or the 00120 * member function compute(const MatrixType&, bool) has been called before 00121 * to compute the Schur decomposition of a matrix, and \p computeU was set 00122 * to true (the default value). 00123 * 00124 * \sa RealSchur(const MatrixType&, bool) for an example 00125 */ 00126 const MatrixType& matrixU() const 00127 { 00128 eigen_assert(m_isInitialized && "RealSchur is not initialized."); 00129 eigen_assert(m_matUisUptodate && "The matrix U has not been computed during the RealSchur decomposition."); 00130 return m_matU; 00131 } 00132 00133 /** \brief Returns the quasi-triangular matrix in the Schur decomposition. 00134 * 00135 * \returns A const reference to the matrix T. 00136 * 00137 * \pre Either the constructor RealSchur(const MatrixType&, bool) or the 00138 * member function compute(const MatrixType&, bool) has been called before 00139 * to compute the Schur decomposition of a matrix. 00140 * 00141 * \sa RealSchur(const MatrixType&, bool) for an example 00142 */ 00143 const MatrixType& matrixT() const 00144 { 00145 eigen_assert(m_isInitialized && "RealSchur is not initialized."); 00146 return m_matT; 00147 } 00148 00149 /** \brief Computes Schur decomposition of given matrix. 00150 * 00151 * \param[in] matrix Square matrix whose Schur decomposition is to be computed. 00152 * \param[in] computeU If true, both T and U are computed; if false, only T is computed. 00153 * \returns Reference to \c *this 00154 * 00155 * The Schur decomposition is computed by first reducing the matrix to 00156 * Hessenberg form using the class HessenbergDecomposition. The Hessenberg 00157 * matrix is then reduced to triangular form by performing Francis QR 00158 * iterations with implicit double shift. The cost of computing the Schur 00159 * decomposition depends on the number of iterations; as a rough guide, it 00160 * may be taken to be \f$25n^3\f$ flops if \a computeU is true and 00161 * \f$10n^3\f$ flops if \a computeU is false. 00162 * 00163 * Example: \include RealSchur_compute.cpp 00164 * Output: \verbinclude RealSchur_compute.out 00165 * 00166 * \sa compute(const MatrixType&, bool, Index) 00167 */ 00168 RealSchur & compute(const MatrixType& matrix, bool computeU = true); 00169 00170 /** \brief Computes Schur decomposition of a Hessenberg matrix H = Z T Z^T 00171 * \param[in] matrixH Matrix in Hessenberg form H 00172 * \param[in] matrixQ orthogonal matrix Q that transform a matrix A to H : A = Q H Q^T 00173 * \param computeU Computes the matriX U of the Schur vectors 00174 * \return Reference to \c *this 00175 * 00176 * This routine assumes that the matrix is already reduced in Hessenberg form matrixH 00177 * using either the class HessenbergDecomposition or another mean. 00178 * It computes the upper quasi-triangular matrix T of the Schur decomposition of H 00179 * When computeU is true, this routine computes the matrix U such that 00180 * A = U T U^T = (QZ) T (QZ)^T = Q H Q^T where A is the initial matrix 00181 * 00182 * NOTE Q is referenced if computeU is true; so, if the initial orthogonal matrix 00183 * is not available, the user should give an identity matrix (Q.setIdentity()) 00184 * 00185 * \sa compute(const MatrixType&, bool) 00186 */ 00187 template<typename HessMatrixType, typename OrthMatrixType> 00188 RealSchur & computeFromHessenberg(const HessMatrixType& matrixH, const OrthMatrixType& matrixQ, bool computeU); 00189 /** \brief Reports whether previous computation was successful. 00190 * 00191 * \returns \c Success if computation was succesful, \c NoConvergence otherwise. 00192 */ 00193 ComputationInfo info() const 00194 { 00195 eigen_assert(m_isInitialized && "RealSchur is not initialized."); 00196 return m_info; 00197 } 00198 00199 /** \brief Sets the maximum number of iterations allowed. 00200 * 00201 * If not specified by the user, the maximum number of iterations is m_maxIterationsPerRow times the size 00202 * of the matrix. 00203 */ 00204 RealSchur & setMaxIterations(Index maxIters) 00205 { 00206 m_maxIters = maxIters; 00207 return *this; 00208 } 00209 00210 /** \brief Returns the maximum number of iterations. */ 00211 Index getMaxIterations() 00212 { 00213 return m_maxIters; 00214 } 00215 00216 /** \brief Maximum number of iterations per row. 00217 * 00218 * If not otherwise specified, the maximum number of iterations is this number times the size of the 00219 * matrix. It is currently set to 40. 00220 */ 00221 static const int m_maxIterationsPerRow = 40; 00222 00223 private: 00224 00225 MatrixType m_matT; 00226 MatrixType m_matU; 00227 ColumnVectorType m_workspaceVector; 00228 HessenbergDecomposition<MatrixType> m_hess; 00229 ComputationInfo m_info; 00230 bool m_isInitialized; 00231 bool m_matUisUptodate; 00232 Index m_maxIters; 00233 00234 typedef Matrix<Scalar,3,1> Vector3s ; 00235 00236 Scalar computeNormOfT(); 00237 Index findSmallSubdiagEntry(Index iu); 00238 void splitOffTwoRows(Index iu, bool computeU, const Scalar& exshift); 00239 void computeShift(Index iu, Index iter, Scalar& exshift, Vector3s & shiftInfo); 00240 void initFrancisQRStep(Index il, Index iu, const Vector3s & shiftInfo, Index& im, Vector3s & firstHouseholderVector); 00241 void performFrancisQRStep(Index il, Index im, Index iu, bool computeU, const Vector3s & firstHouseholderVector, Scalar* workspace); 00242 }; 00243 00244 00245 template<typename MatrixType> 00246 RealSchur<MatrixType> & RealSchur<MatrixType>::compute(const MatrixType& matrix, bool computeU) 00247 { 00248 eigen_assert(matrix.cols() == matrix.rows()); 00249 Index maxIters = m_maxIters; 00250 if (maxIters == -1) 00251 maxIters = m_maxIterationsPerRow * matrix.rows(); 00252 00253 // Step 1. Reduce to Hessenberg form 00254 m_hess.compute(matrix); 00255 00256 // Step 2. Reduce to real Schur form 00257 computeFromHessenberg(m_hess.matrixH(), m_hess.matrixQ(), computeU); 00258 00259 return *this; 00260 } 00261 template<typename MatrixType> 00262 template<typename HessMatrixType, typename OrthMatrixType> 00263 RealSchur<MatrixType> & RealSchur<MatrixType>::computeFromHessenberg (const HessMatrixType& matrixH, const OrthMatrixType& matrixQ, bool computeU) 00264 { 00265 m_matT = matrixH; 00266 if(computeU) 00267 m_matU = matrixQ; 00268 00269 Index maxIters = m_maxIters; 00270 if (maxIters == -1) 00271 maxIters = m_maxIterationsPerRow * matrixH.rows(); 00272 m_workspaceVector.resize(m_matT.cols()); 00273 Scalar* workspace = &m_workspaceVector.coeffRef(0); 00274 00275 // The matrix m_matT is divided in three parts. 00276 // Rows 0,...,il-1 are decoupled from the rest because m_matT(il,il-1) is zero. 00277 // Rows il,...,iu is the part we are working on (the active window). 00278 // Rows iu+1,...,end are already brought in triangular form. 00279 Index iu = m_matT.cols() - 1; 00280 Index iter = 0; // iteration count for current eigenvalue 00281 Index totalIter = 0; // iteration count for whole matrix 00282 Scalar exshift(0); // sum of exceptional shifts 00283 Scalar norm = computeNormOfT(); 00284 00285 if(norm!=0) 00286 { 00287 while (iu >= 0) 00288 { 00289 Index il = findSmallSubdiagEntry(iu); 00290 00291 // Check for convergence 00292 if (il == iu) // One root found 00293 { 00294 m_matT.coeffRef(iu,iu) = m_matT.coeff(iu,iu) + exshift; 00295 if (iu > 0) 00296 m_matT.coeffRef(iu, iu-1) = Scalar(0); 00297 iu--; 00298 iter = 0; 00299 } 00300 else if (il == iu-1) // Two roots found 00301 { 00302 splitOffTwoRows(iu, computeU, exshift); 00303 iu -= 2; 00304 iter = 0; 00305 } 00306 else // No convergence yet 00307 { 00308 // The firstHouseholderVector vector has to be initialized to something to get rid of a silly GCC warning (-O1 -Wall -DNDEBUG ) 00309 Vector3s firstHouseholderVector(0,0,0), shiftInfo; 00310 computeShift(iu, iter, exshift, shiftInfo); 00311 iter = iter + 1; 00312 totalIter = totalIter + 1; 00313 if (totalIter > maxIters) break; 00314 Index im; 00315 initFrancisQRStep(il, iu, shiftInfo, im, firstHouseholderVector); 00316 performFrancisQRStep(il, im, iu, computeU, firstHouseholderVector, workspace); 00317 } 00318 } 00319 } 00320 if(totalIter <= maxIters) 00321 m_info = Success; 00322 else 00323 m_info = NoConvergence; 00324 00325 m_isInitialized = true; 00326 m_matUisUptodate = computeU; 00327 return *this; 00328 } 00329 00330 /** \internal Computes and returns vector L1 norm of T */ 00331 template<typename MatrixType> 00332 inline typename MatrixType::Scalar RealSchur<MatrixType>::computeNormOfT () 00333 { 00334 const Index size = m_matT.cols(); 00335 // FIXME to be efficient the following would requires a triangular reduxion code 00336 // Scalar norm = m_matT.upper().cwiseAbs().sum() 00337 // + m_matT.bottomLeftCorner(size-1,size-1).diagonal().cwiseAbs().sum(); 00338 Scalar norm(0); 00339 for (Index j = 0; j < size; ++j) 00340 norm += m_matT.col(j).segment(0, (std::min)(size,j+2)).cwiseAbs().sum(); 00341 return norm; 00342 } 00343 00344 /** \internal Look for single small sub-diagonal element and returns its index */ 00345 template<typename MatrixType> 00346 inline typename MatrixType::Index RealSchur<MatrixType>::findSmallSubdiagEntry(Index iu) 00347 { 00348 using std::abs; 00349 Index res = iu; 00350 while (res > 0) 00351 { 00352 Scalar s = abs(m_matT.coeff(res-1,res-1)) + abs(m_matT.coeff(res,res)); 00353 if (abs(m_matT.coeff(res,res-1)) <= NumTraits<Scalar>::epsilon() * s) 00354 break; 00355 res--; 00356 } 00357 return res; 00358 } 00359 00360 /** \internal Update T given that rows iu-1 and iu decouple from the rest. */ 00361 template<typename MatrixType> 00362 inline void RealSchur<MatrixType>::splitOffTwoRows(Index iu, bool computeU, const Scalar& exshift) 00363 { 00364 using std::sqrt; 00365 using std::abs; 00366 const Index size = m_matT.cols(); 00367 00368 // The eigenvalues of the 2x2 matrix [a b; c d] are 00369 // trace +/- sqrt(discr/4) where discr = tr^2 - 4*det, tr = a + d, det = ad - bc 00370 Scalar p = Scalar(0.5) * (m_matT.coeff(iu-1,iu-1) - m_matT.coeff(iu,iu)); 00371 Scalar q = p * p + m_matT.coeff(iu,iu-1) * m_matT.coeff(iu-1,iu); // q = tr^2 / 4 - det = discr/4 00372 m_matT.coeffRef(iu,iu) += exshift; 00373 m_matT.coeffRef(iu-1,iu-1) += exshift; 00374 00375 if (q >= Scalar(0)) // Two real eigenvalues 00376 { 00377 Scalar z = sqrt(abs(q)); 00378 JacobiRotation<Scalar> rot; 00379 if (p >= Scalar(0)) 00380 rot.makeGivens(p + z, m_matT.coeff(iu, iu-1)); 00381 else 00382 rot.makeGivens(p - z, m_matT.coeff(iu, iu-1)); 00383 00384 m_matT.rightCols(size-iu+1).applyOnTheLeft(iu-1, iu, rot.adjoint()); 00385 m_matT.topRows(iu+1).applyOnTheRight(iu-1, iu, rot); 00386 m_matT.coeffRef(iu, iu-1) = Scalar(0); 00387 if (computeU) 00388 m_matU.applyOnTheRight(iu-1, iu, rot); 00389 } 00390 00391 if (iu > 1) 00392 m_matT.coeffRef(iu-1, iu-2) = Scalar(0); 00393 } 00394 00395 /** \internal Form shift in shiftInfo, and update exshift if an exceptional shift is performed. */ 00396 template<typename MatrixType> 00397 inline void RealSchur<MatrixType>::computeShift(Index iu, Index iter, Scalar& exshift, Vector3s& shiftInfo) 00398 { 00399 using std::sqrt; 00400 using std::abs; 00401 shiftInfo.coeffRef(0) = m_matT.coeff(iu,iu); 00402 shiftInfo.coeffRef(1) = m_matT.coeff(iu-1,iu-1); 00403 shiftInfo.coeffRef(2) = m_matT.coeff(iu,iu-1) * m_matT.coeff(iu-1,iu); 00404 00405 // Wilkinson's original ad hoc shift 00406 if (iter == 10) 00407 { 00408 exshift += shiftInfo.coeff(0); 00409 for (Index i = 0; i <= iu; ++i) 00410 m_matT.coeffRef(i,i) -= shiftInfo.coeff(0); 00411 Scalar s = abs(m_matT.coeff(iu,iu-1)) + abs(m_matT.coeff(iu-1,iu-2)); 00412 shiftInfo.coeffRef(0) = Scalar(0.75) * s; 00413 shiftInfo.coeffRef(1) = Scalar(0.75) * s; 00414 shiftInfo.coeffRef(2) = Scalar(-0.4375) * s * s; 00415 } 00416 00417 // MATLAB's new ad hoc shift 00418 if (iter == 30) 00419 { 00420 Scalar s = (shiftInfo.coeff(1) - shiftInfo.coeff(0)) / Scalar(2.0); 00421 s = s * s + shiftInfo.coeff(2); 00422 if (s > Scalar(0)) 00423 { 00424 s = sqrt(s); 00425 if (shiftInfo.coeff(1) < shiftInfo.coeff(0)) 00426 s = -s; 00427 s = s + (shiftInfo.coeff(1) - shiftInfo.coeff(0)) / Scalar(2.0); 00428 s = shiftInfo.coeff(0) - shiftInfo.coeff(2) / s; 00429 exshift += s; 00430 for (Index i = 0; i <= iu; ++i) 00431 m_matT.coeffRef(i,i) -= s; 00432 shiftInfo.setConstant(Scalar(0.964)); 00433 } 00434 } 00435 } 00436 00437 /** \internal Compute index im at which Francis QR step starts and the first Householder vector. */ 00438 template<typename MatrixType> 00439 inline void RealSchur<MatrixType>::initFrancisQRStep(Index il, Index iu, const Vector3s& shiftInfo, Index& im, Vector3s& firstHouseholderVector) 00440 { 00441 using std::abs; 00442 Vector3s& v = firstHouseholderVector; // alias to save typing 00443 00444 for (im = iu-2; im >= il; --im) 00445 { 00446 const Scalar Tmm = m_matT.coeff(im,im); 00447 const Scalar r = shiftInfo.coeff(0) - Tmm; 00448 const Scalar s = shiftInfo.coeff(1) - Tmm; 00449 v.coeffRef(0) = (r * s - shiftInfo.coeff(2)) / m_matT.coeff(im+1,im) + m_matT.coeff(im,im+1); 00450 v.coeffRef(1) = m_matT.coeff(im+1,im+1) - Tmm - r - s; 00451 v.coeffRef(2) = m_matT.coeff(im+2,im+1); 00452 if (im == il) { 00453 break; 00454 } 00455 const Scalar lhs = m_matT.coeff(im,im-1) * (abs(v.coeff(1)) + abs(v.coeff(2))); 00456 const Scalar rhs = v.coeff(0) * (abs(m_matT.coeff(im-1,im-1)) + abs(Tmm) + abs(m_matT.coeff(im+1,im+1))); 00457 if (abs(lhs) < NumTraits<Scalar>::epsilon() * rhs) 00458 break; 00459 } 00460 } 00461 00462 /** \internal Perform a Francis QR step involving rows il:iu and columns im:iu. */ 00463 template<typename MatrixType> 00464 inline void RealSchur<MatrixType>::performFrancisQRStep(Index il, Index im, Index iu, bool computeU, const Vector3s& firstHouseholderVector, Scalar* workspace) 00465 { 00466 eigen_assert(im >= il); 00467 eigen_assert(im <= iu-2); 00468 00469 const Index size = m_matT.cols(); 00470 00471 for (Index k = im; k <= iu-2; ++k) 00472 { 00473 bool firstIteration = (k == im); 00474 00475 Vector3s v; 00476 if (firstIteration) 00477 v = firstHouseholderVector; 00478 else 00479 v = m_matT.template block<3,1>(k,k-1); 00480 00481 Scalar tau, beta; 00482 Matrix<Scalar, 2, 1> ess; 00483 v.makeHouseholder(ess, tau, beta); 00484 00485 if (beta != Scalar(0)) // if v is not zero 00486 { 00487 if (firstIteration && k > il) 00488 m_matT.coeffRef(k,k-1) = -m_matT.coeff(k,k-1); 00489 else if (!firstIteration) 00490 m_matT.coeffRef(k,k-1) = beta; 00491 00492 // These Householder transformations form the O(n^3) part of the algorithm 00493 m_matT.block(k, k, 3, size-k).applyHouseholderOnTheLeft(ess, tau, workspace); 00494 m_matT.block(0, k, (std::min)(iu,k+3) + 1, 3).applyHouseholderOnTheRight(ess, tau, workspace); 00495 if (computeU) 00496 m_matU.block(0, k, size, 3).applyHouseholderOnTheRight(ess, tau, workspace); 00497 } 00498 } 00499 00500 Matrix<Scalar, 2, 1> v = m_matT.template block<2,1>(iu-1, iu-2); 00501 Scalar tau, beta; 00502 Matrix<Scalar, 1, 1> ess; 00503 v.makeHouseholder(ess, tau, beta); 00504 00505 if (beta != Scalar(0)) // if v is not zero 00506 { 00507 m_matT.coeffRef(iu-1, iu-2) = beta; 00508 m_matT.block(iu-1, iu-1, 2, size-iu+1).applyHouseholderOnTheLeft(ess, tau, workspace); 00509 m_matT.block(0, iu-1, iu+1, 2).applyHouseholderOnTheRight(ess, tau, workspace); 00510 if (computeU) 00511 m_matU.block(0, iu-1, size, 2).applyHouseholderOnTheRight(ess, tau, workspace); 00512 } 00513 00514 // clean up pollution due to round-off errors 00515 for (Index i = im+2; i <= iu; ++i) 00516 { 00517 m_matT.coeffRef(i,i-2) = Scalar(0); 00518 if (i > im+2) 00519 m_matT.coeffRef(i,i-3) = Scalar(0); 00520 } 00521 } 00522 00523 } // end namespace Eigen 00524 00525 #endif // EIGEN_REAL_SCHUR_H
Generated on Tue Jul 12 2022 17:46:59 by 1.7.2