Eigen::RealSchur< MatrixType_ > Class Template Reference

Performs a real Schur decomposition of a square matrix. More...

Public Types

enum  {
  RowsAtCompileTime ,
  ColsAtCompileTime ,
  Options ,
  MaxRowsAtCompileTime ,
  MaxColsAtCompileTime
}
 
typedef Matrix< Scalar, ColsAtCompileTime, 1, Options &~RowMajor, MaxColsAtCompileTime, 1 > ColumnVectorType
 
typedef std::complex< typename NumTraits< Scalar >::Real > ComplexScalar
 
typedef Matrix< ComplexScalar, ColsAtCompileTime, 1, Options &~RowMajor, MaxColsAtCompileTime, 1 > EigenvalueType
 
typedef Eigen::Index Index
 
typedef MatrixType_ MatrixType
 
typedef MatrixType::Scalar Scalar
 

Public Member Functions

template<typename InputType >
RealSchur< MatrixType > & compute (const EigenBase< InputType > &matrix, bool computeU)
 
template<typename InputType >
RealSchurcompute (const EigenBase< InputType > &matrix, bool computeU=true)
 Computes Schur decomposition of given matrix. More...
 
template<typename HessMatrixType , typename OrthMatrixType >
RealSchurcomputeFromHessenberg (const HessMatrixType &matrixH, const OrthMatrixType &matrixQ, bool computeU)
 Computes Schur decomposition of a Hessenberg matrix H = Z T Z^T. More...
 
template<typename HessMatrixType , typename OrthMatrixType >
RealSchur< MatrixType > & computeFromHessenberg (const HessMatrixType &matrixH, const OrthMatrixType &matrixQ, bool computeU)
 
Index getMaxIterations ()
 Returns the maximum number of iterations. More...
 
ComputationInfo info () const
 Reports whether previous computation was successful. More...
 
const MatrixTypematrixT () const
 Returns the quasi-triangular matrix in the Schur decomposition. More...
 
const MatrixTypematrixU () const
 Returns the orthogonal matrix in the Schur decomposition. More...
 
template<typename InputType >
 RealSchur (const EigenBase< InputType > &matrix, bool computeU=true)
 Constructor; computes real Schur decomposition of given matrix. More...
 
 RealSchur (Index size=RowsAtCompileTime==Dynamic ? 1 :RowsAtCompileTime)
 Default constructor. More...
 
RealSchursetMaxIterations (Index maxIters)
 Sets the maximum number of iterations allowed. More...
 

Static Public Attributes

static const int m_maxIterationsPerRow
 Maximum number of iterations per row. More...
 

Private Types

typedef Matrix< Scalar, 3, 1 > Vector3s
 

Private Member Functions

Scalar computeNormOfT ()
 
void computeShift (Index iu, Index iter, Scalar &exshift, Vector3s &shiftInfo)
 
Index findSmallSubdiagEntry (Index iu, const Scalar &considerAsZero)
 
void initFrancisQRStep (Index il, Index iu, const Vector3s &shiftInfo, Index &im, Vector3s &firstHouseholderVector)
 
void performFrancisQRStep (Index il, Index im, Index iu, bool computeU, const Vector3s &firstHouseholderVector, Scalar *workspace)
 
void splitOffTwoRows (Index iu, bool computeU, const Scalar &exshift)
 

Private Attributes

HessenbergDecomposition< MatrixTypem_hess
 
ComputationInfo m_info
 
bool m_isInitialized
 
MatrixType m_matT
 
MatrixType m_matU
 
bool m_matUisUptodate
 
Index m_maxIters
 
ColumnVectorType m_workspaceVector
 

Detailed Description

template<typename MatrixType_>
class Eigen::RealSchur< MatrixType_ >

Performs a real Schur decomposition of a square matrix.

This is defined in the Eigenvalues module.

Template Parameters
MatrixType_the type of the matrix of which we are computing the real Schur decomposition; this is expected to be an instantiation of the Matrix class template.

Given a real square matrix A, this class computes the real Schur decomposition: \( A = U T U^T \) where U is a real orthogonal matrix and T is a real quasi-triangular matrix. An orthogonal matrix is a matrix whose inverse is equal to its transpose, \( U^{-1} = U^T \). A quasi-triangular matrix is a block-triangular matrix whose diagonal consists of 1-by-1 blocks and 2-by-2 blocks with complex eigenvalues. The eigenvalues of the blocks on the diagonal of T are the same as the eigenvalues of the matrix A, and thus the real Schur decomposition is used in EigenSolver to compute the eigendecomposition of a matrix.

Call the function compute() to compute the real Schur decomposition of a given matrix. Alternatively, you can use the RealSchur(const MatrixType&, bool) constructor which computes the real Schur decomposition at construction time. Once the decomposition is computed, you can use the matrixU() and matrixT() functions to retrieve the matrices U and T in the decomposition.

The documentation of RealSchur(const MatrixType&, bool) contains an example of the typical use of this class.

Note
The implementation is adapted from JAMA (public domain). Their code is based on EISPACK.
See also
class ComplexSchur, class EigenSolver, class ComplexEigenSolver

Definition at line 56 of file RealSchur.h.

Member Typedef Documentation

◆ ColumnVectorType

template<typename MatrixType_ >
typedef Matrix<Scalar, ColsAtCompileTime, 1, Options & ~RowMajor, MaxColsAtCompileTime, 1> Eigen::RealSchur< MatrixType_ >::ColumnVectorType

Definition at line 72 of file RealSchur.h.

◆ ComplexScalar

template<typename MatrixType_ >
typedef std::complex<typename NumTraits<Scalar>::Real> Eigen::RealSchur< MatrixType_ >::ComplexScalar

Definition at line 68 of file RealSchur.h.

◆ EigenvalueType

template<typename MatrixType_ >
typedef Matrix<ComplexScalar, ColsAtCompileTime, 1, Options & ~RowMajor, MaxColsAtCompileTime, 1> Eigen::RealSchur< MatrixType_ >::EigenvalueType

Definition at line 71 of file RealSchur.h.

◆ Index

template<typename MatrixType_ >
typedef Eigen::Index Eigen::RealSchur< MatrixType_ >::Index
Deprecated:
since Eigen 3.3

Definition at line 69 of file RealSchur.h.

◆ MatrixType

template<typename MatrixType_ >
typedef MatrixType_ Eigen::RealSchur< MatrixType_ >::MatrixType

Definition at line 59 of file RealSchur.h.

◆ Scalar

template<typename MatrixType_ >
typedef MatrixType::Scalar Eigen::RealSchur< MatrixType_ >::Scalar

Definition at line 67 of file RealSchur.h.

◆ Vector3s

template<typename MatrixType_ >
typedef Matrix<Scalar,3,1> Eigen::RealSchur< MatrixType_ >::Vector3s
private

Definition at line 238 of file RealSchur.h.

Member Enumeration Documentation

◆ anonymous enum

template<typename MatrixType_ >
anonymous enum
Enumerator
RowsAtCompileTime 
ColsAtCompileTime 
Options 
MaxRowsAtCompileTime 
MaxColsAtCompileTime 

Definition at line 60 of file RealSchur.h.

60  {
61  RowsAtCompileTime = MatrixType::RowsAtCompileTime,
62  ColsAtCompileTime = MatrixType::ColsAtCompileTime,
63  Options = MatrixType::Options,
64  MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
65  MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
66  };

Constructor & Destructor Documentation

◆ RealSchur() [1/2]

template<typename MatrixType_ >
Eigen::RealSchur< MatrixType_ >::RealSchur ( Index  size = RowsAtCompileTime==Dynamic ? 1 : RowsAtCompileTime)
inlineexplicit

Default constructor.

Parameters
[in]sizePositive integer, size of the matrix whose Schur decomposition will be computed.

The default constructor is useful in cases in which the user intends to perform decompositions via compute(). The size parameter is only used as a hint. It is not an error to give a wrong size, but it may impair performance.

See also
compute() for an example.

Definition at line 85 of file RealSchur.h.

86  : m_matT(size, size),
87  m_matU(size, size),
89  m_hess(size),
90  m_isInitialized(false),
91  m_matUisUptodate(false),
92  m_maxIters(-1)
93  { }
bool m_matUisUptodate
Definition: RealSchur.h:235
MatrixType m_matU
Definition: RealSchur.h:230
ColumnVectorType m_workspaceVector
Definition: RealSchur.h:231
MatrixType m_matT
Definition: RealSchur.h:229
HessenbergDecomposition< MatrixType > m_hess
Definition: RealSchur.h:232
bool m_isInitialized
Definition: RealSchur.h:234

◆ RealSchur() [2/2]

template<typename MatrixType_ >
template<typename InputType >
Eigen::RealSchur< MatrixType_ >::RealSchur ( const EigenBase< InputType > &  matrix,
bool  computeU = true 
)
inlineexplicit

Constructor; computes real Schur decomposition of given matrix.

Parameters
[in]matrixSquare matrix whose Schur decomposition is to be computed.
[in]computeUIf true, both T and U are computed; if false, only T is computed.

This constructor calls compute() to compute the Schur decomposition.

Example:

cout << "Here is a random 6x6 matrix, A:" << endl << A << endl << endl;
RealSchur<MatrixXd> schur(A);
cout << "The orthogonal matrix U is:" << endl << schur.matrixU() << endl;
cout << "The quasi-triangular matrix T is:" << endl << schur.matrixT() << endl << endl;
MatrixXd U = schur.matrixU();
MatrixXd T = schur.matrixT();
cout << "U * T * U^T = " << endl << U * T * U.transpose() << endl;
MatrixXcf A
ComplexSchur< MatrixXcf > schur(4)
static const RandomReturnType Random()
Definition: Random.h:114
Matrix< double, Dynamic, Dynamic > MatrixXd
Dynamic×Dynamic matrix of type double.
Definition: Matrix.h:502

Output:

Here is a random 6x6 matrix, A:
   0.68   -0.33   -0.27  -0.717  -0.687  0.0259
 -0.211   0.536  0.0268   0.214  -0.198   0.678
  0.566  -0.444   0.904  -0.967   -0.74   0.225
  0.597   0.108   0.832  -0.514  -0.782  -0.408
  0.823 -0.0452   0.271  -0.726   0.998   0.275
 -0.605   0.258   0.435   0.608  -0.563  0.0486

The orthogonal matrix U is:
  0.348  -0.754 0.00435  -0.351  0.0146   0.432
  -0.16  -0.266  -0.747   0.457  -0.366  0.0571
  0.505  -0.157  0.0746   0.644   0.518  -0.177
  0.703   0.324  -0.409  -0.349  -0.187  -0.275
  0.296   0.372    0.24   0.324  -0.379   0.684
 -0.126   0.305   -0.46  -0.161   0.647   0.485
The quasi-triangular matrix T is:
   -0.2   -1.83   0.864   0.271    1.09   0.139
  0.647   0.298 -0.0536   0.676  -0.288  0.0231
      0       0   0.967  -0.201  -0.429   0.847
      0       0       0   0.353   0.603   0.694
      0       0       0       0   0.572   -1.03
      0       0       0       0  0.0184   0.664

U * T * U^T = 
   0.68   -0.33   -0.27  -0.717  -0.687  0.0259
 -0.211   0.536  0.0268   0.214  -0.198   0.678
  0.566  -0.444   0.904  -0.967   -0.74   0.225
  0.597   0.108   0.832  -0.514  -0.782  -0.408
  0.823 -0.0452   0.271  -0.726   0.998   0.275
 -0.605   0.258   0.435   0.608  -0.563  0.0486

Definition at line 106 of file RealSchur.h.

107  : m_matT(matrix.rows(),matrix.cols()),
108  m_matU(matrix.rows(),matrix.cols()),
109  m_workspaceVector(matrix.rows()),
110  m_hess(matrix.rows()),
111  m_isInitialized(false),
112  m_matUisUptodate(false),
113  m_maxIters(-1)
114  {
115  compute(matrix.derived(), computeU);
116  }
RealSchur & compute(const EigenBase< InputType > &matrix, bool computeU=true)
Computes Schur decomposition of given matrix.

Member Function Documentation

◆ compute() [1/2]

template<typename MatrixType_ >
template<typename InputType >
RealSchur<MatrixType>& Eigen::RealSchur< MatrixType_ >::compute ( const EigenBase< InputType > &  matrix,
bool  computeU 
)

Definition at line 251 of file RealSchur.h.

252 {
253  const Scalar considerAsZero = (std::numeric_limits<Scalar>::min)();
254 
255  eigen_assert(matrix.cols() == matrix.rows());
256  Index maxIters = m_maxIters;
257  if (maxIters == -1)
258  maxIters = m_maxIterationsPerRow * matrix.rows();
259 
260  Scalar scale = matrix.derived().cwiseAbs().maxCoeff();
261  if(scale<considerAsZero)
262  {
263  m_matT.setZero(matrix.rows(),matrix.cols());
264  if(computeU)
265  m_matU.setIdentity(matrix.rows(),matrix.cols());
266  m_info = Success;
267  m_isInitialized = true;
268  m_matUisUptodate = computeU;
269  return *this;
270  }
271 
272  // Step 1. Reduce to Hessenberg form
273  m_hess.compute(matrix.derived()/scale);
274 
275  // Step 2. Reduce to real Schur form
276  // Note: we copy m_hess.matrixQ() into m_matU here and not in computeFromHessenberg
277  // to be able to pass our working-space buffer for the Householder to Dense evaluation.
278  m_workspaceVector.resize(matrix.cols());
279  if(computeU)
282 
283  m_matT *= scale;
284 
285  return *this;
286 }
#define eigen_assert(x)
Definition: Macros.h:902
HessenbergDecomposition & compute(const EigenBase< InputType > &matrix)
Computes Hessenberg decomposition of given matrix.
MatrixHReturnType matrixH() const
Constructs the Hessenberg matrix H in the decomposition.
HouseholderSequenceType matrixQ() const
Reconstructs the orthogonal matrix Q in the decomposition.
void evalTo(DestType &dst) const
constexpr void resize(Index rows, Index cols)
RealSchur & computeFromHessenberg(const HessMatrixType &matrixH, const OrthMatrixType &matrixQ, bool computeU)
Computes Schur decomposition of a Hessenberg matrix H = Z T Z^T.
MatrixType::Scalar Scalar
Definition: RealSchur.h:67
static const int m_maxIterationsPerRow
Maximum number of iterations per row.
Definition: RealSchur.h:225
ComputationInfo m_info
Definition: RealSchur.h:233
Eigen::Index Index
Definition: RealSchur.h:69
@ Success
Definition: Constants.h:446
bfloat16() min(const bfloat16 &a, const bfloat16 &b)
Definition: BFloat16.h:684

◆ compute() [2/2]

template<typename MatrixType_ >
template<typename InputType >
RealSchur& Eigen::RealSchur< MatrixType_ >::compute ( const EigenBase< InputType > &  matrix,
bool  computeU = true 
)

Computes Schur decomposition of given matrix.

Parameters
[in]matrixSquare matrix whose Schur decomposition is to be computed.
[in]computeUIf true, both T and U are computed; if false, only T is computed.
Returns
Reference to *this

The Schur decomposition is computed by first reducing the matrix to Hessenberg form using the class HessenbergDecomposition. The Hessenberg matrix is then reduced to triangular form by performing Francis QR iterations with implicit double shift. The cost of computing the Schur decomposition depends on the number of iterations; as a rough guide, it may be taken to be \(25n^3\) flops if computeU is true and \(10n^3\) flops if computeU is false.

Example:

RealSchur<MatrixXf> schur(4);
schur.compute(A, /* computeU = */ false);
cout << "The matrix T in the decomposition of A is:" << endl << schur.matrixT() << endl;
schur.compute(A.inverse(), /* computeU = */ false);
cout << "The matrix T in the decomposition of A^(-1) is:" << endl << schur.matrixT() << endl;
Matrix< float, Dynamic, Dynamic > MatrixXf
Dynamic×Dynamic matrix of type float.
Definition: Matrix.h:501

Output:

The matrix T in the decomposition of A is:
 0.523 -0.698  0.148  0.742
 0.475  0.986 -0.793  0.721
     0      0  -0.28  -0.77
     0      0 0.0145 -0.367
The matrix T in the decomposition of A^(-1) is:
-3.06 -4.57 -5.97  5.48
0.168 -2.62 -3.27   3.9
    0     0 0.427 0.573
    0     0 -1.05  1.35
See also
compute(const MatrixType&, bool, Index)

◆ computeFromHessenberg() [1/2]

template<typename MatrixType_ >
template<typename HessMatrixType , typename OrthMatrixType >
RealSchur& Eigen::RealSchur< MatrixType_ >::computeFromHessenberg ( const HessMatrixType &  matrixH,
const OrthMatrixType &  matrixQ,
bool  computeU 
)

Computes Schur decomposition of a Hessenberg matrix H = Z T Z^T.

Parameters
[in]matrixHMatrix in Hessenberg form H
[in]matrixQorthogonal matrix Q that transform a matrix A to H : A = Q H Q^T
computeUComputes the matriX U of the Schur vectors
Returns
Reference to *this

This routine assumes that the matrix is already reduced in Hessenberg form matrixH using either the class HessenbergDecomposition or another mean. It computes the upper quasi-triangular matrix T of the Schur decomposition of H When computeU is true, this routine computes the matrix U such that A = U T U^T = (QZ) T (QZ)^T = Q H Q^T where A is the initial matrix

NOTE Q is referenced if computeU is true; so, if the initial orthogonal matrix is not available, the user should give an identity matrix (Q.setIdentity())

See also
compute(const MatrixType&, bool)

◆ computeFromHessenberg() [2/2]

template<typename MatrixType_ >
template<typename HessMatrixType , typename OrthMatrixType >
RealSchur<MatrixType>& Eigen::RealSchur< MatrixType_ >::computeFromHessenberg ( const HessMatrixType &  matrixH,
const OrthMatrixType &  matrixQ,
bool  computeU 
)

Definition at line 289 of file RealSchur.h.

290 {
291  using std::abs;
292 
293  m_matT = matrixH;
295  if(computeU && !internal::is_same_dense(m_matU,matrixQ))
296  m_matU = matrixQ;
297 
298  Index maxIters = m_maxIters;
299  if (maxIters == -1)
300  maxIters = m_maxIterationsPerRow * matrixH.rows();
301  Scalar* workspace = &m_workspaceVector.coeffRef(0);
302 
303  // The matrix m_matT is divided in three parts.
304  // Rows 0,...,il-1 are decoupled from the rest because m_matT(il,il-1) is zero.
305  // Rows il,...,iu is the part we are working on (the active window).
306  // Rows iu+1,...,end are already brought in triangular form.
307  Index iu = m_matT.cols() - 1;
308  Index iter = 0; // iteration count for current eigenvalue
309  Index totalIter = 0; // iteration count for whole matrix
310  Scalar exshift(0); // sum of exceptional shifts
311  Scalar norm = computeNormOfT();
312  // sub-diagonal entries smaller than considerAsZero will be treated as zero.
313  // We use eps^2 to enable more precision in small eigenvalues.
314  Scalar considerAsZero = numext::maxi<Scalar>( norm * numext::abs2(NumTraits<Scalar>::epsilon()),
316 
317  if(!numext::is_exactly_zero(norm))
318  {
319  while (iu >= 0)
320  {
321  Index il = findSmallSubdiagEntry(iu,considerAsZero);
322 
323  // Check for convergence
324  if (il == iu) // One root found
325  {
326  m_matT.coeffRef(iu,iu) = m_matT.coeff(iu,iu) + exshift;
327  if (iu > 0)
328  m_matT.coeffRef(iu, iu-1) = Scalar(0);
329  iu--;
330  iter = 0;
331  }
332  else if (il == iu-1) // Two roots found
333  {
334  splitOffTwoRows(iu, computeU, exshift);
335  iu -= 2;
336  iter = 0;
337  }
338  else // No convergence yet
339  {
340  // The firstHouseholderVector vector has to be initialized to something to get rid of a silly GCC warning (-O1 -Wall -DNDEBUG )
341  Vector3s firstHouseholderVector = Vector3s::Zero(), shiftInfo;
342  computeShift(iu, iter, exshift, shiftInfo);
343  iter = iter + 1;
344  totalIter = totalIter + 1;
345  if (totalIter > maxIters) break;
346  Index im;
347  initFrancisQRStep(il, iu, shiftInfo, im, firstHouseholderVector);
348  performFrancisQRStep(il, im, iu, computeU, firstHouseholderVector, workspace);
349  }
350  }
351  }
352  if(totalIter <= maxIters)
353  m_info = Success;
354  else
356 
357  m_isInitialized = true;
358  m_matUisUptodate = computeU;
359  return *this;
360 }
const AbsReturnType abs() const
constexpr Scalar & coeffRef(Index rowId, Index colId)
Matrix< Scalar, 3, 1 > Vector3s
Definition: RealSchur.h:238
void initFrancisQRStep(Index il, Index iu, const Vector3s &shiftInfo, Index &im, Vector3s &firstHouseholderVector)
Definition: RealSchur.h:474
void computeShift(Index iu, Index iter, Scalar &exshift, Vector3s &shiftInfo)
Definition: RealSchur.h:432
void splitOffTwoRows(Index iu, bool computeU, const Scalar &exshift)
Definition: RealSchur.h:397
Scalar computeNormOfT()
Definition: RealSchur.h:364
Index findSmallSubdiagEntry(Index iu, const Scalar &considerAsZero)
Definition: RealSchur.h:378
void performFrancisQRStep(Index il, Index im, Index iu, bool computeU, const Vector3s &firstHouseholderVector, Scalar *workspace)
Definition: RealSchur.h:499
@ NoConvergence
Definition: Constants.h:450
bool is_same_dense(const T1 &mat1, const T2 &mat2, std::enable_if_t< possibly_same_dense< T1, T2 >::value > *=0)
Definition: XprHelper.h:745
bool abs2(bool x)
bool is_exactly_zero(const X &x)
Definition: Meta.h:475

◆ computeNormOfT()

template<typename MatrixType >
MatrixType::Scalar Eigen::RealSchur< MatrixType >::computeNormOfT
inlineprivate

Definition at line 364 of file RealSchur.h.

365 {
366  const Index size = m_matT.cols();
367  // FIXME to be efficient the following would requires a triangular reduxion code
368  // Scalar norm = m_matT.upper().cwiseAbs().sum()
369  // + m_matT.bottomLeftCorner(size-1,size-1).diagonal().cwiseAbs().sum();
370  Scalar norm(0);
371  for (Index j = 0; j < size; ++j)
372  norm += m_matT.col(j).segment(0, (std::min)(size,j+2)).cwiseAbs().sum();
373  return norm;
374 }
std::ptrdiff_t j

◆ computeShift()

template<typename MatrixType >
void Eigen::RealSchur< MatrixType >::computeShift ( Index  iu,
Index  iter,
Scalar exshift,
Vector3s shiftInfo 
)
inlineprivate

Definition at line 432 of file RealSchur.h.

433 {
434  using std::sqrt;
435  using std::abs;
436  shiftInfo.coeffRef(0) = m_matT.coeff(iu,iu);
437  shiftInfo.coeffRef(1) = m_matT.coeff(iu-1,iu-1);
438  shiftInfo.coeffRef(2) = m_matT.coeff(iu,iu-1) * m_matT.coeff(iu-1,iu);
439 
440  // Wilkinson's original ad hoc shift
441  if (iter == 10)
442  {
443  exshift += shiftInfo.coeff(0);
444  for (Index i = 0; i <= iu; ++i)
445  m_matT.coeffRef(i,i) -= shiftInfo.coeff(0);
446  Scalar s = abs(m_matT.coeff(iu,iu-1)) + abs(m_matT.coeff(iu-1,iu-2));
447  shiftInfo.coeffRef(0) = Scalar(0.75) * s;
448  shiftInfo.coeffRef(1) = Scalar(0.75) * s;
449  shiftInfo.coeffRef(2) = Scalar(-0.4375) * s * s;
450  }
451 
452  // MATLAB's new ad hoc shift
453  if (iter == 30)
454  {
455  Scalar s = (shiftInfo.coeff(1) - shiftInfo.coeff(0)) / Scalar(2.0);
456  s = s * s + shiftInfo.coeff(2);
457  if (s > Scalar(0))
458  {
459  s = sqrt(s);
460  if (shiftInfo.coeff(1) < shiftInfo.coeff(0))
461  s = -s;
462  s = s + (shiftInfo.coeff(1) - shiftInfo.coeff(0)) / Scalar(2.0);
463  s = shiftInfo.coeff(0) - shiftInfo.coeff(2) / s;
464  exshift += s;
465  for (Index i = 0; i <= iu; ++i)
466  m_matT.coeffRef(i,i) -= s;
467  shiftInfo.setConstant(Scalar(0.964));
468  }
469  }
470 }
const SqrtReturnType sqrt() const
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_abs_op< typename Derived::Scalar >, const Derived > abs(const Eigen::ArrayBase< Derived > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_sqrt_op< typename Derived::Scalar >, const Derived > sqrt(const Eigen::ArrayBase< Derived > &x)

◆ findSmallSubdiagEntry()

template<typename MatrixType >
Index Eigen::RealSchur< MatrixType >::findSmallSubdiagEntry ( Index  iu,
const Scalar considerAsZero 
)
inlineprivate

Definition at line 378 of file RealSchur.h.

379 {
380  using std::abs;
381  Index res = iu;
382  while (res > 0)
383  {
384  Scalar s = abs(m_matT.coeff(res-1,res-1)) + abs(m_matT.coeff(res,res));
385 
386  s = numext::maxi<Scalar>(s * NumTraits<Scalar>::epsilon(), considerAsZero);
387 
388  if (abs(m_matT.coeff(res,res-1)) <= s)
389  break;
390  res--;
391  }
392  return res;
393 }
cout<< "Here is the matrix m:"<< endl<< m<< endl;Matrix< ptrdiff_t, 3, 1 > res

◆ getMaxIterations()

template<typename MatrixType_ >
Index Eigen::RealSchur< MatrixType_ >::getMaxIterations ( )
inline

Returns the maximum number of iterations.

Definition at line 215 of file RealSchur.h.

216  {
217  return m_maxIters;
218  }

◆ info()

template<typename MatrixType_ >
ComputationInfo Eigen::RealSchur< MatrixType_ >::info ( ) const
inline

Reports whether previous computation was successful.

Returns
Success if computation was successful, NoConvergence otherwise.

Definition at line 197 of file RealSchur.h.

198  {
199  eigen_assert(m_isInitialized && "RealSchur is not initialized.");
200  return m_info;
201  }

◆ initFrancisQRStep()

template<typename MatrixType >
void Eigen::RealSchur< MatrixType >::initFrancisQRStep ( Index  il,
Index  iu,
const Vector3s shiftInfo,
Index im,
Vector3s firstHouseholderVector 
)
inlineprivate

Definition at line 474 of file RealSchur.h.

475 {
476  using std::abs;
477  Vector3s& v = firstHouseholderVector; // alias to save typing
478 
479  for (im = iu-2; im >= il; --im)
480  {
481  const Scalar Tmm = m_matT.coeff(im,im);
482  const Scalar r = shiftInfo.coeff(0) - Tmm;
483  const Scalar s = shiftInfo.coeff(1) - Tmm;
484  v.coeffRef(0) = (r * s - shiftInfo.coeff(2)) / m_matT.coeff(im+1,im) + m_matT.coeff(im,im+1);
485  v.coeffRef(1) = m_matT.coeff(im+1,im+1) - Tmm - r - s;
486  v.coeffRef(2) = m_matT.coeff(im+2,im+1);
487  if (im == il) {
488  break;
489  }
490  const Scalar lhs = m_matT.coeff(im,im-1) * (abs(v.coeff(1)) + abs(v.coeff(2)));
491  const Scalar rhs = v.coeff(0) * (abs(m_matT.coeff(im-1,im-1)) + abs(Tmm) + abs(m_matT.coeff(im+1,im+1)));
492  if (abs(lhs) < NumTraits<Scalar>::epsilon() * rhs)
493  break;
494  }
495 }
Array< int, Dynamic, 1 > v

◆ matrixT()

template<typename MatrixType_ >
const MatrixType& Eigen::RealSchur< MatrixType_ >::matrixT ( ) const
inline

Returns the quasi-triangular matrix in the Schur decomposition.

Returns
A const reference to the matrix T.
Precondition
Either the constructor RealSchur(const MatrixType&, bool) or the member function compute(const MatrixType&, bool) has been called before to compute the Schur decomposition of a matrix.
See also
RealSchur(const MatrixType&, bool) for an example

Definition at line 146 of file RealSchur.h.

147  {
148  eigen_assert(m_isInitialized && "RealSchur is not initialized.");
149  return m_matT;
150  }

◆ matrixU()

template<typename MatrixType_ >
const MatrixType& Eigen::RealSchur< MatrixType_ >::matrixU ( ) const
inline

Returns the orthogonal matrix in the Schur decomposition.

Returns
A const reference to the matrix U.
Precondition
Either the constructor RealSchur(const MatrixType&, bool) or the member function compute(const MatrixType&, bool) has been called before to compute the Schur decomposition of a matrix, and computeU was set to true (the default value).
See also
RealSchur(const MatrixType&, bool) for an example

Definition at line 129 of file RealSchur.h.

130  {
131  eigen_assert(m_isInitialized && "RealSchur is not initialized.");
132  eigen_assert(m_matUisUptodate && "The matrix U has not been computed during the RealSchur decomposition.");
133  return m_matU;
134  }

◆ performFrancisQRStep()

template<typename MatrixType >
void Eigen::RealSchur< MatrixType >::performFrancisQRStep ( Index  il,
Index  im,
Index  iu,
bool  computeU,
const Vector3s firstHouseholderVector,
Scalar workspace 
)
inlineprivate

Definition at line 499 of file RealSchur.h.

500 {
501  eigen_assert(im >= il);
502  eigen_assert(im <= iu-2);
503 
504  const Index size = m_matT.cols();
505 
506  for (Index k = im; k <= iu-2; ++k)
507  {
508  bool firstIteration = (k == im);
509 
510  Vector3s v;
511  if (firstIteration)
512  v = firstHouseholderVector;
513  else
514  v = m_matT.template block<3,1>(k,k-1);
515 
516  Scalar tau, beta;
517  Matrix<Scalar, 2, 1> ess;
518  v.makeHouseholder(ess, tau, beta);
519 
520  if (!numext::is_exactly_zero(beta)) // if v is not zero
521  {
522  if (firstIteration && k > il)
523  m_matT.coeffRef(k,k-1) = -m_matT.coeff(k,k-1);
524  else if (!firstIteration)
525  m_matT.coeffRef(k,k-1) = beta;
526 
527  // These Householder transformations form the O(n^3) part of the algorithm
528  m_matT.block(k, k, 3, size-k).applyHouseholderOnTheLeft(ess, tau, workspace);
529  m_matT.block(0, k, (std::min)(iu,k+3) + 1, 3).applyHouseholderOnTheRight(ess, tau, workspace);
530  if (computeU)
531  m_matU.block(0, k, size, 3).applyHouseholderOnTheRight(ess, tau, workspace);
532  }
533  }
534 
535  Matrix<Scalar, 2, 1> v = m_matT.template block<2,1>(iu-1, iu-2);
536  Scalar tau, beta;
537  Matrix<Scalar, 1, 1> ess;
538  v.makeHouseholder(ess, tau, beta);
539 
540  if (!numext::is_exactly_zero(beta)) // if v is not zero
541  {
542  m_matT.coeffRef(iu-1, iu-2) = beta;
543  m_matT.block(iu-1, iu-1, 2, size-iu+1).applyHouseholderOnTheLeft(ess, tau, workspace);
544  m_matT.block(0, iu-1, iu+1, 2).applyHouseholderOnTheRight(ess, tau, workspace);
545  if (computeU)
546  m_matU.block(0, iu-1, size, 2).applyHouseholderOnTheRight(ess, tau, workspace);
547  }
548 
549  // clean up pollution due to round-off errors
550  for (Index i = im+2; i <= iu; ++i)
551  {
552  m_matT.coeffRef(i,i-2) = Scalar(0);
553  if (i > im+2)
554  m_matT.coeffRef(i,i-3) = Scalar(0);
555  }
556 }

◆ setMaxIterations()

template<typename MatrixType_ >
RealSchur& Eigen::RealSchur< MatrixType_ >::setMaxIterations ( Index  maxIters)
inline

Sets the maximum number of iterations allowed.

If not specified by the user, the maximum number of iterations is m_maxIterationsPerRow times the size of the matrix.

Definition at line 208 of file RealSchur.h.

209  {
210  m_maxIters = maxIters;
211  return *this;
212  }

◆ splitOffTwoRows()

template<typename MatrixType >
void Eigen::RealSchur< MatrixType >::splitOffTwoRows ( Index  iu,
bool  computeU,
const Scalar exshift 
)
inlineprivate

Definition at line 397 of file RealSchur.h.

398 {
399  using std::sqrt;
400  using std::abs;
401  const Index size = m_matT.cols();
402 
403  // The eigenvalues of the 2x2 matrix [a b; c d] are
404  // trace +/- sqrt(discr/4) where discr = tr^2 - 4*det, tr = a + d, det = ad - bc
405  Scalar p = Scalar(0.5) * (m_matT.coeff(iu-1,iu-1) - m_matT.coeff(iu,iu));
406  Scalar q = p * p + m_matT.coeff(iu,iu-1) * m_matT.coeff(iu-1,iu); // q = tr^2 / 4 - det = discr/4
407  m_matT.coeffRef(iu,iu) += exshift;
408  m_matT.coeffRef(iu-1,iu-1) += exshift;
409 
410  if (q >= Scalar(0)) // Two real eigenvalues
411  {
412  Scalar z = sqrt(abs(q));
413  JacobiRotation<Scalar> rot;
414  if (p >= Scalar(0))
415  rot.makeGivens(p + z, m_matT.coeff(iu, iu-1));
416  else
417  rot.makeGivens(p - z, m_matT.coeff(iu, iu-1));
418 
419  m_matT.rightCols(size-iu+1).applyOnTheLeft(iu-1, iu, rot.adjoint());
420  m_matT.topRows(iu+1).applyOnTheRight(iu-1, iu, rot);
421  m_matT.coeffRef(iu, iu-1) = Scalar(0);
422  if (computeU)
423  m_matU.applyOnTheRight(iu-1, iu, rot);
424  }
425 
426  if (iu > 1)
427  m_matT.coeffRef(iu-1, iu-2) = Scalar(0);
428 }
float * p

Member Data Documentation

◆ m_hess

template<typename MatrixType_ >
HessenbergDecomposition<MatrixType> Eigen::RealSchur< MatrixType_ >::m_hess
private

Definition at line 232 of file RealSchur.h.

◆ m_info

template<typename MatrixType_ >
ComputationInfo Eigen::RealSchur< MatrixType_ >::m_info
private

Definition at line 233 of file RealSchur.h.

◆ m_isInitialized

template<typename MatrixType_ >
bool Eigen::RealSchur< MatrixType_ >::m_isInitialized
private

Definition at line 234 of file RealSchur.h.

◆ m_matT

template<typename MatrixType_ >
MatrixType Eigen::RealSchur< MatrixType_ >::m_matT
private

Definition at line 229 of file RealSchur.h.

◆ m_matU

template<typename MatrixType_ >
MatrixType Eigen::RealSchur< MatrixType_ >::m_matU
private

Definition at line 230 of file RealSchur.h.

◆ m_matUisUptodate

template<typename MatrixType_ >
bool Eigen::RealSchur< MatrixType_ >::m_matUisUptodate
private

Definition at line 235 of file RealSchur.h.

◆ m_maxIterationsPerRow

template<typename MatrixType_ >
const int Eigen::RealSchur< MatrixType_ >::m_maxIterationsPerRow
static

Maximum number of iterations per row.

If not otherwise specified, the maximum number of iterations is this number times the size of the matrix. It is currently set to 40.

Definition at line 225 of file RealSchur.h.

◆ m_maxIters

template<typename MatrixType_ >
Index Eigen::RealSchur< MatrixType_ >::m_maxIters
private

Definition at line 236 of file RealSchur.h.

◆ m_workspaceVector

template<typename MatrixType_ >
ColumnVectorType Eigen::RealSchur< MatrixType_ >::m_workspaceVector
private

Definition at line 231 of file RealSchur.h.


The documentation for this class was generated from the following file: