File:  [local] / rpl / lapack / lapack / zlagtm.f
Revision 1.6: download - view: text, annotated - select for diffs - revision graph
Fri Aug 13 21:04:08 2010 UTC (13 years, 9 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_0_19, rpl-4_0_18, HEAD
Patches pour OS/2

    1:       SUBROUTINE ZLAGTM( TRANS, N, NRHS, ALPHA, DL, D, DU, X, LDX, BETA,
    2:      $                   B, LDB )
    3: *
    4: *  -- LAPACK auxiliary routine (version 3.2) --
    5: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
    6: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
    7: *     November 2006
    8: *
    9: *     .. Scalar Arguments ..
   10:       CHARACTER          TRANS
   11:       INTEGER            LDB, LDX, N, NRHS
   12:       DOUBLE PRECISION   ALPHA, BETA
   13: *     ..
   14: *     .. Array Arguments ..
   15:       COMPLEX*16         B( LDB, * ), D( * ), DL( * ), DU( * ),
   16:      $                   X( LDX, * )
   17: *     ..
   18: *
   19: *  Purpose
   20: *  =======
   21: *
   22: *  ZLAGTM performs a matrix-vector product of the form
   23: *
   24: *     B := alpha * A * X + beta * B
   25: *
   26: *  where A is a tridiagonal matrix of order N, B and X are N by NRHS
   27: *  matrices, and alpha and beta are real scalars, each of which may be
   28: *  0., 1., or -1.
   29: *
   30: *  Arguments
   31: *  =========
   32: *
   33: *  TRANS   (input) CHARACTER*1
   34: *          Specifies the operation applied to A.
   35: *          = 'N':  No transpose, B := alpha * A * X + beta * B
   36: *          = 'T':  Transpose,    B := alpha * A**T * X + beta * B
   37: *          = 'C':  Conjugate transpose, B := alpha * A**H * X + beta * B
   38: *
   39: *  N       (input) INTEGER
   40: *          The order of the matrix A.  N >= 0.
   41: *
   42: *  NRHS    (input) INTEGER
   43: *          The number of right hand sides, i.e., the number of columns
   44: *          of the matrices X and B.
   45: *
   46: *  ALPHA   (input) DOUBLE PRECISION
   47: *          The scalar alpha.  ALPHA must be 0., 1., or -1.; otherwise,
   48: *          it is assumed to be 0.
   49: *
   50: *  DL      (input) COMPLEX*16 array, dimension (N-1)
   51: *          The (n-1) sub-diagonal elements of T.
   52: *
   53: *  D       (input) COMPLEX*16 array, dimension (N)
   54: *          The diagonal elements of T.
   55: *
   56: *  DU      (input) COMPLEX*16 array, dimension (N-1)
   57: *          The (n-1) super-diagonal elements of T.
   58: *
   59: *  X       (input) COMPLEX*16 array, dimension (LDX,NRHS)
   60: *          The N by NRHS matrix X.
   61: *  LDX     (input) INTEGER
   62: *          The leading dimension of the array X.  LDX >= max(N,1).
   63: *
   64: *  BETA    (input) DOUBLE PRECISION
   65: *          The scalar beta.  BETA must be 0., 1., or -1.; otherwise,
   66: *          it is assumed to be 1.
   67: *
   68: *  B       (input/output) COMPLEX*16 array, dimension (LDB,NRHS)
   69: *          On entry, the N by NRHS matrix B.
   70: *          On exit, B is overwritten by the matrix expression
   71: *          B := alpha * A * X + beta * B.
   72: *
   73: *  LDB     (input) INTEGER
   74: *          The leading dimension of the array B.  LDB >= max(N,1).
   75: *
   76: *  =====================================================================
   77: *
   78: *     .. Parameters ..
   79:       DOUBLE PRECISION   ONE, ZERO
   80:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
   81: *     ..
   82: *     .. Local Scalars ..
   83:       INTEGER            I, J
   84: *     ..
   85: *     .. External Functions ..
   86:       LOGICAL            LSAME
   87:       EXTERNAL           LSAME
   88: *     ..
   89: *     .. Intrinsic Functions ..
   90:       INTRINSIC          DCONJG
   91: *     ..
   92: *     .. Executable Statements ..
   93: *
   94:       IF( N.EQ.0 )
   95:      $   RETURN
   96: *
   97: *     Multiply B by BETA if BETA.NE.1.
   98: *
   99:       IF( BETA.EQ.ZERO ) THEN
  100:          DO 20 J = 1, NRHS
  101:             DO 10 I = 1, N
  102:                B( I, J ) = ZERO
  103:    10       CONTINUE
  104:    20    CONTINUE
  105:       ELSE IF( BETA.EQ.-ONE ) THEN
  106:          DO 40 J = 1, NRHS
  107:             DO 30 I = 1, N
  108:                B( I, J ) = -B( I, J )
  109:    30       CONTINUE
  110:    40    CONTINUE
  111:       END IF
  112: *
  113:       IF( ALPHA.EQ.ONE ) THEN
  114:          IF( LSAME( TRANS, 'N' ) ) THEN
  115: *
  116: *           Compute B := B + A*X
  117: *
  118:             DO 60 J = 1, NRHS
  119:                IF( N.EQ.1 ) THEN
  120:                   B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J )
  121:                ELSE
  122:                   B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J ) +
  123:      $                        DU( 1 )*X( 2, J )
  124:                   B( N, J ) = B( N, J ) + DL( N-1 )*X( N-1, J ) +
  125:      $                        D( N )*X( N, J )
  126:                   DO 50 I = 2, N - 1
  127:                      B( I, J ) = B( I, J ) + DL( I-1 )*X( I-1, J ) +
  128:      $                           D( I )*X( I, J ) + DU( I )*X( I+1, J )
  129:    50             CONTINUE
  130:                END IF
  131:    60       CONTINUE
  132:          ELSE IF( LSAME( TRANS, 'T' ) ) THEN
  133: *
  134: *           Compute B := B + A**T * X
  135: *
  136:             DO 80 J = 1, NRHS
  137:                IF( N.EQ.1 ) THEN
  138:                   B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J )
  139:                ELSE
  140:                   B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J ) +
  141:      $                        DL( 1 )*X( 2, J )
  142:                   B( N, J ) = B( N, J ) + DU( N-1 )*X( N-1, J ) +
  143:      $                        D( N )*X( N, J )
  144:                   DO 70 I = 2, N - 1
  145:                      B( I, J ) = B( I, J ) + DU( I-1 )*X( I-1, J ) +
  146:      $                           D( I )*X( I, J ) + DL( I )*X( I+1, J )
  147:    70             CONTINUE
  148:                END IF
  149:    80       CONTINUE
  150:          ELSE IF( LSAME( TRANS, 'C' ) ) THEN
  151: *
  152: *           Compute B := B + A**H * X
  153: *
  154:             DO 100 J = 1, NRHS
  155:                IF( N.EQ.1 ) THEN
  156:                   B( 1, J ) = B( 1, J ) + DCONJG( D( 1 ) )*X( 1, J )
  157:                ELSE
  158:                   B( 1, J ) = B( 1, J ) + DCONJG( D( 1 ) )*X( 1, J ) +
  159:      $                        DCONJG( DL( 1 ) )*X( 2, J )
  160:                   B( N, J ) = B( N, J ) + DCONJG( DU( N-1 ) )*
  161:      $                        X( N-1, J ) + DCONJG( D( N ) )*X( N, J )
  162:                   DO 90 I = 2, N - 1
  163:                      B( I, J ) = B( I, J ) + DCONJG( DU( I-1 ) )*
  164:      $                           X( I-1, J ) + DCONJG( D( I ) )*
  165:      $                           X( I, J ) + DCONJG( DL( I ) )*
  166:      $                           X( I+1, J )
  167:    90             CONTINUE
  168:                END IF
  169:   100       CONTINUE
  170:          END IF
  171:       ELSE IF( ALPHA.EQ.-ONE ) THEN
  172:          IF( LSAME( TRANS, 'N' ) ) THEN
  173: *
  174: *           Compute B := B - A*X
  175: *
  176:             DO 120 J = 1, NRHS
  177:                IF( N.EQ.1 ) THEN
  178:                   B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J )
  179:                ELSE
  180:                   B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J ) -
  181:      $                        DU( 1 )*X( 2, J )
  182:                   B( N, J ) = B( N, J ) - DL( N-1 )*X( N-1, J ) -
  183:      $                        D( N )*X( N, J )
  184:                   DO 110 I = 2, N - 1
  185:                      B( I, J ) = B( I, J ) - DL( I-1 )*X( I-1, J ) -
  186:      $                           D( I )*X( I, J ) - DU( I )*X( I+1, J )
  187:   110             CONTINUE
  188:                END IF
  189:   120       CONTINUE
  190:          ELSE IF( LSAME( TRANS, 'T' ) ) THEN
  191: *
  192: *           Compute B := B - A'*X
  193: *
  194:             DO 140 J = 1, NRHS
  195:                IF( N.EQ.1 ) THEN
  196:                   B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J )
  197:                ELSE
  198:                   B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J ) -
  199:      $                        DL( 1 )*X( 2, J )
  200:                   B( N, J ) = B( N, J ) - DU( N-1 )*X( N-1, J ) -
  201:      $                        D( N )*X( N, J )
  202:                   DO 130 I = 2, N - 1
  203:                      B( I, J ) = B( I, J ) - DU( I-1 )*X( I-1, J ) -
  204:      $                           D( I )*X( I, J ) - DL( I )*X( I+1, J )
  205:   130             CONTINUE
  206:                END IF
  207:   140       CONTINUE
  208:          ELSE IF( LSAME( TRANS, 'C' ) ) THEN
  209: *
  210: *           Compute B := B - A'*X
  211: *
  212:             DO 160 J = 1, NRHS
  213:                IF( N.EQ.1 ) THEN
  214:                   B( 1, J ) = B( 1, J ) - DCONJG( D( 1 ) )*X( 1, J )
  215:                ELSE
  216:                   B( 1, J ) = B( 1, J ) - DCONJG( D( 1 ) )*X( 1, J ) -
  217:      $                        DCONJG( DL( 1 ) )*X( 2, J )
  218:                   B( N, J ) = B( N, J ) - DCONJG( DU( N-1 ) )*
  219:      $                        X( N-1, J ) - DCONJG( D( N ) )*X( N, J )
  220:                   DO 150 I = 2, N - 1
  221:                      B( I, J ) = B( I, J ) - DCONJG( DU( I-1 ) )*
  222:      $                           X( I-1, J ) - DCONJG( D( I ) )*
  223:      $                           X( I, J ) - DCONJG( DL( I ) )*
  224:      $                           X( I+1, J )
  225:   150             CONTINUE
  226:                END IF
  227:   160       CONTINUE
  228:          END IF
  229:       END IF
  230:       RETURN
  231: *
  232: *     End of ZLAGTM
  233: *
  234:       END

CVSweb interface <joel.bertrand@systella.fr>