version 1.7, 2010/12/21 13:53:27
|
version 1.8, 2011/07/22 07:38:05
|
Line 2
|
Line 2
|
$ TOLA, TOLB, K, L, U, LDU, V, LDV, Q, LDQ, |
$ TOLA, TOLB, K, L, U, LDU, V, LDV, Q, LDQ, |
$ IWORK, TAU, WORK, INFO ) |
$ IWORK, TAU, WORK, INFO ) |
* |
* |
* -- LAPACK routine (version 3.2) -- |
* -- LAPACK 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 .. |
CHARACTER JOBQ, JOBU, JOBV |
CHARACTER JOBQ, JOBU, JOBV |
Line 23
|
Line 23
|
* |
* |
* DGGSVP computes orthogonal matrices U, V and Q such that |
* DGGSVP computes orthogonal matrices U, V and Q such that |
* |
* |
* N-K-L K L |
* N-K-L K L |
* U'*A*Q = K ( 0 A12 A13 ) if M-K-L >= 0; |
* U**T*A*Q = K ( 0 A12 A13 ) if M-K-L >= 0; |
* L ( 0 0 A23 ) |
* L ( 0 0 A23 ) |
* M-K-L ( 0 0 0 ) |
* M-K-L ( 0 0 0 ) |
* |
* |
* N-K-L K L |
* N-K-L K L |
* = K ( 0 A12 A13 ) if M-K-L < 0; |
* = K ( 0 A12 A13 ) if M-K-L < 0; |
* M-K ( 0 0 A23 ) |
* M-K ( 0 0 A23 ) |
* |
* |
* N-K-L K L |
* N-K-L K L |
* V'*B*Q = L ( 0 0 B13 ) |
* V**T*B*Q = L ( 0 0 B13 ) |
* P-L ( 0 0 0 ) |
* P-L ( 0 0 0 ) |
* |
* |
* where the K-by-K matrix A12 and L-by-L matrix B13 are nonsingular |
* where the K-by-K matrix A12 and L-by-L matrix B13 are nonsingular |
* upper triangular; A23 is L-by-L upper triangular if M-K-L >= 0, |
* upper triangular; A23 is L-by-L upper triangular if M-K-L >= 0, |
* otherwise A23 is (M-K)-by-L upper trapezoidal. K+L = the effective |
* otherwise A23 is (M-K)-by-L upper trapezoidal. K+L = the effective |
* numerical rank of the (M+P)-by-N matrix (A',B')'. Z' denotes the |
* numerical rank of the (M+P)-by-N matrix (A**T,B**T)**T. |
* transpose of Z. |
|
* |
* |
* This decomposition is the preprocessing step for computing the |
* This decomposition is the preprocessing step for computing the |
* Generalized Singular Value Decomposition (GSVD), see subroutine |
* Generalized Singular Value Decomposition (GSVD), see subroutine |
Line 99
|
Line 98
|
* K (output) INTEGER |
* K (output) INTEGER |
* L (output) INTEGER |
* L (output) INTEGER |
* On exit, K and L specify the dimension of the subblocks |
* On exit, K and L specify the dimension of the subblocks |
* described in Purpose. |
* described in Purpose section. |
* K + L = effective numerical rank of (A',B')'. |
* K + L = effective numerical rank of (A**T,B**T)**T. |
* |
* |
* U (output) DOUBLE PRECISION array, dimension (LDU,M) |
* U (output) DOUBLE PRECISION array, dimension (LDU,M) |
* If JOBU = 'U', U contains the orthogonal matrix U. |
* If JOBU = 'U', U contains the orthogonal matrix U. |
Line 258
|
Line 257
|
* |
* |
CALL DGERQ2( L, N, B, LDB, TAU, WORK, INFO ) |
CALL DGERQ2( L, N, B, LDB, TAU, WORK, INFO ) |
* |
* |
* Update A := A*Z' |
* Update A := A*Z**T |
* |
* |
CALL DORMR2( 'Right', 'Transpose', M, N, L, B, LDB, TAU, A, |
CALL DORMR2( 'Right', 'Transpose', M, N, L, B, LDB, TAU, A, |
$ LDA, WORK, INFO ) |
$ LDA, WORK, INFO ) |
* |
* |
IF( WANTQ ) THEN |
IF( WANTQ ) THEN |
* |
* |
* Update Q := Q*Z' |
* Update Q := Q*Z**T |
* |
* |
CALL DORMR2( 'Right', 'Transpose', N, N, L, B, LDB, TAU, Q, |
CALL DORMR2( 'Right', 'Transpose', N, N, L, B, LDB, TAU, Q, |
$ LDQ, WORK, INFO ) |
$ LDQ, WORK, INFO ) |
Line 287
|
Line 286
|
* |
* |
* then the following does the complete QR decomposition of A11: |
* then the following does the complete QR decomposition of A11: |
* |
* |
* A11 = U*( 0 T12 )*P1' |
* A11 = U*( 0 T12 )*P1**T |
* ( 0 0 ) |
* ( 0 0 ) |
* |
* |
DO 70 I = 1, N - L |
DO 70 I = 1, N - L |
Line 303
|
Line 302
|
$ K = K + 1 |
$ K = K + 1 |
80 CONTINUE |
80 CONTINUE |
* |
* |
* Update A12 := U'*A12, where A12 = A( 1:M, N-L+1:N ) |
* Update A12 := U**T*A12, where A12 = A( 1:M, N-L+1:N ) |
* |
* |
CALL DORM2R( 'Left', 'Transpose', M, L, MIN( M, N-L ), A, LDA, |
CALL DORM2R( 'Left', 'Transpose', M, L, MIN( M, N-L ), A, LDA, |
$ TAU, A( 1, N-L+1 ), LDA, WORK, INFO ) |
$ TAU, A( 1, N-L+1 ), LDA, WORK, INFO ) |
Line 345
|
Line 344
|
* |
* |
IF( WANTQ ) THEN |
IF( WANTQ ) THEN |
* |
* |
* Update Q( 1:N,1:N-L ) = Q( 1:N,1:N-L )*Z1' |
* Update Q( 1:N,1:N-L ) = Q( 1:N,1:N-L )*Z1**T |
* |
* |
CALL DORMR2( 'Right', 'Transpose', N, N-L, K, A, LDA, TAU, |
CALL DORMR2( 'Right', 'Transpose', N, N-L, K, A, LDA, TAU, |
$ Q, LDQ, WORK, INFO ) |
$ Q, LDQ, WORK, INFO ) |