1 SUBROUTINE ZSPTRI( UPLO, N, AP, IPIV, WORK, INFO )
2 *
3 * -- LAPACK routine (version 3.3.1) --
4 * -- LAPACK is a software package provided by Univ. of Tennessee, --
5 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
6 * -- April 2011 --
7 *
8 * .. Scalar Arguments ..
9 CHARACTER UPLO
10 INTEGER INFO, N
11 * ..
12 * .. Array Arguments ..
13 INTEGER IPIV( * )
14 COMPLEX*16 AP( * ), WORK( * )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * ZSPTRI computes the inverse of a complex symmetric indefinite matrix
21 * A in packed storage using the factorization A = U*D*U**T or
22 * A = L*D*L**T computed by ZSPTRF.
23 *
24 * Arguments
25 * =========
26 *
27 * UPLO (input) CHARACTER*1
28 * Specifies whether the details of the factorization are stored
29 * as an upper or lower triangular matrix.
30 * = 'U': Upper triangular, form is A = U*D*U**T;
31 * = 'L': Lower triangular, form is A = L*D*L**T.
32 *
33 * N (input) INTEGER
34 * The order of the matrix A. N >= 0.
35 *
36 * AP (input/output) COMPLEX*16 array, dimension (N*(N+1)/2)
37 * On entry, the block diagonal matrix D and the multipliers
38 * used to obtain the factor U or L as computed by ZSPTRF,
39 * stored as a packed triangular matrix.
40 *
41 * On exit, if INFO = 0, the (symmetric) inverse of the original
42 * matrix, stored as a packed triangular matrix. The j-th column
43 * of inv(A) is stored in the array AP as follows:
44 * if UPLO = 'U', AP(i + (j-1)*j/2) = inv(A)(i,j) for 1<=i<=j;
45 * if UPLO = 'L',
46 * AP(i + (j-1)*(2n-j)/2) = inv(A)(i,j) for j<=i<=n.
47 *
48 * IPIV (input) INTEGER array, dimension (N)
49 * Details of the interchanges and the block structure of D
50 * as determined by ZSPTRF.
51 *
52 * WORK (workspace) COMPLEX*16 array, dimension (N)
53 *
54 * INFO (output) INTEGER
55 * = 0: successful exit
56 * < 0: if INFO = -i, the i-th argument had an illegal value
57 * > 0: if INFO = i, D(i,i) = 0; the matrix is singular and its
58 * inverse could not be computed.
59 *
60 * =====================================================================
61 *
62 * .. Parameters ..
63 COMPLEX*16 ONE, ZERO
64 PARAMETER ( ONE = ( 1.0D+0, 0.0D+0 ),
65 $ ZERO = ( 0.0D+0, 0.0D+0 ) )
66 * ..
67 * .. Local Scalars ..
68 LOGICAL UPPER
69 INTEGER J, K, KC, KCNEXT, KP, KPC, KSTEP, KX, NPP
70 COMPLEX*16 AK, AKKP1, AKP1, D, T, TEMP
71 * ..
72 * .. External Functions ..
73 LOGICAL LSAME
74 COMPLEX*16 ZDOTU
75 EXTERNAL LSAME, ZDOTU
76 * ..
77 * .. External Subroutines ..
78 EXTERNAL XERBLA, ZCOPY, ZSPMV, ZSWAP
79 * ..
80 * .. Intrinsic Functions ..
81 INTRINSIC ABS
82 * ..
83 * .. Executable Statements ..
84 *
85 * Test the input parameters.
86 *
87 INFO = 0
88 UPPER = LSAME( UPLO, 'U' )
89 IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
90 INFO = -1
91 ELSE IF( N.LT.0 ) THEN
92 INFO = -2
93 END IF
94 IF( INFO.NE.0 ) THEN
95 CALL XERBLA( 'ZSPTRI', -INFO )
96 RETURN
97 END IF
98 *
99 * Quick return if possible
100 *
101 IF( N.EQ.0 )
102 $ RETURN
103 *
104 * Check that the diagonal matrix D is nonsingular.
105 *
106 IF( UPPER ) THEN
107 *
108 * Upper triangular storage: examine D from bottom to top
109 *
110 KP = N*( N+1 ) / 2
111 DO 10 INFO = N, 1, -1
112 IF( IPIV( INFO ).GT.0 .AND. AP( KP ).EQ.ZERO )
113 $ RETURN
114 KP = KP - INFO
115 10 CONTINUE
116 ELSE
117 *
118 * Lower triangular storage: examine D from top to bottom.
119 *
120 KP = 1
121 DO 20 INFO = 1, N
122 IF( IPIV( INFO ).GT.0 .AND. AP( KP ).EQ.ZERO )
123 $ RETURN
124 KP = KP + N - INFO + 1
125 20 CONTINUE
126 END IF
127 INFO = 0
128 *
129 IF( UPPER ) THEN
130 *
131 * Compute inv(A) from the factorization A = U*D*U**T.
132 *
133 * K is the main loop index, increasing from 1 to N in steps of
134 * 1 or 2, depending on the size of the diagonal blocks.
135 *
136 K = 1
137 KC = 1
138 30 CONTINUE
139 *
140 * If K > N, exit from loop.
141 *
142 IF( K.GT.N )
143 $ GO TO 50
144 *
145 KCNEXT = KC + K
146 IF( IPIV( K ).GT.0 ) THEN
147 *
148 * 1 x 1 diagonal block
149 *
150 * Invert the diagonal block.
151 *
152 AP( KC+K-1 ) = ONE / AP( KC+K-1 )
153 *
154 * Compute column K of the inverse.
155 *
156 IF( K.GT.1 ) THEN
157 CALL ZCOPY( K-1, AP( KC ), 1, WORK, 1 )
158 CALL ZSPMV( UPLO, K-1, -ONE, AP, WORK, 1, ZERO, AP( KC ),
159 $ 1 )
160 AP( KC+K-1 ) = AP( KC+K-1 ) -
161 $ ZDOTU( K-1, WORK, 1, AP( KC ), 1 )
162 END IF
163 KSTEP = 1
164 ELSE
165 *
166 * 2 x 2 diagonal block
167 *
168 * Invert the diagonal block.
169 *
170 T = AP( KCNEXT+K-1 )
171 AK = AP( KC+K-1 ) / T
172 AKP1 = AP( KCNEXT+K ) / T
173 AKKP1 = AP( KCNEXT+K-1 ) / T
174 D = T*( AK*AKP1-ONE )
175 AP( KC+K-1 ) = AKP1 / D
176 AP( KCNEXT+K ) = AK / D
177 AP( KCNEXT+K-1 ) = -AKKP1 / D
178 *
179 * Compute columns K and K+1 of the inverse.
180 *
181 IF( K.GT.1 ) THEN
182 CALL ZCOPY( K-1, AP( KC ), 1, WORK, 1 )
183 CALL ZSPMV( UPLO, K-1, -ONE, AP, WORK, 1, ZERO, AP( KC ),
184 $ 1 )
185 AP( KC+K-1 ) = AP( KC+K-1 ) -
186 $ ZDOTU( K-1, WORK, 1, AP( KC ), 1 )
187 AP( KCNEXT+K-1 ) = AP( KCNEXT+K-1 ) -
188 $ ZDOTU( K-1, AP( KC ), 1, AP( KCNEXT ),
189 $ 1 )
190 CALL ZCOPY( K-1, AP( KCNEXT ), 1, WORK, 1 )
191 CALL ZSPMV( UPLO, K-1, -ONE, AP, WORK, 1, ZERO,
192 $ AP( KCNEXT ), 1 )
193 AP( KCNEXT+K ) = AP( KCNEXT+K ) -
194 $ ZDOTU( K-1, WORK, 1, AP( KCNEXT ), 1 )
195 END IF
196 KSTEP = 2
197 KCNEXT = KCNEXT + K + 1
198 END IF
199 *
200 KP = ABS( IPIV( K ) )
201 IF( KP.NE.K ) THEN
202 *
203 * Interchange rows and columns K and KP in the leading
204 * submatrix A(1:k+1,1:k+1)
205 *
206 KPC = ( KP-1 )*KP / 2 + 1
207 CALL ZSWAP( KP-1, AP( KC ), 1, AP( KPC ), 1 )
208 KX = KPC + KP - 1
209 DO 40 J = KP + 1, K - 1
210 KX = KX + J - 1
211 TEMP = AP( KC+J-1 )
212 AP( KC+J-1 ) = AP( KX )
213 AP( KX ) = TEMP
214 40 CONTINUE
215 TEMP = AP( KC+K-1 )
216 AP( KC+K-1 ) = AP( KPC+KP-1 )
217 AP( KPC+KP-1 ) = TEMP
218 IF( KSTEP.EQ.2 ) THEN
219 TEMP = AP( KC+K+K-1 )
220 AP( KC+K+K-1 ) = AP( KC+K+KP-1 )
221 AP( KC+K+KP-1 ) = TEMP
222 END IF
223 END IF
224 *
225 K = K + KSTEP
226 KC = KCNEXT
227 GO TO 30
228 50 CONTINUE
229 *
230 ELSE
231 *
232 * Compute inv(A) from the factorization A = L*D*L**T.
233 *
234 * K is the main loop index, increasing from 1 to N in steps of
235 * 1 or 2, depending on the size of the diagonal blocks.
236 *
237 NPP = N*( N+1 ) / 2
238 K = N
239 KC = NPP
240 60 CONTINUE
241 *
242 * If K < 1, exit from loop.
243 *
244 IF( K.LT.1 )
245 $ GO TO 80
246 *
247 KCNEXT = KC - ( N-K+2 )
248 IF( IPIV( K ).GT.0 ) THEN
249 *
250 * 1 x 1 diagonal block
251 *
252 * Invert the diagonal block.
253 *
254 AP( KC ) = ONE / AP( KC )
255 *
256 * Compute column K of the inverse.
257 *
258 IF( K.LT.N ) THEN
259 CALL ZCOPY( N-K, AP( KC+1 ), 1, WORK, 1 )
260 CALL ZSPMV( UPLO, N-K, -ONE, AP( KC+N-K+1 ), WORK, 1,
261 $ ZERO, AP( KC+1 ), 1 )
262 AP( KC ) = AP( KC ) - ZDOTU( N-K, WORK, 1, AP( KC+1 ),
263 $ 1 )
264 END IF
265 KSTEP = 1
266 ELSE
267 *
268 * 2 x 2 diagonal block
269 *
270 * Invert the diagonal block.
271 *
272 T = AP( KCNEXT+1 )
273 AK = AP( KCNEXT ) / T
274 AKP1 = AP( KC ) / T
275 AKKP1 = AP( KCNEXT+1 ) / T
276 D = T*( AK*AKP1-ONE )
277 AP( KCNEXT ) = AKP1 / D
278 AP( KC ) = AK / D
279 AP( KCNEXT+1 ) = -AKKP1 / D
280 *
281 * Compute columns K-1 and K of the inverse.
282 *
283 IF( K.LT.N ) THEN
284 CALL ZCOPY( N-K, AP( KC+1 ), 1, WORK, 1 )
285 CALL ZSPMV( UPLO, N-K, -ONE, AP( KC+( N-K+1 ) ), WORK, 1,
286 $ ZERO, AP( KC+1 ), 1 )
287 AP( KC ) = AP( KC ) - ZDOTU( N-K, WORK, 1, AP( KC+1 ),
288 $ 1 )
289 AP( KCNEXT+1 ) = AP( KCNEXT+1 ) -
290 $ ZDOTU( N-K, AP( KC+1 ), 1,
291 $ AP( KCNEXT+2 ), 1 )
292 CALL ZCOPY( N-K, AP( KCNEXT+2 ), 1, WORK, 1 )
293 CALL ZSPMV( UPLO, N-K, -ONE, AP( KC+( N-K+1 ) ), WORK, 1,
294 $ ZERO, AP( KCNEXT+2 ), 1 )
295 AP( KCNEXT ) = AP( KCNEXT ) -
296 $ ZDOTU( N-K, WORK, 1, AP( KCNEXT+2 ), 1 )
297 END IF
298 KSTEP = 2
299 KCNEXT = KCNEXT - ( N-K+3 )
300 END IF
301 *
302 KP = ABS( IPIV( K ) )
303 IF( KP.NE.K ) THEN
304 *
305 * Interchange rows and columns K and KP in the trailing
306 * submatrix A(k-1:n,k-1:n)
307 *
308 KPC = NPP - ( N-KP+1 )*( N-KP+2 ) / 2 + 1
309 IF( KP.LT.N )
310 $ CALL ZSWAP( N-KP, AP( KC+KP-K+1 ), 1, AP( KPC+1 ), 1 )
311 KX = KC + KP - K
312 DO 70 J = K + 1, KP - 1
313 KX = KX + N - J + 1
314 TEMP = AP( KC+J-K )
315 AP( KC+J-K ) = AP( KX )
316 AP( KX ) = TEMP
317 70 CONTINUE
318 TEMP = AP( KC )
319 AP( KC ) = AP( KPC )
320 AP( KPC ) = TEMP
321 IF( KSTEP.EQ.2 ) THEN
322 TEMP = AP( KC-N+K-1 )
323 AP( KC-N+K-1 ) = AP( KC-N+KP-1 )
324 AP( KC-N+KP-1 ) = TEMP
325 END IF
326 END IF
327 *
328 K = K - KSTEP
329 KC = KCNEXT
330 GO TO 60
331 80 CONTINUE
332 END IF
333 *
334 RETURN
335 *
336 * End of ZSPTRI
337 *
338 END
2 *
3 * -- LAPACK routine (version 3.3.1) --
4 * -- LAPACK is a software package provided by Univ. of Tennessee, --
5 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
6 * -- April 2011 --
7 *
8 * .. Scalar Arguments ..
9 CHARACTER UPLO
10 INTEGER INFO, N
11 * ..
12 * .. Array Arguments ..
13 INTEGER IPIV( * )
14 COMPLEX*16 AP( * ), WORK( * )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * ZSPTRI computes the inverse of a complex symmetric indefinite matrix
21 * A in packed storage using the factorization A = U*D*U**T or
22 * A = L*D*L**T computed by ZSPTRF.
23 *
24 * Arguments
25 * =========
26 *
27 * UPLO (input) CHARACTER*1
28 * Specifies whether the details of the factorization are stored
29 * as an upper or lower triangular matrix.
30 * = 'U': Upper triangular, form is A = U*D*U**T;
31 * = 'L': Lower triangular, form is A = L*D*L**T.
32 *
33 * N (input) INTEGER
34 * The order of the matrix A. N >= 0.
35 *
36 * AP (input/output) COMPLEX*16 array, dimension (N*(N+1)/2)
37 * On entry, the block diagonal matrix D and the multipliers
38 * used to obtain the factor U or L as computed by ZSPTRF,
39 * stored as a packed triangular matrix.
40 *
41 * On exit, if INFO = 0, the (symmetric) inverse of the original
42 * matrix, stored as a packed triangular matrix. The j-th column
43 * of inv(A) is stored in the array AP as follows:
44 * if UPLO = 'U', AP(i + (j-1)*j/2) = inv(A)(i,j) for 1<=i<=j;
45 * if UPLO = 'L',
46 * AP(i + (j-1)*(2n-j)/2) = inv(A)(i,j) for j<=i<=n.
47 *
48 * IPIV (input) INTEGER array, dimension (N)
49 * Details of the interchanges and the block structure of D
50 * as determined by ZSPTRF.
51 *
52 * WORK (workspace) COMPLEX*16 array, dimension (N)
53 *
54 * INFO (output) INTEGER
55 * = 0: successful exit
56 * < 0: if INFO = -i, the i-th argument had an illegal value
57 * > 0: if INFO = i, D(i,i) = 0; the matrix is singular and its
58 * inverse could not be computed.
59 *
60 * =====================================================================
61 *
62 * .. Parameters ..
63 COMPLEX*16 ONE, ZERO
64 PARAMETER ( ONE = ( 1.0D+0, 0.0D+0 ),
65 $ ZERO = ( 0.0D+0, 0.0D+0 ) )
66 * ..
67 * .. Local Scalars ..
68 LOGICAL UPPER
69 INTEGER J, K, KC, KCNEXT, KP, KPC, KSTEP, KX, NPP
70 COMPLEX*16 AK, AKKP1, AKP1, D, T, TEMP
71 * ..
72 * .. External Functions ..
73 LOGICAL LSAME
74 COMPLEX*16 ZDOTU
75 EXTERNAL LSAME, ZDOTU
76 * ..
77 * .. External Subroutines ..
78 EXTERNAL XERBLA, ZCOPY, ZSPMV, ZSWAP
79 * ..
80 * .. Intrinsic Functions ..
81 INTRINSIC ABS
82 * ..
83 * .. Executable Statements ..
84 *
85 * Test the input parameters.
86 *
87 INFO = 0
88 UPPER = LSAME( UPLO, 'U' )
89 IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
90 INFO = -1
91 ELSE IF( N.LT.0 ) THEN
92 INFO = -2
93 END IF
94 IF( INFO.NE.0 ) THEN
95 CALL XERBLA( 'ZSPTRI', -INFO )
96 RETURN
97 END IF
98 *
99 * Quick return if possible
100 *
101 IF( N.EQ.0 )
102 $ RETURN
103 *
104 * Check that the diagonal matrix D is nonsingular.
105 *
106 IF( UPPER ) THEN
107 *
108 * Upper triangular storage: examine D from bottom to top
109 *
110 KP = N*( N+1 ) / 2
111 DO 10 INFO = N, 1, -1
112 IF( IPIV( INFO ).GT.0 .AND. AP( KP ).EQ.ZERO )
113 $ RETURN
114 KP = KP - INFO
115 10 CONTINUE
116 ELSE
117 *
118 * Lower triangular storage: examine D from top to bottom.
119 *
120 KP = 1
121 DO 20 INFO = 1, N
122 IF( IPIV( INFO ).GT.0 .AND. AP( KP ).EQ.ZERO )
123 $ RETURN
124 KP = KP + N - INFO + 1
125 20 CONTINUE
126 END IF
127 INFO = 0
128 *
129 IF( UPPER ) THEN
130 *
131 * Compute inv(A) from the factorization A = U*D*U**T.
132 *
133 * K is the main loop index, increasing from 1 to N in steps of
134 * 1 or 2, depending on the size of the diagonal blocks.
135 *
136 K = 1
137 KC = 1
138 30 CONTINUE
139 *
140 * If K > N, exit from loop.
141 *
142 IF( K.GT.N )
143 $ GO TO 50
144 *
145 KCNEXT = KC + K
146 IF( IPIV( K ).GT.0 ) THEN
147 *
148 * 1 x 1 diagonal block
149 *
150 * Invert the diagonal block.
151 *
152 AP( KC+K-1 ) = ONE / AP( KC+K-1 )
153 *
154 * Compute column K of the inverse.
155 *
156 IF( K.GT.1 ) THEN
157 CALL ZCOPY( K-1, AP( KC ), 1, WORK, 1 )
158 CALL ZSPMV( UPLO, K-1, -ONE, AP, WORK, 1, ZERO, AP( KC ),
159 $ 1 )
160 AP( KC+K-1 ) = AP( KC+K-1 ) -
161 $ ZDOTU( K-1, WORK, 1, AP( KC ), 1 )
162 END IF
163 KSTEP = 1
164 ELSE
165 *
166 * 2 x 2 diagonal block
167 *
168 * Invert the diagonal block.
169 *
170 T = AP( KCNEXT+K-1 )
171 AK = AP( KC+K-1 ) / T
172 AKP1 = AP( KCNEXT+K ) / T
173 AKKP1 = AP( KCNEXT+K-1 ) / T
174 D = T*( AK*AKP1-ONE )
175 AP( KC+K-1 ) = AKP1 / D
176 AP( KCNEXT+K ) = AK / D
177 AP( KCNEXT+K-1 ) = -AKKP1 / D
178 *
179 * Compute columns K and K+1 of the inverse.
180 *
181 IF( K.GT.1 ) THEN
182 CALL ZCOPY( K-1, AP( KC ), 1, WORK, 1 )
183 CALL ZSPMV( UPLO, K-1, -ONE, AP, WORK, 1, ZERO, AP( KC ),
184 $ 1 )
185 AP( KC+K-1 ) = AP( KC+K-1 ) -
186 $ ZDOTU( K-1, WORK, 1, AP( KC ), 1 )
187 AP( KCNEXT+K-1 ) = AP( KCNEXT+K-1 ) -
188 $ ZDOTU( K-1, AP( KC ), 1, AP( KCNEXT ),
189 $ 1 )
190 CALL ZCOPY( K-1, AP( KCNEXT ), 1, WORK, 1 )
191 CALL ZSPMV( UPLO, K-1, -ONE, AP, WORK, 1, ZERO,
192 $ AP( KCNEXT ), 1 )
193 AP( KCNEXT+K ) = AP( KCNEXT+K ) -
194 $ ZDOTU( K-1, WORK, 1, AP( KCNEXT ), 1 )
195 END IF
196 KSTEP = 2
197 KCNEXT = KCNEXT + K + 1
198 END IF
199 *
200 KP = ABS( IPIV( K ) )
201 IF( KP.NE.K ) THEN
202 *
203 * Interchange rows and columns K and KP in the leading
204 * submatrix A(1:k+1,1:k+1)
205 *
206 KPC = ( KP-1 )*KP / 2 + 1
207 CALL ZSWAP( KP-1, AP( KC ), 1, AP( KPC ), 1 )
208 KX = KPC + KP - 1
209 DO 40 J = KP + 1, K - 1
210 KX = KX + J - 1
211 TEMP = AP( KC+J-1 )
212 AP( KC+J-1 ) = AP( KX )
213 AP( KX ) = TEMP
214 40 CONTINUE
215 TEMP = AP( KC+K-1 )
216 AP( KC+K-1 ) = AP( KPC+KP-1 )
217 AP( KPC+KP-1 ) = TEMP
218 IF( KSTEP.EQ.2 ) THEN
219 TEMP = AP( KC+K+K-1 )
220 AP( KC+K+K-1 ) = AP( KC+K+KP-1 )
221 AP( KC+K+KP-1 ) = TEMP
222 END IF
223 END IF
224 *
225 K = K + KSTEP
226 KC = KCNEXT
227 GO TO 30
228 50 CONTINUE
229 *
230 ELSE
231 *
232 * Compute inv(A) from the factorization A = L*D*L**T.
233 *
234 * K is the main loop index, increasing from 1 to N in steps of
235 * 1 or 2, depending on the size of the diagonal blocks.
236 *
237 NPP = N*( N+1 ) / 2
238 K = N
239 KC = NPP
240 60 CONTINUE
241 *
242 * If K < 1, exit from loop.
243 *
244 IF( K.LT.1 )
245 $ GO TO 80
246 *
247 KCNEXT = KC - ( N-K+2 )
248 IF( IPIV( K ).GT.0 ) THEN
249 *
250 * 1 x 1 diagonal block
251 *
252 * Invert the diagonal block.
253 *
254 AP( KC ) = ONE / AP( KC )
255 *
256 * Compute column K of the inverse.
257 *
258 IF( K.LT.N ) THEN
259 CALL ZCOPY( N-K, AP( KC+1 ), 1, WORK, 1 )
260 CALL ZSPMV( UPLO, N-K, -ONE, AP( KC+N-K+1 ), WORK, 1,
261 $ ZERO, AP( KC+1 ), 1 )
262 AP( KC ) = AP( KC ) - ZDOTU( N-K, WORK, 1, AP( KC+1 ),
263 $ 1 )
264 END IF
265 KSTEP = 1
266 ELSE
267 *
268 * 2 x 2 diagonal block
269 *
270 * Invert the diagonal block.
271 *
272 T = AP( KCNEXT+1 )
273 AK = AP( KCNEXT ) / T
274 AKP1 = AP( KC ) / T
275 AKKP1 = AP( KCNEXT+1 ) / T
276 D = T*( AK*AKP1-ONE )
277 AP( KCNEXT ) = AKP1 / D
278 AP( KC ) = AK / D
279 AP( KCNEXT+1 ) = -AKKP1 / D
280 *
281 * Compute columns K-1 and K of the inverse.
282 *
283 IF( K.LT.N ) THEN
284 CALL ZCOPY( N-K, AP( KC+1 ), 1, WORK, 1 )
285 CALL ZSPMV( UPLO, N-K, -ONE, AP( KC+( N-K+1 ) ), WORK, 1,
286 $ ZERO, AP( KC+1 ), 1 )
287 AP( KC ) = AP( KC ) - ZDOTU( N-K, WORK, 1, AP( KC+1 ),
288 $ 1 )
289 AP( KCNEXT+1 ) = AP( KCNEXT+1 ) -
290 $ ZDOTU( N-K, AP( KC+1 ), 1,
291 $ AP( KCNEXT+2 ), 1 )
292 CALL ZCOPY( N-K, AP( KCNEXT+2 ), 1, WORK, 1 )
293 CALL ZSPMV( UPLO, N-K, -ONE, AP( KC+( N-K+1 ) ), WORK, 1,
294 $ ZERO, AP( KCNEXT+2 ), 1 )
295 AP( KCNEXT ) = AP( KCNEXT ) -
296 $ ZDOTU( N-K, WORK, 1, AP( KCNEXT+2 ), 1 )
297 END IF
298 KSTEP = 2
299 KCNEXT = KCNEXT - ( N-K+3 )
300 END IF
301 *
302 KP = ABS( IPIV( K ) )
303 IF( KP.NE.K ) THEN
304 *
305 * Interchange rows and columns K and KP in the trailing
306 * submatrix A(k-1:n,k-1:n)
307 *
308 KPC = NPP - ( N-KP+1 )*( N-KP+2 ) / 2 + 1
309 IF( KP.LT.N )
310 $ CALL ZSWAP( N-KP, AP( KC+KP-K+1 ), 1, AP( KPC+1 ), 1 )
311 KX = KC + KP - K
312 DO 70 J = K + 1, KP - 1
313 KX = KX + N - J + 1
314 TEMP = AP( KC+J-K )
315 AP( KC+J-K ) = AP( KX )
316 AP( KX ) = TEMP
317 70 CONTINUE
318 TEMP = AP( KC )
319 AP( KC ) = AP( KPC )
320 AP( KPC ) = TEMP
321 IF( KSTEP.EQ.2 ) THEN
322 TEMP = AP( KC-N+K-1 )
323 AP( KC-N+K-1 ) = AP( KC-N+KP-1 )
324 AP( KC-N+KP-1 ) = TEMP
325 END IF
326 END IF
327 *
328 K = K - KSTEP
329 KC = KCNEXT
330 GO TO 60
331 80 CONTINUE
332 END IF
333 *
334 RETURN
335 *
336 * End of ZSPTRI
337 *
338 END