294 SUBROUTINE ctgsyl( TRANS, IJOB, M, N, A, LDA, B, LDB, C, LDC, D,
295 $ ldd, e,
lde, f, ldf, scale, dif, work, lwork,
305 INTEGER ijob, info, lda, ldb, ldc, ldd,
lde, ldf,
311 COMPLEX a( lda, * ),
b( ldb, * ), c( ldc, * ),
312 $ d( ldd, * ), e(
lde, * ), f( ldf, * ),
322 parameter( zero = 0.0e+0, one = 1.0e+0 )
324 parameter( czero = (0.0e+0, 0.0e+0) )
327 LOGICAL lquery, notran
328 INTEGER i, ie, ifunc, iround, is, isolve,
j, je,
js, k,
329 $ linfo, lwmin, mb, nb, p, pq, q
330 REAL dscale, dsum, scale2, scaloc
341 INTRINSIC cmplx, max,
REAL, sqrt
348 notran =
lsame( trans,
'N' )
349 lquery = ( lwork.EQ.-1 )
351 IF( .NOT.notran .AND. .NOT.
lsame( trans,
'C' ) )
THEN
353 ELSE IF( notran )
THEN
354 IF( ( ijob.LT.0 ) .OR. ( ijob.GT.4 ) )
THEN
361 ELSE IF( n.LE.0 )
THEN
363 ELSE IF( lda.LT.max( 1, m ) )
THEN
365 ELSE IF( ldb.LT.max( 1, n ) )
THEN
367 ELSE IF( ldc.LT.max( 1, m ) )
THEN
369 ELSE IF( ldd.LT.max( 1, m ) )
THEN
371 ELSE IF(
lde.LT.max( 1, n ) )
THEN
373 ELSE IF( ldf.LT.max( 1, m ) )
THEN
380 IF( ijob.EQ.1 .OR. ijob.EQ.2 )
THEN
381 lwmin = max( 1, 2*m*n )
390 IF( lwork.LT.lwmin .AND. .NOT.lquery )
THEN
396 CALL
xerbla(
'CTGSYL', -info )
398 ELSE IF( lquery )
THEN
404 IF( m.EQ.0 .OR. n.EQ.0 )
THEN
416 mb =
ilaenv( 2,
'CTGSYL', trans, m, n, -1, -1 )
417 nb =
ilaenv( 5,
'CTGSYL', trans, m, n, -1, -1 )
424 CALL
claset(
'F', m, n, czero, czero, c, ldc )
425 CALL
claset(
'F', m, n, czero, czero, f, ldf )
426 ELSE IF( ijob.GE.1 .AND. notran )
THEN
431 IF( ( mb.LE.1 .AND. nb.LE.1 ) .OR. ( mb.GE.m .AND. nb.GE.n ) )
436 DO 30 iround = 1, isolve
442 CALL
ctgsy2( trans, ifunc, m, n, a, lda,
b, ldb, c, ldc, d,
443 $ ldd, e,
lde, f, ldf, scale, dsum, dscale,
445 IF( dscale.NE.zero )
THEN
446 IF( ijob.EQ.1 .OR. ijob.EQ.3 )
THEN
447 dif = sqrt(
REAL( 2*M*N ) ) / ( dscale*sqrt( dsum ) )
449 dif = sqrt(
REAL( PQ ) ) / ( dscale*sqrt( dsum ) )
452 IF( isolve.EQ.2 .AND. iround.EQ.1 )
THEN
457 CALL
clacpy(
'F', m, n, c, ldc, work, m )
458 CALL
clacpy(
'F', m, n, f, ldf, work( m*n+1 ), m )
459 CALL
claset(
'F', m, n, czero, czero, c, ldc )
460 CALL
claset(
'F', m, n, czero, czero, f, ldf )
461 ELSE IF( isolve.EQ.2 .AND. iround.EQ.2 )
THEN
462 CALL
clacpy(
'F', m, n, work, m, c, ldc )
463 CALL
clacpy(
'F', m, n, work( m*n+1 ), m, f, ldf )
487 IF( iwork( p ).EQ.iwork( p+1 ) )
507 IF( iwork( q ).EQ.iwork( q+1 ) )
511 DO 150 iround = 1, isolve
524 je = iwork(
j+1 ) - 1
528 ie = iwork( i+1 ) - 1
530 CALL
ctgsy2( trans, ifunc, mb, nb, a( is, is ), lda,
531 $
b(
js,
js ), ldb, c( is,
js ), ldc,
532 $ d( is, is ), ldd, e(
js,
js ),
lde,
533 $ f( is,
js ), ldf, scaloc, dsum, dscale,
538 IF( scaloc.NE.one )
THEN
540 CALL
cscal( m, cmplx( scaloc, zero ), c( 1, k ),
542 CALL
cscal( m, cmplx( scaloc, zero ), f( 1, k ),
546 CALL
cscal( is-1, cmplx( scaloc, zero ),
548 CALL
cscal( is-1, cmplx( scaloc, zero ),
552 CALL
cscal( m-ie, cmplx( scaloc, zero ),
554 CALL
cscal( m-ie, cmplx( scaloc, zero ),
558 CALL
cscal( m, cmplx( scaloc, zero ), c( 1, k ),
560 CALL
cscal( m, cmplx( scaloc, zero ), f( 1, k ),
569 CALL
cgemm(
'N',
'N', is-1, nb, mb,
570 $ cmplx( -one, zero ), a( 1, is ), lda,
571 $ c( is,
js ), ldc, cmplx( one, zero ),
573 CALL
cgemm(
'N',
'N', is-1, nb, mb,
574 $ cmplx( -one, zero ), d( 1, is ), ldd,
575 $ c( is,
js ), ldc, cmplx( one, zero ),
579 CALL
cgemm(
'N',
'N', mb, n-je, nb,
580 $ cmplx( one, zero ), f( is,
js ), ldf,
581 $
b(
js, je+1 ), ldb, cmplx( one, zero ),
582 $ c( is, je+1 ), ldc )
583 CALL
cgemm(
'N',
'N', mb, n-je, nb,
584 $ cmplx( one, zero ), f( is,
js ), ldf,
585 $ e(
js, je+1 ),
lde, cmplx( one, zero ),
586 $ f( is, je+1 ), ldf )
590 IF( dscale.NE.zero )
THEN
591 IF( ijob.EQ.1 .OR. ijob.EQ.3 )
THEN
592 dif = sqrt(
REAL( 2*M*N ) ) / ( dscale*sqrt( dsum ) )
594 dif = sqrt(
REAL( PQ ) ) / ( dscale*sqrt( dsum ) )
597 IF( isolve.EQ.2 .AND. iround.EQ.1 )
THEN
602 CALL
clacpy(
'F', m, n, c, ldc, work, m )
603 CALL
clacpy(
'F', m, n, f, ldf, work( m*n+1 ), m )
604 CALL
claset(
'F', m, n, czero, czero, c, ldc )
605 CALL
claset(
'F', m, n, czero, czero, f, ldf )
606 ELSE IF( isolve.EQ.2 .AND. iround.EQ.2 )
THEN
607 CALL
clacpy(
'F', m, n, work, m, c, ldc )
608 CALL
clacpy(
'F', m, n, work( m*n+1 ), m, f, ldf )
622 ie = iwork( i+1 ) - 1
624 DO 200
j = q, p + 2, -1
626 je = iwork(
j+1 ) - 1
628 CALL
ctgsy2( trans, ifunc, mb, nb, a( is, is ), lda,
629 $
b(
js,
js ), ldb, c( is,
js ), ldc,
630 $ d( is, is ), ldd, e(
js,
js ),
lde,
631 $ f( is,
js ), ldf, scaloc, dsum, dscale,
635 IF( scaloc.NE.one )
THEN
637 CALL
cscal( m, cmplx( scaloc, zero ), c( 1, k ),
639 CALL
cscal( m, cmplx( scaloc, zero ), f( 1, k ),
643 CALL
cscal( is-1, cmplx( scaloc, zero ), c( 1, k ),
645 CALL
cscal( is-1, cmplx( scaloc, zero ), f( 1, k ),
649 CALL
cscal( m-ie, cmplx( scaloc, zero ),
651 CALL
cscal( m-ie, cmplx( scaloc, zero ),
655 CALL
cscal( m, cmplx( scaloc, zero ), c( 1, k ),
657 CALL
cscal( m, cmplx( scaloc, zero ), f( 1, k ),
666 CALL
cgemm(
'N',
'C', mb,
js-1, nb,
667 $ cmplx( one, zero ), c( is,
js ), ldc,
668 $
b( 1,
js ), ldb, cmplx( one, zero ),
670 CALL
cgemm(
'N',
'C', mb,
js-1, nb,
671 $ cmplx( one, zero ), f( is,
js ), ldf,
672 $ e( 1,
js ),
lde, cmplx( one, zero ),
676 CALL
cgemm(
'C',
'N', m-ie, nb, mb,
677 $ cmplx( -one, zero ), a( is, ie+1 ), lda,
678 $ c( is,
js ), ldc, cmplx( one, zero ),
679 $ c( ie+1,
js ), ldc )
680 CALL
cgemm(
'C',
'N', m-ie, nb, mb,
681 $ cmplx( -one, zero ), d( is, ie+1 ), ldd,
682 $ f( is,
js ), ldf, cmplx( one, zero ),
683 $ c( ie+1,
js ), ldc )
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...
LOGICAL function lsame(CA, CB)
LSAME
subroutine cscal(N, CA, CX, INCX)
CSCAL
LOGICAL function lde(RI, RJ, LR)
subroutine xerbla(SRNAME, INFO)
XERBLA
for consistency with libblas dev *Do not ship a copy of jquery js
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real b(3) integer i
subroutine ctgsy2(TRANS, IJOB, M, N, A, LDA, B, LDB, C, LDC, D, LDD, E, LDE, F, LDF, SCALE, RDSUM, RDSCAL, INFO)
CTGSY2 solves the generalized Sylvester equation (unblocked algorithm).
INTEGER function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
subroutine ctgsyl(TRANS, IJOB, M, N, A, LDA, B, LDB, C, LDC, D, LDD, E, LDE, F, LDF, SCALE, DIF, WORK, LWORK, IWORK, INFO)
CTGSYL
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM