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

SSYTRI_ROOK

Purpose:
``` SSYTRI_ROOK computes the inverse of a real symmetric
matrix A using the factorization A = U*D*U**T or A = L*D*L**T
computed by SSYTRF_ROOK.```
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 REAL 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 SSYTRF_ROOK. 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 SSYTRF_ROOK.``` [out] WORK ` WORK is REAL 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.```
Date
April 2012
Contributors:
```   April 2012, Igor Kozachenko,
Computer Science Division,
University of California, Berkeley

September 2007, Sven Hammarling, Nicholas J. Higham, Craig Lucas,
School of Mathematics,
University of Manchester```

Definition at line 131 of file ssytri_rook.f.

131 *
132 * -- LAPACK computational routine (version 3.4.1) --
133 * -- LAPACK is a software package provided by Univ. of Tennessee, --
134 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
135 * April 2012
136 *
137 * .. Scalar Arguments ..
138  CHARACTER uplo
139  INTEGER info, lda, n
140 * ..
141 * .. Array Arguments ..
142  INTEGER ipiv( * )
143  REAL a( lda, * ), work( * )
144 * ..
145 *
146 * =====================================================================
147 *
148 * .. Parameters ..
149  REAL one, zero
150  parameter ( one = 1.0e+0, zero = 0.0e+0 )
151 * ..
152 * .. Local Scalars ..
153  LOGICAL upper
154  INTEGER k, kp, kstep
155  REAL ak, akkp1, akp1, d, t, temp
156 * ..
157 * .. External Functions ..
158  LOGICAL lsame
159  REAL sdot
160  EXTERNAL lsame, sdot
161 * ..
162 * .. External Subroutines ..
163  EXTERNAL scopy, sswap, ssymv, xerbla
164 * ..
165 * .. Intrinsic Functions ..
166  INTRINSIC abs, max
167 * ..
168 * .. Executable Statements ..
169 *
170 * Test the input parameters.
171 *
172  info = 0
173  upper = lsame( uplo, 'U' )
174  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
175  info = -1
176  ELSE IF( n.LT.0 ) THEN
177  info = -2
178  ELSE IF( lda.LT.max( 1, n ) ) THEN
179  info = -4
180  END IF
181  IF( info.NE.0 ) THEN
182  CALL xerbla( 'SSYTRI_ROOK', -info )
183  RETURN
184  END IF
185 *
186 * Quick return if possible
187 *
188  IF( n.EQ.0 )
189  \$ RETURN
190 *
191 * Check that the diagonal matrix D is nonsingular.
192 *
193  IF( upper ) THEN
194 *
195 * Upper triangular storage: examine D from bottom to top
196 *
197  DO 10 info = n, 1, -1
198  IF( ipiv( info ).GT.0 .AND. a( info, info ).EQ.zero )
199  \$ RETURN
200  10 CONTINUE
201  ELSE
202 *
203 * Lower triangular storage: examine D from top to bottom.
204 *
205  DO 20 info = 1, n
206  IF( ipiv( info ).GT.0 .AND. a( info, info ).EQ.zero )
207  \$ RETURN
208  20 CONTINUE
209  END IF
210  info = 0
211 *
212  IF( upper ) THEN
213 *
214 * Compute inv(A) from the factorization A = U*D*U**T.
215 *
216 * K is the main loop index, increasing from 1 to N in steps of
217 * 1 or 2, depending on the size of the diagonal blocks.
218 *
219  k = 1
220  30 CONTINUE
221 *
222 * If K > N, exit from loop.
223 *
224  IF( k.GT.n )
225  \$ GO TO 40
226 *
227  IF( ipiv( k ).GT.0 ) THEN
228 *
229 * 1 x 1 diagonal block
230 *
231 * Invert the diagonal block.
232 *
233  a( k, k ) = one / a( k, k )
234 *
235 * Compute column K of the inverse.
236 *
237  IF( k.GT.1 ) THEN
238  CALL scopy( k-1, a( 1, k ), 1, work, 1 )
239  CALL ssymv( uplo, k-1, -one, a, lda, work, 1, zero,
240  \$ a( 1, k ), 1 )
241  a( k, k ) = a( k, k ) - sdot( k-1, work, 1, a( 1, k ),
242  \$ 1 )
243  END IF
244  kstep = 1
245  ELSE
246 *
247 * 2 x 2 diagonal block
248 *
249 * Invert the diagonal block.
250 *
251  t = abs( a( k, k+1 ) )
252  ak = a( k, k ) / t
253  akp1 = a( k+1, k+1 ) / t
254  akkp1 = a( k, k+1 ) / t
255  d = t*( ak*akp1-one )
256  a( k, k ) = akp1 / d
257  a( k+1, k+1 ) = ak / d
258  a( k, k+1 ) = -akkp1 / d
259 *
260 * Compute columns K and K+1 of the inverse.
261 *
262  IF( k.GT.1 ) THEN
263  CALL scopy( k-1, a( 1, k ), 1, work, 1 )
264  CALL ssymv( uplo, k-1, -one, a, lda, work, 1, zero,
265  \$ a( 1, k ), 1 )
266  a( k, k ) = a( k, k ) - sdot( k-1, work, 1, a( 1, k ),
267  \$ 1 )
268  a( k, k+1 ) = a( k, k+1 ) -
269  \$ sdot( k-1, a( 1, k ), 1, a( 1, k+1 ), 1 )
270  CALL scopy( k-1, a( 1, k+1 ), 1, work, 1 )
271  CALL ssymv( uplo, k-1, -one, a, lda, work, 1, zero,
272  \$ a( 1, k+1 ), 1 )
273  a( k+1, k+1 ) = a( k+1, k+1 ) -
274  \$ sdot( k-1, work, 1, a( 1, k+1 ), 1 )
275  END IF
276  kstep = 2
277  END IF
278 *
279  IF( kstep.EQ.1 ) THEN
280 *
281 * Interchange rows and columns K and IPIV(K) in the leading
282 * submatrix A(1:k+1,1:k+1)
283 *
284  kp = ipiv( k )
285  IF( kp.NE.k ) THEN
286  IF( kp.GT.1 )
287  \$ CALL sswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
288  CALL sswap( k-kp-1, a( kp+1, k ), 1, a( kp, kp+1 ), lda )
289  temp = a( k, k )
290  a( k, k ) = a( kp, kp )
291  a( kp, kp ) = temp
292  END IF
293  ELSE
294 *
295 * Interchange rows and columns K and K+1 with -IPIV(K) and
296 * -IPIV(K+1)in the leading submatrix A(1:k+1,1:k+1)
297 *
298  kp = -ipiv( k )
299  IF( kp.NE.k ) THEN
300  IF( kp.GT.1 )
301  \$ CALL sswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
302  CALL sswap( k-kp-1, a( kp+1, k ), 1, a( kp, kp+1 ), lda )
303 *
304  temp = a( k, k )
305  a( k, k ) = a( kp, kp )
306  a( kp, kp ) = temp
307  temp = a( k, k+1 )
308  a( k, k+1 ) = a( kp, k+1 )
309  a( kp, k+1 ) = temp
310  END IF
311 *
312  k = k + 1
313  kp = -ipiv( k )
314  IF( kp.NE.k ) THEN
315  IF( kp.GT.1 )
316  \$ CALL sswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
317  CALL sswap( k-kp-1, a( kp+1, k ), 1, a( kp, kp+1 ), lda )
318  temp = a( k, k )
319  a( k, k ) = a( kp, kp )
320  a( kp, kp ) = temp
321  END IF
322  END IF
323 *
324  k = k + 1
325  GO TO 30
326  40 CONTINUE
327 *
328  ELSE
329 *
330 * Compute inv(A) from the factorization A = L*D*L**T.
331 *
332 * K is the main loop index, increasing from 1 to N in steps of
333 * 1 or 2, depending on the size of the diagonal blocks.
334 *
335  k = n
336  50 CONTINUE
337 *
338 * If K < 1, exit from loop.
339 *
340  IF( k.LT.1 )
341  \$ GO TO 60
342 *
343  IF( ipiv( k ).GT.0 ) THEN
344 *
345 * 1 x 1 diagonal block
346 *
347 * Invert the diagonal block.
348 *
349  a( k, k ) = one / a( k, k )
350 *
351 * Compute column K of the inverse.
352 *
353  IF( k.LT.n ) THEN
354  CALL scopy( n-k, a( k+1, k ), 1, work, 1 )
355  CALL ssymv( uplo, n-k, -one, a( k+1, k+1 ), lda, work, 1,
356  \$ zero, a( k+1, k ), 1 )
357  a( k, k ) = a( k, k ) - sdot( n-k, work, 1, a( k+1, k ),
358  \$ 1 )
359  END IF
360  kstep = 1
361  ELSE
362 *
363 * 2 x 2 diagonal block
364 *
365 * Invert the diagonal block.
366 *
367  t = abs( a( k, k-1 ) )
368  ak = a( k-1, k-1 ) / t
369  akp1 = a( k, k ) / t
370  akkp1 = a( k, k-1 ) / t
371  d = t*( ak*akp1-one )
372  a( k-1, k-1 ) = akp1 / d
373  a( k, k ) = ak / d
374  a( k, k-1 ) = -akkp1 / d
375 *
376 * Compute columns K-1 and K of the inverse.
377 *
378  IF( k.LT.n ) THEN
379  CALL scopy( n-k, a( k+1, k ), 1, work, 1 )
380  CALL ssymv( uplo, n-k, -one, a( k+1, k+1 ), lda, work, 1,
381  \$ zero, a( k+1, k ), 1 )
382  a( k, k ) = a( k, k ) - sdot( n-k, work, 1, a( k+1, k ),
383  \$ 1 )
384  a( k, k-1 ) = a( k, k-1 ) -
385  \$ sdot( n-k, a( k+1, k ), 1, a( k+1, k-1 ),
386  \$ 1 )
387  CALL scopy( n-k, a( k+1, k-1 ), 1, work, 1 )
388  CALL ssymv( uplo, n-k, -one, a( k+1, k+1 ), lda, work, 1,
389  \$ zero, a( k+1, k-1 ), 1 )
390  a( k-1, k-1 ) = a( k-1, k-1 ) -
391  \$ sdot( n-k, work, 1, a( k+1, k-1 ), 1 )
392  END IF
393  kstep = 2
394  END IF
395 *
396  IF( kstep.EQ.1 ) THEN
397 *
398 * Interchange rows and columns K and IPIV(K) in the trailing
399 * submatrix A(k-1:n,k-1:n)
400 *
401  kp = ipiv( k )
402  IF( kp.NE.k ) THEN
403  IF( kp.LT.n )
404  \$ CALL sswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
405  CALL sswap( kp-k-1, a( k+1, k ), 1, a( kp, k+1 ), lda )
406  temp = a( k, k )
407  a( k, k ) = a( kp, kp )
408  a( kp, kp ) = temp
409  END IF
410  ELSE
411 *
412 * Interchange rows and columns K and K-1 with -IPIV(K) and
413 * -IPIV(K-1) in the trailing submatrix A(k-1:n,k-1:n)
414 *
415  kp = -ipiv( k )
416  IF( kp.NE.k ) THEN
417  IF( kp.LT.n )
418  \$ CALL sswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
419  CALL sswap( kp-k-1, a( k+1, k ), 1, a( kp, k+1 ), lda )
420 *
421  temp = a( k, k )
422  a( k, k ) = a( kp, kp )
423  a( kp, kp ) = temp
424  temp = a( k, k-1 )
425  a( k, k-1 ) = a( kp, k-1 )
426  a( kp, k-1 ) = temp
427  END IF
428 *
429  k = k - 1
430  kp = -ipiv( k )
431  IF( kp.NE.k ) THEN
432  IF( kp.LT.n )
433  \$ CALL sswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
434  CALL sswap( kp-k-1, a( k+1, k ), 1, a( kp, k+1 ), lda )
435  temp = a( k, k )
436  a( k, k ) = a( kp, kp )
437  a( kp, kp ) = temp
438  END IF
439  END IF
440 *
441  k = k - 1
442  GO TO 50
443  60 CONTINUE
444  END IF
445 *
446  RETURN
447 *
448 * End of SSYTRI_ROOK
449 *
real function sdot(N, SX, INCX, SY, INCY)
SDOT
Definition: sdot.f:53
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine sswap(N, SX, INCX, SY, INCY)
SSWAP
Definition: sswap.f:53
subroutine ssymv(UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
SSYMV
Definition: ssymv.f:154
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:53

Here is the call graph for this function:

Here is the caller graph for this function: