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

ZHETRI

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

Purpose: 
 ZHETRI computes the inverse of a complex Hermitian indefinite matrix
 A using the factorization A = U*D*U**H or A = L*D*L**H computed by
 ZHETRF.
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**H;
          = 'L':  Lower triangular, form is A = L*D*L**H.
[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 ZHETRF.

          On exit, if INFO = 0, the (Hermitian) 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 ZHETRF.
[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 2011

Definition at line 116 of file zhetri.f.

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