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

CHETRI

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

Purpose: 
 CHETRI 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
 CHETRF.
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 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 CHETRF.

          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 CHETRF.
[out]WORK
          WORK is COMPLEX 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 chetri.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 a( lda, * ), work( * )
129 * ..
130 *
131 * =====================================================================
132 *
133 * .. Parameters ..
134  REAL one
135  COMPLEX cone, zero
136  parameter ( one = 1.0e+0, cone = ( 1.0e+0, 0.0e+0 ),
137  $ zero = ( 0.0e+0, 0.0e+0 ) )
138 * ..
139 * .. Local Scalars ..
140  LOGICAL upper
141  INTEGER j, k, kp, kstep
142  REAL ak, akp1, d, t
143  COMPLEX akkp1, temp
144 * ..
145 * .. External Functions ..
146  LOGICAL lsame
147  COMPLEX cdotc
148  EXTERNAL lsame, cdotc
149 * ..
150 * .. External Subroutines ..
151  EXTERNAL ccopy, chemv, cswap, xerbla
152 * ..
153 * .. Intrinsic Functions ..
154  INTRINSIC abs, conjg, max, real
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( 'CHETRI', -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 / REAL( A( K, K ) )
222 *
223 * Compute column K of the inverse.
224 *
225  IF( k.GT.1 ) THEN
226  CALL ccopy( k-1, a( 1, k ), 1, work, 1 )
227  CALL chemv( uplo, k-1, -cone, a, lda, work, 1, zero,
228  $ a( 1, k ), 1 )
229  a( k, k ) = a( k, k ) - REAL( CDOTC( K-1, WORK, 1, A( 1, $ K ), 1 ) )
230  END IF
231  kstep = 1
232  ELSE
233 *
234 * 2 x 2 diagonal block
235 *
236 * Invert the diagonal block.
237 *
238  t = abs( a( k, k+1 ) )
239  ak = REAL( A( K, K ) ) / t
240  akp1 = REAL( A( K+1, K+1 ) ) / t
241  akkp1 = a( k, k+1 ) / t
242  d = t*( ak*akp1-one )
243  a( k, k ) = akp1 / d
244  a( k+1, k+1 ) = ak / d
245  a( k, k+1 ) = -akkp1 / d
246 *
247 * Compute columns K and K+1 of the inverse.
248 *
249  IF( k.GT.1 ) THEN
250  CALL ccopy( k-1, a( 1, k ), 1, work, 1 )
251  CALL chemv( uplo, k-1, -cone, a, lda, work, 1, zero,
252  $ a( 1, k ), 1 )
253  a( k, k ) = a( k, k ) - REAL( CDOTC( K-1, WORK, 1, A( 1, $ K ), 1 ) )
254  a( k, k+1 ) = a( k, k+1 ) -
255  $ cdotc( k-1, a( 1, k ), 1, a( 1, k+1 ), 1 )
256  CALL ccopy( k-1, a( 1, k+1 ), 1, work, 1 )
257  CALL chemv( uplo, k-1, -cone, a, lda, work, 1, zero,
258  $ a( 1, k+1 ), 1 )
259  a( k+1, k+1 ) = a( k+1, k+1 ) -
260  $ REAL( CDOTC( K-1, WORK, 1, A( 1, K+1 ), $ 1 ) )
261  END IF
262  kstep = 2
263  END IF
264 *
265  kp = abs( ipiv( k ) )
266  IF( kp.NE.k ) THEN
267 *
268 * Interchange rows and columns K and KP in the leading
269 * submatrix A(1:k+1,1:k+1)
270 *
271  CALL cswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
272  DO 40 j = kp + 1, k - 1
273  temp = conjg( a( j, k ) )
274  a( j, k ) = conjg( a( kp, j ) )
275  a( kp, j ) = temp
276  40 CONTINUE
277  a( kp, k ) = conjg( a( kp, k ) )
278  temp = a( k, k )
279  a( k, k ) = a( kp, kp )
280  a( kp, kp ) = temp
281  IF( kstep.EQ.2 ) THEN
282  temp = a( k, k+1 )
283  a( k, k+1 ) = a( kp, k+1 )
284  a( kp, k+1 ) = temp
285  END IF
286  END IF
287 *
288  k = k + kstep
289  GO TO 30
290  50 CONTINUE
291 *
292  ELSE
293 *
294 * Compute inv(A) from the factorization A = L*D*L**H.
295 *
296 * K is the main loop index, increasing from 1 to N in steps of
297 * 1 or 2, depending on the size of the diagonal blocks.
298 *
299  k = n
300  60 CONTINUE
301 *
302 * If K < 1, exit from loop.
303 *
304  IF( k.LT.1 )
305  $ GO TO 80
306 *
307  IF( ipiv( k ).GT.0 ) THEN
308 *
309 * 1 x 1 diagonal block
310 *
311 * Invert the diagonal block.
312 *
313  a( k, k ) = one / REAL( A( K, K ) )
314 *
315 * Compute column K of the inverse.
316 *
317  IF( k.LT.n ) THEN
318  CALL ccopy( n-k, a( k+1, k ), 1, work, 1 )
319  CALL chemv( uplo, n-k, -cone, a( k+1, k+1 ), lda, work,
320  $ 1, zero, a( k+1, k ), 1 )
321  a( k, k ) = a( k, k ) - REAL( CDOTC( N-K, WORK, 1, $ A( K+1, K ), 1 ) )
322  END IF
323  kstep = 1
324  ELSE
325 *
326 * 2 x 2 diagonal block
327 *
328 * Invert the diagonal block.
329 *
330  t = abs( a( k, k-1 ) )
331  ak = REAL( A( K-1, K-1 ) ) / t
332  akp1 = REAL( A( K, K ) ) / t
333  akkp1 = a( k, k-1 ) / t
334  d = t*( ak*akp1-one )
335  a( k-1, k-1 ) = akp1 / d
336  a( k, k ) = ak / d
337  a( k, k-1 ) = -akkp1 / d
338 *
339 * Compute columns K-1 and K of the inverse.
340 *
341  IF( k.LT.n ) THEN
342  CALL ccopy( n-k, a( k+1, k ), 1, work, 1 )
343  CALL chemv( uplo, n-k, -cone, a( k+1, k+1 ), lda, work,
344  $ 1, zero, a( k+1, k ), 1 )
345  a( k, k ) = a( k, k ) - REAL( CDOTC( N-K, WORK, 1, $ A( K+1, K ), 1 ) )
346  a( k, k-1 ) = a( k, k-1 ) -
347  $ cdotc( n-k, a( k+1, k ), 1, a( k+1, k-1 ),
348  $ 1 )
349  CALL ccopy( n-k, a( k+1, k-1 ), 1, work, 1 )
350  CALL chemv( uplo, n-k, -cone, a( k+1, k+1 ), lda, work,
351  $ 1, zero, a( k+1, k-1 ), 1 )
352  a( k-1, k-1 ) = a( k-1, k-1 ) -
353  $ REAL( CDOTC( N-K, WORK, 1, A( K+1, K-1 ), $ 1 ) )
354  END IF
355  kstep = 2
356  END IF
357 *
358  kp = abs( ipiv( k ) )
359  IF( kp.NE.k ) THEN
360 *
361 * Interchange rows and columns K and KP in the trailing
362 * submatrix A(k-1:n,k-1:n)
363 *
364  IF( kp.LT.n )
365  $ CALL cswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
366  DO 70 j = k + 1, kp - 1
367  temp = conjg( a( j, k ) )
368  a( j, k ) = conjg( a( kp, j ) )
369  a( kp, j ) = temp
370  70 CONTINUE
371  a( kp, k ) = conjg( a( kp, k ) )
372  temp = a( k, k )
373  a( k, k ) = a( kp, kp )
374  a( kp, kp ) = temp
375  IF( kstep.EQ.2 ) THEN
376  temp = a( k, k-1 )
377  a( k, k-1 ) = a( kp, k-1 )
378  a( kp, k-1 ) = temp
379  END IF
380  END IF
381 *
382  k = k - kstep
383  GO TO 60
384  80 CONTINUE
385  END IF
386 *
387  RETURN
388 *
389 * End of CHETRI
390 *
391 
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
cdotc
complex function cdotc(N, CX, INCX, CY, INCY)
CDOTC
Definition: cdotc.f:54
chemv
subroutine chemv(UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
CHEMV
Definition: chemv.f:156
ccopy
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
Definition: ccopy.f:52
cswap
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
Definition: cswap.f:52
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: