LAPACK  3.8.0
LAPACK: Linear Algebra PACKage

◆ dsytrf()

subroutine dsytrf ( character  UPLO,
integer  N,
double precision, dimension( lda, * )  A,
integer  LDA,
integer, dimension( * )  IPIV,
double precision, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

DSYTRF

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

Purpose:
 DSYTRF computes the factorization of a real symmetric matrix A using
 the Bunch-Kaufman diagonal pivoting method.  The form of the
 factorization is

    A = U*D*U**T  or  A = L*D*L**T

 where U (or L) is a product of permutation and unit upper (lower)
 triangular matrices, and D is symmetric and block diagonal with
 1-by-1 and 2-by-2 diagonal blocks.

 This is the blocked version of the algorithm, calling Level 3 BLAS.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is DOUBLE PRECISION array, dimension (LDA,N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the leading
          N-by-N upper triangular part of A contains the upper
          triangular part of the matrix A, and the strictly lower
          triangular part of A is not referenced.  If UPLO = 'L', the
          leading N-by-N lower triangular part of A contains the lower
          triangular part of the matrix A, and the strictly upper
          triangular part of A is not referenced.

          On exit, the block diagonal matrix D and the multipliers used
          to obtain the factor U or L (see below for further details).
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]IPIV
          IPIV is INTEGER array, dimension (N)
          Details of the interchanges and the block structure of D.
          If IPIV(k) > 0, then rows and columns k and IPIV(k) were
          interchanged and D(k,k) is a 1-by-1 diagonal block.
          If UPLO = 'U' and IPIV(k) = IPIV(k-1) < 0, then rows and
          columns k-1 and -IPIV(k) were interchanged and D(k-1:k,k-1:k)
          is a 2-by-2 diagonal block.  If UPLO = 'L' and IPIV(k) =
          IPIV(k+1) < 0, then rows and columns k+1 and -IPIV(k) were
          interchanged and D(k:k+1,k:k+1) is a 2-by-2 diagonal block.
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The length of WORK.  LWORK >=1.  For best performance
          LWORK >= N*NB, where NB is the block size 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, D(i,i) is exactly zero.  The factorization
                has been completed, but the block diagonal matrix D is
                exactly singular, and division by zero will occur if it
                is used to solve a system of equations.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
December 2016
Further Details:
  If UPLO = 'U', then A = U*D*U**T, where
     U = P(n)*U(n)* ... *P(k)U(k)* ...,
  i.e., U is a product of terms P(k)*U(k), where k decreases from n to
  1 in steps of 1 or 2, and D is a block diagonal matrix with 1-by-1
  and 2-by-2 diagonal blocks D(k).  P(k) is a permutation matrix as
  defined by IPIV(k), and U(k) is a unit upper triangular matrix, such
  that if the diagonal block D(k) is of order s (s = 1 or 2), then

             (   I    v    0   )   k-s
     U(k) =  (   0    I    0   )   s
             (   0    0    I   )   n-k
                k-s   s   n-k

  If s = 1, D(k) overwrites A(k,k), and v overwrites A(1:k-1,k).
  If s = 2, the upper triangle of D(k) overwrites A(k-1,k-1), A(k-1,k),
  and A(k,k), and v overwrites A(1:k-2,k-1:k).

  If UPLO = 'L', then A = L*D*L**T, where
     L = P(1)*L(1)* ... *P(k)*L(k)* ...,
  i.e., L is a product of terms P(k)*L(k), where k increases from 1 to
  n in steps of 1 or 2, and D is a block diagonal matrix with 1-by-1
  and 2-by-2 diagonal blocks D(k).  P(k) is a permutation matrix as
  defined by IPIV(k), and L(k) is a unit lower triangular matrix, such
  that if the diagonal block D(k) is of order s (s = 1 or 2), then

             (   I    0     0   )  k-1
     L(k) =  (   0    I     0   )  s
             (   0    v     I   )  n-k-s+1
                k-1   s  n-k-s+1

  If s = 1, D(k) overwrites A(k,k), and v overwrites A(k+1:n,k).
  If s = 2, the lower triangle of D(k) overwrites A(k,k), A(k+1,k),
  and A(k+1,k+1), and v overwrites A(k+2:n,k:k+1).

Definition at line 184 of file dsytrf.f.

184 *
185 * -- LAPACK computational routine (version 3.7.0) --
186 * -- LAPACK is a software package provided by Univ. of Tennessee, --
187 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
188 * December 2016
189 *
190 * .. Scalar Arguments ..
191  CHARACTER uplo
192  INTEGER info, lda, lwork, n
193 * ..
194 * .. Array Arguments ..
195  INTEGER ipiv( * )
196  DOUBLE PRECISION a( lda, * ), work( * )
197 * ..
198 *
199 * =====================================================================
200 *
201 * .. Local Scalars ..
202  LOGICAL lquery, upper
203  INTEGER iinfo, iws, j, k, kb, ldwork, lwkopt, nb, nbmin
204 * ..
205 * .. External Functions ..
206  LOGICAL lsame
207  INTEGER ilaenv
208  EXTERNAL lsame, ilaenv
209 * ..
210 * .. External Subroutines ..
211  EXTERNAL dlasyf, dsytf2, xerbla
212 * ..
213 * .. Intrinsic Functions ..
214  INTRINSIC max
215 * ..
216 * .. Executable Statements ..
217 *
218 * Test the input parameters.
219 *
220  info = 0
221  upper = lsame( uplo, 'U' )
222  lquery = ( lwork.EQ.-1 )
223  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
224  info = -1
225  ELSE IF( n.LT.0 ) THEN
226  info = -2
227  ELSE IF( lda.LT.max( 1, n ) ) THEN
228  info = -4
229  ELSE IF( lwork.LT.1 .AND. .NOT.lquery ) THEN
230  info = -7
231  END IF
232 *
233  IF( info.EQ.0 ) THEN
234 *
235 * Determine the block size
236 *
237  nb = ilaenv( 1, 'DSYTRF', uplo, n, -1, -1, -1 )
238  lwkopt = n*nb
239  work( 1 ) = lwkopt
240  END IF
241 *
242  IF( info.NE.0 ) THEN
243  CALL xerbla( 'DSYTRF', -info )
244  RETURN
245  ELSE IF( lquery ) THEN
246  RETURN
247  END IF
248 *
249  nbmin = 2
250  ldwork = n
251  IF( nb.GT.1 .AND. nb.LT.n ) THEN
252  iws = ldwork*nb
253  IF( lwork.LT.iws ) THEN
254  nb = max( lwork / ldwork, 1 )
255  nbmin = max( 2, ilaenv( 2, 'DSYTRF', uplo, n, -1, -1, -1 ) )
256  END IF
257  ELSE
258  iws = 1
259  END IF
260  IF( nb.LT.nbmin )
261  $ nb = n
262 *
263  IF( upper ) THEN
264 *
265 * Factorize A as U*D*U**T using the upper triangle of A
266 *
267 * K is the main loop index, decreasing from N to 1 in steps of
268 * KB, where KB is the number of columns factorized by DLASYF;
269 * KB is either NB or NB-1, or K for the last block
270 *
271  k = n
272  10 CONTINUE
273 *
274 * If K < 1, exit from loop
275 *
276  IF( k.LT.1 )
277  $ GO TO 40
278 *
279  IF( k.GT.nb ) THEN
280 *
281 * Factorize columns k-kb+1:k of A and use blocked code to
282 * update columns 1:k-kb
283 *
284  CALL dlasyf( uplo, k, nb, kb, a, lda, ipiv, work, ldwork,
285  $ iinfo )
286  ELSE
287 *
288 * Use unblocked code to factorize columns 1:k of A
289 *
290  CALL dsytf2( uplo, k, a, lda, ipiv, iinfo )
291  kb = k
292  END IF
293 *
294 * Set INFO on the first occurrence of a zero pivot
295 *
296  IF( info.EQ.0 .AND. iinfo.GT.0 )
297  $ info = iinfo
298 *
299 * Decrease K and return to the start of the main loop
300 *
301  k = k - kb
302  GO TO 10
303 *
304  ELSE
305 *
306 * Factorize A as L*D*L**T using the lower triangle of A
307 *
308 * K is the main loop index, increasing from 1 to N in steps of
309 * KB, where KB is the number of columns factorized by DLASYF;
310 * KB is either NB or NB-1, or N-K+1 for the last block
311 *
312  k = 1
313  20 CONTINUE
314 *
315 * If K > N, exit from loop
316 *
317  IF( k.GT.n )
318  $ GO TO 40
319 *
320  IF( k.LE.n-nb ) THEN
321 *
322 * Factorize columns k:k+kb-1 of A and use blocked code to
323 * update columns k+kb:n
324 *
325  CALL dlasyf( uplo, n-k+1, nb, kb, a( k, k ), lda, ipiv( k ),
326  $ work, ldwork, iinfo )
327  ELSE
328 *
329 * Use unblocked code to factorize columns k:n of A
330 *
331  CALL dsytf2( uplo, n-k+1, a( k, k ), lda, ipiv( k ), iinfo )
332  kb = n - k + 1
333  END IF
334 *
335 * Set INFO on the first occurrence of a zero pivot
336 *
337  IF( info.EQ.0 .AND. iinfo.GT.0 )
338  $ info = iinfo + k - 1
339 *
340 * Adjust IPIV
341 *
342  DO 30 j = k, k + kb - 1
343  IF( ipiv( j ).GT.0 ) THEN
344  ipiv( j ) = ipiv( j ) + k - 1
345  ELSE
346  ipiv( j ) = ipiv( j ) - k + 1
347  END IF
348  30 CONTINUE
349 *
350 * Increase K and return to the start of the main loop
351 *
352  k = k + kb
353  GO TO 20
354 *
355  END IF
356 *
357  40 CONTINUE
358  work( 1 ) = lwkopt
359  RETURN
360 *
361 * End of DSYTRF
362 *
subroutine dlasyf(UPLO, N, NB, KB, A, LDA, IPIV, W, LDW, INFO)
DLASYF computes a partial factorization of a real symmetric matrix using the Bunch-Kaufman diagonal p...
Definition: dlasyf.f:178
subroutine dsytf2(UPLO, N, A, LDA, IPIV, INFO)
DSYTF2 computes the factorization of a real symmetric indefinite matrix, using the diagonal pivoting ...
Definition: dsytf2.f:196
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: tstiee.f:83
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
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: