1 *> \brief \b ZLAGTM performs a matrix-matrix product of the form C = αAB+βC, where A is a tridiagonal matrix, B and C are rectangular matrices, and α and β are scalars, which may be 0, 1, or -1.
3 * =========== DOCUMENTATION ===========
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
9 *> Download ZLAGTM + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlagtm.f">
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlagtm.f">
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlagtm.f">
21 * SUBROUTINE ZLAGTM( TRANS, N, NRHS, ALPHA, DL, D, DU, X, LDX, BETA,
24 * .. Scalar Arguments ..
26 * INTEGER LDB, LDX, N, NRHS
27 * DOUBLE PRECISION ALPHA, BETA
29 * .. Array Arguments ..
30 * COMPLEX*16 B( LDB, * ), D( * ), DL( * ), DU( * ),
40 *> ZLAGTM performs a matrix-vector product of the form
42 *> B := alpha * A * X + beta * B
44 *> where A is a tridiagonal matrix of order N, B and X are N by NRHS
45 *> matrices, and alpha and beta are real scalars, each of which may be
54 *> TRANS is CHARACTER*1
55 *> Specifies the operation applied to A.
56 *> = 'N': No transpose, B := alpha * A * X + beta * B
57 *> = 'T': Transpose, B := alpha * A**T * X + beta * B
58 *> = 'C': Conjugate transpose, B := alpha * A**H * X + beta * B
64 *> The order of the matrix A. N >= 0.
70 *> The number of right hand sides, i.e., the number of columns
71 *> of the matrices X and B.
76 *> ALPHA is DOUBLE PRECISION
77 *> The scalar alpha. ALPHA must be 0., 1., or -1.; otherwise,
78 *> it is assumed to be 0.
83 *> DL is COMPLEX*16 array, dimension (N-1)
84 *> The (n-1) sub-diagonal elements of T.
89 *> D is COMPLEX*16 array, dimension (N)
90 *> The diagonal elements of T.
95 *> DU is COMPLEX*16 array, dimension (N-1)
96 *> The (n-1) super-diagonal elements of T.
101 *> X is COMPLEX*16 array, dimension (LDX,NRHS)
102 *> The N by NRHS matrix X.
108 *> The leading dimension of the array X. LDX >= max(N,1).
113 *> BETA is DOUBLE PRECISION
114 *> The scalar beta. BETA must be 0., 1., or -1.; otherwise,
115 *> it is assumed to be 1.
120 *> B is COMPLEX*16 array, dimension (LDB,NRHS)
121 *> On entry, the N by NRHS matrix B.
122 *> On exit, B is overwritten by the matrix expression
123 *> B := alpha * A * X + beta * B.
129 *> The leading dimension of the array B. LDB >= max(N,1).
135 *> \author Univ. of Tennessee
136 *> \author Univ. of California Berkeley
137 *> \author Univ. of Colorado Denver
140 *> \date September 2012
142 *> \ingroup complex16OTHERauxiliary
144 * =====================================================================
145 SUBROUTINE ZLAGTM( TRANS, N, NRHS, ALPHA, DL, D, DU, X, LDX, BETA,
148 * -- LAPACK auxiliary routine (version 3.4.2) --
149 * -- LAPACK is a software package provided by Univ. of Tennessee, --
150 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
153 * .. Scalar Arguments ..
155 INTEGER LDB, LDX, N, NRHS
156 DOUBLE PRECISION ALPHA, BETA
158 * .. Array Arguments ..
159 COMPLEX*16 B( LDB, * ), D( * ), DL( * ), DU( * ),
163 * =====================================================================
166 DOUBLE PRECISION ONE, ZERO
167 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
169 * .. Local Scalars ..
172 * .. External Functions ..
176 * .. Intrinsic Functions ..
179 * .. Executable Statements ..
184 * Multiply B by BETA if BETA.NE.1.
186 IF( BETA.EQ.ZERO ) THEN
192 ELSE IF( BETA.EQ.-ONE ) THEN
195 B( I, J ) = -B( I, J )
200 IF( ALPHA.EQ.ONE ) THEN
201 IF( LSAME( TRANS, 'N' ) ) THEN
203 * Compute B := B + A*X
207 B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J )
209 B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J ) +
211 B( N, J ) = B( N, J ) + DL( N-1 )*X( N-1, J ) +
214 B( I, J ) = B( I, J ) + DL( I-1 )*X( I-1, J ) +
215 $ D( I )*X( I, J ) + DU( I )*X( I+1, J )
219 ELSE IF( LSAME( TRANS, 'T' ) ) THEN
221 * Compute B := B + A**T * X
225 B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J )
227 B( 1, J ) = B( 1, J ) + D( 1 )*X( 1, J ) +
229 B( N, J ) = B( N, J ) + DU( N-1 )*X( N-1, J ) +
232 B( I, J ) = B( I, J ) + DU( I-1 )*X( I-1, J ) +
233 $ D( I )*X( I, J ) + DL( I )*X( I+1, J )
237 ELSE IF( LSAME( TRANS, 'C' ) ) THEN
239 * Compute B := B + A**H * X
243 B( 1, J ) = B( 1, J ) + DCONJG( D( 1 ) )*X( 1, J )
245 B( 1, J ) = B( 1, J ) + DCONJG( D( 1 ) )*X( 1, J ) +
246 $ DCONJG( DL( 1 ) )*X( 2, J )
247 B( N, J ) = B( N, J ) + DCONJG( DU( N-1 ) )*
248 $ X( N-1, J ) + DCONJG( D( N ) )*X( N, J )
250 B( I, J ) = B( I, J ) + DCONJG( DU( I-1 ) )*
251 $ X( I-1, J ) + DCONJG( D( I ) )*
252 $ X( I, J ) + DCONJG( DL( I ) )*
258 ELSE IF( ALPHA.EQ.-ONE ) THEN
259 IF( LSAME( TRANS, 'N' ) ) THEN
261 * Compute B := B - A*X
265 B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J )
267 B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J ) -
269 B( N, J ) = B( N, J ) - DL( N-1 )*X( N-1, J ) -
272 B( I, J ) = B( I, J ) - DL( I-1 )*X( I-1, J ) -
273 $ D( I )*X( I, J ) - DU( I )*X( I+1, J )
277 ELSE IF( LSAME( TRANS, 'T' ) ) THEN
279 * Compute B := B - A**T *X
283 B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J )
285 B( 1, J ) = B( 1, J ) - D( 1 )*X( 1, J ) -
287 B( N, J ) = B( N, J ) - DU( N-1 )*X( N-1, J ) -
290 B( I, J ) = B( I, J ) - DU( I-1 )*X( I-1, J ) -
291 $ D( I )*X( I, J ) - DL( I )*X( I+1, J )
295 ELSE IF( LSAME( TRANS, 'C' ) ) THEN
297 * Compute B := B - A**H *X
301 B( 1, J ) = B( 1, J ) - DCONJG( D( 1 ) )*X( 1, J )
303 B( 1, J ) = B( 1, J ) - DCONJG( D( 1 ) )*X( 1, J ) -
304 $ DCONJG( DL( 1 ) )*X( 2, J )
305 B( N, J ) = B( N, J ) - DCONJG( DU( N-1 ) )*
306 $ X( N-1, J ) - DCONJG( D( N ) )*X( N, J )
308 B( I, J ) = B( I, J ) - DCONJG( DU( I-1 ) )*
309 $ X( I-1, J ) - DCONJG( D( I ) )*
310 $ X( I, J ) - DCONJG( DL( I ) )*