Eigen::ComplexSchur< MatrixType_ > Class Template Reference

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

Public Types

enum  {
  RowsAtCompileTime ,
  ColsAtCompileTime ,
  Options ,
  MaxRowsAtCompileTime ,
  MaxColsAtCompileTime
}
 
typedef Matrix< ComplexScalar, RowsAtCompileTime, ColsAtCompileTime, Options, MaxRowsAtCompileTime, MaxColsAtCompileTimeComplexMatrixType
 Type for the matrices in the Schur decomposition. More...
 
typedef std::complex< RealScalarComplexScalar
 Complex scalar type for MatrixType_. More...
 
typedef Eigen::Index Index
 
typedef MatrixType_ MatrixType
 
typedef NumTraits< Scalar >::Real RealScalar
 
typedef MatrixType::Scalar Scalar
 Scalar type for matrices of type MatrixType_. More...
 

Public Member Functions

template<typename InputType >
 ComplexSchur (const EigenBase< InputType > &matrix, bool computeU=true)
 Constructor; computes Schur decomposition of given matrix. More...
 
 ComplexSchur (Index size=RowsAtCompileTime==Dynamic ? 1 :RowsAtCompileTime)
 Default constructor. More...
 
template<typename InputType >
ComplexSchur< MatrixType > & compute (const EigenBase< InputType > &matrix, bool computeU)
 
template<typename InputType >
ComplexSchurcompute (const EigenBase< InputType > &matrix, bool computeU=true)
 Computes Schur decomposition of given matrix. More...
 
template<typename HessMatrixType , typename OrthMatrixType >
ComplexSchur< MatrixType > & computeFromHessenberg (const HessMatrixType &matrixH, const OrthMatrixType &matrixQ, bool computeU)
 
template<typename HessMatrixType , typename OrthMatrixType >
ComplexSchurcomputeFromHessenberg (const HessMatrixType &matrixH, const OrthMatrixType &matrixQ, bool computeU=true)
 Compute Schur decomposition from a given Hessenberg matrix. More...
 
Index getMaxIterations ()
 Returns the maximum number of iterations. More...
 
ComputationInfo info () const
 Reports whether previous computation was successful. More...
 
const ComplexMatrixTypematrixT () const
 Returns the triangular matrix in the Schur decomposition. More...
 
const ComplexMatrixTypematrixU () const
 Returns the unitary matrix in the Schur decomposition. More...
 
ComplexSchursetMaxIterations (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...
 

Protected Attributes

HessenbergDecomposition< MatrixTypem_hess
 
ComputationInfo m_info
 
bool m_isInitialized
 
ComplexMatrixType m_matT
 
ComplexMatrixType m_matU
 
bool m_matUisUptodate
 
Index m_maxIters
 

Private Member Functions

ComplexScalar computeShift (Index iu, Index iter)
 
void reduceToTriangularForm (bool computeU)
 
bool subdiagonalEntryIsNeglegible (Index i)
 

Detailed Description

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

Performs a complex Schur decomposition of a real or complex square matrix.

This is defined in the Eigenvalues module.

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

Given a real or complex square matrix A, this class computes the Schur decomposition: \( A = U T U^*\) where U is a unitary complex matrix, and T is a complex upper triangular matrix. The diagonal of the matrix T corresponds to the eigenvalues of the matrix A.

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

Note
This code is inspired from Jampack
See also
class RealSchur, class EigenSolver, class ComplexEigenSolver

Definition at line 53 of file ComplexSchur.h.

Member Typedef Documentation

◆ ComplexMatrixType

Type for the matrices in the Schur decomposition.

This is a square matrix with entries of type ComplexScalar. The size is the same as the size of MatrixType_.

Definition at line 83 of file ComplexSchur.h.

◆ ComplexScalar

template<typename MatrixType_ >
typedef std::complex<RealScalar> Eigen::ComplexSchur< MatrixType_ >::ComplexScalar

Complex scalar type for MatrixType_.

This is std::complex<Scalar> if Scalar is real (e.g., float or double) and just Scalar if Scalar is complex.

Definition at line 76 of file ComplexSchur.h.

◆ Index

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

Definition at line 68 of file ComplexSchur.h.

◆ MatrixType

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

Definition at line 56 of file ComplexSchur.h.

◆ RealScalar

template<typename MatrixType_ >
typedef NumTraits<Scalar>::Real Eigen::ComplexSchur< MatrixType_ >::RealScalar

Definition at line 67 of file ComplexSchur.h.

◆ Scalar

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

Scalar type for matrices of type MatrixType_.

Definition at line 66 of file ComplexSchur.h.

Member Enumeration Documentation

◆ anonymous enum

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

Definition at line 57 of file ComplexSchur.h.

57  {
58  RowsAtCompileTime = MatrixType::RowsAtCompileTime,
59  ColsAtCompileTime = MatrixType::ColsAtCompileTime,
60  Options = MatrixType::Options,
61  MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
62  MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
63  };

Constructor & Destructor Documentation

◆ ComplexSchur() [1/2]

template<typename MatrixType_ >
Eigen::ComplexSchur< MatrixType_ >::ComplexSchur ( 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 96 of file ComplexSchur.h.

97  : m_matT(size,size),
98  m_matU(size,size),
99  m_hess(size),
100  m_isInitialized(false),
101  m_matUisUptodate(false),
102  m_maxIters(-1)
103  {}
HessenbergDecomposition< MatrixType > m_hess
Definition: ComplexSchur.h:251
ComplexMatrixType m_matU
Definition: ComplexSchur.h:250
ComplexMatrixType m_matT
Definition: ComplexSchur.h:250

◆ ComplexSchur() [2/2]

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

Constructor; 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.

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

See also
matrixT() and matrixU() for examples.

Definition at line 115 of file ComplexSchur.h.

116  : m_matT(matrix.rows(),matrix.cols()),
117  m_matU(matrix.rows(),matrix.cols()),
118  m_hess(matrix.rows()),
119  m_isInitialized(false),
120  m_matUisUptodate(false),
121  m_maxIters(-1)
122  {
123  compute(matrix.derived(), computeU);
124  }
ComplexSchur & 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 >
ComplexSchur<MatrixType>& Eigen::ComplexSchur< MatrixType_ >::compute ( const EigenBase< InputType > &  matrix,
bool  computeU 
)

Definition at line 324 of file ComplexSchur.h.

325 {
326  m_matUisUptodate = false;
327  eigen_assert(matrix.cols() == matrix.rows());
328 
329  if(matrix.cols() == 1)
330  {
331  m_matT = matrix.derived().template cast<ComplexScalar>();
332  if(computeU) m_matU = ComplexMatrixType::Identity(1,1);
333  m_info = Success;
334  m_isInitialized = true;
335  m_matUisUptodate = computeU;
336  return *this;
337  }
338 
339  internal::complex_schur_reduce_to_hessenberg<MatrixType, NumTraits<Scalar>::IsComplex>::run(*this, matrix.derived(), computeU);
340  computeFromHessenberg(m_matT, m_matU, computeU);
341  return *this;
342 }
#define eigen_assert(x)
Definition: Macros.h:902
ComplexSchur & computeFromHessenberg(const HessMatrixType &matrixH, const OrthMatrixType &matrixQ, bool computeU=true)
Compute Schur decomposition from a given Hessenberg matrix.
ComputationInfo m_info
Definition: ComplexSchur.h:252
@ Success
Definition: Constants.h:446

◆ compute() [2/2]

template<typename MatrixType_ >
template<typename InputType >
ComplexSchur& Eigen::ComplexSchur< 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 QR iterations with a single shift. The cost of computing the Schur decomposition depends on the number of iterations; as a rough guide, it may be taken on the number of iterations; as a rough guide, it may be taken to be \(25n^3\) complex flops, or \(10n^3\) complex flops if computeU is false.

Example:

ComplexSchur<MatrixXcf> schur(4);
schur.compute(A);
cout << "The matrix T in the decomposition of A is:" << endl << schur.matrixT() << endl;
schur.compute(A.inverse());
cout << "The matrix T in the decomposition of A^(-1) is:" << endl << schur.matrixT() << endl;
MatrixXcf A
ComplexSchur< MatrixXcf > schur(4)
static const RandomReturnType Random()
Definition: Random.h:114
Matrix< std::complex< float >, Dynamic, Dynamic > MatrixXcf
Dynamic×Dynamic matrix of type std::complex<float>.
Definition: Matrix.h:503

Output:

The matrix T in the decomposition of A is:
 (-0.691,-1.63)  (0.763,-0.144) (-0.104,-0.836) (-0.462,-0.378)
          (0,0)   (-0.758,1.22)  (-0.65,-0.772)  (-0.244,0.113)
          (0,0)           (0,0)   (0.137,0.505) (0.0687,-0.404)
          (0,0)           (0,0)           (0,0)   (1.52,-0.402)
The matrix T in the decomposition of A^(-1) is:
    (0.501,-1.84)    (-1.01,-0.984)       (0.636,1.3)    (-0.676,0.352)
            (0,0)   (-0.369,-0.593)     (0.0733,0.18) (-0.0658,-0.0263)
            (0,0)             (0,0)    (-0.222,0.521)    (-0.191,0.121)
            (0,0)             (0,0)             (0,0)     (0.614,0.162)
See also
compute(const MatrixType&, bool, Index)

◆ computeFromHessenberg() [1/2]

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

Definition at line 346 of file ComplexSchur.h.

347 {
348  m_matT = matrixH;
349  if(computeU)
350  m_matU = matrixQ;
351  reduceToTriangularForm(computeU);
352  return *this;
353 }
void reduceToTriangularForm(bool computeU)
Definition: ComplexSchur.h:392

◆ computeFromHessenberg() [2/2]

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

Compute Schur decomposition from a given Hessenberg matrix.

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)

◆ computeShift()

template<typename MatrixType >
ComplexSchur< MatrixType >::ComplexScalar Eigen::ComplexSchur< MatrixType >::computeShift ( Index  iu,
Index  iter 
)
private

Compute the shift in the current QR iteration.

Definition at line 283 of file ComplexSchur.h.

284 {
285  using std::abs;
286  if (iter == 10 || iter == 20)
287  {
288  // exceptional shift, taken from http://www.netlib.org/eispack/comqr.f
289  return abs(numext::real(m_matT.coeff(iu,iu-1))) + abs(numext::real(m_matT.coeff(iu-1,iu-2)));
290  }
291 
292  // compute the shift as one of the eigenvalues of t, the 2x2
293  // diagonal block on the bottom of the active submatrix
294  Matrix<ComplexScalar,2,2> t = m_matT.template block<2,2>(iu-1,iu-1);
295  RealScalar normt = t.cwiseAbs().sum();
296  t /= normt; // the normalization by sf is to avoid under/overflow
297 
298  ComplexScalar b = t.coeff(0,1) * t.coeff(1,0);
299  ComplexScalar c = t.coeff(0,0) - t.coeff(1,1);
300  ComplexScalar disc = sqrt(c*c + RealScalar(4)*b);
301  ComplexScalar det = t.coeff(0,0) * t.coeff(1,1) - b;
302  ComplexScalar trace = t.coeff(0,0) + t.coeff(1,1);
303  ComplexScalar eival1 = (trace + disc) / RealScalar(2);
304  ComplexScalar eival2 = (trace - disc) / RealScalar(2);
305  RealScalar eival1_norm = numext::norm1(eival1);
306  RealScalar eival2_norm = numext::norm1(eival2);
307  // A division by zero can only occur if eival1==eival2==0.
308  // In this case, det==0, and all we have to do is checking that eival2_norm!=0
309  if(eival1_norm > eival2_norm)
310  eival2 = det / eival1;
311  else if(!numext::is_exactly_zero(eival2_norm))
312  eival1 = det / eival2;
313 
314  // choose the eigenvalue closest to the bottom entry of the diagonal
315  if(numext::norm1(eival1-t.coeff(1,1)) < numext::norm1(eival2-t.coeff(1,1)))
316  return normt * eival1;
317  else
318  return normt * eival2;
319 }
const AbsReturnType abs() const
Array< int, 3, 1 > b
RealReturnType real() const
Array33i c
NumTraits< Scalar >::Real RealScalar
Definition: ComplexSchur.h:67
std::complex< RealScalar > ComplexScalar
Complex scalar type for MatrixType_.
Definition: ComplexSchur.h:76
constexpr const Scalar & coeff(Index rowId, Index colId) const
bool is_exactly_zero(const X &x)
Definition: Meta.h:475
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)

◆ getMaxIterations()

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

Returns the maximum number of iterations.

Definition at line 237 of file ComplexSchur.h.

238  {
239  return m_maxIters;
240  }

◆ info()

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

Reports whether previous computation was successful.

Returns
Success if computation was successful, NoConvergence otherwise.

Definition at line 219 of file ComplexSchur.h.

220  {
221  eigen_assert(m_isInitialized && "ComplexSchur is not initialized.");
222  return m_info;
223  }

◆ matrixT()

template<typename MatrixType_ >
const ComplexMatrixType& Eigen::ComplexSchur< MatrixType_ >::matrixT ( ) const
inline

Returns the triangular matrix in the Schur decomposition.

Returns
A const reference to the matrix T.

It is assumed that either the constructor ComplexSchur(const MatrixType& matrix, bool computeU) or the member function compute(const MatrixType& matrix, bool computeU) has been called before to compute the Schur decomposition of a matrix.

Note that this function returns a plain square matrix. If you want to reference only the upper triangular part, use:

schur.matrixT().triangularView<Upper>()
@ Upper
Definition: Constants.h:213

Example:

cout << "Here is a random 4x4 matrix, A:" << endl << A << endl << endl;
ComplexSchur<MatrixXcf> schurOfA(A, false); // false means do not compute U
cout << "The triangular matrix T is:" << endl << schurOfA.matrixT() << endl;
cout<< "Here is a random 4x4 matrix, A:"<< endl<< A<< endl<< endl;ComplexSchur< MatrixXcf > schurOfA(A, false)

Output:

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

The triangular matrix T is:
 (-0.691,-1.63)  (0.763,-0.144) (-0.104,-0.836) (-0.462,-0.378)
          (0,0)   (-0.758,1.22)  (-0.65,-0.772)  (-0.244,0.113)
          (0,0)           (0,0)   (0.137,0.505) (0.0687,-0.404)
          (0,0)           (0,0)           (0,0)   (1.52,-0.402)

Definition at line 164 of file ComplexSchur.h.

165  {
166  eigen_assert(m_isInitialized && "ComplexSchur is not initialized.");
167  return m_matT;
168  }

◆ matrixU()

template<typename MatrixType_ >
const ComplexMatrixType& Eigen::ComplexSchur< MatrixType_ >::matrixU ( ) const
inline

Returns the unitary matrix in the Schur decomposition.

Returns
A const reference to the matrix U.

It is assumed that either the constructor ComplexSchur(const MatrixType& matrix, bool computeU) or the member function compute(const MatrixType& matrix, bool computeU) has been called before to compute the Schur decomposition of a matrix, and that computeU was set to true (the default value).

Example:

cout << "Here is a random 4x4 matrix, A:" << endl << A << endl << endl;
ComplexSchur<MatrixXcf> schurOfA(A);
cout << "The unitary matrix U is:" << endl << schurOfA.matrixU() << endl;

Output:

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

The unitary matrix U is:
 (-0.122,0.271)   (0.354,0.255)    (-0.7,0.321) (0.0909,-0.346)
   (0.247,0.23)  (0.435,-0.395)   (0.184,-0.38)  (0.492,-0.347)
(0.859,-0.0877)  (0.00469,0.21) (-0.256,0.0163)   (0.133,0.355)
 (-0.116,0.195) (-0.484,-0.432)  (-0.183,0.359)   (0.559,0.231)

Definition at line 140 of file ComplexSchur.h.

141  {
142  eigen_assert(m_isInitialized && "ComplexSchur is not initialized.");
143  eigen_assert(m_matUisUptodate && "The matrix U has not been computed during the ComplexSchur decomposition.");
144  return m_matU;
145  }

◆ reduceToTriangularForm()

template<typename MatrixType >
void Eigen::ComplexSchur< MatrixType >::reduceToTriangularForm ( bool  computeU)
private

Definition at line 392 of file ComplexSchur.h.

393 {
394  Index maxIters = m_maxIters;
395  if (maxIters == -1)
396  maxIters = m_maxIterationsPerRow * m_matT.rows();
397 
398  // The matrix m_matT is divided in three parts.
399  // Rows 0,...,il-1 are decoupled from the rest because m_matT(il,il-1) is zero.
400  // Rows il,...,iu is the part we are working on (the active submatrix).
401  // Rows iu+1,...,end are already brought in triangular form.
402  Index iu = m_matT.cols() - 1;
403  Index il;
404  Index iter = 0; // number of iterations we are working on the (iu,iu) element
405  Index totalIter = 0; // number of iterations for whole matrix
406 
407  while(true)
408  {
409  // find iu, the bottom row of the active submatrix
410  while(iu > 0)
411  {
412  if(!subdiagonalEntryIsNeglegible(iu-1)) break;
413  iter = 0;
414  --iu;
415  }
416 
417  // if iu is zero then we are done; the whole matrix is triangularized
418  if(iu==0) break;
419 
420  // if we spent too many iterations, we give up
421  iter++;
422  totalIter++;
423  if(totalIter > maxIters) break;
424 
425  // find il, the top row of the active submatrix
426  il = iu-1;
427  while(il > 0 && !subdiagonalEntryIsNeglegible(il-1))
428  {
429  --il;
430  }
431 
432  /* perform the QR step using Givens rotations. The first rotation
433  creates a bulge; the (il+2,il) element becomes nonzero. This
434  bulge is chased down to the bottom of the active submatrix. */
435 
436  ComplexScalar shift = computeShift(iu, iter);
437  JacobiRotation<ComplexScalar> rot;
438  rot.makeGivens(m_matT.coeff(il,il) - shift, m_matT.coeff(il+1,il));
439  m_matT.rightCols(m_matT.cols()-il).applyOnTheLeft(il, il+1, rot.adjoint());
440  m_matT.topRows((std::min)(il+2,iu)+1).applyOnTheRight(il, il+1, rot);
441  if(computeU) m_matU.applyOnTheRight(il, il+1, rot);
442 
443  for(Index i=il+1 ; i<iu ; i++)
444  {
445  rot.makeGivens(m_matT.coeffRef(i,i-1), m_matT.coeffRef(i+1,i-1), &m_matT.coeffRef(i,i-1));
446  m_matT.coeffRef(i+1,i-1) = ComplexScalar(0);
447  m_matT.rightCols(m_matT.cols()-i).applyOnTheLeft(i, i+1, rot.adjoint());
448  m_matT.topRows((std::min)(i+2,iu)+1).applyOnTheRight(i, i+1, rot);
449  if(computeU) m_matU.applyOnTheRight(i, i+1, rot);
450  }
451  }
452 
453  if(totalIter <= maxIters)
454  m_info = Success;
455  else
457 
458  m_isInitialized = true;
459  m_matUisUptodate = computeU;
460 }
static const int m_maxIterationsPerRow
Maximum number of iterations per row.
Definition: ComplexSchur.h:247
Eigen::Index Index
Definition: ComplexSchur.h:68
bool subdiagonalEntryIsNeglegible(Index i)
Definition: ComplexSchur.h:268
ComplexScalar computeShift(Index iu, Index iter)
Definition: ComplexSchur.h:283
EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
constexpr Scalar & coeffRef(Index rowId, Index colId)
@ NoConvergence
Definition: Constants.h:450
bfloat16() min(const bfloat16 &a, const bfloat16 &b)
Definition: BFloat16.h:684

◆ setMaxIterations()

template<typename MatrixType_ >
ComplexSchur& Eigen::ComplexSchur< 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 230 of file ComplexSchur.h.

231  {
232  m_maxIters = maxIters;
233  return *this;
234  }

◆ subdiagonalEntryIsNeglegible()

template<typename MatrixType >
bool Eigen::ComplexSchur< MatrixType >::subdiagonalEntryIsNeglegible ( Index  i)
inlineprivate

If m_matT(i+1,i) is negligible in floating point arithmetic compared to m_matT(i,i) and m_matT(j,j), then set it to zero and return true, else return false.

Definition at line 268 of file ComplexSchur.h.

269 {
270  RealScalar d = numext::norm1(m_matT.coeff(i,i)) + numext::norm1(m_matT.coeff(i+1,i+1));
271  RealScalar sd = numext::norm1(m_matT.coeff(i+1,i));
272  if (internal::isMuchSmallerThan(sd, d, NumTraits<RealScalar>::epsilon()))
273  {
274  m_matT.coeffRef(i+1,i) = ComplexScalar(0);
275  return true;
276  }
277  return false;
278 }
bool isMuchSmallerThan(const Scalar &x, const OtherScalar &y, const typename NumTraits< Scalar >::Real &precision=NumTraits< Scalar >::dummy_precision())

Member Data Documentation

◆ m_hess

template<typename MatrixType_ >
HessenbergDecomposition<MatrixType> Eigen::ComplexSchur< MatrixType_ >::m_hess
protected

Definition at line 251 of file ComplexSchur.h.

◆ m_info

template<typename MatrixType_ >
ComputationInfo Eigen::ComplexSchur< MatrixType_ >::m_info
protected

Definition at line 252 of file ComplexSchur.h.

◆ m_isInitialized

template<typename MatrixType_ >
bool Eigen::ComplexSchur< MatrixType_ >::m_isInitialized
protected

Definition at line 253 of file ComplexSchur.h.

◆ m_matT

template<typename MatrixType_ >
ComplexMatrixType Eigen::ComplexSchur< MatrixType_ >::m_matT
protected

Definition at line 250 of file ComplexSchur.h.

◆ m_matU

template<typename MatrixType_ >
ComplexMatrixType Eigen::ComplexSchur< MatrixType_ >::m_matU
protected

Definition at line 250 of file ComplexSchur.h.

◆ m_matUisUptodate

template<typename MatrixType_ >
bool Eigen::ComplexSchur< MatrixType_ >::m_matUisUptodate
protected

Definition at line 254 of file ComplexSchur.h.

◆ m_maxIterationsPerRow

template<typename MatrixType_ >
const int Eigen::ComplexSchur< 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 30.

Definition at line 247 of file ComplexSchur.h.

◆ m_maxIters

template<typename MatrixType_ >
Index Eigen::ComplexSchur< MatrixType_ >::m_maxIters
protected

Definition at line 255 of file ComplexSchur.h.


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