LAPACK  3.9.0
LAPACK: Linear Algebra PACKage

◆ chetrf_aa_2stage()

subroutine chetrf_aa_2stage ( character  UPLO,
integer  N,
complex, dimension( lda, * )  A,
integer  LDA,
complex, dimension( * )  TB,
integer  LTB,
integer, dimension( * )  IPIV,
integer, dimension( * )  IPIV2,
complex, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

CHETRF_AA_2STAGE

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

Purpose:
 CHETRF_AA_2STAGE computes the factorization of a real hermitian matrix A
 using the Aasen's algorithm.  The form of the factorization is

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

 where U (or L) is a product of permutation and unit upper (lower)
 triangular matrices, and T is a hermitian band matrix with the
 bandwidth of NB (NB is internally selected and stored in TB( 1 ), and T is 
 LU factorized with partial pivoting).

 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, L is stored below (or above) the subdiaonal blocks,
          when UPLO  is 'L' (or 'U').
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]TB
          TB is COMPLEX array, dimension (LTB)
          On exit, details of the LU factorization of the band matrix.
[in]LTB
          LTB is INTEGER
          The size of the array TB. LTB >= 4*N, internally
          used to select NB such that LTB >= (3*NB+1)*N.

          If LTB = -1, then a workspace query is assumed; the
          routine only calculates the optimal size of LTB, 
          returns this value as the first entry of TB, and
          no error message related to LTB is issued by XERBLA.
[out]IPIV
          IPIV is INTEGER array, dimension (N)
          On exit, it contains the details of the interchanges, i.e.,
          the row and column k of A were interchanged with the
          row and column IPIV(k).
[out]IPIV2
          IPIV2 is INTEGER array, dimension (N)
          On exit, it contains the details of the interchanges, i.e.,
          the row and column k of T were interchanged with the
          row and column IPIV(k).
[out]WORK
          WORK is COMPLEX workspace of size LWORK
[in]LWORK
          LWORK is INTEGER
          The size of WORK. LWORK >= N, internally used to select NB
          such that LWORK >= N*NB.

          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, band LU factorization failed on i-th column
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2017

Definition at line 162 of file chetrf_aa_2stage.f.

162 *
163 * -- LAPACK computational routine (version 3.8.0) --
164 * -- LAPACK is a software package provided by Univ. of Tennessee, --
165 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
166 * November 2017
167 *
168  IMPLICIT NONE
169 *
170 * .. Scalar Arguments ..
171  CHARACTER UPLO
172  INTEGER N, LDA, LTB, LWORK, INFO
173 * ..
174 * .. Array Arguments ..
175  INTEGER IPIV( * ), IPIV2( * )
176  COMPLEX A( LDA, * ), TB( * ), WORK( * )
177 * ..
178 *
179 * =====================================================================
180 * .. Parameters ..
181  COMPLEX ZERO, ONE
182  parameter( zero = ( 0.0e+0, 0.0e+0 ),
183  $ one = ( 1.0e+0, 0.0e+0 ) )
184 *
185 * .. Local Scalars ..
186  LOGICAL UPPER, TQUERY, WQUERY
187  INTEGER I, J, K, I1, I2, TD
188  INTEGER LDTB, NB, KB, JB, NT, IINFO
189  COMPLEX PIV
190 * ..
191 * .. External Functions ..
192  LOGICAL LSAME
193  INTEGER ILAENV
194  EXTERNAL lsame, ilaenv
195 
196 * ..
197 * .. External Subroutines ..
198  EXTERNAL xerbla, ccopy, clacgv, clacpy,
199  $ claset, cgbtrf, cgemm, cgetrf,
200  $ chegst, cswap, ctrsm
201 * ..
202 * .. Intrinsic Functions ..
203  INTRINSIC conjg, min, max
204 * ..
205 * .. Executable Statements ..
206 *
207 * Test the input parameters.
208 *
209  info = 0
210  upper = lsame( uplo, 'U' )
211  wquery = ( lwork.EQ.-1 )
212  tquery = ( ltb.EQ.-1 )
213  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
214  info = -1
215  ELSE IF( n.LT.0 ) THEN
216  info = -2
217  ELSE IF( lda.LT.max( 1, n ) ) THEN
218  info = -4
219  ELSE IF ( ltb .LT. 4*n .AND. .NOT.tquery ) THEN
220  info = -6
221  ELSE IF ( lwork .LT. n .AND. .NOT.wquery ) THEN
222  info = -10
223  END IF
224 *
225  IF( info.NE.0 ) THEN
226  CALL xerbla( 'CHETRF_AA_2STAGE', -info )
227  RETURN
228  END IF
229 *
230 * Answer the query
231 *
232  nb = ilaenv( 1, 'CHETRF_AA_2STAGE', uplo, n, -1, -1, -1 )
233  IF( info.EQ.0 ) THEN
234  IF( tquery ) THEN
235  tb( 1 ) = (3*nb+1)*n
236  END IF
237  IF( wquery ) THEN
238  work( 1 ) = n*nb
239  END IF
240  END IF
241  IF( tquery .OR. wquery ) THEN
242  RETURN
243  END IF
244 *
245 * Quick return
246 *
247  IF ( n.EQ.0 ) THEN
248  RETURN
249  ENDIF
250 *
251 * Determine the number of the block size
252 *
253  ldtb = ltb/n
254  IF( ldtb .LT. 3*nb+1 ) THEN
255  nb = (ldtb-1)/3
256  END IF
257  IF( lwork .LT. nb*n ) THEN
258  nb = lwork/n
259  END IF
260 *
261 * Determine the number of the block columns
262 *
263  nt = (n+nb-1)/nb
264  td = 2*nb
265  kb = min(nb, n)
266 *
267 * Initialize vectors/matrices
268 *
269  DO j = 1, kb
270  ipiv( j ) = j
271  END DO
272 *
273 * Save NB
274 *
275  tb( 1 ) = nb
276 *
277  IF( upper ) THEN
278 *
279 * .....................................................
280 * Factorize A as U**T*D*U using the upper triangle of A
281 * .....................................................
282 *
283  DO j = 0, nt-1
284 *
285 * Generate Jth column of W and H
286 *
287  kb = min(nb, n-j*nb)
288  DO i = 1, j-1
289  IF( i.EQ.1 ) THEN
290 * H(I,J) = T(I,I)*U(I,J) + T(I+1,I)*U(I+1,J)
291  IF( i .EQ. (j-1) ) THEN
292  jb = nb+kb
293  ELSE
294  jb = 2*nb
295  END IF
296  CALL cgemm( 'NoTranspose', 'NoTranspose',
297  $ nb, kb, jb,
298  $ one, tb( td+1 + (i*nb)*ldtb ), ldtb-1,
299  $ a( (i-1)*nb+1, j*nb+1 ), lda,
300  $ zero, work( i*nb+1 ), n )
301  ELSE
302 * H(I,J) = T(I,I-1)*U(I-1,J) + T(I,I)*U(I,J) + T(I,I+1)*U(I+1,J)
303  IF( i .EQ. (j-1) ) THEN
304  jb = 2*nb+kb
305  ELSE
306  jb = 3*nb
307  END IF
308  CALL cgemm( 'NoTranspose', 'NoTranspose',
309  $ nb, kb, jb,
310  $ one, tb( td+nb+1 + ((i-1)*nb)*ldtb ),
311  $ ldtb-1,
312  $ a( (i-2)*nb+1, j*nb+1 ), lda,
313  $ zero, work( i*nb+1 ), n )
314  END IF
315  END DO
316 *
317 * Compute T(J,J)
318 *
319  CALL clacpy( 'Upper', kb, kb, a( j*nb+1, j*nb+1 ), lda,
320  $ tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
321  IF( j.GT.1 ) THEN
322 * T(J,J) = U(1:J,J)'*H(1:J)
323  CALL cgemm( 'Conjugate transpose', 'NoTranspose',
324  $ kb, kb, (j-1)*nb,
325  $ -one, a( 1, j*nb+1 ), lda,
326  $ work( nb+1 ), n,
327  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
328 * T(J,J) += U(J,J)'*T(J,J-1)*U(J-1,J)
329  CALL cgemm( 'Conjugate transpose', 'NoTranspose',
330  $ kb, nb, kb,
331  $ one, a( (j-1)*nb+1, j*nb+1 ), lda,
332  $ tb( td+nb+1 + ((j-1)*nb)*ldtb ), ldtb-1,
333  $ zero, work( 1 ), n )
334  CALL cgemm( 'NoTranspose', 'NoTranspose',
335  $ kb, kb, nb,
336  $ -one, work( 1 ), n,
337  $ a( (j-2)*nb+1, j*nb+1 ), lda,
338  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
339  END IF
340  IF( j.GT.0 ) THEN
341  CALL chegst( 1, 'Upper', kb,
342  $ tb( td+1 + (j*nb)*ldtb ), ldtb-1,
343  $ a( (j-1)*nb+1, j*nb+1 ), lda, iinfo )
344  END IF
345 *
346 * Expand T(J,J) into full format
347 *
348  DO i = 1, kb
349  tb( td+1 + (j*nb+i-1)*ldtb )
350  $ = real( tb( td+1 + (j*nb+i-1)*ldtb ) )
351  DO k = i+1, kb
352  tb( td+(k-i)+1 + (j*nb+i-1)*ldtb )
353  $ = conjg( tb( td-(k-(i+1)) + (j*nb+k-1)*ldtb ) )
354  END DO
355  END DO
356 *
357  IF( j.LT.nt-1 ) THEN
358  IF( j.GT.0 ) THEN
359 *
360 * Compute H(J,J)
361 *
362  IF( j.EQ.1 ) THEN
363  CALL cgemm( 'NoTranspose', 'NoTranspose',
364  $ kb, kb, kb,
365  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1,
366  $ a( (j-1)*nb+1, j*nb+1 ), lda,
367  $ zero, work( j*nb+1 ), n )
368  ELSE
369  CALL cgemm( 'NoTranspose', 'NoTranspose',
370  $ kb, kb, nb+kb,
371  $ one, tb( td+nb+1 + ((j-1)*nb)*ldtb ),
372  $ ldtb-1,
373  $ a( (j-2)*nb+1, j*nb+1 ), lda,
374  $ zero, work( j*nb+1 ), n )
375  END IF
376 *
377 * Update with the previous column
378 *
379  CALL cgemm( 'Conjugate transpose', 'NoTranspose',
380  $ nb, n-(j+1)*nb, j*nb,
381  $ -one, work( nb+1 ), n,
382  $ a( 1, (j+1)*nb+1 ), lda,
383  $ one, a( j*nb+1, (j+1)*nb+1 ), lda )
384  END IF
385 *
386 * Copy panel to workspace to call CGETRF
387 *
388  DO k = 1, nb
389  CALL ccopy( n-(j+1)*nb,
390  $ a( j*nb+k, (j+1)*nb+1 ), lda,
391  $ work( 1+(k-1)*n ), 1 )
392  END DO
393 *
394 * Factorize panel
395 *
396  CALL cgetrf( n-(j+1)*nb, nb,
397  $ work, n,
398  $ ipiv( (j+1)*nb+1 ), iinfo )
399 c IF (IINFO.NE.0 .AND. INFO.EQ.0) THEN
400 c INFO = IINFO+(J+1)*NB
401 c END IF
402 *
403 * Copy panel back
404 *
405  DO k = 1, nb
406 *
407 * Copy only L-factor
408 *
409  CALL ccopy( n-k-(j+1)*nb,
410  $ work( k+1+(k-1)*n ), 1,
411  $ a( j*nb+k, (j+1)*nb+k+1 ), lda )
412 *
413 * Transpose U-factor to be copied back into T(J+1, J)
414 *
415  CALL clacgv( k, work( 1+(k-1)*n ), 1 )
416  END DO
417 *
418 * Compute T(J+1, J), zero out for GEMM update
419 *
420  kb = min(nb, n-(j+1)*nb)
421  CALL claset( 'Full', kb, nb, zero, zero,
422  $ tb( td+nb+1 + (j*nb)*ldtb), ldtb-1 )
423  CALL clacpy( 'Upper', kb, nb,
424  $ work, n,
425  $ tb( td+nb+1 + (j*nb)*ldtb ), ldtb-1 )
426  IF( j.GT.0 ) THEN
427  CALL ctrsm( 'R', 'U', 'N', 'U', kb, nb, one,
428  $ a( (j-1)*nb+1, j*nb+1 ), lda,
429  $ tb( td+nb+1 + (j*nb)*ldtb ), ldtb-1 )
430  END IF
431 *
432 * Copy T(J,J+1) into T(J+1, J), both upper/lower for GEMM
433 * updates
434 *
435  DO k = 1, nb
436  DO i = 1, kb
437  tb( td-nb+k-i+1 + (j*nb+nb+i-1)*ldtb )
438  $ = conjg( tb( td+nb+i-k+1 + (j*nb+k-1)*ldtb ) )
439  END DO
440  END DO
441  CALL claset( 'Lower', kb, nb, zero, one,
442  $ a( j*nb+1, (j+1)*nb+1), lda )
443 *
444 * Apply pivots to trailing submatrix of A
445 *
446  DO k = 1, kb
447 * > Adjust ipiv
448  ipiv( (j+1)*nb+k ) = ipiv( (j+1)*nb+k ) + (j+1)*nb
449 *
450  i1 = (j+1)*nb+k
451  i2 = ipiv( (j+1)*nb+k )
452  IF( i1.NE.i2 ) THEN
453 * > Apply pivots to previous columns of L
454  CALL cswap( k-1, a( (j+1)*nb+1, i1 ), 1,
455  $ a( (j+1)*nb+1, i2 ), 1 )
456 * > Swap A(I1+1:M, I1) with A(I2, I1+1:M)
457  IF( i2.GT.(i1+1) ) THEN
458  CALL cswap( i2-i1-1, a( i1, i1+1 ), lda,
459  $ a( i1+1, i2 ), 1 )
460  CALL clacgv( i2-i1-1, a( i1+1, i2 ), 1 )
461  END IF
462  CALL clacgv( i2-i1, a( i1, i1+1 ), lda )
463 * > Swap A(I2+1:M, I1) with A(I2+1:M, I2)
464  IF( i2.LT.n )
465  $ CALL cswap( n-i2, a( i1, i2+1 ), lda,
466  $ a( i2, i2+1 ), lda )
467 * > Swap A(I1, I1) with A(I2, I2)
468  piv = a( i1, i1 )
469  a( i1, i1 ) = a( i2, i2 )
470  a( i2, i2 ) = piv
471 * > Apply pivots to previous columns of L
472  IF( j.GT.0 ) THEN
473  CALL cswap( j*nb, a( 1, i1 ), 1,
474  $ a( 1, i2 ), 1 )
475  END IF
476  ENDIF
477  END DO
478  END IF
479  END DO
480  ELSE
481 *
482 * .....................................................
483 * Factorize A as L*D*L**T using the lower triangle of A
484 * .....................................................
485 *
486  DO j = 0, nt-1
487 *
488 * Generate Jth column of W and H
489 *
490  kb = min(nb, n-j*nb)
491  DO i = 1, j-1
492  IF( i.EQ.1 ) THEN
493 * H(I,J) = T(I,I)*L(J,I)' + T(I+1,I)'*L(J,I+1)'
494  IF( i .EQ. (j-1) ) THEN
495  jb = nb+kb
496  ELSE
497  jb = 2*nb
498  END IF
499  CALL cgemm( 'NoTranspose', 'Conjugate transpose',
500  $ nb, kb, jb,
501  $ one, tb( td+1 + (i*nb)*ldtb ), ldtb-1,
502  $ a( j*nb+1, (i-1)*nb+1 ), lda,
503  $ zero, work( i*nb+1 ), n )
504  ELSE
505 * H(I,J) = T(I,I-1)*L(J,I-1)' + T(I,I)*L(J,I)' + T(I,I+1)*L(J,I+1)'
506  IF( i .EQ. (j-1) ) THEN
507  jb = 2*nb+kb
508  ELSE
509  jb = 3*nb
510  END IF
511  CALL cgemm( 'NoTranspose', 'Conjugate transpose',
512  $ nb, kb, jb,
513  $ one, tb( td+nb+1 + ((i-1)*nb)*ldtb ),
514  $ ldtb-1,
515  $ a( j*nb+1, (i-2)*nb+1 ), lda,
516  $ zero, work( i*nb+1 ), n )
517  END IF
518  END DO
519 *
520 * Compute T(J,J)
521 *
522  CALL clacpy( 'Lower', kb, kb, a( j*nb+1, j*nb+1 ), lda,
523  $ tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
524  IF( j.GT.1 ) THEN
525 * T(J,J) = L(J,1:J)*H(1:J)
526  CALL cgemm( 'NoTranspose', 'NoTranspose',
527  $ kb, kb, (j-1)*nb,
528  $ -one, a( j*nb+1, 1 ), lda,
529  $ work( nb+1 ), n,
530  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
531 * T(J,J) += L(J,J)*T(J,J-1)*L(J,J-1)'
532  CALL cgemm( 'NoTranspose', 'NoTranspose',
533  $ kb, nb, kb,
534  $ one, a( j*nb+1, (j-1)*nb+1 ), lda,
535  $ tb( td+nb+1 + ((j-1)*nb)*ldtb ), ldtb-1,
536  $ zero, work( 1 ), n )
537  CALL cgemm( 'NoTranspose', 'Conjugate transpose',
538  $ kb, kb, nb,
539  $ -one, work( 1 ), n,
540  $ a( j*nb+1, (j-2)*nb+1 ), lda,
541  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1 )
542  END IF
543  IF( j.GT.0 ) THEN
544  CALL chegst( 1, 'Lower', kb,
545  $ tb( td+1 + (j*nb)*ldtb ), ldtb-1,
546  $ a( j*nb+1, (j-1)*nb+1 ), lda, iinfo )
547  END IF
548 *
549 * Expand T(J,J) into full format
550 *
551  DO i = 1, kb
552  tb( td+1 + (j*nb+i-1)*ldtb )
553  $ = real( tb( td+1 + (j*nb+i-1)*ldtb ) )
554  DO k = i+1, kb
555  tb( td-(k-(i+1)) + (j*nb+k-1)*ldtb )
556  $ = conjg( tb( td+(k-i)+1 + (j*nb+i-1)*ldtb ) )
557  END DO
558  END DO
559 *
560  IF( j.LT.nt-1 ) THEN
561  IF( j.GT.0 ) THEN
562 *
563 * Compute H(J,J)
564 *
565  IF( j.EQ.1 ) THEN
566  CALL cgemm( 'NoTranspose', 'Conjugate transpose',
567  $ kb, kb, kb,
568  $ one, tb( td+1 + (j*nb)*ldtb ), ldtb-1,
569  $ a( j*nb+1, (j-1)*nb+1 ), lda,
570  $ zero, work( j*nb+1 ), n )
571  ELSE
572  CALL cgemm( 'NoTranspose', 'Conjugate transpose',
573  $ kb, kb, nb+kb,
574  $ one, tb( td+nb+1 + ((j-1)*nb)*ldtb ),
575  $ ldtb-1,
576  $ a( j*nb+1, (j-2)*nb+1 ), lda,
577  $ zero, work( j*nb+1 ), n )
578  END IF
579 *
580 * Update with the previous column
581 *
582  CALL cgemm( 'NoTranspose', 'NoTranspose',
583  $ n-(j+1)*nb, nb, j*nb,
584  $ -one, a( (j+1)*nb+1, 1 ), lda,
585  $ work( nb+1 ), n,
586  $ one, a( (j+1)*nb+1, j*nb+1 ), lda )
587  END IF
588 *
589 * Factorize panel
590 *
591  CALL cgetrf( n-(j+1)*nb, nb,
592  $ a( (j+1)*nb+1, j*nb+1 ), lda,
593  $ ipiv( (j+1)*nb+1 ), iinfo )
594 c IF (IINFO.NE.0 .AND. INFO.EQ.0) THEN
595 c INFO = IINFO+(J+1)*NB
596 c END IF
597 *
598 * Compute T(J+1, J), zero out for GEMM update
599 *
600  kb = min(nb, n-(j+1)*nb)
601  CALL claset( 'Full', kb, nb, zero, zero,
602  $ tb( td+nb+1 + (j*nb)*ldtb), ldtb-1 )
603  CALL clacpy( 'Upper', kb, nb,
604  $ a( (j+1)*nb+1, j*nb+1 ), lda,
605  $ tb( td+nb+1 + (j*nb)*ldtb ), ldtb-1 )
606  IF( j.GT.0 ) THEN
607  CALL ctrsm( 'R', 'L', 'C', 'U', kb, nb, one,
608  $ a( j*nb+1, (j-1)*nb+1 ), lda,
609  $ tb( td+nb+1 + (j*nb)*ldtb ), ldtb-1 )
610  END IF
611 *
612 * Copy T(J+1,J) into T(J, J+1), both upper/lower for GEMM
613 * updates
614 *
615  DO k = 1, nb
616  DO i = 1, kb
617  tb( td-nb+k-i+1 + (j*nb+nb+i-1)*ldtb )
618  $ = conjg( tb( td+nb+i-k+1 + (j*nb+k-1)*ldtb ) )
619  END DO
620  END DO
621  CALL claset( 'Upper', kb, nb, zero, one,
622  $ a( (j+1)*nb+1, j*nb+1), lda )
623 *
624 * Apply pivots to trailing submatrix of A
625 *
626  DO k = 1, kb
627 * > Adjust ipiv
628  ipiv( (j+1)*nb+k ) = ipiv( (j+1)*nb+k ) + (j+1)*nb
629 *
630  i1 = (j+1)*nb+k
631  i2 = ipiv( (j+1)*nb+k )
632  IF( i1.NE.i2 ) THEN
633 * > Apply pivots to previous columns of L
634  CALL cswap( k-1, a( i1, (j+1)*nb+1 ), lda,
635  $ a( i2, (j+1)*nb+1 ), lda )
636 * > Swap A(I1+1:M, I1) with A(I2, I1+1:M)
637  IF( i2.GT.(i1+1) ) THEN
638  CALL cswap( i2-i1-1, a( i1+1, i1 ), 1,
639  $ a( i2, i1+1 ), lda )
640  CALL clacgv( i2-i1-1, a( i2, i1+1 ), lda )
641  END IF
642  CALL clacgv( i2-i1, a( i1+1, i1 ), 1 )
643 * > Swap A(I2+1:M, I1) with A(I2+1:M, I2)
644  IF( i2.LT.n )
645  $ CALL cswap( n-i2, a( i2+1, i1 ), 1,
646  $ a( i2+1, i2 ), 1 )
647 * > Swap A(I1, I1) with A(I2, I2)
648  piv = a( i1, i1 )
649  a( i1, i1 ) = a( i2, i2 )
650  a( i2, i2 ) = piv
651 * > Apply pivots to previous columns of L
652  IF( j.GT.0 ) THEN
653  CALL cswap( j*nb, a( i1, 1 ), lda,
654  $ a( i2, 1 ), lda )
655  END IF
656  ENDIF
657  END DO
658 *
659 * Apply pivots to previous columns of L
660 *
661 c CALL CLASWP( J*NB, A( 1, 1 ), LDA,
662 c $ (J+1)*NB+1, (J+1)*NB+KB, IPIV, 1 )
663  END IF
664  END DO
665  END IF
666 *
667 * Factor the band matrix
668  CALL cgbtrf( n, n, nb, nb, tb, ldtb, ipiv2, info )
669 *
670  RETURN
671 *
672 * End of CHETRF_AA_2STAGE
673 *
Here is the call graph for this function:
Here is the caller graph for this function:
cgetrf
subroutine cgetrf(M, N, A, LDA, IPIV, INFO)
CGETRF
Definition: cgetrf.f:110
cgemm
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
Definition: cgemm.f:189
ctrsm
subroutine ctrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
CTRSM
Definition: ctrsm.f:182
clacgv
subroutine clacgv(N, X, INCX)
CLACGV conjugates a complex vector.
Definition: clacgv.f:76
clacpy
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
chegst
subroutine chegst(ITYPE, UPLO, N, A, LDA, B, LDB, INFO)
CHEGST
Definition: chegst.f:130
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
cswap
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
Definition: cswap.f:83
lsame
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
claset
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: claset.f:108
ccopy
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
Definition: ccopy.f:83
ilaenv
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: tstiee.f:83
cgbtrf
subroutine cgbtrf(M, N, KL, KU, AB, LDAB, IPIV, INFO)
CGBTRF
Definition: cgbtrf.f:146