110 SUBROUTINE zspt03( UPLO, N, A, AINV, WORK, LDW, RWORK, RCOND,
121 DOUBLE PRECISION rcond, resid
124 DOUBLE PRECISION rwork( * )
125 COMPLEX*16 a( * ), ainv( * ), work( ldw, * )
131 DOUBLE PRECISION zero, one
132 parameter( zero = 0.0d+0, one = 1.0d+0 )
135 INTEGER i, icol,
j, jcol, k, kcol, nall
136 DOUBLE PRECISION ainvnm, anorm, eps
161 anorm =
zlansp(
'1', uplo, n, a, rwork )
162 ainvnm =
zlansp(
'1', uplo, n, ainv, rwork )
163 IF( anorm.LE.zero .OR. ainvnm.LE.zero )
THEN
168 rcond = ( one / anorm ) / ainvnm
174 IF(
lsame( uplo,
'U' ) )
THEN
176 icol = ( ( i-1 )*i ) / 2 + 1
181 jcol = ( (
j-1 )*
j ) / 2 + 1
182 t =
zdotu(
j, a( icol ), 1, ainv( jcol ), 1 )
183 jcol = jcol + 2*
j - 1
186 t = t + a( kcol+k )*ainv( jcol )
191 t = t + a( kcol )*ainv( jcol )
201 jcol = ( (
j-1 )*
j ) / 2 + 1
202 t =
zdotu( i, a( icol ), 1, ainv( jcol ), 1 )
204 kcol = icol + 2*i - 1
206 t = t + a( kcol )*ainv( jcol+k )
211 t = t + a( kcol )*ainv( jcol )
222 nall = ( n*( n+1 ) ) / 2
227 icol = nall - ( ( n-i+1 )*( n-i+2 ) ) / 2 + 1
229 jcol = nall - ( ( n-
j )*( n-
j+1 ) ) / 2 - ( n-i )
230 t =
zdotu( n-i+1, a( icol ), 1, ainv( jcol ), 1 )
234 t = t + a( kcol )*ainv( jcol )
240 t = t + a( kcol )*ainv( jcol+k )
248 icol = nall - ( ( n-i )*( n-i+1 ) ) / 2
250 jcol = nall - ( ( n-
j+1 )*( n-
j+2 ) ) / 2 + 1
251 t =
zdotu( n-
j+1, a( icol-n+
j ), 1, ainv( jcol ), 1 )
255 t = t + a( kcol )*ainv( jcol )
261 t = t + a( kcol+k )*ainv( jcol )
272 work( i, i ) = work( i, i ) + one
277 resid =
zlange(
'1', n, n, work, ldw, rwork )
279 resid = ( ( resid*rcond ) / eps ) / dble( n )
subroutine zspt03(UPLO, N, A, AINV, WORK, LDW, RWORK, RCOND, RESID)
ZSPT03
logical function lsame(CA, CB)
LSAME
double precision function zlange(NORM, M, N, A, LDA, WORK)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
double precision function dlamch(CMACH)
DLAMCH
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
complex *16 function zdotu(N, ZX, INCX, ZY, INCY)
ZDOTU
double precision function zlansp(NORM, UPLO, N, AP, WORK)
ZLANSP returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a symmetric matrix supplied in packed form.