KTH framework for Nek5000 toolboxes; testing version  0.0.1
dtrsm.f
Go to the documentation of this file.
1  SUBROUTINE dtrsm ( SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA,
2  $ B, LDB )
3 * .. Scalar Arguments ..
4  CHARACTER*1 SIDE, UPLO, TRANSA, DIAG
5  INTEGER M, N, LDA, LDB
6  DOUBLE PRECISION ALPHA
7 * .. Array Arguments ..
8  DOUBLE PRECISION A( LDA, * ), B( LDB, * )
9 * ..
10 *
11 * Purpose
12 * =======
13 *
14 * DTRSM 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'.
22 *
23 * The matrix X is overwritten on B.
24 *
25 * Parameters
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 ) = 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 - DOUBLE PRECISION.
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 - DOUBLE PRECISION 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 - DOUBLE PRECISION 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 *
120 * Level 3 Blas routine.
121 *
122 *
123 * -- Written on 8-February-1989.
124 * Jack Dongarra, Argonne National Laboratory.
125 * Iain Duff, AERE Harwell.
126 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
127 * Sven Hammarling, Numerical Algorithms Group Ltd.
128 *
129 *
130 * .. External Functions ..
131  LOGICAL LSAME
132  EXTERNAL lsame
133 * .. External Subroutines ..
134  EXTERNAL xerbla
135 * .. Intrinsic Functions ..
136  INTRINSIC max
137 * .. Local Scalars ..
138  LOGICAL LSIDE, NOUNIT, UPPER
139  INTEGER I, INFO, J, K, NROWA
140  DOUBLE PRECISION TEMP
141 * .. Parameters ..
142  DOUBLE PRECISION ONE , ZERO
143  parameter( one = 1.0d+0, zero = 0.0d+0 )
144 * ..
145 * .. Executable Statements ..
146 *
147 * Test the input parameters.
148 *
149  lside = lsame( side , 'L' )
150  IF( lside )THEN
151  nrowa = m
152  ELSE
153  nrowa = n
154  END IF
155  nounit = lsame( diag , 'N' )
156  upper = lsame( uplo , 'U' )
157 *
158  info = 0
159  IF( ( .NOT.lside ).AND.
160  $ ( .NOT.lsame( side , 'R' ) ) )THEN
161  info = 1
162  ELSE IF( ( .NOT.upper ).AND.
163  $ ( .NOT.lsame( uplo , 'L' ) ) )THEN
164  info = 2
165  ELSE IF( ( .NOT.lsame( transa, 'N' ) ).AND.
166  $ ( .NOT.lsame( transa, 'T' ) ).AND.
167  $ ( .NOT.lsame( transa, 'C' ) ) )THEN
168  info = 3
169  ELSE IF( ( .NOT.lsame( diag , 'U' ) ).AND.
170  $ ( .NOT.lsame( diag , 'N' ) ) )THEN
171  info = 4
172  ELSE IF( m .LT.0 )THEN
173  info = 5
174  ELSE IF( n .LT.0 )THEN
175  info = 6
176  ELSE IF( lda.LT.max( 1, nrowa ) )THEN
177  info = 9
178  ELSE IF( ldb.LT.max( 1, m ) )THEN
179  info = 11
180  END IF
181  IF( info.NE.0 )THEN
182  CALL xerbla( 'DTRSM ', info )
183  RETURN
184  END IF
185 *
186 * Quick return if possible.
187 *
188  IF( n.EQ.0 )
189  $ RETURN
190 *
191 * And when alpha.eq.zero.
192 *
193  IF( alpha.EQ.zero )THEN
194  DO 20, j = 1, n
195  DO 10, i = 1, m
196  b( i, j ) = zero
197  10 CONTINUE
198  20 CONTINUE
199  RETURN
200  END IF
201 *
202 * Start the operations.
203 *
204  IF( lside )THEN
205  IF( lsame( transa, 'N' ) )THEN
206 *
207 * Form B := alpha*inv( A )*B.
208 *
209  IF( upper )THEN
210  DO 60, j = 1, n
211  IF( alpha.NE.one )THEN
212  DO 30, i = 1, m
213  b( i, j ) = alpha*b( i, j )
214  30 CONTINUE
215  END IF
216  DO 50, k = m, 1, -1
217  IF( b( k, j ).NE.zero )THEN
218  IF( nounit )
219  $ b( k, j ) = b( k, j )/a( k, k )
220  DO 40, i = 1, k - 1
221  b( i, j ) = b( i, j ) - b( k, j )*a( i, k )
222  40 CONTINUE
223  END IF
224  50 CONTINUE
225  60 CONTINUE
226  ELSE
227  DO 100, j = 1, n
228  IF( alpha.NE.one )THEN
229  DO 70, i = 1, m
230  b( i, j ) = alpha*b( i, j )
231  70 CONTINUE
232  END IF
233  DO 90 k = 1, m
234  IF( b( k, j ).NE.zero )THEN
235  IF( nounit )
236  $ b( k, j ) = b( k, j )/a( k, k )
237  DO 80, i = k + 1, m
238  b( i, j ) = b( i, j ) - b( k, j )*a( i, k )
239  80 CONTINUE
240  END IF
241  90 CONTINUE
242  100 CONTINUE
243  END IF
244  ELSE
245 *
246 * Form B := alpha*inv( A' )*B.
247 *
248  IF( upper )THEN
249  DO 130, j = 1, n
250  DO 120, i = 1, m
251  temp = alpha*b( i, j )
252  DO 110, k = 1, i - 1
253  temp = temp - a( k, i )*b( k, j )
254  110 CONTINUE
255  IF( nounit )
256  $ temp = temp/a( i, i )
257  b( i, j ) = temp
258  120 CONTINUE
259  130 CONTINUE
260  ELSE
261  DO 160, j = 1, n
262  DO 150, i = m, 1, -1
263  temp = alpha*b( i, j )
264  DO 140, k = i + 1, m
265  temp = temp - a( k, i )*b( k, j )
266  140 CONTINUE
267  IF( nounit )
268  $ temp = temp/a( i, i )
269  b( i, j ) = temp
270  150 CONTINUE
271  160 CONTINUE
272  END IF
273  END IF
274  ELSE
275  IF( lsame( transa, 'N' ) )THEN
276 *
277 * Form B := alpha*B*inv( A ).
278 *
279  IF( upper )THEN
280  DO 210, j = 1, n
281  IF( alpha.NE.one )THEN
282  DO 170, i = 1, m
283  b( i, j ) = alpha*b( i, j )
284  170 CONTINUE
285  END IF
286  DO 190, k = 1, j - 1
287  IF( a( k, j ).NE.zero )THEN
288  DO 180, i = 1, m
289  b( i, j ) = b( i, j ) - a( k, j )*b( i, k )
290  180 CONTINUE
291  END IF
292  190 CONTINUE
293  IF( nounit )THEN
294  temp = one/a( j, j )
295  DO 200, i = 1, m
296  b( i, j ) = temp*b( i, j )
297  200 CONTINUE
298  END IF
299  210 CONTINUE
300  ELSE
301  DO 260, j = n, 1, -1
302  IF( alpha.NE.one )THEN
303  DO 220, i = 1, m
304  b( i, j ) = alpha*b( i, j )
305  220 CONTINUE
306  END IF
307  DO 240, k = j + 1, n
308  IF( a( k, j ).NE.zero )THEN
309  DO 230, i = 1, m
310  b( i, j ) = b( i, j ) - a( k, j )*b( i, k )
311  230 CONTINUE
312  END IF
313  240 CONTINUE
314  IF( nounit )THEN
315  temp = one/a( j, j )
316  DO 250, i = 1, m
317  b( i, j ) = temp*b( i, j )
318  250 CONTINUE
319  END IF
320  260 CONTINUE
321  END IF
322  ELSE
323 *
324 * Form B := alpha*B*inv( A' ).
325 *
326  IF( upper )THEN
327  DO 310, k = n, 1, -1
328  IF( nounit )THEN
329  temp = one/a( k, k )
330  DO 270, i = 1, m
331  b( i, k ) = temp*b( i, k )
332  270 CONTINUE
333  END IF
334  DO 290, j = 1, k - 1
335  IF( a( j, k ).NE.zero )THEN
336  temp = a( j, k )
337  DO 280, i = 1, m
338  b( i, j ) = b( i, j ) - temp*b( i, k )
339  280 CONTINUE
340  END IF
341  290 CONTINUE
342  IF( alpha.NE.one )THEN
343  DO 300, i = 1, m
344  b( i, k ) = alpha*b( i, k )
345  300 CONTINUE
346  END IF
347  310 CONTINUE
348  ELSE
349  DO 360, k = 1, n
350  IF( nounit )THEN
351  temp = one/a( k, k )
352  DO 320, i = 1, m
353  b( i, k ) = temp*b( i, k )
354  320 CONTINUE
355  END IF
356  DO 340, j = k + 1, n
357  IF( a( j, k ).NE.zero )THEN
358  temp = a( j, k )
359  DO 330, i = 1, m
360  b( i, j ) = b( i, j ) - temp*b( i, k )
361  330 CONTINUE
362  END IF
363  340 CONTINUE
364  IF( alpha.NE.one )THEN
365  DO 350, i = 1, m
366  b( i, k ) = alpha*b( i, k )
367  350 CONTINUE
368  END IF
369  360 CONTINUE
370  END IF
371  END IF
372  END IF
373 *
374  RETURN
375 *
376 * End of DTRSM .
377 *
378  END
dtrsm
subroutine dtrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
Definition: dtrsm.f:3
xerbla
subroutine xerbla(SRNAME, INFO)
Definition: xerbla.f:2