Annotation of rpl/lapack/lapack/zggqrf.f, revision 1.2
1.1 bertrand 1: SUBROUTINE ZGGQRF( N, M, P, A, LDA, TAUA, B, LDB, TAUB, WORK,
2: $ LWORK, INFO )
3: *
4: * -- LAPACK routine (version 3.2) --
5: * -- LAPACK is a software package provided by Univ. of Tennessee, --
6: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7: * November 2006
8: *
9: * .. Scalar Arguments ..
10: INTEGER INFO, LDA, LDB, LWORK, M, N, P
11: * ..
12: * .. Array Arguments ..
13: COMPLEX*16 A( LDA, * ), B( LDB, * ), TAUA( * ), TAUB( * ),
14: $ WORK( * )
15: * ..
16: *
17: * Purpose
18: * =======
19: *
20: * ZGGQRF computes a generalized QR factorization of an N-by-M matrix A
21: * and an N-by-P matrix B:
22: *
23: * A = Q*R, B = Q*T*Z,
24: *
25: * where Q is an N-by-N unitary matrix, Z is a P-by-P unitary matrix,
26: * and R and T assume one of the forms:
27: *
28: * if N >= M, R = ( R11 ) M , or if N < M, R = ( R11 R12 ) N,
29: * ( 0 ) N-M N M-N
30: * M
31: *
32: * where R11 is upper triangular, and
33: *
34: * if N <= P, T = ( 0 T12 ) N, or if N > P, T = ( T11 ) N-P,
35: * P-N N ( T21 ) P
36: * P
37: *
38: * where T12 or T21 is upper triangular.
39: *
40: * In particular, if B is square and nonsingular, the GQR factorization
41: * of A and B implicitly gives the QR factorization of inv(B)*A:
42: *
43: * inv(B)*A = Z'*(inv(T)*R)
44: *
45: * where inv(B) denotes the inverse of the matrix B, and Z' denotes the
46: * conjugate transpose of matrix Z.
47: *
48: * Arguments
49: * =========
50: *
51: * N (input) INTEGER
52: * The number of rows of the matrices A and B. N >= 0.
53: *
54: * M (input) INTEGER
55: * The number of columns of the matrix A. M >= 0.
56: *
57: * P (input) INTEGER
58: * The number of columns of the matrix B. P >= 0.
59: *
60: * A (input/output) COMPLEX*16 array, dimension (LDA,M)
61: * On entry, the N-by-M matrix A.
62: * On exit, the elements on and above the diagonal of the array
63: * contain the min(N,M)-by-M upper trapezoidal matrix R (R is
64: * upper triangular if N >= M); the elements below the diagonal,
65: * with the array TAUA, represent the unitary matrix Q as a
66: * product of min(N,M) elementary reflectors (see Further
67: * Details).
68: *
69: * LDA (input) INTEGER
70: * The leading dimension of the array A. LDA >= max(1,N).
71: *
72: * TAUA (output) COMPLEX*16 array, dimension (min(N,M))
73: * The scalar factors of the elementary reflectors which
74: * represent the unitary matrix Q (see Further Details).
75: *
76: * B (input/output) COMPLEX*16 array, dimension (LDB,P)
77: * On entry, the N-by-P matrix B.
78: * On exit, if N <= P, the upper triangle of the subarray
79: * B(1:N,P-N+1:P) contains the N-by-N upper triangular matrix T;
80: * if N > P, the elements on and above the (N-P)-th subdiagonal
81: * contain the N-by-P upper trapezoidal matrix T; the remaining
82: * elements, with the array TAUB, represent the unitary
83: * matrix Z as a product of elementary reflectors (see Further
84: * Details).
85: *
86: * LDB (input) INTEGER
87: * The leading dimension of the array B. LDB >= max(1,N).
88: *
89: * TAUB (output) COMPLEX*16 array, dimension (min(N,P))
90: * The scalar factors of the elementary reflectors which
91: * represent the unitary matrix Z (see Further Details).
92: *
93: * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
94: * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
95: *
96: * LWORK (input) INTEGER
97: * The dimension of the array WORK. LWORK >= max(1,N,M,P).
98: * For optimum performance LWORK >= max(N,M,P)*max(NB1,NB2,NB3),
99: * where NB1 is the optimal blocksize for the QR factorization
100: * of an N-by-M matrix, NB2 is the optimal blocksize for the
101: * RQ factorization of an N-by-P matrix, and NB3 is the optimal
102: * blocksize for a call of ZUNMQR.
103: *
104: * If LWORK = -1, then a workspace query is assumed; the routine
105: * only calculates the optimal size of the WORK array, returns
106: * this value as the first entry of the WORK array, and no error
107: * message related to LWORK is issued by XERBLA.
108: *
109: * INFO (output) INTEGER
110: * = 0: successful exit
111: * < 0: if INFO = -i, the i-th argument had an illegal value.
112: *
113: * Further Details
114: * ===============
115: *
116: * The matrix Q is represented as a product of elementary reflectors
117: *
118: * Q = H(1) H(2) . . . H(k), where k = min(n,m).
119: *
120: * Each H(i) has the form
121: *
122: * H(i) = I - taua * v * v'
123: *
124: * where taua is a complex scalar, and v is a complex vector with
125: * v(1:i-1) = 0 and v(i) = 1; v(i+1:n) is stored on exit in A(i+1:n,i),
126: * and taua in TAUA(i).
127: * To form Q explicitly, use LAPACK subroutine ZUNGQR.
128: * To use Q to update another matrix, use LAPACK subroutine ZUNMQR.
129: *
130: * The matrix Z is represented as a product of elementary reflectors
131: *
132: * Z = H(1) H(2) . . . H(k), where k = min(n,p).
133: *
134: * Each H(i) has the form
135: *
136: * H(i) = I - taub * v * v'
137: *
138: * where taub is a complex scalar, and v is a complex vector with
139: * v(p-k+i+1:p) = 0 and v(p-k+i) = 1; v(1:p-k+i-1) is stored on exit in
140: * B(n-k+i,1:p-k+i-1), and taub in TAUB(i).
141: * To form Z explicitly, use LAPACK subroutine ZUNGRQ.
142: * To use Z to update another matrix, use LAPACK subroutine ZUNMRQ.
143: *
144: * =====================================================================
145: *
146: * .. Local Scalars ..
147: LOGICAL LQUERY
148: INTEGER LOPT, LWKOPT, NB, NB1, NB2, NB3
149: * ..
150: * .. External Subroutines ..
151: EXTERNAL XERBLA, ZGEQRF, ZGERQF, ZUNMQR
152: * ..
153: * .. External Functions ..
154: INTEGER ILAENV
155: EXTERNAL ILAENV
156: * ..
157: * .. Intrinsic Functions ..
158: INTRINSIC INT, MAX, MIN
159: * ..
160: * .. Executable Statements ..
161: *
162: * Test the input parameters
163: *
164: INFO = 0
165: NB1 = ILAENV( 1, 'ZGEQRF', ' ', N, M, -1, -1 )
166: NB2 = ILAENV( 1, 'ZGERQF', ' ', N, P, -1, -1 )
167: NB3 = ILAENV( 1, 'ZUNMQR', ' ', N, M, P, -1 )
168: NB = MAX( NB1, NB2, NB3 )
169: LWKOPT = MAX( N, M, P )*NB
170: WORK( 1 ) = LWKOPT
171: LQUERY = ( LWORK.EQ.-1 )
172: IF( N.LT.0 ) THEN
173: INFO = -1
174: ELSE IF( M.LT.0 ) THEN
175: INFO = -2
176: ELSE IF( P.LT.0 ) THEN
177: INFO = -3
178: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
179: INFO = -5
180: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
181: INFO = -8
182: ELSE IF( LWORK.LT.MAX( 1, N, M, P ) .AND. .NOT.LQUERY ) THEN
183: INFO = -11
184: END IF
185: IF( INFO.NE.0 ) THEN
186: CALL XERBLA( 'ZGGQRF', -INFO )
187: RETURN
188: ELSE IF( LQUERY ) THEN
189: RETURN
190: END IF
191: *
192: * QR factorization of N-by-M matrix A: A = Q*R
193: *
194: CALL ZGEQRF( N, M, A, LDA, TAUA, WORK, LWORK, INFO )
195: LOPT = WORK( 1 )
196: *
197: * Update B := Q'*B.
198: *
199: CALL ZUNMQR( 'Left', 'Conjugate Transpose', N, P, MIN( N, M ), A,
200: $ LDA, TAUA, B, LDB, WORK, LWORK, INFO )
201: LOPT = MAX( LOPT, INT( WORK( 1 ) ) )
202: *
203: * RQ factorization of N-by-P matrix B: B = T*Z.
204: *
205: CALL ZGERQF( N, P, B, LDB, TAUB, WORK, LWORK, INFO )
206: WORK( 1 ) = MAX( LOPT, INT( WORK( 1 ) ) )
207: *
208: RETURN
209: *
210: * End of ZGGQRF
211: *
212: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to zggqrf.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack