EPSDenseSchur

Computes the upper (quasi-)triangular form of a dense upper Hessenberg matrix.

Synopsis

#include "slepceps.h" 
PetscErrorCode EPSDenseSchur(int n,int k,PetscScalar *H,int ldh,PetscScalar *Z,PetscScalar *wr,PetscScalar *wi)
Not Collective

Input Parameters

n - dimension of the matrix
k - first active column
ldh - leading dimension of H

Input/Output Parameters

H - on entry, the upper Hessenber matrix; on exit, the upper (quasi-)triangular matrix (T)
Z - on entry, initial transformation matrix; on exit, orthogonal matrix of Schur vectors

Output Parameters

wr - pointer to the array to store the computed eigenvalues
wi - imaginary part of the eigenvalues (only when using real numbers)

Notes

This function computes the (real) Schur decomposition of an upper Hessenberg matrix H: H*Z = Z*T, where T is an upper (quasi-)triangular matrix (returned in H), and Z is the orthogonal matrix of Schur vectors. Eigenvalues are extracted from the diagonal blocks of T and returned in wr,wi. Transformations are accumulated in Z so that on entry it can contain the transformation matrix associated to the Hessenberg reduction.

Only active columns (from k to n) are computed.

Both H and Z are overwritten.

This routine uses LAPACK routines xHSEQR.

See Also

EPSDenseHessenberg(), EPSSortDenseSchur(), EPSDenseTridiagonal()

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