LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ 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: