EPSQRDecomposition

Compute the QR factorization of a set of vectors.

Synopsis

#include "slepceps.h" 
PetscErrorCode EPSQRDecomposition(EPS eps,Vec *V,int m,int n,PetscScalar *R,int ldr)
Collective on EPS

Input Parameters

eps - the eigenproblem solver context
V - set of vectors
m - starting column
n - ending column
ldr - leading dimension of R

Output Parameter

R - triangular matrix of QR factorization

Notes

This routine orthonormalizes the columns of V so that V'*V=I up to working precision. It assumes that the first m columns (from 0 to m-1) are already orthonormal. The coefficients of the orthogonalization are stored in R so that V*R is equal to the original V.

In some cases, this routine makes V B-orthonormal, that is, V'*B*V=I.

If one of the vectors is linearly dependent wrt the rest (at working precision) then it is replaced by a random vector.

See Also

EPSOrthogonalize(), STNorm(), STInnerProduct().

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