54
   55
   56
   57
   58
   59
   60      INTEGER            NIN, NOUT
   61
   62
   63
   64
   65
   66      INTEGER            LDA, LDB, LDVL, LDVR
   67      parameter( lda = 50, ldb = 50, ldvl = 50, ldvr = 50 )
   68      INTEGER            LDE, LDF, LDWORK, LRWORK
   69      parameter( 
lde = 50, ldf = 50, ldwork = 50,
 
   70     $                   lrwork = 6*50 )
   71      REAL               ZERO
   72      parameter( zero = 0.0e+0 )
   73      COMPLEX            CZERO, CONE
   74      parameter( czero = ( 0.0e+0, 0.0e+0 ),
   75     $                   cone = ( 1.0e+0, 0.0e+0 ) )
   76
   77
   78      INTEGER            I, IHI, ILO, INFO, J, KNT, M, N, NINFO
   79      REAL               ANORM, BNORM, EPS, RMAX, VMAX
   80      COMPLEX            CDUM
   81
   82
   83      INTEGER            LMAX( 4 )
   84      REAL               LSCALE( LDA ), RSCALE( LDA ), RWORK( LRWORK )
   85      COMPLEX            A( LDA, LDA ), AF( LDA, LDA ), B( LDB, LDB ),
   86     $                   BF( LDB, LDB ), E( LDE, LDE ), F( LDF, LDF ),
   87     $                   VL( LDVL, LDVL ), VLF( LDVL, LDVL ),
   88     $                   VR( LDVR, LDVR ), VRF( LDVR, LDVR ),
   89     $                   WORK( LDWORK, LDWORK )
   90
   91
   92      REAL               CLANGE, SLAMCH
   94
   95
   97
   98
   99      INTRINSIC          abs, aimag, max, real
  100
  101
  102      REAL               CABS1
  103
  104
  105      cabs1( cdum ) = abs( real( cdum ) ) + abs( aimag( cdum ) )
  106
  107
  108
  109      lmax( 1 ) = 0
  110      lmax( 2 ) = 0
  111      lmax( 3 ) = 0
  112      lmax( 4 ) = 0
  113      ninfo = 0
  114      knt = 0
  115      rmax = zero
  116
  117      eps = 
slamch( 
'Precision' )
 
  118
  119   10 CONTINUE
  120      READ( nin, fmt = * )n, m
  121      IF( n.EQ.0 )
  122     $   GO TO 100
  123
  124      DO 20 i = 1, n
  125         READ( nin, fmt = * )( a( i, j ), j = 1, n )
  126   20 CONTINUE
  127
  128      DO 30 i = 1, n
  129         READ( nin, fmt = * )( b( i, j ), j = 1, n )
  130   30 CONTINUE
  131
  132      DO 40 i = 1, n
  133         READ( nin, fmt = * )( vl( i, j ), j = 1, m )
  134   40 CONTINUE
  135
  136      DO 50 i = 1, n
  137         READ( nin, fmt = * )( vr( i, j ), j = 1, m )
  138   50 CONTINUE
  139
  140      knt = knt + 1
  141
  142      anorm = 
clange( 
'M', n, n, a, lda, rwork )
 
  143      bnorm = 
clange( 
'M', n, n, b, ldb, rwork )
 
  144
  145      CALL clacpy( 
'FULL', n, n, a, lda, af, lda )
 
  146      CALL clacpy( 
'FULL', n, n, b, ldb, bf, ldb )
 
  147
  148      CALL cggbal( 
'B', n, a, lda, b, ldb, ilo, ihi, lscale, rscale,
 
  149     $             rwork, info )
  150      IF( info.NE.0 ) THEN
  151         ninfo = ninfo + 1
  152         lmax( 1 ) = knt
  153      END IF
  154
  155      CALL clacpy( 
'FULL', n, m, vl, ldvl, vlf, ldvl )
 
  156      CALL clacpy( 
'FULL', n, m, vr, ldvr, vrf, ldvr )
 
  157
  158      CALL cggbak( 
'B', 
'L', n, ilo, ihi, lscale, rscale, m, vl, ldvl,
 
  159     $             info )
  160      IF( info.NE.0 ) THEN
  161         ninfo = ninfo + 1
  162         lmax( 2 ) = knt
  163      END IF
  164
  165      CALL cggbak( 
'B', 
'R', n, ilo, ihi, lscale, rscale, m, vr, ldvr,
 
  166     $             info )
  167      IF( info.NE.0 ) THEN
  168         ninfo = ninfo + 1
  169         lmax( 3 ) = knt
  170      END IF
  171
  172
  173
  174
  175
  176
  177      CALL cgemm( 
'N', 
'N', n, m, n, cone, af, lda, vr, ldvr, czero,
 
  178     $            work, ldwork )
  179      CALL cgemm( 
'C', 
'N', m, m, n, cone, vl, ldvl, work, ldwork,
 
  181
  182      CALL cgemm( 
'N', 
'N', n, m, n, cone, a, lda, vrf, ldvr, czero,
 
  183     $            work, ldwork )
  184      CALL cgemm( 
'C', 
'N', m, m, n, cone, vlf, ldvl, work, ldwork,
 
  185     $            czero, f, ldf )
  186
  187      vmax = zero
  188      DO 70 j = 1, m
  189         DO 60 i = 1, m
  190            vmax = max( vmax, cabs1( e( i, j )-f( i, j ) ) )
  191   60    CONTINUE
  192   70 CONTINUE
  193      vmax = vmax / ( eps*max( anorm, bnorm ) )
  194      IF( vmax.GT.rmax ) THEN
  195         lmax( 4 ) = knt
  196         rmax = vmax
  197      END IF
  198
  199
  200
  201      CALL cgemm( 
'N', 
'N', n, m, n, cone, bf, ldb, vr, ldvr, czero,
 
  202     $            work, ldwork )
  203      CALL cgemm( 
'C', 
'N', m, m, n, cone, vl, ldvl, work, ldwork,
 
  205
  206      CALL cgemm( 
'n', 
'n', n, m, n, cone, b, ldb, vrf, ldvr, czero,
 
  207     $            work, ldwork )
  208      CALL cgemm( 
'C', 
'N', m, m, n, cone, vlf, ldvl, work, ldwork,
 
  209     $            czero, f, ldf )
  210
  211      vmax = zero
  212      DO 90 j = 1, m
  213         DO 80 i = 1, m
  214            vmax = max( vmax, cabs1( e( i, j )-f( i, j ) ) )
  215   80    CONTINUE
  216   90 CONTINUE
  217      vmax = vmax / ( eps*max( anorm, bnorm ) )
  218      IF( vmax.GT.rmax ) THEN
  219         lmax( 4 ) = knt
  220         rmax = vmax
  221      END IF
  222
  223      GO TO 10
  224
  225  100 CONTINUE
  226
  227      WRITE( nout, fmt = 9999 )
  228 9999 FORMAT( 1x, '.. test output of CGGBAK .. ' )
  229
  230      WRITE( nout, fmt = 9998 )rmax
  231 9998 FORMAT( ' value of largest test error                  =', e12.3 )
  232      WRITE( nout, fmt = 9997 )lmax( 1 )
  233 9997 FORMAT( ' example number where CGGBAL info is not 0    =', i4 )
  234      WRITE( nout, fmt = 9996 )lmax( 2 )
  235 9996 FORMAT( ' example number where CGGBAK(L) info is not 0 =', i4 )
  236      WRITE( nout, fmt = 9995 )lmax( 3 )
  237 9995 FORMAT( ' example number where CGGBAK(R) info is not 0 =', i4 )
  238      WRITE( nout, fmt = 9994 )lmax( 4 )
  239 9994 FORMAT( ' example number having largest error          =', i4 )
  240      WRITE( nout, fmt = 9992 )ninfo
  241 9992 FORMAT( ' number of examples where info is not 0       =', i4 )
  242      WRITE( nout, fmt = 9991 )knt
  243 9991 FORMAT( ' total number of examples tested              =', i4 )
  244
  245      RETURN
  246
  247
  248
logical function lde(ri, rj, lr)
subroutine cgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
CGEMM
subroutine cggbak(job, side, n, ilo, ihi, lscale, rscale, m, v, ldv, info)
CGGBAK
subroutine cggbal(job, n, a, lda, b, ldb, ilo, ihi, lscale, rscale, work, info)
CGGBAL
subroutine clacpy(uplo, m, n, a, lda, b, ldb)
CLACPY copies all or part of one two-dimensional array to another.
real function slamch(cmach)
SLAMCH
real function clange(norm, m, n, a, lda, work)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...