--- rpl/lapack/lapack/dlasd6.f 2011/07/22 07:38:08 1.10 +++ rpl/lapack/lapack/dlasd6.f 2011/11/21 20:42:58 1.11 @@ -1,12 +1,323 @@ +*> \brief \b DLASD6 +* +* =========== DOCUMENTATION =========== +* +* Online html documentation available at +* http://www.netlib.org/lapack/explore-html/ +* +*> \htmlonly +*> Download DLASD6 + dependencies +*> +*> [TGZ] +*> +*> [ZIP] +*> +*> [TXT] +*> \endhtmlonly +* +* Definition: +* =========== +* +* SUBROUTINE DLASD6( ICOMPQ, NL, NR, SQRE, D, VF, VL, ALPHA, BETA, +* IDXQ, PERM, GIVPTR, GIVCOL, LDGCOL, GIVNUM, +* LDGNUM, POLES, DIFL, DIFR, Z, K, C, S, WORK, +* IWORK, INFO ) +* +* .. Scalar Arguments .. +* INTEGER GIVPTR, ICOMPQ, INFO, K, LDGCOL, LDGNUM, NL, +* $ NR, SQRE +* DOUBLE PRECISION ALPHA, BETA, C, S +* .. +* .. Array Arguments .. +* INTEGER GIVCOL( LDGCOL, * ), IDXQ( * ), IWORK( * ), +* $ PERM( * ) +* DOUBLE PRECISION D( * ), DIFL( * ), DIFR( * ), +* $ GIVNUM( LDGNUM, * ), POLES( LDGNUM, * ), +* $ VF( * ), VL( * ), WORK( * ), Z( * ) +* .. +* +* +*> \par Purpose: +* ============= +*> +*> \verbatim +*> +*> DLASD6 computes the SVD of an updated upper bidiagonal matrix B +*> obtained by merging two smaller ones by appending a row. This +*> routine is used only for the problem which requires all singular +*> values and optionally singular vector matrices in factored form. +*> B is an N-by-M matrix with N = NL + NR + 1 and M = N + SQRE. +*> A related subroutine, DLASD1, handles the case in which all singular +*> values and singular vectors of the bidiagonal matrix are desired. +*> +*> DLASD6 computes the SVD as follows: +*> +*> ( D1(in) 0 0 0 ) +*> B = U(in) * ( Z1**T a Z2**T b ) * VT(in) +*> ( 0 0 D2(in) 0 ) +*> +*> = U(out) * ( D(out) 0) * VT(out) +*> +*> where Z**T = (Z1**T a Z2**T b) = u**T VT**T, and u is a vector of dimension M +*> with ALPHA and BETA in the NL+1 and NL+2 th entries and zeros +*> elsewhere; and the entry b is empty if SQRE = 0. +*> +*> The singular values of B can be computed using D1, D2, the first +*> components of all the right singular vectors of the lower block, and +*> the last components of all the right singular vectors of the upper +*> block. These components are stored and updated in VF and VL, +*> respectively, in DLASD6. Hence U and VT are not explicitly +*> referenced. +*> +*> The singular values are stored in D. The algorithm consists of two +*> stages: +*> +*> The first stage consists of deflating the size of the problem +*> when there are multiple singular values or if there is a zero +*> in the Z vector. For each such occurence the dimension of the +*> secular equation problem is reduced by one. This stage is +*> performed by the routine DLASD7. +*> +*> The second stage consists of calculating the updated +*> singular values. This is done by finding the roots of the +*> secular equation via the routine DLASD4 (as called by DLASD8). +*> This routine also updates VF and VL and computes the distances +*> between the updated singular values and the old singular +*> values. +*> +*> DLASD6 is called from DLASDA. +*> \endverbatim +* +* Arguments: +* ========== +* +*> \param[in] ICOMPQ +*> \verbatim +*> ICOMPQ is INTEGER +*> Specifies whether singular vectors are to be computed in +*> factored form: +*> = 0: Compute singular values only. +*> = 1: Compute singular vectors in factored form as well. +*> \endverbatim +*> +*> \param[in] NL +*> \verbatim +*> NL is INTEGER +*> The row dimension of the upper block. NL >= 1. +*> \endverbatim +*> +*> \param[in] NR +*> \verbatim +*> NR is INTEGER +*> The row dimension of the lower block. NR >= 1. +*> \endverbatim +*> +*> \param[in] SQRE +*> \verbatim +*> SQRE is INTEGER +*> = 0: the lower block is an NR-by-NR square matrix. +*> = 1: the lower block is an NR-by-(NR+1) rectangular matrix. +*> +*> The bidiagonal matrix has row dimension N = NL + NR + 1, +*> and column dimension M = N + SQRE. +*> \endverbatim +*> +*> \param[in,out] D +*> \verbatim +*> D is DOUBLE PRECISION array, dimension ( NL+NR+1 ). +*> On entry D(1:NL,1:NL) contains the singular values of the +*> upper block, and D(NL+2:N) contains the singular values +*> of the lower block. On exit D(1:N) contains the singular +*> values of the modified matrix. +*> \endverbatim +*> +*> \param[in,out] VF +*> \verbatim +*> VF is DOUBLE PRECISION array, dimension ( M ) +*> On entry, VF(1:NL+1) contains the first components of all +*> right singular vectors of the upper block; and VF(NL+2:M) +*> contains the first components of all right singular vectors +*> of the lower block. On exit, VF contains the first components +*> of all right singular vectors of the bidiagonal matrix. +*> \endverbatim +*> +*> \param[in,out] VL +*> \verbatim +*> VL is DOUBLE PRECISION array, dimension ( M ) +*> On entry, VL(1:NL+1) contains the last components of all +*> right singular vectors of the upper block; and VL(NL+2:M) +*> contains the last components of all right singular vectors of +*> the lower block. On exit, VL contains the last components of +*> all right singular vectors of the bidiagonal matrix. +*> \endverbatim +*> +*> \param[in,out] ALPHA +*> \verbatim +*> ALPHA is DOUBLE PRECISION +*> Contains the diagonal element associated with the added row. +*> \endverbatim +*> +*> \param[in,out] BETA +*> \verbatim +*> BETA is DOUBLE PRECISION +*> Contains the off-diagonal element associated with the added +*> row. +*> \endverbatim +*> +*> \param[out] IDXQ +*> \verbatim +*> IDXQ is INTEGER array, dimension ( N ) +*> This contains the permutation which will reintegrate the +*> subproblem just solved back into sorted order, i.e. +*> D( IDXQ( I = 1, N ) ) will be in ascending order. +*> \endverbatim +*> +*> \param[out] PERM +*> \verbatim +*> PERM is INTEGER array, dimension ( N ) +*> The permutations (from deflation and sorting) to be applied +*> to each block. Not referenced if ICOMPQ = 0. +*> \endverbatim +*> +*> \param[out] GIVPTR +*> \verbatim +*> GIVPTR is INTEGER +*> The number of Givens rotations which took place in this +*> subproblem. Not referenced if ICOMPQ = 0. +*> \endverbatim +*> +*> \param[out] GIVCOL +*> \verbatim +*> GIVCOL is INTEGER array, dimension ( LDGCOL, 2 ) +*> Each pair of numbers indicates a pair of columns to take place +*> in a Givens rotation. Not referenced if ICOMPQ = 0. +*> \endverbatim +*> +*> \param[in] LDGCOL +*> \verbatim +*> LDGCOL is INTEGER +*> leading dimension of GIVCOL, must be at least N. +*> \endverbatim +*> +*> \param[out] GIVNUM +*> \verbatim +*> GIVNUM is DOUBLE PRECISION array, dimension ( LDGNUM, 2 ) +*> Each number indicates the C or S value to be used in the +*> corresponding Givens rotation. Not referenced if ICOMPQ = 0. +*> \endverbatim +*> +*> \param[in] LDGNUM +*> \verbatim +*> LDGNUM is INTEGER +*> The leading dimension of GIVNUM and POLES, must be at least N. +*> \endverbatim +*> +*> \param[out] POLES +*> \verbatim +*> POLES is DOUBLE PRECISION array, dimension ( LDGNUM, 2 ) +*> On exit, POLES(1,*) is an array containing the new singular +*> values obtained from solving the secular equation, and +*> POLES(2,*) is an array containing the poles in the secular +*> equation. Not referenced if ICOMPQ = 0. +*> \endverbatim +*> +*> \param[out] DIFL +*> \verbatim +*> DIFL is DOUBLE PRECISION array, dimension ( N ) +*> On exit, DIFL(I) is the distance between I-th updated +*> (undeflated) singular value and the I-th (undeflated) old +*> singular value. +*> \endverbatim +*> +*> \param[out] DIFR +*> \verbatim +*> DIFR is DOUBLE PRECISION array, +*> dimension ( LDGNUM, 2 ) if ICOMPQ = 1 and +*> dimension ( N ) if ICOMPQ = 0. +*> On exit, DIFR(I, 1) is the distance between I-th updated +*> (undeflated) singular value and the I+1-th (undeflated) old +*> singular value. +*> +*> If ICOMPQ = 1, DIFR(1:K,2) is an array containing the +*> normalizing factors for the right singular vector matrix. +*> +*> See DLASD8 for details on DIFL and DIFR. +*> \endverbatim +*> +*> \param[out] Z +*> \verbatim +*> Z is DOUBLE PRECISION array, dimension ( M ) +*> The first elements of this array contain the components +*> of the deflation-adjusted updating row vector. +*> \endverbatim +*> +*> \param[out] K +*> \verbatim +*> K is INTEGER +*> Contains the dimension of the non-deflated matrix, +*> This is the order of the related secular equation. 1 <= K <=N. +*> \endverbatim +*> +*> \param[out] C +*> \verbatim +*> C is DOUBLE PRECISION +*> C contains garbage if SQRE =0 and the C-value of a Givens +*> rotation related to the right null space if SQRE = 1. +*> \endverbatim +*> +*> \param[out] S +*> \verbatim +*> S is DOUBLE PRECISION +*> S contains garbage if SQRE =0 and the S-value of a Givens +*> rotation related to the right null space if SQRE = 1. +*> \endverbatim +*> +*> \param[out] WORK +*> \verbatim +*> WORK is DOUBLE PRECISION array, dimension ( 4 * M ) +*> \endverbatim +*> +*> \param[out] IWORK +*> \verbatim +*> IWORK is INTEGER array, dimension ( 3 * N ) +*> \endverbatim +*> +*> \param[out] INFO +*> \verbatim +*> INFO is INTEGER +*> = 0: successful exit. +*> < 0: if INFO = -i, the i-th argument had an illegal value. +*> > 0: if INFO = 1, a singular value did not converge +*> \endverbatim +* +* Authors: +* ======== +* +*> \author Univ. of Tennessee +*> \author Univ. of California Berkeley +*> \author Univ. of Colorado Denver +*> \author NAG Ltd. +* +*> \date November 2011 +* +*> \ingroup auxOTHERauxiliary +* +*> \par Contributors: +* ================== +*> +*> Ming Gu and Huan Ren, Computer Science Division, University of +*> California at Berkeley, USA +*> +* ===================================================================== SUBROUTINE DLASD6( ICOMPQ, NL, NR, SQRE, D, VF, VL, ALPHA, BETA, $ IDXQ, PERM, GIVPTR, GIVCOL, LDGCOL, GIVNUM, $ LDGNUM, POLES, DIFL, DIFR, Z, K, C, S, WORK, $ IWORK, INFO ) * -* -- LAPACK auxiliary routine (version 3.3.0) -- +* -- LAPACK auxiliary routine (version 3.4.0) -- * -- LAPACK is a software package provided by Univ. of Tennessee, -- * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- -* November 2010 +* November 2011 * * .. Scalar Arguments .. INTEGER GIVPTR, ICOMPQ, INFO, K, LDGCOL, LDGNUM, NL, @@ -21,185 +332,6 @@ $ VF( * ), VL( * ), WORK( * ), Z( * ) * .. * -* Purpose -* ======= -* -* DLASD6 computes the SVD of an updated upper bidiagonal matrix B -* obtained by merging two smaller ones by appending a row. This -* routine is used only for the problem which requires all singular -* values and optionally singular vector matrices in factored form. -* B is an N-by-M matrix with N = NL + NR + 1 and M = N + SQRE. -* A related subroutine, DLASD1, handles the case in which all singular -* values and singular vectors of the bidiagonal matrix are desired. -* -* DLASD6 computes the SVD as follows: -* -* ( D1(in) 0 0 0 ) -* B = U(in) * ( Z1**T a Z2**T b ) * VT(in) -* ( 0 0 D2(in) 0 ) -* -* = U(out) * ( D(out) 0) * VT(out) -* -* where Z**T = (Z1**T a Z2**T b) = u**T VT**T, and u is a vector of dimension M -* with ALPHA and BETA in the NL+1 and NL+2 th entries and zeros -* elsewhere; and the entry b is empty if SQRE = 0. -* -* The singular values of B can be computed using D1, D2, the first -* components of all the right singular vectors of the lower block, and -* the last components of all the right singular vectors of the upper -* block. These components are stored and updated in VF and VL, -* respectively, in DLASD6. Hence U and VT are not explicitly -* referenced. -* -* The singular values are stored in D. The algorithm consists of two -* stages: -* -* The first stage consists of deflating the size of the problem -* when there are multiple singular values or if there is a zero -* in the Z vector. For each such occurence the dimension of the -* secular equation problem is reduced by one. This stage is -* performed by the routine DLASD7. -* -* The second stage consists of calculating the updated -* singular values. This is done by finding the roots of the -* secular equation via the routine DLASD4 (as called by DLASD8). -* This routine also updates VF and VL and computes the distances -* between the updated singular values and the old singular -* values. -* -* DLASD6 is called from DLASDA. -* -* Arguments -* ========= -* -* ICOMPQ (input) INTEGER -* Specifies whether singular vectors are to be computed in -* factored form: -* = 0: Compute singular values only. -* = 1: Compute singular vectors in factored form as well. -* -* NL (input) INTEGER -* The row dimension of the upper block. NL >= 1. -* -* NR (input) INTEGER -* The row dimension of the lower block. NR >= 1. -* -* SQRE (input) INTEGER -* = 0: the lower block is an NR-by-NR square matrix. -* = 1: the lower block is an NR-by-(NR+1) rectangular matrix. -* -* The bidiagonal matrix has row dimension N = NL + NR + 1, -* and column dimension M = N + SQRE. -* -* D (input/output) DOUBLE PRECISION array, dimension ( NL+NR+1 ). -* On entry D(1:NL,1:NL) contains the singular values of the -* upper block, and D(NL+2:N) contains the singular values -* of the lower block. On exit D(1:N) contains the singular -* values of the modified matrix. -* -* VF (input/output) DOUBLE PRECISION array, dimension ( M ) -* On entry, VF(1:NL+1) contains the first components of all -* right singular vectors of the upper block; and VF(NL+2:M) -* contains the first components of all right singular vectors -* of the lower block. On exit, VF contains the first components -* of all right singular vectors of the bidiagonal matrix. -* -* VL (input/output) DOUBLE PRECISION array, dimension ( M ) -* On entry, VL(1:NL+1) contains the last components of all -* right singular vectors of the upper block; and VL(NL+2:M) -* contains the last components of all right singular vectors of -* the lower block. On exit, VL contains the last components of -* all right singular vectors of the bidiagonal matrix. -* -* ALPHA (input/output) DOUBLE PRECISION -* Contains the diagonal element associated with the added row. -* -* BETA (input/output) DOUBLE PRECISION -* Contains the off-diagonal element associated with the added -* row. -* -* IDXQ (output) INTEGER array, dimension ( N ) -* This contains the permutation which will reintegrate the -* subproblem just solved back into sorted order, i.e. -* D( IDXQ( I = 1, N ) ) will be in ascending order. -* -* PERM (output) INTEGER array, dimension ( N ) -* The permutations (from deflation and sorting) to be applied -* to each block. Not referenced if ICOMPQ = 0. -* -* GIVPTR (output) INTEGER -* The number of Givens rotations which took place in this -* subproblem. Not referenced if ICOMPQ = 0. -* -* GIVCOL (output) INTEGER array, dimension ( LDGCOL, 2 ) -* Each pair of numbers indicates a pair of columns to take place -* in a Givens rotation. Not referenced if ICOMPQ = 0. -* -* LDGCOL (input) INTEGER -* leading dimension of GIVCOL, must be at least N. -* -* GIVNUM (output) DOUBLE PRECISION array, dimension ( LDGNUM, 2 ) -* Each number indicates the C or S value to be used in the -* corresponding Givens rotation. Not referenced if ICOMPQ = 0. -* -* LDGNUM (input) INTEGER -* The leading dimension of GIVNUM and POLES, must be at least N. -* -* POLES (output) DOUBLE PRECISION array, dimension ( LDGNUM, 2 ) -* On exit, POLES(1,*) is an array containing the new singular -* values obtained from solving the secular equation, and -* POLES(2,*) is an array containing the poles in the secular -* equation. Not referenced if ICOMPQ = 0. -* -* DIFL (output) DOUBLE PRECISION array, dimension ( N ) -* On exit, DIFL(I) is the distance between I-th updated -* (undeflated) singular value and the I-th (undeflated) old -* singular value. -* -* DIFR (output) DOUBLE PRECISION array, -* dimension ( LDGNUM, 2 ) if ICOMPQ = 1 and -* dimension ( N ) if ICOMPQ = 0. -* On exit, DIFR(I, 1) is the distance between I-th updated -* (undeflated) singular value and the I+1-th (undeflated) old -* singular value. -* -* If ICOMPQ = 1, DIFR(1:K,2) is an array containing the -* normalizing factors for the right singular vector matrix. -* -* See DLASD8 for details on DIFL and DIFR. -* -* Z (output) DOUBLE PRECISION array, dimension ( M ) -* The first elements of this array contain the components -* of the deflation-adjusted updating row vector. -* -* K (output) INTEGER -* Contains the dimension of the non-deflated matrix, -* This is the order of the related secular equation. 1 <= K <=N. -* -* C (output) DOUBLE PRECISION -* C contains garbage if SQRE =0 and the C-value of a Givens -* rotation related to the right null space if SQRE = 1. -* -* S (output) DOUBLE PRECISION -* S contains garbage if SQRE =0 and the S-value of a Givens -* rotation related to the right null space if SQRE = 1. -* -* WORK (workspace) DOUBLE PRECISION array, dimension ( 4 * M ) -* -* IWORK (workspace) INTEGER array, dimension ( 3 * N ) -* -* INFO (output) INTEGER -* = 0: successful exit. -* < 0: if INFO = -i, the i-th argument had an illegal value. -* > 0: if INFO = 1, a singular value did not converge -* -* Further Details -* =============== -* -* Based on contributions by -* Ming Gu and Huan Ren, Computer Science Division, University of -* California at Berkeley, USA -* * ===================================================================== * * .. Parameters ..