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

◆ zhetrf_rook()

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

ZHETRF_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 ZHETRF_ROOK + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!> ZHETRF_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*16 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*16 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 zhetrf_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*16 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 EXTERNAL lsame, ilaenv
235* ..
236* .. External Subroutines ..
238* ..
239* .. Intrinsic Functions ..
240 INTRINSIC max
241* ..
242* .. Executable Statements ..
243*
244* Test the input parameters.
245*
246 info = 0
247 upper = lsame( uplo, 'U' )
248 lquery = ( lwork.EQ.-1 )
249 IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
250 info = -1
251 ELSE IF( n.LT.0 ) THEN
252 info = -2
253 ELSE IF( lda.LT.max( 1, n ) ) THEN
254 info = -4
255 ELSE IF( lwork.LT.1 .AND. .NOT.lquery ) THEN
256 info = -7
257 END IF
258*
259 IF( info.EQ.0 ) THEN
260*
261* Determine the block size
262*
263 nb = ilaenv( 1, 'ZHETRF_ROOK', uplo, n, -1, -1, -1 )
264 lwkopt = max( 1, n*nb )
265 work( 1 ) = lwkopt
266 END IF
267*
268 IF( info.NE.0 ) THEN
269 CALL xerbla( 'ZHETRF_ROOK', -info )
270 RETURN
271 ELSE IF( lquery ) THEN
272 RETURN
273 END IF
274*
275 nbmin = 2
276 ldwork = n
277 IF( nb.GT.1 .AND. nb.LT.n ) THEN
278 iws = ldwork*nb
279 IF( lwork.LT.iws ) THEN
280 nb = max( lwork / ldwork, 1 )
281 nbmin = max( 2, ilaenv( 2, 'ZHETRF_ROOK',
282 $ uplo, n, -1, -1, -1 ) )
283 END IF
284 ELSE
285 iws = 1
286 END IF
287 IF( nb.LT.nbmin )
288 $ nb = n
289*
290 IF( upper ) THEN
291*
292* Factorize A as U*D*U**T using the upper triangle of A
293*
294* K is the main loop index, decreasing from N to 1 in steps of
295* KB, where KB is the number of columns factorized by ZLAHEF_ROOK;
296* KB is either NB or NB-1, or K for the last block
297*
298 k = n
299 10 CONTINUE
300*
301* If K < 1, exit from loop
302*
303 IF( k.LT.1 )
304 $ GO TO 40
305*
306 IF( k.GT.nb ) THEN
307*
308* Factorize columns k-kb+1:k of A and use blocked code to
309* update columns 1:k-kb
310*
311 CALL zlahef_rook( uplo, k, nb, kb, a, lda,
312 $ ipiv, work, ldwork, iinfo )
313 ELSE
314*
315* Use unblocked code to factorize columns 1:k of A
316*
317 CALL zhetf2_rook( uplo, k, a, lda, ipiv, iinfo )
318 kb = k
319 END IF
320*
321* Set INFO on the first occurrence of a zero pivot
322*
323 IF( info.EQ.0 .AND. iinfo.GT.0 )
324 $ info = iinfo
325*
326* No need to adjust IPIV
327*
328* Decrease K and return to the start of the main loop
329*
330 k = k - kb
331 GO TO 10
332*
333 ELSE
334*
335* Factorize A as L*D*L**T using the lower triangle of A
336*
337* K is the main loop index, increasing from 1 to N in steps of
338* KB, where KB is the number of columns factorized by ZLAHEF_ROOK;
339* KB is either NB or NB-1, or N-K+1 for the last block
340*
341 k = 1
342 20 CONTINUE
343*
344* If K > N, exit from loop
345*
346 IF( k.GT.n )
347 $ GO TO 40
348*
349 IF( k.LE.n-nb ) THEN
350*
351* Factorize columns k:k+kb-1 of A and use blocked code to
352* update columns k+kb:n
353*
354 CALL zlahef_rook( uplo, n-k+1, nb, kb, a( k, k ), lda,
355 $ ipiv( k ), work, ldwork, iinfo )
356 ELSE
357*
358* Use unblocked code to factorize columns k:n of A
359*
360 CALL zhetf2_rook( uplo, n-k+1, a( k, k ), lda, ipiv( k ),
361 $ iinfo )
362 kb = n - k + 1
363 END IF
364*
365* Set INFO on the first occurrence of a zero pivot
366*
367 IF( info.EQ.0 .AND. iinfo.GT.0 )
368 $ info = iinfo + k - 1
369*
370* Adjust IPIV
371*
372 DO 30 j = k, k + kb - 1
373 IF( ipiv( j ).GT.0 ) THEN
374 ipiv( j ) = ipiv( j ) + k - 1
375 ELSE
376 ipiv( j ) = ipiv( j ) - k + 1
377 END IF
378 30 CONTINUE
379*
380* Increase K and return to the start of the main loop
381*
382 k = k + kb
383 GO TO 20
384*
385 END IF
386*
387 40 CONTINUE
388 work( 1 ) = lwkopt
389 RETURN
390*
391* End of ZHETRF_ROOK
392*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine zhetf2_rook(uplo, n, a, lda, ipiv, info)
ZHETF2_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 zlahef_rook(uplo, n, nb, kb, a, lda, ipiv, w, ldw, info)
Download ZLAHEF_ROOK + dependencies [TGZ] [ZIP] [TXT]
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: