195 SUBROUTINE dlarfb( SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV,
196 $ t, ldt, c, ldc, work, ldwork )
204 CHARACTER direct, side, storev, trans
205 INTEGER k, ldc, ldt, ldv, ldwork, m, n
208 DOUBLE PRECISION c( ldc, * ), t( ldt, * ), v( ldv, * ),
216 parameter( one = 1.0d+0 )
233 IF( m.LE.0 .OR. n.LE.0 )
236 IF(
lsame( trans,
'N' ) )
THEN
242 IF(
lsame( storev,
'C' ) )
THEN
244 IF(
lsame( direct,
'F' ) )
THEN
250 IF(
lsame( side,
'L' ) )
THEN
260 CALL
dcopy( n, c(
j, 1 ), ldc, work( 1,
j ), 1 )
265 CALL
dtrmm(
'Right',
'Lower',
'No transpose',
'Unit', n,
266 $ k, one, v, ldv, work, ldwork )
271 CALL
dgemm(
'Transpose',
'No transpose', n, k, m-k,
272 $ one, c( k+1, 1 ), ldc, v( k+1, 1 ), ldv,
273 $ one, work, ldwork )
278 CALL
dtrmm(
'Right',
'Upper', transt,
'Non-unit', n, k,
279 $ one, t, ldt, work, ldwork )
287 CALL
dgemm(
'No transpose',
'Transpose', m-k, n, k,
288 $ -one, v( k+1, 1 ), ldv, work, ldwork, one,
294 CALL
dtrmm(
'Right',
'Lower',
'Transpose',
'Unit', n, k,
295 $ one, v, ldv, work, ldwork )
301 c(
j, i ) = c(
j, i ) - work( i,
j )
305 ELSE IF(
lsame( side,
'R' ) )
THEN
314 CALL
dcopy( m, c( 1,
j ), 1, work( 1,
j ), 1 )
319 CALL
dtrmm(
'Right',
'Lower',
'No transpose',
'Unit', m,
320 $ k, one, v, ldv, work, ldwork )
325 CALL
dgemm(
'No transpose',
'No transpose', m, k, n-k,
326 $ one, c( 1, k+1 ), ldc, v( k+1, 1 ), ldv,
327 $ one, work, ldwork )
332 CALL
dtrmm(
'Right',
'Upper', trans,
'Non-unit', m, k,
333 $ one, t, ldt, work, ldwork )
341 CALL
dgemm(
'No transpose',
'Transpose', m, n-k, k,
342 $ -one, work, ldwork, v( k+1, 1 ), ldv, one,
348 CALL
dtrmm(
'Right',
'Lower',
'Transpose',
'Unit', m, k,
349 $ one, v, ldv, work, ldwork )
355 c( i,
j ) = c( i,
j ) - work( i,
j )
366 IF(
lsame( side,
'L' ) )
THEN
376 CALL
dcopy( n, c( m-k+
j, 1 ), ldc, work( 1,
j ), 1 )
381 CALL
dtrmm(
'Right',
'Upper',
'No transpose',
'Unit', n,
382 $ k, one, v( m-k+1, 1 ), ldv, work, ldwork )
387 CALL
dgemm(
'Transpose',
'No transpose', n, k, m-k,
388 $ one, c, ldc, v, ldv, one, work, ldwork )
393 CALL
dtrmm(
'Right',
'Lower', transt,
'Non-unit', n, k,
394 $ one, t, ldt, work, ldwork )
402 CALL
dgemm(
'No transpose',
'Transpose', m-k, n, k,
403 $ -one, v, ldv, work, ldwork, one, c, ldc )
408 CALL
dtrmm(
'Right',
'Upper',
'Transpose',
'Unit', n, k,
409 $ one, v( m-k+1, 1 ), ldv, work, ldwork )
415 c( m-k+
j, i ) = c( m-k+
j, i ) - work( i,
j )
419 ELSE IF(
lsame( side,
'R' ) )
THEN
428 CALL
dcopy( m, c( 1, n-k+
j ), 1, work( 1,
j ), 1 )
433 CALL
dtrmm(
'Right',
'Upper',
'No transpose',
'Unit', m,
434 $ k, one, v( n-k+1, 1 ), ldv, work, ldwork )
439 CALL
dgemm(
'No transpose',
'No transpose', m, k, n-k,
440 $ one, c, ldc, v, ldv, one, work, ldwork )
445 CALL
dtrmm(
'Right',
'Lower', trans,
'Non-unit', m, k,
446 $ one, t, ldt, work, ldwork )
454 CALL
dgemm(
'No transpose',
'Transpose', m, n-k, k,
455 $ -one, work, ldwork, v, ldv, one, c, ldc )
460 CALL
dtrmm(
'Right',
'Upper',
'Transpose',
'Unit', m, k,
461 $ one, v( n-k+1, 1 ), ldv, work, ldwork )
467 c( i, n-k+
j ) = c( i, n-k+
j ) - work( i,
j )
473 ELSE IF(
lsame( storev,
'R' ) )
THEN
475 IF(
lsame( direct,
'F' ) )
THEN
480 IF(
lsame( side,
'L' ) )
THEN
490 CALL
dcopy( n, c(
j, 1 ), ldc, work( 1,
j ), 1 )
495 CALL
dtrmm(
'Right',
'Upper',
'Transpose',
'Unit', n, k,
496 $ one, v, ldv, work, ldwork )
501 CALL
dgemm(
'Transpose',
'Transpose', n, k, m-k, one,
502 $ c( k+1, 1 ), ldc, v( 1, k+1 ), ldv, one,
508 CALL
dtrmm(
'Right',
'Upper', transt,
'Non-unit', n, k,
509 $ one, t, ldt, work, ldwork )
517 CALL
dgemm(
'Transpose',
'Transpose', m-k, n, k, -one,
518 $ v( 1, k+1 ), ldv, work, ldwork, one,
524 CALL
dtrmm(
'Right',
'Upper',
'No transpose',
'Unit', n,
525 $ k, one, v, ldv, work, ldwork )
531 c(
j, i ) = c(
j, i ) - work( i,
j )
535 ELSE IF(
lsame( side,
'R' ) )
THEN
544 CALL
dcopy( m, c( 1,
j ), 1, work( 1,
j ), 1 )
549 CALL
dtrmm(
'Right',
'Upper',
'Transpose',
'Unit', m, k,
550 $ one, v, ldv, work, ldwork )
555 CALL
dgemm(
'No transpose',
'Transpose', m, k, n-k,
556 $ one, c( 1, k+1 ), ldc, v( 1, k+1 ), ldv,
557 $ one, work, ldwork )
562 CALL
dtrmm(
'Right',
'Upper', trans,
'Non-unit', m, k,
563 $ one, t, ldt, work, ldwork )
571 CALL
dgemm(
'No transpose',
'No transpose', m, n-k, k,
572 $ -one, work, ldwork, v( 1, k+1 ), ldv, one,
578 CALL
dtrmm(
'Right',
'Upper',
'No transpose',
'Unit', m,
579 $ k, one, v, ldv, work, ldwork )
585 c( i,
j ) = c( i,
j ) - work( i,
j )
596 IF(
lsame( side,
'L' ) )
THEN
606 CALL
dcopy( n, c( m-k+
j, 1 ), ldc, work( 1,
j ), 1 )
611 CALL
dtrmm(
'Right',
'Lower',
'Transpose',
'Unit', n, k,
612 $ one, v( 1, m-k+1 ), ldv, work, ldwork )
617 CALL
dgemm(
'Transpose',
'Transpose', n, k, m-k, one,
618 $ c, ldc, v, ldv, one, work, ldwork )
623 CALL
dtrmm(
'Right',
'Lower', transt,
'Non-unit', n, k,
624 $ one, t, ldt, work, ldwork )
632 CALL
dgemm(
'Transpose',
'Transpose', m-k, n, k, -one,
633 $ v, ldv, work, ldwork, one, c, ldc )
638 CALL
dtrmm(
'Right',
'Lower',
'No transpose',
'Unit', n,
639 $ k, one, v( 1, m-k+1 ), ldv, work, ldwork )
645 c( m-k+
j, i ) = c( m-k+
j, i ) - work( i,
j )
649 ELSE IF(
lsame( side,
'R' ) )
THEN
658 CALL
dcopy( m, c( 1, n-k+
j ), 1, work( 1,
j ), 1 )
663 CALL
dtrmm(
'Right',
'Lower',
'Transpose',
'Unit', m, k,
664 $ one, v( 1, n-k+1 ), ldv, work, ldwork )
669 CALL
dgemm(
'No transpose',
'Transpose', m, k, n-k,
670 $ one, c, ldc, v, ldv, one, work, ldwork )
675 CALL
dtrmm(
'Right',
'Lower', trans,
'Non-unit', m, k,
676 $ one, t, ldt, work, ldwork )
684 CALL
dgemm(
'No transpose',
'No transpose', m, n-k, k,
685 $ -one, work, ldwork, v, ldv, one, c, ldc )
690 CALL
dtrmm(
'Right',
'Lower',
'No transpose',
'Unit', m,
691 $ k, one, v( 1, n-k+1 ), ldv, work, ldwork )
697 c( i, n-k+
j ) = c( i, n-k+
j ) - work( i,
j )
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
subroutine dlarfb(SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV, T, LDT, C, LDC, WORK, LDWORK)
DLARFB applies a block reflector or its transpose to a general rectangular matrix.
logical function lsame(CA, CB)
LSAME
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine dtrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRMM