1 SUBROUTINE SCHKGK( NIN, NOUT )
2 *
3 * -- LAPACK test routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8 INTEGER NIN, NOUT
9 * ..
10 *
11 * Purpose
12 * =======
13 *
14 * SCHKGK tests SGGBAK, a routine for backward balancing of
15 * a matrix pair (A, B).
16 *
17 * Arguments
18 * =========
19 *
20 * NIN (input) INTEGER
21 * The logical unit number for input. NIN > 0.
22 *
23 * NOUT (input) INTEGER
24 * The logical unit number for output. NOUT > 0.
25 *
26 * =====================================================================
27 *
28 * .. Parameters ..
29 INTEGER LDA, LDB, LDVL, LDVR
30 PARAMETER ( LDA = 50, LDB = 50, LDVL = 50, LDVR = 50 )
31 INTEGER LDE, LDF, LDWORK
32 PARAMETER ( LDE = 50, LDF = 50, LDWORK = 50 )
33 REAL ZERO, ONE
34 PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 )
35 * ..
36 * .. Local Scalars ..
37 INTEGER I, IHI, ILO, INFO, J, KNT, M, N, NINFO
38 REAL ANORM, BNORM, EPS, RMAX, VMAX
39 * ..
40 * .. Local Arrays ..
41 INTEGER LMAX( 4 )
42 REAL A( LDA, LDA ), AF( LDA, LDA ), B( LDB, LDB ),
43 $ BF( LDB, LDB ), E( LDE, LDE ), F( LDF, LDF ),
44 $ LSCALE( LDA ), RSCALE( LDA ), VL( LDVL, LDVL ),
45 $ VLF( LDVL, LDVL ), VR( LDVR, LDVR ),
46 $ VRF( LDVR, LDVR ), WORK( LDWORK, LDWORK )
47 * ..
48 * .. External Functions ..
49 REAL SLAMCH, SLANGE
50 EXTERNAL SLAMCH, SLANGE
51 * ..
52 * .. External Subroutines ..
53 EXTERNAL SGEMM, SGGBAK, SGGBAL, SLACPY
54 * ..
55 * .. Intrinsic Functions ..
56 INTRINSIC ABS, MAX
57 * ..
58 * .. Executable Statements ..
59 *
60 * Initialization
61 *
62 LMAX( 1 ) = 0
63 LMAX( 2 ) = 0
64 LMAX( 3 ) = 0
65 LMAX( 4 ) = 0
66 NINFO = 0
67 KNT = 0
68 RMAX = ZERO
69 *
70 EPS = SLAMCH( 'Precision' )
71 *
72 10 CONTINUE
73 READ( NIN, FMT = * )N, M
74 IF( N.EQ.0 )
75 $ GO TO 100
76 *
77 DO 20 I = 1, N
78 READ( NIN, FMT = * )( A( I, J ), J = 1, N )
79 20 CONTINUE
80 *
81 DO 30 I = 1, N
82 READ( NIN, FMT = * )( B( I, J ), J = 1, N )
83 30 CONTINUE
84 *
85 DO 40 I = 1, N
86 READ( NIN, FMT = * )( VL( I, J ), J = 1, M )
87 40 CONTINUE
88 *
89 DO 50 I = 1, N
90 READ( NIN, FMT = * )( VR( I, J ), J = 1, M )
91 50 CONTINUE
92 *
93 KNT = KNT + 1
94 *
95 ANORM = SLANGE( 'M', N, N, A, LDA, WORK )
96 BNORM = SLANGE( 'M', N, N, B, LDB, WORK )
97 *
98 CALL SLACPY( 'FULL', N, N, A, LDA, AF, LDA )
99 CALL SLACPY( 'FULL', N, N, B, LDB, BF, LDB )
100 *
101 CALL SGGBAL( 'B', N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE,
102 $ WORK, INFO )
103 IF( INFO.NE.0 ) THEN
104 NINFO = NINFO + 1
105 LMAX( 1 ) = KNT
106 END IF
107 *
108 CALL SLACPY( 'FULL', N, M, VL, LDVL, VLF, LDVL )
109 CALL SLACPY( 'FULL', N, M, VR, LDVR, VRF, LDVR )
110 *
111 CALL SGGBAK( 'B', 'L', N, ILO, IHI, LSCALE, RSCALE, M, VL, LDVL,
112 $ INFO )
113 IF( INFO.NE.0 ) THEN
114 NINFO = NINFO + 1
115 LMAX( 2 ) = KNT
116 END IF
117 *
118 CALL SGGBAK( 'B', 'R', N, ILO, IHI, LSCALE, RSCALE, M, VR, LDVR,
119 $ INFO )
120 IF( INFO.NE.0 ) THEN
121 NINFO = NINFO + 1
122 LMAX( 3 ) = KNT
123 END IF
124 *
125 * Test of SGGBAK
126 *
127 * Check tilde(VL)'*A*tilde(VR) - VL'*tilde(A)*VR
128 * where tilde(A) denotes the transformed matrix.
129 *
130 CALL SGEMM( 'N', 'N', N, M, N, ONE, AF, LDA, VR, LDVR, ZERO, WORK,
131 $ LDWORK )
132 CALL SGEMM( 'T', 'N', M, M, N, ONE, VL, LDVL, WORK, LDWORK, ZERO,
133 $ E, LDE )
134 *
135 CALL SGEMM( 'N', 'N', N, M, N, ONE, A, LDA, VRF, LDVR, ZERO, WORK,
136 $ LDWORK )
137 CALL SGEMM( 'T', 'N', M, M, N, ONE, VLF, LDVL, WORK, LDWORK, ZERO,
138 $ F, LDF )
139 *
140 VMAX = ZERO
141 DO 70 J = 1, M
142 DO 60 I = 1, M
143 VMAX = MAX( VMAX, ABS( E( I, J )-F( I, J ) ) )
144 60 CONTINUE
145 70 CONTINUE
146 VMAX = VMAX / ( EPS*MAX( ANORM, BNORM ) )
147 IF( VMAX.GT.RMAX ) THEN
148 LMAX( 4 ) = KNT
149 RMAX = VMAX
150 END IF
151 *
152 * Check tilde(VL)'*B*tilde(VR) - VL'*tilde(B)*VR
153 *
154 CALL SGEMM( 'N', 'N', N, M, N, ONE, BF, LDB, VR, LDVR, ZERO, WORK,
155 $ LDWORK )
156 CALL SGEMM( 'T', 'N', M, M, N, ONE, VL, LDVL, WORK, LDWORK, ZERO,
157 $ E, LDE )
158 *
159 CALL SGEMM( 'N', 'N', N, M, N, ONE, B, LDB, VRF, LDVR, ZERO, WORK,
160 $ LDWORK )
161 CALL SGEMM( 'T', 'N', M, M, N, ONE, VLF, LDVL, WORK, LDWORK, ZERO,
162 $ F, LDF )
163 *
164 VMAX = ZERO
165 DO 90 J = 1, M
166 DO 80 I = 1, M
167 VMAX = MAX( VMAX, ABS( E( I, J )-F( I, J ) ) )
168 80 CONTINUE
169 90 CONTINUE
170 VMAX = VMAX / ( EPS*MAX( ANORM, BNORM ) )
171 IF( VMAX.GT.RMAX ) THEN
172 LMAX( 4 ) = KNT
173 RMAX = VMAX
174 END IF
175 *
176 GO TO 10
177 *
178 100 CONTINUE
179 *
180 WRITE( NOUT, FMT = 9999 )
181 9999 FORMAT( 1X, '.. test output of SGGBAK .. ' )
182 *
183 WRITE( NOUT, FMT = 9998 )RMAX
184 9998 FORMAT( ' value of largest test error =', E12.3 )
185 WRITE( NOUT, FMT = 9997 )LMAX( 1 )
186 9997 FORMAT( ' example number where SGGBAL info is not 0 =', I4 )
187 WRITE( NOUT, FMT = 9996 )LMAX( 2 )
188 9996 FORMAT( ' example number where SGGBAK(L) info is not 0 =', I4 )
189 WRITE( NOUT, FMT = 9995 )LMAX( 3 )
190 9995 FORMAT( ' example number where SGGBAK(R) info is not 0 =', I4 )
191 WRITE( NOUT, FMT = 9994 )LMAX( 4 )
192 9994 FORMAT( ' example number having largest error =', I4 )
193 WRITE( NOUT, FMT = 9992 )NINFO
194 9992 FORMAT( ' number of examples where info is not 0 =', I4 )
195 WRITE( NOUT, FMT = 9991 )KNT
196 9991 FORMAT( ' total number of examples tested =', I4 )
197 *
198 RETURN
199 *
200 * End of SCHKGK
201 *
202 END
2 *
3 * -- LAPACK test routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8 INTEGER NIN, NOUT
9 * ..
10 *
11 * Purpose
12 * =======
13 *
14 * SCHKGK tests SGGBAK, a routine for backward balancing of
15 * a matrix pair (A, B).
16 *
17 * Arguments
18 * =========
19 *
20 * NIN (input) INTEGER
21 * The logical unit number for input. NIN > 0.
22 *
23 * NOUT (input) INTEGER
24 * The logical unit number for output. NOUT > 0.
25 *
26 * =====================================================================
27 *
28 * .. Parameters ..
29 INTEGER LDA, LDB, LDVL, LDVR
30 PARAMETER ( LDA = 50, LDB = 50, LDVL = 50, LDVR = 50 )
31 INTEGER LDE, LDF, LDWORK
32 PARAMETER ( LDE = 50, LDF = 50, LDWORK = 50 )
33 REAL ZERO, ONE
34 PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 )
35 * ..
36 * .. Local Scalars ..
37 INTEGER I, IHI, ILO, INFO, J, KNT, M, N, NINFO
38 REAL ANORM, BNORM, EPS, RMAX, VMAX
39 * ..
40 * .. Local Arrays ..
41 INTEGER LMAX( 4 )
42 REAL A( LDA, LDA ), AF( LDA, LDA ), B( LDB, LDB ),
43 $ BF( LDB, LDB ), E( LDE, LDE ), F( LDF, LDF ),
44 $ LSCALE( LDA ), RSCALE( LDA ), VL( LDVL, LDVL ),
45 $ VLF( LDVL, LDVL ), VR( LDVR, LDVR ),
46 $ VRF( LDVR, LDVR ), WORK( LDWORK, LDWORK )
47 * ..
48 * .. External Functions ..
49 REAL SLAMCH, SLANGE
50 EXTERNAL SLAMCH, SLANGE
51 * ..
52 * .. External Subroutines ..
53 EXTERNAL SGEMM, SGGBAK, SGGBAL, SLACPY
54 * ..
55 * .. Intrinsic Functions ..
56 INTRINSIC ABS, MAX
57 * ..
58 * .. Executable Statements ..
59 *
60 * Initialization
61 *
62 LMAX( 1 ) = 0
63 LMAX( 2 ) = 0
64 LMAX( 3 ) = 0
65 LMAX( 4 ) = 0
66 NINFO = 0
67 KNT = 0
68 RMAX = ZERO
69 *
70 EPS = SLAMCH( 'Precision' )
71 *
72 10 CONTINUE
73 READ( NIN, FMT = * )N, M
74 IF( N.EQ.0 )
75 $ GO TO 100
76 *
77 DO 20 I = 1, N
78 READ( NIN, FMT = * )( A( I, J ), J = 1, N )
79 20 CONTINUE
80 *
81 DO 30 I = 1, N
82 READ( NIN, FMT = * )( B( I, J ), J = 1, N )
83 30 CONTINUE
84 *
85 DO 40 I = 1, N
86 READ( NIN, FMT = * )( VL( I, J ), J = 1, M )
87 40 CONTINUE
88 *
89 DO 50 I = 1, N
90 READ( NIN, FMT = * )( VR( I, J ), J = 1, M )
91 50 CONTINUE
92 *
93 KNT = KNT + 1
94 *
95 ANORM = SLANGE( 'M', N, N, A, LDA, WORK )
96 BNORM = SLANGE( 'M', N, N, B, LDB, WORK )
97 *
98 CALL SLACPY( 'FULL', N, N, A, LDA, AF, LDA )
99 CALL SLACPY( 'FULL', N, N, B, LDB, BF, LDB )
100 *
101 CALL SGGBAL( 'B', N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE,
102 $ WORK, INFO )
103 IF( INFO.NE.0 ) THEN
104 NINFO = NINFO + 1
105 LMAX( 1 ) = KNT
106 END IF
107 *
108 CALL SLACPY( 'FULL', N, M, VL, LDVL, VLF, LDVL )
109 CALL SLACPY( 'FULL', N, M, VR, LDVR, VRF, LDVR )
110 *
111 CALL SGGBAK( 'B', 'L', N, ILO, IHI, LSCALE, RSCALE, M, VL, LDVL,
112 $ INFO )
113 IF( INFO.NE.0 ) THEN
114 NINFO = NINFO + 1
115 LMAX( 2 ) = KNT
116 END IF
117 *
118 CALL SGGBAK( 'B', 'R', N, ILO, IHI, LSCALE, RSCALE, M, VR, LDVR,
119 $ INFO )
120 IF( INFO.NE.0 ) THEN
121 NINFO = NINFO + 1
122 LMAX( 3 ) = KNT
123 END IF
124 *
125 * Test of SGGBAK
126 *
127 * Check tilde(VL)'*A*tilde(VR) - VL'*tilde(A)*VR
128 * where tilde(A) denotes the transformed matrix.
129 *
130 CALL SGEMM( 'N', 'N', N, M, N, ONE, AF, LDA, VR, LDVR, ZERO, WORK,
131 $ LDWORK )
132 CALL SGEMM( 'T', 'N', M, M, N, ONE, VL, LDVL, WORK, LDWORK, ZERO,
133 $ E, LDE )
134 *
135 CALL SGEMM( 'N', 'N', N, M, N, ONE, A, LDA, VRF, LDVR, ZERO, WORK,
136 $ LDWORK )
137 CALL SGEMM( 'T', 'N', M, M, N, ONE, VLF, LDVL, WORK, LDWORK, ZERO,
138 $ F, LDF )
139 *
140 VMAX = ZERO
141 DO 70 J = 1, M
142 DO 60 I = 1, M
143 VMAX = MAX( VMAX, ABS( E( I, J )-F( I, J ) ) )
144 60 CONTINUE
145 70 CONTINUE
146 VMAX = VMAX / ( EPS*MAX( ANORM, BNORM ) )
147 IF( VMAX.GT.RMAX ) THEN
148 LMAX( 4 ) = KNT
149 RMAX = VMAX
150 END IF
151 *
152 * Check tilde(VL)'*B*tilde(VR) - VL'*tilde(B)*VR
153 *
154 CALL SGEMM( 'N', 'N', N, M, N, ONE, BF, LDB, VR, LDVR, ZERO, WORK,
155 $ LDWORK )
156 CALL SGEMM( 'T', 'N', M, M, N, ONE, VL, LDVL, WORK, LDWORK, ZERO,
157 $ E, LDE )
158 *
159 CALL SGEMM( 'N', 'N', N, M, N, ONE, B, LDB, VRF, LDVR, ZERO, WORK,
160 $ LDWORK )
161 CALL SGEMM( 'T', 'N', M, M, N, ONE, VLF, LDVL, WORK, LDWORK, ZERO,
162 $ F, LDF )
163 *
164 VMAX = ZERO
165 DO 90 J = 1, M
166 DO 80 I = 1, M
167 VMAX = MAX( VMAX, ABS( E( I, J )-F( I, J ) ) )
168 80 CONTINUE
169 90 CONTINUE
170 VMAX = VMAX / ( EPS*MAX( ANORM, BNORM ) )
171 IF( VMAX.GT.RMAX ) THEN
172 LMAX( 4 ) = KNT
173 RMAX = VMAX
174 END IF
175 *
176 GO TO 10
177 *
178 100 CONTINUE
179 *
180 WRITE( NOUT, FMT = 9999 )
181 9999 FORMAT( 1X, '.. test output of SGGBAK .. ' )
182 *
183 WRITE( NOUT, FMT = 9998 )RMAX
184 9998 FORMAT( ' value of largest test error =', E12.3 )
185 WRITE( NOUT, FMT = 9997 )LMAX( 1 )
186 9997 FORMAT( ' example number where SGGBAL info is not 0 =', I4 )
187 WRITE( NOUT, FMT = 9996 )LMAX( 2 )
188 9996 FORMAT( ' example number where SGGBAK(L) info is not 0 =', I4 )
189 WRITE( NOUT, FMT = 9995 )LMAX( 3 )
190 9995 FORMAT( ' example number where SGGBAK(R) info is not 0 =', I4 )
191 WRITE( NOUT, FMT = 9994 )LMAX( 4 )
192 9994 FORMAT( ' example number having largest error =', I4 )
193 WRITE( NOUT, FMT = 9992 )NINFO
194 9992 FORMAT( ' number of examples where info is not 0 =', I4 )
195 WRITE( NOUT, FMT = 9991 )KNT
196 9991 FORMAT( ' total number of examples tested =', I4 )
197 *
198 RETURN
199 *
200 * End of SCHKGK
201 *
202 END