Annotation of rpl/lapack/lapack/dlasd0.f, revision 1.1
1.1 ! bertrand 1: SUBROUTINE DLASD0( N, SQRE, D, E, U, LDU, VT, LDVT, SMLSIZ, IWORK,
! 2: $ WORK, INFO )
! 3: *
! 4: * -- LAPACK auxiliary routine (version 3.2) --
! 5: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 6: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 7: * November 2006
! 8: *
! 9: * .. Scalar Arguments ..
! 10: INTEGER INFO, LDU, LDVT, N, SMLSIZ, SQRE
! 11: * ..
! 12: * .. Array Arguments ..
! 13: INTEGER IWORK( * )
! 14: DOUBLE PRECISION D( * ), E( * ), U( LDU, * ), VT( LDVT, * ),
! 15: $ WORK( * )
! 16: * ..
! 17: *
! 18: * Purpose
! 19: * =======
! 20: *
! 21: * Using a divide and conquer approach, DLASD0 computes the singular
! 22: * value decomposition (SVD) of a real upper bidiagonal N-by-M
! 23: * matrix B with diagonal D and offdiagonal E, where M = N + SQRE.
! 24: * The algorithm computes orthogonal matrices U and VT such that
! 25: * B = U * S * VT. The singular values S are overwritten on D.
! 26: *
! 27: * A related subroutine, DLASDA, computes only the singular values,
! 28: * and optionally, the singular vectors in compact form.
! 29: *
! 30: * Arguments
! 31: * =========
! 32: *
! 33: * N (input) INTEGER
! 34: * On entry, the row dimension of the upper bidiagonal matrix.
! 35: * This is also the dimension of the main diagonal array D.
! 36: *
! 37: * SQRE (input) INTEGER
! 38: * Specifies the column dimension of the bidiagonal matrix.
! 39: * = 0: The bidiagonal matrix has column dimension M = N;
! 40: * = 1: The bidiagonal matrix has column dimension M = N+1;
! 41: *
! 42: * D (input/output) DOUBLE PRECISION array, dimension (N)
! 43: * On entry D contains the main diagonal of the bidiagonal
! 44: * matrix.
! 45: * On exit D, if INFO = 0, contains its singular values.
! 46: *
! 47: * E (input) DOUBLE PRECISION array, dimension (M-1)
! 48: * Contains the subdiagonal entries of the bidiagonal matrix.
! 49: * On exit, E has been destroyed.
! 50: *
! 51: * U (output) DOUBLE PRECISION array, dimension at least (LDQ, N)
! 52: * On exit, U contains the left singular vectors.
! 53: *
! 54: * LDU (input) INTEGER
! 55: * On entry, leading dimension of U.
! 56: *
! 57: * VT (output) DOUBLE PRECISION array, dimension at least (LDVT, M)
! 58: * On exit, VT' contains the right singular vectors.
! 59: *
! 60: * LDVT (input) INTEGER
! 61: * On entry, leading dimension of VT.
! 62: *
! 63: * SMLSIZ (input) INTEGER
! 64: * On entry, maximum size of the subproblems at the
! 65: * bottom of the computation tree.
! 66: *
! 67: * IWORK (workspace) INTEGER work array.
! 68: * Dimension must be at least (8 * N)
! 69: *
! 70: * WORK (workspace) DOUBLE PRECISION work array.
! 71: * Dimension must be at least (3 * M**2 + 2 * M)
! 72: *
! 73: * INFO (output) INTEGER
! 74: * = 0: successful exit.
! 75: * < 0: if INFO = -i, the i-th argument had an illegal value.
! 76: * > 0: if INFO = 1, an singular value did not converge
! 77: *
! 78: * Further Details
! 79: * ===============
! 80: *
! 81: * Based on contributions by
! 82: * Ming Gu and Huan Ren, Computer Science Division, University of
! 83: * California at Berkeley, USA
! 84: *
! 85: * =====================================================================
! 86: *
! 87: * .. Local Scalars ..
! 88: INTEGER I, I1, IC, IDXQ, IDXQC, IM1, INODE, ITEMP, IWK,
! 89: $ J, LF, LL, LVL, M, NCC, ND, NDB1, NDIML, NDIMR,
! 90: $ NL, NLF, NLP1, NLVL, NR, NRF, NRP1, SQREI
! 91: DOUBLE PRECISION ALPHA, BETA
! 92: * ..
! 93: * .. External Subroutines ..
! 94: EXTERNAL DLASD1, DLASDQ, DLASDT, XERBLA
! 95: * ..
! 96: * .. Executable Statements ..
! 97: *
! 98: * Test the input parameters.
! 99: *
! 100: INFO = 0
! 101: *
! 102: IF( N.LT.0 ) THEN
! 103: INFO = -1
! 104: ELSE IF( ( SQRE.LT.0 ) .OR. ( SQRE.GT.1 ) ) THEN
! 105: INFO = -2
! 106: END IF
! 107: *
! 108: M = N + SQRE
! 109: *
! 110: IF( LDU.LT.N ) THEN
! 111: INFO = -6
! 112: ELSE IF( LDVT.LT.M ) THEN
! 113: INFO = -8
! 114: ELSE IF( SMLSIZ.LT.3 ) THEN
! 115: INFO = -9
! 116: END IF
! 117: IF( INFO.NE.0 ) THEN
! 118: CALL XERBLA( 'DLASD0', -INFO )
! 119: RETURN
! 120: END IF
! 121: *
! 122: * If the input matrix is too small, call DLASDQ to find the SVD.
! 123: *
! 124: IF( N.LE.SMLSIZ ) THEN
! 125: CALL DLASDQ( 'U', SQRE, N, M, N, 0, D, E, VT, LDVT, U, LDU, U,
! 126: $ LDU, WORK, INFO )
! 127: RETURN
! 128: END IF
! 129: *
! 130: * Set up the computation tree.
! 131: *
! 132: INODE = 1
! 133: NDIML = INODE + N
! 134: NDIMR = NDIML + N
! 135: IDXQ = NDIMR + N
! 136: IWK = IDXQ + N
! 137: CALL DLASDT( N, NLVL, ND, IWORK( INODE ), IWORK( NDIML ),
! 138: $ IWORK( NDIMR ), SMLSIZ )
! 139: *
! 140: * For the nodes on bottom level of the tree, solve
! 141: * their subproblems by DLASDQ.
! 142: *
! 143: NDB1 = ( ND+1 ) / 2
! 144: NCC = 0
! 145: DO 30 I = NDB1, ND
! 146: *
! 147: * IC : center row of each node
! 148: * NL : number of rows of left subproblem
! 149: * NR : number of rows of right subproblem
! 150: * NLF: starting row of the left subproblem
! 151: * NRF: starting row of the right subproblem
! 152: *
! 153: I1 = I - 1
! 154: IC = IWORK( INODE+I1 )
! 155: NL = IWORK( NDIML+I1 )
! 156: NLP1 = NL + 1
! 157: NR = IWORK( NDIMR+I1 )
! 158: NRP1 = NR + 1
! 159: NLF = IC - NL
! 160: NRF = IC + 1
! 161: SQREI = 1
! 162: CALL DLASDQ( 'U', SQREI, NL, NLP1, NL, NCC, D( NLF ), E( NLF ),
! 163: $ VT( NLF, NLF ), LDVT, U( NLF, NLF ), LDU,
! 164: $ U( NLF, NLF ), LDU, WORK, INFO )
! 165: IF( INFO.NE.0 ) THEN
! 166: RETURN
! 167: END IF
! 168: ITEMP = IDXQ + NLF - 2
! 169: DO 10 J = 1, NL
! 170: IWORK( ITEMP+J ) = J
! 171: 10 CONTINUE
! 172: IF( I.EQ.ND ) THEN
! 173: SQREI = SQRE
! 174: ELSE
! 175: SQREI = 1
! 176: END IF
! 177: NRP1 = NR + SQREI
! 178: CALL DLASDQ( 'U', SQREI, NR, NRP1, NR, NCC, D( NRF ), E( NRF ),
! 179: $ VT( NRF, NRF ), LDVT, U( NRF, NRF ), LDU,
! 180: $ U( NRF, NRF ), LDU, WORK, INFO )
! 181: IF( INFO.NE.0 ) THEN
! 182: RETURN
! 183: END IF
! 184: ITEMP = IDXQ + IC
! 185: DO 20 J = 1, NR
! 186: IWORK( ITEMP+J-1 ) = J
! 187: 20 CONTINUE
! 188: 30 CONTINUE
! 189: *
! 190: * Now conquer each subproblem bottom-up.
! 191: *
! 192: DO 50 LVL = NLVL, 1, -1
! 193: *
! 194: * Find the first node LF and last node LL on the
! 195: * current level LVL.
! 196: *
! 197: IF( LVL.EQ.1 ) THEN
! 198: LF = 1
! 199: LL = 1
! 200: ELSE
! 201: LF = 2**( LVL-1 )
! 202: LL = 2*LF - 1
! 203: END IF
! 204: DO 40 I = LF, LL
! 205: IM1 = I - 1
! 206: IC = IWORK( INODE+IM1 )
! 207: NL = IWORK( NDIML+IM1 )
! 208: NR = IWORK( NDIMR+IM1 )
! 209: NLF = IC - NL
! 210: IF( ( SQRE.EQ.0 ) .AND. ( I.EQ.LL ) ) THEN
! 211: SQREI = SQRE
! 212: ELSE
! 213: SQREI = 1
! 214: END IF
! 215: IDXQC = IDXQ + NLF - 1
! 216: ALPHA = D( IC )
! 217: BETA = E( IC )
! 218: CALL DLASD1( NL, NR, SQREI, D( NLF ), ALPHA, BETA,
! 219: $ U( NLF, NLF ), LDU, VT( NLF, NLF ), LDVT,
! 220: $ IWORK( IDXQC ), IWORK( IWK ), WORK, INFO )
! 221: IF( INFO.NE.0 ) THEN
! 222: RETURN
! 223: END IF
! 224: 40 CONTINUE
! 225: 50 CONTINUE
! 226: *
! 227: RETURN
! 228: *
! 229: * End of DLASD0
! 230: *
! 231: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dlasd0.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack