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>