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

◆ dsytrf_rook()

subroutine dsytrf_rook ( character uplo,
integer n,
double precision, dimension( lda, * ) a,
integer lda,
integer, dimension( * ) ipiv,
double precision, dimension( * ) work,
integer lwork,
integer info )

DSYTRF_ROOK

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

Purpose:
!>
!> DSYTRF_ROOK computes the factorization of a real symmetric matrix A
!> using the bounded 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 UPLO = 'U':
!>               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 IPIV(k) < 0 and IPIV(k-1) < 0, then rows and
!>               columns k and -IPIV(k) were interchanged and rows and
!>               columns k-1 and -IPIV(k-1) were inerchaged,
!>               D(k-1:k,k-1:k) is a 2-by-2 diagonal block.
!>
!>          If UPLO = 'L':
!>               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 IPIV(k) < 0 and IPIV(k+1) < 0, then rows and
!>               columns k and -IPIV(k) were interchanged and rows and
!>               columns k+1 and -IPIV(k+1) were inerchaged,
!>               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.
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).
!>
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 205 of file dsytrf_rook.f.

207*
208* -- LAPACK computational routine --
209* -- LAPACK is a software package provided by Univ. of Tennessee, --
210* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
211*
212* .. Scalar Arguments ..
213 CHARACTER UPLO
214 INTEGER INFO, LDA, LWORK, N
215* ..
216* .. Array Arguments ..
217 INTEGER IPIV( * )
218 DOUBLE PRECISION A( LDA, * ), WORK( * )
219* ..
220*
221* =====================================================================
222*
223* .. Local Scalars ..
224 LOGICAL LQUERY, UPPER
225 INTEGER IINFO, IWS, J, K, KB, LDWORK, LWKOPT, NB, NBMIN
226* ..
227* .. External Functions ..
228 LOGICAL LSAME
229 INTEGER ILAENV
230 EXTERNAL lsame, ilaenv
231* ..
232* .. External Subroutines ..
233 EXTERNAL dlasyf_rook, dsytf2_rook, xerbla
234* ..
235* .. Intrinsic Functions ..
236 INTRINSIC max
237* ..
238* .. Executable Statements ..
239*
240* Test the input parameters.
241*
242 info = 0
243 upper = lsame( uplo, 'U' )
244 lquery = ( lwork.EQ.-1 )
245 IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
246 info = -1
247 ELSE IF( n.LT.0 ) THEN
248 info = -2
249 ELSE IF( lda.LT.max( 1, n ) ) THEN
250 info = -4
251 ELSE IF( lwork.LT.1 .AND. .NOT.lquery ) THEN
252 info = -7
253 END IF
254*
255 IF( info.EQ.0 ) THEN
256*
257* Determine the block size
258*
259 nb = ilaenv( 1, 'DSYTRF_ROOK', uplo, n, -1, -1, -1 )
260 lwkopt = max( 1, n*nb )
261 work( 1 ) = lwkopt
262 END IF
263*
264 IF( info.NE.0 ) THEN
265 CALL xerbla( 'DSYTRF_ROOK', -info )
266 RETURN
267 ELSE IF( lquery ) THEN
268 RETURN
269 END IF
270*
271 nbmin = 2
272 ldwork = n
273 IF( nb.GT.1 .AND. nb.LT.n ) THEN
274 iws = ldwork*nb
275 IF( lwork.LT.iws ) THEN
276 nb = max( lwork / ldwork, 1 )
277 nbmin = max( 2, ilaenv( 2, 'DSYTRF_ROOK',
278 $ uplo, n, -1, -1, -1 ) )
279 END IF
280 ELSE
281 iws = 1
282 END IF
283 IF( nb.LT.nbmin )
284 $ nb = n
285*
286 IF( upper ) THEN
287*
288* Factorize A as U*D*U**T using the upper triangle of A
289*
290* K is the main loop index, decreasing from N to 1 in steps of
291* KB, where KB is the number of columns factorized by DLASYF_ROOK;
292* KB is either NB or NB-1, or K for the last block
293*
294 k = n
295 10 CONTINUE
296*
297* If K < 1, exit from loop
298*
299 IF( k.LT.1 )
300 $ GO TO 40
301*
302 IF( k.GT.nb ) THEN
303*
304* Factorize columns k-kb+1:k of A and use blocked code to
305* update columns 1:k-kb
306*
307 CALL dlasyf_rook( uplo, k, nb, kb, a, lda,
308 $ ipiv, work, ldwork, iinfo )
309 ELSE
310*
311* Use unblocked code to factorize columns 1:k of A
312*
313 CALL dsytf2_rook( uplo, k, a, lda, ipiv, iinfo )
314 kb = k
315 END IF
316*
317* Set INFO on the first occurrence of a zero pivot
318*
319 IF( info.EQ.0 .AND. iinfo.GT.0 )
320 $ info = iinfo
321*
322* No need to adjust IPIV
323*
324* Decrease K and return to the start of the main loop
325*
326 k = k - kb
327 GO TO 10
328*
329 ELSE
330*
331* Factorize A as L*D*L**T using the lower triangle of A
332*
333* K is the main loop index, increasing from 1 to N in steps of
334* KB, where KB is the number of columns factorized by DLASYF_ROOK;
335* KB is either NB or NB-1, or N-K+1 for the last block
336*
337 k = 1
338 20 CONTINUE
339*
340* If K > N, exit from loop
341*
342 IF( k.GT.n )
343 $ GO TO 40
344*
345 IF( k.LE.n-nb ) THEN
346*
347* Factorize columns k:k+kb-1 of A and use blocked code to
348* update columns k+kb:n
349*
350 CALL dlasyf_rook( uplo, n-k+1, nb, kb, a( k, k ), lda,
351 $ ipiv( k ), work, ldwork, iinfo )
352 ELSE
353*
354* Use unblocked code to factorize columns k:n of A
355*
356 CALL dsytf2_rook( uplo, n-k+1, a( k, k ), lda, ipiv( k ),
357 $ iinfo )
358 kb = n - k + 1
359 END IF
360*
361* Set INFO on the first occurrence of a zero pivot
362*
363 IF( info.EQ.0 .AND. iinfo.GT.0 )
364 $ info = iinfo + k - 1
365*
366* Adjust IPIV
367*
368 DO 30 j = k, k + kb - 1
369 IF( ipiv( j ).GT.0 ) THEN
370 ipiv( j ) = ipiv( j ) + k - 1
371 ELSE
372 ipiv( j ) = ipiv( j ) - k + 1
373 END IF
374 30 CONTINUE
375*
376* Increase K and return to the start of the main loop
377*
378 k = k + kb
379 GO TO 20
380*
381 END IF
382*
383 40 CONTINUE
384 work( 1 ) = lwkopt
385 RETURN
386*
387* End of DSYTRF_ROOK
388*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
dsytf2_rook
subroutine dsytf2_rook(uplo, n, a, lda, ipiv, info)
DSYTF2_ROOK computes the factorization of a real symmetric indefinite matrix using the bounded Bunch-...
Definition dsytf2_rook.f:192
ilaenv
integer function ilaenv(ispec, name, opts, n1, n2, n3, n4)
ILAENV
Definition ilaenv.f:160
dlasyf_rook
subroutine dlasyf_rook(uplo, n, nb, kb, a, lda, ipiv, w, ldw, info)
DLASYF_ROOK *> DLASYF_ROOK computes a partial factorization of a real symmetric matrix using the boun...
Definition dlasyf_rook.f:182
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
Here is the call graph for this function:
Here is the caller graph for this function: