LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ dtfsm()

subroutine dtfsm ( character transr,
character side,
character uplo,
character trans,
character diag,
integer m,
integer n,
double precision alpha,
double precision, dimension( 0: * ) a,
double precision, dimension( 0: ldb-1, 0: * ) b,
integer ldb )

DTFSM solves a matrix equation (one operand is a triangular matrix in RFP format).

Download DTFSM + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!> Level 3 BLAS like routine for A in RFP Format.
!>
!> DTFSM  solves the matrix equation
!>
!>    op( A )*X = alpha*B  or  X*op( A ) = alpha*B
!>
!> where alpha is a scalar, X and B are m by n matrices, A is a unit, or
!> non-unit,  upper or lower triangular matrix  and  op( A )  is one  of
!>
!>    op( A ) = A   or   op( A ) = A**T.
!>
!> A is in Rectangular Full Packed (RFP) Format.
!>
!> The matrix X is overwritten on B.
!>
Parameters
[in]TRANSR
!>          TRANSR is CHARACTER*1
!>          = 'N':  The Normal Form of RFP A is stored;
!>          = 'T':  The Transpose Form of RFP A is stored.
!>
[in]SIDE
!>          SIDE is CHARACTER*1
!>           On entry, SIDE specifies whether op( A ) appears on the left
!>           or right of X as follows:
!>
!>              SIDE = 'L' or 'l'   op( A )*X = alpha*B.
!>
!>              SIDE = 'R' or 'r'   X*op( A ) = alpha*B.
!>
!>           Unchanged on exit.
!>
[in]UPLO
!>          UPLO is CHARACTER*1
!>           On entry, UPLO specifies whether the RFP matrix A came from
!>           an upper or lower triangular matrix as follows:
!>           UPLO = 'U' or 'u' RFP A came from an upper triangular matrix
!>           UPLO = 'L' or 'l' RFP A came from a  lower triangular matrix
!>
!>           Unchanged on exit.
!>
[in]TRANS
!>          TRANS is CHARACTER*1
!>           On entry, TRANS  specifies the form of op( A ) to be used
!>           in the matrix multiplication as follows:
!>
!>              TRANS  = 'N' or 'n'   op( A ) = A.
!>
!>              TRANS  = 'T' or 't'   op( A ) = A'.
!>
!>           Unchanged on exit.
!>
[in]DIAG
!>          DIAG is CHARACTER*1
!>           On entry, DIAG specifies whether or not RFP A is unit
!>           triangular as follows:
!>
!>              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
!>
!>              DIAG = 'N' or 'n'   A is not assumed to be unit
!>                                  triangular.
!>
!>           Unchanged on exit.
!>
[in]M
!>          M is INTEGER
!>           On entry, M specifies the number of rows of B. M must be at
!>           least zero.
!>           Unchanged on exit.
!>
[in]N
!>          N is INTEGER
!>           On entry, N specifies the number of columns of B.  N must be
!>           at least zero.
!>           Unchanged on exit.
!>
[in]ALPHA
!>          ALPHA is DOUBLE PRECISION
!>           On entry,  ALPHA specifies the scalar  alpha. When  alpha is
!>           zero then  A is not referenced and  B need not be set before
!>           entry.
!>           Unchanged on exit.
!>
[in]A
!>          A is DOUBLE PRECISION array, dimension (NT)
!>           NT = N*(N+1)/2 if SIDE='R' and NT = M*(M+1)/2 otherwise.
!>           On entry, the matrix A in RFP Format.
!>           RFP Format is described by TRANSR, UPLO and N as follows:
!>           If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even;
!>           K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If
!>           TRANSR = 'T' then RFP is the transpose of RFP A as
!>           defined when TRANSR = 'N'. The contents of RFP A are defined
!>           by UPLO as follows: If UPLO = 'U' the RFP A contains the NT
!>           elements of upper packed A either in normal or
!>           transpose Format. If UPLO = 'L' the RFP A contains
!>           the NT elements of lower packed A either in normal or
!>           transpose Format. The LDA of RFP A is (N+1)/2 when
!>           TRANSR = 'T'. When TRANSR is 'N' the LDA is N+1 when N is
!>           even and is N when is odd.
!>           See the Note below for more details. Unchanged on exit.
!>
[in,out]B
!>          B is DOUBLE PRECISION array, dimension (LDB,N)
!>           Before entry,  the leading  m by n part of the array  B must
!>           contain  the  right-hand  side  matrix  B,  and  on exit  is
!>           overwritten by the solution matrix  X.
!>
[in]LDB
!>          LDB is INTEGER
!>           On entry, LDB specifies the first dimension of B as declared
!>           in  the  calling  (sub)  program.   LDB  must  be  at  least
!>           max( 1, m ).
!>           Unchanged on exit.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  We first consider Rectangular Full Packed (RFP) Format when N is
!>  even. We give an example where N = 6.
!>
!>      AP is Upper             AP is Lower
!>
!>   00 01 02 03 04 05       00
!>      11 12 13 14 15       10 11
!>         22 23 24 25       20 21 22
!>            33 34 35       30 31 32 33
!>               44 45       40 41 42 43 44
!>                  55       50 51 52 53 54 55
!>
!>
!>  Let TRANSR = 'N'. RFP holds AP as follows:
!>  For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
!>  three columns of AP upper. The lower triangle A(4:6,0:2) consists of
!>  the transpose of the first three columns of AP upper.
!>  For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
!>  three columns of AP lower. The upper triangle A(0:2,0:2) consists of
!>  the transpose of the last three columns of AP lower.
!>  This covers the case N even and TRANSR = 'N'.
!>
!>         RFP A                   RFP A
!>
!>        03 04 05                33 43 53
!>        13 14 15                00 44 54
!>        23 24 25                10 11 55
!>        33 34 35                20 21 22
!>        00 44 45                30 31 32
!>        01 11 55                40 41 42
!>        02 12 22                50 51 52
!>
!>  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
!>  transpose of RFP A above. One therefore gets:
!>
!>
!>           RFP A                   RFP A
!>
!>     03 13 23 33 00 01 02    33 00 10 20 30 40 50
!>     04 14 24 34 44 11 12    43 44 11 21 31 41 51
!>     05 15 25 35 45 55 22    53 54 55 22 32 42 52
!>
!>
!>  We then consider Rectangular Full Packed (RFP) Format when N is
!>  odd. We give an example where N = 5.
!>
!>     AP is Upper                 AP is Lower
!>
!>   00 01 02 03 04              00
!>      11 12 13 14              10 11
!>         22 23 24              20 21 22
!>            33 34              30 31 32 33
!>               44              40 41 42 43 44
!>
!>
!>  Let TRANSR = 'N'. RFP holds AP as follows:
!>  For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
!>  three columns of AP upper. The lower triangle A(3:4,0:1) consists of
!>  the transpose of the first two columns of AP upper.
!>  For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
!>  three columns of AP lower. The upper triangle A(0:1,1:2) consists of
!>  the transpose of the last two columns of AP lower.
!>  This covers the case N odd and TRANSR = 'N'.
!>
!>         RFP A                   RFP A
!>
!>        02 03 04                00 33 43
!>        12 13 14                10 11 44
!>        22 23 24                20 21 22
!>        00 33 34                30 31 32
!>        01 11 44                40 41 42
!>
!>  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
!>  transpose of RFP A above. One therefore gets:
!>
!>           RFP A                   RFP A
!>
!>     02 12 22 00 01             00 10 20 30 40 50
!>     03 13 23 33 11             33 11 21 31 41 51
!>     04 14 24 34 44             43 44 22 32 42 52
!>

Definition at line 274 of file dtfsm.f.

277*
278* -- LAPACK computational routine --
279* -- LAPACK is a software package provided by Univ. of Tennessee, --
280* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
281*
282* .. Scalar Arguments ..
283 CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO
284 INTEGER LDB, M, N
285 DOUBLE PRECISION ALPHA
286* ..
287* .. Array Arguments ..
288 DOUBLE PRECISION A( 0: * ), B( 0: LDB-1, 0: * )
289* ..
290*
291* =====================================================================
292*
293* ..
294* .. Parameters ..
295 DOUBLE PRECISION ONE, ZERO
296 parameter( one = 1.0d+0, zero = 0.0d+0 )
297* ..
298* .. Local Scalars ..
299 LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR,
300 $ NOTRANS
301 INTEGER M1, M2, N1, N2, K, INFO, I, J
302* ..
303* .. External Functions ..
304 LOGICAL LSAME
305 EXTERNAL lsame
306* ..
307* .. External Subroutines ..
308 EXTERNAL xerbla, dgemm, dtrsm
309* ..
310* .. Intrinsic Functions ..
311 INTRINSIC max, mod
312* ..
313* .. Executable Statements ..
314*
315* Test the input parameters.
316*
317 info = 0
318 normaltransr = lsame( transr, 'N' )
319 lside = lsame( side, 'L' )
320 lower = lsame( uplo, 'L' )
321 notrans = lsame( trans, 'N' )
322 IF( .NOT.normaltransr .AND. .NOT.lsame( transr, 'T' ) ) THEN
323 info = -1
324 ELSE IF( .NOT.lside .AND. .NOT.lsame( side, 'R' ) ) THEN
325 info = -2
326 ELSE IF( .NOT.lower .AND. .NOT.lsame( uplo, 'U' ) ) THEN
327 info = -3
328 ELSE IF( .NOT.notrans .AND. .NOT.lsame( trans, 'T' ) ) THEN
329 info = -4
330 ELSE IF( .NOT.lsame( diag, 'N' ) .AND.
331 $ .NOT.lsame( diag, 'U' ) )
332 $ THEN
333 info = -5
334 ELSE IF( m.LT.0 ) THEN
335 info = -6
336 ELSE IF( n.LT.0 ) THEN
337 info = -7
338 ELSE IF( ldb.LT.max( 1, m ) ) THEN
339 info = -11
340 END IF
341 IF( info.NE.0 ) THEN
342 CALL xerbla( 'DTFSM ', -info )
343 RETURN
344 END IF
345*
346* Quick return when ( (N.EQ.0).OR.(M.EQ.0) )
347*
348 IF( ( m.EQ.0 ) .OR. ( n.EQ.0 ) )
349 $ RETURN
350*
351* Quick return when ALPHA.EQ.(0D+0)
352*
353 IF( alpha.EQ.zero ) THEN
354 DO 20 j = 0, n - 1
355 DO 10 i = 0, m - 1
356 b( i, j ) = zero
357 10 CONTINUE
358 20 CONTINUE
359 RETURN
360 END IF
361*
362 IF( lside ) THEN
363*
364* SIDE = 'L'
365*
366* A is M-by-M.
367* If M is odd, set NISODD = .TRUE., and M1 and M2.
368* If M is even, NISODD = .FALSE., and M.
369*
370 IF( mod( m, 2 ).EQ.0 ) THEN
371 misodd = .false.
372 k = m / 2
373 ELSE
374 misodd = .true.
375 IF( lower ) THEN
376 m2 = m / 2
377 m1 = m - m2
378 ELSE
379 m1 = m / 2
380 m2 = m - m1
381 END IF
382 END IF
383*
384*
385 IF( misodd ) THEN
386*
387* SIDE = 'L' and N is odd
388*
389 IF( normaltransr ) THEN
390*
391* SIDE = 'L', N is odd, and TRANSR = 'N'
392*
393 IF( lower ) THEN
394*
395* SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'L'
396*
397 IF( notrans ) THEN
398*
399* SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
400* TRANS = 'N'
401*
402 IF( m.EQ.1 ) THEN
403 CALL dtrsm( 'L', 'L', 'N', diag, m1, n,
404 $ alpha,
405 $ a, m, b, ldb )
406 ELSE
407 CALL dtrsm( 'L', 'L', 'N', diag, m1, n,
408 $ alpha,
409 $ a( 0 ), m, b, ldb )
410 CALL dgemm( 'N', 'N', m2, n, m1, -one,
411 $ a( m1 ),
412 $ m, b, ldb, alpha, b( m1, 0 ), ldb )
413 CALL dtrsm( 'L', 'U', 'T', diag, m2, n, one,
414 $ a( m ), m, b( m1, 0 ), ldb )
415 END IF
416*
417 ELSE
418*
419* SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
420* TRANS = 'T'
421*
422 IF( m.EQ.1 ) THEN
423 CALL dtrsm( 'L', 'L', 'T', diag, m1, n,
424 $ alpha,
425 $ a( 0 ), m, b, ldb )
426 ELSE
427 CALL dtrsm( 'L', 'U', 'N', diag, m2, n,
428 $ alpha,
429 $ a( m ), m, b( m1, 0 ), ldb )
430 CALL dgemm( 'T', 'N', m1, n, m2, -one,
431 $ a( m1 ),
432 $ m, b( m1, 0 ), ldb, alpha, b, ldb )
433 CALL dtrsm( 'L', 'L', 'T', diag, m1, n, one,
434 $ a( 0 ), m, b, ldb )
435 END IF
436*
437 END IF
438*
439 ELSE
440*
441* SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'U'
442*
443 IF( .NOT.notrans ) THEN
444*
445* SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
446* TRANS = 'N'
447*
448 CALL dtrsm( 'L', 'L', 'N', diag, m1, n, alpha,
449 $ a( m2 ), m, b, ldb )
450 CALL dgemm( 'T', 'N', m2, n, m1, -one, a( 0 ),
451 $ m,
452 $ b, ldb, alpha, b( m1, 0 ), ldb )
453 CALL dtrsm( 'L', 'U', 'T', diag, m2, n, one,
454 $ a( m1 ), m, b( m1, 0 ), ldb )
455*
456 ELSE
457*
458* SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
459* TRANS = 'T'
460*
461 CALL dtrsm( 'L', 'U', 'N', diag, m2, n, alpha,
462 $ a( m1 ), m, b( m1, 0 ), ldb )
463 CALL dgemm( 'N', 'N', m1, n, m2, -one, a( 0 ),
464 $ m,
465 $ b( m1, 0 ), ldb, alpha, b, ldb )
466 CALL dtrsm( 'L', 'L', 'T', diag, m1, n, one,
467 $ a( m2 ), m, b, ldb )
468*
469 END IF
470*
471 END IF
472*
473 ELSE
474*
475* SIDE = 'L', N is odd, and TRANSR = 'T'
476*
477 IF( lower ) THEN
478*
479* SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'L'
480*
481 IF( notrans ) THEN
482*
483* SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
484* TRANS = 'N'
485*
486 IF( m.EQ.1 ) THEN
487 CALL dtrsm( 'L', 'U', 'T', diag, m1, n,
488 $ alpha,
489 $ a( 0 ), m1, b, ldb )
490 ELSE
491 CALL dtrsm( 'L', 'U', 'T', diag, m1, n,
492 $ alpha,
493 $ a( 0 ), m1, b, ldb )
494 CALL dgemm( 'T', 'N', m2, n, m1, -one,
495 $ a( m1*m1 ), m1, b, ldb, alpha,
496 $ b( m1, 0 ), ldb )
497 CALL dtrsm( 'L', 'L', 'N', diag, m2, n, one,
498 $ a( 1 ), m1, b( m1, 0 ), ldb )
499 END IF
500*
501 ELSE
502*
503* SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and
504* TRANS = 'T'
505*
506 IF( m.EQ.1 ) THEN
507 CALL dtrsm( 'L', 'U', 'N', diag, m1, n,
508 $ alpha,
509 $ a( 0 ), m1, b, ldb )
510 ELSE
511 CALL dtrsm( 'L', 'L', 'T', diag, m2, n,
512 $ alpha,
513 $ a( 1 ), m1, b( m1, 0 ), ldb )
514 CALL dgemm( 'N', 'N', m1, n, m2, -one,
515 $ a( m1*m1 ), m1, b( m1, 0 ), ldb,
516 $ alpha, b, ldb )
517 CALL dtrsm( 'L', 'U', 'N', diag, m1, n, one,
518 $ a( 0 ), m1, b, ldb )
519 END IF
520*
521 END IF
522*
523 ELSE
524*
525* SIDE ='L', N is odd, TRANSR = 'T', and UPLO = 'U'
526*
527 IF( .NOT.notrans ) THEN
528*
529* SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
530* TRANS = 'N'
531*
532 CALL dtrsm( 'L', 'U', 'T', diag, m1, n, alpha,
533 $ a( m2*m2 ), m2, b, ldb )
534 CALL dgemm( 'N', 'N', m2, n, m1, -one, a( 0 ),
535 $ m2,
536 $ b, ldb, alpha, b( m1, 0 ), ldb )
537 CALL dtrsm( 'L', 'L', 'N', diag, m2, n, one,
538 $ a( m1*m2 ), m2, b( m1, 0 ), ldb )
539*
540 ELSE
541*
542* SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'U', and
543* TRANS = 'T'
544*
545 CALL dtrsm( 'L', 'L', 'T', diag, m2, n, alpha,
546 $ a( m1*m2 ), m2, b( m1, 0 ), ldb )
547 CALL dgemm( 'T', 'N', m1, n, m2, -one, a( 0 ),
548 $ m2,
549 $ b( m1, 0 ), ldb, alpha, b, ldb )
550 CALL dtrsm( 'L', 'U', 'N', diag, m1, n, one,
551 $ a( m2*m2 ), m2, b, ldb )
552*
553 END IF
554*
555 END IF
556*
557 END IF
558*
559 ELSE
560*
561* SIDE = 'L' and N is even
562*
563 IF( normaltransr ) THEN
564*
565* SIDE = 'L', N is even, and TRANSR = 'N'
566*
567 IF( lower ) THEN
568*
569* SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'L'
570*
571 IF( notrans ) THEN
572*
573* SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
574* and TRANS = 'N'
575*
576 CALL dtrsm( 'L', 'L', 'N', diag, k, n, alpha,
577 $ a( 1 ), m+1, b, ldb )
578 CALL dgemm( 'N', 'N', k, n, k, -one, a( k+1 ),
579 $ m+1, b, ldb, alpha, b( k, 0 ), ldb )
580 CALL dtrsm( 'L', 'U', 'T', diag, k, n, one,
581 $ a( 0 ), m+1, b( k, 0 ), ldb )
582*
583 ELSE
584*
585* SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
586* and TRANS = 'T'
587*
588 CALL dtrsm( 'L', 'U', 'N', diag, k, n, alpha,
589 $ a( 0 ), m+1, b( k, 0 ), ldb )
590 CALL dgemm( 'T', 'N', k, n, k, -one, a( k+1 ),
591 $ m+1, b( k, 0 ), ldb, alpha, b, ldb )
592 CALL dtrsm( 'L', 'L', 'T', diag, k, n, one,
593 $ a( 1 ), m+1, b, ldb )
594*
595 END IF
596*
597 ELSE
598*
599* SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'U'
600*
601 IF( .NOT.notrans ) THEN
602*
603* SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
604* and TRANS = 'N'
605*
606 CALL dtrsm( 'L', 'L', 'N', diag, k, n, alpha,
607 $ a( k+1 ), m+1, b, ldb )
608 CALL dgemm( 'T', 'N', k, n, k, -one, a( 0 ),
609 $ m+1,
610 $ b, ldb, alpha, b( k, 0 ), ldb )
611 CALL dtrsm( 'L', 'U', 'T', diag, k, n, one,
612 $ a( k ), m+1, b( k, 0 ), ldb )
613*
614 ELSE
615*
616* SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
617* and TRANS = 'T'
618 CALL dtrsm( 'L', 'U', 'N', diag, k, n, alpha,
619 $ a( k ), m+1, b( k, 0 ), ldb )
620 CALL dgemm( 'N', 'N', k, n, k, -one, a( 0 ),
621 $ m+1,
622 $ b( k, 0 ), ldb, alpha, b, ldb )
623 CALL dtrsm( 'L', 'L', 'T', diag, k, n, one,
624 $ a( k+1 ), m+1, b, ldb )
625*
626 END IF
627*
628 END IF
629*
630 ELSE
631*
632* SIDE = 'L', N is even, and TRANSR = 'T'
633*
634 IF( lower ) THEN
635*
636* SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'L'
637*
638 IF( notrans ) THEN
639*
640* SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L',
641* and TRANS = 'N'
642*
643 CALL dtrsm( 'L', 'U', 'T', diag, k, n, alpha,
644 $ a( k ), k, b, ldb )
645 CALL dgemm( 'T', 'N', k, n, k, -one,
646 $ a( k*( k+1 ) ), k, b, ldb, alpha,
647 $ b( k, 0 ), ldb )
648 CALL dtrsm( 'L', 'L', 'N', diag, k, n, one,
649 $ a( 0 ), k, b( k, 0 ), ldb )
650*
651 ELSE
652*
653* SIDE ='L', N is even, TRANSR = 'T', UPLO = 'L',
654* and TRANS = 'T'
655*
656 CALL dtrsm( 'L', 'L', 'T', diag, k, n, alpha,
657 $ a( 0 ), k, b( k, 0 ), ldb )
658 CALL dgemm( 'N', 'N', k, n, k, -one,
659 $ a( k*( k+1 ) ), k, b( k, 0 ), ldb,
660 $ alpha, b, ldb )
661 CALL dtrsm( 'L', 'U', 'N', diag, k, n, one,
662 $ a( k ), k, b, ldb )
663*
664 END IF
665*
666 ELSE
667*
668* SIDE ='L', N is even, TRANSR = 'T', and UPLO = 'U'
669*
670 IF( .NOT.notrans ) THEN
671*
672* SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U',
673* and TRANS = 'N'
674*
675 CALL dtrsm( 'L', 'U', 'T', diag, k, n, alpha,
676 $ a( k*( k+1 ) ), k, b, ldb )
677 CALL dgemm( 'N', 'N', k, n, k, -one, a( 0 ), k,
678 $ b,
679 $ ldb, alpha, b( k, 0 ), ldb )
680 CALL dtrsm( 'L', 'L', 'N', diag, k, n, one,
681 $ a( k*k ), k, b( k, 0 ), ldb )
682*
683 ELSE
684*
685* SIDE ='L', N is even, TRANSR = 'T', UPLO = 'U',
686* and TRANS = 'T'
687*
688 CALL dtrsm( 'L', 'L', 'T', diag, k, n, alpha,
689 $ a( k*k ), k, b( k, 0 ), ldb )
690 CALL dgemm( 'T', 'N', k, n, k, -one, a( 0 ), k,
691 $ b( k, 0 ), ldb, alpha, b, ldb )
692 CALL dtrsm( 'L', 'U', 'N', diag, k, n, one,
693 $ a( k*( k+1 ) ), k, b, ldb )
694*
695 END IF
696*
697 END IF
698*
699 END IF
700*
701 END IF
702*
703 ELSE
704*
705* SIDE = 'R'
706*
707* A is N-by-N.
708* If N is odd, set NISODD = .TRUE., and N1 and N2.
709* If N is even, NISODD = .FALSE., and K.
710*
711 IF( mod( n, 2 ).EQ.0 ) THEN
712 nisodd = .false.
713 k = n / 2
714 ELSE
715 nisodd = .true.
716 IF( lower ) THEN
717 n2 = n / 2
718 n1 = n - n2
719 ELSE
720 n1 = n / 2
721 n2 = n - n1
722 END IF
723 END IF
724*
725 IF( nisodd ) THEN
726*
727* SIDE = 'R' and N is odd
728*
729 IF( normaltransr ) THEN
730*
731* SIDE = 'R', N is odd, and TRANSR = 'N'
732*
733 IF( lower ) THEN
734*
735* SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'L'
736*
737 IF( notrans ) THEN
738*
739* SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
740* TRANS = 'N'
741*
742 CALL dtrsm( 'R', 'U', 'T', diag, m, n2, alpha,
743 $ a( n ), n, b( 0, n1 ), ldb )
744 CALL dgemm( 'N', 'N', m, n1, n2, -one, b( 0,
745 $ n1 ),
746 $ ldb, a( n1 ), n, alpha, b( 0, 0 ),
747 $ ldb )
748 CALL dtrsm( 'R', 'L', 'N', diag, m, n1, one,
749 $ a( 0 ), n, b( 0, 0 ), ldb )
750*
751 ELSE
752*
753* SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
754* TRANS = 'T'
755*
756 CALL dtrsm( 'R', 'L', 'T', diag, m, n1, alpha,
757 $ a( 0 ), n, b( 0, 0 ), ldb )
758 CALL dgemm( 'N', 'T', m, n2, n1, -one, b( 0,
759 $ 0 ),
760 $ ldb, a( n1 ), n, alpha, b( 0, n1 ),
761 $ ldb )
762 CALL dtrsm( 'R', 'U', 'N', diag, m, n2, one,
763 $ a( n ), n, b( 0, n1 ), ldb )
764*
765 END IF
766*
767 ELSE
768*
769* SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'U'
770*
771 IF( notrans ) THEN
772*
773* SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
774* TRANS = 'N'
775*
776 CALL dtrsm( 'R', 'L', 'T', diag, m, n1, alpha,
777 $ a( n2 ), n, b( 0, 0 ), ldb )
778 CALL dgemm( 'N', 'N', m, n2, n1, -one, b( 0,
779 $ 0 ),
780 $ ldb, a( 0 ), n, alpha, b( 0, n1 ),
781 $ ldb )
782 CALL dtrsm( 'R', 'U', 'N', diag, m, n2, one,
783 $ a( n1 ), n, b( 0, n1 ), ldb )
784*
785 ELSE
786*
787* SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
788* TRANS = 'T'
789*
790 CALL dtrsm( 'R', 'U', 'T', diag, m, n2, alpha,
791 $ a( n1 ), n, b( 0, n1 ), ldb )
792 CALL dgemm( 'N', 'T', m, n1, n2, -one, b( 0,
793 $ n1 ),
794 $ ldb, a( 0 ), n, alpha, b( 0, 0 ), ldb )
795 CALL dtrsm( 'R', 'L', 'N', diag, m, n1, one,
796 $ a( n2 ), n, b( 0, 0 ), ldb )
797*
798 END IF
799*
800 END IF
801*
802 ELSE
803*
804* SIDE = 'R', N is odd, and TRANSR = 'T'
805*
806 IF( lower ) THEN
807*
808* SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'L'
809*
810 IF( notrans ) THEN
811*
812* SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
813* TRANS = 'N'
814*
815 CALL dtrsm( 'R', 'L', 'N', diag, m, n2, alpha,
816 $ a( 1 ), n1, b( 0, n1 ), ldb )
817 CALL dgemm( 'N', 'T', m, n1, n2, -one, b( 0,
818 $ n1 ),
819 $ ldb, a( n1*n1 ), n1, alpha, b( 0, 0 ),
820 $ ldb )
821 CALL dtrsm( 'R', 'U', 'T', diag, m, n1, one,
822 $ a( 0 ), n1, b( 0, 0 ), ldb )
823*
824 ELSE
825*
826* SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'L', and
827* TRANS = 'T'
828*
829 CALL dtrsm( 'R', 'U', 'N', diag, m, n1, alpha,
830 $ a( 0 ), n1, b( 0, 0 ), ldb )
831 CALL dgemm( 'N', 'N', m, n2, n1, -one, b( 0,
832 $ 0 ),
833 $ ldb, a( n1*n1 ), n1, alpha, b( 0, n1 ),
834 $ ldb )
835 CALL dtrsm( 'R', 'L', 'T', diag, m, n2, one,
836 $ a( 1 ), n1, b( 0, n1 ), ldb )
837*
838 END IF
839*
840 ELSE
841*
842* SIDE ='R', N is odd, TRANSR = 'T', and UPLO = 'U'
843*
844 IF( notrans ) THEN
845*
846* SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
847* TRANS = 'N'
848*
849 CALL dtrsm( 'R', 'U', 'N', diag, m, n1, alpha,
850 $ a( n2*n2 ), n2, b( 0, 0 ), ldb )
851 CALL dgemm( 'N', 'T', m, n2, n1, -one, b( 0,
852 $ 0 ),
853 $ ldb, a( 0 ), n2, alpha, b( 0, n1 ),
854 $ ldb )
855 CALL dtrsm( 'R', 'L', 'T', diag, m, n2, one,
856 $ a( n1*n2 ), n2, b( 0, n1 ), ldb )
857*
858 ELSE
859*
860* SIDE ='R', N is odd, TRANSR = 'T', UPLO = 'U', and
861* TRANS = 'T'
862*
863 CALL dtrsm( 'R', 'L', 'N', diag, m, n2, alpha,
864 $ a( n1*n2 ), n2, b( 0, n1 ), ldb )
865 CALL dgemm( 'N', 'N', m, n1, n2, -one, b( 0,
866 $ n1 ),
867 $ ldb, a( 0 ), n2, alpha, b( 0, 0 ),
868 $ ldb )
869 CALL dtrsm( 'R', 'U', 'T', diag, m, n1, one,
870 $ a( n2*n2 ), n2, b( 0, 0 ), ldb )
871*
872 END IF
873*
874 END IF
875*
876 END IF
877*
878 ELSE
879*
880* SIDE = 'R' and N is even
881*
882 IF( normaltransr ) THEN
883*
884* SIDE = 'R', N is even, and TRANSR = 'N'
885*
886 IF( lower ) THEN
887*
888* SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'L'
889*
890 IF( notrans ) THEN
891*
892* SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
893* and TRANS = 'N'
894*
895 CALL dtrsm( 'R', 'U', 'T', diag, m, k, alpha,
896 $ a( 0 ), n+1, b( 0, k ), ldb )
897 CALL dgemm( 'N', 'N', m, k, k, -one, b( 0, k ),
898 $ ldb, a( k+1 ), n+1, alpha, b( 0, 0 ),
899 $ ldb )
900 CALL dtrsm( 'R', 'L', 'N', diag, m, k, one,
901 $ a( 1 ), n+1, b( 0, 0 ), ldb )
902*
903 ELSE
904*
905* SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
906* and TRANS = 'T'
907*
908 CALL dtrsm( 'R', 'L', 'T', diag, m, k, alpha,
909 $ a( 1 ), n+1, b( 0, 0 ), ldb )
910 CALL dgemm( 'N', 'T', m, k, k, -one, b( 0, 0 ),
911 $ ldb, a( k+1 ), n+1, alpha, b( 0, k ),
912 $ ldb )
913 CALL dtrsm( 'R', 'U', 'N', diag, m, k, one,
914 $ a( 0 ), n+1, b( 0, k ), ldb )
915*
916 END IF
917*
918 ELSE
919*
920* SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'U'
921*
922 IF( notrans ) THEN
923*
924* SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
925* and TRANS = 'N'
926*
927 CALL dtrsm( 'R', 'L', 'T', diag, m, k, alpha,
928 $ a( k+1 ), n+1, b( 0, 0 ), ldb )
929 CALL dgemm( 'N', 'N', m, k, k, -one, b( 0, 0 ),
930 $ ldb, a( 0 ), n+1, alpha, b( 0, k ),
931 $ ldb )
932 CALL dtrsm( 'R', 'U', 'N', diag, m, k, one,
933 $ a( k ), n+1, b( 0, k ), ldb )
934*
935 ELSE
936*
937* SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
938* and TRANS = 'T'
939*
940 CALL dtrsm( 'R', 'U', 'T', diag, m, k, alpha,
941 $ a( k ), n+1, b( 0, k ), ldb )
942 CALL dgemm( 'N', 'T', m, k, k, -one, b( 0, k ),
943 $ ldb, a( 0 ), n+1, alpha, b( 0, 0 ),
944 $ ldb )
945 CALL dtrsm( 'R', 'L', 'N', diag, m, k, one,
946 $ a( k+1 ), n+1, b( 0, 0 ), ldb )
947*
948 END IF
949*
950 END IF
951*
952 ELSE
953*
954* SIDE = 'R', N is even, and TRANSR = 'T'
955*
956 IF( lower ) THEN
957*
958* SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'L'
959*
960 IF( notrans ) THEN
961*
962* SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L',
963* and TRANS = 'N'
964*
965 CALL dtrsm( 'R', 'L', 'N', diag, m, k, alpha,
966 $ a( 0 ), k, b( 0, k ), ldb )
967 CALL dgemm( 'N', 'T', m, k, k, -one, b( 0, k ),
968 $ ldb, a( ( k+1 )*k ), k, alpha,
969 $ b( 0, 0 ), ldb )
970 CALL dtrsm( 'R', 'U', 'T', diag, m, k, one,
971 $ a( k ), k, b( 0, 0 ), ldb )
972*
973 ELSE
974*
975* SIDE ='R', N is even, TRANSR = 'T', UPLO = 'L',
976* and TRANS = 'T'
977*
978 CALL dtrsm( 'R', 'U', 'N', diag, m, k, alpha,
979 $ a( k ), k, b( 0, 0 ), ldb )
980 CALL dgemm( 'N', 'N', m, k, k, -one, b( 0, 0 ),
981 $ ldb, a( ( k+1 )*k ), k, alpha,
982 $ b( 0, k ), ldb )
983 CALL dtrsm( 'R', 'L', 'T', diag, m, k, one,
984 $ a( 0 ), k, b( 0, k ), ldb )
985*
986 END IF
987*
988 ELSE
989*
990* SIDE ='R', N is even, TRANSR = 'T', and UPLO = 'U'
991*
992 IF( notrans ) THEN
993*
994* SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U',
995* and TRANS = 'N'
996*
997 CALL dtrsm( 'R', 'U', 'N', diag, m, k, alpha,
998 $ a( ( k+1 )*k ), k, b( 0, 0 ), ldb )
999 CALL dgemm( 'N', 'T', m, k, k, -one, b( 0, 0 ),
1000 $ ldb, a( 0 ), k, alpha, b( 0, k ), ldb )
1001 CALL dtrsm( 'R', 'L', 'T', diag, m, k, one,
1002 $ a( k*k ), k, b( 0, k ), ldb )
1003*
1004 ELSE
1005*
1006* SIDE ='R', N is even, TRANSR = 'T', UPLO = 'U',
1007* and TRANS = 'T'
1008*
1009 CALL dtrsm( 'R', 'L', 'N', diag, m, k, alpha,
1010 $ a( k*k ), k, b( 0, k ), ldb )
1011 CALL dgemm( 'N', 'N', m, k, k, -one, b( 0, k ),
1012 $ ldb, a( 0 ), k, alpha, b( 0, 0 ), ldb )
1013 CALL dtrsm( 'R', 'U', 'T', diag, m, k, one,
1014 $ a( ( k+1 )*k ), k, b( 0, 0 ), ldb )
1015*
1016 END IF
1017*
1018 END IF
1019*
1020 END IF
1021*
1022 END IF
1023 END IF
1024*
1025 RETURN
1026*
1027* End of DTFSM
1028*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
dgemm
subroutine dgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
DGEMM
Definition dgemm.f:188
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
dtrsm
subroutine dtrsm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
DTRSM
Definition dtrsm.f:181
Here is the call graph for this function:
Here is the caller graph for this function: