[BACK]Return to dgejsv.f CVS log [TXT] [DIR] Up to [local] / rpl / lapack / lapack

Diff for /rpl/lapack/lapack/dgejsv.f between versions 1.15 and 1.16

version 1.15, 2016/08/27 15:34:21 version 1.16, 2017/06/17 10:53:48
Line 2 Line 2
 *  *
 *  =========== DOCUMENTATION ===========  *  =========== DOCUMENTATION ===========
 *  *
 * Online html documentation available at   * Online html documentation available at
 *            http://www.netlib.org/lapack/explore-html/   *            http://www.netlib.org/lapack/explore-html/
 *  *
 *> \htmlonly  *> \htmlonly
 *> Download DGEJSV + dependencies   *> Download DGEJSV + dependencies
 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgejsv.f">   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgejsv.f">
 *> [TGZ]</a>   *> [TGZ]</a>
 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgejsv.f">   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgejsv.f">
 *> [ZIP]</a>   *> [ZIP]</a>
 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgejsv.f">   *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgejsv.f">
 *> [TXT]</a>  *> [TXT]</a>
 *> \endhtmlonly   *> \endhtmlonly
 *  *
 *  Definition:  *  Definition:
 *  ===========  *  ===========
Line 21 Line 21
 *       SUBROUTINE DGEJSV( JOBA, JOBU, JOBV, JOBR, JOBT, JOBP,  *       SUBROUTINE DGEJSV( JOBA, JOBU, JOBV, JOBR, JOBT, JOBP,
 *                          M, N, A, LDA, SVA, U, LDU, V, LDV,  *                          M, N, A, LDA, SVA, U, LDU, V, LDV,
 *                          WORK, LWORK, IWORK, INFO )  *                          WORK, LWORK, IWORK, INFO )
 *   *
 *       .. Scalar Arguments ..  *       .. Scalar Arguments ..
 *       IMPLICIT    NONE  *       IMPLICIT    NONE
 *       INTEGER     INFO, LDA, LDU, LDV, LWORK, M, N  *       INTEGER     INFO, LDA, LDU, LDV, LWORK, M, N
Line 32 Line 32
 *       INTEGER     IWORK( * )  *       INTEGER     IWORK( * )
 *       CHARACTER*1 JOBA, JOBP, JOBR, JOBT, JOBU, JOBV  *       CHARACTER*1 JOBA, JOBP, JOBR, JOBT, JOBU, JOBV
 *       ..  *       ..
 *    *
 *  *
 *> \par Purpose:  *> \par Purpose:
 *  =============  *  =============
Line 87 Line 87
 *>             rows, then using this condition number gives too pessimistic  *>             rows, then using this condition number gives too pessimistic
 *>             error bound.  *>             error bound.
 *>       = 'A': Small singular values are the noise and the matrix is treated  *>       = 'A': Small singular values are the noise and the matrix is treated
 *>             as numerically rank defficient. The error in the computed  *>             as numerically rank deficient. The error in the computed
 *>             singular values is bounded by f(m,n)*epsilon*||A||.  *>             singular values is bounded by f(m,n)*epsilon*||A||.
 *>             The computed SVD A = U * S * V^t restores A up to  *>             The computed SVD A = U * S * V^t restores A up to
 *>             f(m,n)*epsilon*||A||.  *>             f(m,n)*epsilon*||A||.
Line 272 Line 272
 *> \param[out] WORK  *> \param[out] WORK
 *> \verbatim  *> \verbatim
 *>          WORK is DOUBLE PRECISION array, dimension at least LWORK.  *>          WORK is DOUBLE PRECISION array, dimension at least LWORK.
 *>          On exit, if N.GT.0 .AND. M.GT.0 (else not referenced),   *>          On exit, if N.GT.0 .AND. M.GT.0 (else not referenced),
 *>          WORK(1) = SCALE = WORK(2) / WORK(1) is the scaling factor such  *>          WORK(1) = SCALE = WORK(2) / WORK(1) is the scaling factor such
 *>                    that SCALE*SVA(1:N) are the computed singular values  *>                    that SCALE*SVA(1:N) are the computed singular values
 *>                    of A. (See the description of SVA().)  *>                    of A. (See the description of SVA().)
Line 319 Line 319
 *>               ->> For optimal performance (blocked code) the optimal value  *>               ->> For optimal performance (blocked code) the optimal value
 *>               is LWORK >= max(2*M+N,3*N+(N+1)*NB,7). Here NB is the optimal  *>               is LWORK >= max(2*M+N,3*N+(N+1)*NB,7). Here NB is the optimal
 *>               block size for DGEQP3 and DGEQRF.  *>               block size for DGEQP3 and DGEQRF.
 *>               In general, optimal LWORK is computed as   *>               In general, optimal LWORK is computed as
 *>               LWORK >= max(2*M+N,N+LWORK(DGEQP3),N+LWORK(DGEQRF), 7).          *>               LWORK >= max(2*M+N,N+LWORK(DGEQP3),N+LWORK(DGEQRF), 7).
 *>            -> .. an estimate of the scaled condition number of A is  *>            -> .. an estimate of the scaled condition number of A is
 *>               required (JOBA='E', 'G'). In this case, LWORK is the maximum  *>               required (JOBA='E', 'G'). In this case, LWORK is the maximum
 *>               of the above and N*N+4*N, i.e. LWORK >= max(2*M+N,N*N+4*N,7).  *>               of the above and N*N+4*N, i.e. LWORK >= max(2*M+N,N*N+4*N,7).
 *>               ->> For optimal performance (blocked code) the optimal value   *>               ->> For optimal performance (blocked code) the optimal value
 *>               is LWORK >= max(2*M+N,3*N+(N+1)*NB, N*N+4*N, 7).  *>               is LWORK >= max(2*M+N,3*N+(N+1)*NB, N*N+4*N, 7).
 *>               In general, the optimal length LWORK is computed as  *>               In general, the optimal length LWORK is computed as
 *>               LWORK >= max(2*M+N,N+LWORK(DGEQP3),N+LWORK(DGEQRF),   *>               LWORK >= max(2*M+N,N+LWORK(DGEQP3),N+LWORK(DGEQRF),
 *>                                                     N+N*N+LWORK(DPOCON),7).  *>                                                     N+N*N+LWORK(DPOCON),7).
 *>  *>
 *>          If SIGMA and the right singular vectors are needed (JOBV.EQ.'V'),  *>          If SIGMA and the right singular vectors are needed (JOBV.EQ.'V'),
Line 335 Line 335
 *>            -> For optimal performance, LWORK >= max(2*M+N,3*N+(N+1)*NB,7),  *>            -> For optimal performance, LWORK >= max(2*M+N,3*N+(N+1)*NB,7),
 *>               where NB is the optimal block size for DGEQP3, DGEQRF, DGELQF,  *>               where NB is the optimal block size for DGEQP3, DGEQRF, DGELQF,
 *>               DORMLQ. In general, the optimal length LWORK is computed as  *>               DORMLQ. In general, the optimal length LWORK is computed as
 *>               LWORK >= max(2*M+N,N+LWORK(DGEQP3), N+LWORK(DPOCON),   *>               LWORK >= max(2*M+N,N+LWORK(DGEQP3), N+LWORK(DPOCON),
 *>                       N+LWORK(DGELQF), 2*N+LWORK(DGEQRF), N+LWORK(DORMLQ)).  *>                       N+LWORK(DGELQF), 2*N+LWORK(DGEQRF), N+LWORK(DORMLQ)).
 *>  *>
 *>          If SIGMA and the left singular vectors are needed  *>          If SIGMA and the left singular vectors are needed
Line 346 Line 346
 *>               where NB is the optimal block size for DGEQP3, DGEQRF, DORMQR.  *>               where NB is the optimal block size for DGEQP3, DGEQRF, DORMQR.
 *>               In general, the optimal length LWORK is computed as  *>               In general, the optimal length LWORK is computed as
 *>               LWORK >= max(2*M+N,N+LWORK(DGEQP3),N+LWORK(DPOCON),  *>               LWORK >= max(2*M+N,N+LWORK(DGEQP3),N+LWORK(DPOCON),
 *>                        2*N+LWORK(DGEQRF), N+LWORK(DORMQR)).   *>                        2*N+LWORK(DGEQRF), N+LWORK(DORMQR)).
 *>               Here LWORK(DORMQR) equals N*NB (for JOBU.EQ.'U') or   *>               Here LWORK(DORMQR) equals N*NB (for JOBU.EQ.'U') or
 *>               M*NB (for JOBU.EQ.'F').  *>               M*NB (for JOBU.EQ.'F').
 *>                 *>
 *>          If the full SVD is needed: (JOBU.EQ.'U' or JOBU.EQ.'F') and   *>          If the full SVD is needed: (JOBU.EQ.'U' or JOBU.EQ.'F') and
 *>            -> if JOBV.EQ.'V'    *>            -> if JOBV.EQ.'V'
 *>               the minimal requirement is LWORK >= max(2*M+N,6*N+2*N*N).   *>               the minimal requirement is LWORK >= max(2*M+N,6*N+2*N*N).
 *>            -> if JOBV.EQ.'J' the minimal requirement is   *>            -> if JOBV.EQ.'J' the minimal requirement is
 *>               LWORK >= max(2*M+N, 4*N+N*N,2*N+N*N+6).  *>               LWORK >= max(2*M+N, 4*N+N*N,2*N+N*N+6).
 *>            -> For optimal performance, LWORK should be additionally  *>            -> For optimal performance, LWORK should be additionally
 *>               larger than N+M*NB, where NB is the optimal block size  *>               larger than N+M*NB, where NB is the optimal block size
Line 378 Line 378
 *> \verbatim  *> \verbatim
 *>          INFO is INTEGER  *>          INFO is INTEGER
 *>           < 0  : if INFO = -i, then the i-th argument had an illegal value.  *>           < 0  : if INFO = -i, then the i-th argument had an illegal value.
 *>           = 0 :  successfull exit;  *>           = 0 :  successful exit;
 *>           > 0 :  DGEJSV  did not converge in the maximal allowed number  *>           > 0 :  DGEJSV  did not converge in the maximal allowed number
 *>                  of sweeps. The computed values may be inaccurate.  *>                  of sweeps. The computed values may be inaccurate.
 *> \endverbatim  *> \endverbatim
Line 386 Line 386
 *  Authors:  *  Authors:
 *  ========  *  ========
 *  *
 *> \author Univ. of Tennessee   *> \author Univ. of Tennessee
 *> \author Univ. of California Berkeley   *> \author Univ. of California Berkeley
 *> \author Univ. of Colorado Denver   *> \author Univ. of Colorado Denver
 *> \author NAG Ltd.   *> \author NAG Ltd.
 *  *
 *> \date June 2016  *> \date June 2016
 *  *
Line 428 Line 428
 *>     The rank revealing QR factorization (in this code: DGEQP3) should be  *>     The rank revealing QR factorization (in this code: DGEQP3) should be
 *>  implemented as in [3]. We have a new version of DGEQP3 under development  *>  implemented as in [3]. We have a new version of DGEQP3 under development
 *>  that is more robust than the current one in LAPACK, with a cleaner cut in  *>  that is more robust than the current one in LAPACK, with a cleaner cut in
 *>  rank defficient cases. It will be available in the SIGMA library [4].  *>  rank deficient cases. It will be available in the SIGMA library [4].
 *>  If M is much larger than N, it is obvious that the inital QRF with  *>  If M is much larger than N, it is obvious that the initial QRF with
 *>  column pivoting can be preprocessed by the QRF without pivoting. That  *>  column pivoting can be preprocessed by the QRF without pivoting. That
 *>  well known trick is not used in DGEJSV because in some cases heavy row  *>  well known trick is not used in DGEJSV because in some cases heavy row
 *>  weighting can be treated with complete pivoting. The overhead in cases  *>  weighting can be treated with complete pivoting. The overhead in cases
Line 476 Line 476
      $                   M, N, A, LDA, SVA, U, LDU, V, LDV,       $                   M, N, A, LDA, SVA, U, LDU, V, LDV,
      $                   WORK, LWORK, IWORK, INFO )       $                   WORK, LWORK, IWORK, INFO )
 *  *
 *  -- LAPACK computational routine (version 3.6.1) --  *  -- LAPACK computational routine (version 3.7.0) --
 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --  *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--  *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 *     June 2016  *     June 2016
Line 562 Line 562
       ELSE IF ( LSVEC .AND. ( LDU .LT. M ) ) THEN        ELSE IF ( LSVEC .AND. ( LDU .LT. M ) ) THEN
          INFO = - 13           INFO = - 13
       ELSE IF ( RSVEC .AND. ( LDV .LT. N ) ) THEN        ELSE IF ( RSVEC .AND. ( LDV .LT. N ) ) THEN
          INFO = - 14           INFO = - 15
       ELSE IF ( (.NOT.(LSVEC .OR. RSVEC .OR. ERREST).AND.        ELSE IF ( (.NOT.(LSVEC .OR. RSVEC .OR. ERREST).AND.
      &                           (LWORK .LT. MAX(7,4*N+1,2*M+N))) .OR.       &                           (LWORK .LT. MAX(7,4*N+1,2*M+N))) .OR.
      & (.NOT.(LSVEC .OR. RSVEC) .AND. ERREST .AND.       & (.NOT.(LSVEC .OR. RSVEC) .AND. ERREST .AND.
Line 571 Line 571
      & .OR.       & .OR.
      & (RSVEC .AND. (.NOT.LSVEC) .AND. (LWORK .LT. MAX(7,2*M+N,4*N+1)))       & (RSVEC .AND. (.NOT.LSVEC) .AND. (LWORK .LT. MAX(7,2*M+N,4*N+1)))
      & .OR.       & .OR.
      & (LSVEC .AND. RSVEC .AND. (.NOT.JRACC) .AND.        & (LSVEC .AND. RSVEC .AND. (.NOT.JRACC) .AND.
      &                          (LWORK.LT.MAX(2*M+N,6*N+2*N*N)))       &                          (LWORK.LT.MAX(2*M+N,6*N+2*N*N)))
      & .OR. (LSVEC .AND. RSVEC .AND. JRACC .AND.       & .OR. (LSVEC .AND. RSVEC .AND. JRACC .AND.
      &                          LWORK.LT.MAX(2*M+N,4*N+N*N,2*N+N*N+6)))       &                          LWORK.LT.MAX(2*M+N,4*N+N*N,2*N+N*N+6)))
Line 967 Line 967
       ELSE IF ( L2RANK ) THEN        ELSE IF ( L2RANK ) THEN
 *        .. similarly as above, only slightly more gentle (less agressive).  *        .. similarly as above, only slightly more gentle (less agressive).
 *        Sudden drop on the diagonal of R1 is used as the criterion for  *        Sudden drop on the diagonal of R1 is used as the criterion for
 *        close-to-rank-defficient.  *        close-to-rank-deficient.
          TEMP1 = DSQRT(SFMIN)           TEMP1 = DSQRT(SFMIN)
          DO 3401 p = 2, N           DO 3401 p = 2, N
             IF ( ( DABS(A(p,p)) .LT. (EPSLN*DABS(A(p-1,p-1))) ) .OR.              IF ( ( DABS(A(p,p)) .LT. (EPSLN*DABS(A(p-1,p-1))) ) .OR.
Removed from v.1.15  
changed lines
  Added in v.1.16

CVSweb interface <joel.bertrand@systella.fr>