General to bidiagonal

Reduces a general fully-populated \(m \times n\) matrix to bidiagonal form through two-sided Householder transformations; when the \(m \ge n\), the result is upper bidiagonal, otherwise it is lower bidiagonal. This routine is most commonly used as a preprocessing step in computing the SVD of a general matrix.

Reduction

Python API

Bidiag(A[, bidiagOnly=False])

C++ API

void Bidiag(Matrix<F> &A, Matrix<F> &tP, Matrix<F> &tQ)
void Bidiag(AbstractDistMatrix<F> &A, AbstractDistMatrix<F> &tP, AbstractDistMatrix<F> &tQ)

Overwrites the main and sub (or super) diagonal of the real matrix A with the resulting bidiagonal matrix and stores the scaled Householder vectors in the remainder of the matrix. The complex case must also store the scalings of the Householder transformations (in tP and tQ) if they are to be applied.

void bidiag::Explicit(Matrix<F> &A, Matrix<F> &P, Matrix<F> &Q)
void bidiag::Explicit(AbstractDistMatrix<F> &A, AbstractDistMatrix<F> &P, AbstractDistMatrix<F> &Q)

Overwrite \(A\) with the bidiagonal matrix, \(B = Q^H A P\), and also return \(P\) and \(Q\).

void bidiag::ExplicitCondensed(Matrix<F> &A)
void bidiag::ExplicitCondensed(AbstractDistMatrix<F> &A)

Returns just the resulting bidiagonal matrix, \(B = Q^H A P\).

C API

Single-precision

ElError ElBidiag_s(ElMatrix_s A, ElMatrix_s tP, ElMatrix_s tQ)
ElError ElBidiagDist_s(ElDistMatrix_s A, ElDistMatrix_s tP, ElDistMatrix_s tQ)

Overwrites the main and sub (or super) diagonal of the real matrix A with the resulting bidiagonal matrix and stores the scaled Householder vectors in the remainder of the matrix. The complex case must also store the scalings of the Householder transformations (in tP and tQ) if they are to be applied.

ElError ElBidiagExplicit_s(ElMatrix_s A, ElMatrix_s P, ElMatrix_s Q)
ElError ElBidiagExplicitDist_s(ElDistMatrix_s A, ElDistMatrix_s P, ElDistMatrix_s Q)

Overwrite \(A\) with the bidiagonal matrix, \(B = Q^H A P\), and also return \(P\) and \(Q\).

ElError ElBidiagExplicitCondensed_s(ElMatrix_s A)
ElError ElBidiagExplicitCondensedDist_s(ElDistMatrix_s A)

Returns just the resulting bidiagonal matrix, \(B = Q^H A P\).

Double-precision

ElError ElBidiag_d(ElMatrix_d A, ElMatrix_d tP, ElMatrix_d tQ)
ElError ElBidiagDist_d(ElDistMatrix_d A, ElDistMatrix_d tP, ElDistMatrix_d tQ)

Overwrites the main and sub (or super) diagonal of the real matrix A with the resulting bidiagonal matrix and stores the scaled Householder vectors in the remainder of the matrix. The complex case must also store the scalings of the Householder transformations (in tP and tQ) if they are to be applied.

ElError ElBidiagExplicit_d(ElMatrix_d A, ElMatrix_d P, ElMatrix_d Q)
ElError ElBidiagExplicitDist_d(ElDistMatrix_d A, ElDistMatrix_d P, ElDistMatrix_d Q)

Overwrite \(A\) with the bidiagonal matrix, \(B = Q^H A P\), and also return \(P\) and \(Q\).

ElError ElBidiagExplicitCondensed_d(ElMatrix_d A)
ElError ElBidiagExplicitCondensedDist_d(ElDistMatrix_d A)

Returns just the resulting bidiagonal matrix, \(B = Q^H A P\).

Single-precision complex

ElError ElBidiag_c(ElMatrix_c A, ElMatrix_c tP, ElMatrix_c tQ)
ElError ElBidiagDist_c(ElDistMatrix_c A, ElDistMatrix_c tP, ElDistMatrix_c tQ)

Overwrites the main and sub (or super) diagonal of the real matrix A with the resulting bidiagonal matrix and stores the scaled Householder vectors in the remainder of the matrix. The complex case must also store the scalings of the Householder transformations (in tP and tQ) if they are to be applied.

ElError ElBidiagExplicit_c(ElMatrix_c A, ElMatrix_c P, ElMatrix_c Q)
ElError ElBidiagExplicitDist_c(ElDistMatrix_c A, ElDistMatrix_c P, ElDistMatrix_c Q)

Overwrite \(A\) with the bidiagonal matrix, \(B = Q^H A P\), and also return \(P\) and \(Q\).

ElError ElBidiagExplicitCondensed_c(ElMatrix_c A)
ElError ElBidiagExplicitCondensedDist_c(ElDistMatrix_c A)

Returns just the resulting bidiagonal matrix, \(B = Q^H A P\).

Double-precision complex

ElError ElBidiag_z(ElMatrix_z A, ElMatrix_z tP, ElMatrix_z tQ)
ElError ElBidiagDist_z(ElDistMatrix_z A, ElDistMatrix_z tP, ElDistMatrix_z tQ)

Overwrites the main and sub (or super) diagonal of the real matrix A with the resulting bidiagonal matrix and stores the scaled Householder vectors in the remainder of the matrix. The complex case must also store the scalings of the Householder transformations (in tP and tQ) if they are to be applied.

ElError ElBidiagExplicit_z(ElMatrix_z A, ElMatrix_z P, ElMatrix_z Q)
ElError ElBidiagExplicitDist_z(ElDistMatrix_z A, ElDistMatrix_z P, ElDistMatrix_z Q)

Overwrite \(A\) with the bidiagonal matrix, \(B = Q^H A P\), and also return \(P\) and \(Q\).

ElError ElBidiagExplicitCondensed_z(ElMatrix_z A)
ElError ElBidiagExplicitCondensedDist_z(ElDistMatrix_z A)

Returns just the resulting bidiagonal matrix, \(B = Q^H A P\).

Applying the changes of basis

Python API

ApplyQAfterBidiag(side, orient, A, t, B)

C++ API

void bidiag::ApplyQ(LeftOrRight side, Orientation orientation, const Matrix<F> &A, const Matrix<F> &t, Matrix<F> &B)
void bidiag::ApplyQ(LeftOrRight side, Orientation orientation, const AbstractDistMatrix<F> &A, const AbstractDistMatrix<F> &t, AbstractDistMatrix<F> &B)
void bidiag::ApplyP(LeftOrRight side, Orientation orientation, const Matrix<F> &A, const Matrix<F> &t, Matrix<F> &B)
void bidiag::ApplyP(LeftOrRight side, Orientation orientation, const AbstractDistMatrix<F> &A, const AbstractDistMatrix<F> &t, AbstractDistMatrix<F> &B)

C API

Single-precision

ElError ElApplyQAfterBidiag_s(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_s A, ElConstMatrix_s t, ElMatrix_s B)
ElError ElApplyQAfterBidiagDist_s(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_s A, ElConstDistMatrix_s t, ElDistMatrix_s B)
ElError ElApplyPAfterBidiag_s(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_s A, ElConstMatrix_s t, ElMatrix_s B)
ElError ElApplyPAfterBidiagDist_s(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_s A, ElConstDistMatrix_s t, ElDistMatrix_s B)

Double-precision

ElError ElApplyQAfterBidiag_d(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_d A, ElConstMatrix_d t, ElMatrix_d B)
ElError ElApplyQAfterBidiagDist_d(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_d A, ElConstDistMatrix_d t, ElDistMatrix_d B)
ElError ElApplyPAfterBidiag_d(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_d A, ElConstMatrix_d t, ElMatrix_d B)
ElError ElApplyPAfterBidiagDist_d(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_d A, ElConstDistMatrix_d t, ElDistMatrix_d B)

Single-precision complex

ElError ElApplyQAfterBidiag_c(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_c A, ElConstMatrix_c t, ElMatrix_c B)
ElError ElApplyQAfterBidiagDist_c(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_c A, ElConstDistMatrix_c t, ElDistMatrix_c B)
ElError ElApplyPAfterBidiag_c(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_c A, ElConstMatrix_c t, ElMatrix_c B)
ElError ElApplyPAfterBidiagDist_c(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_c A, ElConstDistMatrix_c t, ElDistMatrix_c B)

Double-precision complex

ElError ElApplyQAfterBidiag_z(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_z A, ElConstMatrix_z t, ElMatrix_z B)
ElError ElApplyQAfterBidiagDist_z(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_z A, ElConstDistMatrix_z t, ElDistMatrix_z B)
ElError ElApplyPAfterBidiag_z(ElLeftOrRight side, ElOrientation orientation, ElConstMatrix_z A, ElConstMatrix_z t, ElMatrix_z B)
ElError ElApplyPAfterBidiagDist_z(ElLeftOrRight side, ElOrientation orientation, ElConstDistMatrix_z A, ElConstDistMatrix_z t, ElDistMatrix_z B)