LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
dtrsm.f
Go to the documentation of this file.
1 *> \brief \b DTRSM
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * SUBROUTINE DTRSM(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
12 *
13 * .. Scalar Arguments ..
14 * DOUBLE PRECISION ALPHA
15 * INTEGER LDA,LDB,M,N
16 * CHARACTER DIAG,SIDE,TRANSA,UPLO
17 * ..
18 * .. Array Arguments ..
19 * DOUBLE PRECISION A(LDA,*),B(LDB,*)
20 * ..
21 *
22 *
23 *> \par Purpose:
24 * =============
25 *>
26 *> \verbatim
27 *>
28 *> DTRSM solves one of the matrix equations
29 *>
30 *> op( A )*X = alpha*B, or X*op( A ) = alpha*B,
31 *>
32 *> where alpha is a scalar, X and B are m by n matrices, A is a unit, or
33 *> non-unit, upper or lower triangular matrix and op( A ) is one of
34 *>
35 *> op( A ) = A or op( A ) = A**T.
36 *>
37 *> The matrix X is overwritten on B.
38 *> \endverbatim
39 *
40 * Arguments:
41 * ==========
42 *
43 *> \param[in] SIDE
44 *> \verbatim
45 *> SIDE is CHARACTER*1
46 *> On entry, SIDE specifies whether op( A ) appears on the left
47 *> or right of X as follows:
48 *>
49 *> SIDE = 'L' or 'l' op( A )*X = alpha*B.
50 *>
51 *> SIDE = 'R' or 'r' X*op( A ) = alpha*B.
52 *> \endverbatim
53 *>
54 *> \param[in] UPLO
55 *> \verbatim
56 *> UPLO is CHARACTER*1
57 *> On entry, UPLO specifies whether the matrix A is an upper or
58 *> lower triangular matrix as follows:
59 *>
60 *> UPLO = 'U' or 'u' A is an upper triangular matrix.
61 *>
62 *> UPLO = 'L' or 'l' A is a lower triangular matrix.
63 *> \endverbatim
64 *>
65 *> \param[in] TRANSA
66 *> \verbatim
67 *> TRANSA is CHARACTER*1
68 *> On entry, TRANSA specifies the form of op( A ) to be used in
69 *> the matrix multiplication as follows:
70 *>
71 *> TRANSA = 'N' or 'n' op( A ) = A.
72 *>
73 *> TRANSA = 'T' or 't' op( A ) = A**T.
74 *>
75 *> TRANSA = 'C' or 'c' op( A ) = A**T.
76 *> \endverbatim
77 *>
78 *> \param[in] DIAG
79 *> \verbatim
80 *> DIAG is CHARACTER*1
81 *> On entry, DIAG specifies whether or not A is unit triangular
82 *> as follows:
83 *>
84 *> DIAG = 'U' or 'u' A is assumed to be unit triangular.
85 *>
86 *> DIAG = 'N' or 'n' A is not assumed to be unit
87 *> triangular.
88 *> \endverbatim
89 *>
90 *> \param[in] M
91 *> \verbatim
92 *> M is INTEGER
93 *> On entry, M specifies the number of rows of B. M must be at
94 *> least zero.
95 *> \endverbatim
96 *>
97 *> \param[in] N
98 *> \verbatim
99 *> N is INTEGER
100 *> On entry, N specifies the number of columns of B. N must be
101 *> at least zero.
102 *> \endverbatim
103 *>
104 *> \param[in] ALPHA
105 *> \verbatim
106 *> ALPHA is DOUBLE PRECISION.
107 *> On entry, ALPHA specifies the scalar alpha. When alpha is
108 *> zero then A is not referenced and B need not be set before
109 *> entry.
110 *> \endverbatim
111 *>
112 *> \param[in] A
113 *> \verbatim
114 *> A is DOUBLE PRECISION array of DIMENSION ( LDA, k ),
115 *> where k is m when SIDE = 'L' or 'l'
116 *> and k is n when SIDE = 'R' or 'r'.
117 *> Before entry with UPLO = 'U' or 'u', the leading k by k
118 *> upper triangular part of the array A must contain the upper
119 *> triangular matrix and the strictly lower triangular part of
120 *> A is not referenced.
121 *> Before entry with UPLO = 'L' or 'l', the leading k by k
122 *> lower triangular part of the array A must contain the lower
123 *> triangular matrix and the strictly upper triangular part of
124 *> A is not referenced.
125 *> Note that when DIAG = 'U' or 'u', the diagonal elements of
126 *> A are not referenced either, but are assumed to be unity.
127 *> \endverbatim
128 *>
129 *> \param[in] LDA
130 *> \verbatim
131 *> LDA is INTEGER
132 *> On entry, LDA specifies the first dimension of A as declared
133 *> in the calling (sub) program. When SIDE = 'L' or 'l' then
134 *> LDA must be at least max( 1, m ), when SIDE = 'R' or 'r'
135 *> then LDA must be at least max( 1, n ).
136 *> \endverbatim
137 *>
138 *> \param[in,out] B
139 *> \verbatim
140 *> B is DOUBLE PRECISION array of DIMENSION ( LDB, n ).
141 *> Before entry, the leading m by n part of the array B must
142 *> contain the right-hand side matrix B, and on exit is
143 *> overwritten by the solution matrix X.
144 *> \endverbatim
145 *>
146 *> \param[in] LDB
147 *> \verbatim
148 *> LDB is INTEGER
149 *> On entry, LDB specifies the first dimension of B as declared
150 *> in the calling (sub) program. LDB must be at least
151 *> max( 1, m ).
152 *> \endverbatim
153 *
154 * Authors:
155 * ========
156 *
157 *> \author Univ. of Tennessee
158 *> \author Univ. of California Berkeley
159 *> \author Univ. of Colorado Denver
160 *> \author NAG Ltd.
161 *
162 *> \date November 2011
163 *
164 *> \ingroup double_blas_level3
165 *
166 *> \par Further Details:
167 * =====================
168 *>
169 *> \verbatim
170 *>
171 *> Level 3 Blas routine.
172 *>
173 *>
174 *> -- Written on 8-February-1989.
175 *> Jack Dongarra, Argonne National Laboratory.
176 *> Iain Duff, AERE Harwell.
177 *> Jeremy Du Croz, Numerical Algorithms Group Ltd.
178 *> Sven Hammarling, Numerical Algorithms Group Ltd.
179 *> \endverbatim
180 *>
181 * =====================================================================
182  SUBROUTINE dtrsm(SIDE,UPLO,TRANSA,DIAG,M,N,ALPHA,A,LDA,B,LDB)
183 *
184 * -- Reference BLAS level3 routine (version 3.4.0) --
185 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
186 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
187 * November 2011
188 *
189 * .. Scalar Arguments ..
190  DOUBLE PRECISION alpha
191  INTEGER lda,ldb,m,n
192  CHARACTER diag,side,transa,uplo
193 * ..
194 * .. Array Arguments ..
195  DOUBLE PRECISION a(lda,*),b(ldb,*)
196 * ..
197 *
198 * =====================================================================
199 *
200 * .. External Functions ..
201  LOGICAL lsame
202  EXTERNAL lsame
203 * ..
204 * .. External Subroutines ..
205  EXTERNAL xerbla
206 * ..
207 * .. Intrinsic Functions ..
208  INTRINSIC max
209 * ..
210 * .. Local Scalars ..
211  DOUBLE PRECISION temp
212  INTEGER i,info,j,k,nrowa
213  LOGICAL lside,nounit,upper
214 * ..
215 * .. Parameters ..
216  DOUBLE PRECISION one,zero
217  parameter(one=1.0d+0,zero=0.0d+0)
218 * ..
219 *
220 * Test the input parameters.
221 *
222  lside = lsame(side,'L')
223  IF (lside) THEN
224  nrowa = m
225  ELSE
226  nrowa = n
227  END IF
228  nounit = lsame(diag,'N')
229  upper = lsame(uplo,'U')
230 *
231  info = 0
232  IF ((.NOT.lside) .AND. (.NOT.lsame(side,'R'))) THEN
233  info = 1
234  ELSE IF ((.NOT.upper) .AND. (.NOT.lsame(uplo,'L'))) THEN
235  info = 2
236  ELSE IF ((.NOT.lsame(transa,'N')) .AND.
237  + (.NOT.lsame(transa,'T')) .AND.
238  + (.NOT.lsame(transa,'C'))) THEN
239  info = 3
240  ELSE IF ((.NOT.lsame(diag,'U')) .AND. (.NOT.lsame(diag,'N'))) THEN
241  info = 4
242  ELSE IF (m.LT.0) THEN
243  info = 5
244  ELSE IF (n.LT.0) THEN
245  info = 6
246  ELSE IF (lda.LT.max(1,nrowa)) THEN
247  info = 9
248  ELSE IF (ldb.LT.max(1,m)) THEN
249  info = 11
250  END IF
251  IF (info.NE.0) THEN
252  CALL xerbla('DTRSM ',info)
253  return
254  END IF
255 *
256 * Quick return if possible.
257 *
258  IF (m.EQ.0 .OR. n.EQ.0) return
259 *
260 * And when alpha.eq.zero.
261 *
262  IF (alpha.EQ.zero) THEN
263  DO 20 j = 1,n
264  DO 10 i = 1,m
265  b(i,j) = zero
266  10 continue
267  20 continue
268  return
269  END IF
270 *
271 * Start the operations.
272 *
273  IF (lside) THEN
274  IF (lsame(transa,'N')) THEN
275 *
276 * Form B := alpha*inv( A )*B.
277 *
278  IF (upper) THEN
279  DO 60 j = 1,n
280  IF (alpha.NE.one) THEN
281  DO 30 i = 1,m
282  b(i,j) = alpha*b(i,j)
283  30 continue
284  END IF
285  DO 50 k = m,1,-1
286  IF (b(k,j).NE.zero) THEN
287  IF (nounit) b(k,j) = b(k,j)/a(k,k)
288  DO 40 i = 1,k - 1
289  b(i,j) = b(i,j) - b(k,j)*a(i,k)
290  40 continue
291  END IF
292  50 continue
293  60 continue
294  ELSE
295  DO 100 j = 1,n
296  IF (alpha.NE.one) THEN
297  DO 70 i = 1,m
298  b(i,j) = alpha*b(i,j)
299  70 continue
300  END IF
301  DO 90 k = 1,m
302  IF (b(k,j).NE.zero) THEN
303  IF (nounit) b(k,j) = b(k,j)/a(k,k)
304  DO 80 i = k + 1,m
305  b(i,j) = b(i,j) - b(k,j)*a(i,k)
306  80 continue
307  END IF
308  90 continue
309  100 continue
310  END IF
311  ELSE
312 *
313 * Form B := alpha*inv( A**T )*B.
314 *
315  IF (upper) THEN
316  DO 130 j = 1,n
317  DO 120 i = 1,m
318  temp = alpha*b(i,j)
319  DO 110 k = 1,i - 1
320  temp = temp - a(k,i)*b(k,j)
321  110 continue
322  IF (nounit) temp = temp/a(i,i)
323  b(i,j) = temp
324  120 continue
325  130 continue
326  ELSE
327  DO 160 j = 1,n
328  DO 150 i = m,1,-1
329  temp = alpha*b(i,j)
330  DO 140 k = i + 1,m
331  temp = temp - a(k,i)*b(k,j)
332  140 continue
333  IF (nounit) temp = temp/a(i,i)
334  b(i,j) = temp
335  150 continue
336  160 continue
337  END IF
338  END IF
339  ELSE
340  IF (lsame(transa,'N')) THEN
341 *
342 * Form B := alpha*B*inv( A ).
343 *
344  IF (upper) THEN
345  DO 210 j = 1,n
346  IF (alpha.NE.one) THEN
347  DO 170 i = 1,m
348  b(i,j) = alpha*b(i,j)
349  170 continue
350  END IF
351  DO 190 k = 1,j - 1
352  IF (a(k,j).NE.zero) THEN
353  DO 180 i = 1,m
354  b(i,j) = b(i,j) - a(k,j)*b(i,k)
355  180 continue
356  END IF
357  190 continue
358  IF (nounit) THEN
359  temp = one/a(j,j)
360  DO 200 i = 1,m
361  b(i,j) = temp*b(i,j)
362  200 continue
363  END IF
364  210 continue
365  ELSE
366  DO 260 j = n,1,-1
367  IF (alpha.NE.one) THEN
368  DO 220 i = 1,m
369  b(i,j) = alpha*b(i,j)
370  220 continue
371  END IF
372  DO 240 k = j + 1,n
373  IF (a(k,j).NE.zero) THEN
374  DO 230 i = 1,m
375  b(i,j) = b(i,j) - a(k,j)*b(i,k)
376  230 continue
377  END IF
378  240 continue
379  IF (nounit) THEN
380  temp = one/a(j,j)
381  DO 250 i = 1,m
382  b(i,j) = temp*b(i,j)
383  250 continue
384  END IF
385  260 continue
386  END IF
387  ELSE
388 *
389 * Form B := alpha*B*inv( A**T ).
390 *
391  IF (upper) THEN
392  DO 310 k = n,1,-1
393  IF (nounit) THEN
394  temp = one/a(k,k)
395  DO 270 i = 1,m
396  b(i,k) = temp*b(i,k)
397  270 continue
398  END IF
399  DO 290 j = 1,k - 1
400  IF (a(j,k).NE.zero) THEN
401  temp = a(j,k)
402  DO 280 i = 1,m
403  b(i,j) = b(i,j) - temp*b(i,k)
404  280 continue
405  END IF
406  290 continue
407  IF (alpha.NE.one) THEN
408  DO 300 i = 1,m
409  b(i,k) = alpha*b(i,k)
410  300 continue
411  END IF
412  310 continue
413  ELSE
414  DO 360 k = 1,n
415  IF (nounit) THEN
416  temp = one/a(k,k)
417  DO 320 i = 1,m
418  b(i,k) = temp*b(i,k)
419  320 continue
420  END IF
421  DO 340 j = k + 1,n
422  IF (a(j,k).NE.zero) THEN
423  temp = a(j,k)
424  DO 330 i = 1,m
425  b(i,j) = b(i,j) - temp*b(i,k)
426  330 continue
427  END IF
428  340 continue
429  IF (alpha.NE.one) THEN
430  DO 350 i = 1,m
431  b(i,k) = alpha*b(i,k)
432  350 continue
433  END IF
434  360 continue
435  END IF
436  END IF
437  END IF
438 *
439  return
440 *
441 * End of DTRSM .
442 *
443  END