1
SUBROUTINE ZTRSV ( UPLO, TRANS, DIAG, N, A, LDA, X, INCX )
3
* $Id: ztrsv.f 19695 2010-10-29 16:51:02Z d3y133 $
3
* =========== DOCUMENTATION ===========
5
* Online html documentation available at
6
* http://www.netlib.org/lapack/explore-html/
11
* SUBROUTINE ZTRSV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
13
* .. Scalar Arguments ..
15
* CHARACTER DIAG,TRANS,UPLO
17
* .. Array Arguments ..
18
* COMPLEX*16 A(LDA,*),X(*)
27
*> ZTRSV solves one of the systems of equations
29
*> A*x = b, or A**T*x = b, or A**H*x = b,
31
*> where b and x are n element vectors and A is an n by n unit, or
32
*> non-unit, upper or lower triangular matrix.
34
*> No test for singularity or near-singularity is included in this
35
*> routine. Such tests must be performed before calling this routine.
43
*> UPLO is CHARACTER*1
44
*> On entry, UPLO specifies whether the matrix is an upper or
45
*> lower triangular matrix as follows:
47
*> UPLO = 'U' or 'u' A is an upper triangular matrix.
49
*> UPLO = 'L' or 'l' A is a lower triangular matrix.
54
*> TRANS is CHARACTER*1
55
*> On entry, TRANS specifies the equations to be solved as
58
*> TRANS = 'N' or 'n' A*x = b.
60
*> TRANS = 'T' or 't' A**T*x = b.
62
*> TRANS = 'C' or 'c' A**H*x = b.
67
*> DIAG is CHARACTER*1
68
*> On entry, DIAG specifies whether or not A is unit
69
*> triangular as follows:
71
*> DIAG = 'U' or 'u' A is assumed to be unit triangular.
73
*> DIAG = 'N' or 'n' A is not assumed to be unit
80
*> On entry, N specifies the order of the matrix A.
81
*> N must be at least zero.
86
*> A is COMPLEX*16 array of DIMENSION ( LDA, n ).
87
*> Before entry with UPLO = 'U' or 'u', the leading n by n
88
*> upper triangular part of the array A must contain the upper
89
*> triangular matrix and the strictly lower triangular part of
90
*> A is not referenced.
91
*> Before entry with UPLO = 'L' or 'l', the leading n by n
92
*> lower triangular part of the array A must contain the lower
93
*> triangular matrix and the strictly upper triangular part of
94
*> A is not referenced.
95
*> Note that when DIAG = 'U' or 'u', the diagonal elements of
96
*> A are not referenced either, but are assumed to be unity.
102
*> On entry, LDA specifies the first dimension of A as declared
103
*> in the calling (sub) program. LDA must be at least
109
*> X is COMPLEX*16 array of dimension at least
110
*> ( 1 + ( n - 1 )*abs( INCX ) ).
111
*> Before entry, the incremented array X must contain the n
112
*> element right-hand side vector b. On exit, X is overwritten
113
*> with the solution vector x.
119
*> On entry, INCX specifies the increment for the elements of
120
*> X. INCX must not be zero.
126
*> \author Univ. of Tennessee
127
*> \author Univ. of California Berkeley
128
*> \author Univ. of Colorado Denver
131
*> \date November 2011
133
*> \ingroup complex16_blas_level2
135
*> \par Further Details:
136
* =====================
140
*> Level 2 Blas routine.
142
*> -- Written on 22-October-1986.
143
*> Jack Dongarra, Argonne National Lab.
144
*> Jeremy Du Croz, Nag Central Office.
145
*> Sven Hammarling, Nag Central Office.
146
*> Richard Hanson, Sandia National Labs.
149
* =====================================================================
150
SUBROUTINE ZTRSV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
152
* -- Reference BLAS level2 routine (version 3.4.0) --
153
* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
154
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
5
157
* .. Scalar Arguments ..
7
CHARACTER*1 DIAG, TRANS, UPLO
159
CHARACTER DIAG,TRANS,UPLO
8
161
* .. Array Arguments ..
9
COMPLEX*16 A( LDA, * ), X( * )
162
COMPLEX*16 A(LDA,*),X(*)
15
* ZTRSV solves one of the systems of equations
17
* A*x = b, or A'*x = b, or conjg( A' )*x = b,
19
* where b and x are n element vectors and A is an n by n unit, or
20
* non-unit, upper or lower triangular matrix.
22
* No test for singularity or near-singularity is included in this
23
* routine. Such tests must be performed before calling this routine.
29
* On entry, UPLO specifies whether the matrix is an upper or
30
* lower triangular matrix as follows:
32
* UPLO = 'U' or 'u' A is an upper triangular matrix.
34
* UPLO = 'L' or 'l' A is a lower triangular matrix.
38
* TRANS - CHARACTER*1.
39
* On entry, TRANS specifies the equations to be solved as
42
* TRANS = 'N' or 'n' A*x = b.
44
* TRANS = 'T' or 't' A'*x = b.
46
* TRANS = 'C' or 'c' conjg( A' )*x = b.
51
* On entry, DIAG specifies whether or not A is unit
52
* triangular as follows:
54
* DIAG = 'U' or 'u' A is assumed to be unit triangular.
56
* DIAG = 'N' or 'n' A is not assumed to be unit
62
* On entry, N specifies the order of the matrix A.
63
* N must be at least zero.
66
* A - COMPLEX*16 array of DIMENSION ( LDA, n ).
67
* Before entry with UPLO = 'U' or 'u', the leading n by n
68
* upper triangular part of the array A must contain the upper
69
* triangular matrix and the strictly lower triangular part of
70
* A is not referenced.
71
* Before entry with UPLO = 'L' or 'l', the leading n by n
72
* lower triangular part of the array A must contain the lower
73
* triangular matrix and the strictly upper triangular part of
74
* A is not referenced.
75
* Note that when DIAG = 'U' or 'u', the diagonal elements of
76
* A are not referenced either, but are assumed to be unity.
80
* On entry, LDA specifies the first dimension of A as declared
81
* in the calling (sub) program. LDA must be at least
85
* X - COMPLEX*16 array of dimension at least
86
* ( 1 + ( n - 1 )*abs( INCX ) ).
87
* Before entry, the incremented array X must contain the n
88
* element right-hand side vector b. On exit, X is overwritten
89
* with the solution vector x.
92
* On entry, INCX specifies the increment for the elements of
93
* X. INCX must not be zero.
97
* Level 2 Blas routine.
99
* -- Written on 22-October-1986.
100
* Jack Dongarra, Argonne National Lab.
101
* Jeremy Du Croz, Nag Central Office.
102
* Sven Hammarling, Nag Central Office.
103
* Richard Hanson, Sandia National Labs.
165
* =====================================================================
106
167
* .. Parameters ..
108
PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ) )
169
PARAMETER (ZERO= (0.0D+0,0.0D+0))
109
171
* .. Local Scalars ..
111
INTEGER I, INFO, IX, J, JX, KX
112
LOGICAL NOCONJ, NOUNIT
173
INTEGER I,INFO,IX,J,JX,KX
174
LOGICAL NOCONJ,NOUNIT
113
176
* .. External Functions ..
116
180
* .. External Subroutines ..
118
183
* .. Intrinsic Functions ..
119
INTRINSIC DCONJG, MAX
121
* .. Executable Statements ..
123
187
* Test the input parameters.
126
IF ( .NOT.LSAME( UPLO , 'U' ).AND.
127
$ .NOT.LSAME( UPLO , 'L' ) )THEN
129
ELSE IF( .NOT.LSAME( TRANS, 'N' ).AND.
130
$ .NOT.LSAME( TRANS, 'T' ).AND.
131
$ .NOT.LSAME( TRANS, 'C' ) )THEN
133
ELSE IF( .NOT.LSAME( DIAG , 'U' ).AND.
134
$ .NOT.LSAME( DIAG , 'N' ) )THEN
136
ELSE IF( N.LT.0 )THEN
138
ELSE IF( LDA.LT.MAX( 1, N ) )THEN
140
ELSE IF( INCX.EQ.0 )THEN
190
IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
192
ELSE IF (.NOT.LSAME(TRANS,'N') .AND. .NOT.LSAME(TRANS,'T') .AND.
193
+ .NOT.LSAME(TRANS,'C')) THEN
195
ELSE IF (.NOT.LSAME(DIAG,'U') .AND. .NOT.LSAME(DIAG,'N')) THEN
197
ELSE IF (N.LT.0) THEN
199
ELSE IF (LDA.LT.MAX(1,N)) THEN
201
ELSE IF (INCX.EQ.0) THEN
144
CALL XERBLA( 'ZTRSV ', INFO )
205
CALL XERBLA('ZTRSV ',INFO)
148
209
* Quick return if possible.
153
NOCONJ = LSAME( TRANS, 'T' )
154
NOUNIT = LSAME( DIAG , 'N' )
213
NOCONJ = LSAME(TRANS,'T')
214
NOUNIT = LSAME(DIAG,'N')
156
216
* Set up the start point in X if the increment is not unity. This
157
217
* will be ( N - 1 )*INCX too small for descending loops.
158
* The next line is to satisfy compiler warnings.
162
KX = 1 - ( N - 1 )*INCX
163
ELSE IF( INCX.NE.1 )THEN
221
ELSE IF (INCX.NE.1) THEN
167
225
* Start the operations. In this version the elements of A are
168
226
* accessed sequentially with one pass through A.
170
IF( LSAME( TRANS, 'N' ) )THEN
228
IF (LSAME(TRANS,'N')) THEN
172
230
* Form x := inv( A )*x.
174
IF( LSAME( UPLO, 'U' ) )THEN
177
IF( X( J ).NE.ZERO )THEN
179
$ X( J ) = X( J )/A( J, J )
181
DO 10, I = J - 1, 1, -1
182
X( I ) = X( I ) - TEMP*A( I, J )
187
JX = KX + ( N - 1 )*INCX
189
IF( X( JX ).NE.ZERO )THEN
191
$ X( JX ) = X( JX )/A( J, J )
194
DO 30, I = J - 1, 1, -1
196
X( IX ) = X( IX ) - TEMP*A( I, J )
205
IF( X( J ).NE.ZERO )THEN
207
$ X( J ) = X( J )/A( J, J )
210
X( I ) = X( I ) - TEMP*A( I, J )
217
IF( X( JX ).NE.ZERO )THEN
219
$ X( JX ) = X( JX )/A( J, J )
224
X( IX ) = X( IX ) - TEMP*A( I, J )
232
IF (LSAME(UPLO,'U')) THEN
235
IF (X(J).NE.ZERO) THEN
236
IF (NOUNIT) X(J) = X(J)/A(J,J)
239
X(I) = X(I) - TEMP*A(I,J)
246
IF (X(JX).NE.ZERO) THEN
247
IF (NOUNIT) X(JX) = X(JX)/A(J,J)
252
X(IX) = X(IX) - TEMP*A(I,J)
261
IF (X(J).NE.ZERO) THEN
262
IF (NOUNIT) X(J) = X(J)/A(J,J)
265
X(I) = X(I) - TEMP*A(I,J)
272
IF (X(JX).NE.ZERO) THEN
273
IF (NOUNIT) X(JX) = X(JX)/A(J,J)
278
X(IX) = X(IX) - TEMP*A(I,J)
233
* Form x := inv( A' )*x or x := inv( conjg( A' ) )*x.
287
* Form x := inv( A**T )*x or x := inv( A**H )*x.
235
IF( LSAME( UPLO, 'U' ) )THEN
241
TEMP = TEMP - A( I, J )*X( I )
244
$ TEMP = TEMP/A( J, J )
247
TEMP = TEMP - DCONJG( A( I, J ) )*X( I )
250
$ TEMP = TEMP/DCONJG( A( J, J ) )
261
TEMP = TEMP - A( I, J )*X( IX )
265
$ TEMP = TEMP/A( J, J )
268
TEMP = TEMP - DCONJG( A( I, J ) )*X( IX )
272
$ TEMP = TEMP/DCONJG( A( J, J ) )
283
DO 150, I = N, J + 1, -1
284
TEMP = TEMP - A( I, J )*X( I )
287
$ TEMP = TEMP/A( J, J )
289
DO 160, I = N, J + 1, -1
290
TEMP = TEMP - DCONJG( A( I, J ) )*X( I )
293
$ TEMP = TEMP/DCONJG( A( J, J ) )
298
KX = KX + ( N - 1 )*INCX
304
DO 180, I = N, J + 1, -1
305
TEMP = TEMP - A( I, J )*X( IX )
309
$ TEMP = TEMP/A( J, J )
311
DO 190, I = N, J + 1, -1
312
TEMP = TEMP - DCONJG( A( I, J ) )*X( IX )
316
$ TEMP = TEMP/DCONJG( A( J, J ) )
289
IF (LSAME(UPLO,'U')) THEN
295
TEMP = TEMP - A(I,J)*X(I)
297
IF (NOUNIT) TEMP = TEMP/A(J,J)
300
TEMP = TEMP - DCONJG(A(I,J))*X(I)
302
IF (NOUNIT) TEMP = TEMP/DCONJG(A(J,J))
313
TEMP = TEMP - A(I,J)*X(IX)
316
IF (NOUNIT) TEMP = TEMP/A(J,J)
319
TEMP = TEMP - DCONJG(A(I,J))*X(IX)
322
IF (NOUNIT) TEMP = TEMP/DCONJG(A(J,J))
333
DO 150 I = N,J + 1,-1
334
TEMP = TEMP - A(I,J)*X(I)
336
IF (NOUNIT) TEMP = TEMP/A(J,J)
338
DO 160 I = N,J + 1,-1
339
TEMP = TEMP - DCONJG(A(I,J))*X(I)
341
IF (NOUNIT) TEMP = TEMP/DCONJG(A(J,J))
352
DO 180 I = N,J + 1,-1
353
TEMP = TEMP - A(I,J)*X(IX)
356
IF (NOUNIT) TEMP = TEMP/A(J,J)
358
DO 190 I = N,J + 1,-1
359
TEMP = TEMP - DCONJG(A(I,J))*X(IX)
362
IF (NOUNIT) TEMP = TEMP/DCONJG(A(J,J))