81 SUBROUTINE sget34( RMAX, LMAX, NINFO, KNT )
99 parameter( zero = 0.0e0, half = 0.5e0, one = 1.0e0 )
101 parameter( two = 2.0e0, three = 3.0e0 )
103 parameter( lwork = 32 )
106 INTEGER I, IA, IA11, IA12, IA21, IA22, IAM, IB, IC,
107 $ IC11, IC12, IC21, IC22, ICM, INFO, J
108 REAL BIGNUM, EPS, RES, SMLNUM, TNRM
111 REAL Q( 4, 4 ), RESULT( 2 ), T( 4, 4 ), T1( 4, 4 ),
112 $ VAL( 9 ), VM( 2 ), WORK( LWORK )
122 INTRINSIC abs, max, real, sign, sqrt
129 smlnum = slamch(
'S' ) / eps
130 bignum = one / smlnum
131 CALL slabad( smlnum, bignum )
136 val( 2 ) = sqrt( smlnum )
139 val( 5 ) = sqrt( bignum )
140 val( 6 ) = -sqrt( smlnum )
143 val( 9 ) = -sqrt( bignum )
145 vm( 2 ) = one + two*eps
146 CALL scopy( 16, val( 4 ), 0, t( 1, 1 ), 1 )
160 t( 1, 1 ) = val( ia )*vm( iam )
161 t( 2, 2 ) = val( ic )
162 t( 1, 2 ) = val( ib )
164 tnrm = max( abs( t( 1, 1 ) ), abs( t( 2, 2 ) ),
166 CALL scopy( 16, t, 1, t1, 1 )
167 CALL scopy( 16, val( 1 ), 0, q, 1 )
168 CALL scopy( 4, val( 3 ), 0, q, 5 )
169 CALL slaexc( .true., 2, t, 4, q, 4, 1, 1, 1, work,
172 $ ninfo( info ) = ninfo( info ) + 1
173 CALL shst01( 2, 1, 2, t1, 4, t, 4, q, 4, work, lwork,
175 res = result( 1 ) + result( 2 )
177 $ res = res + one / eps
178 IF( t( 1, 1 ).NE.t1( 2, 2 ) )
179 $ res = res + one / eps
180 IF( t( 2, 2 ).NE.t1( 1, 1 ) )
181 $ res = res + one / eps
182 IF( t( 2, 1 ).NE.zero )
183 $ res = res + one / eps
185 IF( res.GT.rmax )
THEN
200 DO 50 ic22 = -1, 1, 2
201 t( 1, 1 ) = val( ia )*vm( iam )
202 t( 1, 2 ) = val( ib )
203 t( 1, 3 ) = -two*val( ib )
205 t( 2, 2 ) = val( ic11 )
206 t( 2, 3 ) = val( ic12 )
208 t( 3, 2 ) = -val( ic21 )
209 t( 3, 3 ) = val( ic11 )*real( ic22 )
210 tnrm = max( abs( t( 1, 1 ) ),
211 $ abs( t( 1, 2 ) ), abs( t( 1, 3 ) ),
212 $ abs( t( 2, 2 ) ), abs( t( 2, 3 ) ),
213 $ abs( t( 3, 2 ) ), abs( t( 3, 3 ) ) )
214 CALL scopy( 16, t, 1, t1, 1 )
215 CALL scopy( 16, val( 1 ), 0, q, 1 )
216 CALL scopy( 4, val( 3 ), 0, q, 5 )
217 CALL slaexc( .true., 3, t, 4, q, 4, 1, 1, 2,
220 $ ninfo( info ) = ninfo( info ) + 1
221 CALL shst01( 3, 1, 3, t1, 4, t, 4, q, 4,
222 $ work, lwork, result )
223 res = result( 1 ) + result( 2 )
225 IF( t1( 1, 1 ).NE.t( 3, 3 ) )
226 $ res = res + one / eps
227 IF( t( 3, 1 ).NE.zero )
228 $ res = res + one / eps
229 IF( t( 3, 2 ).NE.zero )
230 $ res = res + one / eps
231 IF( t( 2, 1 ).NE.0 .AND.
232 $ ( t( 1, 1 ).NE.t( 2,
233 $ 2 ) .OR. sign( one, t( 1,
234 $ 2 ) ).EQ.sign( one, t( 2, 1 ) ) ) )
235 $ res = res + one / eps
238 IF( res.GT.rmax )
THEN
253 DO 150 ia22 = -1, 1, 2
257 t( 1, 1 ) = val( ia11 )
258 t( 1, 2 ) = val( ia12 )
259 t( 1, 3 ) = -two*val( ib )
260 t( 2, 1 ) = -val( ia21 )
261 t( 2, 2 ) = val( ia11 )*real( ia22 )
262 t( 2, 3 ) = val( ib )
265 t( 3, 3 ) = val( ic )*vm( icm )
266 tnrm = max( abs( t( 1, 1 ) ),
267 $ abs( t( 1, 2 ) ), abs( t( 1, 3 ) ),
268 $ abs( t( 2, 2 ) ), abs( t( 2, 3 ) ),
269 $ abs( t( 3, 2 ) ), abs( t( 3, 3 ) ) )
270 CALL scopy( 16, t, 1, t1, 1 )
271 CALL scopy( 16, val( 1 ), 0, q, 1 )
272 CALL scopy( 4, val( 3 ), 0, q, 5 )
273 CALL slaexc( .true., 3, t, 4, q, 4, 1, 2, 1,
276 $ ninfo( info ) = ninfo( info ) + 1
277 CALL shst01( 3, 1, 3, t1, 4, t, 4, q, 4,
278 $ work, lwork, result )
279 res = result( 1 ) + result( 2 )
281 IF( t1( 3, 3 ).NE.t( 1, 1 ) )
282 $ res = res + one / eps
283 IF( t( 2, 1 ).NE.zero )
284 $ res = res + one / eps
285 IF( t( 3, 1 ).NE.zero )
286 $ res = res + one / eps
287 IF( t( 3, 2 ).NE.0 .AND.
288 $ ( t( 2, 2 ).NE.t( 3,
289 $ 3 ) .OR. sign( one, t( 2,
290 $ 3 ) ).EQ.sign( one, t( 3, 2 ) ) ) )
291 $ res = res + one / eps
294 IF( res.GT.rmax )
THEN
309 DO 270 ia22 = -1, 1, 2
314 DO 220 ic22 = -1, 1, 2
317 t( 1, 1 ) = val( ia11 )*vm( iam )
318 t( 1, 2 ) = val( ia12 )*vm( iam )
319 t( 1, 3 ) = -two*val( ib )
320 t( 1, 4 ) = half*val( ib )
321 t( 2, 1 ) = -t( 1, 2 )*val( ia21 )
322 t( 2, 2 ) = val( ia11 )*
323 $ real( ia22 )*vm( iam )
324 t( 2, 3 ) = val( ib )
325 t( 2, 4 ) = three*val( ib )
328 t( 3, 3 ) = val( ic11 )*
330 t( 3, 4 ) = val( ic12 )*
334 t( 4, 3 ) = -t( 3, 4 )*val( ic21 )*
336 t( 4, 4 ) = val( ic11 )*
346 CALL scopy( 16, t, 1, t1, 1 )
347 CALL scopy( 16, val( 1 ), 0, q, 1 )
348 CALL scopy( 4, val( 3 ), 0, q, 5 )
349 CALL slaexc( .true., 4, t, 4, q, 4,
350 $ 1, 2, 2, work, info )
352 $ ninfo( info ) = ninfo( info ) + 1
353 CALL shst01( 4, 1, 4, t1, 4, t, 4,
356 res = result( 1 ) + result( 2 )
358 IF( t( 3, 1 ).NE.zero )
359 $ res = res + one / eps
360 IF( t( 4, 1 ).NE.zero )
361 $ res = res + one / eps
362 IF( t( 3, 2 ).NE.zero )
363 $ res = res + one / eps
364 IF( t( 4, 2 ).NE.zero )
365 $ res = res + one / eps
366 IF( t( 2, 1 ).NE.0 .AND.
367 $ ( t( 1, 1 ).NE.t( 2,
368 $ 2 ) .OR. sign( one, t( 1,
369 $ 2 ) ).EQ.sign( one, t( 2,
370 $ 1 ) ) ) )res = res +
372 IF( t( 4, 3 ).NE.0 .AND.
373 $ ( t( 3, 3 ).NE.t( 4,
374 $ 4 ) .OR. sign( one, t( 3,
375 $ 4 ) ).EQ.sign( one, t( 4,
376 $ 3 ) ) ) )res = res +
380 IF( res.GT.rmax )
THEN
subroutine slabad(SMALL, LARGE)
SLABAD
subroutine slaexc(WANTQ, N, T, LDT, Q, LDQ, J1, N1, N2, WORK, INFO)
SLAEXC swaps adjacent diagonal blocks of a real upper quasi-triangular matrix in Schur canonical form...
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
subroutine sget34(RMAX, LMAX, NINFO, KNT)
SGET34
subroutine shst01(N, ILO, IHI, A, LDA, H, LDH, Q, LDQ, WORK, LWORK, RESULT)
SHST01