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.

dget23.f 28 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875
  1. *> \brief \b DGET23
  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 DGET23( COMP, BALANC, JTYPE, THRESH, ISEED, NOUNIT, N,
  12. * A, LDA, H, WR, WI, WR1, WI1, VL, LDVL, VR,
  13. * LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN,
  14. * RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT,
  15. * WORK, LWORK, IWORK, INFO )
  16. *
  17. * .. Scalar Arguments ..
  18. * LOGICAL COMP
  19. * CHARACTER BALANC
  20. * INTEGER INFO, JTYPE, LDA, LDLRE, LDVL, LDVR, LWORK, N,
  21. * $ NOUNIT
  22. * DOUBLE PRECISION THRESH
  23. * ..
  24. * .. Array Arguments ..
  25. * INTEGER ISEED( 4 ), IWORK( * )
  26. * DOUBLE PRECISION A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
  27. * $ RCDEIN( * ), RCDVIN( * ), RCNDE1( * ),
  28. * $ RCNDV1( * ), RCONDE( * ), RCONDV( * ),
  29. * $ RESULT( 11 ), SCALE( * ), SCALE1( * ),
  30. * $ VL( LDVL, * ), VR( LDVR, * ), WI( * ),
  31. * $ WI1( * ), WORK( * ), WR( * ), WR1( * )
  32. * ..
  33. *
  34. *
  35. *> \par Purpose:
  36. * =============
  37. *>
  38. *> \verbatim
  39. *>
  40. *> DGET23 checks the nonsymmetric eigenvalue problem driver SGEEVX.
  41. *> If COMP = .FALSE., the first 8 of the following tests will be
  42. *> performed on the input matrix A, and also test 9 if LWORK is
  43. *> sufficiently large.
  44. *> if COMP is .TRUE. all 11 tests will be performed.
  45. *>
  46. *> (1) | A * VR - VR * W | / ( n |A| ulp )
  47. *>
  48. *> Here VR is the matrix of unit right eigenvectors.
  49. *> W is a block diagonal matrix, with a 1x1 block for each
  50. *> real eigenvalue and a 2x2 block for each complex conjugate
  51. *> pair. If eigenvalues j and j+1 are a complex conjugate pair,
  52. *> so WR(j) = WR(j+1) = wr and WI(j) = - WI(j+1) = wi, then the
  53. *> 2 x 2 block corresponding to the pair will be:
  54. *>
  55. *> ( wr wi )
  56. *> ( -wi wr )
  57. *>
  58. *> Such a block multiplying an n x 2 matrix ( ur ui ) on the
  59. *> right will be the same as multiplying ur + i*ui by wr + i*wi.
  60. *>
  61. *> (2) | A**H * VL - VL * W**H | / ( n |A| ulp )
  62. *>
  63. *> Here VL is the matrix of unit left eigenvectors, A**H is the
  64. *> conjugate transpose of A, and W is as above.
  65. *>
  66. *> (3) | |VR(i)| - 1 | / ulp and largest component real
  67. *>
  68. *> VR(i) denotes the i-th column of VR.
  69. *>
  70. *> (4) | |VL(i)| - 1 | / ulp and largest component real
  71. *>
  72. *> VL(i) denotes the i-th column of VL.
  73. *>
  74. *> (5) 0 if W(full) = W(partial), 1/ulp otherwise
  75. *>
  76. *> W(full) denotes the eigenvalues computed when VR, VL, RCONDV
  77. *> and RCONDE are also computed, and W(partial) denotes the
  78. *> eigenvalues computed when only some of VR, VL, RCONDV, and
  79. *> RCONDE are computed.
  80. *>
  81. *> (6) 0 if VR(full) = VR(partial), 1/ulp otherwise
  82. *>
  83. *> VR(full) denotes the right eigenvectors computed when VL, RCONDV
  84. *> and RCONDE are computed, and VR(partial) denotes the result
  85. *> when only some of VL and RCONDV are computed.
  86. *>
  87. *> (7) 0 if VL(full) = VL(partial), 1/ulp otherwise
  88. *>
  89. *> VL(full) denotes the left eigenvectors computed when VR, RCONDV
  90. *> and RCONDE are computed, and VL(partial) denotes the result
  91. *> when only some of VR and RCONDV are computed.
  92. *>
  93. *> (8) 0 if SCALE, ILO, IHI, ABNRM (full) =
  94. *> SCALE, ILO, IHI, ABNRM (partial)
  95. *> 1/ulp otherwise
  96. *>
  97. *> SCALE, ILO, IHI and ABNRM describe how the matrix is balanced.
  98. *> (full) is when VR, VL, RCONDE and RCONDV are also computed, and
  99. *> (partial) is when some are not computed.
  100. *>
  101. *> (9) 0 if RCONDV(full) = RCONDV(partial), 1/ulp otherwise
  102. *>
  103. *> RCONDV(full) denotes the reciprocal condition numbers of the
  104. *> right eigenvectors computed when VR, VL and RCONDE are also
  105. *> computed. RCONDV(partial) denotes the reciprocal condition
  106. *> numbers when only some of VR, VL and RCONDE are computed.
  107. *>
  108. *> (10) |RCONDV - RCDVIN| / cond(RCONDV)
  109. *>
  110. *> RCONDV is the reciprocal right eigenvector condition number
  111. *> computed by DGEEVX and RCDVIN (the precomputed true value)
  112. *> is supplied as input. cond(RCONDV) is the condition number of
  113. *> RCONDV, and takes errors in computing RCONDV into account, so
  114. *> that the resulting quantity should be O(ULP). cond(RCONDV) is
  115. *> essentially given by norm(A)/RCONDE.
  116. *>
  117. *> (11) |RCONDE - RCDEIN| / cond(RCONDE)
  118. *>
  119. *> RCONDE is the reciprocal eigenvalue condition number
  120. *> computed by DGEEVX and RCDEIN (the precomputed true value)
  121. *> is supplied as input. cond(RCONDE) is the condition number
  122. *> of RCONDE, and takes errors in computing RCONDE into account,
  123. *> so that the resulting quantity should be O(ULP). cond(RCONDE)
  124. *> is essentially given by norm(A)/RCONDV.
  125. *> \endverbatim
  126. *
  127. * Arguments:
  128. * ==========
  129. *
  130. *> \param[in] COMP
  131. *> \verbatim
  132. *> COMP is LOGICAL
  133. *> COMP describes which input tests to perform:
  134. *> = .FALSE. if the computed condition numbers are not to
  135. *> be tested against RCDVIN and RCDEIN
  136. *> = .TRUE. if they are to be compared
  137. *> \endverbatim
  138. *>
  139. *> \param[in] BALANC
  140. *> \verbatim
  141. *> BALANC is CHARACTER
  142. *> Describes the balancing option to be tested.
  143. *> = 'N' for no permuting or diagonal scaling
  144. *> = 'P' for permuting but no diagonal scaling
  145. *> = 'S' for no permuting but diagonal scaling
  146. *> = 'B' for permuting and diagonal scaling
  147. *> \endverbatim
  148. *>
  149. *> \param[in] JTYPE
  150. *> \verbatim
  151. *> JTYPE is INTEGER
  152. *> Type of input matrix. Used to label output if error occurs.
  153. *> \endverbatim
  154. *>
  155. *> \param[in] THRESH
  156. *> \verbatim
  157. *> THRESH is DOUBLE PRECISION
  158. *> A test will count as "failed" if the "error", computed as
  159. *> described above, exceeds THRESH. Note that the error
  160. *> is scaled to be O(1), so THRESH should be a reasonably
  161. *> small multiple of 1, e.g., 10 or 100. In particular,
  162. *> it should not depend on the precision (single vs. double)
  163. *> or the size of the matrix. It must be at least zero.
  164. *> \endverbatim
  165. *>
  166. *> \param[in] ISEED
  167. *> \verbatim
  168. *> ISEED is INTEGER array, dimension (4)
  169. *> If COMP = .FALSE., the random number generator seed
  170. *> used to produce matrix.
  171. *> If COMP = .TRUE., ISEED(1) = the number of the example.
  172. *> Used to label output if error occurs.
  173. *> \endverbatim
  174. *>
  175. *> \param[in] NOUNIT
  176. *> \verbatim
  177. *> NOUNIT is INTEGER
  178. *> The FORTRAN unit number for printing out error messages
  179. *> (e.g., if a routine returns INFO not equal to 0.)
  180. *> \endverbatim
  181. *>
  182. *> \param[in] N
  183. *> \verbatim
  184. *> N is INTEGER
  185. *> The dimension of A. N must be at least 0.
  186. *> \endverbatim
  187. *>
  188. *> \param[in,out] A
  189. *> \verbatim
  190. *> A is DOUBLE PRECISION array, dimension (LDA,N)
  191. *> Used to hold the matrix whose eigenvalues are to be
  192. *> computed.
  193. *> \endverbatim
  194. *>
  195. *> \param[in] LDA
  196. *> \verbatim
  197. *> LDA is INTEGER
  198. *> The leading dimension of A, and H. LDA must be at
  199. *> least 1 and at least N.
  200. *> \endverbatim
  201. *>
  202. *> \param[out] H
  203. *> \verbatim
  204. *> H is DOUBLE PRECISION array, dimension (LDA,N)
  205. *> Another copy of the test matrix A, modified by DGEEVX.
  206. *> \endverbatim
  207. *>
  208. *> \param[out] WR
  209. *> \verbatim
  210. *> WR is DOUBLE PRECISION array, dimension (N)
  211. *> \endverbatim
  212. *>
  213. *> \param[out] WI
  214. *> \verbatim
  215. *> WI is DOUBLE PRECISION array, dimension (N)
  216. *>
  217. *> The real and imaginary parts of the eigenvalues of A.
  218. *> On exit, WR + WI*i are the eigenvalues of the matrix in A.
  219. *> \endverbatim
  220. *>
  221. *> \param[out] WR1
  222. *> \verbatim
  223. *> WR1 is DOUBLE PRECISION array, dimension (N)
  224. *> \endverbatim
  225. *>
  226. *> \param[out] WI1
  227. *> \verbatim
  228. *> WI1 is DOUBLE PRECISION array, dimension (N)
  229. *>
  230. *> Like WR, WI, these arrays contain the eigenvalues of A,
  231. *> but those computed when DGEEVX only computes a partial
  232. *> eigendecomposition, i.e. not the eigenvalues and left
  233. *> and right eigenvectors.
  234. *> \endverbatim
  235. *>
  236. *> \param[out] VL
  237. *> \verbatim
  238. *> VL is DOUBLE PRECISION array, dimension (LDVL,N)
  239. *> VL holds the computed left eigenvectors.
  240. *> \endverbatim
  241. *>
  242. *> \param[in] LDVL
  243. *> \verbatim
  244. *> LDVL is INTEGER
  245. *> Leading dimension of VL. Must be at least max(1,N).
  246. *> \endverbatim
  247. *>
  248. *> \param[out] VR
  249. *> \verbatim
  250. *> VR is DOUBLE PRECISION array, dimension (LDVR,N)
  251. *> VR holds the computed right eigenvectors.
  252. *> \endverbatim
  253. *>
  254. *> \param[in] LDVR
  255. *> \verbatim
  256. *> LDVR is INTEGER
  257. *> Leading dimension of VR. Must be at least max(1,N).
  258. *> \endverbatim
  259. *>
  260. *> \param[out] LRE
  261. *> \verbatim
  262. *> LRE is DOUBLE PRECISION array, dimension (LDLRE,N)
  263. *> LRE holds the computed right or left eigenvectors.
  264. *> \endverbatim
  265. *>
  266. *> \param[in] LDLRE
  267. *> \verbatim
  268. *> LDLRE is INTEGER
  269. *> Leading dimension of LRE. Must be at least max(1,N).
  270. *> \endverbatim
  271. *>
  272. *> \param[out] RCONDV
  273. *> \verbatim
  274. *> RCONDV is DOUBLE PRECISION array, dimension (N)
  275. *> RCONDV holds the computed reciprocal condition numbers
  276. *> for eigenvectors.
  277. *> \endverbatim
  278. *>
  279. *> \param[out] RCNDV1
  280. *> \verbatim
  281. *> RCNDV1 is DOUBLE PRECISION array, dimension (N)
  282. *> RCNDV1 holds more computed reciprocal condition numbers
  283. *> for eigenvectors.
  284. *> \endverbatim
  285. *>
  286. *> \param[in] RCDVIN
  287. *> \verbatim
  288. *> RCDVIN is DOUBLE PRECISION array, dimension (N)
  289. *> When COMP = .TRUE. RCDVIN holds the precomputed reciprocal
  290. *> condition numbers for eigenvectors to be compared with
  291. *> RCONDV.
  292. *> \endverbatim
  293. *>
  294. *> \param[out] RCONDE
  295. *> \verbatim
  296. *> RCONDE is DOUBLE PRECISION array, dimension (N)
  297. *> RCONDE holds the computed reciprocal condition numbers
  298. *> for eigenvalues.
  299. *> \endverbatim
  300. *>
  301. *> \param[out] RCNDE1
  302. *> \verbatim
  303. *> RCNDE1 is DOUBLE PRECISION array, dimension (N)
  304. *> RCNDE1 holds more computed reciprocal condition numbers
  305. *> for eigenvalues.
  306. *> \endverbatim
  307. *>
  308. *> \param[in] RCDEIN
  309. *> \verbatim
  310. *> RCDEIN is DOUBLE PRECISION array, dimension (N)
  311. *> When COMP = .TRUE. RCDEIN holds the precomputed reciprocal
  312. *> condition numbers for eigenvalues to be compared with
  313. *> RCONDE.
  314. *> \endverbatim
  315. *>
  316. *> \param[out] SCALE
  317. *> \verbatim
  318. *> SCALE is DOUBLE PRECISION array, dimension (N)
  319. *> Holds information describing balancing of matrix.
  320. *> \endverbatim
  321. *>
  322. *> \param[out] SCALE1
  323. *> \verbatim
  324. *> SCALE1 is DOUBLE PRECISION array, dimension (N)
  325. *> Holds information describing balancing of matrix.
  326. *> \endverbatim
  327. *>
  328. *> \param[out] RESULT
  329. *> \verbatim
  330. *> RESULT is DOUBLE PRECISION array, dimension (11)
  331. *> The values computed by the 11 tests described above.
  332. *> The values are currently limited to 1/ulp, to avoid
  333. *> overflow.
  334. *> \endverbatim
  335. *>
  336. *> \param[out] WORK
  337. *> \verbatim
  338. *> WORK is DOUBLE PRECISION array, dimension (LWORK)
  339. *> \endverbatim
  340. *>
  341. *> \param[in] LWORK
  342. *> \verbatim
  343. *> LWORK is INTEGER
  344. *> The number of entries in WORK. This must be at least
  345. *> 3*N, and 6*N+N**2 if tests 9, 10 or 11 are to be performed.
  346. *> \endverbatim
  347. *>
  348. *> \param[out] IWORK
  349. *> \verbatim
  350. *> IWORK is INTEGER array, dimension (2*N)
  351. *> \endverbatim
  352. *>
  353. *> \param[out] INFO
  354. *> \verbatim
  355. *> INFO is INTEGER
  356. *> If 0, successful exit.
  357. *> If <0, input parameter -INFO had an incorrect value.
  358. *> If >0, DGEEVX returned an error code, the absolute
  359. *> value of which is returned.
  360. *> \endverbatim
  361. *
  362. * Authors:
  363. * ========
  364. *
  365. *> \author Univ. of Tennessee
  366. *> \author Univ. of California Berkeley
  367. *> \author Univ. of Colorado Denver
  368. *> \author NAG Ltd.
  369. *
  370. *> \date December 2016
  371. *
  372. *> \ingroup double_eig
  373. *
  374. * =====================================================================
  375. SUBROUTINE DGET23( COMP, BALANC, JTYPE, THRESH, ISEED, NOUNIT, N,
  376. $ A, LDA, H, WR, WI, WR1, WI1, VL, LDVL, VR,
  377. $ LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN,
  378. $ RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT,
  379. $ WORK, LWORK, IWORK, INFO )
  380. *
  381. * -- LAPACK test routine (version 3.7.0) --
  382. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  383. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  384. * December 2016
  385. *
  386. * .. Scalar Arguments ..
  387. LOGICAL COMP
  388. CHARACTER BALANC
  389. INTEGER INFO, JTYPE, LDA, LDLRE, LDVL, LDVR, LWORK, N,
  390. $ NOUNIT
  391. DOUBLE PRECISION THRESH
  392. * ..
  393. * .. Array Arguments ..
  394. INTEGER ISEED( 4 ), IWORK( * )
  395. DOUBLE PRECISION A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
  396. $ RCDEIN( * ), RCDVIN( * ), RCNDE1( * ),
  397. $ RCNDV1( * ), RCONDE( * ), RCONDV( * ),
  398. $ RESULT( 11 ), SCALE( * ), SCALE1( * ),
  399. $ VL( LDVL, * ), VR( LDVR, * ), WI( * ),
  400. $ WI1( * ), WORK( * ), WR( * ), WR1( * )
  401. * ..
  402. *
  403. * =====================================================================
  404. *
  405. *
  406. * .. Parameters ..
  407. DOUBLE PRECISION ZERO, ONE, TWO
  408. PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0, TWO = 2.0D0 )
  409. DOUBLE PRECISION EPSIN
  410. PARAMETER ( EPSIN = 5.9605D-8 )
  411. * ..
  412. * .. Local Scalars ..
  413. LOGICAL BALOK, NOBAL
  414. CHARACTER SENSE
  415. INTEGER I, IHI, IHI1, IINFO, ILO, ILO1, ISENS, ISENSM,
  416. $ J, JJ, KMIN
  417. DOUBLE PRECISION ABNRM, ABNRM1, EPS, SMLNUM, TNRM, TOL, TOLIN,
  418. $ ULP, ULPINV, V, VIMIN, VMAX, VMX, VRMIN, VRMX,
  419. $ VTST
  420. * ..
  421. * .. Local Arrays ..
  422. CHARACTER SENS( 2 )
  423. DOUBLE PRECISION DUM( 1 ), RES( 2 )
  424. * ..
  425. * .. External Functions ..
  426. LOGICAL LSAME
  427. DOUBLE PRECISION DLAMCH, DLAPY2, DNRM2
  428. EXTERNAL LSAME, DLAMCH, DLAPY2, DNRM2
  429. * ..
  430. * .. External Subroutines ..
  431. EXTERNAL DGEEVX, DGET22, DLACPY, XERBLA
  432. * ..
  433. * .. Intrinsic Functions ..
  434. INTRINSIC ABS, DBLE, MAX, MIN
  435. * ..
  436. * .. Data statements ..
  437. DATA SENS / 'N', 'V' /
  438. * ..
  439. * .. Executable Statements ..
  440. *
  441. * Check for errors
  442. *
  443. NOBAL = LSAME( BALANC, 'N' )
  444. BALOK = NOBAL .OR. LSAME( BALANC, 'P' ) .OR.
  445. $ LSAME( BALANC, 'S' ) .OR. LSAME( BALANC, 'B' )
  446. INFO = 0
  447. IF( .NOT.BALOK ) THEN
  448. INFO = -2
  449. ELSE IF( THRESH.LT.ZERO ) THEN
  450. INFO = -4
  451. ELSE IF( NOUNIT.LE.0 ) THEN
  452. INFO = -6
  453. ELSE IF( N.LT.0 ) THEN
  454. INFO = -7
  455. ELSE IF( LDA.LT.1 .OR. LDA.LT.N ) THEN
  456. INFO = -9
  457. ELSE IF( LDVL.LT.1 .OR. LDVL.LT.N ) THEN
  458. INFO = -16
  459. ELSE IF( LDVR.LT.1 .OR. LDVR.LT.N ) THEN
  460. INFO = -18
  461. ELSE IF( LDLRE.LT.1 .OR. LDLRE.LT.N ) THEN
  462. INFO = -20
  463. ELSE IF( LWORK.LT.3*N .OR. ( COMP .AND. LWORK.LT.6*N+N*N ) ) THEN
  464. INFO = -31
  465. END IF
  466. *
  467. IF( INFO.NE.0 ) THEN
  468. CALL XERBLA( 'DGET23', -INFO )
  469. RETURN
  470. END IF
  471. *
  472. * Quick return if nothing to do
  473. *
  474. DO 10 I = 1, 11
  475. RESULT( I ) = -ONE
  476. 10 CONTINUE
  477. *
  478. IF( N.EQ.0 )
  479. $ RETURN
  480. *
  481. * More Important constants
  482. *
  483. ULP = DLAMCH( 'Precision' )
  484. SMLNUM = DLAMCH( 'S' )
  485. ULPINV = ONE / ULP
  486. *
  487. * Compute eigenvalues and eigenvectors, and test them
  488. *
  489. IF( LWORK.GE.6*N+N*N ) THEN
  490. SENSE = 'B'
  491. ISENSM = 2
  492. ELSE
  493. SENSE = 'E'
  494. ISENSM = 1
  495. END IF
  496. CALL DLACPY( 'F', N, N, A, LDA, H, LDA )
  497. CALL DGEEVX( BALANC, 'V', 'V', SENSE, N, H, LDA, WR, WI, VL, LDVL,
  498. $ VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV,
  499. $ WORK, LWORK, IWORK, IINFO )
  500. IF( IINFO.NE.0 ) THEN
  501. RESULT( 1 ) = ULPINV
  502. IF( JTYPE.NE.22 ) THEN
  503. WRITE( NOUNIT, FMT = 9998 )'DGEEVX1', IINFO, N, JTYPE,
  504. $ BALANC, ISEED
  505. ELSE
  506. WRITE( NOUNIT, FMT = 9999 )'DGEEVX1', IINFO, N, ISEED( 1 )
  507. END IF
  508. INFO = ABS( IINFO )
  509. RETURN
  510. END IF
  511. *
  512. * Do Test (1)
  513. *
  514. CALL DGET22( 'N', 'N', 'N', N, A, LDA, VR, LDVR, WR, WI, WORK,
  515. $ RES )
  516. RESULT( 1 ) = RES( 1 )
  517. *
  518. * Do Test (2)
  519. *
  520. CALL DGET22( 'T', 'N', 'T', N, A, LDA, VL, LDVL, WR, WI, WORK,
  521. $ RES )
  522. RESULT( 2 ) = RES( 1 )
  523. *
  524. * Do Test (3)
  525. *
  526. DO 30 J = 1, N
  527. TNRM = ONE
  528. IF( WI( J ).EQ.ZERO ) THEN
  529. TNRM = DNRM2( N, VR( 1, J ), 1 )
  530. ELSE IF( WI( J ).GT.ZERO ) THEN
  531. TNRM = DLAPY2( DNRM2( N, VR( 1, J ), 1 ),
  532. $ DNRM2( N, VR( 1, J+1 ), 1 ) )
  533. END IF
  534. RESULT( 3 ) = MAX( RESULT( 3 ),
  535. $ MIN( ULPINV, ABS( TNRM-ONE ) / ULP ) )
  536. IF( WI( J ).GT.ZERO ) THEN
  537. VMX = ZERO
  538. VRMX = ZERO
  539. DO 20 JJ = 1, N
  540. VTST = DLAPY2( VR( JJ, J ), VR( JJ, J+1 ) )
  541. IF( VTST.GT.VMX )
  542. $ VMX = VTST
  543. IF( VR( JJ, J+1 ).EQ.ZERO .AND. ABS( VR( JJ, J ) ).GT.
  544. $ VRMX )VRMX = ABS( VR( JJ, J ) )
  545. 20 CONTINUE
  546. IF( VRMX / VMX.LT.ONE-TWO*ULP )
  547. $ RESULT( 3 ) = ULPINV
  548. END IF
  549. 30 CONTINUE
  550. *
  551. * Do Test (4)
  552. *
  553. DO 50 J = 1, N
  554. TNRM = ONE
  555. IF( WI( J ).EQ.ZERO ) THEN
  556. TNRM = DNRM2( N, VL( 1, J ), 1 )
  557. ELSE IF( WI( J ).GT.ZERO ) THEN
  558. TNRM = DLAPY2( DNRM2( N, VL( 1, J ), 1 ),
  559. $ DNRM2( N, VL( 1, J+1 ), 1 ) )
  560. END IF
  561. RESULT( 4 ) = MAX( RESULT( 4 ),
  562. $ MIN( ULPINV, ABS( TNRM-ONE ) / ULP ) )
  563. IF( WI( J ).GT.ZERO ) THEN
  564. VMX = ZERO
  565. VRMX = ZERO
  566. DO 40 JJ = 1, N
  567. VTST = DLAPY2( VL( JJ, J ), VL( JJ, J+1 ) )
  568. IF( VTST.GT.VMX )
  569. $ VMX = VTST
  570. IF( VL( JJ, J+1 ).EQ.ZERO .AND. ABS( VL( JJ, J ) ).GT.
  571. $ VRMX )VRMX = ABS( VL( JJ, J ) )
  572. 40 CONTINUE
  573. IF( VRMX / VMX.LT.ONE-TWO*ULP )
  574. $ RESULT( 4 ) = ULPINV
  575. END IF
  576. 50 CONTINUE
  577. *
  578. * Test for all options of computing condition numbers
  579. *
  580. DO 200 ISENS = 1, ISENSM
  581. *
  582. SENSE = SENS( ISENS )
  583. *
  584. * Compute eigenvalues only, and test them
  585. *
  586. CALL DLACPY( 'F', N, N, A, LDA, H, LDA )
  587. CALL DGEEVX( BALANC, 'N', 'N', SENSE, N, H, LDA, WR1, WI1, DUM,
  588. $ 1, DUM, 1, ILO1, IHI1, SCALE1, ABNRM1, RCNDE1,
  589. $ RCNDV1, WORK, LWORK, IWORK, IINFO )
  590. IF( IINFO.NE.0 ) THEN
  591. RESULT( 1 ) = ULPINV
  592. IF( JTYPE.NE.22 ) THEN
  593. WRITE( NOUNIT, FMT = 9998 )'DGEEVX2', IINFO, N, JTYPE,
  594. $ BALANC, ISEED
  595. ELSE
  596. WRITE( NOUNIT, FMT = 9999 )'DGEEVX2', IINFO, N,
  597. $ ISEED( 1 )
  598. END IF
  599. INFO = ABS( IINFO )
  600. GO TO 190
  601. END IF
  602. *
  603. * Do Test (5)
  604. *
  605. DO 60 J = 1, N
  606. IF( WR( J ).NE.WR1( J ) .OR. WI( J ).NE.WI1( J ) )
  607. $ RESULT( 5 ) = ULPINV
  608. 60 CONTINUE
  609. *
  610. * Do Test (8)
  611. *
  612. IF( .NOT.NOBAL ) THEN
  613. DO 70 J = 1, N
  614. IF( SCALE( J ).NE.SCALE1( J ) )
  615. $ RESULT( 8 ) = ULPINV
  616. 70 CONTINUE
  617. IF( ILO.NE.ILO1 )
  618. $ RESULT( 8 ) = ULPINV
  619. IF( IHI.NE.IHI1 )
  620. $ RESULT( 8 ) = ULPINV
  621. IF( ABNRM.NE.ABNRM1 )
  622. $ RESULT( 8 ) = ULPINV
  623. END IF
  624. *
  625. * Do Test (9)
  626. *
  627. IF( ISENS.EQ.2 .AND. N.GT.1 ) THEN
  628. DO 80 J = 1, N
  629. IF( RCONDV( J ).NE.RCNDV1( J ) )
  630. $ RESULT( 9 ) = ULPINV
  631. 80 CONTINUE
  632. END IF
  633. *
  634. * Compute eigenvalues and right eigenvectors, and test them
  635. *
  636. CALL DLACPY( 'F', N, N, A, LDA, H, LDA )
  637. CALL DGEEVX( BALANC, 'N', 'V', SENSE, N, H, LDA, WR1, WI1, DUM,
  638. $ 1, LRE, LDLRE, ILO1, IHI1, SCALE1, ABNRM1, RCNDE1,
  639. $ RCNDV1, WORK, LWORK, IWORK, IINFO )
  640. IF( IINFO.NE.0 ) THEN
  641. RESULT( 1 ) = ULPINV
  642. IF( JTYPE.NE.22 ) THEN
  643. WRITE( NOUNIT, FMT = 9998 )'DGEEVX3', IINFO, N, JTYPE,
  644. $ BALANC, ISEED
  645. ELSE
  646. WRITE( NOUNIT, FMT = 9999 )'DGEEVX3', IINFO, N,
  647. $ ISEED( 1 )
  648. END IF
  649. INFO = ABS( IINFO )
  650. GO TO 190
  651. END IF
  652. *
  653. * Do Test (5) again
  654. *
  655. DO 90 J = 1, N
  656. IF( WR( J ).NE.WR1( J ) .OR. WI( J ).NE.WI1( J ) )
  657. $ RESULT( 5 ) = ULPINV
  658. 90 CONTINUE
  659. *
  660. * Do Test (6)
  661. *
  662. DO 110 J = 1, N
  663. DO 100 JJ = 1, N
  664. IF( VR( J, JJ ).NE.LRE( J, JJ ) )
  665. $ RESULT( 6 ) = ULPINV
  666. 100 CONTINUE
  667. 110 CONTINUE
  668. *
  669. * Do Test (8) again
  670. *
  671. IF( .NOT.NOBAL ) THEN
  672. DO 120 J = 1, N
  673. IF( SCALE( J ).NE.SCALE1( J ) )
  674. $ RESULT( 8 ) = ULPINV
  675. 120 CONTINUE
  676. IF( ILO.NE.ILO1 )
  677. $ RESULT( 8 ) = ULPINV
  678. IF( IHI.NE.IHI1 )
  679. $ RESULT( 8 ) = ULPINV
  680. IF( ABNRM.NE.ABNRM1 )
  681. $ RESULT( 8 ) = ULPINV
  682. END IF
  683. *
  684. * Do Test (9) again
  685. *
  686. IF( ISENS.EQ.2 .AND. N.GT.1 ) THEN
  687. DO 130 J = 1, N
  688. IF( RCONDV( J ).NE.RCNDV1( J ) )
  689. $ RESULT( 9 ) = ULPINV
  690. 130 CONTINUE
  691. END IF
  692. *
  693. * Compute eigenvalues and left eigenvectors, and test them
  694. *
  695. CALL DLACPY( 'F', N, N, A, LDA, H, LDA )
  696. CALL DGEEVX( BALANC, 'V', 'N', SENSE, N, H, LDA, WR1, WI1, LRE,
  697. $ LDLRE, DUM, 1, ILO1, IHI1, SCALE1, ABNRM1, RCNDE1,
  698. $ RCNDV1, WORK, LWORK, IWORK, IINFO )
  699. IF( IINFO.NE.0 ) THEN
  700. RESULT( 1 ) = ULPINV
  701. IF( JTYPE.NE.22 ) THEN
  702. WRITE( NOUNIT, FMT = 9998 )'DGEEVX4', IINFO, N, JTYPE,
  703. $ BALANC, ISEED
  704. ELSE
  705. WRITE( NOUNIT, FMT = 9999 )'DGEEVX4', IINFO, N,
  706. $ ISEED( 1 )
  707. END IF
  708. INFO = ABS( IINFO )
  709. GO TO 190
  710. END IF
  711. *
  712. * Do Test (5) again
  713. *
  714. DO 140 J = 1, N
  715. IF( WR( J ).NE.WR1( J ) .OR. WI( J ).NE.WI1( J ) )
  716. $ RESULT( 5 ) = ULPINV
  717. 140 CONTINUE
  718. *
  719. * Do Test (7)
  720. *
  721. DO 160 J = 1, N
  722. DO 150 JJ = 1, N
  723. IF( VL( J, JJ ).NE.LRE( J, JJ ) )
  724. $ RESULT( 7 ) = ULPINV
  725. 150 CONTINUE
  726. 160 CONTINUE
  727. *
  728. * Do Test (8) again
  729. *
  730. IF( .NOT.NOBAL ) THEN
  731. DO 170 J = 1, N
  732. IF( SCALE( J ).NE.SCALE1( J ) )
  733. $ RESULT( 8 ) = ULPINV
  734. 170 CONTINUE
  735. IF( ILO.NE.ILO1 )
  736. $ RESULT( 8 ) = ULPINV
  737. IF( IHI.NE.IHI1 )
  738. $ RESULT( 8 ) = ULPINV
  739. IF( ABNRM.NE.ABNRM1 )
  740. $ RESULT( 8 ) = ULPINV
  741. END IF
  742. *
  743. * Do Test (9) again
  744. *
  745. IF( ISENS.EQ.2 .AND. N.GT.1 ) THEN
  746. DO 180 J = 1, N
  747. IF( RCONDV( J ).NE.RCNDV1( J ) )
  748. $ RESULT( 9 ) = ULPINV
  749. 180 CONTINUE
  750. END IF
  751. *
  752. 190 CONTINUE
  753. *
  754. 200 CONTINUE
  755. *
  756. * If COMP, compare condition numbers to precomputed ones
  757. *
  758. IF( COMP ) THEN
  759. CALL DLACPY( 'F', N, N, A, LDA, H, LDA )
  760. CALL DGEEVX( 'N', 'V', 'V', 'B', N, H, LDA, WR, WI, VL, LDVL,
  761. $ VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV,
  762. $ WORK, LWORK, IWORK, IINFO )
  763. IF( IINFO.NE.0 ) THEN
  764. RESULT( 1 ) = ULPINV
  765. WRITE( NOUNIT, FMT = 9999 )'DGEEVX5', IINFO, N, ISEED( 1 )
  766. INFO = ABS( IINFO )
  767. GO TO 250
  768. END IF
  769. *
  770. * Sort eigenvalues and condition numbers lexicographically
  771. * to compare with inputs
  772. *
  773. DO 220 I = 1, N - 1
  774. KMIN = I
  775. VRMIN = WR( I )
  776. VIMIN = WI( I )
  777. DO 210 J = I + 1, N
  778. IF( WR( J ).LT.VRMIN ) THEN
  779. KMIN = J
  780. VRMIN = WR( J )
  781. VIMIN = WI( J )
  782. END IF
  783. 210 CONTINUE
  784. WR( KMIN ) = WR( I )
  785. WI( KMIN ) = WI( I )
  786. WR( I ) = VRMIN
  787. WI( I ) = VIMIN
  788. VRMIN = RCONDE( KMIN )
  789. RCONDE( KMIN ) = RCONDE( I )
  790. RCONDE( I ) = VRMIN
  791. VRMIN = RCONDV( KMIN )
  792. RCONDV( KMIN ) = RCONDV( I )
  793. RCONDV( I ) = VRMIN
  794. 220 CONTINUE
  795. *
  796. * Compare condition numbers for eigenvectors
  797. * taking their condition numbers into account
  798. *
  799. RESULT( 10 ) = ZERO
  800. EPS = MAX( EPSIN, ULP )
  801. V = MAX( DBLE( N )*EPS*ABNRM, SMLNUM )
  802. IF( ABNRM.EQ.ZERO )
  803. $ V = ONE
  804. DO 230 I = 1, N
  805. IF( V.GT.RCONDV( I )*RCONDE( I ) ) THEN
  806. TOL = RCONDV( I )
  807. ELSE
  808. TOL = V / RCONDE( I )
  809. END IF
  810. IF( V.GT.RCDVIN( I )*RCDEIN( I ) ) THEN
  811. TOLIN = RCDVIN( I )
  812. ELSE
  813. TOLIN = V / RCDEIN( I )
  814. END IF
  815. TOL = MAX( TOL, SMLNUM / EPS )
  816. TOLIN = MAX( TOLIN, SMLNUM / EPS )
  817. IF( EPS*( RCDVIN( I )-TOLIN ).GT.RCONDV( I )+TOL ) THEN
  818. VMAX = ONE / EPS
  819. ELSE IF( RCDVIN( I )-TOLIN.GT.RCONDV( I )+TOL ) THEN
  820. VMAX = ( RCDVIN( I )-TOLIN ) / ( RCONDV( I )+TOL )
  821. ELSE IF( RCDVIN( I )+TOLIN.LT.EPS*( RCONDV( I )-TOL ) ) THEN
  822. VMAX = ONE / EPS
  823. ELSE IF( RCDVIN( I )+TOLIN.LT.RCONDV( I )-TOL ) THEN
  824. VMAX = ( RCONDV( I )-TOL ) / ( RCDVIN( I )+TOLIN )
  825. ELSE
  826. VMAX = ONE
  827. END IF
  828. RESULT( 10 ) = MAX( RESULT( 10 ), VMAX )
  829. 230 CONTINUE
  830. *
  831. * Compare condition numbers for eigenvalues
  832. * taking their condition numbers into account
  833. *
  834. RESULT( 11 ) = ZERO
  835. DO 240 I = 1, N
  836. IF( V.GT.RCONDV( I ) ) THEN
  837. TOL = ONE
  838. ELSE
  839. TOL = V / RCONDV( I )
  840. END IF
  841. IF( V.GT.RCDVIN( I ) ) THEN
  842. TOLIN = ONE
  843. ELSE
  844. TOLIN = V / RCDVIN( I )
  845. END IF
  846. TOL = MAX( TOL, SMLNUM / EPS )
  847. TOLIN = MAX( TOLIN, SMLNUM / EPS )
  848. IF( EPS*( RCDEIN( I )-TOLIN ).GT.RCONDE( I )+TOL ) THEN
  849. VMAX = ONE / EPS
  850. ELSE IF( RCDEIN( I )-TOLIN.GT.RCONDE( I )+TOL ) THEN
  851. VMAX = ( RCDEIN( I )-TOLIN ) / ( RCONDE( I )+TOL )
  852. ELSE IF( RCDEIN( I )+TOLIN.LT.EPS*( RCONDE( I )-TOL ) ) THEN
  853. VMAX = ONE / EPS
  854. ELSE IF( RCDEIN( I )+TOLIN.LT.RCONDE( I )-TOL ) THEN
  855. VMAX = ( RCONDE( I )-TOL ) / ( RCDEIN( I )+TOLIN )
  856. ELSE
  857. VMAX = ONE
  858. END IF
  859. RESULT( 11 ) = MAX( RESULT( 11 ), VMAX )
  860. 240 CONTINUE
  861. 250 CONTINUE
  862. *
  863. END IF
  864. *
  865. 9999 FORMAT( ' DGET23: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
  866. $ I6, ', INPUT EXAMPLE NUMBER = ', I4 )
  867. 9998 FORMAT( ' DGET23: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
  868. $ I6, ', JTYPE=', I6, ', BALANC = ', A, ', ISEED=(',
  869. $ 3( I5, ',' ), I5, ')' )
  870. *
  871. RETURN
  872. *
  873. * End of DGET23
  874. *
  875. END