145 SUBROUTINE cgbtrf( M, N, KL, KU, AB, LDAB, IPIV, INFO )
153 INTEGER INFO, KL, KU, LDAB, M, N
157 COMPLEX AB( ldab, * )
164 parameter( one = ( 1.0e+0, 0.0e+0 ),
165 $ zero = ( 0.0e+0, 0.0e+0 ) )
166 INTEGER NBMAX, LDWORK
167 parameter( nbmax = 64, ldwork = nbmax+1 )
170 INTEGER I, I2, I3, II, IP, J, J2, J3, JB, JJ, JM, JP,
171 $ ju, k2, km, kv, nb, nw
175 COMPLEX WORK13( ldwork, nbmax ),
176 $ work31( ldwork, nbmax )
179 INTEGER ICAMAX, ILAENV
180 EXTERNAL icamax, ilaenv
201 ELSE IF( n.LT.0 )
THEN
203 ELSE IF( kl.LT.0 )
THEN
205 ELSE IF( ku.LT.0 )
THEN
207 ELSE IF( ldab.LT.kl+kv+1 )
THEN
211 CALL
xerbla(
'CGBTRF', -info )
217 IF( m.EQ.0 .OR. n.EQ.0 )
222 nb = ilaenv( 1,
'CGBTRF',
' ', m, n, kl, ku )
227 nb = min( nb, nbmax )
229 IF( nb.LE.1 .OR. nb.GT.kl )
THEN
233 CALL
cgbtf2( m, n, kl, ku, ab, ldab, ipiv, info )
242 work13( i, j ) = zero
250 work31( i, j ) = zero
258 DO 60 j = ku + 2, min( kv, n )
259 DO 50 i = kv - j + 2, kl
269 DO 180 j = 1, min( m, n ), nb
270 jb = min( nb, min( m, n )-j+1 )
284 i2 = min( kl-jb, m-j-jb+1 )
285 i3 = min( jb, m-j-kl+1 )
291 DO 80 jj = j, j + jb - 1
295 IF( jj+kv.LE.n )
THEN
297 ab( i, jj+kv ) = zero
305 jp = icamax( km+1, ab( kv+1, jj ), 1 )
306 ipiv( jj ) = jp + jj - j
307 IF( ab( kv+jp, jj ).NE.zero )
THEN
308 ju = max( ju, min( jj+ku+jp-1, n ) )
313 IF( jp+jj-1.LT.j+kl )
THEN
315 CALL
cswap( jb, ab( kv+1+jj-j, j ), ldab-1,
316 $ ab( kv+jp+jj-j, j ), ldab-1 )
322 CALL
cswap( jj-j, ab( kv+1+jj-j, j ), ldab-1,
323 $ work31( jp+jj-j-kl, 1 ), ldwork )
324 CALL
cswap( j+jb-jj, ab( kv+1, jj ), ldab-1,
325 $ ab( kv+jp, jj ), ldab-1 )
331 CALL
cscal( km, one / ab( kv+1, jj ), ab( kv+2, jj ),
338 jm = min( ju, j+jb-1 )
340 $ CALL
cgeru( km, jm-jj, -one, ab( kv+2, jj ), 1,
341 $ ab( kv, jj+1 ), ldab-1,
342 $ ab( kv+1, jj+1 ), ldab-1 )
354 nw = min( jj-j+1, i3 )
356 $ CALL
ccopy( nw, ab( kv+kl+1-jj+j, jj ), 1,
357 $ work31( 1, jj-j+1 ), 1 )
363 j2 = min( ju-j+1, kv ) - jb
364 j3 = max( 0, ju-j-kv+1 )
369 CALL
claswp( j2, ab( kv+1-jb, j+jb ), ldab-1, 1, jb,
374 DO 90 i = j, j + jb - 1
375 ipiv( i ) = ipiv( i ) + j - 1
384 DO 100 ii = j + i - 1, j + jb - 1
387 temp = ab( kv+1+ii-jj, jj )
388 ab( kv+1+ii-jj, jj ) = ab( kv+1+ip-jj, jj )
389 ab( kv+1+ip-jj, jj ) = temp
400 CALL
ctrsm(
'Left',
'Lower',
'No transpose',
'Unit',
401 $ jb, j2, one, ab( kv+1, j ), ldab-1,
402 $ ab( kv+1-jb, j+jb ), ldab-1 )
408 CALL
cgemm(
'No transpose',
'No transpose', i2, j2,
409 $ jb, -one, ab( kv+1+jb, j ), ldab-1,
410 $ ab( kv+1-jb, j+jb ), ldab-1, one,
411 $ ab( kv+1, j+jb ), ldab-1 )
418 CALL
cgemm(
'No transpose',
'No transpose', i3, j2,
419 $ jb, -one, work31, ldwork,
420 $ ab( kv+1-jb, j+jb ), ldab-1, one,
421 $ ab( kv+kl+1-jb, j+jb ), ldab-1 )
432 work13( ii, jj ) = ab( ii-jj+1, jj+j+kv-1 )
438 CALL
ctrsm(
'Left',
'Lower',
'No transpose',
'Unit',
439 $ jb, j3, one, ab( kv+1, j ), ldab-1,
446 CALL
cgemm(
'No transpose',
'No transpose', i2, j3,
447 $ jb, -one, ab( kv+1+jb, j ), ldab-1,
448 $ work13, ldwork, one, ab( 1+jb, j+kv ),
456 CALL
cgemm(
'No transpose',
'No transpose', i3, j3,
457 $ jb, -one, work31, ldwork, work13,
458 $ ldwork, one, ab( 1+kl, j+kv ), ldab-1 )
465 ab( ii-jj+1, jj+j+kv-1 ) = work13( ii, jj )
473 DO 160 i = j, j + jb - 1
474 ipiv( i ) = ipiv( i ) + j - 1
482 DO 170 jj = j + jb - 1, j, -1
483 jp = ipiv( jj ) - jj + 1
488 IF( jp+jj-1.LT.j+kl )
THEN
492 CALL
cswap( jj-j, ab( kv+1+jj-j, j ), ldab-1,
493 $ ab( kv+jp+jj-j, j ), ldab-1 )
498 CALL
cswap( jj-j, ab( kv+1+jj-j, j ), ldab-1,
499 $ work31( jp+jj-j-kl, 1 ), ldwork )
505 nw = min( i3, jj-j+1 )
507 $ CALL
ccopy( nw, work31( 1, jj-j+1 ), 1,
508 $ ab( kv+kl+1-jj+j, jj ), 1 )
subroutine claswp(N, A, LDA, K1, K2, IPIV, INCX)
CLASWP performs a series of row interchanges on a general rectangular matrix.
subroutine ctrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
CTRSM
subroutine cscal(N, CA, CX, INCX)
CSCAL
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
subroutine cgbtf2(M, N, KL, KU, AB, LDAB, IPIV, INFO)
CGBTF2 computes the LU factorization of a general band matrix using the unblocked version of the algo...
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
subroutine cgbtrf(M, N, KL, KU, AB, LDAB, IPIV, INFO)
CGBTRF
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
subroutine cgeru(M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
CGERU