LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ zsytri()

subroutine zsytri ( character  UPLO,
integer  N,
complex*16, dimension( lda, * )  A,
integer  LDA,
integer, dimension( * )  IPIV,
complex*16, dimension( * )  WORK,
integer  INFO 
)

ZSYTRI

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

Purpose:
 ZSYTRI computes the inverse of a complex symmetric indefinite matrix
 A using the factorization A = U*D*U**T or A = L*D*L**T computed by
 ZSYTRF.
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.

          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.
[out]WORK
          WORK is COMPLEX*16 array, dimension (2*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.

Definition at line 113 of file zsytri.f.

114 *
115 * -- LAPACK computational routine --
116 * -- LAPACK is a software package provided by Univ. of Tennessee, --
117 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
118 *
119 * .. Scalar Arguments ..
120  CHARACTER UPLO
121  INTEGER INFO, LDA, N
122 * ..
123 * .. Array Arguments ..
124  INTEGER IPIV( * )
125  COMPLEX*16 A( LDA, * ), WORK( * )
126 * ..
127 *
128 * =====================================================================
129 *
130 * .. Parameters ..
131  COMPLEX*16 ONE, ZERO
132  parameter( one = ( 1.0d+0, 0.0d+0 ),
133  $ zero = ( 0.0d+0, 0.0d+0 ) )
134 * ..
135 * .. Local Scalars ..
136  LOGICAL UPPER
137  INTEGER K, KP, KSTEP
138  COMPLEX*16 AK, AKKP1, AKP1, D, T, TEMP
139 * ..
140 * .. External Functions ..
141  LOGICAL LSAME
142  COMPLEX*16 ZDOTU
143  EXTERNAL lsame, zdotu
144 * ..
145 * .. External Subroutines ..
146  EXTERNAL xerbla, zcopy, zswap, zsymv
147 * ..
148 * .. Intrinsic Functions ..
149  INTRINSIC abs, max
150 * ..
151 * .. Executable Statements ..
152 *
153 * Test the input parameters.
154 *
155  info = 0
156  upper = lsame( uplo, 'U' )
157  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
158  info = -1
159  ELSE IF( n.LT.0 ) THEN
160  info = -2
161  ELSE IF( lda.LT.max( 1, n ) ) THEN
162  info = -4
163  END IF
164  IF( info.NE.0 ) THEN
165  CALL xerbla( 'ZSYTRI', -info )
166  RETURN
167  END IF
168 *
169 * Quick return if possible
170 *
171  IF( n.EQ.0 )
172  $ RETURN
173 *
174 * Check that the diagonal matrix D is nonsingular.
175 *
176  IF( upper ) THEN
177 *
178 * Upper triangular storage: examine D from bottom to top
179 *
180  DO 10 info = n, 1, -1
181  IF( ipiv( info ).GT.0 .AND. a( info, info ).EQ.zero )
182  $ RETURN
183  10 CONTINUE
184  ELSE
185 *
186 * Lower triangular storage: examine D from top to bottom.
187 *
188  DO 20 info = 1, n
189  IF( ipiv( info ).GT.0 .AND. a( info, info ).EQ.zero )
190  $ RETURN
191  20 CONTINUE
192  END IF
193  info = 0
194 *
195  IF( upper ) THEN
196 *
197 * Compute inv(A) from the factorization A = U*D*U**T.
198 *
199 * K is the main loop index, increasing from 1 to N in steps of
200 * 1 or 2, depending on the size of the diagonal blocks.
201 *
202  k = 1
203  30 CONTINUE
204 *
205 * If K > N, exit from loop.
206 *
207  IF( k.GT.n )
208  $ GO TO 40
209 *
210  IF( ipiv( k ).GT.0 ) THEN
211 *
212 * 1 x 1 diagonal block
213 *
214 * Invert the diagonal block.
215 *
216  a( k, k ) = one / a( k, k )
217 *
218 * Compute column K of the inverse.
219 *
220  IF( k.GT.1 ) THEN
221  CALL zcopy( k-1, a( 1, k ), 1, work, 1 )
222  CALL zsymv( uplo, k-1, -one, a, lda, work, 1, zero,
223  $ a( 1, k ), 1 )
224  a( k, k ) = a( k, k ) - zdotu( k-1, work, 1, a( 1, k ),
225  $ 1 )
226  END IF
227  kstep = 1
228  ELSE
229 *
230 * 2 x 2 diagonal block
231 *
232 * Invert the diagonal block.
233 *
234  t = a( k, k+1 )
235  ak = a( k, k ) / t
236  akp1 = a( k+1, k+1 ) / t
237  akkp1 = a( k, k+1 ) / t
238  d = t*( ak*akp1-one )
239  a( k, k ) = akp1 / d
240  a( k+1, k+1 ) = ak / d
241  a( k, k+1 ) = -akkp1 / d
242 *
243 * Compute columns K and K+1 of the inverse.
244 *
245  IF( k.GT.1 ) THEN
246  CALL zcopy( k-1, a( 1, k ), 1, work, 1 )
247  CALL zsymv( uplo, k-1, -one, a, lda, work, 1, zero,
248  $ a( 1, k ), 1 )
249  a( k, k ) = a( k, k ) - zdotu( k-1, work, 1, a( 1, k ),
250  $ 1 )
251  a( k, k+1 ) = a( k, k+1 ) -
252  $ zdotu( k-1, a( 1, k ), 1, a( 1, k+1 ), 1 )
253  CALL zcopy( k-1, a( 1, k+1 ), 1, work, 1 )
254  CALL zsymv( uplo, k-1, -one, a, lda, work, 1, zero,
255  $ a( 1, k+1 ), 1 )
256  a( k+1, k+1 ) = a( k+1, k+1 ) -
257  $ zdotu( k-1, work, 1, a( 1, k+1 ), 1 )
258  END IF
259  kstep = 2
260  END IF
261 *
262  kp = abs( ipiv( k ) )
263  IF( kp.NE.k ) THEN
264 *
265 * Interchange rows and columns K and KP in the leading
266 * submatrix A(1:k+1,1:k+1)
267 *
268  CALL zswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
269  CALL zswap( k-kp-1, a( kp+1, k ), 1, a( kp, kp+1 ), lda )
270  temp = a( k, k )
271  a( k, k ) = a( kp, kp )
272  a( kp, kp ) = temp
273  IF( kstep.EQ.2 ) THEN
274  temp = a( k, k+1 )
275  a( k, k+1 ) = a( kp, k+1 )
276  a( kp, k+1 ) = temp
277  END IF
278  END IF
279 *
280  k = k + kstep
281  GO TO 30
282  40 CONTINUE
283 *
284  ELSE
285 *
286 * Compute inv(A) from the factorization A = L*D*L**T.
287 *
288 * K is the main loop index, increasing from 1 to N in steps of
289 * 1 or 2, depending on the size of the diagonal blocks.
290 *
291  k = n
292  50 CONTINUE
293 *
294 * If K < 1, exit from loop.
295 *
296  IF( k.LT.1 )
297  $ GO TO 60
298 *
299  IF( ipiv( k ).GT.0 ) THEN
300 *
301 * 1 x 1 diagonal block
302 *
303 * Invert the diagonal block.
304 *
305  a( k, k ) = one / a( k, k )
306 *
307 * Compute column K of the inverse.
308 *
309  IF( k.LT.n ) THEN
310  CALL zcopy( n-k, a( k+1, k ), 1, work, 1 )
311  CALL zsymv( uplo, n-k, -one, a( k+1, k+1 ), lda, work, 1,
312  $ zero, a( k+1, k ), 1 )
313  a( k, k ) = a( k, k ) - zdotu( n-k, work, 1, a( k+1, k ),
314  $ 1 )
315  END IF
316  kstep = 1
317  ELSE
318 *
319 * 2 x 2 diagonal block
320 *
321 * Invert the diagonal block.
322 *
323  t = a( k, k-1 )
324  ak = a( k-1, k-1 ) / t
325  akp1 = a( k, k ) / t
326  akkp1 = a( k, k-1 ) / t
327  d = t*( ak*akp1-one )
328  a( k-1, k-1 ) = akp1 / d
329  a( k, k ) = ak / d
330  a( k, k-1 ) = -akkp1 / d
331 *
332 * Compute columns K-1 and K of the inverse.
333 *
334  IF( k.LT.n ) THEN
335  CALL zcopy( n-k, a( k+1, k ), 1, work, 1 )
336  CALL zsymv( uplo, n-k, -one, a( k+1, k+1 ), lda, work, 1,
337  $ zero, a( k+1, k ), 1 )
338  a( k, k ) = a( k, k ) - zdotu( n-k, work, 1, a( k+1, k ),
339  $ 1 )
340  a( k, k-1 ) = a( k, k-1 ) -
341  $ zdotu( n-k, a( k+1, k ), 1, a( k+1, k-1 ),
342  $ 1 )
343  CALL zcopy( n-k, a( k+1, k-1 ), 1, work, 1 )
344  CALL zsymv( uplo, n-k, -one, a( k+1, k+1 ), lda, work, 1,
345  $ zero, a( k+1, k-1 ), 1 )
346  a( k-1, k-1 ) = a( k-1, k-1 ) -
347  $ zdotu( n-k, work, 1, a( k+1, k-1 ), 1 )
348  END IF
349  kstep = 2
350  END IF
351 *
352  kp = abs( ipiv( k ) )
353  IF( kp.NE.k ) THEN
354 *
355 * Interchange rows and columns K and KP in the trailing
356 * submatrix A(k-1:n,k-1:n)
357 *
358  IF( kp.LT.n )
359  $ CALL zswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
360  CALL zswap( kp-k-1, a( k+1, k ), 1, a( kp, k+1 ), lda )
361  temp = a( k, k )
362  a( k, k ) = a( kp, kp )
363  a( kp, kp ) = temp
364  IF( kstep.EQ.2 ) THEN
365  temp = a( k, k-1 )
366  a( k, k-1 ) = a( kp, k-1 )
367  a( kp, k-1 ) = temp
368  END IF
369  END IF
370 *
371  k = k - kstep
372  GO TO 50
373  60 CONTINUE
374  END IF
375 *
376  RETURN
377 *
378 * End of ZSYTRI
379 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine zswap(N, ZX, INCX, ZY, INCY)
ZSWAP
Definition: zswap.f:81
complex *16 function zdotu(N, ZX, INCX, ZY, INCY)
ZDOTU
Definition: zdotu.f:83
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:81
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:157
Here is the call graph for this function:
Here is the caller graph for this function: