CGESVJ(3)
๐ NAME ๐
okay so this is CGESVJ(3) โ aka cgesvj.f(3) from the Library Functions Manual โ and the file it points at is /home/abuild/rpmbuild/BUILD/lapack-static-3.12.0-build/lapack-3.12.0/SRC/cgesvj.f. Weโre doing this in two parts, bestie, and this is part 1. Itโs the CGESVJ routine page, no cap.
โจ SYNOPSIS โจ
๐ง Functions/Subroutines ๐ง
subroutine CGESVJ (joba, jobu, jobv, m, n, a, lda, sva, mv, v, ldv,
cwork, lwork, rwork, lrwork, info)
CGESVJ
๐ฆ Function/Subroutine Documentation ๐ฆ
subroutine CGESVJ (character*1 joba, character*1 jobu, character*1 jobv,
integer m, integer n, complex, dimension( lda, * ) a, integer lda, real,
dimension( n ) sva, integer mv, complex, dimension( ldv, * ) v, integer
ldv, complex, dimension( lwork ) cwork, integer lwork, real, dimension(
lrwork ) rwork, integer lrwork, integer info)
CGESVJ โ yes bestie, the whole page is about this one. ๐
๐ธ Purpose ๐ธ
okay so CGESVJ computes the singular value decomposition (SVD) of a complex M-by-N matrix A, where M >= N. The SVD of A is written as
[++] [xx] [x0] [xx]
A = U * SIGMA * V^*, [++] = [xx] * [ox] * [xx]
[++] [xx]
where SIGMA is an N-by-N diagonal matrix, U is an M-by-N orthonormal matrix, and V is an N-by-N unitary matrix. The diagonal elements of SIGMA are the singular values of A. The columns of U and V are the left and the right singular vectors of A, respectively. Thatโs the whole main-character plot, periodt. ๐
๐ Parameters ๐
๐ JOBA
JOBA is CHARACTER*1. Specifies the structure of A.
'L'โ The input matrixAis lower triangular;'U'โ The input matrixAis upper triangular;'G'โ The input matrixAis general M-by-N matrix,M >= N.
๐ฅ JOBU
JOBU is CHARACTER*1. Specifies whether to compute the left singular vectors (columns of U):
'U'or'F'โ The left singular vectors corresponding to the nonzero singular values are computed and returned in the leading columns ofA. See more details in the description ofA. The default numerical orthogonality threshold is set to approximatelyTOL=CTOL*EPS,CTOL=SQRT(M),EPS=SLAMCH('E').'C'โ Analogous toJOBU='U', except that user can control the level of numerical orthogonality of the computed left singular vectors.TOLcan be set toTOL = CTOL*EPS, whereCTOLis given on input in the arrayWORK. NoCTOLsmaller thanONEis allowed.CTOLgreater than1 / EPSis meaningless. The option'C'can be used ifM*EPSis satisfactory orthogonality of the computed left singular vectors, soCTOL=Mcould save few sweeps of Jacobi rotations. See the descriptions ofAandWORK(1).'N'โ The matrixUis not computed. However, see the description ofA.
๐ฆ JOBV
JOBV is CHARACTER*1. Specifies whether to compute the right singular vectors, that is, the matrix V:
'V'or'J'โ the matrixVis computed and returned in the arrayV'A'โ the Jacobi rotations are applied to the MV-by-N arrayV. In other words, the right singular vector matrixVis not computed explicitly; instead it is applied to an MV-by-N matrix initially stored in the first MV rows ofV.'N'โ the matrixVis not computed and the arrayVis not referenced
โญ M
M is INTEGER. The number of rows of the input matrix A. 1/SLAMCH('E') > M >= 0.
โญ N
N is INTEGER. The number of columns of the input matrix A. M >= N >= 0.
๐ A
A is COMPLEX array, dimension (LDA,N). On entry, the M-by-N matrix A.
On exit,
If JOBU = 'U' .OR. JOBU = 'C':
- If
INFO = 0:RANKAorthonormal columns ofUare returned in the leadingRANKAcolumns of the arrayA. HereRANKA <= Nis the number of computed singular values ofAthat are above the underflow thresholdSLAMCH('S'). The singular vectors corresponding to underflowed or zero singular values are not computed. The value ofRANKAis returned in the arrayRWORKasRANKA=NINT(RWORK(2)). Also see the descriptions ofSVAandRWORK. The computed columns ofUare mutually numerically orthogonal up to approximatelyTOL=SQRT(M)*EPS(default); orTOL=CTOL*EPS(JOBU = 'C'), see the description ofJOBU. - If
INFO > 0, the procedureCGESVJdid not converge in the given number of iterations (sweeps). In that case, the computed columns ofUmay not be orthogonal up toTOL. The outputU(stored inA),SIGMA(given by the computed singular values inSVA(1:N)) andVis still a decomposition of the input matrixAin the sense that the residual|| A - SCALE * U * SIGMA * V^* ||_2 / ||A||_2is small.
If JOBU = 'N':
- If
INFO = 0: Note that the left singular vectors are 'for free' in the one-sided Jacobi SVD algorithm. However, if only the singular values are needed, the level of numerical orthogonality ofUis not an issue and iterations are stopped when the columns of the iterated matrix are numerically orthogonal up to approximatelyM*EPS. Thus, on exit,Acontains the columns ofUscaled with the corresponding singular values. - If
INFO > 0: the procedureCGESVJdid not converge in the given number of iterations (sweeps).
โญ LDA
LDA is INTEGER. The leading dimension of the array A. LDA >= max(1,M).
๐ SVA
SVA is REAL array, dimension (N).
On exit,
If INFO = 0:
depending on the value SCALE = RWORK(1), we have:
- If
SCALE = ONE:SVA(1:N)contains the computed singular values ofA. During the computationSVAcontains the Euclidean column norms of the iterated matrices in the arrayA. - If
SCALE .NE. ONE: The singular values ofAareSCALE*SVA(1:N), and this factored representation is due to the fact that some of the singular values ofAmight underflow or overflow.
If INFO > 0:
the procedure CGESVJ did not converge in the given number of iterations (sweeps) and SCALE*SVA(1:N) may not be accurate.
๐ MV
MV is INTEGER. If JOBV = 'A', then the product of Jacobi rotations in CGESVJ is applied to the first MV rows of V. See the description of JOBV.
๐ธ V
V is COMPLEX array, dimension (LDV,N).
- If
JOBV = 'V', thenVcontains on exit the N-by-N matrix of the right singular vectors; - If
JOBV = 'A', thenVcontains the product of the computed right singular vector matrix and the initial matrix in the arrayV. - If
JOBV = 'N', thenVis not referenced.
๐ LDV
LDV is INTEGER. The leading dimension of the array V, LDV >= 1. If JOBV = 'V', then LDV >= max(1,N). If JOBV = 'A', then LDV >= max(1,MV) .
๐ซถ CWORK
CWORK is COMPLEX array, dimension (max(1,LWORK)). Used as workspace. If on entry LWORK = -1, then a workspace query is assumed and no computation is done; CWORK(1) is set to the minial (and optimal) length of CWORK.
๐ญ LWORK
LWORK is INTEGER. Length of CWORK, LWORK >= M+N.
๐ RWORK
Part 2 picks up right here with RWORK, bestie.
okay so back to RWORK โ it's a REAL array, dimension (max(6,LRWORK)), and bestie it is literally doing double duty here, so pay attention to whether you're looking at it on entry or on exit ๐
On entry: if JOBU = 'C':
RWORK(1) = CTOL, whereCTOLdefines the threshold for convergence. The process stops if all columns ofAare mutually orthogonal up toCTOL*EPS,EPS=SLAMCH('E'). It is required thatCTOL >= ONE, i.e. it is not allowed to force the routine to obtain orthogonality belowEPSILON๐ซ
On exit:
RWORK(1) = SCALEis the scaling factor such thatSCALE*SVA(1:N)are the computed singular values ofA. (See description ofSVA().) โจRWORK(2) = NINT(RWORK(2))is the number of the computed nonzero singular values.RWORK(3) = NINT(RWORK(3))is the number of the computed singular values that are larger than the underflow threshold.RWORK(4) = NINT(RWORK(4))is the number of sweeps of Jacobi rotations needed for numerical convergence.RWORK(5) = max_{i.NE.j} |COS(A(:,i),A(:,j))|in the last sweep. This is useful information in cases whenCGESVJdid not converge, as it can be used to estimate whether the output is still useful and for post festum analysis ๐RWORK(6)= the largest absolute value over all sines of the Jacobi rotation angles in the last sweep. It can be useful for a post festum analysis.
If on entryLRWORK = -1, then a workspace query is assumed and no computation is done;RWORK(1)is set to theminial(and optimal) length ofRWORK. (Yes, "minial". That's the source's spelling, not mine, I would never ๐ญ๐)
๐ธ LRWORK ๐ธ
LRWORK is INTEGER, the length of RWORK, and it must satisfy LRWORK >= MAX(6,N).
๐
INFO ๐
INFO is INTEGER:
= 0: successful exit ๐< 0: ifINFO = -i, then thei-th argument had an illegal value ๐> 0:CGESVJdid not converge in the maximal allowed number (NSWEEP=30) of sweeps. The output may still be useful. See the description ofRWORKโ ๏ธ
๐ Author ๐
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
๐ง Further Details: ๐ง
The orthogonal N-by-N matrix V is obtained as a product of Jacobi plane rotations. In the case of underflow of the tangent of the Jacobi angle, a modified Jacobi transformation of Drmac [3] is used. Pivot strategy uses column interchanges of de Rijk [1]. The relative accuracy of the computed singular values and the accuracy of the computed singular vectors (in angle metric) is as guaranteed by the theory of Demmel and Veselic [2]. The condition number that determines the accuracy in the full rank case is essentially min_{D=diag} kappa(A*D), where kappa(.) is the spectral condition number. The best performance of this Jacobi SVD procedure is achieved if used in an accelerated version of Drmac and Veselic [4,5], and it is the kernel routine in the SIGMA library [6]. Some tuning parameters (marked with [TP]) are available for the implementer.
The computational range for the nonzero singular values is the machine number interval ( UNDERFLOW , OVERFLOW ). In extreme cases, even denormalized singular values can be computed with the corresponding gradual loss of accurate digits ๐ฎโ๐จ
๐ฆ Contributor: ๐ฆ
============
Zlatko Drmac (Zagreb, Croatia)
๐ References: ๐
[1]P. P. M. De Rijk: A one-sided Jacobi algorithm for computing the singular value decomposition on a vector computer. SIAM J. Sci. Stat. Comp., Vol. 10 (1998), pp. 359-371.[2]J. Demmel and K. Veselic: Jacobi method is more accurate than QR.[3]Z. Drmac: Implementation of Jacobi rotations for accurate singular value computation in floating point arithmetic. SIAM J. Sci. Comp., Vol. 18 (1997), pp. 1200-1222.[4]Z. Drmac and K. Veselic: New fast and accurate Jacobi SVD algorithm I. SIAM J. Matrix Anal. Appl. Vol. 35, No. 2 (2008), pp. 1322-1342. LAPACK Working note 169.[5]Z. Drmac and K. Veselic: New fast and accurate Jacobi SVD algorithm II. SIAM J. Matrix Anal. Appl. Vol. 35, No. 2 (2008), pp. 1343-1362. LAPACK Working note 170.[6]Z. Drmac: SIGMA - mathematical software library for accurate SVD, PSV, QSVD, (H,K)-SVD computations. Department of Mathematics, University of Zagreb, 2008, 2015.
๐ Bugs, examples and comments: ๐
===========================
Please report all bugs and send interesting test examples and comments to drmac@math.hr. Thank you ๐ซถ
Definition at line 349 of file cgesvj.f.
๐ค Author (the robot edition) ๐ค
Generated automatically by Doxygen for LAPACK from the source code.
okay bestie that's the whole page, go forth and diagonalize responsibly ๐โจ This was LAPACK, Version 3.12.0, straight out of /home/abui...RC/cgesvj.f(3) โ mwah, slay, periodt ๐๐
๐ฆ where this page comes from
- opensuse-leap 16.0: lapack-man 3.12.0-160000.2.2 ยท BSD-3-Clause
- opensuse-leap 16.0: lapack-man 3.12.0-160000.2.2 ยท BSD-3-Clause
- opensuse-leap 16.0: lapack-man 3.12.0-160000.2.2 ยท BSD-3-Clause
- opensuse-leap 16.0: lapack-man 3.12.0-160000.2.2 ยท BSD-3-Clause
- ๐ girlified by deepseek-ai/DeepSeek-V4.1-Flash-1017 on : a language model's rewrite of the original, under its terms