![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b SCHKGK 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * SUBROUTINE SCHKGK( NIN, NOUT ) 00012 * 00013 * .. Scalar Arguments .. 00014 * INTEGER NIN, NOUT 00015 * .. 00016 * 00017 * 00018 *> \par Purpose: 00019 * ============= 00020 *> 00021 *> \verbatim 00022 *> 00023 *> SCHKGK tests SGGBAK, a routine for backward balancing of 00024 *> a matrix pair (A, B). 00025 *> \endverbatim 00026 * 00027 * Arguments: 00028 * ========== 00029 * 00030 *> \param[in] NIN 00031 *> \verbatim 00032 *> NIN is INTEGER 00033 *> The logical unit number for input. NIN > 0. 00034 *> \endverbatim 00035 *> 00036 *> \param[in] NOUT 00037 *> \verbatim 00038 *> NOUT is INTEGER 00039 *> The logical unit number for output. NOUT > 0. 00040 *> \endverbatim 00041 * 00042 * Authors: 00043 * ======== 00044 * 00045 *> \author Univ. of Tennessee 00046 *> \author Univ. of California Berkeley 00047 *> \author Univ. of Colorado Denver 00048 *> \author NAG Ltd. 00049 * 00050 *> \date November 2011 00051 * 00052 *> \ingroup single_eig 00053 * 00054 * ===================================================================== 00055 SUBROUTINE SCHKGK( NIN, NOUT ) 00056 * 00057 * -- LAPACK test routine (version 3.4.0) -- 00058 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00059 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00060 * November 2011 00061 * 00062 * .. Scalar Arguments .. 00063 INTEGER NIN, NOUT 00064 * .. 00065 * 00066 * ===================================================================== 00067 * 00068 * .. Parameters .. 00069 INTEGER LDA, LDB, LDVL, LDVR 00070 PARAMETER ( LDA = 50, LDB = 50, LDVL = 50, LDVR = 50 ) 00071 INTEGER LDE, LDF, LDWORK 00072 PARAMETER ( LDE = 50, LDF = 50, LDWORK = 50 ) 00073 REAL ZERO, ONE 00074 PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 ) 00075 * .. 00076 * .. Local Scalars .. 00077 INTEGER I, IHI, ILO, INFO, J, KNT, M, N, NINFO 00078 REAL ANORM, BNORM, EPS, RMAX, VMAX 00079 * .. 00080 * .. Local Arrays .. 00081 INTEGER LMAX( 4 ) 00082 REAL A( LDA, LDA ), AF( LDA, LDA ), B( LDB, LDB ), 00083 $ BF( LDB, LDB ), E( LDE, LDE ), F( LDF, LDF ), 00084 $ LSCALE( LDA ), RSCALE( LDA ), VL( LDVL, LDVL ), 00085 $ VLF( LDVL, LDVL ), VR( LDVR, LDVR ), 00086 $ VRF( LDVR, LDVR ), WORK( LDWORK, LDWORK ) 00087 * .. 00088 * .. External Functions .. 00089 REAL SLAMCH, SLANGE 00090 EXTERNAL SLAMCH, SLANGE 00091 * .. 00092 * .. External Subroutines .. 00093 EXTERNAL SGEMM, SGGBAK, SGGBAL, SLACPY 00094 * .. 00095 * .. Intrinsic Functions .. 00096 INTRINSIC ABS, MAX 00097 * .. 00098 * .. Executable Statements .. 00099 * 00100 * Initialization 00101 * 00102 LMAX( 1 ) = 0 00103 LMAX( 2 ) = 0 00104 LMAX( 3 ) = 0 00105 LMAX( 4 ) = 0 00106 NINFO = 0 00107 KNT = 0 00108 RMAX = ZERO 00109 * 00110 EPS = SLAMCH( 'Precision' ) 00111 * 00112 10 CONTINUE 00113 READ( NIN, FMT = * )N, M 00114 IF( N.EQ.0 ) 00115 $ GO TO 100 00116 * 00117 DO 20 I = 1, N 00118 READ( NIN, FMT = * )( A( I, J ), J = 1, N ) 00119 20 CONTINUE 00120 * 00121 DO 30 I = 1, N 00122 READ( NIN, FMT = * )( B( I, J ), J = 1, N ) 00123 30 CONTINUE 00124 * 00125 DO 40 I = 1, N 00126 READ( NIN, FMT = * )( VL( I, J ), J = 1, M ) 00127 40 CONTINUE 00128 * 00129 DO 50 I = 1, N 00130 READ( NIN, FMT = * )( VR( I, J ), J = 1, M ) 00131 50 CONTINUE 00132 * 00133 KNT = KNT + 1 00134 * 00135 ANORM = SLANGE( 'M', N, N, A, LDA, WORK ) 00136 BNORM = SLANGE( 'M', N, N, B, LDB, WORK ) 00137 * 00138 CALL SLACPY( 'FULL', N, N, A, LDA, AF, LDA ) 00139 CALL SLACPY( 'FULL', N, N, B, LDB, BF, LDB ) 00140 * 00141 CALL SGGBAL( 'B', N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, 00142 $ WORK, INFO ) 00143 IF( INFO.NE.0 ) THEN 00144 NINFO = NINFO + 1 00145 LMAX( 1 ) = KNT 00146 END IF 00147 * 00148 CALL SLACPY( 'FULL', N, M, VL, LDVL, VLF, LDVL ) 00149 CALL SLACPY( 'FULL', N, M, VR, LDVR, VRF, LDVR ) 00150 * 00151 CALL SGGBAK( 'B', 'L', N, ILO, IHI, LSCALE, RSCALE, M, VL, LDVL, 00152 $ INFO ) 00153 IF( INFO.NE.0 ) THEN 00154 NINFO = NINFO + 1 00155 LMAX( 2 ) = KNT 00156 END IF 00157 * 00158 CALL SGGBAK( 'B', 'R', N, ILO, IHI, LSCALE, RSCALE, M, VR, LDVR, 00159 $ INFO ) 00160 IF( INFO.NE.0 ) THEN 00161 NINFO = NINFO + 1 00162 LMAX( 3 ) = KNT 00163 END IF 00164 * 00165 * Test of SGGBAK 00166 * 00167 * Check tilde(VL)'*A*tilde(VR) - VL'*tilde(A)*VR 00168 * where tilde(A) denotes the transformed matrix. 00169 * 00170 CALL SGEMM( 'N', 'N', N, M, N, ONE, AF, LDA, VR, LDVR, ZERO, WORK, 00171 $ LDWORK ) 00172 CALL SGEMM( 'T', 'N', M, M, N, ONE, VL, LDVL, WORK, LDWORK, ZERO, 00173 $ E, LDE ) 00174 * 00175 CALL SGEMM( 'N', 'N', N, M, N, ONE, A, LDA, VRF, LDVR, ZERO, WORK, 00176 $ LDWORK ) 00177 CALL SGEMM( 'T', 'N', M, M, N, ONE, VLF, LDVL, WORK, LDWORK, ZERO, 00178 $ F, LDF ) 00179 * 00180 VMAX = ZERO 00181 DO 70 J = 1, M 00182 DO 60 I = 1, M 00183 VMAX = MAX( VMAX, ABS( E( I, J )-F( I, J ) ) ) 00184 60 CONTINUE 00185 70 CONTINUE 00186 VMAX = VMAX / ( EPS*MAX( ANORM, BNORM ) ) 00187 IF( VMAX.GT.RMAX ) THEN 00188 LMAX( 4 ) = KNT 00189 RMAX = VMAX 00190 END IF 00191 * 00192 * Check tilde(VL)'*B*tilde(VR) - VL'*tilde(B)*VR 00193 * 00194 CALL SGEMM( 'N', 'N', N, M, N, ONE, BF, LDB, VR, LDVR, ZERO, WORK, 00195 $ LDWORK ) 00196 CALL SGEMM( 'T', 'N', M, M, N, ONE, VL, LDVL, WORK, LDWORK, ZERO, 00197 $ E, LDE ) 00198 * 00199 CALL SGEMM( 'N', 'N', N, M, N, ONE, B, LDB, VRF, LDVR, ZERO, WORK, 00200 $ LDWORK ) 00201 CALL SGEMM( 'T', 'N', M, M, N, ONE, VLF, LDVL, WORK, LDWORK, ZERO, 00202 $ F, LDF ) 00203 * 00204 VMAX = ZERO 00205 DO 90 J = 1, M 00206 DO 80 I = 1, M 00207 VMAX = MAX( VMAX, ABS( E( I, J )-F( I, J ) ) ) 00208 80 CONTINUE 00209 90 CONTINUE 00210 VMAX = VMAX / ( EPS*MAX( ANORM, BNORM ) ) 00211 IF( VMAX.GT.RMAX ) THEN 00212 LMAX( 4 ) = KNT 00213 RMAX = VMAX 00214 END IF 00215 * 00216 GO TO 10 00217 * 00218 100 CONTINUE 00219 * 00220 WRITE( NOUT, FMT = 9999 ) 00221 9999 FORMAT( 1X, '.. test output of SGGBAK .. ' ) 00222 * 00223 WRITE( NOUT, FMT = 9998 )RMAX 00224 9998 FORMAT( ' value of largest test error =', E12.3 ) 00225 WRITE( NOUT, FMT = 9997 )LMAX( 1 ) 00226 9997 FORMAT( ' example number where SGGBAL info is not 0 =', I4 ) 00227 WRITE( NOUT, FMT = 9996 )LMAX( 2 ) 00228 9996 FORMAT( ' example number where SGGBAK(L) info is not 0 =', I4 ) 00229 WRITE( NOUT, FMT = 9995 )LMAX( 3 ) 00230 9995 FORMAT( ' example number where SGGBAK(R) info is not 0 =', I4 ) 00231 WRITE( NOUT, FMT = 9994 )LMAX( 4 ) 00232 9994 FORMAT( ' example number having largest error =', I4 ) 00233 WRITE( NOUT, FMT = 9992 )NINFO 00234 9992 FORMAT( ' number of examples where info is not 0 =', I4 ) 00235 WRITE( NOUT, FMT = 9991 )KNT 00236 9991 FORMAT( ' total number of examples tested =', I4 ) 00237 * 00238 RETURN 00239 * 00240 * End of SCHKGK 00241 * 00242 END