QEPSortDenseSchur

Reorders the Schur decomposition computed by EPSDenseSchur().

Synopsis

#include "slepcqep.h" 
PetscErrorCode QEPSortDenseSchur(QEP qep,PetscInt n_,PetscInt k,PetscScalar *T,PetscInt ldt_,PetscScalar *Q,PetscScalar *wr,PetscScalar *wi)
Not Collective

Input Parameters

qep - the qep solver context
n - dimension of the matrix
k - first active column
ldt - leading dimension of T

Input/Output Parameters

T - the upper (quasi-)triangular matrix
Q - the orthogonal matrix of Schur vectors
wr - pointer to the array to store the computed eigenvalues
wi - imaginary part of the eigenvalues (only when using real numbers)

Notes

This function reorders the eigenvalues in wr,wi located in positions k to n according to the sort order specified in EPSetWhicheigenpairs. The Schur decomposition Z*T*Z^T, is also reordered by means of rotations so that eigenvalues in the diagonal blocks of T follow the same order.

Both T and Q are overwritten.

This routine uses LAPACK routines xTREXC.

See Also

EPSDenseHessenberg(), EPSDenseSchur(), EPSDenseTridiagonal()

Location: src/qep/interface/qepdense.c
Index of all QEP routines
Table of Contents for all manual pages
Index of all manual pages