![]() ![]() | ![]() |
Mise à jour de lapack vers la version 3.3.0.
1: SUBROUTINE ZUPMTR( SIDE, UPLO, TRANS, M, N, AP, TAU, C, LDC, WORK, 2: $ 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: CHARACTER SIDE, TRANS, UPLO 11: INTEGER INFO, LDC, M, N 12: * .. 13: * .. Array Arguments .. 14: COMPLEX*16 AP( * ), C( LDC, * ), TAU( * ), WORK( * ) 15: * .. 16: * 17: * Purpose 18: * ======= 19: * 20: * ZUPMTR overwrites the general complex M-by-N matrix C with 21: * 22: * SIDE = 'L' SIDE = 'R' 23: * TRANS = 'N': Q * C C * Q 24: * TRANS = 'C': Q**H * C C * Q**H 25: * 26: * where Q is a complex unitary matrix of order nq, with nq = m if 27: * SIDE = 'L' and nq = n if SIDE = 'R'. Q is defined as the product of 28: * nq-1 elementary reflectors, as returned by ZHPTRD using packed 29: * storage: 30: * 31: * if UPLO = 'U', Q = H(nq-1) . . . H(2) H(1); 32: * 33: * if UPLO = 'L', Q = H(1) H(2) . . . H(nq-1). 34: * 35: * Arguments 36: * ========= 37: * 38: * SIDE (input) CHARACTER*1 39: * = 'L': apply Q or Q**H from the Left; 40: * = 'R': apply Q or Q**H from the Right. 41: * 42: * UPLO (input) CHARACTER*1 43: * = 'U': Upper triangular packed storage used in previous 44: * call to ZHPTRD; 45: * = 'L': Lower triangular packed storage used in previous 46: * call to ZHPTRD. 47: * 48: * TRANS (input) CHARACTER*1 49: * = 'N': No transpose, apply Q; 50: * = 'C': Conjugate transpose, apply Q**H. 51: * 52: * M (input) INTEGER 53: * The number of rows of the matrix C. M >= 0. 54: * 55: * N (input) INTEGER 56: * The number of columns of the matrix C. N >= 0. 57: * 58: * AP (input) COMPLEX*16 array, dimension 59: * (M*(M+1)/2) if SIDE = 'L' 60: * (N*(N+1)/2) if SIDE = 'R' 61: * The vectors which define the elementary reflectors, as 62: * returned by ZHPTRD. AP is modified by the routine but 63: * restored on exit. 64: * 65: * TAU (input) COMPLEX*16 array, dimension (M-1) if SIDE = 'L' 66: * or (N-1) if SIDE = 'R' 67: * TAU(i) must contain the scalar factor of the elementary 68: * reflector H(i), as returned by ZHPTRD. 69: * 70: * C (input/output) COMPLEX*16 array, dimension (LDC,N) 71: * On entry, the M-by-N matrix C. 72: * On exit, C is overwritten by Q*C or Q**H*C or C*Q**H or C*Q. 73: * 74: * LDC (input) INTEGER 75: * The leading dimension of the array C. LDC >= max(1,M). 76: * 77: * WORK (workspace) COMPLEX*16 array, dimension 78: * (N) if SIDE = 'L' 79: * (M) if SIDE = 'R' 80: * 81: * INFO (output) INTEGER 82: * = 0: successful exit 83: * < 0: if INFO = -i, the i-th argument had an illegal value 84: * 85: * ===================================================================== 86: * 87: * .. Parameters .. 88: COMPLEX*16 ONE 89: PARAMETER ( ONE = ( 1.0D+0, 0.0D+0 ) ) 90: * .. 91: * .. Local Scalars .. 92: LOGICAL FORWRD, LEFT, NOTRAN, UPPER 93: INTEGER I, I1, I2, I3, IC, II, JC, MI, NI, NQ 94: COMPLEX*16 AII, TAUI 95: * .. 96: * .. External Functions .. 97: LOGICAL LSAME 98: EXTERNAL LSAME 99: * .. 100: * .. External Subroutines .. 101: EXTERNAL XERBLA, ZLARF 102: * .. 103: * .. Intrinsic Functions .. 104: INTRINSIC DCONJG, MAX 105: * .. 106: * .. Executable Statements .. 107: * 108: * Test the input arguments 109: * 110: INFO = 0 111: LEFT = LSAME( SIDE, 'L' ) 112: NOTRAN = LSAME( TRANS, 'N' ) 113: UPPER = LSAME( UPLO, 'U' ) 114: * 115: * NQ is the order of Q 116: * 117: IF( LEFT ) THEN 118: NQ = M 119: ELSE 120: NQ = N 121: END IF 122: IF( .NOT.LEFT .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN 123: INFO = -1 124: ELSE IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN 125: INFO = -2 126: ELSE IF( .NOT.NOTRAN .AND. .NOT.LSAME( TRANS, 'C' ) ) THEN 127: INFO = -3 128: ELSE IF( M.LT.0 ) THEN 129: INFO = -4 130: ELSE IF( N.LT.0 ) THEN 131: INFO = -5 132: ELSE IF( LDC.LT.MAX( 1, M ) ) THEN 133: INFO = -9 134: END IF 135: IF( INFO.NE.0 ) THEN 136: CALL XERBLA( 'ZUPMTR', -INFO ) 137: RETURN 138: END IF 139: * 140: * Quick return if possible 141: * 142: IF( M.EQ.0 .OR. N.EQ.0 ) 143: $ RETURN 144: * 145: IF( UPPER ) THEN 146: * 147: * Q was determined by a call to ZHPTRD with UPLO = 'U' 148: * 149: FORWRD = ( LEFT .AND. NOTRAN ) .OR. 150: $ ( .NOT.LEFT .AND. .NOT.NOTRAN ) 151: * 152: IF( FORWRD ) THEN 153: I1 = 1 154: I2 = NQ - 1 155: I3 = 1 156: II = 2 157: ELSE 158: I1 = NQ - 1 159: I2 = 1 160: I3 = -1 161: II = NQ*( NQ+1 ) / 2 - 1 162: END IF 163: * 164: IF( LEFT ) THEN 165: NI = N 166: ELSE 167: MI = M 168: END IF 169: * 170: DO 10 I = I1, I2, I3 171: IF( LEFT ) THEN 172: * 173: * H(i) or H(i)' is applied to C(1:i,1:n) 174: * 175: MI = I 176: ELSE 177: * 178: * H(i) or H(i)' is applied to C(1:m,1:i) 179: * 180: NI = I 181: END IF 182: * 183: * Apply H(i) or H(i)' 184: * 185: IF( NOTRAN ) THEN 186: TAUI = TAU( I ) 187: ELSE 188: TAUI = DCONJG( TAU( I ) ) 189: END IF 190: AII = AP( II ) 191: AP( II ) = ONE 192: CALL ZLARF( SIDE, MI, NI, AP( II-I+1 ), 1, TAUI, C, LDC, 193: $ WORK ) 194: AP( II ) = AII 195: * 196: IF( FORWRD ) THEN 197: II = II + I + 2 198: ELSE 199: II = II - I - 1 200: END IF 201: 10 CONTINUE 202: ELSE 203: * 204: * Q was determined by a call to ZHPTRD with UPLO = 'L'. 205: * 206: FORWRD = ( LEFT .AND. .NOT.NOTRAN ) .OR. 207: $ ( .NOT.LEFT .AND. NOTRAN ) 208: * 209: IF( FORWRD ) THEN 210: I1 = 1 211: I2 = NQ - 1 212: I3 = 1 213: II = 2 214: ELSE 215: I1 = NQ - 1 216: I2 = 1 217: I3 = -1 218: II = NQ*( NQ+1 ) / 2 - 1 219: END IF 220: * 221: IF( LEFT ) THEN 222: NI = N 223: JC = 1 224: ELSE 225: MI = M 226: IC = 1 227: END IF 228: * 229: DO 20 I = I1, I2, I3 230: AII = AP( II ) 231: AP( II ) = ONE 232: IF( LEFT ) THEN 233: * 234: * H(i) or H(i)' is applied to C(i+1:m,1:n) 235: * 236: MI = M - I 237: IC = I + 1 238: ELSE 239: * 240: * H(i) or H(i)' is applied to C(1:m,i+1:n) 241: * 242: NI = N - I 243: JC = I + 1 244: END IF 245: * 246: * Apply H(i) or H(i)' 247: * 248: IF( NOTRAN ) THEN 249: TAUI = TAU( I ) 250: ELSE 251: TAUI = DCONJG( TAU( I ) ) 252: END IF 253: CALL ZLARF( SIDE, MI, NI, AP( II ), 1, TAUI, C( IC, JC ), 254: $ LDC, WORK ) 255: AP( II ) = AII 256: * 257: IF( FORWRD ) THEN 258: II = II + NQ - I + 1 259: ELSE 260: II = II - NQ + I - 2 261: END IF 262: 20 CONTINUE 263: END IF 264: RETURN 265: * 266: * End of ZUPMTR 267: * 268: END