![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to ztrsm.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / blas|
Return to ztrsm.f CVS log | Up to [local] / rpl / lapack / blas |
1.1 bertrand 1: SUBROUTINE ZTRSM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
2: * .. Scalar Arguments ..
3: DOUBLE COMPLEX ALPHA
4: INTEGER LDA,LDB,M,N
5: CHARACTER DIAG,SIDE,TRANSA,UPLO
6: * ..
7: * .. Array Arguments ..
8: DOUBLE COMPLEX A(LDA,*),B(LDB,*)
9: * ..
10: *
11: * Purpose
12: * =======
13: *
14: * ZTRSM solves one of the matrix equations
15: *
16: * op( A )*X = alpha*B, or X*op( A ) = alpha*B,
17: *
18: * where alpha is a scalar, X and B are m by n matrices, A is a unit, or
19: * non-unit, upper or lower triangular matrix and op( A ) is one of
20: *
21: * op( A ) = A or op( A ) = A' or op( A ) = conjg( A' ).
22: *
23: * The matrix X is overwritten on B.
24: *
25: * Arguments
26: * ==========
27: *
28: * SIDE - CHARACTER*1.
29: * On entry, SIDE specifies whether op( A ) appears on the left
30: * or right of X as follows:
31: *
32: * SIDE = 'L' or 'l' op( A )*X = alpha*B.
33: *
34: * SIDE = 'R' or 'r' X*op( A ) = alpha*B.
35: *
36: * Unchanged on exit.
37: *
38: * UPLO - CHARACTER*1.
39: * On entry, UPLO specifies whether the matrix A is an upper or
40: * lower triangular matrix as follows:
41: *
42: * UPLO = 'U' or 'u' A is an upper triangular matrix.
43: *
44: * UPLO = 'L' or 'l' A is a lower triangular matrix.
45: *
46: * Unchanged on exit.
47: *
48: * TRANSA - CHARACTER*1.
49: * On entry, TRANSA specifies the form of op( A ) to be used in
50: * the matrix multiplication as follows:
51: *
52: * TRANSA = 'N' or 'n' op( A ) = A.
53: *
54: * TRANSA = 'T' or 't' op( A ) = A'.
55: *
56: * TRANSA = 'C' or 'c' op( A ) = conjg( A' ).
57: *
58: * Unchanged on exit.
59: *
60: * DIAG - CHARACTER*1.
61: * On entry, DIAG specifies whether or not A is unit triangular
62: * as follows:
63: *
64: * DIAG = 'U' or 'u' A is assumed to be unit triangular.
65: *
66: * DIAG = 'N' or 'n' A is not assumed to be unit
67: * triangular.
68: *
69: * Unchanged on exit.
70: *
71: * M - INTEGER.
72: * On entry, M specifies the number of rows of B. M must be at
73: * least zero.
74: * Unchanged on exit.
75: *
76: * N - INTEGER.
77: * On entry, N specifies the number of columns of B. N must be
78: * at least zero.
79: * Unchanged on exit.
80: *
81: * ALPHA - COMPLEX*16 .
82: * On entry, ALPHA specifies the scalar alpha. When alpha is
83: * zero then A is not referenced and B need not be set before
84: * entry.
85: * Unchanged on exit.
86: *
87: * A - COMPLEX*16 array of DIMENSION ( LDA, k ), where k is m
88: * when SIDE = 'L' or 'l' and is n when SIDE = 'R' or 'r'.
89: * Before entry with UPLO = 'U' or 'u', the leading k by k
90: * upper triangular part of the array A must contain the upper
91: * triangular matrix and the strictly lower triangular part of
92: * A is not referenced.
93: * Before entry with UPLO = 'L' or 'l', the leading k by k
94: * lower triangular part of the array A must contain the lower
95: * triangular matrix and the strictly upper triangular part of
96: * A is not referenced.
97: * Note that when DIAG = 'U' or 'u', the diagonal elements of
98: * A are not referenced either, but are assumed to be unity.
99: * Unchanged on exit.
100: *
101: * LDA - INTEGER.
102: * On entry, LDA specifies the first dimension of A as declared
103: * in the calling (sub) program. When SIDE = 'L' or 'l' then
104: * LDA must be at least max( 1, m ), when SIDE = 'R' or 'r'
105: * then LDA must be at least max( 1, n ).
106: * Unchanged on exit.
107: *
108: * B - COMPLEX*16 array of DIMENSION ( LDB, n ).
109: * Before entry, the leading m by n part of the array B must
110: * contain the right-hand side matrix B, and on exit is
111: * overwritten by the solution matrix X.
112: *
113: * LDB - INTEGER.
114: * On entry, LDB specifies the first dimension of B as declared
115: * in the calling (sub) program. LDB must be at least
116: * max( 1, m ).
117: * Unchanged on exit.
118: *
119: * Further Details
120: * ===============
121: *
122: * Level 3 Blas routine.
123: *
124: * -- Written on 8-February-1989.
125: * Jack Dongarra, Argonne National Laboratory.
126: * Iain Duff, AERE Harwell.
127: * Jeremy Du Croz, Numerical Algorithms Group Ltd.
128: * Sven Hammarling, Numerical Algorithms Group Ltd.
129: *
130: * =====================================================================
131: *
132: * .. External Functions ..
133: LOGICAL LSAME
134: EXTERNAL LSAME
135: * ..
136: * .. External Subroutines ..
137: EXTERNAL XERBLA
138: * ..
139: * .. Intrinsic Functions ..
140: INTRINSIC DCONJG,MAX
141: * ..
142: * .. Local Scalars ..
143: DOUBLE COMPLEX TEMP
144: INTEGER I,INFO,J,K,NROWA
145: LOGICAL LSIDE,NOCONJ,NOUNIT,UPPER
146: * ..
147: * .. Parameters ..
148: DOUBLE COMPLEX ONE
149: PARAMETER (ONE= (1.0D+0,0.0D+0))
150: DOUBLE COMPLEX ZERO
151: PARAMETER (ZERO= (0.0D+0,0.0D+0))
152: * ..
153: *
154: * Test the input parameters.
155: *
156: LSIDE = LSAME(SIDE,'L')
157: IF (LSIDE) THEN
158: NROWA = M
159: ELSE
160: NROWA = N
161: END IF
162: NOCONJ = LSAME(TRANSA,'T')
163: NOUNIT = LSAME(DIAG,'N')
164: UPPER = LSAME(UPLO,'U')
165: *
166: INFO = 0
167: IF ((.NOT.LSIDE) .AND. (.NOT.LSAME(SIDE,'R'))) THEN
168: INFO = 1
169: ELSE IF ((.NOT.UPPER) .AND. (.NOT.LSAME(UPLO,'L'))) THEN
170: INFO = 2
171: ELSE IF ((.NOT.LSAME(TRANSA,'N')) .AND.
172: + (.NOT.LSAME(TRANSA,'T')) .AND.
173: + (.NOT.LSAME(TRANSA,'C'))) THEN
174: INFO = 3
175: ELSE IF ((.NOT.LSAME(DIAG,'U')) .AND. (.NOT.LSAME(DIAG,'N'))) THEN
176: INFO = 4
177: ELSE IF (M.LT.0) THEN
178: INFO = 5
179: ELSE IF (N.LT.0) THEN
180: INFO = 6
181: ELSE IF (LDA.LT.MAX(1,NROWA)) THEN
182: INFO = 9
183: ELSE IF (LDB.LT.MAX(1,M)) THEN
184: INFO = 11
185: END IF
186: IF (INFO.NE.0) THEN
187: CALL XERBLA('ZTRSM ',INFO)
188: RETURN
189: END IF
190: *
191: * Quick return if possible.
192: *
193: IF (M.EQ.0 .OR. N.EQ.0) RETURN
194: *
195: * And when alpha.eq.zero.
196: *
197: IF (ALPHA.EQ.ZERO) THEN
198: DO 20 J = 1,N
199: DO 10 I = 1,M
200: B(I,J) = ZERO
201: 10 CONTINUE
202: 20 CONTINUE
203: RETURN
204: END IF
205: *
206: * Start the operations.
207: *
208: IF (LSIDE) THEN
209: IF (LSAME(TRANSA,'N')) THEN
210: *
211: * Form B := alpha*inv( A )*B.
212: *
213: IF (UPPER) THEN
214: DO 60 J = 1,N
215: IF (ALPHA.NE.ONE) THEN
216: DO 30 I = 1,M
217: B(I,J) = ALPHA*B(I,J)
218: 30 CONTINUE
219: END IF
220: DO 50 K = M,1,-1
221: IF (B(K,J).NE.ZERO) THEN
222: IF (NOUNIT) B(K,J) = B(K,J)/A(K,K)
223: DO 40 I = 1,K - 1
224: B(I,J) = B(I,J) - B(K,J)*A(I,K)
225: 40 CONTINUE
226: END IF
227: 50 CONTINUE
228: 60 CONTINUE
229: ELSE
230: DO 100 J = 1,N
231: IF (ALPHA.NE.ONE) THEN
232: DO 70 I = 1,M
233: B(I,J) = ALPHA*B(I,J)
234: 70 CONTINUE
235: END IF
236: DO 90 K = 1,M
237: IF (B(K,J).NE.ZERO) THEN
238: IF (NOUNIT) B(K,J) = B(K,J)/A(K,K)
239: DO 80 I = K + 1,M
240: B(I,J) = B(I,J) - B(K,J)*A(I,K)
241: 80 CONTINUE
242: END IF
243: 90 CONTINUE
244: 100 CONTINUE
245: END IF
246: ELSE
247: *
248: * Form B := alpha*inv( A' )*B
249: * or B := alpha*inv( conjg( A' ) )*B.
250: *
251: IF (UPPER) THEN
252: DO 140 J = 1,N
253: DO 130 I = 1,M
254: TEMP = ALPHA*B(I,J)
255: IF (NOCONJ) THEN
256: DO 110 K = 1,I - 1
257: TEMP = TEMP - A(K,I)*B(K,J)
258: 110 CONTINUE
259: IF (NOUNIT) TEMP = TEMP/A(I,I)
260: ELSE
261: DO 120 K = 1,I - 1
262: TEMP = TEMP - DCONJG(A(K,I))*B(K,J)
263: 120 CONTINUE
264: IF (NOUNIT) TEMP = TEMP/DCONJG(A(I,I))
265: END IF
266: B(I,J) = TEMP
267: 130 CONTINUE
268: 140 CONTINUE
269: ELSE
270: DO 180 J = 1,N
271: DO 170 I = M,1,-1
272: TEMP = ALPHA*B(I,J)
273: IF (NOCONJ) THEN
274: DO 150 K = I + 1,M
275: TEMP = TEMP - A(K,I)*B(K,J)
276: 150 CONTINUE
277: IF (NOUNIT) TEMP = TEMP/A(I,I)
278: ELSE
279: DO 160 K = I + 1,M
280: TEMP = TEMP - DCONJG(A(K,I))*B(K,J)
281: 160 CONTINUE
282: IF (NOUNIT) TEMP = TEMP/DCONJG(A(I,I))
283: END IF
284: B(I,J) = TEMP
285: 170 CONTINUE
286: 180 CONTINUE
287: END IF
288: END IF
289: ELSE
290: IF (LSAME(TRANSA,'N')) THEN
291: *
292: * Form B := alpha*B*inv( A ).
293: *
294: IF (UPPER) THEN
295: DO 230 J = 1,N
296: IF (ALPHA.NE.ONE) THEN
297: DO 190 I = 1,M
298: B(I,J) = ALPHA*B(I,J)
299: 190 CONTINUE
300: END IF
301: DO 210 K = 1,J - 1
302: IF (A(K,J).NE.ZERO) THEN
303: DO 200 I = 1,M
304: B(I,J) = B(I,J) - A(K,J)*B(I,K)
305: 200 CONTINUE
306: END IF
307: 210 CONTINUE
308: IF (NOUNIT) THEN
309: TEMP = ONE/A(J,J)
310: DO 220 I = 1,M
311: B(I,J) = TEMP*B(I,J)
312: 220 CONTINUE
313: END IF
314: 230 CONTINUE
315: ELSE
316: DO 280 J = N,1,-1
317: IF (ALPHA.NE.ONE) THEN
318: DO 240 I = 1,M
319: B(I,J) = ALPHA*B(I,J)
320: 240 CONTINUE
321: END IF
322: DO 260 K = J + 1,N
323: IF (A(K,J).NE.ZERO) THEN
324: DO 250 I = 1,M
325: B(I,J) = B(I,J) - A(K,J)*B(I,K)
326: 250 CONTINUE
327: END IF
328: 260 CONTINUE
329: IF (NOUNIT) THEN
330: TEMP = ONE/A(J,J)
331: DO 270 I = 1,M
332: B(I,J) = TEMP*B(I,J)
333: 270 CONTINUE
334: END IF
335: 280 CONTINUE
336: END IF
337: ELSE
338: *
339: * Form B := alpha*B*inv( A' )
340: * or B := alpha*B*inv( conjg( A' ) ).
341: *
342: IF (UPPER) THEN
343: DO 330 K = N,1,-1
344: IF (NOUNIT) THEN
345: IF (NOCONJ) THEN
346: TEMP = ONE/A(K,K)
347: ELSE
348: TEMP = ONE/DCONJG(A(K,K))
349: END IF
350: DO 290 I = 1,M
351: B(I,K) = TEMP*B(I,K)
352: 290 CONTINUE
353: END IF
354: DO 310 J = 1,K - 1
355: IF (A(J,K).NE.ZERO) THEN
356: IF (NOCONJ) THEN
357: TEMP = A(J,K)
358: ELSE
359: TEMP = DCONJG(A(J,K))
360: END IF
361: DO 300 I = 1,M
362: B(I,J) = B(I,J) - TEMP*B(I,K)
363: 300 CONTINUE
364: END IF
365: 310 CONTINUE
366: IF (ALPHA.NE.ONE) THEN
367: DO 320 I = 1,M
368: B(I,K) = ALPHA*B(I,K)
369: 320 CONTINUE
370: END IF
371: 330 CONTINUE
372: ELSE
373: DO 380 K = 1,N
374: IF (NOUNIT) THEN
375: IF (NOCONJ) THEN
376: TEMP = ONE/A(K,K)
377: ELSE
378: TEMP = ONE/DCONJG(A(K,K))
379: END IF
380: DO 340 I = 1,M
381: B(I,K) = TEMP*B(I,K)
382: 340 CONTINUE
383: END IF
384: DO 360 J = K + 1,N
385: IF (A(J,K).NE.ZERO) THEN
386: IF (NOCONJ) THEN
387: TEMP = A(J,K)
388: ELSE
389: TEMP = DCONJG(A(J,K))
390: END IF
391: DO 350 I = 1,M
392: B(I,J) = B(I,J) - TEMP*B(I,K)
393: 350 CONTINUE
394: END IF
395: 360 CONTINUE
396: IF (ALPHA.NE.ONE) THEN
397: DO 370 I = 1,M
398: B(I,K) = ALPHA*B(I,K)
399: 370 CONTINUE
400: END IF
401: 380 CONTINUE
402: END IF
403: END IF
404: END IF
405: *
406: RETURN
407: *
408: * End of ZTRSM .
409: *
410: END