Diff for /rpl/lapack/lapack/dtgex2.f between versions 1.1 and 1.9

version 1.1, 2010/01/26 15:22:45 version 1.9, 2011/07/22 07:38:12
Line 1 Line 1
       SUBROUTINE DTGEX2( WANTQ, WANTZ, N, A, LDA, B, LDB, Q, LDQ, Z,        SUBROUTINE DTGEX2( WANTQ, WANTZ, N, A, LDA, B, LDB, Q, LDQ, Z,
      $                   LDZ, J1, N1, N2, WORK, LWORK, INFO )       $                   LDZ, J1, N1, N2, WORK, LWORK, INFO )
 *  *
 *  -- LAPACK auxiliary routine (version 3.2) --  *  -- LAPACK auxiliary routine (version 3.3.1) --
 *  -- 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..--
 *     November 2006  *  -- April 2011                                                      --
 *  *
 *     .. Scalar Arguments ..  *     .. Scalar Arguments ..
       LOGICAL            WANTQ, WANTZ        LOGICAL            WANTQ, WANTZ
Line 29 Line 29
 *  Optionally, the matrices Q and Z of generalized Schur vectors are  *  Optionally, the matrices Q and Z of generalized Schur vectors are
 *  updated.  *  updated.
 *  *
 *         Q(in) * A(in) * Z(in)' = Q(out) * A(out) * Z(out)'  *         Q(in) * A(in) * Z(in)**T = Q(out) * A(out) * Z(out)**T
 *         Q(in) * B(in) * Z(in)' = Q(out) * B(out) * Z(out)'  *         Q(in) * B(in) * Z(in)**T = Q(out) * B(out) * Z(out)**T
 *  *
 *  *
 *  Arguments  *  Arguments
Line 47 Line 47
 *  N       (input) INTEGER  *  N       (input) INTEGER
 *          The order of the matrices A and B. N >= 0.  *          The order of the matrices A and B. N >= 0.
 *  *
 *  A      (input/output) DOUBLE PRECISION arrays, dimensions (LDA,N)  *  A       (input/output) DOUBLE PRECISION array, dimensions (LDA,N)
 *          On entry, the matrix A in the pair (A, B).  *          On entry, the matrix A in the pair (A, B).
 *          On exit, the updated matrix A.  *          On exit, the updated matrix A.
 *  *
 *  LDA     (input)  INTEGER  *  LDA     (input) INTEGER
 *          The leading dimension of the array A. LDA >= max(1,N).  *          The leading dimension of the array A. LDA >= max(1,N).
 *  *
 *  B      (input/output) DOUBLE PRECISION arrays, dimensions (LDB,N)  *  B       (input/output) DOUBLE PRECISION array, dimensions (LDB,N)
 *          On entry, the matrix B in the pair (A, B).  *          On entry, the matrix B in the pair (A, B).
 *          On exit, the updated matrix B.  *          On exit, the updated matrix B.
 *  *
 *  LDB     (input)  INTEGER  *  LDB     (input) INTEGER
 *          The leading dimension of the array B. LDB >= max(1,N).  *          The leading dimension of the array B. LDB >= max(1,N).
 *  *
 *  Q       (input/output) DOUBLE PRECISION array, dimension (LDZ,N)  *  Q       (input/output) DOUBLE PRECISION array, dimension (LDQ,N)
 *          On entry, if WANTQ = .TRUE., the orthogonal matrix Q.  *          On entry, if WANTQ = .TRUE., the orthogonal matrix Q.
 *          On exit, the updated matrix Q.  *          On exit, the updated matrix Q.
 *          Not referenced if WANTQ = .FALSE..  *          Not referenced if WANTQ = .FALSE..
Line 134 Line 134
 *     .. Parameters ..  *     .. Parameters ..
       DOUBLE PRECISION   ZERO, ONE        DOUBLE PRECISION   ZERO, ONE
       PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )        PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )
       DOUBLE PRECISION   TEN        DOUBLE PRECISION   TWENTY
       PARAMETER          ( TEN = 1.0D+01 )        PARAMETER          ( TWENTY = 2.0D+01 )
       INTEGER            LDST        INTEGER            LDST
       PARAMETER          ( LDST = 4 )        PARAMETER          ( LDST = 4 )
       LOGICAL            WANDS        LOGICAL            WANDS
Line 205 Line 205
       CALL DLACPY( 'Full', M, M, T, LDST, WORK, M )        CALL DLACPY( 'Full', M, M, T, LDST, WORK, M )
       CALL DLASSQ( M*M, WORK, 1, DSCALE, DSUM )        CALL DLASSQ( M*M, WORK, 1, DSCALE, DSUM )
       DNORM = DSCALE*SQRT( DSUM )        DNORM = DSCALE*SQRT( DSUM )
       THRESH = MAX( TEN*EPS*DNORM, SMLNUM )  *
   *     THRES has been changed from 
   *        THRESH = MAX( TEN*EPS*SA, SMLNUM )
   *     to
   *        THRESH = MAX( TWENTY*EPS*SA, SMLNUM )
   *     on 04/01/10.
   *     "Bug" reported by Ondra Kamenik, confirmed by Julie Langou, fixed by
   *     Jim Demmel and Guillaume Revy. See forum post 1783.
   *
         THRESH = MAX( TWENTY*EPS*DNORM, SMLNUM )
 *  *
       IF( M.EQ.2 ) THEN        IF( M.EQ.2 ) THEN
 *  *
Line 250 Line 259
          IF( WANDS ) THEN           IF( WANDS ) THEN
 *  *
 *           Strong stability test:  *           Strong stability test:
 *             F-norm((A-QL'*S*QR, B-QL'*T*QR)) <= O(EPS*F-norm((A,B)))  *             F-norm((A-QL**T*S*QR, B-QL**T*T*QR)) <= O(EPS*F-norm((A,B)))
 *  *
             CALL DLACPY( 'Full', M, M, A( J1, J1 ), LDA, WORK( M*M+1 ),              CALL DLACPY( 'Full', M, M, A( J1, J1 ), LDA, WORK( M*M+1 ),
      $                   M )       $                   M )
Line 325 Line 334
 *  *
 *        Compute orthogonal matrix QL:  *        Compute orthogonal matrix QL:
 *  *
 *                    QL' * LI = [ TL ]  *                    QL**T * LI = [ TL ]
 *                               [ 0  ]  *                                 [ 0  ]
 *        where  *        where
 *                    LI =  [      -L              ]  *                    LI =  [      -L              ]
 *                          [ SCALE * identity(N2) ]  *                          [ SCALE * identity(N2) ]
Line 344 Line 353
 *  *
 *        Compute orthogonal matrix RQ:  *        Compute orthogonal matrix RQ:
 *  *
 *                    IR * RQ' =   [ 0  TR],  *                    IR * RQ**T =   [ 0  TR],
 *  *
 *         where IR = [ SCALE * identity(N1), R ]  *         where IR = [ SCALE * identity(N1), R ]
 *  *
Line 439 Line 448
          IF( WANDS ) THEN           IF( WANDS ) THEN
 *  *
 *           Strong stability test:  *           Strong stability test:
 *              F-norm((A-QL*S*QR', B-QL*T*QR')) <= O(EPS*F-norm((A,B)))  *              F-norm((A-QL*S*QR**T, B-QL*T*QR**T)) <= O(EPS*F-norm((A,B)))
 *  *
             CALL DLACPY( 'Full', M, M, A( J1, J1 ), LDA, WORK( M*M+1 ),              CALL DLACPY( 'Full', M, M, A( J1, J1 ), LDA, WORK( M*M+1 ),
      $                   M )       $                   M )

Removed from v.1.1  
changed lines
  Added in v.1.9


CVSweb interface <joel.bertrand@systella.fr>