00001 SUBROUTINE ZUNMBR( VECT, SIDE, TRANS, M, N, K, A, LDA, TAU, C, 00002 $ LDC, WORK, LWORK, INFO ) 00003 * 00004 * -- LAPACK routine (version 3.2) -- 00005 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00006 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00007 * November 2006 00008 * 00009 * .. Scalar Arguments .. 00010 CHARACTER SIDE, TRANS, VECT 00011 INTEGER INFO, K, LDA, LDC, LWORK, M, N 00012 * .. 00013 * .. Array Arguments .. 00014 COMPLEX*16 A( LDA, * ), C( LDC, * ), TAU( * ), WORK( * ) 00015 * .. 00016 * 00017 * Purpose 00018 * ======= 00019 * 00020 * If VECT = 'Q', ZUNMBR overwrites the general complex M-by-N matrix C 00021 * with 00022 * SIDE = 'L' SIDE = 'R' 00023 * TRANS = 'N': Q * C C * Q 00024 * TRANS = 'C': Q**H * C C * Q**H 00025 * 00026 * If VECT = 'P', ZUNMBR overwrites the general complex M-by-N matrix C 00027 * with 00028 * SIDE = 'L' SIDE = 'R' 00029 * TRANS = 'N': P * C C * P 00030 * TRANS = 'C': P**H * C C * P**H 00031 * 00032 * Here Q and P**H are the unitary matrices determined by ZGEBRD when 00033 * reducing a complex matrix A to bidiagonal form: A = Q * B * P**H. Q 00034 * and P**H are defined as products of elementary reflectors H(i) and 00035 * G(i) respectively. 00036 * 00037 * Let nq = m if SIDE = 'L' and nq = n if SIDE = 'R'. Thus nq is the 00038 * order of the unitary matrix Q or P**H that is applied. 00039 * 00040 * If VECT = 'Q', A is assumed to have been an NQ-by-K matrix: 00041 * if nq >= k, Q = H(1) H(2) . . . H(k); 00042 * if nq < k, Q = H(1) H(2) . . . H(nq-1). 00043 * 00044 * If VECT = 'P', A is assumed to have been a K-by-NQ matrix: 00045 * if k < nq, P = G(1) G(2) . . . G(k); 00046 * if k >= nq, P = G(1) G(2) . . . G(nq-1). 00047 * 00048 * Arguments 00049 * ========= 00050 * 00051 * VECT (input) CHARACTER*1 00052 * = 'Q': apply Q or Q**H; 00053 * = 'P': apply P or P**H. 00054 * 00055 * SIDE (input) CHARACTER*1 00056 * = 'L': apply Q, Q**H, P or P**H from the Left; 00057 * = 'R': apply Q, Q**H, P or P**H from the Right. 00058 * 00059 * TRANS (input) CHARACTER*1 00060 * = 'N': No transpose, apply Q or P; 00061 * = 'C': Conjugate transpose, apply Q**H or P**H. 00062 * 00063 * M (input) INTEGER 00064 * The number of rows of the matrix C. M >= 0. 00065 * 00066 * N (input) INTEGER 00067 * The number of columns of the matrix C. N >= 0. 00068 * 00069 * K (input) INTEGER 00070 * If VECT = 'Q', the number of columns in the original 00071 * matrix reduced by ZGEBRD. 00072 * If VECT = 'P', the number of rows in the original 00073 * matrix reduced by ZGEBRD. 00074 * K >= 0. 00075 * 00076 * A (input) COMPLEX*16 array, dimension 00077 * (LDA,min(nq,K)) if VECT = 'Q' 00078 * (LDA,nq) if VECT = 'P' 00079 * The vectors which define the elementary reflectors H(i) and 00080 * G(i), whose products determine the matrices Q and P, as 00081 * returned by ZGEBRD. 00082 * 00083 * LDA (input) INTEGER 00084 * The leading dimension of the array A. 00085 * If VECT = 'Q', LDA >= max(1,nq); 00086 * if VECT = 'P', LDA >= max(1,min(nq,K)). 00087 * 00088 * TAU (input) COMPLEX*16 array, dimension (min(nq,K)) 00089 * TAU(i) must contain the scalar factor of the elementary 00090 * reflector H(i) or G(i) which determines Q or P, as returned 00091 * by ZGEBRD in the array argument TAUQ or TAUP. 00092 * 00093 * C (input/output) COMPLEX*16 array, dimension (LDC,N) 00094 * On entry, the M-by-N matrix C. 00095 * On exit, C is overwritten by Q*C or Q**H*C or C*Q**H or C*Q 00096 * or P*C or P**H*C or C*P or C*P**H. 00097 * 00098 * LDC (input) INTEGER 00099 * The leading dimension of the array C. LDC >= max(1,M). 00100 * 00101 * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK)) 00102 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK. 00103 * 00104 * LWORK (input) INTEGER 00105 * The dimension of the array WORK. 00106 * If SIDE = 'L', LWORK >= max(1,N); 00107 * if SIDE = 'R', LWORK >= max(1,M); 00108 * if N = 0 or M = 0, LWORK >= 1. 00109 * For optimum performance LWORK >= max(1,N*NB) if SIDE = 'L', 00110 * and LWORK >= max(1,M*NB) if SIDE = 'R', where NB is the 00111 * optimal blocksize. (NB = 0 if M = 0 or N = 0.) 00112 * 00113 * If LWORK = -1, then a workspace query is assumed; the routine 00114 * only calculates the optimal size of the WORK array, returns 00115 * this value as the first entry of the WORK array, and no error 00116 * message related to LWORK is issued by XERBLA. 00117 * 00118 * INFO (output) INTEGER 00119 * = 0: successful exit 00120 * < 0: if INFO = -i, the i-th argument had an illegal value 00121 * 00122 * ===================================================================== 00123 * 00124 * .. Local Scalars .. 00125 LOGICAL APPLYQ, LEFT, LQUERY, NOTRAN 00126 CHARACTER TRANST 00127 INTEGER I1, I2, IINFO, LWKOPT, MI, NB, NI, NQ, NW 00128 * .. 00129 * .. External Functions .. 00130 LOGICAL LSAME 00131 INTEGER ILAENV 00132 EXTERNAL LSAME, ILAENV 00133 * .. 00134 * .. External Subroutines .. 00135 EXTERNAL XERBLA, ZUNMLQ, ZUNMQR 00136 * .. 00137 * .. Intrinsic Functions .. 00138 INTRINSIC MAX, MIN 00139 * .. 00140 * .. Executable Statements .. 00141 * 00142 * Test the input arguments 00143 * 00144 INFO = 0 00145 APPLYQ = LSAME( VECT, 'Q' ) 00146 LEFT = LSAME( SIDE, 'L' ) 00147 NOTRAN = LSAME( TRANS, 'N' ) 00148 LQUERY = ( LWORK.EQ.-1 ) 00149 * 00150 * NQ is the order of Q or P and NW is the minimum dimension of WORK 00151 * 00152 IF( LEFT ) THEN 00153 NQ = M 00154 NW = N 00155 ELSE 00156 NQ = N 00157 NW = M 00158 END IF 00159 IF( M.EQ.0 .OR. N.EQ.0 ) THEN 00160 NW = 0 00161 END IF 00162 IF( .NOT.APPLYQ .AND. .NOT.LSAME( VECT, 'P' ) ) THEN 00163 INFO = -1 00164 ELSE IF( .NOT.LEFT .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN 00165 INFO = -2 00166 ELSE IF( .NOT.NOTRAN .AND. .NOT.LSAME( TRANS, 'C' ) ) THEN 00167 INFO = -3 00168 ELSE IF( M.LT.0 ) THEN 00169 INFO = -4 00170 ELSE IF( N.LT.0 ) THEN 00171 INFO = -5 00172 ELSE IF( K.LT.0 ) THEN 00173 INFO = -6 00174 ELSE IF( ( APPLYQ .AND. LDA.LT.MAX( 1, NQ ) ) .OR. 00175 $ ( .NOT.APPLYQ .AND. LDA.LT.MAX( 1, MIN( NQ, K ) ) ) ) 00176 $ THEN 00177 INFO = -8 00178 ELSE IF( LDC.LT.MAX( 1, M ) ) THEN 00179 INFO = -11 00180 ELSE IF( LWORK.LT.MAX( 1, NW ) .AND. .NOT.LQUERY ) THEN 00181 INFO = -13 00182 END IF 00183 * 00184 IF( INFO.EQ.0 ) THEN 00185 IF( NW.GT.0 ) THEN 00186 IF( APPLYQ ) THEN 00187 IF( LEFT ) THEN 00188 NB = ILAENV( 1, 'ZUNMQR', SIDE // TRANS, M-1, N, M-1, 00189 $ -1 ) 00190 ELSE 00191 NB = ILAENV( 1, 'ZUNMQR', SIDE // TRANS, M, N-1, N-1, 00192 $ -1 ) 00193 END IF 00194 ELSE 00195 IF( LEFT ) THEN 00196 NB = ILAENV( 1, 'ZUNMLQ', SIDE // TRANS, M-1, N, M-1, 00197 $ -1 ) 00198 ELSE 00199 NB = ILAENV( 1, 'ZUNMLQ', SIDE // TRANS, M, N-1, N-1, 00200 $ -1 ) 00201 END IF 00202 END IF 00203 LWKOPT = MAX( 1, NW*NB ) 00204 ELSE 00205 LWKOPT = 1 00206 END IF 00207 WORK( 1 ) = LWKOPT 00208 END IF 00209 * 00210 IF( INFO.NE.0 ) THEN 00211 CALL XERBLA( 'ZUNMBR', -INFO ) 00212 RETURN 00213 ELSE IF( LQUERY ) THEN 00214 RETURN 00215 END IF 00216 * 00217 * Quick return if possible 00218 * 00219 IF( M.EQ.0 .OR. N.EQ.0 ) 00220 $ RETURN 00221 * 00222 IF( APPLYQ ) THEN 00223 * 00224 * Apply Q 00225 * 00226 IF( NQ.GE.K ) THEN 00227 * 00228 * Q was determined by a call to ZGEBRD with nq >= k 00229 * 00230 CALL ZUNMQR( SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, 00231 $ WORK, LWORK, IINFO ) 00232 ELSE IF( NQ.GT.1 ) THEN 00233 * 00234 * Q was determined by a call to ZGEBRD with nq < k 00235 * 00236 IF( LEFT ) THEN 00237 MI = M - 1 00238 NI = N 00239 I1 = 2 00240 I2 = 1 00241 ELSE 00242 MI = M 00243 NI = N - 1 00244 I1 = 1 00245 I2 = 2 00246 END IF 00247 CALL ZUNMQR( SIDE, TRANS, MI, NI, NQ-1, A( 2, 1 ), LDA, TAU, 00248 $ C( I1, I2 ), LDC, WORK, LWORK, IINFO ) 00249 END IF 00250 ELSE 00251 * 00252 * Apply P 00253 * 00254 IF( NOTRAN ) THEN 00255 TRANST = 'C' 00256 ELSE 00257 TRANST = 'N' 00258 END IF 00259 IF( NQ.GT.K ) THEN 00260 * 00261 * P was determined by a call to ZGEBRD with nq > k 00262 * 00263 CALL ZUNMLQ( SIDE, TRANST, M, N, K, A, LDA, TAU, C, LDC, 00264 $ WORK, LWORK, IINFO ) 00265 ELSE IF( NQ.GT.1 ) THEN 00266 * 00267 * P was determined by a call to ZGEBRD with nq <= k 00268 * 00269 IF( LEFT ) THEN 00270 MI = M - 1 00271 NI = N 00272 I1 = 2 00273 I2 = 1 00274 ELSE 00275 MI = M 00276 NI = N - 1 00277 I1 = 1 00278 I2 = 2 00279 END IF 00280 CALL ZUNMLQ( SIDE, TRANST, MI, NI, NQ-1, A( 1, 2 ), LDA, 00281 $ TAU, C( I1, I2 ), LDC, WORK, LWORK, IINFO ) 00282 END IF 00283 END IF 00284 WORK( 1 ) = LWKOPT 00285 RETURN 00286 * 00287 * End of ZUNMBR 00288 * 00289 END
1.7.2