1 SUBROUTINE DTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A,
2 $ B, LDB )
3 *
4 * -- LAPACK routine (version 3.3.1) --
5 *
6 * -- Contributed by Fred Gustavson of the IBM Watson Research Center --
7 * -- April 2011 --
8 *
9 * -- LAPACK is a software package provided by Univ. of Tennessee, --
10 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
11 *
12 * ..
13 * .. Scalar Arguments ..
14 CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO
15 INTEGER LDB, M, N
16 DOUBLE PRECISION ALPHA
17 * ..
18 * .. Array Arguments ..
19 DOUBLE PRECISION A( 0: * ), B( 0: LDB-1, 0: * )
20 * ..
21 *
22 * Purpose
23 * =======
24 *
25 * Level 3 BLAS like routine for A in RFP Format.
26 *
27 * DTFSM solves the matrix equation
28 *
29 * op( A )*X = alpha*B or X*op( A ) = alpha*B
30 *
31 * where alpha is a scalar, X and B are m by n matrices, A is a unit, or
32 * non-unit, upper or lower triangular matrix and op( A ) is one of
33 *
34 * op( A ) = A or op( A ) = A**T.
35 *
36 * A is in Rectangular Full Packed (RFP) Format.
37 *
38 * The matrix X is overwritten on B.
39 *
40 * Arguments
41 * ==========
42 *
43 * TRANSR (input) CHARACTER*1
44 * = 'N': The Normal Form of RFP A is stored;
45 * = 'T': The Transpose Form of RFP A is stored.
46 *
47 * SIDE (input) CHARACTER*1
48 * On entry, SIDE specifies whether op( A ) appears on the left
49 * or right of X as follows:
50 *
51 * SIDE = 'L' or 'l' op( A )*X = alpha*B.
52 *
53 * SIDE = 'R' or 'r' X*op( A ) = alpha*B.
54 *
55 * Unchanged on exit.
56 *
57 * UPLO (input) CHARACTER*1
58 * On entry, UPLO specifies whether the RFP matrix A came from
59 * an upper or lower triangular matrix as follows:
60 * UPLO = 'U' or 'u' RFP A came from an upper triangular matrix
61 * UPLO = 'L' or 'l' RFP A came from a lower triangular matrix
62 *
63 * Unchanged on exit.
64 *
65 * TRANS (input) CHARACTER*1
66 * On entry, TRANS specifies the form of op( A ) to be used
67 * in the matrix multiplication as follows:
68 *
69 * TRANS = 'N' or 'n' op( A ) = A.
70 *
71 * TRANS = 'T' or 't' op( A ) = A'.
72 *
73 * Unchanged on exit.
74 *
75 * DIAG (input) CHARACTER*1
76 * On entry, DIAG specifies whether or not RFP A is unit
77 * triangular as follows:
78 *
79 * DIAG = 'U' or 'u' A is assumed to be unit triangular.
80 *
81 * DIAG = 'N' or 'n' A is not assumed to be unit
82 * triangular.
83 *
84 * Unchanged on exit.
85 *
86 * M (input) INTEGER
87 * On entry, M specifies the number of rows of B. M must be at
88 * least zero.
89 * Unchanged on exit.
90 *
91 * N (input) INTEGER
92 * On entry, N specifies the number of columns of B. N must be
93 * at least zero.
94 * Unchanged on exit.
95 *
96 * ALPHA (input) DOUBLE PRECISION
97 * On entry, ALPHA specifies the scalar alpha. When alpha is
98 * zero then A is not referenced and B need not be set before
99 * entry.
100 * Unchanged on exit.
101 *
102 * A (input) DOUBLE PRECISION array, dimension (NT)
103 * NT = N*(N+1)/2. On entry, the matrix A in RFP Format.
104 * RFP Format is described by TRANSR, UPLO and N as follows:
105 * If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even;
106 * K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If
107 * TRANSR = 'T' then RFP is the transpose of RFP A as
108 * defined when TRANSR = 'N'. The contents of RFP A are defined
109 * by UPLO as follows: If UPLO = 'U' the RFP A contains the NT
110 * elements of upper packed A either in normal or
111 * transpose Format. If UPLO = 'L' the RFP A contains
112 * the NT elements of lower packed A either in normal or
113 * transpose Format. The LDA of RFP A is (N+1)/2 when
114 * TRANSR = 'T'. When TRANSR is 'N' the LDA is N+1 when N is
115 * even and is N when is odd.
116 * See the Note below for more details. Unchanged on exit.
117 *
118 * B (input/output) DOUBLE PRECISION array, dimension (LDB,N)
119 * Before entry, the leading m by n part of the array B must
120 * contain the right-hand side matrix B, and on exit is
121 * overwritten by the solution matrix X.
122 *
123 * LDB (input) INTEGER
124 * On entry, LDB specifies the first dimension of B as declared
125 * in the calling (sub) program. LDB must be at least
126 * max( 1, m ).
127 * Unchanged on exit.
128 *
129 * Further Details
130 * ===============
131 *
132 * We first consider Rectangular Full Packed (RFP) Format when N is
133 * even. We give an example where N = 6.
134 *
135 * AP is Upper AP is Lower
136 *
137 * 00 01 02 03 04 05 00
138 * 11 12 13 14 15 10 11
139 * 22 23 24 25 20 21 22
140 * 33 34 35 30 31 32 33
141 * 44 45 40 41 42 43 44
142 * 55 50 51 52 53 54 55
143 *
144 *
145 * Let TRANSR = 'N'. RFP holds AP as follows:
146 * For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
147 * three columns of AP upper. The lower triangle A(4:6,0:2) consists of
148 * the transpose of the first three columns of AP upper.
149 * For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
150 * three columns of AP lower. The upper triangle A(0:2,0:2) consists of
151 * the transpose of the last three columns of AP lower.
152 * This covers the case N even and TRANSR = 'N'.
153 *
154 * RFP A RFP A
155 *
156 * 03 04 05 33 43 53
157 * 13 14 15 00 44 54
158 * 23 24 25 10 11 55
159 * 33 34 35 20 21 22
160 * 00 44 45 30 31 32
161 * 01 11 55 40 41 42
162 * 02 12 22 50 51 52
163 *
164 * Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
165 * transpose of RFP A above. One therefore gets:
166 *
167 *
168 * RFP A RFP A
169 *
170 * 03 13 23 33 00 01 02 33 00 10 20 30 40 50
171 * 04 14 24 34 44 11 12 43 44 11 21 31 41 51
172 * 05 15 25 35 45 55 22 53 54 55 22 32 42 52
173 *
174 *
175 * We then consider Rectangular Full Packed (RFP) Format when N is
176 * odd. We give an example where N = 5.
177 *
178 * AP is Upper AP is Lower
179 *
180 * 00 01 02 03 04 00
181 * 11 12 13 14 10 11
182 * 22 23 24 20 21 22
183 * 33 34 30 31 32 33
184 * 44 40 41 42 43 44
185 *
186 *
187 * Let TRANSR = 'N'. RFP holds AP as follows:
188 * For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
189 * three columns of AP upper. The lower triangle A(3:4,0:1) consists of
190 * the transpose of the first two columns of AP upper.
191 * For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
192 * three columns of AP lower. The upper triangle A(0:1,1:2) consists of
193 * the transpose of the last two columns of AP lower.
194 * This covers the case N odd and TRANSR = 'N'.
195 *
196 * RFP A RFP A
197 *
198 * 02 03 04 00 33 43
199 * 12 13 14 10 11 44
200 * 22 23 24 20 21 22
201 * 00 33 34 30 31 32
202 * 01 11 44 40 41 42
203 *
204 * Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
205 * transpose of RFP A above. One therefore gets:
206 *
207 * RFP A RFP A
208 *
209 * 02 12 22 00 01 00 10 20 30 40 50
210 * 03 13 23 33 11 33 11 21 31 41 51
211 * 04 14 24 34 44 43 44 22 32 42 52
212 *
213 * Reference
214 * =========
215 *
216 * =====================================================================
217 *
218 * ..
219 * .. Parameters ..
220 DOUBLE PRECISION ONE, ZERO
221 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
222 * ..
223 * .. Local Scalars ..
224 LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR,
225 $ NOTRANS
226 INTEGER M1, M2, N1, N2, K, INFO, I, J
227 * ..
228 * .. External Functions ..
229 LOGICAL LSAME
230 EXTERNAL LSAME
231 * ..
232 * .. External Subroutines ..
233 EXTERNAL XERBLA, DGEMM, DTRSM
234 * ..
235 * .. Intrinsic Functions ..
236 INTRINSIC MAX, MOD
237 * ..
238 * .. Executable Statements ..
239 *
240 * Test the input parameters.
241 *
242 INFO = 0
243 NORMALTRANSR = LSAME( TRANSR, 'N' )
244 LSIDE = LSAME( SIDE, 'L' )
245 LOWER = LSAME( UPLO, 'L' )
246 NOTRANS = LSAME( TRANS, 'N' )
247 IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'T' ) ) THEN
248 INFO = -1
249 ELSE IF( .NOT.LSIDE .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
250 INFO = -2
251 ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
252 INFO = -3
253 ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'T' ) ) THEN
254 INFO = -4
255 ELSE IF( .NOT.LSAME( DIAG, 'N' ) .AND. .NOT.LSAME( DIAG, 'U' ) )
256 $ THEN
257 INFO = -5
258 ELSE IF( M.LT.0 ) THEN
259 INFO = -6
260 ELSE IF( N.LT.0 ) THEN
261 INFO = -7
262 ELSE IF( LDB.LT.MAX( 1, M ) ) THEN
263 INFO = -11
264 END IF
265 IF( INFO.NE.0 ) THEN
266 CALL XERBLA( 'DTFSM ', -INFO )
267 RETURN
268 END IF
269 *
270 * Quick return when ( (N.EQ.0).OR.(M.EQ.0) )
271 *
272 IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) )
273 $ RETURN
274 *
275 * Quick return when ALPHA.EQ.(0D+0)
276 *
277 IF( ALPHA.EQ.ZERO ) THEN
278 DO 20 J = 0, N - 1
279 DO 10 I = 0, M - 1
280 B( I, J ) = ZERO
281 10 CONTINUE
282 20 CONTINUE
283 RETURN
284 END IF
285 *
286 IF( LSIDE ) THEN
287 *
288 * SIDE = 'L'
289 *
290 * A is M-by-M.
291 * If M is odd, set NISODD = .TRUE., and M1 and M2.
292 * If M is even, NISODD = .FALSE., and M.
293 *
294 IF( MOD( M, 2 ).EQ.0 ) THEN
295 MISODD = .FALSE.
296 K = M / 2
297 ELSE
298 MISODD = .TRUE.
299 IF( LOWER ) THEN
300 M2 = M / 2
301 M1 = M - M2
302 ELSE
303 M1 = M / 2
304 M2 = M - M1
305 END IF
306 END IF
307 *
308 *
309 IF( MISODD ) THEN
310 *
311 * SIDE = 'L' and N is odd
312 *
313 IF( NORMALTRANSR ) THEN
314 *
315 * SIDE = 'L', N is odd, and TRANSR = 'N'
316 *
317 IF( LOWER ) THEN
318 *
319 * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'L'
320 *
321 IF( NOTRANS ) THEN
322 *
323 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
324 * TRANS = 'N'
325 *
326 IF( M.EQ.1 ) THEN
327 CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
328 $ A, M, B, LDB )
329 ELSE
330 CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
331 $ A( 0 ), M, B, LDB )
332 CALL DGEMM( 'N', 'N', M2, N, M1, -ONE, A( M1 ),
333 $ M, B, LDB, ALPHA, B( M1, 0 ), LDB )
334 CALL DTRSM( 'L', 'U', 'T', DIAG, M2, N, ONE,
335 $ A( M ), M, B( M1, 0 ), LDB )
336 END IF
337 *
338 ELSE
339 *
340 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
341 * TRANS = 'T'
342 *
343 IF( M.EQ.1 ) THEN
344 CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ALPHA,
345 $ A( 0 ), M, B, LDB )
346 ELSE
347 CALL DTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
348 $ A( M ), M, B( M1, 0 ), LDB )
349 CALL DGEMM( 'T', 'N', M1, N, M2, -ONE, A( M1 ),
350 $ M, B( M1, 0 ), LDB, ALPHA, B, LDB )
351 CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ONE,
352 $ A( 0 ), M, B, LDB )
353 END IF
354 *
355 END IF
356 *
357 ELSE
358 *
359 * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'U'
360 *
361 IF( .NOT.NOTRANS ) THEN
362 *
363 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
364 * TRANS = 'N'
365 *
366 CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
367 $ A( M2 ), M, B, LDB )
368 CALL DGEMM( 'T', 'N', M2, N, M1, -ONE, A( 0 ), M,
369 $ B, LDB, ALPHA, B( M1, 0 ), LDB )
370 CALL DTRSM( 'L', 'U', 'T', DIAG, M2, N, ONE,
371 $ A( M1 ), M, B( M1, 0 ), LDB )
372 *
373 ELSE
374 *
375 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
376 * TRANS = 'T'
377 *
378 CALL DTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
379 $ A( M1 ), M, B( M1, 0 ), LDB )
380 CALL DGEMM( 'N', 'N', M1, N, M2, -ONE, A( 0 ), M,
381 $ B( M1, 0 ), LDB, ALPHA, B, LDB )
382 CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ONE,
383 $ A( M2 ), M, B, LDB )
384 *
385 END IF
386 *
387 END IF
388 *
389 ELSE
390 *
391 * SIDE = 'L', N is odd, and TRANSR = 'T'
392 *
393 IF( LOWER ) THEN
394 *
395 * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'L'
396 *
397 IF( NOTRANS ) THEN
398 *
399 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
400 * TRANS = 'N'
401 *
402 IF( M.EQ.1 ) THEN
403 CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA,
404 $ A( 0 ), M1, B, LDB )
405 ELSE
406 CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA,
407 $ A( 0 ), M1, B, LDB )
408 CALL DGEMM( 'T', 'N', M2, N, M1, -ONE,
409 $ A( M1*M1 ), M1, B, LDB, ALPHA,
410 $ B( M1, 0 ), LDB )
411 CALL DTRSM( 'L', 'L', 'N', DIAG, M2, N, ONE,
412 $ A( 1 ), M1, B( M1, 0 ), LDB )
413 END IF
414 *
415 ELSE
416 *
417 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
418 * TRANS = 'T'
419 *
420 IF( M.EQ.1 ) THEN
421 CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ALPHA,
422 $ A( 0 ), M1, B, LDB )
423 ELSE
424 CALL DTRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA,
425 $ A( 1 ), M1, B( M1, 0 ), LDB )
426 CALL DGEMM( 'N', 'N', M1, N, M2, -ONE,
427 $ A( M1*M1 ), M1, B( M1, 0 ), LDB,
428 $ ALPHA, B, LDB )
429 CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ONE,
430 $ A( 0 ), M1, B, LDB )
431 END IF
432 *
433 END IF
434 *
435 ELSE
436 *
437 * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'U'
438 *
439 IF( .NOT.NOTRANS ) THEN
440 *
441 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
442 * TRANS = 'N'
443 *
444 CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA,
445 $ A( M2*M2 ), M2, B, LDB )
446 CALL DGEMM( 'N', 'N', M2, N, M1, -ONE, A( 0 ), M2,
447 $ B, LDB, ALPHA, B( M1, 0 ), LDB )
448 CALL DTRSM( 'L', 'L', 'N', DIAG, M2, N, ONE,
449 $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
450 *
451 ELSE
452 *
453 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
454 * TRANS = 'T'
455 *
456 CALL DTRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA,
457 $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
458 CALL DGEMM( 'T', 'N', M1, N, M2, -ONE, A( 0 ), M2,
459 $ B( M1, 0 ), LDB, ALPHA, B, LDB )
460 CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ONE,
461 $ A( M2*M2 ), M2, B, LDB )
462 *
463 END IF
464 *
465 END IF
466 *
467 END IF
468 *
469 ELSE
470 *
471 * SIDE = 'L' and N is even
472 *
473 IF( NORMALTRANSR ) THEN
474 *
475 * SIDE = 'L', N is even, and TRANSR = 'N'
476 *
477 IF( LOWER ) THEN
478 *
479 * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'L'
480 *
481 IF( NOTRANS ) THEN
482 *
483 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
484 * and TRANS = 'N'
485 *
486 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
487 $ A( 1 ), M+1, B, LDB )
488 CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( K+1 ),
489 $ M+1, B, LDB, ALPHA, B( K, 0 ), LDB )
490 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ONE,
491 $ A( 0 ), M+1, B( K, 0 ), LDB )
492 *
493 ELSE
494 *
495 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
496 * and TRANS = 'T'
497 *
498 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
499 $ A( 0 ), M+1, B( K, 0 ), LDB )
500 CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( K+1 ),
501 $ M+1, B( K, 0 ), LDB, ALPHA, B, LDB )
502 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ONE,
503 $ A( 1 ), M+1, B, LDB )
504 *
505 END IF
506 *
507 ELSE
508 *
509 * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'U'
510 *
511 IF( .NOT.NOTRANS ) THEN
512 *
513 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
514 * and TRANS = 'N'
515 *
516 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
517 $ A( K+1 ), M+1, B, LDB )
518 CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), M+1,
519 $ B, LDB, ALPHA, B( K, 0 ), LDB )
520 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ONE,
521 $ A( K ), M+1, B( K, 0 ), LDB )
522 *
523 ELSE
524 *
525 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
526 * and TRANS = 'T'
527 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
528 $ A( K ), M+1, B( K, 0 ), LDB )
529 CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), M+1,
530 $ B( K, 0 ), LDB, ALPHA, B, LDB )
531 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ONE,
532 $ A( K+1 ), M+1, B, LDB )
533 *
534 END IF
535 *
536 END IF
537 *
538 ELSE
539 *
540 * SIDE = 'L', N is even, and TRANSR = 'T'
541 *
542 IF( LOWER ) THEN
543 *
544 * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'L'
545 *
546 IF( NOTRANS ) THEN
547 *
548 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L',
549 * and TRANS = 'N'
550 *
551 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA,
552 $ A( K ), K, B, LDB )
553 CALL DGEMM( 'T', 'N', K, N, K, -ONE,
554 $ A( K*( K+1 ) ), K, B, LDB, ALPHA,
555 $ B( K, 0 ), LDB )
556 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ONE,
557 $ A( 0 ), K, B( K, 0 ), LDB )
558 *
559 ELSE
560 *
561 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L',
562 * and TRANS = 'T'
563 *
564 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA,
565 $ A( 0 ), K, B( K, 0 ), LDB )
566 CALL DGEMM( 'N', 'N', K, N, K, -ONE,
567 $ A( K*( K+1 ) ), K, B( K, 0 ), LDB,
568 $ ALPHA, B, LDB )
569 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ONE,
570 $ A( K ), K, B, LDB )
571 *
572 END IF
573 *
574 ELSE
575 *
576 * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'U'
577 *
578 IF( .NOT.NOTRANS ) THEN
579 *
580 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U',
581 * and TRANS = 'N'
582 *
583 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA,
584 $ A( K*( K+1 ) ), K, B, LDB )
585 CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), K, B,
586 $ LDB, ALPHA, B( K, 0 ), LDB )
587 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ONE,
588 $ A( K*K ), K, B( K, 0 ), LDB )
589 *
590 ELSE
591 *
592 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U',
593 * and TRANS = 'T'
594 *
595 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA,
596 $ A( K*K ), K, B( K, 0 ), LDB )
597 CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), K,
598 $ B( K, 0 ), LDB, ALPHA, B, LDB )
599 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ONE,
600 $ A( K*( K+1 ) ), K, B, LDB )
601 *
602 END IF
603 *
604 END IF
605 *
606 END IF
607 *
608 END IF
609 *
610 ELSE
611 *
612 * SIDE = 'R'
613 *
614 * A is N-by-N.
615 * If N is odd, set NISODD = .TRUE., and N1 and N2.
616 * If N is even, NISODD = .FALSE., and K.
617 *
618 IF( MOD( N, 2 ).EQ.0 ) THEN
619 NISODD = .FALSE.
620 K = N / 2
621 ELSE
622 NISODD = .TRUE.
623 IF( LOWER ) THEN
624 N2 = N / 2
625 N1 = N - N2
626 ELSE
627 N1 = N / 2
628 N2 = N - N1
629 END IF
630 END IF
631 *
632 IF( NISODD ) THEN
633 *
634 * SIDE = 'R' and N is odd
635 *
636 IF( NORMALTRANSR ) THEN
637 *
638 * SIDE = 'R', N is odd, and TRANSR = 'N'
639 *
640 IF( LOWER ) THEN
641 *
642 * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'L'
643 *
644 IF( NOTRANS ) THEN
645 *
646 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
647 * TRANS = 'N'
648 *
649 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA,
650 $ A( N ), N, B( 0, N1 ), LDB )
651 CALL DGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ),
652 $ LDB, A( N1 ), N, ALPHA, B( 0, 0 ),
653 $ LDB )
654 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N1, ONE,
655 $ A( 0 ), N, B( 0, 0 ), LDB )
656 *
657 ELSE
658 *
659 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
660 * TRANS = 'T'
661 *
662 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA,
663 $ A( 0 ), N, B( 0, 0 ), LDB )
664 CALL DGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ),
665 $ LDB, A( N1 ), N, ALPHA, B( 0, N1 ),
666 $ LDB )
667 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N2, ONE,
668 $ A( N ), N, B( 0, N1 ), LDB )
669 *
670 END IF
671 *
672 ELSE
673 *
674 * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'U'
675 *
676 IF( NOTRANS ) THEN
677 *
678 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
679 * TRANS = 'N'
680 *
681 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA,
682 $ A( N2 ), N, B( 0, 0 ), LDB )
683 CALL DGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ),
684 $ LDB, A( 0 ), N, ALPHA, B( 0, N1 ),
685 $ LDB )
686 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N2, ONE,
687 $ A( N1 ), N, B( 0, N1 ), LDB )
688 *
689 ELSE
690 *
691 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
692 * TRANS = 'T'
693 *
694 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA,
695 $ A( N1 ), N, B( 0, N1 ), LDB )
696 CALL DGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ),
697 $ LDB, A( 0 ), N, ALPHA, B( 0, 0 ), LDB )
698 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N1, ONE,
699 $ A( N2 ), N, B( 0, 0 ), LDB )
700 *
701 END IF
702 *
703 END IF
704 *
705 ELSE
706 *
707 * SIDE = 'R', N is odd, and TRANSR = 'T'
708 *
709 IF( LOWER ) THEN
710 *
711 * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'L'
712 *
713 IF( NOTRANS ) THEN
714 *
715 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
716 * TRANS = 'N'
717 *
718 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
719 $ A( 1 ), N1, B( 0, N1 ), LDB )
720 CALL DGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ),
721 $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, 0 ),
722 $ LDB )
723 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N1, ONE,
724 $ A( 0 ), N1, B( 0, 0 ), LDB )
725 *
726 ELSE
727 *
728 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
729 * TRANS = 'T'
730 *
731 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
732 $ A( 0 ), N1, B( 0, 0 ), LDB )
733 CALL DGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ),
734 $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, N1 ),
735 $ LDB )
736 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N2, ONE,
737 $ A( 1 ), N1, B( 0, N1 ), LDB )
738 *
739 END IF
740 *
741 ELSE
742 *
743 * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'U'
744 *
745 IF( NOTRANS ) THEN
746 *
747 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
748 * TRANS = 'N'
749 *
750 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
751 $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
752 CALL DGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ),
753 $ LDB, A( 0 ), N2, ALPHA, B( 0, N1 ),
754 $ LDB )
755 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N2, ONE,
756 $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
757 *
758 ELSE
759 *
760 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
761 * TRANS = 'T'
762 *
763 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
764 $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
765 CALL DGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ),
766 $ LDB, A( 0 ), N2, ALPHA, B( 0, 0 ),
767 $ LDB )
768 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N1, ONE,
769 $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
770 *
771 END IF
772 *
773 END IF
774 *
775 END IF
776 *
777 ELSE
778 *
779 * SIDE = 'R' and N is even
780 *
781 IF( NORMALTRANSR ) THEN
782 *
783 * SIDE = 'R', N is even, and TRANSR = 'N'
784 *
785 IF( LOWER ) THEN
786 *
787 * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'L'
788 *
789 IF( NOTRANS ) THEN
790 *
791 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
792 * and TRANS = 'N'
793 *
794 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA,
795 $ A( 0 ), N+1, B( 0, K ), LDB )
796 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ),
797 $ LDB, A( K+1 ), N+1, ALPHA, B( 0, 0 ),
798 $ LDB )
799 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ONE,
800 $ A( 1 ), N+1, B( 0, 0 ), LDB )
801 *
802 ELSE
803 *
804 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
805 * and TRANS = 'T'
806 *
807 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA,
808 $ A( 1 ), N+1, B( 0, 0 ), LDB )
809 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ),
810 $ LDB, A( K+1 ), N+1, ALPHA, B( 0, K ),
811 $ LDB )
812 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ONE,
813 $ A( 0 ), N+1, B( 0, K ), LDB )
814 *
815 END IF
816 *
817 ELSE
818 *
819 * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'U'
820 *
821 IF( NOTRANS ) THEN
822 *
823 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
824 * and TRANS = 'N'
825 *
826 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA,
827 $ A( K+1 ), N+1, B( 0, 0 ), LDB )
828 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ),
829 $ LDB, A( 0 ), N+1, ALPHA, B( 0, K ),
830 $ LDB )
831 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ONE,
832 $ A( K ), N+1, B( 0, K ), LDB )
833 *
834 ELSE
835 *
836 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
837 * and TRANS = 'T'
838 *
839 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA,
840 $ A( K ), N+1, B( 0, K ), LDB )
841 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ),
842 $ LDB, A( 0 ), N+1, ALPHA, B( 0, 0 ),
843 $ LDB )
844 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ONE,
845 $ A( K+1 ), N+1, B( 0, 0 ), LDB )
846 *
847 END IF
848 *
849 END IF
850 *
851 ELSE
852 *
853 * SIDE = 'R', N is even, and TRANSR = 'T'
854 *
855 IF( LOWER ) THEN
856 *
857 * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'L'
858 *
859 IF( NOTRANS ) THEN
860 *
861 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L',
862 * and TRANS = 'N'
863 *
864 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
865 $ A( 0 ), K, B( 0, K ), LDB )
866 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ),
867 $ LDB, A( ( K+1 )*K ), K, ALPHA,
868 $ B( 0, 0 ), LDB )
869 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ONE,
870 $ A( K ), K, B( 0, 0 ), LDB )
871 *
872 ELSE
873 *
874 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L',
875 * and TRANS = 'T'
876 *
877 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
878 $ A( K ), K, B( 0, 0 ), LDB )
879 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ),
880 $ LDB, A( ( K+1 )*K ), K, ALPHA,
881 $ B( 0, K ), LDB )
882 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ONE,
883 $ A( 0 ), K, B( 0, K ), LDB )
884 *
885 END IF
886 *
887 ELSE
888 *
889 * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'U'
890 *
891 IF( NOTRANS ) THEN
892 *
893 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U',
894 * and TRANS = 'N'
895 *
896 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
897 $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
898 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ),
899 $ LDB, A( 0 ), K, ALPHA, B( 0, K ), LDB )
900 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ONE,
901 $ A( K*K ), K, B( 0, K ), LDB )
902 *
903 ELSE
904 *
905 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U',
906 * and TRANS = 'T'
907 *
908 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
909 $ A( K*K ), K, B( 0, K ), LDB )
910 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ),
911 $ LDB, A( 0 ), K, ALPHA, B( 0, 0 ), LDB )
912 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ONE,
913 $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
914 *
915 END IF
916 *
917 END IF
918 *
919 END IF
920 *
921 END IF
922 END IF
923 *
924 RETURN
925 *
926 * End of DTFSM
927 *
928 END
2 $ B, LDB )
3 *
4 * -- LAPACK routine (version 3.3.1) --
5 *
6 * -- Contributed by Fred Gustavson of the IBM Watson Research Center --
7 * -- April 2011 --
8 *
9 * -- LAPACK is a software package provided by Univ. of Tennessee, --
10 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
11 *
12 * ..
13 * .. Scalar Arguments ..
14 CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO
15 INTEGER LDB, M, N
16 DOUBLE PRECISION ALPHA
17 * ..
18 * .. Array Arguments ..
19 DOUBLE PRECISION A( 0: * ), B( 0: LDB-1, 0: * )
20 * ..
21 *
22 * Purpose
23 * =======
24 *
25 * Level 3 BLAS like routine for A in RFP Format.
26 *
27 * DTFSM solves the matrix equation
28 *
29 * op( A )*X = alpha*B or X*op( A ) = alpha*B
30 *
31 * where alpha is a scalar, X and B are m by n matrices, A is a unit, or
32 * non-unit, upper or lower triangular matrix and op( A ) is one of
33 *
34 * op( A ) = A or op( A ) = A**T.
35 *
36 * A is in Rectangular Full Packed (RFP) Format.
37 *
38 * The matrix X is overwritten on B.
39 *
40 * Arguments
41 * ==========
42 *
43 * TRANSR (input) CHARACTER*1
44 * = 'N': The Normal Form of RFP A is stored;
45 * = 'T': The Transpose Form of RFP A is stored.
46 *
47 * SIDE (input) CHARACTER*1
48 * On entry, SIDE specifies whether op( A ) appears on the left
49 * or right of X as follows:
50 *
51 * SIDE = 'L' or 'l' op( A )*X = alpha*B.
52 *
53 * SIDE = 'R' or 'r' X*op( A ) = alpha*B.
54 *
55 * Unchanged on exit.
56 *
57 * UPLO (input) CHARACTER*1
58 * On entry, UPLO specifies whether the RFP matrix A came from
59 * an upper or lower triangular matrix as follows:
60 * UPLO = 'U' or 'u' RFP A came from an upper triangular matrix
61 * UPLO = 'L' or 'l' RFP A came from a lower triangular matrix
62 *
63 * Unchanged on exit.
64 *
65 * TRANS (input) CHARACTER*1
66 * On entry, TRANS specifies the form of op( A ) to be used
67 * in the matrix multiplication as follows:
68 *
69 * TRANS = 'N' or 'n' op( A ) = A.
70 *
71 * TRANS = 'T' or 't' op( A ) = A'.
72 *
73 * Unchanged on exit.
74 *
75 * DIAG (input) CHARACTER*1
76 * On entry, DIAG specifies whether or not RFP A is unit
77 * triangular as follows:
78 *
79 * DIAG = 'U' or 'u' A is assumed to be unit triangular.
80 *
81 * DIAG = 'N' or 'n' A is not assumed to be unit
82 * triangular.
83 *
84 * Unchanged on exit.
85 *
86 * M (input) INTEGER
87 * On entry, M specifies the number of rows of B. M must be at
88 * least zero.
89 * Unchanged on exit.
90 *
91 * N (input) INTEGER
92 * On entry, N specifies the number of columns of B. N must be
93 * at least zero.
94 * Unchanged on exit.
95 *
96 * ALPHA (input) DOUBLE PRECISION
97 * On entry, ALPHA specifies the scalar alpha. When alpha is
98 * zero then A is not referenced and B need not be set before
99 * entry.
100 * Unchanged on exit.
101 *
102 * A (input) DOUBLE PRECISION array, dimension (NT)
103 * NT = N*(N+1)/2. On entry, the matrix A in RFP Format.
104 * RFP Format is described by TRANSR, UPLO and N as follows:
105 * If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even;
106 * K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If
107 * TRANSR = 'T' then RFP is the transpose of RFP A as
108 * defined when TRANSR = 'N'. The contents of RFP A are defined
109 * by UPLO as follows: If UPLO = 'U' the RFP A contains the NT
110 * elements of upper packed A either in normal or
111 * transpose Format. If UPLO = 'L' the RFP A contains
112 * the NT elements of lower packed A either in normal or
113 * transpose Format. The LDA of RFP A is (N+1)/2 when
114 * TRANSR = 'T'. When TRANSR is 'N' the LDA is N+1 when N is
115 * even and is N when is odd.
116 * See the Note below for more details. Unchanged on exit.
117 *
118 * B (input/output) DOUBLE PRECISION array, dimension (LDB,N)
119 * Before entry, the leading m by n part of the array B must
120 * contain the right-hand side matrix B, and on exit is
121 * overwritten by the solution matrix X.
122 *
123 * LDB (input) INTEGER
124 * On entry, LDB specifies the first dimension of B as declared
125 * in the calling (sub) program. LDB must be at least
126 * max( 1, m ).
127 * Unchanged on exit.
128 *
129 * Further Details
130 * ===============
131 *
132 * We first consider Rectangular Full Packed (RFP) Format when N is
133 * even. We give an example where N = 6.
134 *
135 * AP is Upper AP is Lower
136 *
137 * 00 01 02 03 04 05 00
138 * 11 12 13 14 15 10 11
139 * 22 23 24 25 20 21 22
140 * 33 34 35 30 31 32 33
141 * 44 45 40 41 42 43 44
142 * 55 50 51 52 53 54 55
143 *
144 *
145 * Let TRANSR = 'N'. RFP holds AP as follows:
146 * For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
147 * three columns of AP upper. The lower triangle A(4:6,0:2) consists of
148 * the transpose of the first three columns of AP upper.
149 * For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
150 * three columns of AP lower. The upper triangle A(0:2,0:2) consists of
151 * the transpose of the last three columns of AP lower.
152 * This covers the case N even and TRANSR = 'N'.
153 *
154 * RFP A RFP A
155 *
156 * 03 04 05 33 43 53
157 * 13 14 15 00 44 54
158 * 23 24 25 10 11 55
159 * 33 34 35 20 21 22
160 * 00 44 45 30 31 32
161 * 01 11 55 40 41 42
162 * 02 12 22 50 51 52
163 *
164 * Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
165 * transpose of RFP A above. One therefore gets:
166 *
167 *
168 * RFP A RFP A
169 *
170 * 03 13 23 33 00 01 02 33 00 10 20 30 40 50
171 * 04 14 24 34 44 11 12 43 44 11 21 31 41 51
172 * 05 15 25 35 45 55 22 53 54 55 22 32 42 52
173 *
174 *
175 * We then consider Rectangular Full Packed (RFP) Format when N is
176 * odd. We give an example where N = 5.
177 *
178 * AP is Upper AP is Lower
179 *
180 * 00 01 02 03 04 00
181 * 11 12 13 14 10 11
182 * 22 23 24 20 21 22
183 * 33 34 30 31 32 33
184 * 44 40 41 42 43 44
185 *
186 *
187 * Let TRANSR = 'N'. RFP holds AP as follows:
188 * For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
189 * three columns of AP upper. The lower triangle A(3:4,0:1) consists of
190 * the transpose of the first two columns of AP upper.
191 * For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
192 * three columns of AP lower. The upper triangle A(0:1,1:2) consists of
193 * the transpose of the last two columns of AP lower.
194 * This covers the case N odd and TRANSR = 'N'.
195 *
196 * RFP A RFP A
197 *
198 * 02 03 04 00 33 43
199 * 12 13 14 10 11 44
200 * 22 23 24 20 21 22
201 * 00 33 34 30 31 32
202 * 01 11 44 40 41 42
203 *
204 * Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
205 * transpose of RFP A above. One therefore gets:
206 *
207 * RFP A RFP A
208 *
209 * 02 12 22 00 01 00 10 20 30 40 50
210 * 03 13 23 33 11 33 11 21 31 41 51
211 * 04 14 24 34 44 43 44 22 32 42 52
212 *
213 * Reference
214 * =========
215 *
216 * =====================================================================
217 *
218 * ..
219 * .. Parameters ..
220 DOUBLE PRECISION ONE, ZERO
221 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
222 * ..
223 * .. Local Scalars ..
224 LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR,
225 $ NOTRANS
226 INTEGER M1, M2, N1, N2, K, INFO, I, J
227 * ..
228 * .. External Functions ..
229 LOGICAL LSAME
230 EXTERNAL LSAME
231 * ..
232 * .. External Subroutines ..
233 EXTERNAL XERBLA, DGEMM, DTRSM
234 * ..
235 * .. Intrinsic Functions ..
236 INTRINSIC MAX, MOD
237 * ..
238 * .. Executable Statements ..
239 *
240 * Test the input parameters.
241 *
242 INFO = 0
243 NORMALTRANSR = LSAME( TRANSR, 'N' )
244 LSIDE = LSAME( SIDE, 'L' )
245 LOWER = LSAME( UPLO, 'L' )
246 NOTRANS = LSAME( TRANS, 'N' )
247 IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'T' ) ) THEN
248 INFO = -1
249 ELSE IF( .NOT.LSIDE .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
250 INFO = -2
251 ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
252 INFO = -3
253 ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'T' ) ) THEN
254 INFO = -4
255 ELSE IF( .NOT.LSAME( DIAG, 'N' ) .AND. .NOT.LSAME( DIAG, 'U' ) )
256 $ THEN
257 INFO = -5
258 ELSE IF( M.LT.0 ) THEN
259 INFO = -6
260 ELSE IF( N.LT.0 ) THEN
261 INFO = -7
262 ELSE IF( LDB.LT.MAX( 1, M ) ) THEN
263 INFO = -11
264 END IF
265 IF( INFO.NE.0 ) THEN
266 CALL XERBLA( 'DTFSM ', -INFO )
267 RETURN
268 END IF
269 *
270 * Quick return when ( (N.EQ.0).OR.(M.EQ.0) )
271 *
272 IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) )
273 $ RETURN
274 *
275 * Quick return when ALPHA.EQ.(0D+0)
276 *
277 IF( ALPHA.EQ.ZERO ) THEN
278 DO 20 J = 0, N - 1
279 DO 10 I = 0, M - 1
280 B( I, J ) = ZERO
281 10 CONTINUE
282 20 CONTINUE
283 RETURN
284 END IF
285 *
286 IF( LSIDE ) THEN
287 *
288 * SIDE = 'L'
289 *
290 * A is M-by-M.
291 * If M is odd, set NISODD = .TRUE., and M1 and M2.
292 * If M is even, NISODD = .FALSE., and M.
293 *
294 IF( MOD( M, 2 ).EQ.0 ) THEN
295 MISODD = .FALSE.
296 K = M / 2
297 ELSE
298 MISODD = .TRUE.
299 IF( LOWER ) THEN
300 M2 = M / 2
301 M1 = M - M2
302 ELSE
303 M1 = M / 2
304 M2 = M - M1
305 END IF
306 END IF
307 *
308 *
309 IF( MISODD ) THEN
310 *
311 * SIDE = 'L' and N is odd
312 *
313 IF( NORMALTRANSR ) THEN
314 *
315 * SIDE = 'L', N is odd, and TRANSR = 'N'
316 *
317 IF( LOWER ) THEN
318 *
319 * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'L'
320 *
321 IF( NOTRANS ) THEN
322 *
323 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
324 * TRANS = 'N'
325 *
326 IF( M.EQ.1 ) THEN
327 CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
328 $ A, M, B, LDB )
329 ELSE
330 CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
331 $ A( 0 ), M, B, LDB )
332 CALL DGEMM( 'N', 'N', M2, N, M1, -ONE, A( M1 ),
333 $ M, B, LDB, ALPHA, B( M1, 0 ), LDB )
334 CALL DTRSM( 'L', 'U', 'T', DIAG, M2, N, ONE,
335 $ A( M ), M, B( M1, 0 ), LDB )
336 END IF
337 *
338 ELSE
339 *
340 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
341 * TRANS = 'T'
342 *
343 IF( M.EQ.1 ) THEN
344 CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ALPHA,
345 $ A( 0 ), M, B, LDB )
346 ELSE
347 CALL DTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
348 $ A( M ), M, B( M1, 0 ), LDB )
349 CALL DGEMM( 'T', 'N', M1, N, M2, -ONE, A( M1 ),
350 $ M, B( M1, 0 ), LDB, ALPHA, B, LDB )
351 CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ONE,
352 $ A( 0 ), M, B, LDB )
353 END IF
354 *
355 END IF
356 *
357 ELSE
358 *
359 * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'U'
360 *
361 IF( .NOT.NOTRANS ) THEN
362 *
363 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
364 * TRANS = 'N'
365 *
366 CALL DTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
367 $ A( M2 ), M, B, LDB )
368 CALL DGEMM( 'T', 'N', M2, N, M1, -ONE, A( 0 ), M,
369 $ B, LDB, ALPHA, B( M1, 0 ), LDB )
370 CALL DTRSM( 'L', 'U', 'T', DIAG, M2, N, ONE,
371 $ A( M1 ), M, B( M1, 0 ), LDB )
372 *
373 ELSE
374 *
375 * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
376 * TRANS = 'T'
377 *
378 CALL DTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
379 $ A( M1 ), M, B( M1, 0 ), LDB )
380 CALL DGEMM( 'N', 'N', M1, N, M2, -ONE, A( 0 ), M,
381 $ B( M1, 0 ), LDB, ALPHA, B, LDB )
382 CALL DTRSM( 'L', 'L', 'T', DIAG, M1, N, ONE,
383 $ A( M2 ), M, B, LDB )
384 *
385 END IF
386 *
387 END IF
388 *
389 ELSE
390 *
391 * SIDE = 'L', N is odd, and TRANSR = 'T'
392 *
393 IF( LOWER ) THEN
394 *
395 * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'L'
396 *
397 IF( NOTRANS ) THEN
398 *
399 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
400 * TRANS = 'N'
401 *
402 IF( M.EQ.1 ) THEN
403 CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA,
404 $ A( 0 ), M1, B, LDB )
405 ELSE
406 CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA,
407 $ A( 0 ), M1, B, LDB )
408 CALL DGEMM( 'T', 'N', M2, N, M1, -ONE,
409 $ A( M1*M1 ), M1, B, LDB, ALPHA,
410 $ B( M1, 0 ), LDB )
411 CALL DTRSM( 'L', 'L', 'N', DIAG, M2, N, ONE,
412 $ A( 1 ), M1, B( M1, 0 ), LDB )
413 END IF
414 *
415 ELSE
416 *
417 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
418 * TRANS = 'T'
419 *
420 IF( M.EQ.1 ) THEN
421 CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ALPHA,
422 $ A( 0 ), M1, B, LDB )
423 ELSE
424 CALL DTRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA,
425 $ A( 1 ), M1, B( M1, 0 ), LDB )
426 CALL DGEMM( 'N', 'N', M1, N, M2, -ONE,
427 $ A( M1*M1 ), M1, B( M1, 0 ), LDB,
428 $ ALPHA, B, LDB )
429 CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ONE,
430 $ A( 0 ), M1, B, LDB )
431 END IF
432 *
433 END IF
434 *
435 ELSE
436 *
437 * SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'U'
438 *
439 IF( .NOT.NOTRANS ) THEN
440 *
441 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
442 * TRANS = 'N'
443 *
444 CALL DTRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA,
445 $ A( M2*M2 ), M2, B, LDB )
446 CALL DGEMM( 'N', 'N', M2, N, M1, -ONE, A( 0 ), M2,
447 $ B, LDB, ALPHA, B( M1, 0 ), LDB )
448 CALL DTRSM( 'L', 'L', 'N', DIAG, M2, N, ONE,
449 $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
450 *
451 ELSE
452 *
453 * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
454 * TRANS = 'T'
455 *
456 CALL DTRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA,
457 $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
458 CALL DGEMM( 'T', 'N', M1, N, M2, -ONE, A( 0 ), M2,
459 $ B( M1, 0 ), LDB, ALPHA, B, LDB )
460 CALL DTRSM( 'L', 'U', 'N', DIAG, M1, N, ONE,
461 $ A( M2*M2 ), M2, B, LDB )
462 *
463 END IF
464 *
465 END IF
466 *
467 END IF
468 *
469 ELSE
470 *
471 * SIDE = 'L' and N is even
472 *
473 IF( NORMALTRANSR ) THEN
474 *
475 * SIDE = 'L', N is even, and TRANSR = 'N'
476 *
477 IF( LOWER ) THEN
478 *
479 * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'L'
480 *
481 IF( NOTRANS ) THEN
482 *
483 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
484 * and TRANS = 'N'
485 *
486 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
487 $ A( 1 ), M+1, B, LDB )
488 CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( K+1 ),
489 $ M+1, B, LDB, ALPHA, B( K, 0 ), LDB )
490 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ONE,
491 $ A( 0 ), M+1, B( K, 0 ), LDB )
492 *
493 ELSE
494 *
495 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
496 * and TRANS = 'T'
497 *
498 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
499 $ A( 0 ), M+1, B( K, 0 ), LDB )
500 CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( K+1 ),
501 $ M+1, B( K, 0 ), LDB, ALPHA, B, LDB )
502 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ONE,
503 $ A( 1 ), M+1, B, LDB )
504 *
505 END IF
506 *
507 ELSE
508 *
509 * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'U'
510 *
511 IF( .NOT.NOTRANS ) THEN
512 *
513 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
514 * and TRANS = 'N'
515 *
516 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
517 $ A( K+1 ), M+1, B, LDB )
518 CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), M+1,
519 $ B, LDB, ALPHA, B( K, 0 ), LDB )
520 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ONE,
521 $ A( K ), M+1, B( K, 0 ), LDB )
522 *
523 ELSE
524 *
525 * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
526 * and TRANS = 'T'
527 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
528 $ A( K ), M+1, B( K, 0 ), LDB )
529 CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), M+1,
530 $ B( K, 0 ), LDB, ALPHA, B, LDB )
531 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ONE,
532 $ A( K+1 ), M+1, B, LDB )
533 *
534 END IF
535 *
536 END IF
537 *
538 ELSE
539 *
540 * SIDE = 'L', N is even, and TRANSR = 'T'
541 *
542 IF( LOWER ) THEN
543 *
544 * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'L'
545 *
546 IF( NOTRANS ) THEN
547 *
548 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L',
549 * and TRANS = 'N'
550 *
551 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA,
552 $ A( K ), K, B, LDB )
553 CALL DGEMM( 'T', 'N', K, N, K, -ONE,
554 $ A( K*( K+1 ) ), K, B, LDB, ALPHA,
555 $ B( K, 0 ), LDB )
556 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ONE,
557 $ A( 0 ), K, B( K, 0 ), LDB )
558 *
559 ELSE
560 *
561 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L',
562 * and TRANS = 'T'
563 *
564 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA,
565 $ A( 0 ), K, B( K, 0 ), LDB )
566 CALL DGEMM( 'N', 'N', K, N, K, -ONE,
567 $ A( K*( K+1 ) ), K, B( K, 0 ), LDB,
568 $ ALPHA, B, LDB )
569 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ONE,
570 $ A( K ), K, B, LDB )
571 *
572 END IF
573 *
574 ELSE
575 *
576 * SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'U'
577 *
578 IF( .NOT.NOTRANS ) THEN
579 *
580 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U',
581 * and TRANS = 'N'
582 *
583 CALL DTRSM( 'L', 'U', 'T', DIAG, K, N, ALPHA,
584 $ A( K*( K+1 ) ), K, B, LDB )
585 CALL DGEMM( 'N', 'N', K, N, K, -ONE, A( 0 ), K, B,
586 $ LDB, ALPHA, B( K, 0 ), LDB )
587 CALL DTRSM( 'L', 'L', 'N', DIAG, K, N, ONE,
588 $ A( K*K ), K, B( K, 0 ), LDB )
589 *
590 ELSE
591 *
592 * SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U',
593 * and TRANS = 'T'
594 *
595 CALL DTRSM( 'L', 'L', 'T', DIAG, K, N, ALPHA,
596 $ A( K*K ), K, B( K, 0 ), LDB )
597 CALL DGEMM( 'T', 'N', K, N, K, -ONE, A( 0 ), K,
598 $ B( K, 0 ), LDB, ALPHA, B, LDB )
599 CALL DTRSM( 'L', 'U', 'N', DIAG, K, N, ONE,
600 $ A( K*( K+1 ) ), K, B, LDB )
601 *
602 END IF
603 *
604 END IF
605 *
606 END IF
607 *
608 END IF
609 *
610 ELSE
611 *
612 * SIDE = 'R'
613 *
614 * A is N-by-N.
615 * If N is odd, set NISODD = .TRUE., and N1 and N2.
616 * If N is even, NISODD = .FALSE., and K.
617 *
618 IF( MOD( N, 2 ).EQ.0 ) THEN
619 NISODD = .FALSE.
620 K = N / 2
621 ELSE
622 NISODD = .TRUE.
623 IF( LOWER ) THEN
624 N2 = N / 2
625 N1 = N - N2
626 ELSE
627 N1 = N / 2
628 N2 = N - N1
629 END IF
630 END IF
631 *
632 IF( NISODD ) THEN
633 *
634 * SIDE = 'R' and N is odd
635 *
636 IF( NORMALTRANSR ) THEN
637 *
638 * SIDE = 'R', N is odd, and TRANSR = 'N'
639 *
640 IF( LOWER ) THEN
641 *
642 * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'L'
643 *
644 IF( NOTRANS ) THEN
645 *
646 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
647 * TRANS = 'N'
648 *
649 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA,
650 $ A( N ), N, B( 0, N1 ), LDB )
651 CALL DGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ),
652 $ LDB, A( N1 ), N, ALPHA, B( 0, 0 ),
653 $ LDB )
654 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N1, ONE,
655 $ A( 0 ), N, B( 0, 0 ), LDB )
656 *
657 ELSE
658 *
659 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
660 * TRANS = 'T'
661 *
662 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA,
663 $ A( 0 ), N, B( 0, 0 ), LDB )
664 CALL DGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ),
665 $ LDB, A( N1 ), N, ALPHA, B( 0, N1 ),
666 $ LDB )
667 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N2, ONE,
668 $ A( N ), N, B( 0, N1 ), LDB )
669 *
670 END IF
671 *
672 ELSE
673 *
674 * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'U'
675 *
676 IF( NOTRANS ) THEN
677 *
678 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
679 * TRANS = 'N'
680 *
681 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N1, ALPHA,
682 $ A( N2 ), N, B( 0, 0 ), LDB )
683 CALL DGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ),
684 $ LDB, A( 0 ), N, ALPHA, B( 0, N1 ),
685 $ LDB )
686 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N2, ONE,
687 $ A( N1 ), N, B( 0, N1 ), LDB )
688 *
689 ELSE
690 *
691 * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
692 * TRANS = 'T'
693 *
694 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N2, ALPHA,
695 $ A( N1 ), N, B( 0, N1 ), LDB )
696 CALL DGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ),
697 $ LDB, A( 0 ), N, ALPHA, B( 0, 0 ), LDB )
698 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N1, ONE,
699 $ A( N2 ), N, B( 0, 0 ), LDB )
700 *
701 END IF
702 *
703 END IF
704 *
705 ELSE
706 *
707 * SIDE = 'R', N is odd, and TRANSR = 'T'
708 *
709 IF( LOWER ) THEN
710 *
711 * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'L'
712 *
713 IF( NOTRANS ) THEN
714 *
715 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
716 * TRANS = 'N'
717 *
718 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
719 $ A( 1 ), N1, B( 0, N1 ), LDB )
720 CALL DGEMM( 'N', 'T', M, N1, N2, -ONE, B( 0, N1 ),
721 $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, 0 ),
722 $ LDB )
723 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N1, ONE,
724 $ A( 0 ), N1, B( 0, 0 ), LDB )
725 *
726 ELSE
727 *
728 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
729 * TRANS = 'T'
730 *
731 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
732 $ A( 0 ), N1, B( 0, 0 ), LDB )
733 CALL DGEMM( 'N', 'N', M, N2, N1, -ONE, B( 0, 0 ),
734 $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, N1 ),
735 $ LDB )
736 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N2, ONE,
737 $ A( 1 ), N1, B( 0, N1 ), LDB )
738 *
739 END IF
740 *
741 ELSE
742 *
743 * SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'U'
744 *
745 IF( NOTRANS ) THEN
746 *
747 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
748 * TRANS = 'N'
749 *
750 CALL DTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
751 $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
752 CALL DGEMM( 'N', 'T', M, N2, N1, -ONE, B( 0, 0 ),
753 $ LDB, A( 0 ), N2, ALPHA, B( 0, N1 ),
754 $ LDB )
755 CALL DTRSM( 'R', 'L', 'T', DIAG, M, N2, ONE,
756 $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
757 *
758 ELSE
759 *
760 * SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
761 * TRANS = 'T'
762 *
763 CALL DTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
764 $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
765 CALL DGEMM( 'N', 'N', M, N1, N2, -ONE, B( 0, N1 ),
766 $ LDB, A( 0 ), N2, ALPHA, B( 0, 0 ),
767 $ LDB )
768 CALL DTRSM( 'R', 'U', 'T', DIAG, M, N1, ONE,
769 $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
770 *
771 END IF
772 *
773 END IF
774 *
775 END IF
776 *
777 ELSE
778 *
779 * SIDE = 'R' and N is even
780 *
781 IF( NORMALTRANSR ) THEN
782 *
783 * SIDE = 'R', N is even, and TRANSR = 'N'
784 *
785 IF( LOWER ) THEN
786 *
787 * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'L'
788 *
789 IF( NOTRANS ) THEN
790 *
791 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
792 * and TRANS = 'N'
793 *
794 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA,
795 $ A( 0 ), N+1, B( 0, K ), LDB )
796 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ),
797 $ LDB, A( K+1 ), N+1, ALPHA, B( 0, 0 ),
798 $ LDB )
799 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ONE,
800 $ A( 1 ), N+1, B( 0, 0 ), LDB )
801 *
802 ELSE
803 *
804 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
805 * and TRANS = 'T'
806 *
807 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA,
808 $ A( 1 ), N+1, B( 0, 0 ), LDB )
809 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ),
810 $ LDB, A( K+1 ), N+1, ALPHA, B( 0, K ),
811 $ LDB )
812 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ONE,
813 $ A( 0 ), N+1, B( 0, K ), LDB )
814 *
815 END IF
816 *
817 ELSE
818 *
819 * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'U'
820 *
821 IF( NOTRANS ) THEN
822 *
823 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
824 * and TRANS = 'N'
825 *
826 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ALPHA,
827 $ A( K+1 ), N+1, B( 0, 0 ), LDB )
828 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ),
829 $ LDB, A( 0 ), N+1, ALPHA, B( 0, K ),
830 $ LDB )
831 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ONE,
832 $ A( K ), N+1, B( 0, K ), LDB )
833 *
834 ELSE
835 *
836 * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
837 * and TRANS = 'T'
838 *
839 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ALPHA,
840 $ A( K ), N+1, B( 0, K ), LDB )
841 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ),
842 $ LDB, A( 0 ), N+1, ALPHA, B( 0, 0 ),
843 $ LDB )
844 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ONE,
845 $ A( K+1 ), N+1, B( 0, 0 ), LDB )
846 *
847 END IF
848 *
849 END IF
850 *
851 ELSE
852 *
853 * SIDE = 'R', N is even, and TRANSR = 'T'
854 *
855 IF( LOWER ) THEN
856 *
857 * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'L'
858 *
859 IF( NOTRANS ) THEN
860 *
861 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L',
862 * and TRANS = 'N'
863 *
864 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
865 $ A( 0 ), K, B( 0, K ), LDB )
866 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, K ),
867 $ LDB, A( ( K+1 )*K ), K, ALPHA,
868 $ B( 0, 0 ), LDB )
869 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ONE,
870 $ A( K ), K, B( 0, 0 ), LDB )
871 *
872 ELSE
873 *
874 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L',
875 * and TRANS = 'T'
876 *
877 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
878 $ A( K ), K, B( 0, 0 ), LDB )
879 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, 0 ),
880 $ LDB, A( ( K+1 )*K ), K, ALPHA,
881 $ B( 0, K ), LDB )
882 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ONE,
883 $ A( 0 ), K, B( 0, K ), LDB )
884 *
885 END IF
886 *
887 ELSE
888 *
889 * SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'U'
890 *
891 IF( NOTRANS ) THEN
892 *
893 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U',
894 * and TRANS = 'N'
895 *
896 CALL DTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
897 $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
898 CALL DGEMM( 'N', 'T', M, K, K, -ONE, B( 0, 0 ),
899 $ LDB, A( 0 ), K, ALPHA, B( 0, K ), LDB )
900 CALL DTRSM( 'R', 'L', 'T', DIAG, M, K, ONE,
901 $ A( K*K ), K, B( 0, K ), LDB )
902 *
903 ELSE
904 *
905 * SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U',
906 * and TRANS = 'T'
907 *
908 CALL DTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
909 $ A( K*K ), K, B( 0, K ), LDB )
910 CALL DGEMM( 'N', 'N', M, K, K, -ONE, B( 0, K ),
911 $ LDB, A( 0 ), K, ALPHA, B( 0, 0 ), LDB )
912 CALL DTRSM( 'R', 'U', 'T', DIAG, M, K, ONE,
913 $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
914 *
915 END IF
916 *
917 END IF
918 *
919 END IF
920 *
921 END IF
922 END IF
923 *
924 RETURN
925 *
926 * End of DTFSM
927 *
928 END