LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine zsytri_rook ( character  UPLO,
integer  N,
complex*16, dimension( lda, * )  A,
integer  LDA,
integer, dimension( * )  IPIV,
complex*16, dimension( * )  WORK,
integer  INFO 
)

ZSYTRI_ROOK

Download ZSYTRI_ROOK + dependencies [TGZ] [ZIP] [TXT]

Purpose: 
 ZSYTRI_ROOK computes the inverse of a complex symmetric
 matrix A using the factorization A = U*D*U**T or A = L*D*L**T
 computed by ZSYTRF_ROOK.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          Specifies whether the details of the factorization are stored
          as an upper or lower triangular matrix.
          = 'U':  Upper triangular, form is A = U*D*U**T;
          = 'L':  Lower triangular, form is A = L*D*L**T.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is COMPLEX*16 array, dimension (LDA,N)
          On entry, the block diagonal matrix D and the multipliers
          used to obtain the factor U or L as computed by ZSYTRF_ROOK.

          On exit, if INFO = 0, the (symmetric) inverse of the original
          matrix.  If UPLO = 'U', the upper triangular part of the
          inverse is formed and the part of A below the diagonal is not
          referenced; if UPLO = 'L' the lower triangular part of the
          inverse is formed and the part of A above the diagonal is
          not referenced.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          Details of the interchanges and the block structure of D
          as determined by ZSYTRF_ROOK.
[out]WORK
          WORK is COMPLEX*16 array, dimension (N)
[out]INFO
          INFO is INTEGER
          = 0: successful exit
          < 0: if INFO = -i, the i-th argument had an illegal value
          > 0: if INFO = i, D(i,i) = 0; the matrix is singular and its
               inverse could not be computed.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015
Contributors: 
   November 2015, Igor Kozachenko,
                  Computer Science Division,
                  University of California, Berkeley

  September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas,
                  School of Mathematics,
                  University of Manchester

Definition at line 131 of file zsytri_rook.f.

131 *
132 * -- LAPACK computational routine (version 3.6.0) --
133 * -- LAPACK is a software package provided by Univ. of Tennessee, --
134 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
135 * November 2015
136 *
137 * .. Scalar Arguments ..
138  CHARACTER uplo
139  INTEGER info, lda, n
140 * ..
141 * .. Array Arguments ..
142  INTEGER ipiv( * )
143  COMPLEX*16 a( lda, * ), work( * )
144 * ..
145 *
146 * =====================================================================
147 *
148 * .. Parameters ..
149  COMPLEX*16 cone, czero
150  parameter ( cone = ( 1.0d+0, 0.0d+0 ),
151  $ czero = ( 0.0d+0, 0.0d+0 ) )
152 * ..
153 * .. Local Scalars ..
154  LOGICAL upper
155  INTEGER k, kp, kstep
156  COMPLEX*16 ak, akkp1, akp1, d, t, temp
157 * ..
158 * .. External Functions ..
159  LOGICAL lsame
160  COMPLEX*16 zdotu
161  EXTERNAL lsame, zdotu
162 * ..
163 * .. External Subroutines ..
164  EXTERNAL zcopy, zswap, zsymv, xerbla
165 * ..
166 * .. Intrinsic Functions ..
167  INTRINSIC max
168 * ..
169 * .. Executable Statements ..
170 *
171 * Test the input parameters.
172 *
173  info = 0
174  upper = lsame( uplo, 'U' )
175  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
176  info = -1
177  ELSE IF( n.LT.0 ) THEN
178  info = -2
179  ELSE IF( lda.LT.max( 1, n ) ) THEN
180  info = -4
181  END IF
182  IF( info.NE.0 ) THEN
183  CALL xerbla( 'ZSYTRI_ROOK', -info )
184  RETURN
185  END IF
186 *
187 * Quick return if possible
188 *
189  IF( n.EQ.0 )
190  $ RETURN
191 *
192 * Check that the diagonal matrix D is nonsingular.
193 *
194  IF( upper ) THEN
195 *
196 * Upper triangular storage: examine D from bottom to top
197 *
198  DO 10 info = n, 1, -1
199  IF( ipiv( info ).GT.0 .AND. a( info, info ).EQ.czero )
200  $ RETURN
201  10 CONTINUE
202  ELSE
203 *
204 * Lower triangular storage: examine D from top to bottom.
205 *
206  DO 20 info = 1, n
207  IF( ipiv( info ).GT.0 .AND. a( info, info ).EQ.czero )
208  $ RETURN
209  20 CONTINUE
210  END IF
211  info = 0
212 *
213  IF( upper ) THEN
214 *
215 * Compute inv(A) from the factorization A = U*D*U**T.
216 *
217 * K is the main loop index, increasing from 1 to N in steps of
218 * 1 or 2, depending on the size of the diagonal blocks.
219 *
220  k = 1
221  30 CONTINUE
222 *
223 * If K > N, exit from loop.
224 *
225  IF( k.GT.n )
226  $ GO TO 40
227 *
228  IF( ipiv( k ).GT.0 ) THEN
229 *
230 * 1 x 1 diagonal block
231 *
232 * Invert the diagonal block.
233 *
234  a( k, k ) = cone / a( k, k )
235 *
236 * Compute column K of the inverse.
237 *
238  IF( k.GT.1 ) THEN
239  CALL zcopy( k-1, a( 1, k ), 1, work, 1 )
240  CALL zsymv( uplo, k-1, -cone, a, lda, work, 1, czero,
241  $ a( 1, k ), 1 )
242  a( k, k ) = a( k, k ) - zdotu( k-1, work, 1, a( 1, k ),
243  $ 1 )
244  END IF
245  kstep = 1
246  ELSE
247 *
248 * 2 x 2 diagonal block
249 *
250 * Invert the diagonal block.
251 *
252  t = a( k, k+1 )
253  ak = a( k, k ) / t
254  akp1 = a( k+1, k+1 ) / t
255  akkp1 = a( k, k+1 ) / t
256  d = t*( ak*akp1-cone )
257  a( k, k ) = akp1 / d
258  a( k+1, k+1 ) = ak / d
259  a( k, k+1 ) = -akkp1 / d
260 *
261 * Compute columns K and K+1 of the inverse.
262 *
263  IF( k.GT.1 ) THEN
264  CALL zcopy( k-1, a( 1, k ), 1, work, 1 )
265  CALL zsymv( uplo, k-1, -cone, a, lda, work, 1, czero,
266  $ a( 1, k ), 1 )
267  a( k, k ) = a( k, k ) - zdotu( k-1, work, 1, a( 1, k ),
268  $ 1 )
269  a( k, k+1 ) = a( k, k+1 ) -
270  $ zdotu( k-1, a( 1, k ), 1, a( 1, k+1 ), 1 )
271  CALL zcopy( k-1, a( 1, k+1 ), 1, work, 1 )
272  CALL zsymv( uplo, k-1, -cone, a, lda, work, 1, czero,
273  $ a( 1, k+1 ), 1 )
274  a( k+1, k+1 ) = a( k+1, k+1 ) -
275  $ zdotu( k-1, work, 1, a( 1, k+1 ), 1 )
276  END IF
277  kstep = 2
278  END IF
279 *
280  IF( kstep.EQ.1 ) THEN
281 *
282 * Interchange rows and columns K and IPIV(K) in the leading
283 * submatrix A(1:k+1,1:k+1)
284 *
285  kp = ipiv( k )
286  IF( kp.NE.k ) THEN
287  IF( kp.GT.1 )
288  $ CALL zswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
289  CALL zswap( k-kp-1, a( kp+1, k ), 1, a( kp, kp+1 ), lda )
290  temp = a( k, k )
291  a( k, k ) = a( kp, kp )
292  a( kp, kp ) = temp
293  END IF
294  ELSE
295 *
296 * Interchange rows and columns K and K+1 with -IPIV(K) and
297 * -IPIV(K+1)in the leading submatrix A(1:k+1,1:k+1)
298 *
299  kp = -ipiv( k )
300  IF( kp.NE.k ) THEN
301  IF( kp.GT.1 )
302  $ CALL zswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
303  CALL zswap( k-kp-1, a( kp+1, k ), 1, a( kp, kp+1 ), lda )
304 *
305  temp = a( k, k )
306  a( k, k ) = a( kp, kp )
307  a( kp, kp ) = temp
308  temp = a( k, k+1 )
309  a( k, k+1 ) = a( kp, k+1 )
310  a( kp, k+1 ) = temp
311  END IF
312 *
313  k = k + 1
314  kp = -ipiv( k )
315  IF( kp.NE.k ) THEN
316  IF( kp.GT.1 )
317  $ CALL zswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
318  CALL zswap( k-kp-1, a( kp+1, k ), 1, a( kp, kp+1 ), lda )
319  temp = a( k, k )
320  a( k, k ) = a( kp, kp )
321  a( kp, kp ) = temp
322  END IF
323  END IF
324 *
325  k = k + 1
326  GO TO 30
327  40 CONTINUE
328 *
329  ELSE
330 *
331 * Compute inv(A) from the factorization A = L*D*L**T.
332 *
333 * K is the main loop index, increasing from 1 to N in steps of
334 * 1 or 2, depending on the size of the diagonal blocks.
335 *
336  k = n
337  50 CONTINUE
338 *
339 * If K < 1, exit from loop.
340 *
341  IF( k.LT.1 )
342  $ GO TO 60
343 *
344  IF( ipiv( k ).GT.0 ) THEN
345 *
346 * 1 x 1 diagonal block
347 *
348 * Invert the diagonal block.
349 *
350  a( k, k ) = cone / a( k, k )
351 *
352 * Compute column K of the inverse.
353 *
354  IF( k.LT.n ) THEN
355  CALL zcopy( n-k, a( k+1, k ), 1, work, 1 )
356  CALL zsymv( uplo, n-k,-cone, a( k+1, k+1 ), lda, work, 1,
357  $ czero, a( k+1, k ), 1 )
358  a( k, k ) = a( k, k ) - zdotu( n-k, work, 1, a( k+1, k ),
359  $ 1 )
360  END IF
361  kstep = 1
362  ELSE
363 *
364 * 2 x 2 diagonal block
365 *
366 * Invert the diagonal block.
367 *
368  t = a( k, k-1 )
369  ak = a( k-1, k-1 ) / t
370  akp1 = a( k, k ) / t
371  akkp1 = a( k, k-1 ) / t
372  d = t*( ak*akp1-cone )
373  a( k-1, k-1 ) = akp1 / d
374  a( k, k ) = ak / d
375  a( k, k-1 ) = -akkp1 / d
376 *
377 * Compute columns K-1 and K of the inverse.
378 *
379  IF( k.LT.n ) THEN
380  CALL zcopy( n-k, a( k+1, k ), 1, work, 1 )
381  CALL zsymv( uplo, n-k,-cone, a( k+1, k+1 ), lda, work, 1,
382  $ czero, a( k+1, k ), 1 )
383  a( k, k ) = a( k, k ) - zdotu( n-k, work, 1, a( k+1, k ),
384  $ 1 )
385  a( k, k-1 ) = a( k, k-1 ) -
386  $ zdotu( n-k, a( k+1, k ), 1, a( k+1, k-1 ),
387  $ 1 )
388  CALL zcopy( n-k, a( k+1, k-1 ), 1, work, 1 )
389  CALL zsymv( uplo, n-k,-cone, a( k+1, k+1 ), lda, work, 1,
390  $ czero, a( k+1, k-1 ), 1 )
391  a( k-1, k-1 ) = a( k-1, k-1 ) -
392  $ zdotu( n-k, work, 1, a( k+1, k-1 ), 1 )
393  END IF
394  kstep = 2
395  END IF
396 *
397  IF( kstep.EQ.1 ) THEN
398 *
399 * Interchange rows and columns K and IPIV(K) in the trailing
400 * submatrix A(k-1:n,k-1:n)
401 *
402  kp = ipiv( k )
403  IF( kp.NE.k ) THEN
404  IF( kp.LT.n )
405  $ CALL zswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
406  CALL zswap( kp-k-1, a( k+1, k ), 1, a( kp, k+1 ), lda )
407  temp = a( k, k )
408  a( k, k ) = a( kp, kp )
409  a( kp, kp ) = temp
410  END IF
411  ELSE
412 *
413 * Interchange rows and columns K and K-1 with -IPIV(K) and
414 * -IPIV(K-1) in the trailing submatrix A(k-1:n,k-1:n)
415 *
416  kp = -ipiv( k )
417  IF( kp.NE.k ) THEN
418  IF( kp.LT.n )
419  $ CALL zswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
420  CALL zswap( kp-k-1, a( k+1, k ), 1, a( kp, k+1 ), lda )
421 *
422  temp = a( k, k )
423  a( k, k ) = a( kp, kp )
424  a( kp, kp ) = temp
425  temp = a( k, k-1 )
426  a( k, k-1 ) = a( kp, k-1 )
427  a( kp, k-1 ) = temp
428  END IF
429 *
430  k = k - 1
431  kp = -ipiv( k )
432  IF( kp.NE.k ) THEN
433  IF( kp.LT.n )
434  $ CALL zswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
435  CALL zswap( kp-k-1, a( k+1, k ), 1, a( kp, k+1 ), lda )
436  temp = a( k, k )
437  a( k, k ) = a( kp, kp )
438  a( kp, kp ) = temp
439  END IF
440  END IF
441 *
442  k = k - 1
443  GO TO 50
444  60 CONTINUE
445  END IF
446 *
447  RETURN
448 *
449 * End of ZSYTRI_ROOK
450 *
zcopy
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:52
zswap
subroutine zswap(N, ZX, INCX, ZY, INCY)
ZSWAP
Definition: zswap.f:52
zdotu
complex *16 function zdotu(N, ZX, INCX, ZY, INCY)
ZDOTU
Definition: zdotu.f:54
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
zsymv
subroutine zsymv(UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
ZSYMV computes a matrix-vector product for a complex symmetric matrix.
Definition: zsymv.f:159
lsame
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function: