File:  [local] / rpl / lapack / lapack / dlagtm.f
Revision 1.3: download - view: text, annotated - select for diffs - revision graph
Fri Aug 6 15:28:39 2010 UTC (13 years, 9 months ago) by bertrand
Branches: MAIN
CVS tags: HEAD
Cohérence

    1:       SUBROUTINE DLAGTM( 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:       DOUBLE PRECISION   B( LDB, * ), D( * ), DL( * ), DU( * ),
   16:      $                   X( LDX, * )
   17: *     ..
   18: *
   19: *  Purpose
   20: *  =======
   21: *
   22: *  DLAGTM 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'* X + beta * B
   37: *          = 'C':  Conjugate transpose = Transpose
   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) DOUBLE PRECISION array, dimension (N-1)
   51: *          The (n-1) sub-diagonal elements of T.
   52: *
   53: *  D       (input) DOUBLE PRECISION array, dimension (N)
   54: *          The diagonal elements of T.
   55: *
   56: *  DU      (input) DOUBLE PRECISION array, dimension (N-1)
   57: *          The (n-1) super-diagonal elements of T.
   58: *
   59: *  X       (input) DOUBLE PRECISION 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) DOUBLE PRECISION 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: *     .. Executable Statements ..
   90: *
   91:       IF( N.EQ.0 )
   92:      $   RETURN
   93: *
   94: *     Multiply B by BETA if BETA.NE.1.
   95: *
   96:       IF( BETA.EQ.ZERO ) THEN
   97:          DO 20 J = 1, NRHS
   98:             DO 10 I = 1, N
   99:                B( I, J ) = ZERO
  100:    10       CONTINUE
  101:    20    CONTINUE
  102:       ELSE IF( BETA.EQ.-ONE ) THEN
  103:          DO 40 J = 1, NRHS
  104:             DO 30 I = 1, N
  105:                B( I, J ) = -B( I, J )
  106:    30       CONTINUE
  107:    40    CONTINUE
  108:       END IF
  109: *
  110:       IF( ALPHA.EQ.ONE ) THEN
  111:          IF( LSAME( TRANS, 'N' ) ) THEN
  112: *
  113: *           Compute B := B + A*X
  114: *
  115:             DO 60 J = 1, NRHS
  116:                IF( N.EQ.1 ) THEN
  117:                   B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J )
  118:                ELSE
  119:                   B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J ) +
  120:      $                        DU( 1 )*X( 2, J )
  121:                   B( N, J ) = B( N, J ) + DL( N-1 )*X( N-1, J ) +
  122:      $                        D( N )*X( N, J )
  123:                   DO 50 I = 2, N - 1
  124:                      B( I, J ) = B( I, J ) + DL( I-1 )*X( I-1, J ) +
  125:      $                           D( I )*X( I, J ) + DU( I )*X( I+1, J )
  126:    50             CONTINUE
  127:                END IF
  128:    60       CONTINUE
  129:          ELSE
  130: *
  131: *           Compute B := B + A'*X
  132: *
  133:             DO 80 J = 1, NRHS
  134:                IF( N.EQ.1 ) THEN
  135:                   B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J )
  136:                ELSE
  137:                   B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J ) +
  138:      $                        DL( 1 )*X( 2, J )
  139:                   B( N, J ) = B( N, J ) + DU( N-1 )*X( N-1, J ) +
  140:      $                        D( N )*X( N, J )
  141:                   DO 70 I = 2, N - 1
  142:                      B( I, J ) = B( I, J ) + DU( I-1 )*X( I-1, J ) +
  143:      $                           D( I )*X( I, J ) + DL( I )*X( I+1, J )
  144:    70             CONTINUE
  145:                END IF
  146:    80       CONTINUE
  147:          END IF
  148:       ELSE IF( ALPHA.EQ.-ONE ) THEN
  149:          IF( LSAME( TRANS, 'N' ) ) THEN
  150: *
  151: *           Compute B := B - A*X
  152: *
  153:             DO 100 J = 1, NRHS
  154:                IF( N.EQ.1 ) THEN
  155:                   B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J )
  156:                ELSE
  157:                   B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J ) -
  158:      $                        DU( 1 )*X( 2, J )
  159:                   B( N, J ) = B( N, J ) - DL( N-1 )*X( N-1, J ) -
  160:      $                        D( N )*X( N, J )
  161:                   DO 90 I = 2, N - 1
  162:                      B( I, J ) = B( I, J ) - DL( I-1 )*X( I-1, J ) -
  163:      $                           D( I )*X( I, J ) - DU( I )*X( I+1, J )
  164:    90             CONTINUE
  165:                END IF
  166:   100       CONTINUE
  167:          ELSE
  168: *
  169: *           Compute B := B - A'*X
  170: *
  171:             DO 120 J = 1, NRHS
  172:                IF( N.EQ.1 ) THEN
  173:                   B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J )
  174:                ELSE
  175:                   B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J ) -
  176:      $                        DL( 1 )*X( 2, J )
  177:                   B( N, J ) = B( N, J ) - DU( N-1 )*X( N-1, J ) -
  178:      $                        D( N )*X( N, J )
  179:                   DO 110 I = 2, N - 1
  180:                      B( I, J ) = B( I, J ) - DU( I-1 )*X( I-1, J ) -
  181:      $                           D( I )*X( I, J ) - DL( I )*X( I+1, J )
  182:   110             CONTINUE
  183:                END IF
  184:   120       CONTINUE
  185:          END IF
  186:       END IF
  187:       RETURN
  188: *
  189: *     End of DLAGTM
  190: *
  191:       END

CVSweb interface <joel.bertrand@systella.fr>