LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ sgetri()

subroutine sgetri ( integer  N,
real, dimension( lda, * )  A,
integer  LDA,
integer, dimension( * )  IPIV,
real, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

SGETRI

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

Purpose:
 SGETRI computes the inverse of a matrix using the LU factorization
 computed by SGETRF.

 This method inverts U and then computes inv(A) by solving the system
 inv(A)*L = inv(U) for inv(A).
Parameters
[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 factors L and U from the factorization
          A = P*L*U as computed by SGETRF.
          On exit, if INFO = 0, the inverse of the original matrix A.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[in]IPIV
          IPIV is INTEGER array, dimension (N)
          The pivot indices from SGETRF; for 1<=i<=N, row i of the
          matrix was interchanged with row IPIV(i).
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO=0, then WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.  LWORK >= max(1,N).
          For optimal performance LWORK >= N*NB, where NB is
          the optimal blocksize returned by ILAENV.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, U(i,i) is exactly zero; 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 sgetri.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  INTEGER INFO, LDA, LWORK, N
121 * ..
122 * .. Array Arguments ..
123  INTEGER IPIV( * )
124  REAL A( LDA, * ), WORK( * )
125 * ..
126 *
127 * =====================================================================
128 *
129 * .. Parameters ..
130  REAL ZERO, ONE
131  parameter( zero = 0.0e+0, one = 1.0e+0 )
132 * ..
133 * .. Local Scalars ..
134  LOGICAL LQUERY
135  INTEGER I, IWS, J, JB, JJ, JP, LDWORK, LWKOPT, NB,
136  $ NBMIN, NN
137 * ..
138 * .. External Functions ..
139  INTEGER ILAENV
140  EXTERNAL ilaenv
141 * ..
142 * .. External Subroutines ..
143  EXTERNAL sgemm, sgemv, sswap, strsm, strtri, xerbla
144 * ..
145 * .. Intrinsic Functions ..
146  INTRINSIC max, min
147 * ..
148 * .. Executable Statements ..
149 *
150 * Test the input parameters.
151 *
152  info = 0
153  nb = ilaenv( 1, 'SGETRI', ' ', n, -1, -1, -1 )
154  lwkopt = n*nb
155  work( 1 ) = lwkopt
156  lquery = ( lwork.EQ.-1 )
157  IF( n.LT.0 ) THEN
158  info = -1
159  ELSE IF( lda.LT.max( 1, n ) ) THEN
160  info = -3
161  ELSE IF( lwork.LT.max( 1, n ) .AND. .NOT.lquery ) THEN
162  info = -6
163  END IF
164  IF( info.NE.0 ) THEN
165  CALL xerbla( 'SGETRI', -info )
166  RETURN
167  ELSE IF( lquery ) THEN
168  RETURN
169  END IF
170 *
171 * Quick return if possible
172 *
173  IF( n.EQ.0 )
174  $ RETURN
175 *
176 * Form inv(U). If INFO > 0 from STRTRI, then U is singular,
177 * and the inverse is not computed.
178 *
179  CALL strtri( 'Upper', 'Non-unit', n, a, lda, info )
180  IF( info.GT.0 )
181  $ RETURN
182 *
183  nbmin = 2
184  ldwork = n
185  IF( nb.GT.1 .AND. nb.LT.n ) THEN
186  iws = max( ldwork*nb, 1 )
187  IF( lwork.LT.iws ) THEN
188  nb = lwork / ldwork
189  nbmin = max( 2, ilaenv( 2, 'SGETRI', ' ', n, -1, -1, -1 ) )
190  END IF
191  ELSE
192  iws = n
193  END IF
194 *
195 * Solve the equation inv(A)*L = inv(U) for inv(A).
196 *
197  IF( nb.LT.nbmin .OR. nb.GE.n ) THEN
198 *
199 * Use unblocked code.
200 *
201  DO 20 j = n, 1, -1
202 *
203 * Copy current column of L to WORK and replace with zeros.
204 *
205  DO 10 i = j + 1, n
206  work( i ) = a( i, j )
207  a( i, j ) = zero
208  10 CONTINUE
209 *
210 * Compute current column of inv(A).
211 *
212  IF( j.LT.n )
213  $ CALL sgemv( 'No transpose', n, n-j, -one, a( 1, j+1 ),
214  $ lda, work( j+1 ), 1, one, a( 1, j ), 1 )
215  20 CONTINUE
216  ELSE
217 *
218 * Use blocked code.
219 *
220  nn = ( ( n-1 ) / nb )*nb + 1
221  DO 50 j = nn, 1, -nb
222  jb = min( nb, n-j+1 )
223 *
224 * Copy current block column of L to WORK and replace with
225 * zeros.
226 *
227  DO 40 jj = j, j + jb - 1
228  DO 30 i = jj + 1, n
229  work( i+( jj-j )*ldwork ) = a( i, jj )
230  a( i, jj ) = zero
231  30 CONTINUE
232  40 CONTINUE
233 *
234 * Compute current block column of inv(A).
235 *
236  IF( j+jb.LE.n )
237  $ CALL sgemm( 'No transpose', 'No transpose', n, jb,
238  $ n-j-jb+1, -one, a( 1, j+jb ), lda,
239  $ work( j+jb ), ldwork, one, a( 1, j ), lda )
240  CALL strsm( 'Right', 'Lower', 'No transpose', 'Unit', n, jb,
241  $ one, work( j ), ldwork, a( 1, j ), lda )
242  50 CONTINUE
243  END IF
244 *
245 * Apply column interchanges.
246 *
247  DO 60 j = n - 1, 1, -1
248  jp = ipiv( j )
249  IF( jp.NE.j )
250  $ CALL sswap( n, a( 1, j ), 1, a( 1, jp ), 1 )
251  60 CONTINUE
252 *
253  work( 1 ) = iws
254  RETURN
255 *
256 * End of SGETRI
257 *
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine strtri(UPLO, DIAG, N, A, LDA, INFO)
STRTRI
Definition: strtri.f:109
subroutine sswap(N, SX, INCX, SY, INCY)
SSWAP
Definition: sswap.f:82
subroutine sgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
SGEMV
Definition: sgemv.f:156
subroutine strsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRSM
Definition: strsm.f:181
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
Definition: sgemm.f:187
Here is the call graph for this function:
Here is the caller graph for this function: