You can not select more than 25 topics Topics must start with a chinese character,a letter or number, can include dashes ('-') and can be up to 35 characters long.

zchkhb.f 25 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731
  1. *> \brief \b ZCHKHB
  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 ZCHKHB( NSIZES, NN, NWDTHS, KK, NTYPES, DOTYPE, ISEED,
  12. * THRESH, NOUNIT, A, LDA, SD, SE, U, LDU, WORK,
  13. * LWORK, RWORK, RESULT, INFO )
  14. *
  15. * .. Scalar Arguments ..
  16. * INTEGER INFO, LDA, LDU, LWORK, NOUNIT, NSIZES, NTYPES,
  17. * $ NWDTHS
  18. * DOUBLE PRECISION THRESH
  19. * ..
  20. * .. Array Arguments ..
  21. * LOGICAL DOTYPE( * )
  22. * INTEGER ISEED( 4 ), KK( * ), NN( * )
  23. * DOUBLE PRECISION RESULT( * ), RWORK( * ), SD( * ), SE( * )
  24. * COMPLEX*16 A( LDA, * ), U( LDU, * ), WORK( * )
  25. * ..
  26. *
  27. *
  28. *> \par Purpose:
  29. * =============
  30. *>
  31. *> \verbatim
  32. *>
  33. *> ZCHKHB tests the reduction of a Hermitian band matrix to tridiagonal
  34. *> from, used with the Hermitian eigenvalue problem.
  35. *>
  36. *> ZHBTRD factors a Hermitian band matrix A as U S U* , where * means
  37. *> conjugate transpose, S is symmetric tridiagonal, and U is unitary.
  38. *> ZHBTRD can use either just the lower or just the upper triangle
  39. *> of A; ZCHKHB checks both cases.
  40. *>
  41. *> When ZCHKHB is called, a number of matrix "sizes" ("n's"), a number
  42. *> of bandwidths ("k's"), and a number of matrix "types" are
  43. *> specified. For each size ("n"), each bandwidth ("k") less than or
  44. *> equal to "n", and each type of matrix, one matrix will be generated
  45. *> and used to test the hermitian banded reduction routine. For each
  46. *> matrix, a number of tests will be performed:
  47. *>
  48. *> (1) | A - V S V* | / ( |A| n ulp ) computed by ZHBTRD with
  49. *> UPLO='U'
  50. *>
  51. *> (2) | I - UU* | / ( n ulp )
  52. *>
  53. *> (3) | A - V S V* | / ( |A| n ulp ) computed by ZHBTRD with
  54. *> UPLO='L'
  55. *>
  56. *> (4) | I - UU* | / ( n ulp )
  57. *>
  58. *> The "sizes" are specified by an array NN(1:NSIZES); the value of
  59. *> each element NN(j) specifies one size.
  60. *> The "types" are specified by a logical array DOTYPE( 1:NTYPES );
  61. *> if DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
  62. *> Currently, the list of possible types is:
  63. *>
  64. *> (1) The zero matrix.
  65. *> (2) The identity matrix.
  66. *>
  67. *> (3) A diagonal matrix with evenly spaced entries
  68. *> 1, ..., ULP and random signs.
  69. *> (ULP = (first number larger than 1) - 1 )
  70. *> (4) A diagonal matrix with geometrically spaced entries
  71. *> 1, ..., ULP and random signs.
  72. *> (5) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
  73. *> and random signs.
  74. *>
  75. *> (6) Same as (4), but multiplied by SQRT( overflow threshold )
  76. *> (7) Same as (4), but multiplied by SQRT( underflow threshold )
  77. *>
  78. *> (8) A matrix of the form U* D U, where U is unitary and
  79. *> D has evenly spaced entries 1, ..., ULP with random signs
  80. *> on the diagonal.
  81. *>
  82. *> (9) A matrix of the form U* D U, where U is unitary and
  83. *> D has geometrically spaced entries 1, ..., ULP with random
  84. *> signs on the diagonal.
  85. *>
  86. *> (10) A matrix of the form U* D U, where U is unitary and
  87. *> D has "clustered" entries 1, ULP,..., ULP with random
  88. *> signs on the diagonal.
  89. *>
  90. *> (11) Same as (8), but multiplied by SQRT( overflow threshold )
  91. *> (12) Same as (8), but multiplied by SQRT( underflow threshold )
  92. *>
  93. *> (13) Hermitian matrix with random entries chosen from (-1,1).
  94. *> (14) Same as (13), but multiplied by SQRT( overflow threshold )
  95. *> (15) Same as (13), but multiplied by SQRT( underflow threshold )
  96. *> \endverbatim
  97. *
  98. * Arguments:
  99. * ==========
  100. *
  101. *> \param[in] NSIZES
  102. *> \verbatim
  103. *> NSIZES is INTEGER
  104. *> The number of sizes of matrices to use. If it is zero,
  105. *> ZCHKHB does nothing. It must be at least zero.
  106. *> \endverbatim
  107. *>
  108. *> \param[in] NN
  109. *> \verbatim
  110. *> NN is INTEGER array, dimension (NSIZES)
  111. *> An array containing the sizes to be used for the matrices.
  112. *> Zero values will be skipped. The values must be at least
  113. *> zero.
  114. *> \endverbatim
  115. *>
  116. *> \param[in] NWDTHS
  117. *> \verbatim
  118. *> NWDTHS is INTEGER
  119. *> The number of bandwidths to use. If it is zero,
  120. *> ZCHKHB does nothing. It must be at least zero.
  121. *> \endverbatim
  122. *>
  123. *> \param[in] KK
  124. *> \verbatim
  125. *> KK is INTEGER array, dimension (NWDTHS)
  126. *> An array containing the bandwidths to be used for the band
  127. *> matrices. The values must be at least zero.
  128. *> \endverbatim
  129. *>
  130. *> \param[in] NTYPES
  131. *> \verbatim
  132. *> NTYPES is INTEGER
  133. *> The number of elements in DOTYPE. If it is zero, ZCHKHB
  134. *> does nothing. It must be at least zero. If it is MAXTYP+1
  135. *> and NSIZES is 1, then an additional type, MAXTYP+1 is
  136. *> defined, which is to use whatever matrix is in A. This
  137. *> is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
  138. *> DOTYPE(MAXTYP+1) is .TRUE. .
  139. *> \endverbatim
  140. *>
  141. *> \param[in] DOTYPE
  142. *> \verbatim
  143. *> DOTYPE is LOGICAL array, dimension (NTYPES)
  144. *> If DOTYPE(j) is .TRUE., then for each size in NN a
  145. *> matrix of that size and of type j will be generated.
  146. *> If NTYPES is smaller than the maximum number of types
  147. *> defined (PARAMETER MAXTYP), then types NTYPES+1 through
  148. *> MAXTYP will not be generated. If NTYPES is larger
  149. *> than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
  150. *> will be ignored.
  151. *> \endverbatim
  152. *>
  153. *> \param[in,out] ISEED
  154. *> \verbatim
  155. *> ISEED is INTEGER array, dimension (4)
  156. *> On entry ISEED specifies the seed of the random number
  157. *> generator. The array elements should be between 0 and 4095;
  158. *> if not they will be reduced mod 4096. Also, ISEED(4) must
  159. *> be odd. The random number generator uses a linear
  160. *> congruential sequence limited to small integers, and so
  161. *> should produce machine independent random numbers. The
  162. *> values of ISEED are changed on exit, and can be used in the
  163. *> next call to ZCHKHB to continue the same random number
  164. *> sequence.
  165. *> \endverbatim
  166. *>
  167. *> \param[in] THRESH
  168. *> \verbatim
  169. *> THRESH is DOUBLE PRECISION
  170. *> A test will count as "failed" if the "error", computed as
  171. *> described above, exceeds THRESH. Note that the error
  172. *> is scaled to be O(1), so THRESH should be a reasonably
  173. *> small multiple of 1, e.g., 10 or 100. In particular,
  174. *> it should not depend on the precision (single vs. double)
  175. *> or the size of the matrix. It must be at least zero.
  176. *> \endverbatim
  177. *>
  178. *> \param[in] NOUNIT
  179. *> \verbatim
  180. *> NOUNIT is INTEGER
  181. *> The FORTRAN unit number for printing out error messages
  182. *> (e.g., if a routine returns IINFO not equal to 0.)
  183. *> \endverbatim
  184. *>
  185. *> \param[in,out] A
  186. *> \verbatim
  187. *> A is COMPLEX*16 array, dimension
  188. *> (LDA, max(NN))
  189. *> Used to hold the matrix whose eigenvalues are to be
  190. *> computed.
  191. *> \endverbatim
  192. *>
  193. *> \param[in] LDA
  194. *> \verbatim
  195. *> LDA is INTEGER
  196. *> The leading dimension of A. It must be at least 2 (not 1!)
  197. *> and at least max( KK )+1.
  198. *> \endverbatim
  199. *>
  200. *> \param[out] SD
  201. *> \verbatim
  202. *> SD is DOUBLE PRECISION array, dimension (max(NN))
  203. *> Used to hold the diagonal of the tridiagonal matrix computed
  204. *> by ZHBTRD.
  205. *> \endverbatim
  206. *>
  207. *> \param[out] SE
  208. *> \verbatim
  209. *> SE is DOUBLE PRECISION array, dimension (max(NN))
  210. *> Used to hold the off-diagonal of the tridiagonal matrix
  211. *> computed by ZHBTRD.
  212. *> \endverbatim
  213. *>
  214. *> \param[out] U
  215. *> \verbatim
  216. *> U is COMPLEX*16 array, dimension (LDU, max(NN))
  217. *> Used to hold the unitary matrix computed by ZHBTRD.
  218. *> \endverbatim
  219. *>
  220. *> \param[in] LDU
  221. *> \verbatim
  222. *> LDU is INTEGER
  223. *> The leading dimension of U. It must be at least 1
  224. *> and at least max( NN ).
  225. *> \endverbatim
  226. *>
  227. *> \param[out] WORK
  228. *> \verbatim
  229. *> WORK is COMPLEX*16 array, dimension (LWORK)
  230. *> \endverbatim
  231. *>
  232. *> \param[in] LWORK
  233. *> \verbatim
  234. *> LWORK is INTEGER
  235. *> The number of entries in WORK. This must be at least
  236. *> max( LDA+1, max(NN)+1 )*max(NN).
  237. *> \endverbatim
  238. *>
  239. *> \param[out] RWORK
  240. *> \verbatim
  241. *> RWORK is DOUBLE PRECISION array
  242. *> \endverbatim
  243. *>
  244. *> \param[out] RESULT
  245. *> \verbatim
  246. *> RESULT is DOUBLE PRECISION array, dimension (4)
  247. *> The values computed by the tests described above.
  248. *> The values are currently limited to 1/ulp, to avoid
  249. *> overflow.
  250. *> \endverbatim
  251. *>
  252. *> \param[out] INFO
  253. *> \verbatim
  254. *> INFO is INTEGER
  255. *> If 0, then everything ran OK.
  256. *>
  257. *>-----------------------------------------------------------------------
  258. *>
  259. *> Some Local Variables and Parameters:
  260. *> ---- ----- --------- --- ----------
  261. *> ZERO, ONE Real 0 and 1.
  262. *> MAXTYP The number of types defined.
  263. *> NTEST The number of tests performed, or which can
  264. *> be performed so far, for the current matrix.
  265. *> NTESTT The total number of tests performed so far.
  266. *> NMAX Largest value in NN.
  267. *> NMATS The number of matrices generated so far.
  268. *> NERRS The number of tests which have exceeded THRESH
  269. *> so far.
  270. *> COND, IMODE Values to be passed to the matrix generators.
  271. *> ANORM Norm of A; passed to matrix generators.
  272. *>
  273. *> OVFL, UNFL Overflow and underflow thresholds.
  274. *> ULP, ULPINV Finest relative precision and its inverse.
  275. *> RTOVFL, RTUNFL Square roots of the previous 2 values.
  276. *> The following four arrays decode JTYPE:
  277. *> KTYPE(j) The general type (1-10) for type "j".
  278. *> KMODE(j) The MODE value to be passed to the matrix
  279. *> generator for type "j".
  280. *> KMAGN(j) The order of magnitude ( O(1),
  281. *> O(overflow^(1/2) ), O(underflow^(1/2) )
  282. *> \endverbatim
  283. *
  284. * Authors:
  285. * ========
  286. *
  287. *> \author Univ. of Tennessee
  288. *> \author Univ. of California Berkeley
  289. *> \author Univ. of Colorado Denver
  290. *> \author NAG Ltd.
  291. *
  292. *> \ingroup complex16_eig
  293. *
  294. * =====================================================================
  295. SUBROUTINE ZCHKHB( NSIZES, NN, NWDTHS, KK, NTYPES, DOTYPE, ISEED,
  296. $ THRESH, NOUNIT, A, LDA, SD, SE, U, LDU, WORK,
  297. $ LWORK, RWORK, RESULT, INFO )
  298. *
  299. * -- LAPACK test routine --
  300. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  301. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  302. *
  303. * .. Scalar Arguments ..
  304. INTEGER INFO, LDA, LDU, LWORK, NOUNIT, NSIZES, NTYPES,
  305. $ NWDTHS
  306. DOUBLE PRECISION THRESH
  307. * ..
  308. * .. Array Arguments ..
  309. LOGICAL DOTYPE( * )
  310. INTEGER ISEED( 4 ), KK( * ), NN( * )
  311. DOUBLE PRECISION RESULT( * ), RWORK( * ), SD( * ), SE( * )
  312. COMPLEX*16 A( LDA, * ), U( LDU, * ), WORK( * )
  313. * ..
  314. *
  315. * =====================================================================
  316. *
  317. * .. Parameters ..
  318. COMPLEX*16 CZERO, CONE
  319. PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ),
  320. $ CONE = ( 1.0D+0, 0.0D+0 ) )
  321. DOUBLE PRECISION ZERO, ONE, TWO, TEN
  322. PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0, TWO = 2.0D+0,
  323. $ TEN = 10.0D+0 )
  324. DOUBLE PRECISION HALF
  325. PARAMETER ( HALF = ONE / TWO )
  326. INTEGER MAXTYP
  327. PARAMETER ( MAXTYP = 15 )
  328. * ..
  329. * .. Local Scalars ..
  330. LOGICAL BADNN, BADNNB
  331. INTEGER I, IINFO, IMODE, ITYPE, J, JC, JCOL, JR, JSIZE,
  332. $ JTYPE, JWIDTH, K, KMAX, MTYPES, N, NERRS,
  333. $ NMATS, NMAX, NTEST, NTESTT
  334. DOUBLE PRECISION ANINV, ANORM, COND, OVFL, RTOVFL, RTUNFL,
  335. $ TEMP1, ULP, ULPINV, UNFL
  336. * ..
  337. * .. Local Arrays ..
  338. INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KMAGN( MAXTYP ),
  339. $ KMODE( MAXTYP ), KTYPE( MAXTYP )
  340. * ..
  341. * .. External Functions ..
  342. DOUBLE PRECISION DLAMCH
  343. EXTERNAL DLAMCH
  344. * ..
  345. * .. External Subroutines ..
  346. EXTERNAL DLASUM, XERBLA, ZHBT21, ZHBTRD, ZLACPY, ZLASET,
  347. $ ZLATMR, ZLATMS
  348. * ..
  349. * .. Intrinsic Functions ..
  350. INTRINSIC ABS, DBLE, DCONJG, MAX, MIN, SQRT
  351. * ..
  352. * .. Data statements ..
  353. DATA KTYPE / 1, 2, 5*4, 5*5, 3*8 /
  354. DATA KMAGN / 2*1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1,
  355. $ 2, 3 /
  356. DATA KMODE / 2*0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0,
  357. $ 0, 0 /
  358. * ..
  359. * .. Executable Statements ..
  360. *
  361. * Check for errors
  362. *
  363. NTESTT = 0
  364. INFO = 0
  365. *
  366. * Important constants
  367. *
  368. BADNN = .FALSE.
  369. NMAX = 1
  370. DO 10 J = 1, NSIZES
  371. NMAX = MAX( NMAX, NN( J ) )
  372. IF( NN( J ).LT.0 )
  373. $ BADNN = .TRUE.
  374. 10 CONTINUE
  375. *
  376. BADNNB = .FALSE.
  377. KMAX = 0
  378. DO 20 J = 1, NSIZES
  379. KMAX = MAX( KMAX, KK( J ) )
  380. IF( KK( J ).LT.0 )
  381. $ BADNNB = .TRUE.
  382. 20 CONTINUE
  383. KMAX = MIN( NMAX-1, KMAX )
  384. *
  385. * Check for errors
  386. *
  387. IF( NSIZES.LT.0 ) THEN
  388. INFO = -1
  389. ELSE IF( BADNN ) THEN
  390. INFO = -2
  391. ELSE IF( NWDTHS.LT.0 ) THEN
  392. INFO = -3
  393. ELSE IF( BADNNB ) THEN
  394. INFO = -4
  395. ELSE IF( NTYPES.LT.0 ) THEN
  396. INFO = -5
  397. ELSE IF( LDA.LT.KMAX+1 ) THEN
  398. INFO = -11
  399. ELSE IF( LDU.LT.NMAX ) THEN
  400. INFO = -15
  401. ELSE IF( ( MAX( LDA, NMAX )+1 )*NMAX.GT.LWORK ) THEN
  402. INFO = -17
  403. END IF
  404. *
  405. IF( INFO.NE.0 ) THEN
  406. CALL XERBLA( 'ZCHKHB', -INFO )
  407. RETURN
  408. END IF
  409. *
  410. * Quick return if possible
  411. *
  412. IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 .OR. NWDTHS.EQ.0 )
  413. $ RETURN
  414. *
  415. * More Important constants
  416. *
  417. UNFL = DLAMCH( 'Safe minimum' )
  418. OVFL = ONE / UNFL
  419. ULP = DLAMCH( 'Epsilon' )*DLAMCH( 'Base' )
  420. ULPINV = ONE / ULP
  421. RTUNFL = SQRT( UNFL )
  422. RTOVFL = SQRT( OVFL )
  423. *
  424. * Loop over sizes, types
  425. *
  426. NERRS = 0
  427. NMATS = 0
  428. *
  429. DO 190 JSIZE = 1, NSIZES
  430. N = NN( JSIZE )
  431. ANINV = ONE / DBLE( MAX( 1, N ) )
  432. *
  433. DO 180 JWIDTH = 1, NWDTHS
  434. K = KK( JWIDTH )
  435. IF( K.GT.N )
  436. $ GO TO 180
  437. K = MAX( 0, MIN( N-1, K ) )
  438. *
  439. IF( NSIZES.NE.1 ) THEN
  440. MTYPES = MIN( MAXTYP, NTYPES )
  441. ELSE
  442. MTYPES = MIN( MAXTYP+1, NTYPES )
  443. END IF
  444. *
  445. DO 170 JTYPE = 1, MTYPES
  446. IF( .NOT.DOTYPE( JTYPE ) )
  447. $ GO TO 170
  448. NMATS = NMATS + 1
  449. NTEST = 0
  450. *
  451. DO 30 J = 1, 4
  452. IOLDSD( J ) = ISEED( J )
  453. 30 CONTINUE
  454. *
  455. * Compute "A".
  456. * Store as "Upper"; later, we will copy to other format.
  457. *
  458. * Control parameters:
  459. *
  460. * KMAGN KMODE KTYPE
  461. * =1 O(1) clustered 1 zero
  462. * =2 large clustered 2 identity
  463. * =3 small exponential (none)
  464. * =4 arithmetic diagonal, (w/ eigenvalues)
  465. * =5 random log hermitian, w/ eigenvalues
  466. * =6 random (none)
  467. * =7 random diagonal
  468. * =8 random hermitian
  469. * =9 positive definite
  470. * =10 diagonally dominant tridiagonal
  471. *
  472. IF( MTYPES.GT.MAXTYP )
  473. $ GO TO 100
  474. *
  475. ITYPE = KTYPE( JTYPE )
  476. IMODE = KMODE( JTYPE )
  477. *
  478. * Compute norm
  479. *
  480. GO TO ( 40, 50, 60 )KMAGN( JTYPE )
  481. *
  482. 40 CONTINUE
  483. ANORM = ONE
  484. GO TO 70
  485. *
  486. 50 CONTINUE
  487. ANORM = ( RTOVFL*ULP )*ANINV
  488. GO TO 70
  489. *
  490. 60 CONTINUE
  491. ANORM = RTUNFL*N*ULPINV
  492. GO TO 70
  493. *
  494. 70 CONTINUE
  495. *
  496. CALL ZLASET( 'Full', LDA, N, CZERO, CZERO, A, LDA )
  497. IINFO = 0
  498. IF( JTYPE.LE.15 ) THEN
  499. COND = ULPINV
  500. ELSE
  501. COND = ULPINV*ANINV / TEN
  502. END IF
  503. *
  504. * Special Matrices -- Identity & Jordan block
  505. *
  506. * Zero
  507. *
  508. IF( ITYPE.EQ.1 ) THEN
  509. IINFO = 0
  510. *
  511. ELSE IF( ITYPE.EQ.2 ) THEN
  512. *
  513. * Identity
  514. *
  515. DO 80 JCOL = 1, N
  516. A( K+1, JCOL ) = ANORM
  517. 80 CONTINUE
  518. *
  519. ELSE IF( ITYPE.EQ.4 ) THEN
  520. *
  521. * Diagonal Matrix, [Eigen]values Specified
  522. *
  523. CALL ZLATMS( N, N, 'S', ISEED, 'H', RWORK, IMODE,
  524. $ COND, ANORM, 0, 0, 'Q', A( K+1, 1 ), LDA,
  525. $ WORK, IINFO )
  526. *
  527. ELSE IF( ITYPE.EQ.5 ) THEN
  528. *
  529. * Hermitian, eigenvalues specified
  530. *
  531. CALL ZLATMS( N, N, 'S', ISEED, 'H', RWORK, IMODE,
  532. $ COND, ANORM, K, K, 'Q', A, LDA, WORK,
  533. $ IINFO )
  534. *
  535. ELSE IF( ITYPE.EQ.7 ) THEN
  536. *
  537. * Diagonal, random eigenvalues
  538. *
  539. CALL ZLATMR( N, N, 'S', ISEED, 'H', WORK, 6, ONE,
  540. $ CONE, 'T', 'N', WORK( N+1 ), 1, ONE,
  541. $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0,
  542. $ ZERO, ANORM, 'Q', A( K+1, 1 ), LDA,
  543. $ IDUMMA, IINFO )
  544. *
  545. ELSE IF( ITYPE.EQ.8 ) THEN
  546. *
  547. * Hermitian, random eigenvalues
  548. *
  549. CALL ZLATMR( N, N, 'S', ISEED, 'H', WORK, 6, ONE,
  550. $ CONE, 'T', 'N', WORK( N+1 ), 1, ONE,
  551. $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, K, K,
  552. $ ZERO, ANORM, 'Q', A, LDA, IDUMMA, IINFO )
  553. *
  554. ELSE IF( ITYPE.EQ.9 ) THEN
  555. *
  556. * Positive definite, eigenvalues specified.
  557. *
  558. CALL ZLATMS( N, N, 'S', ISEED, 'P', RWORK, IMODE,
  559. $ COND, ANORM, K, K, 'Q', A, LDA,
  560. $ WORK( N+1 ), IINFO )
  561. *
  562. ELSE IF( ITYPE.EQ.10 ) THEN
  563. *
  564. * Positive definite tridiagonal, eigenvalues specified.
  565. *
  566. IF( N.GT.1 )
  567. $ K = MAX( 1, K )
  568. CALL ZLATMS( N, N, 'S', ISEED, 'P', RWORK, IMODE,
  569. $ COND, ANORM, 1, 1, 'Q', A( K, 1 ), LDA,
  570. $ WORK, IINFO )
  571. DO 90 I = 2, N
  572. TEMP1 = ABS( A( K, I ) ) /
  573. $ SQRT( ABS( A( K+1, I-1 )*A( K+1, I ) ) )
  574. IF( TEMP1.GT.HALF ) THEN
  575. A( K, I ) = HALF*SQRT( ABS( A( K+1,
  576. $ I-1 )*A( K+1, I ) ) )
  577. END IF
  578. 90 CONTINUE
  579. *
  580. ELSE
  581. *
  582. IINFO = 1
  583. END IF
  584. *
  585. IF( IINFO.NE.0 ) THEN
  586. WRITE( NOUNIT, FMT = 9999 )'Generator', IINFO, N,
  587. $ JTYPE, IOLDSD
  588. INFO = ABS( IINFO )
  589. RETURN
  590. END IF
  591. *
  592. 100 CONTINUE
  593. *
  594. * Call ZHBTRD to compute S and U from upper triangle.
  595. *
  596. CALL ZLACPY( ' ', K+1, N, A, LDA, WORK, LDA )
  597. *
  598. NTEST = 1
  599. CALL ZHBTRD( 'V', 'U', N, K, WORK, LDA, SD, SE, U, LDU,
  600. $ WORK( LDA*N+1 ), IINFO )
  601. *
  602. IF( IINFO.NE.0 ) THEN
  603. WRITE( NOUNIT, FMT = 9999 )'ZHBTRD(U)', IINFO, N,
  604. $ JTYPE, IOLDSD
  605. INFO = ABS( IINFO )
  606. IF( IINFO.LT.0 ) THEN
  607. RETURN
  608. ELSE
  609. RESULT( 1 ) = ULPINV
  610. GO TO 150
  611. END IF
  612. END IF
  613. *
  614. * Do tests 1 and 2
  615. *
  616. CALL ZHBT21( 'Upper', N, K, 1, A, LDA, SD, SE, U, LDU,
  617. $ WORK, RWORK, RESULT( 1 ) )
  618. *
  619. * Convert A from Upper-Triangle-Only storage to
  620. * Lower-Triangle-Only storage.
  621. *
  622. DO 120 JC = 1, N
  623. DO 110 JR = 0, MIN( K, N-JC )
  624. A( JR+1, JC ) = DCONJG( A( K+1-JR, JC+JR ) )
  625. 110 CONTINUE
  626. 120 CONTINUE
  627. DO 140 JC = N + 1 - K, N
  628. DO 130 JR = MIN( K, N-JC ) + 1, K
  629. A( JR+1, JC ) = ZERO
  630. 130 CONTINUE
  631. 140 CONTINUE
  632. *
  633. * Call ZHBTRD to compute S and U from lower triangle
  634. *
  635. CALL ZLACPY( ' ', K+1, N, A, LDA, WORK, LDA )
  636. *
  637. NTEST = 3
  638. CALL ZHBTRD( 'V', 'L', N, K, WORK, LDA, SD, SE, U, LDU,
  639. $ WORK( LDA*N+1 ), IINFO )
  640. *
  641. IF( IINFO.NE.0 ) THEN
  642. WRITE( NOUNIT, FMT = 9999 )'ZHBTRD(L)', IINFO, N,
  643. $ JTYPE, IOLDSD
  644. INFO = ABS( IINFO )
  645. IF( IINFO.LT.0 ) THEN
  646. RETURN
  647. ELSE
  648. RESULT( 3 ) = ULPINV
  649. GO TO 150
  650. END IF
  651. END IF
  652. NTEST = 4
  653. *
  654. * Do tests 3 and 4
  655. *
  656. CALL ZHBT21( 'Lower', N, K, 1, A, LDA, SD, SE, U, LDU,
  657. $ WORK, RWORK, RESULT( 3 ) )
  658. *
  659. * End of Loop -- Check for RESULT(j) > THRESH
  660. *
  661. 150 CONTINUE
  662. NTESTT = NTESTT + NTEST
  663. *
  664. * Print out tests which fail.
  665. *
  666. DO 160 JR = 1, NTEST
  667. IF( RESULT( JR ).GE.THRESH ) THEN
  668. *
  669. * If this is the first test to fail,
  670. * print a header to the data file.
  671. *
  672. IF( NERRS.EQ.0 ) THEN
  673. WRITE( NOUNIT, FMT = 9998 )'ZHB'
  674. WRITE( NOUNIT, FMT = 9997 )
  675. WRITE( NOUNIT, FMT = 9996 )
  676. WRITE( NOUNIT, FMT = 9995 )'Hermitian'
  677. WRITE( NOUNIT, FMT = 9994 )'unitary', '*',
  678. $ 'conjugate transpose', ( '*', J = 1, 4 )
  679. END IF
  680. NERRS = NERRS + 1
  681. WRITE( NOUNIT, FMT = 9993 )N, K, IOLDSD, JTYPE,
  682. $ JR, RESULT( JR )
  683. END IF
  684. 160 CONTINUE
  685. *
  686. 170 CONTINUE
  687. 180 CONTINUE
  688. 190 CONTINUE
  689. *
  690. * Summary
  691. *
  692. CALL DLASUM( 'ZHB', NOUNIT, NERRS, NTESTT )
  693. RETURN
  694. *
  695. 9999 FORMAT( ' ZCHKHB: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
  696. $ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
  697. 9998 FORMAT( / 1X, A3,
  698. $ ' -- Complex Hermitian Banded Tridiagonal Reduction Routines'
  699. $ )
  700. 9997 FORMAT( ' Matrix types (see DCHK23 for details): ' )
  701. *
  702. 9996 FORMAT( / ' Special Matrices:',
  703. $ / ' 1=Zero matrix. ',
  704. $ ' 5=Diagonal: clustered entries.',
  705. $ / ' 2=Identity matrix. ',
  706. $ ' 6=Diagonal: large, evenly spaced.',
  707. $ / ' 3=Diagonal: evenly spaced entries. ',
  708. $ ' 7=Diagonal: small, evenly spaced.',
  709. $ / ' 4=Diagonal: geometr. spaced entries.' )
  710. 9995 FORMAT( ' Dense ', A, ' Banded Matrices:',
  711. $ / ' 8=Evenly spaced eigenvals. ',
  712. $ ' 12=Small, evenly spaced eigenvals.',
  713. $ / ' 9=Geometrically spaced eigenvals. ',
  714. $ ' 13=Matrix with random O(1) entries.',
  715. $ / ' 10=Clustered eigenvalues. ',
  716. $ ' 14=Matrix with large random entries.',
  717. $ / ' 11=Large, evenly spaced eigenvals. ',
  718. $ ' 15=Matrix with small random entries.' )
  719. *
  720. 9994 FORMAT( / ' Tests performed: (S is Tridiag, U is ', A, ',',
  721. $ / 20X, A, ' means ', A, '.', / ' UPLO=''U'':',
  722. $ / ' 1= | A - U S U', A1, ' | / ( |A| n ulp ) ',
  723. $ ' 2= | I - U U', A1, ' | / ( n ulp )', / ' UPLO=''L'':',
  724. $ / ' 3= | A - U S U', A1, ' | / ( |A| n ulp ) ',
  725. $ ' 4= | I - U U', A1, ' | / ( n ulp )' )
  726. 9993 FORMAT( ' N=', I5, ', K=', I4, ', seed=', 4( I4, ',' ), ' type ',
  727. $ I2, ', test(', I2, ')=', G10.3 )
  728. *
  729. * End of ZCHKHB
  730. *
  731. END