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

◆ chetrf_rook()

subroutine chetrf_rook ( character uplo,
integer n,
complex, dimension( lda, * ) a,
integer lda,
integer, dimension( * ) ipiv,
complex, dimension( * ) work,
integer lwork,
integer info )

CHETRF_ROOK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch-Kaufman ("rook") diagonal pivoting method (blocked algorithm, calling Level 3 BLAS).

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

Purpose:
!>
!> CHETRF_ROOK computes the factorization of a complex Hermitian 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 Hermitian 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 COMPLEX array, dimension (LDA,N)
!>          On entry, the Hermitian 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':
!>             Only the last KB elements of IPIV are set.
!>
!>             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':
!>             Only the first KB elements of IPIV are set.
!>
!>             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 COMPLEX 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:
!>
!>  June 2016,  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 209 of file chetrf_rook.f.

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