LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ csytri()

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

CSYTRI

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

Purpose:
 CSYTRI 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
 CSYTRF.
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 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 CSYTRF.

          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 CSYTRF.
[out]WORK
          WORK is COMPLEX 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 csytri.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 A( LDA, * ), WORK( * )
126 * ..
127 *
128 * =====================================================================
129 *
130 * .. Parameters ..
131  COMPLEX ONE, ZERO
132  parameter( one = ( 1.0e+0, 0.0e+0 ),
133  $ zero = ( 0.0e+0, 0.0e+0 ) )
134 * ..
135 * .. Local Scalars ..
136  LOGICAL UPPER
137  INTEGER K, KP, KSTEP
138  COMPLEX AK, AKKP1, AKP1, D, T, TEMP
139 * ..
140 * .. External Functions ..
141  LOGICAL LSAME
142  COMPLEX CDOTU
143  EXTERNAL lsame, cdotu
144 * ..
145 * .. External Subroutines ..
146  EXTERNAL ccopy, cswap, csymv, xerbla
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( 'CSYTRI', -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 ccopy( k-1, a( 1, k ), 1, work, 1 )
222  CALL csymv( uplo, k-1, -one, a, lda, work, 1, zero,
223  $ a( 1, k ), 1 )
224  a( k, k ) = a( k, k ) - cdotu( 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 ccopy( k-1, a( 1, k ), 1, work, 1 )
247  CALL csymv( uplo, k-1, -one, a, lda, work, 1, zero,
248  $ a( 1, k ), 1 )
249  a( k, k ) = a( k, k ) - cdotu( k-1, work, 1, a( 1, k ),
250  $ 1 )
251  a( k, k+1 ) = a( k, k+1 ) -
252  $ cdotu( k-1, a( 1, k ), 1, a( 1, k+1 ), 1 )
253  CALL ccopy( k-1, a( 1, k+1 ), 1, work, 1 )
254  CALL csymv( 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  $ cdotu( 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 cswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
269  CALL cswap( 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 ccopy( n-k, a( k+1, k ), 1, work, 1 )
311  CALL csymv( 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 ) - cdotu( 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 ccopy( n-k, a( k+1, k ), 1, work, 1 )
336  CALL csymv( 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 ) - cdotu( n-k, work, 1, a( k+1, k ),
339  $ 1 )
340  a( k, k-1 ) = a( k, k-1 ) -
341  $ cdotu( n-k, a( k+1, k ), 1, a( k+1, k-1 ),
342  $ 1 )
343  CALL ccopy( n-k, a( k+1, k-1 ), 1, work, 1 )
344  CALL csymv( 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  $ cdotu( 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 cswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
360  CALL cswap( 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 CSYTRI
379 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
Definition: ccopy.f:81
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
Definition: cswap.f:81
complex function cdotu(N, CX, INCX, CY, INCY)
CDOTU
Definition: cdotu.f:83
subroutine csymv(UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
CSYMV computes a matrix-vector product for a complex symmetric matrix.
Definition: csymv.f:157
Here is the call graph for this function:
Here is the caller graph for this function: