110 SUBROUTINE zhptri( UPLO, N, AP, IPIV, WORK, INFO )
123 COMPLEX*16 AP( * ), WORK( * )
130 COMPLEX*16 CONE, ZERO
131 parameter ( one = 1.0d+0, cone = ( 1.0d+0, 0.0d+0 ),
132 $ zero = ( 0.0d+0, 0.0d+0 ) )
136 INTEGER J, K, KC, KCNEXT, KP, KPC, KSTEP, KX, NPP
137 DOUBLE PRECISION AK, AKP1, D, T
138 COMPLEX*16 AKKP1, TEMP
143 EXTERNAL lsame, zdotc
149 INTRINSIC abs, dble, dconjg
156 upper = lsame( uplo,
'U' )
157 IF( .NOT.upper .AND. .NOT.lsame( uplo,
'L' ) )
THEN
159 ELSE IF( n.LT.0 )
THEN
163 CALL xerbla(
'ZHPTRI', -info )
179 DO 10 info = n, 1, -1
180 IF( ipiv( info ).GT.0 .AND. ap( kp ).EQ.zero )
190 IF( ipiv( info ).GT.0 .AND. ap( kp ).EQ.zero )
192 kp = kp + n - info + 1
214 IF( ipiv( k ).GT.0 )
THEN
220 ap( kc+k-1 ) = one / dble( ap( kc+k-1 ) )
225 CALL zcopy( k-1, ap( kc ), 1, work, 1 )
226 CALL zhpmv( uplo, k-1, -cone, ap, work, 1, zero,
228 ap( kc+k-1 ) = ap( kc+k-1 ) -
229 $ dble( zdotc( k-1, work, 1, ap( kc ), 1 ) )
238 t = abs( ap( kcnext+k-1 ) )
239 ak = dble( ap( kc+k-1 ) ) / t
240 akp1 = dble( ap( kcnext+k ) ) / t
241 akkp1 = ap( kcnext+k-1 ) / t
242 d = t*( ak*akp1-one )
243 ap( kc+k-1 ) = akp1 / d
244 ap( kcnext+k ) = ak / d
245 ap( kcnext+k-1 ) = -akkp1 / d
250 CALL zcopy( k-1, ap( kc ), 1, work, 1 )
251 CALL zhpmv( uplo, k-1, -cone, ap, work, 1, zero,
253 ap( kc+k-1 ) = ap( kc+k-1 ) -
254 $ dble( zdotc( k-1, work, 1, ap( kc ), 1 ) )
255 ap( kcnext+k-1 ) = ap( kcnext+k-1 ) -
256 $ zdotc( k-1, ap( kc ), 1, ap( kcnext ),
258 CALL zcopy( k-1, ap( kcnext ), 1, work, 1 )
259 CALL zhpmv( uplo, k-1, -cone, ap, work, 1, zero,
261 ap( kcnext+k ) = ap( kcnext+k ) -
262 $ dble( zdotc( k-1, work, 1, ap( kcnext ),
266 kcnext = kcnext + k + 1
269 kp = abs( ipiv( k ) )
275 kpc = ( kp-1 )*kp / 2 + 1
276 CALL zswap( kp-1, ap( kc ), 1, ap( kpc ), 1 )
278 DO 40 j = kp + 1, k - 1
280 temp = dconjg( ap( kc+j-1 ) )
281 ap( kc+j-1 ) = dconjg( ap( kx ) )
284 ap( kc+kp-1 ) = dconjg( ap( kc+kp-1 ) )
286 ap( kc+k-1 ) = ap( kpc+kp-1 )
287 ap( kpc+kp-1 ) = temp
288 IF( kstep.EQ.2 )
THEN
289 temp = ap( kc+k+k-1 )
290 ap( kc+k+k-1 ) = ap( kc+k+kp-1 )
291 ap( kc+k+kp-1 ) = temp
317 kcnext = kc - ( n-k+2 )
318 IF( ipiv( k ).GT.0 )
THEN
324 ap( kc ) = one / dble( ap( kc ) )
329 CALL zcopy( n-k, ap( kc+1 ), 1, work, 1 )
330 CALL zhpmv( uplo, n-k, -cone, ap( kc+n-k+1 ), work, 1,
331 $ zero, ap( kc+1 ), 1 )
332 ap( kc ) = ap( kc ) - dble( zdotc( n-k, work, 1,
342 t = abs( ap( kcnext+1 ) )
343 ak = dble( ap( kcnext ) ) / t
344 akp1 = dble( ap( kc ) ) / t
345 akkp1 = ap( kcnext+1 ) / t
346 d = t*( ak*akp1-one )
347 ap( kcnext ) = akp1 / d
349 ap( kcnext+1 ) = -akkp1 / d
354 CALL zcopy( n-k, ap( kc+1 ), 1, work, 1 )
355 CALL zhpmv( uplo, n-k, -cone, ap( kc+( n-k+1 ) ), work,
356 $ 1, zero, ap( kc+1 ), 1 )
357 ap( kc ) = ap( kc ) - dble( zdotc( n-k, work, 1,
359 ap( kcnext+1 ) = ap( kcnext+1 ) -
360 $ zdotc( n-k, ap( kc+1 ), 1,
361 $ ap( kcnext+2 ), 1 )
362 CALL zcopy( n-k, ap( kcnext+2 ), 1, work, 1 )
363 CALL zhpmv( uplo, n-k, -cone, ap( kc+( n-k+1 ) ), work,
364 $ 1, zero, ap( kcnext+2 ), 1 )
365 ap( kcnext ) = ap( kcnext ) -
366 $ dble( zdotc( n-k, work, 1, ap( kcnext+2 ),
370 kcnext = kcnext - ( n-k+3 )
373 kp = abs( ipiv( k ) )
379 kpc = npp - ( n-kp+1 )*( n-kp+2 ) / 2 + 1
381 $
CALL zswap( n-kp, ap( kc+kp-k+1 ), 1, ap( kpc+1 ), 1 )
383 DO 70 j = k + 1, kp - 1
385 temp = dconjg( ap( kc+j-k ) )
386 ap( kc+j-k ) = dconjg( ap( kx ) )
389 ap( kc+kp-k ) = dconjg( ap( kc+kp-k ) )
393 IF( kstep.EQ.2 )
THEN
394 temp = ap( kc-n+k-1 )
395 ap( kc-n+k-1 ) = ap( kc-n+kp-1 )
396 ap( kc-n+kp-1 ) = temp
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
subroutine zswap(N, ZX, INCX, ZY, INCY)
ZSWAP
subroutine zhptri(UPLO, N, AP, IPIV, WORK, INFO)
ZHPTRI
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine zhpmv(UPLO, N, ALPHA, AP, X, INCX, BETA, Y, INCY)
ZHPMV