1*> \brief \b DCHKGK 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8* Definition: 9* =========== 10* 11* SUBROUTINE DCHKGK( NIN, NOUT ) 12* 13* .. Scalar Arguments .. 14* INTEGER NIN, NOUT 15* .. 16* 17* 18*> \par Purpose: 19* ============= 20*> 21*> \verbatim 22*> 23*> DCHKGK tests DGGBAK, a routine for backward balancing of 24*> a matrix pair (A, B). 25*> \endverbatim 26* 27* Arguments: 28* ========== 29* 30*> \param[in] NIN 31*> \verbatim 32*> NIN is INTEGER 33*> The logical unit number for input. NIN > 0. 34*> \endverbatim 35*> 36*> \param[in] NOUT 37*> \verbatim 38*> NOUT is INTEGER 39*> The logical unit number for output. NOUT > 0. 40*> \endverbatim 41* 42* Authors: 43* ======== 44* 45*> \author Univ. of Tennessee 46*> \author Univ. of California Berkeley 47*> \author Univ. of Colorado Denver 48*> \author NAG Ltd. 49* 50*> \date December 2016 51* 52*> \ingroup double_eig 53* 54* ===================================================================== 55 SUBROUTINE DCHKGK( NIN, NOUT ) 56* 57* -- LAPACK test routine (version 3.7.0) -- 58* -- LAPACK is a software package provided by Univ. of Tennessee, -- 59* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 60* December 2016 61* 62* .. Scalar Arguments .. 63 INTEGER NIN, NOUT 64* .. 65* 66* ===================================================================== 67* 68* .. Parameters .. 69 INTEGER LDA, LDB, LDVL, LDVR 70 PARAMETER ( LDA = 50, LDB = 50, LDVL = 50, LDVR = 50 ) 71 INTEGER LDE, LDF, LDWORK 72 PARAMETER ( LDE = 50, LDF = 50, LDWORK = 50 ) 73 DOUBLE PRECISION ZERO, ONE 74 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 75* .. 76* .. Local Scalars .. 77 INTEGER I, IHI, ILO, INFO, J, KNT, M, N, NINFO 78 DOUBLE PRECISION ANORM, BNORM, EPS, RMAX, VMAX 79* .. 80* .. Local Arrays .. 81 INTEGER LMAX( 4 ) 82 DOUBLE PRECISION A( LDA, LDA ), AF( LDA, LDA ), B( LDB, LDB ), 83 $ BF( LDB, LDB ), E( LDE, LDE ), F( LDF, LDF ), 84 $ LSCALE( LDA ), RSCALE( LDA ), VL( LDVL, LDVL ), 85 $ VLF( LDVL, LDVL ), VR( LDVR, LDVR ), 86 $ VRF( LDVR, LDVR ), WORK( LDWORK, LDWORK ) 87* .. 88* .. External Functions .. 89 DOUBLE PRECISION DLAMCH, DLANGE 90 EXTERNAL DLAMCH, DLANGE 91* .. 92* .. External Subroutines .. 93 EXTERNAL DGEMM, DGGBAK, DGGBAL, DLACPY 94* .. 95* .. Intrinsic Functions .. 96 INTRINSIC ABS, MAX 97* .. 98* .. Executable Statements .. 99* 100* Initialization 101* 102 LMAX( 1 ) = 0 103 LMAX( 2 ) = 0 104 LMAX( 3 ) = 0 105 LMAX( 4 ) = 0 106 NINFO = 0 107 KNT = 0 108 RMAX = ZERO 109* 110 EPS = DLAMCH( 'Precision' ) 111* 112 10 CONTINUE 113 READ( NIN, FMT = * )N, M 114 IF( N.EQ.0 ) 115 $ GO TO 100 116* 117 DO 20 I = 1, N 118 READ( NIN, FMT = * )( A( I, J ), J = 1, N ) 119 20 CONTINUE 120* 121 DO 30 I = 1, N 122 READ( NIN, FMT = * )( B( I, J ), J = 1, N ) 123 30 CONTINUE 124* 125 DO 40 I = 1, N 126 READ( NIN, FMT = * )( VL( I, J ), J = 1, M ) 127 40 CONTINUE 128* 129 DO 50 I = 1, N 130 READ( NIN, FMT = * )( VR( I, J ), J = 1, M ) 131 50 CONTINUE 132* 133 KNT = KNT + 1 134* 135 ANORM = DLANGE( 'M', N, N, A, LDA, WORK ) 136 BNORM = DLANGE( 'M', N, N, B, LDB, WORK ) 137* 138 CALL DLACPY( 'FULL', N, N, A, LDA, AF, LDA ) 139 CALL DLACPY( 'FULL', N, N, B, LDB, BF, LDB ) 140* 141 CALL DGGBAL( 'B', N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, 142 $ WORK, INFO ) 143 IF( INFO.NE.0 ) THEN 144 NINFO = NINFO + 1 145 LMAX( 1 ) = KNT 146 END IF 147* 148 CALL DLACPY( 'FULL', N, M, VL, LDVL, VLF, LDVL ) 149 CALL DLACPY( 'FULL', N, M, VR, LDVR, VRF, LDVR ) 150* 151 CALL DGGBAK( 'B', 'L', N, ILO, IHI, LSCALE, RSCALE, M, VL, LDVL, 152 $ INFO ) 153 IF( INFO.NE.0 ) THEN 154 NINFO = NINFO + 1 155 LMAX( 2 ) = KNT 156 END IF 157* 158 CALL DGGBAK( 'B', 'R', N, ILO, IHI, LSCALE, RSCALE, M, VR, LDVR, 159 $ INFO ) 160 IF( INFO.NE.0 ) THEN 161 NINFO = NINFO + 1 162 LMAX( 3 ) = KNT 163 END IF 164* 165* Test of DGGBAK 166* 167* Check tilde(VL)'*A*tilde(VR) - VL'*tilde(A)*VR 168* where tilde(A) denotes the transformed matrix. 169* 170 CALL DGEMM( 'N', 'N', N, M, N, ONE, AF, LDA, VR, LDVR, ZERO, WORK, 171 $ LDWORK ) 172 CALL DGEMM( 'T', 'N', M, M, N, ONE, VL, LDVL, WORK, LDWORK, ZERO, 173 $ E, LDE ) 174* 175 CALL DGEMM( 'N', 'N', N, M, N, ONE, A, LDA, VRF, LDVR, ZERO, WORK, 176 $ LDWORK ) 177 CALL DGEMM( 'T', 'N', M, M, N, ONE, VLF, LDVL, WORK, LDWORK, ZERO, 178 $ F, LDF ) 179* 180 VMAX = ZERO 181 DO 70 J = 1, M 182 DO 60 I = 1, M 183 VMAX = MAX( VMAX, ABS( E( I, J )-F( I, J ) ) ) 184 60 CONTINUE 185 70 CONTINUE 186 VMAX = VMAX / ( EPS*MAX( ANORM, BNORM ) ) 187 IF( VMAX.GT.RMAX ) THEN 188 LMAX( 4 ) = KNT 189 RMAX = VMAX 190 END IF 191* 192* Check tilde(VL)'*B*tilde(VR) - VL'*tilde(B)*VR 193* 194 CALL DGEMM( 'N', 'N', N, M, N, ONE, BF, LDB, VR, LDVR, ZERO, WORK, 195 $ LDWORK ) 196 CALL DGEMM( 'T', 'N', M, M, N, ONE, VL, LDVL, WORK, LDWORK, ZERO, 197 $ E, LDE ) 198* 199 CALL DGEMM( 'N', 'N', N, M, N, ONE, B, LDB, VRF, LDVR, ZERO, WORK, 200 $ LDWORK ) 201 CALL DGEMM( 'T', 'N', M, M, N, ONE, VLF, LDVL, WORK, LDWORK, ZERO, 202 $ F, LDF ) 203* 204 VMAX = ZERO 205 DO 90 J = 1, M 206 DO 80 I = 1, M 207 VMAX = MAX( VMAX, ABS( E( I, J )-F( I, J ) ) ) 208 80 CONTINUE 209 90 CONTINUE 210 VMAX = VMAX / ( EPS*MAX( ANORM, BNORM ) ) 211 IF( VMAX.GT.RMAX ) THEN 212 LMAX( 4 ) = KNT 213 RMAX = VMAX 214 END IF 215* 216 GO TO 10 217* 218 100 CONTINUE 219* 220 WRITE( NOUT, FMT = 9999 ) 221 9999 FORMAT( 1X, '.. test output of DGGBAK .. ' ) 222* 223 WRITE( NOUT, FMT = 9998 )RMAX 224 9998 FORMAT( ' value of largest test error =', D12.3 ) 225 WRITE( NOUT, FMT = 9997 )LMAX( 1 ) 226 9997 FORMAT( ' example number where DGGBAL info is not 0 =', I4 ) 227 WRITE( NOUT, FMT = 9996 )LMAX( 2 ) 228 9996 FORMAT( ' example number where DGGBAK(L) info is not 0 =', I4 ) 229 WRITE( NOUT, FMT = 9995 )LMAX( 3 ) 230 9995 FORMAT( ' example number where DGGBAK(R) info is not 0 =', I4 ) 231 WRITE( NOUT, FMT = 9994 )LMAX( 4 ) 232 9994 FORMAT( ' example number having largest error =', I4 ) 233 WRITE( NOUT, FMT = 9993 )NINFO 234 9993 FORMAT( ' number of examples where info is not 0 =', I4 ) 235 WRITE( NOUT, FMT = 9992 )KNT 236 9992 FORMAT( ' total number of examples tested =', I4 ) 237* 238 RETURN 239* 240* End of DCHKGK 241* 242 END 243