StarNEig User's Guide  master branch
A task-based library for solving dense nonsymmetric eigenvalue problems
Generalized eigenvalue problem

Given a matrix pair \((A, B)\), where \(A, B \in \mathbb{R}^{n \times n}\), the generalized eigenvalue problem (GEP) consists of finding generalized eigenvalues \(\lambda_i \in \mathbb C\) and associated generalized eigenvectors \(0 \neq v_i \in \mathbb C^{n}\) such that

\[ A v_{i} = \lambda_{i} B v_{i}, \text{ for } i = 1, 2, \dots, n. \]

The eigenvalues are the \(n\) (potentially complex) roots of the polynomial \(\text{det}(A - \lambda B) = 0\) of degree \(n\). There is often a full set of \(n\) linearly independent generalized eigenvectors, but if there are multiple eigenvalues (i.e., if \(\lambda_{i} = \lambda_{j}\) for some \(i \neq j\)) then there might not be a full set of independent eigenvectors.

At least in principle, a GEP can be transformed into a SEP provided that \(B\) is invertible, since

\[ A v = \lambda B v \Leftrightarrow (B^{-1} A) v = \lambda v. \]

However, in finite precision arithmetic this practice is not recommended.

The library provides 12 interface functions for the generalized case.

Hessenberg-triangular reduction

Given a general matrix pair \((A,B)\), the starneig_GEP_SM_HessenbergTriangular() and starneig_GEP_DM_HessenbergTriangular() interface functions compute a Hessenberg-triangular decomposition

\[ (A,B) = U_1 * (H,R) * U_2^T, \]

where \(H\) is upper Hessenberg, \(R\) is upper triangular, and \(U_1\) and \(U_2\) are orthogonal. On exit, \(A\) is overwritten by \(H\), \(B\) is overwritten by \(R\), and \(Q\) and \(Z\) (which are orthogonal matrices on entry) are overwritten by

\[ Q \gets Q * U_1 \text{ and } Z \gets Z * U_2. \]

This is done in order to greatly accelerate the subsequent computation of a generalized Schur decomposition.

Generalized Schur reduction

Given a Hessenberg-triangular decomposition

\[ (A,B) = Q * (H,R) * Z^T \]

of a general matrix pair \((A,B)\), the starneig_GEP_SM_Schur() and starneig_GEP_DM_Schur() interface functions function compute a generalized Schur decomposition

\[ (A,B) = Q * ( U_1 * (S,T) * U_2^T ) * Z^T, \]

where \(S\) is upper quasi-triangular with \(1 \times 1\) and \(2 \times 2\) blocks on the diagonal, \(T\) is a upper triangular matrix, and \(U_1\) and \(U_2\) are orthogonal. The \(1 \times 1\) blocks on the diagonal of \((S,T)\) correspond to the real generalized eigenvalues and each \(2 \times 2\) block corresponds to a pair of complex conjugate generalized eigenvalues.

On exit, \(H\) is overwritten by \(S\), \(R\) is overwritten by \(T\), and \(Q\) and \(Z\) are overwritten by

\[ Q \gets Q * U_1 \text{ and } Z \gets Z * U_2. \]

The computed generalized eigenvalues are returned as a pair of values \((\alpha,\beta)\) such that \(\alpha/\beta\) gives the actual generalized eigenvalue. The quantity \(\alpha/\beta\) may overflow.

Generalized eigenvalue reordering

Given a generalized Schur decomposition

\[ (A,B) = Q * (S,T) * Z^T \]

of a general matrix pair \((A,B)\) and a selection of generalized eigenvalues, the starneig_GEP_SM_ReorderSchur() and starneig_GEP_DM_ReorderSchur() interface functions attempt to compute an updated generalized Schur decomposition

\[ (A,B) = Q * \left( U_1 * \left ( \begin{bmatrix} \hat S_{11} & \hat S_{12} \\ 0 & \hat S_{22} \end{bmatrix}, \begin{bmatrix} \hat T_{11} & \hat T_{12} \\ 0 & \hat T_{22} \end{bmatrix} \right) * U_2^T \right) * Z^T. \]

where the selected generalized eigenvalues appear in \((\hat S_{11},\hat T_{11})\).

On exit, \(S\) is overwritten by \(\hat{S}\), \(T\) is overwritten by \(\hat{T}\), and \(Q\) and \(Z\) are overwritten by

\[ Q \gets Q * U_1 \text{ and } Z \gets Z * U_2. \]

Reordering may in rare cases fail. In such cases the output is guaranteed to be a Schur decomposition and all (if any) selected generalized eigenvalues that are correctly placed are marked in the selection array on exit.

Reordering may perturb the generalized eigenvalues and the generalized eigenvalues after reordering are returned. The computed generalized eigenvalues are returned as a pair of values \((\alpha,\beta)\) such that \(\alpha/\beta\) gives the actual generalized eigenvalue. The quantity \(\alpha/\beta\) may overflow.

Combined reduction to generalized Schur form and eigenvalue reordering

Given a general matrix pair \((A,B)\), the starneig_GEP_SM_Reduce() and starneig_GEP_DM_Reduce() interface functions compute a (reordered) generalized Schur decomposition

\[ (A,B) = U_1 * (S,T) * U_2^T. \]

Optionally, the interface functions attempt to reorder selected generalized eigenvalues to the top left corner of the generalized Schur decomposition.

On exit, \(A\) is overwritten by \(S\), \(B\) is overwritten by \(T\), and \(Q\) and \(Z\) (which are orthogonal matrices on entry) are overwritten by

\[ Q \gets Q * U_1 \text{ and } Z \gets Z * U_2. \]

The computed generalized eigenvalues are returned as a pair of values \((\alpha,\beta)\) such that \(\alpha/\beta\) gives the actual generalized eigenvalue. The quantity \(\alpha/\beta\) may overflow.

Reordering may in rare cases fail. In such cases the output is guaranteed to be a Schur-triangular decomposition and all (if any) selected generalized eigenvalues that are correctly placed are marked in the selection array on exit.

Generalized eigenvectors

Given a subset consisting of \(m \leq n\) of the eigenvalues \(\lambda_{i}\) for \(i = 1, 2, \ldots, m\) and a generalized Schur decomposition \((A, B) = Q (S, T) Z^{H}\), we can compute for each \(\lambda_{i}\) a generalized eigenvector \(v_{i} \neq 0\) such that \(A v_{i} = \lambda_{i} B v_{i}\) by first computing a generalized eigenvector \(w_{i}\) of \(S - \lambda_{i} T\) and then transform it back to the original basis by pre-multiplication with \(Z\).

Given a generalized Schur decomposition

\[ (A,B) = Q * (S,T) * Z^T \]

of a general matrix pair \((A,B)\) and a selection of generalized eigenvalues, the starneig_GEP_SM_Eigenvectors() and starneig_GEP_DM_Eigenvectors() interface functions compute and return a generalized eigenvector for each of the selected generalized eigenvalues.

The generalized eigenvectors are stored as columns in the output matrix \(X\) in the same order as their corresponding generalized eigenvalues appear in the selection array. A real generalized eigenvector is stored as a single column. The real and imaginary parts of a complex generalized eigenvector are stored as consecutive columns.

For a selected pair of complex conjugate generalized eigenvalues, a generalized eigenvector is computed only for the generalized eigenvalue with positive imaginary part. Thus, every selected generalized eigenvalue contributes one column to the output matrix and thus the number of selected generalized eigenvalues is equal to the number of columns of \(X\).

Eigenvalue selection helper

Given a Schur-triangular matrix pair \((S,T)\) and a predicate function, the starneig_GEP_SM_Select() and starneig_GEP_DM_Select() interface functions conveniently generate a correct selection array and count the number of selected generalized eigenvalues. The count is useful when allocating storage for the generalized eigenvector matrix computed by starneig_GEP_DM_Eigenvectors().

// a predicate function that selects all finite generalized eigenvalues that
// have a real part that is larger than a given value
static int predicate(double real, double imag, double beta, void *arg)
{
double value = * (double *) arg;
if (beta != 0.0 && value < real/beta)
return 1;
return 0;
}
void func(...)
{
...
double value = 0.5;
int num_selected, *selected = malloc(n*sizeof(int));
n, S, ldS, T, ldT, &predicate, &value, selected, &num_selected);
...
}

See modules Shared Memory / Generalized EVP and Distributed Memory / Generalized EVP for further information. See also examples gep_sm_full_chain.c, gep_dm_full_chain.c and gep_sm_eigenvectors.c.

starneig_GEP_SM_Select
starneig_error_t starneig_GEP_SM_Select(int n, double S[], int ldS, double T[], int ldT, int(*predicate)(double real, double imag, double beta, void *arg), void *arg, int selected[], int *num_selected)
Generates a selection array for a Schur-triangular matrix pencil using a user-supplied predicate func...