LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ zget37()

subroutine zget37 ( double precision, dimension( 3 )  RMAX,
integer, dimension( 3 )  LMAX,
integer, dimension( 3 )  NINFO,
integer  KNT,
integer  NIN 
)

ZGET37

Purpose:
 ZGET37 tests ZTRSNA, a routine for estimating condition numbers of
 eigenvalues and/or right eigenvectors of a matrix.

 The test matrices are read from a file with logical unit number NIN.
Parameters
[out]RMAX
          RMAX is DOUBLE PRECISION array, dimension (3)
          Value of the largest test ratio.
          RMAX(1) = largest ratio comparing different calls to ZTRSNA
          RMAX(2) = largest error in reciprocal condition
                    numbers taking their conditioning into account
          RMAX(3) = largest error in reciprocal condition
                    numbers not taking their conditioning into
                    account (may be larger than RMAX(2))
[out]LMAX
          LMAX is INTEGER array, dimension (3)
          LMAX(i) is example number where largest test ratio
          RMAX(i) is achieved. Also:
          If ZGEHRD returns INFO nonzero on example i, LMAX(1)=i
          If ZHSEQR returns INFO nonzero on example i, LMAX(2)=i
          If ZTRSNA returns INFO nonzero on example i, LMAX(3)=i
[out]NINFO
          NINFO is INTEGER array, dimension (3)
          NINFO(1) = No. of times ZGEHRD returned INFO nonzero
          NINFO(2) = No. of times ZHSEQR returned INFO nonzero
          NINFO(3) = No. of times ZTRSNA returned INFO nonzero
[out]KNT
          KNT is INTEGER
          Total number of examples tested.
[in]NIN
          NIN is INTEGER
          Input logical unit number
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 89 of file zget37.f.

90*
91* -- LAPACK test routine --
92* -- LAPACK is a software package provided by Univ. of Tennessee, --
93* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
94*
95* .. Scalar Arguments ..
96 INTEGER KNT, NIN
97* ..
98* .. Array Arguments ..
99 INTEGER LMAX( 3 ), NINFO( 3 )
100 DOUBLE PRECISION RMAX( 3 )
101* ..
102*
103* =====================================================================
104*
105* .. Parameters ..
106 DOUBLE PRECISION ZERO, ONE, TWO
107 parameter( zero = 0.0d0, one = 1.0d0, two = 2.0d0 )
108 DOUBLE PRECISION EPSIN
109 parameter( epsin = 5.9605d-8 )
110 INTEGER LDT, LWORK
111 parameter( ldt = 20, lwork = 2*ldt*( 10+ldt ) )
112* ..
113* .. Local Scalars ..
114 INTEGER I, ICMP, INFO, ISCL, ISRT, J, KMIN, M, N
115 DOUBLE PRECISION BIGNUM, EPS, SMLNUM, TNRM, TOL, TOLIN, V,
116 $ VCMIN, VMAX, VMIN, VMUL
117* ..
118* .. Local Arrays ..
119 LOGICAL SELECT( LDT )
120 INTEGER LCMP( 3 )
121 DOUBLE PRECISION DUM( 1 ), RWORK( 2*LDT ), S( LDT ), SEP( LDT ),
122 $ SEPIN( LDT ), SEPTMP( LDT ), SIN( LDT ),
123 $ STMP( LDT ), VAL( 3 ), WIIN( LDT ),
124 $ WRIN( LDT ), WSRT( LDT )
125 COMPLEX*16 CDUM( 1 ), LE( LDT, LDT ), RE( LDT, LDT ),
126 $ T( LDT, LDT ), TMP( LDT, LDT ), W( LDT ),
127 $ WORK( LWORK ), WTMP( LDT )
128* ..
129* .. External Functions ..
130 DOUBLE PRECISION DLAMCH, ZLANGE
131 EXTERNAL dlamch, zlange
132* ..
133* .. External Subroutines ..
134 EXTERNAL dcopy, dlabad, dscal, zcopy, zdscal, zgehrd,
135 $ zhseqr, zlacpy, ztrevc, ztrsna
136* ..
137* .. Intrinsic Functions ..
138 INTRINSIC dble, dimag, max, sqrt
139* ..
140* .. Executable Statements ..
141*
142 eps = dlamch( 'P' )
143 smlnum = dlamch( 'S' ) / eps
144 bignum = one / smlnum
145 CALL dlabad( smlnum, bignum )
146*
147* EPSIN = 2**(-24) = precision to which input data computed
148*
149 eps = max( eps, epsin )
150 rmax( 1 ) = zero
151 rmax( 2 ) = zero
152 rmax( 3 ) = zero
153 lmax( 1 ) = 0
154 lmax( 2 ) = 0
155 lmax( 3 ) = 0
156 knt = 0
157 ninfo( 1 ) = 0
158 ninfo( 2 ) = 0
159 ninfo( 3 ) = 0
160 val( 1 ) = sqrt( smlnum )
161 val( 2 ) = one
162 val( 3 ) = sqrt( bignum )
163*
164* Read input data until N=0. Assume input eigenvalues are sorted
165* lexicographically (increasing by real part if ISRT = 0,
166* increasing by imaginary part if ISRT = 1)
167*
168 10 CONTINUE
169 READ( nin, fmt = * )n, isrt
170 IF( n.EQ.0 )
171 $ RETURN
172 DO 20 i = 1, n
173 READ( nin, fmt = * )( tmp( i, j ), j = 1, n )
174 20 CONTINUE
175 DO 30 i = 1, n
176 READ( nin, fmt = * )wrin( i ), wiin( i ), sin( i ), sepin( i )
177 30 CONTINUE
178 tnrm = zlange( 'M', n, n, tmp, ldt, rwork )
179 DO 260 iscl = 1, 3
180*
181* Scale input matrix
182*
183 knt = knt + 1
184 CALL zlacpy( 'F', n, n, tmp, ldt, t, ldt )
185 vmul = val( iscl )
186 DO 40 i = 1, n
187 CALL zdscal( n, vmul, t( 1, i ), 1 )
188 40 CONTINUE
189 IF( tnrm.EQ.zero )
190 $ vmul = one
191*
192* Compute eigenvalues and eigenvectors
193*
194 CALL zgehrd( n, 1, n, t, ldt, work( 1 ), work( n+1 ), lwork-n,
195 $ info )
196 IF( info.NE.0 ) THEN
197 lmax( 1 ) = knt
198 ninfo( 1 ) = ninfo( 1 ) + 1
199 GO TO 260
200 END IF
201 DO 60 j = 1, n - 2
202 DO 50 i = j + 2, n
203 t( i, j ) = zero
204 50 CONTINUE
205 60 CONTINUE
206*
207* Compute Schur form
208*
209 CALL zhseqr( 'S', 'N', n, 1, n, t, ldt, w, cdum, 1, work,
210 $ lwork, info )
211 IF( info.NE.0 ) THEN
212 lmax( 2 ) = knt
213 ninfo( 2 ) = ninfo( 2 ) + 1
214 GO TO 260
215 END IF
216*
217* Compute eigenvectors
218*
219 DO 70 i = 1, n
220 SELECT( i ) = .true.
221 70 CONTINUE
222 CALL ztrevc( 'B', 'A', SELECT, n, t, ldt, le, ldt, re, ldt, n,
223 $ m, work, rwork, info )
224*
225* Compute condition numbers
226*
227 CALL ztrsna( 'B', 'A', SELECT, n, t, ldt, le, ldt, re, ldt, s,
228 $ sep, n, m, work, n, rwork, info )
229 IF( info.NE.0 ) THEN
230 lmax( 3 ) = knt
231 ninfo( 3 ) = ninfo( 3 ) + 1
232 GO TO 260
233 END IF
234*
235* Sort eigenvalues and condition numbers lexicographically
236* to compare with inputs
237*
238 CALL zcopy( n, w, 1, wtmp, 1 )
239 IF( isrt.EQ.0 ) THEN
240*
241* Sort by increasing real part
242*
243 DO 80 i = 1, n
244 wsrt( i ) = dble( w( i ) )
245 80 CONTINUE
246 ELSE
247*
248* Sort by increasing imaginary part
249*
250 DO 90 i = 1, n
251 wsrt( i ) = dimag( w( i ) )
252 90 CONTINUE
253 END IF
254 CALL dcopy( n, s, 1, stmp, 1 )
255 CALL dcopy( n, sep, 1, septmp, 1 )
256 CALL dscal( n, one / vmul, septmp, 1 )
257 DO 110 i = 1, n - 1
258 kmin = i
259 vmin = wsrt( i )
260 DO 100 j = i + 1, n
261 IF( wsrt( j ).LT.vmin ) THEN
262 kmin = j
263 vmin = wsrt( j )
264 END IF
265 100 CONTINUE
266 wsrt( kmin ) = wsrt( i )
267 wsrt( i ) = vmin
268 vcmin = dble( wtmp( i ) )
269 wtmp( i ) = w( kmin )
270 wtmp( kmin ) = vcmin
271 vmin = stmp( kmin )
272 stmp( kmin ) = stmp( i )
273 stmp( i ) = vmin
274 vmin = septmp( kmin )
275 septmp( kmin ) = septmp( i )
276 septmp( i ) = vmin
277 110 CONTINUE
278*
279* Compare condition numbers for eigenvalues
280* taking their condition numbers into account
281*
282 v = max( two*dble( n )*eps*tnrm, smlnum )
283 IF( tnrm.EQ.zero )
284 $ v = one
285 DO 120 i = 1, n
286 IF( v.GT.septmp( i ) ) THEN
287 tol = one
288 ELSE
289 tol = v / septmp( i )
290 END IF
291 IF( v.GT.sepin( i ) ) THEN
292 tolin = one
293 ELSE
294 tolin = v / sepin( i )
295 END IF
296 tol = max( tol, smlnum / eps )
297 tolin = max( tolin, smlnum / eps )
298 IF( eps*( sin( i )-tolin ).GT.stmp( i )+tol ) THEN
299 vmax = one / eps
300 ELSE IF( sin( i )-tolin.GT.stmp( i )+tol ) THEN
301 vmax = ( sin( i )-tolin ) / ( stmp( i )+tol )
302 ELSE IF( sin( i )+tolin.LT.eps*( stmp( i )-tol ) ) THEN
303 vmax = one / eps
304 ELSE IF( sin( i )+tolin.LT.stmp( i )-tol ) THEN
305 vmax = ( stmp( i )-tol ) / ( sin( i )+tolin )
306 ELSE
307 vmax = one
308 END IF
309 IF( vmax.GT.rmax( 2 ) ) THEN
310 rmax( 2 ) = vmax
311 IF( ninfo( 2 ).EQ.0 )
312 $ lmax( 2 ) = knt
313 END IF
314 120 CONTINUE
315*
316* Compare condition numbers for eigenvectors
317* taking their condition numbers into account
318*
319 DO 130 i = 1, n
320 IF( v.GT.septmp( i )*stmp( i ) ) THEN
321 tol = septmp( i )
322 ELSE
323 tol = v / stmp( i )
324 END IF
325 IF( v.GT.sepin( i )*sin( i ) ) THEN
326 tolin = sepin( i )
327 ELSE
328 tolin = v / sin( i )
329 END IF
330 tol = max( tol, smlnum / eps )
331 tolin = max( tolin, smlnum / eps )
332 IF( eps*( sepin( i )-tolin ).GT.septmp( i )+tol ) THEN
333 vmax = one / eps
334 ELSE IF( sepin( i )-tolin.GT.septmp( i )+tol ) THEN
335 vmax = ( sepin( i )-tolin ) / ( septmp( i )+tol )
336 ELSE IF( sepin( i )+tolin.LT.eps*( septmp( i )-tol ) ) THEN
337 vmax = one / eps
338 ELSE IF( sepin( i )+tolin.LT.septmp( i )-tol ) THEN
339 vmax = ( septmp( i )-tol ) / ( sepin( i )+tolin )
340 ELSE
341 vmax = one
342 END IF
343 IF( vmax.GT.rmax( 2 ) ) THEN
344 rmax( 2 ) = vmax
345 IF( ninfo( 2 ).EQ.0 )
346 $ lmax( 2 ) = knt
347 END IF
348 130 CONTINUE
349*
350* Compare condition numbers for eigenvalues
351* without taking their condition numbers into account
352*
353 DO 140 i = 1, n
354 IF( sin( i ).LE.dble( 2*n )*eps .AND. stmp( i ).LE.
355 $ dble( 2*n )*eps ) THEN
356 vmax = one
357 ELSE IF( eps*sin( i ).GT.stmp( i ) ) THEN
358 vmax = one / eps
359 ELSE IF( sin( i ).GT.stmp( i ) ) THEN
360 vmax = sin( i ) / stmp( i )
361 ELSE IF( sin( i ).LT.eps*stmp( i ) ) THEN
362 vmax = one / eps
363 ELSE IF( sin( i ).LT.stmp( i ) ) THEN
364 vmax = stmp( i ) / sin( i )
365 ELSE
366 vmax = one
367 END IF
368 IF( vmax.GT.rmax( 3 ) ) THEN
369 rmax( 3 ) = vmax
370 IF( ninfo( 3 ).EQ.0 )
371 $ lmax( 3 ) = knt
372 END IF
373 140 CONTINUE
374*
375* Compare condition numbers for eigenvectors
376* without taking their condition numbers into account
377*
378 DO 150 i = 1, n
379 IF( sepin( i ).LE.v .AND. septmp( i ).LE.v ) THEN
380 vmax = one
381 ELSE IF( eps*sepin( i ).GT.septmp( i ) ) THEN
382 vmax = one / eps
383 ELSE IF( sepin( i ).GT.septmp( i ) ) THEN
384 vmax = sepin( i ) / septmp( i )
385 ELSE IF( sepin( i ).LT.eps*septmp( i ) ) THEN
386 vmax = one / eps
387 ELSE IF( sepin( i ).LT.septmp( i ) ) THEN
388 vmax = septmp( i ) / sepin( i )
389 ELSE
390 vmax = one
391 END IF
392 IF( vmax.GT.rmax( 3 ) ) THEN
393 rmax( 3 ) = vmax
394 IF( ninfo( 3 ).EQ.0 )
395 $ lmax( 3 ) = knt
396 END IF
397 150 CONTINUE
398*
399* Compute eigenvalue condition numbers only and compare
400*
401 vmax = zero
402 dum( 1 ) = -one
403 CALL dcopy( n, dum, 0, stmp, 1 )
404 CALL dcopy( n, dum, 0, septmp, 1 )
405 CALL ztrsna( 'E', 'A', SELECT, n, t, ldt, le, ldt, re, ldt,
406 $ stmp, septmp, n, m, work, n, rwork, info )
407 IF( info.NE.0 ) THEN
408 lmax( 3 ) = knt
409 ninfo( 3 ) = ninfo( 3 ) + 1
410 GO TO 260
411 END IF
412 DO 160 i = 1, n
413 IF( stmp( i ).NE.s( i ) )
414 $ vmax = one / eps
415 IF( septmp( i ).NE.dum( 1 ) )
416 $ vmax = one / eps
417 160 CONTINUE
418*
419* Compute eigenvector condition numbers only and compare
420*
421 CALL dcopy( n, dum, 0, stmp, 1 )
422 CALL dcopy( n, dum, 0, septmp, 1 )
423 CALL ztrsna( 'V', 'A', SELECT, n, t, ldt, le, ldt, re, ldt,
424 $ stmp, septmp, n, m, work, n, rwork, info )
425 IF( info.NE.0 ) THEN
426 lmax( 3 ) = knt
427 ninfo( 3 ) = ninfo( 3 ) + 1
428 GO TO 260
429 END IF
430 DO 170 i = 1, n
431 IF( stmp( i ).NE.dum( 1 ) )
432 $ vmax = one / eps
433 IF( septmp( i ).NE.sep( i ) )
434 $ vmax = one / eps
435 170 CONTINUE
436*
437* Compute all condition numbers using SELECT and compare
438*
439 DO 180 i = 1, n
440 SELECT( i ) = .true.
441 180 CONTINUE
442 CALL dcopy( n, dum, 0, stmp, 1 )
443 CALL dcopy( n, dum, 0, septmp, 1 )
444 CALL ztrsna( 'B', 'S', SELECT, n, t, ldt, le, ldt, re, ldt,
445 $ stmp, septmp, n, m, work, n, rwork, info )
446 IF( info.NE.0 ) THEN
447 lmax( 3 ) = knt
448 ninfo( 3 ) = ninfo( 3 ) + 1
449 GO TO 260
450 END IF
451 DO 190 i = 1, n
452 IF( septmp( i ).NE.sep( i ) )
453 $ vmax = one / eps
454 IF( stmp( i ).NE.s( i ) )
455 $ vmax = one / eps
456 190 CONTINUE
457*
458* Compute eigenvalue condition numbers using SELECT and compare
459*
460 CALL dcopy( n, dum, 0, stmp, 1 )
461 CALL dcopy( n, dum, 0, septmp, 1 )
462 CALL ztrsna( 'E', 'S', SELECT, n, t, ldt, le, ldt, re, ldt,
463 $ stmp, septmp, n, m, work, n, rwork, info )
464 IF( info.NE.0 ) THEN
465 lmax( 3 ) = knt
466 ninfo( 3 ) = ninfo( 3 ) + 1
467 GO TO 260
468 END IF
469 DO 200 i = 1, n
470 IF( stmp( i ).NE.s( i ) )
471 $ vmax = one / eps
472 IF( septmp( i ).NE.dum( 1 ) )
473 $ vmax = one / eps
474 200 CONTINUE
475*
476* Compute eigenvector condition numbers using SELECT and compare
477*
478 CALL dcopy( n, dum, 0, stmp, 1 )
479 CALL dcopy( n, dum, 0, septmp, 1 )
480 CALL ztrsna( 'V', 'S', SELECT, n, t, ldt, le, ldt, re, ldt,
481 $ stmp, septmp, n, m, work, n, rwork, info )
482 IF( info.NE.0 ) THEN
483 lmax( 3 ) = knt
484 ninfo( 3 ) = ninfo( 3 ) + 1
485 GO TO 260
486 END IF
487 DO 210 i = 1, n
488 IF( stmp( i ).NE.dum( 1 ) )
489 $ vmax = one / eps
490 IF( septmp( i ).NE.sep( i ) )
491 $ vmax = one / eps
492 210 CONTINUE
493 IF( vmax.GT.rmax( 1 ) ) THEN
494 rmax( 1 ) = vmax
495 IF( ninfo( 1 ).EQ.0 )
496 $ lmax( 1 ) = knt
497 END IF
498*
499* Select second and next to last eigenvalues
500*
501 DO 220 i = 1, n
502 SELECT( i ) = .false.
503 220 CONTINUE
504 icmp = 0
505 IF( n.GT.1 ) THEN
506 icmp = 1
507 lcmp( 1 ) = 2
508 SELECT( 2 ) = .true.
509 CALL zcopy( n, re( 1, 2 ), 1, re( 1, 1 ), 1 )
510 CALL zcopy( n, le( 1, 2 ), 1, le( 1, 1 ), 1 )
511 END IF
512 IF( n.GT.3 ) THEN
513 icmp = 2
514 lcmp( 2 ) = n - 1
515 SELECT( n-1 ) = .true.
516 CALL zcopy( n, re( 1, n-1 ), 1, re( 1, 2 ), 1 )
517 CALL zcopy( n, le( 1, n-1 ), 1, le( 1, 2 ), 1 )
518 END IF
519*
520* Compute all selected condition numbers
521*
522 CALL dcopy( icmp, dum, 0, stmp, 1 )
523 CALL dcopy( icmp, dum, 0, septmp, 1 )
524 CALL ztrsna( 'B', 'S', SELECT, n, t, ldt, le, ldt, re, ldt,
525 $ stmp, septmp, n, m, work, n, rwork, info )
526 IF( info.NE.0 ) THEN
527 lmax( 3 ) = knt
528 ninfo( 3 ) = ninfo( 3 ) + 1
529 GO TO 260
530 END IF
531 DO 230 i = 1, icmp
532 j = lcmp( i )
533 IF( septmp( i ).NE.sep( j ) )
534 $ vmax = one / eps
535 IF( stmp( i ).NE.s( j ) )
536 $ vmax = one / eps
537 230 CONTINUE
538*
539* Compute selected eigenvalue condition numbers
540*
541 CALL dcopy( icmp, dum, 0, stmp, 1 )
542 CALL dcopy( icmp, dum, 0, septmp, 1 )
543 CALL ztrsna( 'E', 'S', SELECT, n, t, ldt, le, ldt, re, ldt,
544 $ stmp, septmp, n, m, work, n, rwork, info )
545 IF( info.NE.0 ) THEN
546 lmax( 3 ) = knt
547 ninfo( 3 ) = ninfo( 3 ) + 1
548 GO TO 260
549 END IF
550 DO 240 i = 1, icmp
551 j = lcmp( i )
552 IF( stmp( i ).NE.s( j ) )
553 $ vmax = one / eps
554 IF( septmp( i ).NE.dum( 1 ) )
555 $ vmax = one / eps
556 240 CONTINUE
557*
558* Compute selected eigenvector condition numbers
559*
560 CALL dcopy( icmp, dum, 0, stmp, 1 )
561 CALL dcopy( icmp, dum, 0, septmp, 1 )
562 CALL ztrsna( 'V', 'S', SELECT, n, t, ldt, le, ldt, re, ldt,
563 $ stmp, septmp, n, m, work, n, rwork, info )
564 IF( info.NE.0 ) THEN
565 lmax( 3 ) = knt
566 ninfo( 3 ) = ninfo( 3 ) + 1
567 GO TO 260
568 END IF
569 DO 250 i = 1, icmp
570 j = lcmp( i )
571 IF( stmp( i ).NE.dum( 1 ) )
572 $ vmax = one / eps
573 IF( septmp( i ).NE.sep( j ) )
574 $ vmax = one / eps
575 250 CONTINUE
576 IF( vmax.GT.rmax( 1 ) ) THEN
577 rmax( 1 ) = vmax
578 IF( ninfo( 1 ).EQ.0 )
579 $ lmax( 1 ) = knt
580 END IF
581 260 CONTINUE
582 GO TO 10
583*
584* End of ZGET37
585*
dlamch
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
dlabad
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
zdscal
subroutine zdscal(N, DA, ZX, INCX)
ZDSCAL
Definition: zdscal.f:78
zcopy
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:81
zlange
double precision function zlange(NORM, M, N, A, LDA, WORK)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: zlange.f:115
zgehrd
subroutine zgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
ZGEHRD
Definition: zgehrd.f:167
zlacpy
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:103
zhseqr
subroutine zhseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, W, Z, LDZ, WORK, LWORK, INFO)
ZHSEQR
Definition: zhseqr.f:299
ztrevc
subroutine ztrevc(SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
ZTREVC
Definition: ztrevc.f:218
ztrsna
subroutine ztrsna(JOB, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, S, SEP, MM, M, WORK, LDWORK, RWORK, INFO)
ZTRSNA
Definition: ztrsna.f:249
dcopy
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:82
dscal
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:79
Here is the call graph for this function:
Here is the caller graph for this function: