BVSVDAndRank#

Compute the SVD (left singular vectors only, and singular values) and determine the numerical rank according to a tolerance.

Synopsis#

#include "slepcbv.h" 
PetscErrorCode BVSVDAndRank(BV S,PetscInt m,PetscInt l,PetscReal delta,BVSVDMethod meth,PetscScalar A[],PetscReal sigma[],PetscInt *rank)

Collective

Input Parameters#

  • S - the basis vectors object

  • m - the moment degree

  • l - the block size

  • delta - the tolerance used to determine the rank

  • meth - the method to be used

Output Parameters#

  • A - workspace, on output contains relevant values in the CAA method

  • sigma - computed singular values

  • rank - estimated rank (optional, pass NULL if not needed)

Notes#

This function computes the SVD \(S=U\Sigma V^*\) and replaces \(S\) with \(U\). The current implementation computes this via \(S^*S\), and it may include some kind of iterative refinement to improve accuracy in some cases.

The parameters m and l refer to the moment and block size of contour integral methods. All columns up to m*l are modified, and the active columns are set to 0..m*l.

See BVSVDMethod for available methods.

The A workspace should be m*l*m*l in size.

Once the decomposition is computed, the numerical rank is estimated by counting the number of singular values that are larger than the tolerance delta, relative to the first singular value.

See Also#

BV: Basis Vectors, BVSetActiveColumns()

Level#

developer

Location#

src/sys/classes/bv/interface/bvcontour.c


Index of all BV routines Table of Contents for all manual pages Index of all manual pages