LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
dgsvj0.f
Go to the documentation of this file.
00001 *> \brief \b DGSVJ0
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *> \htmlonly
00009 *> Download DGSVJ0 + dependencies 
00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgsvj0.f"> 
00011 *> [TGZ]</a> 
00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgsvj0.f"> 
00013 *> [ZIP]</a> 
00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgsvj0.f"> 
00015 *> [TXT]</a>
00016 *> \endhtmlonly 
00017 *
00018 *  Definition:
00019 *  ===========
00020 *
00021 *       SUBROUTINE DGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
00022 *                          SFMIN, TOL, NSWEEP, WORK, LWORK, INFO )
00023 * 
00024 *       .. Scalar Arguments ..
00025 *       INTEGER            INFO, LDA, LDV, LWORK, M, MV, N, NSWEEP
00026 *       DOUBLE PRECISION   EPS, SFMIN, TOL
00027 *       CHARACTER*1        JOBV
00028 *       ..
00029 *       .. Array Arguments ..
00030 *       DOUBLE PRECISION   A( LDA, * ), SVA( N ), D( N ), V( LDV, * ),
00031 *      $                   WORK( LWORK )
00032 *       ..
00033 *  
00034 *
00035 *> \par Purpose:
00036 *  =============
00037 *>
00038 *> \verbatim
00039 *>
00040 *> DGSVJ0 is called from DGESVJ as a pre-processor and that is its main
00041 *> purpose. It applies Jacobi rotations in the same way as DGESVJ does, but
00042 *> it does not check convergence (stopping criterion). Few tuning
00043 *> parameters (marked by [TP]) are available for the implementer.
00044 *> \endverbatim
00045 *
00046 *  Arguments:
00047 *  ==========
00048 *
00049 *> \param[in] JOBV
00050 *> \verbatim
00051 *>          JOBV is CHARACTER*1
00052 *>          Specifies whether the output from this procedure is used
00053 *>          to compute the matrix V:
00054 *>          = 'V': the product of the Jacobi rotations is accumulated
00055 *>                 by postmulyiplying the N-by-N array V.
00056 *>                (See the description of V.)
00057 *>          = 'A': the product of the Jacobi rotations is accumulated
00058 *>                 by postmulyiplying the MV-by-N array V.
00059 *>                (See the descriptions of MV and V.)
00060 *>          = 'N': the Jacobi rotations are not accumulated.
00061 *> \endverbatim
00062 *>
00063 *> \param[in] M
00064 *> \verbatim
00065 *>          M is INTEGER
00066 *>          The number of rows of the input matrix A.  M >= 0.
00067 *> \endverbatim
00068 *>
00069 *> \param[in] N
00070 *> \verbatim
00071 *>          N is INTEGER
00072 *>          The number of columns of the input matrix A.
00073 *>          M >= N >= 0.
00074 *> \endverbatim
00075 *>
00076 *> \param[in,out] A
00077 *> \verbatim
00078 *>          A is DOUBLE PRECISION array, dimension (LDA,N)
00079 *>          On entry, M-by-N matrix A, such that A*diag(D) represents
00080 *>          the input matrix.
00081 *>          On exit,
00082 *>          A_onexit * D_onexit represents the input matrix A*diag(D)
00083 *>          post-multiplied by a sequence of Jacobi rotations, where the
00084 *>          rotation threshold and the total number of sweeps are given in
00085 *>          TOL and NSWEEP, respectively.
00086 *>          (See the descriptions of D, TOL and NSWEEP.)
00087 *> \endverbatim
00088 *>
00089 *> \param[in] LDA
00090 *> \verbatim
00091 *>          LDA is INTEGER
00092 *>          The leading dimension of the array A.  LDA >= max(1,M).
00093 *> \endverbatim
00094 *>
00095 *> \param[in,out] D
00096 *> \verbatim
00097 *>          D is DOUBLE PRECISION array, dimension (N)
00098 *>          The array D accumulates the scaling factors from the fast scaled
00099 *>          Jacobi rotations.
00100 *>          On entry, A*diag(D) represents the input matrix.
00101 *>          On exit, A_onexit*diag(D_onexit) represents the input matrix
00102 *>          post-multiplied by a sequence of Jacobi rotations, where the
00103 *>          rotation threshold and the total number of sweeps are given in
00104 *>          TOL and NSWEEP, respectively.
00105 *>          (See the descriptions of A, TOL and NSWEEP.)
00106 *> \endverbatim
00107 *>
00108 *> \param[in,out] SVA
00109 *> \verbatim
00110 *>          SVA is DOUBLE PRECISION array, dimension (N)
00111 *>          On entry, SVA contains the Euclidean norms of the columns of
00112 *>          the matrix A*diag(D).
00113 *>          On exit, SVA contains the Euclidean norms of the columns of
00114 *>          the matrix onexit*diag(D_onexit).
00115 *> \endverbatim
00116 *>
00117 *> \param[in] MV
00118 *> \verbatim
00119 *>          MV is INTEGER
00120 *>          If JOBV .EQ. 'A', then MV rows of V are post-multipled by a
00121 *>                           sequence of Jacobi rotations.
00122 *>          If JOBV = 'N',   then MV is not referenced.
00123 *> \endverbatim
00124 *>
00125 *> \param[in,out] V
00126 *> \verbatim
00127 *>          V is DOUBLE PRECISION array, dimension (LDV,N)
00128 *>          If JOBV .EQ. 'V' then N rows of V are post-multipled by a
00129 *>                           sequence of Jacobi rotations.
00130 *>          If JOBV .EQ. 'A' then MV rows of V are post-multipled by a
00131 *>                           sequence of Jacobi rotations.
00132 *>          If JOBV = 'N',   then V is not referenced.
00133 *> \endverbatim
00134 *>
00135 *> \param[in] LDV
00136 *> \verbatim
00137 *>          LDV is INTEGER
00138 *>          The leading dimension of the array V,  LDV >= 1.
00139 *>          If JOBV = 'V', LDV .GE. N.
00140 *>          If JOBV = 'A', LDV .GE. MV.
00141 *> \endverbatim
00142 *>
00143 *> \param[in] EPS
00144 *> \verbatim
00145 *>          EPS is DOUBLE PRECISION
00146 *>          EPS = DLAMCH('Epsilon')
00147 *> \endverbatim
00148 *>
00149 *> \param[in] SFMIN
00150 *> \verbatim
00151 *>          SFMIN is DOUBLE PRECISION
00152 *>          SFMIN = DLAMCH('Safe Minimum')
00153 *> \endverbatim
00154 *>
00155 *> \param[in] TOL
00156 *> \verbatim
00157 *>          TOL is DOUBLE PRECISION
00158 *>          TOL is the threshold for Jacobi rotations. For a pair
00159 *>          A(:,p), A(:,q) of pivot columns, the Jacobi rotation is
00160 *>          applied only if DABS(COS(angle(A(:,p),A(:,q)))) .GT. TOL.
00161 *> \endverbatim
00162 *>
00163 *> \param[in] NSWEEP
00164 *> \verbatim
00165 *>          NSWEEP is INTEGER
00166 *>          NSWEEP is the number of sweeps of Jacobi rotations to be
00167 *>          performed.
00168 *> \endverbatim
00169 *>
00170 *> \param[out] WORK
00171 *> \verbatim
00172 *>          WORK is DOUBLE PRECISION array, dimension (LWORK)
00173 *> \endverbatim
00174 *>
00175 *> \param[in] LWORK
00176 *> \verbatim
00177 *>          LWORK is INTEGER
00178 *>          LWORK is the dimension of WORK. LWORK .GE. M.
00179 *> \endverbatim
00180 *>
00181 *> \param[out] INFO
00182 *> \verbatim
00183 *>          INFO is INTEGER
00184 *>          = 0 : successful exit.
00185 *>          < 0 : if INFO = -i, then the i-th argument had an illegal value
00186 *> \endverbatim
00187 *
00188 *  Authors:
00189 *  ========
00190 *
00191 *> \author Univ. of Tennessee 
00192 *> \author Univ. of California Berkeley 
00193 *> \author Univ. of Colorado Denver 
00194 *> \author NAG Ltd. 
00195 *
00196 *> \date November 2011
00197 *
00198 *> \ingroup doubleOTHERcomputational
00199 *
00200 *> \par Further Details:
00201 *  =====================
00202 *>
00203 *> DGSVJ0 is used just to enable DGESVJ to call a simplified version of
00204 *> itself to work on a submatrix of the original matrix.
00205 *>
00206 *> \par Contributors:
00207 *  ==================
00208 *>
00209 *> Zlatko Drmac (Zagreb, Croatia) and Kresimir Veselic (Hagen, Germany)
00210 *>
00211 *> \par Bugs, Examples and Comments:
00212 *  =================================
00213 *>
00214 *> Please report all bugs and send interesting test examples and comments to
00215 *> drmac@math.hr. Thank you.
00216 *
00217 *  =====================================================================
00218       SUBROUTINE DGSVJ0( JOBV, M, N, A, LDA, D, SVA, MV, V, LDV, EPS,
00219      $                   SFMIN, TOL, NSWEEP, WORK, LWORK, INFO )
00220 *
00221 *  -- LAPACK computational routine (version 3.4.0) --
00222 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00223 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00224 *     November 2011
00225 *
00226 *     .. Scalar Arguments ..
00227       INTEGER            INFO, LDA, LDV, LWORK, M, MV, N, NSWEEP
00228       DOUBLE PRECISION   EPS, SFMIN, TOL
00229       CHARACTER*1        JOBV
00230 *     ..
00231 *     .. Array Arguments ..
00232       DOUBLE PRECISION   A( LDA, * ), SVA( N ), D( N ), V( LDV, * ),
00233      $                   WORK( LWORK )
00234 *     ..
00235 *
00236 *  =====================================================================
00237 *
00238 *     .. Local Parameters ..
00239       DOUBLE PRECISION   ZERO, HALF, ONE
00240       PARAMETER          ( ZERO = 0.0D0, HALF = 0.5D0, ONE = 1.0D0)
00241 *     ..
00242 *     .. Local Scalars ..
00243       DOUBLE PRECISION   AAPP, AAPP0, AAPQ, AAQQ, APOAQ, AQOAP, BIG,
00244      $                   BIGTHETA, CS, MXAAPQ, MXSINJ, ROOTBIG, ROOTEPS,
00245      $                   ROOTSFMIN, ROOTTOL, SMALL, SN, T, TEMP1, THETA,
00246      $                   THSIGN
00247       INTEGER            BLSKIP, EMPTSW, i, ibr, IERR, igl, IJBLSK, ir1,
00248      $                   ISWROT, jbc, jgl, KBL, LKAHEAD, MVL, NBL,
00249      $                   NOTROT, p, PSKIPPED, q, ROWSKIP, SWBAND
00250       LOGICAL            APPLV, ROTOK, RSVEC
00251 *     ..
00252 *     .. Local Arrays ..
00253       DOUBLE PRECISION   FASTR( 5 )
00254 *     ..
00255 *     .. Intrinsic Functions ..
00256       INTRINSIC          DABS, DMAX1, DBLE, MIN0, DSIGN, DSQRT
00257 *     ..
00258 *     .. External Functions ..
00259       DOUBLE PRECISION   DDOT, DNRM2
00260       INTEGER            IDAMAX
00261       LOGICAL            LSAME
00262       EXTERNAL           IDAMAX, LSAME, DDOT, DNRM2
00263 *     ..
00264 *     .. External Subroutines ..
00265       EXTERNAL           DAXPY, DCOPY, DLASCL, DLASSQ, DROTM, DSWAP
00266 *     ..
00267 *     .. Executable Statements ..
00268 *
00269 *     Test the input parameters.
00270 *
00271       APPLV = LSAME( JOBV, 'A' )
00272       RSVEC = LSAME( JOBV, 'V' )
00273       IF( .NOT.( RSVEC .OR. APPLV .OR. LSAME( JOBV, 'N' ) ) ) THEN
00274          INFO = -1
00275       ELSE IF( M.LT.0 ) THEN
00276          INFO = -2
00277       ELSE IF( ( N.LT.0 ) .OR. ( N.GT.M ) ) THEN
00278          INFO = -3
00279       ELSE IF( LDA.LT.M ) THEN
00280          INFO = -5
00281       ELSE IF( ( RSVEC.OR.APPLV ) .AND. ( MV.LT.0 ) ) THEN
00282          INFO = -8
00283       ELSE IF( ( RSVEC.AND.( LDV.LT.N ) ).OR. 
00284      $         ( APPLV.AND.( LDV.LT.MV ) ) ) THEN
00285          INFO = -10
00286       ELSE IF( TOL.LE.EPS ) THEN
00287          INFO = -13
00288       ELSE IF( NSWEEP.LT.0 ) THEN
00289          INFO = -14
00290       ELSE IF( LWORK.LT.M ) THEN
00291          INFO = -16
00292       ELSE
00293          INFO = 0
00294       END IF
00295 *
00296 *     #:(
00297       IF( INFO.NE.0 ) THEN
00298          CALL XERBLA( 'DGSVJ0', -INFO )
00299          RETURN
00300       END IF
00301 *
00302       IF( RSVEC ) THEN
00303          MVL = N
00304       ELSE IF( APPLV ) THEN
00305          MVL = MV
00306       END IF
00307       RSVEC = RSVEC .OR. APPLV
00308 
00309       ROOTEPS = DSQRT( EPS )
00310       ROOTSFMIN = DSQRT( SFMIN )
00311       SMALL = SFMIN / EPS
00312       BIG = ONE / SFMIN
00313       ROOTBIG = ONE / ROOTSFMIN
00314       BIGTHETA = ONE / ROOTEPS
00315       ROOTTOL = DSQRT( TOL )
00316 *
00317 *     -#- Row-cyclic Jacobi SVD algorithm with column pivoting -#-
00318 *
00319       EMPTSW = ( N*( N-1 ) ) / 2
00320       NOTROT = 0
00321       FASTR( 1 ) = ZERO
00322 *
00323 *     -#- Row-cyclic pivot strategy with de Rijk's pivoting -#-
00324 *
00325 
00326       SWBAND = 0
00327 *[TP] SWBAND is a tuning parameter. It is meaningful and effective
00328 *     if SGESVJ is used as a computational routine in the preconditioned
00329 *     Jacobi SVD algorithm SGESVJ. For sweeps i=1:SWBAND the procedure
00330 *     ......
00331 
00332       KBL = MIN0( 8, N )
00333 *[TP] KBL is a tuning parameter that defines the tile size in the
00334 *     tiling of the p-q loops of pivot pairs. In general, an optimal
00335 *     value of KBL depends on the matrix dimensions and on the
00336 *     parameters of the computer's memory.
00337 *
00338       NBL = N / KBL
00339       IF( ( NBL*KBL ).NE.N )NBL = NBL + 1
00340 
00341       BLSKIP = ( KBL**2 ) + 1
00342 *[TP] BLKSKIP is a tuning parameter that depends on SWBAND and KBL.
00343 
00344       ROWSKIP = MIN0( 5, KBL )
00345 *[TP] ROWSKIP is a tuning parameter.
00346 
00347       LKAHEAD = 1
00348 *[TP] LKAHEAD is a tuning parameter.
00349       SWBAND = 0
00350       PSKIPPED = 0
00351 *
00352       DO 1993 i = 1, NSWEEP
00353 *     .. go go go ...
00354 *
00355          MXAAPQ = ZERO
00356          MXSINJ = ZERO
00357          ISWROT = 0
00358 *
00359          NOTROT = 0
00360          PSKIPPED = 0
00361 *
00362          DO 2000 ibr = 1, NBL
00363 
00364             igl = ( ibr-1 )*KBL + 1
00365 *
00366             DO 1002 ir1 = 0, MIN0( LKAHEAD, NBL-ibr )
00367 *
00368                igl = igl + ir1*KBL
00369 *
00370                DO 2001 p = igl, MIN0( igl+KBL-1, N-1 )
00371 
00372 *     .. de Rijk's pivoting
00373                   q = IDAMAX( N-p+1, SVA( p ), 1 ) + p - 1
00374                   IF( p.NE.q ) THEN
00375                      CALL DSWAP( M, A( 1, p ), 1, A( 1, q ), 1 )
00376                      IF( RSVEC )CALL DSWAP( MVL, V( 1, p ), 1,
00377      $                                      V( 1, q ), 1 )
00378                      TEMP1 = SVA( p )
00379                      SVA( p ) = SVA( q )
00380                      SVA( q ) = TEMP1
00381                      TEMP1 = D( p )
00382                      D( p ) = D( q )
00383                      D( q ) = TEMP1
00384                   END IF
00385 *
00386                   IF( ir1.EQ.0 ) THEN
00387 *
00388 *        Column norms are periodically updated by explicit
00389 *        norm computation.
00390 *        Caveat:
00391 *        Some BLAS implementations compute DNRM2(M,A(1,p),1)
00392 *        as DSQRT(DDOT(M,A(1,p),1,A(1,p),1)), which may result in
00393 *        overflow for ||A(:,p)||_2 > DSQRT(overflow_threshold), and
00394 *        undeflow for ||A(:,p)||_2 < DSQRT(underflow_threshold).
00395 *        Hence, DNRM2 cannot be trusted, not even in the case when
00396 *        the true norm is far from the under(over)flow boundaries.
00397 *        If properly implemented DNRM2 is available, the IF-THEN-ELSE
00398 *        below should read "AAPP = DNRM2( M, A(1,p), 1 ) * D(p)".
00399 *
00400                      IF( ( SVA( p ).LT.ROOTBIG ) .AND.
00401      $                   ( SVA( p ).GT.ROOTSFMIN ) ) THEN
00402                         SVA( p ) = DNRM2( M, A( 1, p ), 1 )*D( p )
00403                      ELSE
00404                         TEMP1 = ZERO
00405                         AAPP = ONE
00406                         CALL DLASSQ( M, A( 1, p ), 1, TEMP1, AAPP )
00407                         SVA( p ) = TEMP1*DSQRT( AAPP )*D( p )
00408                      END IF
00409                      AAPP = SVA( p )
00410                   ELSE
00411                      AAPP = SVA( p )
00412                   END IF
00413 
00414 *
00415                   IF( AAPP.GT.ZERO ) THEN
00416 *
00417                      PSKIPPED = 0
00418 *
00419                      DO 2002 q = p + 1, MIN0( igl+KBL-1, N )
00420 *
00421                         AAQQ = SVA( q )
00422 
00423                         IF( AAQQ.GT.ZERO ) THEN
00424 *
00425                            AAPP0 = AAPP
00426                            IF( AAQQ.GE.ONE ) THEN
00427                               ROTOK = ( SMALL*AAPP ).LE.AAQQ
00428                               IF( AAPP.LT.( BIG / AAQQ ) ) THEN
00429                                  AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1,
00430      $                                  q ), 1 )*D( p )*D( q ) / AAQQ )
00431      $                                  / AAPP
00432                               ELSE
00433                                  CALL DCOPY( M, A( 1, p ), 1, WORK, 1 )
00434                                  CALL DLASCL( 'G', 0, 0, AAPP, D( p ),
00435      $                                        M, 1, WORK, LDA, IERR )
00436                                  AAPQ = DDOT( M, WORK, 1, A( 1, q ),
00437      $                                  1 )*D( q ) / AAQQ
00438                               END IF
00439                            ELSE
00440                               ROTOK = AAPP.LE.( AAQQ / SMALL )
00441                               IF( AAPP.GT.( SMALL / AAQQ ) ) THEN
00442                                  AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1,
00443      $                                  q ), 1 )*D( p )*D( q ) / AAQQ )
00444      $                                  / AAPP
00445                               ELSE
00446                                  CALL DCOPY( M, A( 1, q ), 1, WORK, 1 )
00447                                  CALL DLASCL( 'G', 0, 0, AAQQ, D( q ),
00448      $                                        M, 1, WORK, LDA, IERR )
00449                                  AAPQ = DDOT( M, WORK, 1, A( 1, p ),
00450      $                                  1 )*D( p ) / AAPP
00451                               END IF
00452                            END IF
00453 *
00454                            MXAAPQ = DMAX1( MXAAPQ, DABS( AAPQ ) )
00455 *
00456 *        TO rotate or NOT to rotate, THAT is the question ...
00457 *
00458                            IF( DABS( AAPQ ).GT.TOL ) THEN
00459 *
00460 *           .. rotate
00461 *           ROTATED = ROTATED + ONE
00462 *
00463                               IF( ir1.EQ.0 ) THEN
00464                                  NOTROT = 0
00465                                  PSKIPPED = 0
00466                                  ISWROT = ISWROT + 1
00467                               END IF
00468 *
00469                               IF( ROTOK ) THEN
00470 *
00471                                  AQOAP = AAQQ / AAPP
00472                                  APOAQ = AAPP / AAQQ
00473                                  THETA = -HALF*DABS( AQOAP-APOAQ )/AAPQ
00474 *
00475                                  IF( DABS( THETA ).GT.BIGTHETA ) THEN
00476 *
00477                                     T = HALF / THETA
00478                                     FASTR( 3 ) = T*D( p ) / D( q )
00479                                     FASTR( 4 ) = -T*D( q ) / D( p )
00480                                     CALL DROTM( M, A( 1, p ), 1,
00481      $                                          A( 1, q ), 1, FASTR )
00482                                     IF( RSVEC )CALL DROTM( MVL,
00483      $                                              V( 1, p ), 1,
00484      $                                              V( 1, q ), 1,
00485      $                                              FASTR )
00486                                     SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO,
00487      $                                         ONE+T*APOAQ*AAPQ ) )
00488                                     AAPP = AAPP*DSQRT( DMAX1( ZERO, 
00489      $                                     ONE-T*AQOAP*AAPQ ) )
00490                                     MXSINJ = DMAX1( MXSINJ, DABS( T ) )
00491 *
00492                                  ELSE
00493 *
00494 *                 .. choose correct signum for THETA and rotate
00495 *
00496                                     THSIGN = -DSIGN( ONE, AAPQ )
00497                                     T = ONE / ( THETA+THSIGN*
00498      $                                  DSQRT( ONE+THETA*THETA ) )
00499                                     CS = DSQRT( ONE / ( ONE+T*T ) )
00500                                     SN = T*CS
00501 *
00502                                     MXSINJ = DMAX1( MXSINJ, DABS( SN ) )
00503                                     SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO,
00504      $                                         ONE+T*APOAQ*AAPQ ) )
00505                                     AAPP = AAPP*DSQRT( DMAX1( ZERO,
00506      $                                     ONE-T*AQOAP*AAPQ ) )
00507 *
00508                                     APOAQ = D( p ) / D( q )
00509                                     AQOAP = D( q ) / D( p )
00510                                     IF( D( p ).GE.ONE ) THEN
00511                                        IF( D( q ).GE.ONE ) THEN
00512                                           FASTR( 3 ) = T*APOAQ
00513                                           FASTR( 4 ) = -T*AQOAP
00514                                           D( p ) = D( p )*CS
00515                                           D( q ) = D( q )*CS
00516                                           CALL DROTM( M, A( 1, p ), 1,
00517      $                                                A( 1, q ), 1,
00518      $                                                FASTR )
00519                                           IF( RSVEC )CALL DROTM( MVL,
00520      $                                        V( 1, p ), 1, V( 1, q ),
00521      $                                        1, FASTR )
00522                                        ELSE
00523                                           CALL DAXPY( M, -T*AQOAP,
00524      $                                                A( 1, q ), 1,
00525      $                                                A( 1, p ), 1 )
00526                                           CALL DAXPY( M, CS*SN*APOAQ,
00527      $                                                A( 1, p ), 1,
00528      $                                                A( 1, q ), 1 )
00529                                           D( p ) = D( p )*CS
00530                                           D( q ) = D( q ) / CS
00531                                           IF( RSVEC ) THEN
00532                                              CALL DAXPY( MVL, -T*AQOAP,
00533      $                                                   V( 1, q ), 1,
00534      $                                                   V( 1, p ), 1 )
00535                                              CALL DAXPY( MVL,
00536      $                                                   CS*SN*APOAQ,
00537      $                                                   V( 1, p ), 1,
00538      $                                                   V( 1, q ), 1 )
00539                                           END IF
00540                                        END IF
00541                                     ELSE
00542                                        IF( D( q ).GE.ONE ) THEN
00543                                           CALL DAXPY( M, T*APOAQ,
00544      $                                                A( 1, p ), 1,
00545      $                                                A( 1, q ), 1 )
00546                                           CALL DAXPY( M, -CS*SN*AQOAP,
00547      $                                                A( 1, q ), 1,
00548      $                                                A( 1, p ), 1 )
00549                                           D( p ) = D( p ) / CS
00550                                           D( q ) = D( q )*CS
00551                                           IF( RSVEC ) THEN
00552                                              CALL DAXPY( MVL, T*APOAQ,
00553      $                                                   V( 1, p ), 1,
00554      $                                                   V( 1, q ), 1 )
00555                                              CALL DAXPY( MVL,
00556      $                                                   -CS*SN*AQOAP,
00557      $                                                   V( 1, q ), 1,
00558      $                                                   V( 1, p ), 1 )
00559                                           END IF
00560                                        ELSE
00561                                           IF( D( p ).GE.D( q ) ) THEN
00562                                              CALL DAXPY( M, -T*AQOAP,
00563      $                                                   A( 1, q ), 1,
00564      $                                                   A( 1, p ), 1 )
00565                                              CALL DAXPY( M, CS*SN*APOAQ,
00566      $                                                   A( 1, p ), 1,
00567      $                                                   A( 1, q ), 1 )
00568                                              D( p ) = D( p )*CS
00569                                              D( q ) = D( q ) / CS
00570                                              IF( RSVEC ) THEN
00571                                                 CALL DAXPY( MVL,
00572      $                                               -T*AQOAP,
00573      $                                               V( 1, q ), 1,
00574      $                                               V( 1, p ), 1 )
00575                                                 CALL DAXPY( MVL,
00576      $                                               CS*SN*APOAQ,
00577      $                                               V( 1, p ), 1,
00578      $                                               V( 1, q ), 1 )
00579                                              END IF
00580                                           ELSE
00581                                              CALL DAXPY( M, T*APOAQ,
00582      $                                                   A( 1, p ), 1,
00583      $                                                   A( 1, q ), 1 )
00584                                              CALL DAXPY( M,
00585      $                                                   -CS*SN*AQOAP,
00586      $                                                   A( 1, q ), 1,
00587      $                                                   A( 1, p ), 1 )
00588                                              D( p ) = D( p ) / CS
00589                                              D( q ) = D( q )*CS
00590                                              IF( RSVEC ) THEN
00591                                                 CALL DAXPY( MVL,
00592      $                                               T*APOAQ, V( 1, p ),
00593      $                                               1, V( 1, q ), 1 )
00594                                                 CALL DAXPY( MVL,
00595      $                                               -CS*SN*AQOAP,
00596      $                                               V( 1, q ), 1,
00597      $                                               V( 1, p ), 1 )
00598                                              END IF
00599                                           END IF
00600                                        END IF
00601                                     END IF
00602                                  END IF
00603 *
00604                               ELSE
00605 *              .. have to use modified Gram-Schmidt like transformation
00606                                  CALL DCOPY( M, A( 1, p ), 1, WORK, 1 )
00607                                  CALL DLASCL( 'G', 0, 0, AAPP, ONE, M,
00608      $                                        1, WORK, LDA, IERR )
00609                                  CALL DLASCL( 'G', 0, 0, AAQQ, ONE, M,
00610      $                                        1, A( 1, q ), LDA, IERR )
00611                                  TEMP1 = -AAPQ*D( p ) / D( q )
00612                                  CALL DAXPY( M, TEMP1, WORK, 1,
00613      $                                       A( 1, q ), 1 )
00614                                  CALL DLASCL( 'G', 0, 0, ONE, AAQQ, M,
00615      $                                        1, A( 1, q ), LDA, IERR )
00616                                  SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO,
00617      $                                      ONE-AAPQ*AAPQ ) )
00618                                  MXSINJ = DMAX1( MXSINJ, SFMIN )
00619                               END IF
00620 *           END IF ROTOK THEN ... ELSE
00621 *
00622 *           In the case of cancellation in updating SVA(q), SVA(p)
00623 *           recompute SVA(q), SVA(p).
00624                               IF( ( SVA( q ) / AAQQ )**2.LE.ROOTEPS )
00625      $                            THEN
00626                                  IF( ( AAQQ.LT.ROOTBIG ) .AND.
00627      $                               ( AAQQ.GT.ROOTSFMIN ) ) THEN
00628                                     SVA( q ) = DNRM2( M, A( 1, q ), 1 )*
00629      $                                         D( q )
00630                                  ELSE
00631                                     T = ZERO
00632                                     AAQQ = ONE
00633                                     CALL DLASSQ( M, A( 1, q ), 1, T,
00634      $                                           AAQQ )
00635                                     SVA( q ) = T*DSQRT( AAQQ )*D( q )
00636                                  END IF
00637                               END IF
00638                               IF( ( AAPP / AAPP0 ).LE.ROOTEPS ) THEN
00639                                  IF( ( AAPP.LT.ROOTBIG ) .AND.
00640      $                               ( AAPP.GT.ROOTSFMIN ) ) THEN
00641                                     AAPP = DNRM2( M, A( 1, p ), 1 )*
00642      $                                     D( p )
00643                                  ELSE
00644                                     T = ZERO
00645                                     AAPP = ONE
00646                                     CALL DLASSQ( M, A( 1, p ), 1, T,
00647      $                                           AAPP )
00648                                     AAPP = T*DSQRT( AAPP )*D( p )
00649                                  END IF
00650                                  SVA( p ) = AAPP
00651                               END IF
00652 *
00653                            ELSE
00654 *        A(:,p) and A(:,q) already numerically orthogonal
00655                               IF( ir1.EQ.0 )NOTROT = NOTROT + 1
00656                               PSKIPPED = PSKIPPED + 1
00657                            END IF
00658                         ELSE
00659 *        A(:,q) is zero column
00660                            IF( ir1.EQ.0 )NOTROT = NOTROT + 1
00661                            PSKIPPED = PSKIPPED + 1
00662                         END IF
00663 *
00664                         IF( ( i.LE.SWBAND ) .AND.
00665      $                      ( PSKIPPED.GT.ROWSKIP ) ) THEN
00666                            IF( ir1.EQ.0 )AAPP = -AAPP
00667                            NOTROT = 0
00668                            GO TO 2103
00669                         END IF
00670 *
00671  2002                CONTINUE
00672 *     END q-LOOP
00673 *
00674  2103                CONTINUE
00675 *     bailed out of q-loop
00676 
00677                      SVA( p ) = AAPP
00678 
00679                   ELSE
00680                      SVA( p ) = AAPP
00681                      IF( ( ir1.EQ.0 ) .AND. ( AAPP.EQ.ZERO ) )
00682      $                   NOTROT = NOTROT + MIN0( igl+KBL-1, N ) - p
00683                   END IF
00684 *
00685  2001          CONTINUE
00686 *     end of the p-loop
00687 *     end of doing the block ( ibr, ibr )
00688  1002       CONTINUE
00689 *     end of ir1-loop
00690 *
00691 *........................................................
00692 * ... go to the off diagonal blocks
00693 *
00694             igl = ( ibr-1 )*KBL + 1
00695 *
00696             DO 2010 jbc = ibr + 1, NBL
00697 *
00698                jgl = ( jbc-1 )*KBL + 1
00699 *
00700 *        doing the block at ( ibr, jbc )
00701 *
00702                IJBLSK = 0
00703                DO 2100 p = igl, MIN0( igl+KBL-1, N )
00704 *
00705                   AAPP = SVA( p )
00706 *
00707                   IF( AAPP.GT.ZERO ) THEN
00708 *
00709                      PSKIPPED = 0
00710 *
00711                      DO 2200 q = jgl, MIN0( jgl+KBL-1, N )
00712 *
00713                         AAQQ = SVA( q )
00714 *
00715                         IF( AAQQ.GT.ZERO ) THEN
00716                            AAPP0 = AAPP
00717 *
00718 *     -#- M x 2 Jacobi SVD -#-
00719 *
00720 *        -#- Safe Gram matrix computation -#-
00721 *
00722                            IF( AAQQ.GE.ONE ) THEN
00723                               IF( AAPP.GE.AAQQ ) THEN
00724                                  ROTOK = ( SMALL*AAPP ).LE.AAQQ
00725                               ELSE
00726                                  ROTOK = ( SMALL*AAQQ ).LE.AAPP
00727                               END IF
00728                               IF( AAPP.LT.( BIG / AAQQ ) ) THEN
00729                                  AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1,
00730      $                                  q ), 1 )*D( p )*D( q ) / AAQQ )
00731      $                                  / AAPP
00732                               ELSE
00733                                  CALL DCOPY( M, A( 1, p ), 1, WORK, 1 )
00734                                  CALL DLASCL( 'G', 0, 0, AAPP, D( p ),
00735      $                                        M, 1, WORK, LDA, IERR )
00736                                  AAPQ = DDOT( M, WORK, 1, A( 1, q ),
00737      $                                  1 )*D( q ) / AAQQ
00738                               END IF
00739                            ELSE
00740                               IF( AAPP.GE.AAQQ ) THEN
00741                                  ROTOK = AAPP.LE.( AAQQ / SMALL )
00742                               ELSE
00743                                  ROTOK = AAQQ.LE.( AAPP / SMALL )
00744                               END IF
00745                               IF( AAPP.GT.( SMALL / AAQQ ) ) THEN
00746                                  AAPQ = ( DDOT( M, A( 1, p ), 1, A( 1,
00747      $                                  q ), 1 )*D( p )*D( q ) / AAQQ )
00748      $                                  / AAPP
00749                               ELSE
00750                                  CALL DCOPY( M, A( 1, q ), 1, WORK, 1 )
00751                                  CALL DLASCL( 'G', 0, 0, AAQQ, D( q ),
00752      $                                        M, 1, WORK, LDA, IERR )
00753                                  AAPQ = DDOT( M, WORK, 1, A( 1, p ),
00754      $                                  1 )*D( p ) / AAPP
00755                               END IF
00756                            END IF
00757 *
00758                            MXAAPQ = DMAX1( MXAAPQ, DABS( AAPQ ) )
00759 *
00760 *        TO rotate or NOT to rotate, THAT is the question ...
00761 *
00762                            IF( DABS( AAPQ ).GT.TOL ) THEN
00763                               NOTROT = 0
00764 *           ROTATED  = ROTATED + 1
00765                               PSKIPPED = 0
00766                               ISWROT = ISWROT + 1
00767 *
00768                               IF( ROTOK ) THEN
00769 *
00770                                  AQOAP = AAQQ / AAPP
00771                                  APOAQ = AAPP / AAQQ
00772                                  THETA = -HALF*DABS( AQOAP-APOAQ )/AAPQ
00773                                  IF( AAQQ.GT.AAPP0 )THETA = -THETA
00774 *
00775                                  IF( DABS( THETA ).GT.BIGTHETA ) THEN
00776                                     T = HALF / THETA
00777                                     FASTR( 3 ) = T*D( p ) / D( q )
00778                                     FASTR( 4 ) = -T*D( q ) / D( p )
00779                                     CALL DROTM( M, A( 1, p ), 1,
00780      $                                          A( 1, q ), 1, FASTR )
00781                                     IF( RSVEC )CALL DROTM( MVL,
00782      $                                              V( 1, p ), 1,
00783      $                                              V( 1, q ), 1,
00784      $                                              FASTR )
00785                                     SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO,
00786      $                                         ONE+T*APOAQ*AAPQ ) )
00787                                     AAPP = AAPP*DSQRT( DMAX1( ZERO,
00788      $                                     ONE-T*AQOAP*AAPQ ) )
00789                                     MXSINJ = DMAX1( MXSINJ, DABS( T ) )
00790                                  ELSE
00791 *
00792 *                 .. choose correct signum for THETA and rotate
00793 *
00794                                     THSIGN = -DSIGN( ONE, AAPQ )
00795                                     IF( AAQQ.GT.AAPP0 )THSIGN = -THSIGN
00796                                     T = ONE / ( THETA+THSIGN*
00797      $                                  DSQRT( ONE+THETA*THETA ) )
00798                                     CS = DSQRT( ONE / ( ONE+T*T ) )
00799                                     SN = T*CS
00800                                     MXSINJ = DMAX1( MXSINJ, DABS( SN ) )
00801                                     SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO,
00802      $                                         ONE+T*APOAQ*AAPQ ) )
00803                                     AAPP = AAPP*DSQRT( DMAX1( ZERO, 
00804      $                                     ONE-T*AQOAP*AAPQ ) )
00805 *
00806                                     APOAQ = D( p ) / D( q )
00807                                     AQOAP = D( q ) / D( p )
00808                                     IF( D( p ).GE.ONE ) THEN
00809 *
00810                                        IF( D( q ).GE.ONE ) THEN
00811                                           FASTR( 3 ) = T*APOAQ
00812                                           FASTR( 4 ) = -T*AQOAP
00813                                           D( p ) = D( p )*CS
00814                                           D( q ) = D( q )*CS
00815                                           CALL DROTM( M, A( 1, p ), 1,
00816      $                                                A( 1, q ), 1,
00817      $                                                FASTR )
00818                                           IF( RSVEC )CALL DROTM( MVL,
00819      $                                        V( 1, p ), 1, V( 1, q ),
00820      $                                        1, FASTR )
00821                                        ELSE
00822                                           CALL DAXPY( M, -T*AQOAP,
00823      $                                                A( 1, q ), 1,
00824      $                                                A( 1, p ), 1 )
00825                                           CALL DAXPY( M, CS*SN*APOAQ,
00826      $                                                A( 1, p ), 1,
00827      $                                                A( 1, q ), 1 )
00828                                           IF( RSVEC ) THEN
00829                                              CALL DAXPY( MVL, -T*AQOAP,
00830      $                                                   V( 1, q ), 1,
00831      $                                                   V( 1, p ), 1 )
00832                                              CALL DAXPY( MVL,
00833      $                                                   CS*SN*APOAQ,
00834      $                                                   V( 1, p ), 1,
00835      $                                                   V( 1, q ), 1 )
00836                                           END IF
00837                                           D( p ) = D( p )*CS
00838                                           D( q ) = D( q ) / CS
00839                                        END IF
00840                                     ELSE
00841                                        IF( D( q ).GE.ONE ) THEN
00842                                           CALL DAXPY( M, T*APOAQ,
00843      $                                                A( 1, p ), 1,
00844      $                                                A( 1, q ), 1 )
00845                                           CALL DAXPY( M, -CS*SN*AQOAP,
00846      $                                                A( 1, q ), 1,
00847      $                                                A( 1, p ), 1 )
00848                                           IF( RSVEC ) THEN
00849                                              CALL DAXPY( MVL, T*APOAQ,
00850      $                                                   V( 1, p ), 1,
00851      $                                                   V( 1, q ), 1 )
00852                                              CALL DAXPY( MVL,
00853      $                                                   -CS*SN*AQOAP,
00854      $                                                   V( 1, q ), 1,
00855      $                                                   V( 1, p ), 1 )
00856                                           END IF
00857                                           D( p ) = D( p ) / CS
00858                                           D( q ) = D( q )*CS
00859                                        ELSE
00860                                           IF( D( p ).GE.D( q ) ) THEN
00861                                              CALL DAXPY( M, -T*AQOAP,
00862      $                                                   A( 1, q ), 1,
00863      $                                                   A( 1, p ), 1 )
00864                                              CALL DAXPY( M, CS*SN*APOAQ,
00865      $                                                   A( 1, p ), 1,
00866      $                                                   A( 1, q ), 1 )
00867                                              D( p ) = D( p )*CS
00868                                              D( q ) = D( q ) / CS
00869                                              IF( RSVEC ) THEN
00870                                                 CALL DAXPY( MVL,
00871      $                                               -T*AQOAP,
00872      $                                               V( 1, q ), 1,
00873      $                                               V( 1, p ), 1 )
00874                                                 CALL DAXPY( MVL,
00875      $                                               CS*SN*APOAQ,
00876      $                                               V( 1, p ), 1,
00877      $                                               V( 1, q ), 1 )
00878                                              END IF
00879                                           ELSE
00880                                              CALL DAXPY( M, T*APOAQ,
00881      $                                                   A( 1, p ), 1,
00882      $                                                   A( 1, q ), 1 )
00883                                              CALL DAXPY( M,
00884      $                                                   -CS*SN*AQOAP,
00885      $                                                   A( 1, q ), 1,
00886      $                                                   A( 1, p ), 1 )
00887                                              D( p ) = D( p ) / CS
00888                                              D( q ) = D( q )*CS
00889                                              IF( RSVEC ) THEN
00890                                                 CALL DAXPY( MVL,
00891      $                                               T*APOAQ, V( 1, p ),
00892      $                                               1, V( 1, q ), 1 )
00893                                                 CALL DAXPY( MVL,
00894      $                                               -CS*SN*AQOAP,
00895      $                                               V( 1, q ), 1,
00896      $                                               V( 1, p ), 1 )
00897                                              END IF
00898                                           END IF
00899                                        END IF
00900                                     END IF
00901                                  END IF
00902 *
00903                               ELSE
00904                                  IF( AAPP.GT.AAQQ ) THEN
00905                                     CALL DCOPY( M, A( 1, p ), 1, WORK,
00906      $                                          1 )
00907                                     CALL DLASCL( 'G', 0, 0, AAPP, ONE,
00908      $                                           M, 1, WORK, LDA, IERR )
00909                                     CALL DLASCL( 'G', 0, 0, AAQQ, ONE,
00910      $                                           M, 1, A( 1, q ), LDA,
00911      $                                           IERR )
00912                                     TEMP1 = -AAPQ*D( p ) / D( q )
00913                                     CALL DAXPY( M, TEMP1, WORK, 1,
00914      $                                          A( 1, q ), 1 )
00915                                     CALL DLASCL( 'G', 0, 0, ONE, AAQQ,
00916      $                                           M, 1, A( 1, q ), LDA,
00917      $                                           IERR )
00918                                     SVA( q ) = AAQQ*DSQRT( DMAX1( ZERO,
00919      $                                         ONE-AAPQ*AAPQ ) )
00920                                     MXSINJ = DMAX1( MXSINJ, SFMIN )
00921                                  ELSE
00922                                     CALL DCOPY( M, A( 1, q ), 1, WORK,
00923      $                                          1 )
00924                                     CALL DLASCL( 'G', 0, 0, AAQQ, ONE,
00925      $                                           M, 1, WORK, LDA, IERR )
00926                                     CALL DLASCL( 'G', 0, 0, AAPP, ONE,
00927      $                                           M, 1, A( 1, p ), LDA,
00928      $                                           IERR )
00929                                     TEMP1 = -AAPQ*D( q ) / D( p )
00930                                     CALL DAXPY( M, TEMP1, WORK, 1,
00931      $                                          A( 1, p ), 1 )
00932                                     CALL DLASCL( 'G', 0, 0, ONE, AAPP,
00933      $                                           M, 1, A( 1, p ), LDA,
00934      $                                           IERR )
00935                                     SVA( p ) = AAPP*DSQRT( DMAX1( ZERO,
00936      $                                         ONE-AAPQ*AAPQ ) )
00937                                     MXSINJ = DMAX1( MXSINJ, SFMIN )
00938                                  END IF
00939                               END IF
00940 *           END IF ROTOK THEN ... ELSE
00941 *
00942 *           In the case of cancellation in updating SVA(q)
00943 *           .. recompute SVA(q)
00944                               IF( ( SVA( q ) / AAQQ )**2.LE.ROOTEPS )
00945      $                            THEN
00946                                  IF( ( AAQQ.LT.ROOTBIG ) .AND.
00947      $                               ( AAQQ.GT.ROOTSFMIN ) ) THEN
00948                                     SVA( q ) = DNRM2( M, A( 1, q ), 1 )*
00949      $                                         D( q )
00950                                  ELSE
00951                                     T = ZERO
00952                                     AAQQ = ONE
00953                                     CALL DLASSQ( M, A( 1, q ), 1, T,
00954      $                                           AAQQ )
00955                                     SVA( q ) = T*DSQRT( AAQQ )*D( q )
00956                                  END IF
00957                               END IF
00958                               IF( ( AAPP / AAPP0 )**2.LE.ROOTEPS ) THEN
00959                                  IF( ( AAPP.LT.ROOTBIG ) .AND.
00960      $                               ( AAPP.GT.ROOTSFMIN ) ) THEN
00961                                     AAPP = DNRM2( M, A( 1, p ), 1 )*
00962      $                                     D( p )
00963                                  ELSE
00964                                     T = ZERO
00965                                     AAPP = ONE
00966                                     CALL DLASSQ( M, A( 1, p ), 1, T,
00967      $                                           AAPP )
00968                                     AAPP = T*DSQRT( AAPP )*D( p )
00969                                  END IF
00970                                  SVA( p ) = AAPP
00971                               END IF
00972 *              end of OK rotation
00973                            ELSE
00974                               NOTROT = NOTROT + 1
00975                               PSKIPPED = PSKIPPED + 1
00976                               IJBLSK = IJBLSK + 1
00977                            END IF
00978                         ELSE
00979                            NOTROT = NOTROT + 1
00980                            PSKIPPED = PSKIPPED + 1
00981                            IJBLSK = IJBLSK + 1
00982                         END IF
00983 *
00984                         IF( ( i.LE.SWBAND ) .AND. ( IJBLSK.GE.BLSKIP ) )
00985      $                      THEN
00986                            SVA( p ) = AAPP
00987                            NOTROT = 0
00988                            GO TO 2011
00989                         END IF
00990                         IF( ( i.LE.SWBAND ) .AND.
00991      $                      ( PSKIPPED.GT.ROWSKIP ) ) THEN
00992                            AAPP = -AAPP
00993                            NOTROT = 0
00994                            GO TO 2203
00995                         END IF
00996 *
00997  2200                CONTINUE
00998 *        end of the q-loop
00999  2203                CONTINUE
01000 *
01001                      SVA( p ) = AAPP
01002 *
01003                   ELSE
01004                      IF( AAPP.EQ.ZERO )NOTROT = NOTROT +
01005      $                   MIN0( jgl+KBL-1, N ) - jgl + 1
01006                      IF( AAPP.LT.ZERO )NOTROT = 0
01007                   END IF
01008 
01009  2100          CONTINUE
01010 *     end of the p-loop
01011  2010       CONTINUE
01012 *     end of the jbc-loop
01013  2011       CONTINUE
01014 *2011 bailed out of the jbc-loop
01015             DO 2012 p = igl, MIN0( igl+KBL-1, N )
01016                SVA( p ) = DABS( SVA( p ) )
01017  2012       CONTINUE
01018 *
01019  2000    CONTINUE
01020 *2000 :: end of the ibr-loop
01021 *
01022 *     .. update SVA(N)
01023          IF( ( SVA( N ).LT.ROOTBIG ) .AND. ( SVA( N ).GT.ROOTSFMIN ) )
01024      $       THEN
01025             SVA( N ) = DNRM2( M, A( 1, N ), 1 )*D( N )
01026          ELSE
01027             T = ZERO
01028             AAPP = ONE
01029             CALL DLASSQ( M, A( 1, N ), 1, T, AAPP )
01030             SVA( N ) = T*DSQRT( AAPP )*D( N )
01031          END IF
01032 *
01033 *     Additional steering devices
01034 *
01035          IF( ( i.LT.SWBAND ) .AND. ( ( MXAAPQ.LE.ROOTTOL ) .OR.
01036      $       ( ISWROT.LE.N ) ) )SWBAND = i
01037 *
01038          IF( ( i.GT.SWBAND+1 ) .AND. ( MXAAPQ.LT.DBLE( N )*TOL ) .AND.
01039      $       ( DBLE( N )*MXAAPQ*MXSINJ.LT.TOL ) ) THEN
01040             GO TO 1994
01041          END IF
01042 *
01043          IF( NOTROT.GE.EMPTSW )GO TO 1994
01044 
01045  1993 CONTINUE
01046 *     end i=1:NSWEEP loop
01047 * #:) Reaching this point means that the procedure has comleted the given
01048 *     number of iterations.
01049       INFO = NSWEEP - 1
01050       GO TO 1995
01051  1994 CONTINUE
01052 * #:) Reaching this point means that during the i-th sweep all pivots were
01053 *     below the given tolerance, causing early exit.
01054 *
01055       INFO = 0
01056 * #:) INFO = 0 confirms successful iterations.
01057  1995 CONTINUE
01058 *
01059 *     Sort the vector D.
01060       DO 5991 p = 1, N - 1
01061          q = IDAMAX( N-p+1, SVA( p ), 1 ) + p - 1
01062          IF( p.NE.q ) THEN
01063             TEMP1 = SVA( p )
01064             SVA( p ) = SVA( q )
01065             SVA( q ) = TEMP1
01066             TEMP1 = D( p )
01067             D( p ) = D( q )
01068             D( q ) = TEMP1
01069             CALL DSWAP( M, A( 1, p ), 1, A( 1, q ), 1 )
01070             IF( RSVEC )CALL DSWAP( MVL, V( 1, p ), 1, V( 1, q ), 1 )
01071          END IF
01072  5991 CONTINUE
01073 *
01074       RETURN
01075 *     ..
01076 *     .. END OF DGSVJ0
01077 *     ..
01078       END
 All Files Functions