--- rpl/lapack/lapack/dlasq1.f 2010/01/26 15:22:45 1.1.1.1
+++ rpl/lapack/lapack/dlasq1.f 2017/06/17 10:53:57 1.16
@@ -1,14 +1,117 @@
- SUBROUTINE DLASQ1( N, D, E, WORK, INFO )
+*> \brief \b DLASQ1 computes the singular values of a real square bidiagonal matrix. Used by sbdsqr.
+*
+* =========== DOCUMENTATION ===========
+*
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
+*
+*> \htmlonly
+*> Download DLASQ1 + dependencies
+*>
+*> [TGZ]
+*>
+*> [ZIP]
+*>
+*> [TXT]
+*> \endhtmlonly
+*
+* Definition:
+* ===========
+*
+* SUBROUTINE DLASQ1( N, D, E, WORK, INFO )
+*
+* .. Scalar Arguments ..
+* INTEGER INFO, N
+* ..
+* .. Array Arguments ..
+* DOUBLE PRECISION D( * ), E( * ), WORK( * )
+* ..
+*
+*
+*> \par Purpose:
+* =============
+*>
+*> \verbatim
+*>
+*> DLASQ1 computes the singular values of a real N-by-N bidiagonal
+*> matrix with diagonal D and off-diagonal E. The singular values
+*> are computed to high relative accuracy, in the absence of
+*> denormalization, underflow and overflow. The algorithm was first
+*> presented in
+*>
+*> "Accurate singular values and differential qd algorithms" by K. V.
+*> Fernando and B. N. Parlett, Numer. Math., Vol-67, No. 2, pp. 191-230,
+*> 1994,
+*>
+*> and the present implementation is described in "An implementation of
+*> the dqds Algorithm (Positive Case)", LAPACK Working Note.
+*> \endverbatim
+*
+* Arguments:
+* ==========
+*
+*> \param[in] N
+*> \verbatim
+*> N is INTEGER
+*> The number of rows and columns in the matrix. N >= 0.
+*> \endverbatim
+*>
+*> \param[in,out] D
+*> \verbatim
+*> D is DOUBLE PRECISION array, dimension (N)
+*> On entry, D contains the diagonal elements of the
+*> bidiagonal matrix whose SVD is desired. On normal exit,
+*> D contains the singular values in decreasing order.
+*> \endverbatim
+*>
+*> \param[in,out] E
+*> \verbatim
+*> E is DOUBLE PRECISION array, dimension (N)
+*> On entry, elements E(1:N-1) contain the off-diagonal elements
+*> of the bidiagonal matrix whose SVD is desired.
+*> On exit, E is overwritten.
+*> \endverbatim
+*>
+*> \param[out] WORK
+*> \verbatim
+*> WORK is DOUBLE PRECISION array, dimension (4*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: the algorithm failed
+*> = 1, a split was marked by a positive value in E
+*> = 2, current block of Z not diagonalized after 100*N
+*> iterations (in inner while loop) On exit D and E
+*> represent a matrix with the same singular values
+*> which the calling subroutine could use to finish the
+*> computation, or even feed back into DLASQ1
+*> = 3, termination criterion of outer while loop not met
+*> (program created more than N unreduced blocks)
+*> \endverbatim
+*
+* Authors:
+* ========
+*
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
*
-* -- LAPACK routine (version 3.2) --
+*> \date December 2016
*
-* -- Contributed by Osni Marques of the Lawrence Berkeley National --
-* -- Laboratory and Beresford Parlett of the Univ. of California at --
-* -- Berkeley --
-* -- November 2008 --
+*> \ingroup auxOTHERcomputational
+*
+* =====================================================================
+ SUBROUTINE DLASQ1( N, D, E, WORK, INFO )
*
+* -- LAPACK computational routine (version 3.7.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
+* December 2016
*
* .. Scalar Arguments ..
INTEGER INFO, N
@@ -17,50 +120,6 @@
DOUBLE PRECISION D( * ), E( * ), WORK( * )
* ..
*
-* Purpose
-* =======
-*
-* DLASQ1 computes the singular values of a real N-by-N bidiagonal
-* matrix with diagonal D and off-diagonal E. The singular values
-* are computed to high relative accuracy, in the absence of
-* denormalization, underflow and overflow. The algorithm was first
-* presented in
-*
-* "Accurate singular values and differential qd algorithms" by K. V.
-* Fernando and B. N. Parlett, Numer. Math., Vol-67, No. 2, pp. 191-230,
-* 1994,
-*
-* and the present implementation is described in "An implementation of
-* the dqds Algorithm (Positive Case)", LAPACK Working Note.
-*
-* Arguments
-* =========
-*
-* N (input) INTEGER
-* The number of rows and columns in the matrix. N >= 0.
-*
-* D (input/output) DOUBLE PRECISION array, dimension (N)
-* On entry, D contains the diagonal elements of the
-* bidiagonal matrix whose SVD is desired. On normal exit,
-* D contains the singular values in decreasing order.
-*
-* E (input/output) DOUBLE PRECISION array, dimension (N)
-* On entry, elements E(1:N-1) contain the off-diagonal elements
-* of the bidiagonal matrix whose SVD is desired.
-* On exit, E is overwritten.
-*
-* WORK (workspace) DOUBLE PRECISION array, dimension (4*N)
-*
-* INFO (output) INTEGER
-* = 0: successful exit
-* < 0: if INFO = -i, the i-th argument had an illegal value
-* > 0: the algorithm failed
-* = 1, a split was marked by a positive value in E
-* = 2, current block of Z not diagonalized after 30*N
-* iterations (in inner while loop)
-* = 3, termination criterion of outer while loop not met
-* (program created more than N unreduced blocks)
-*
* =====================================================================
*
* .. Parameters ..
@@ -85,7 +144,7 @@
*
INFO = 0
IF( N.LT.0 ) THEN
- INFO = -2
+ INFO = -1
CALL XERBLA( 'DLASQ1', -INFO )
RETURN
ELSE IF( N.EQ.0 ) THEN
@@ -130,7 +189,7 @@
CALL DCOPY( N-1, E, 1, WORK( 2 ), 2 )
CALL DLASCL( 'G', 0, 0, SIGMX, SCALE, 2*N-1, 1, WORK, 2*N-1,
$ IINFO )
-*
+*
* Compute the q's and e's.
*
DO 30 I = 1, 2*N - 1
@@ -145,6 +204,17 @@
D( I ) = SQRT( WORK( I ) )
40 CONTINUE
CALL DLASCL( 'G', 0, 0, SCALE, SIGMX, N, 1, D, N, IINFO )
+ ELSE IF( INFO.EQ.2 ) THEN
+*
+* Maximum number of iterations exceeded. Move data from WORK
+* into D and E so the calling subroutine can try to finish
+*
+ DO I = 1, N
+ D( I ) = SQRT( WORK( 2*I-1 ) )
+ E( I ) = SQRT( WORK( 2*I ) )
+ END DO
+ CALL DLASCL( 'G', 0, 0, SCALE, SIGMX, N, 1, D, N, IINFO )
+ CALL DLASCL( 'G', 0, 0, SCALE, SIGMX, N, 1, E, N, IINFO )
END IF
*
RETURN