EPSSortDenseSchurGeneralized

Reorders a generalized Schur decomposition.

Synopsis

#include "slepceps.h" 
PetscErrorCode EPSSortDenseSchurGeneralized(EPS eps,PetscInt n_,PetscInt k0,PetscInt k1,PetscScalar *T,PetscScalar *S,PetscInt ldt_,PetscScalar *Q,PetscScalar *Z,PetscScalar *wr,PetscScalar *wi)
Not Collective

Input Parameters

eps - the eigensolver context
n - dimension of the matrix
k0 - first active column
k1 - last column to be ordered
ldt - leading dimension of T

Input/Output Parameters

T,S - the upper (quasi-)triangular matrices
Q,Z - 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 k0 to n according to the sort order specified in EPSetWhicheigenpairs. The selection sort is the method used to sort the eigenvalues, and it stops when the column k1-1 is ordered. 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.

T,S,Q and Z are overwritten.

This routine uses LAPACK routines xTGEXC.

See Also

EPSSortDenseSchur(), EPSDenseHessenberg(), EPSDenseSchur(), EPSDenseTridiagonal()

Location: src/eps/interface/dense.c
Index of all EPS routines
Table of Contents for all manual pages
Index of all manual pages