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.

cdrvbd.f 49 kB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586878889909192939495969798991001011021031041051061071081091101111121131141151161171181191201211221231241251261271281291301311321331341351361371381391401411421431441451461471481491501511521531541551561571581591601611621631641651661671681691701711721731741751761771781791801811821831841851861871881891901911921931941951961971981992002012022032042052062072082092102112122132142152162172182192202212222232242252262272282292302312322332342352362372382392402412422432442452462472482492502512522532542552562572582592602612622632642652662672682692702712722732742752762772782792802812822832842852862872882892902912922932942952962972982993003013023033043053063073083093103113123133143153163173183193203213223233243253263273283293303313323333343353363373383393403413423433443453463473483493503513523533543553563573583593603613623633643653663673683693703713723733743753763773783793803813823833843853863873883893903913923933943953963973983994004014024034044054064074084094104114124134144154164174184194204214224234244254264274284294304314324334344354364374384394404414424434444454464474484494504514524534544554564574584594604614624634644654664674684694704714724734744754764774784794804814824834844854864874884894904914924934944954964974984995005015025035045055065075085095105115125135145155165175185195205215225235245255265275285295305315325335345355365375385395405415425435445455465475485495505515525535545555565575585595605615625635645655665675685695705715725735745755765775785795805815825835845855865875885895905915925935945955965975985996006016026036046056066076086096106116126136146156166176186196206216226236246256266276286296306316326336346356366376386396406416426436446456466476486496506516526536546556566576586596606616626636646656666676686696706716726736746756766776786796806816826836846856866876886896906916926936946956966976986997007017027037047057067077087097107117127137147157167177187197207217227237247257267277287297307317327337347357367377387397407417427437447457467477487497507517527537547557567577587597607617627637647657667677687697707717727737747757767777787797807817827837847857867877887897907917927937947957967977987998008018028038048058068078088098108118128138148158168178188198208218228238248258268278288298308318328338348358368378388398408418428438448458468478488498508518528538548558568578588598608618628638648658668678688698708718728738748758768778788798808818828838848858868878888898908918928938948958968978988999009019029039049059069079089099109119129139149159169179189199209219229239249259269279289299309319329339349359369379389399409419429439449459469479489499509519529539549559569579589599609619629639649659669679689699709719729739749759769779789799809819829839849859869879889899909919929939949959969979989991000100110021003100410051006100710081009101010111012101310141015101610171018101910201021102210231024102510261027102810291030103110321033103410351036103710381039104010411042104310441045104610471048104910501051105210531054105510561057105810591060106110621063106410651066106710681069107010711072107310741075107610771078107910801081108210831084108510861087108810891090109110921093109410951096109710981099110011011102110311041105110611071108110911101111111211131114111511161117111811191120112111221123112411251126112711281129113011311132113311341135113611371138113911401141114211431144114511461147114811491150115111521153115411551156115711581159116011611162116311641165116611671168116911701171117211731174117511761177117811791180118111821183118411851186118711881189119011911192119311941195119611971198119912001201120212031204120512061207120812091210121112121213121412151216121712181219122012211222122312241225122612271228122912301231123212331234123512361237123812391240124112421243124412451246124712481249125012511252125312541255125612571258125912601261126212631264126512661267126812691270127112721273127412751276127712781279128012811282128312841285128612871288128912901291129212931294129512961297129812991300130113021303130413051306130713081309131013111312131313141315131613171318131913201321132213231324132513261327132813291330133113321333133413351336133713381339134013411342134313441345
  1. *> \brief \b CDRVBD
  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 CDRVBD( NSIZES, MM, NN, NTYPES, DOTYPE, ISEED, THRESH,
  12. * A, LDA, U, LDU, VT, LDVT, ASAV, USAV, VTSAV, S,
  13. * SSAV, E, WORK, LWORK, RWORK, IWORK, NOUNIT,
  14. * INFO )
  15. *
  16. * .. Scalar Arguments ..
  17. * INTEGER INFO, LDA, LDU, LDVT, LWORK, NOUNIT, NSIZES,
  18. * $ NTYPES
  19. * REAL THRESH
  20. * ..
  21. * .. Array Arguments ..
  22. * LOGICAL DOTYPE( * )
  23. * INTEGER ISEED( 4 ), IWORK( * ), MM( * ), NN( * )
  24. * REAL E( * ), RWORK( * ), S( * ), SSAV( * )
  25. * COMPLEX A( LDA, * ), ASAV( LDA, * ), U( LDU, * ),
  26. * $ USAV( LDU, * ), VT( LDVT, * ),
  27. * $ VTSAV( LDVT, * ), WORK( * )
  28. * ..
  29. *
  30. *
  31. *> \par Purpose:
  32. * =============
  33. *>
  34. *> \verbatim
  35. *>
  36. *> CDRVBD checks the singular value decomposition (SVD) driver CGESVD,
  37. *> CGESDD, CGESVJ, CGEJSV, CGESVDX, and CGESVDQ.
  38. *>
  39. *> CGESVD and CGESDD factors A = U diag(S) VT, where U and VT are
  40. *> unitary and diag(S) is diagonal with the entries of the array S on
  41. *> its diagonal. The entries of S are the singular values, nonnegative
  42. *> and stored in decreasing order. U and VT can be optionally not
  43. *> computed, overwritten on A, or computed partially.
  44. *>
  45. *> A is M by N. Let MNMIN = min( M, N ). S has dimension MNMIN.
  46. *> U can be M by M or M by MNMIN. VT can be N by N or MNMIN by N.
  47. *>
  48. *> When CDRVBD is called, a number of matrix "sizes" (M's and N's)
  49. *> and a number of matrix "types" are specified. For each size (M,N)
  50. *> and each type of matrix, and for the minimal workspace as well as
  51. *> workspace adequate to permit blocking, an M x N matrix "A" will be
  52. *> generated and used to test the SVD routines. For each matrix, A will
  53. *> be factored as A = U diag(S) VT and the following 12 tests computed:
  54. *>
  55. *> Test for CGESVD:
  56. *>
  57. *> (1) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
  58. *>
  59. *> (2) | I - U'U | / ( M ulp )
  60. *>
  61. *> (3) | I - VT VT' | / ( N ulp )
  62. *>
  63. *> (4) S contains MNMIN nonnegative values in decreasing order.
  64. *> (Return 0 if true, 1/ULP if false.)
  65. *>
  66. *> (5) | U - Upartial | / ( M ulp ) where Upartial is a partially
  67. *> computed U.
  68. *>
  69. *> (6) | VT - VTpartial | / ( N ulp ) where VTpartial is a partially
  70. *> computed VT.
  71. *>
  72. *> (7) | S - Spartial | / ( MNMIN ulp |S| ) where Spartial is the
  73. *> vector of singular values from the partial SVD
  74. *>
  75. *> Test for CGESDD:
  76. *>
  77. *> (8) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
  78. *>
  79. *> (9) | I - U'U | / ( M ulp )
  80. *>
  81. *> (10) | I - VT VT' | / ( N ulp )
  82. *>
  83. *> (11) S contains MNMIN nonnegative values in decreasing order.
  84. *> (Return 0 if true, 1/ULP if false.)
  85. *>
  86. *> (12) | U - Upartial | / ( M ulp ) where Upartial is a partially
  87. *> computed U.
  88. *>
  89. *> (13) | VT - VTpartial | / ( N ulp ) where VTpartial is a partially
  90. *> computed VT.
  91. *>
  92. *> (14) | S - Spartial | / ( MNMIN ulp |S| ) where Spartial is the
  93. *> vector of singular values from the partial SVD
  94. *>
  95. *> Test for CGESVDQ:
  96. *>
  97. *> (36) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
  98. *>
  99. *> (37) | I - U'U | / ( M ulp )
  100. *>
  101. *> (38) | I - VT VT' | / ( N ulp )
  102. *>
  103. *> (39) S contains MNMIN nonnegative values in decreasing order.
  104. *> (Return 0 if true, 1/ULP if false.)
  105. *>
  106. *> Test for CGESVJ:
  107. *>
  108. *> (15) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
  109. *>
  110. *> (16) | I - U'U | / ( M ulp )
  111. *>
  112. *> (17) | I - VT VT' | / ( N ulp )
  113. *>
  114. *> (18) S contains MNMIN nonnegative values in decreasing order.
  115. *> (Return 0 if true, 1/ULP if false.)
  116. *>
  117. *> Test for CGEJSV:
  118. *>
  119. *> (19) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
  120. *>
  121. *> (20) | I - U'U | / ( M ulp )
  122. *>
  123. *> (21) | I - VT VT' | / ( N ulp )
  124. *>
  125. *> (22) S contains MNMIN nonnegative values in decreasing order.
  126. *> (Return 0 if true, 1/ULP if false.)
  127. *>
  128. *> Test for CGESVDX( 'V', 'V', 'A' )/CGESVDX( 'N', 'N', 'A' )
  129. *>
  130. *> (23) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
  131. *>
  132. *> (24) | I - U'U | / ( M ulp )
  133. *>
  134. *> (25) | I - VT VT' | / ( N ulp )
  135. *>
  136. *> (26) S contains MNMIN nonnegative values in decreasing order.
  137. *> (Return 0 if true, 1/ULP if false.)
  138. *>
  139. *> (27) | U - Upartial | / ( M ulp ) where Upartial is a partially
  140. *> computed U.
  141. *>
  142. *> (28) | VT - VTpartial | / ( N ulp ) where VTpartial is a partially
  143. *> computed VT.
  144. *>
  145. *> (29) | S - Spartial | / ( MNMIN ulp |S| ) where Spartial is the
  146. *> vector of singular values from the partial SVD
  147. *>
  148. *> Test for CGESVDX( 'V', 'V', 'I' )
  149. *>
  150. *> (30) | U' A VT''' - diag(S) | / ( |A| max(M,N) ulp )
  151. *>
  152. *> (31) | I - U'U | / ( M ulp )
  153. *>
  154. *> (32) | I - VT VT' | / ( N ulp )
  155. *>
  156. *> Test for CGESVDX( 'V', 'V', 'V' )
  157. *>
  158. *> (33) | U' A VT''' - diag(S) | / ( |A| max(M,N) ulp )
  159. *>
  160. *> (34) | I - U'U | / ( M ulp )
  161. *>
  162. *> (35) | I - VT VT' | / ( N ulp )
  163. *>
  164. *> The "sizes" are specified by the arrays MM(1:NSIZES) and
  165. *> NN(1:NSIZES); the value of each element pair (MM(j),NN(j))
  166. *> specifies one size. The "types" are specified by a logical array
  167. *> DOTYPE( 1:NTYPES ); if DOTYPE(j) is .TRUE., then matrix type "j"
  168. *> will be generated.
  169. *> Currently, the list of possible types is:
  170. *>
  171. *> (1) The zero matrix.
  172. *> (2) The identity matrix.
  173. *> (3) A matrix of the form U D V, where U and V are unitary and
  174. *> D has evenly spaced entries 1, ..., ULP with random signs
  175. *> on the diagonal.
  176. *> (4) Same as (3), but multiplied by the underflow-threshold / ULP.
  177. *> (5) Same as (3), but multiplied by the overflow-threshold * ULP.
  178. *> \endverbatim
  179. *
  180. * Arguments:
  181. * ==========
  182. *
  183. *> \param[in] NSIZES
  184. *> \verbatim
  185. *> NSIZES is INTEGER
  186. *> The number of sizes of matrices to use. If it is zero,
  187. *> CDRVBD does nothing. It must be at least zero.
  188. *> \endverbatim
  189. *>
  190. *> \param[in] MM
  191. *> \verbatim
  192. *> MM is INTEGER array, dimension (NSIZES)
  193. *> An array containing the matrix "heights" to be used. For
  194. *> each j=1,...,NSIZES, if MM(j) is zero, then MM(j) and NN(j)
  195. *> will be ignored. The MM(j) values must be at least zero.
  196. *> \endverbatim
  197. *>
  198. *> \param[in] NN
  199. *> \verbatim
  200. *> NN is INTEGER array, dimension (NSIZES)
  201. *> An array containing the matrix "widths" to be used. For
  202. *> each j=1,...,NSIZES, if NN(j) is zero, then MM(j) and NN(j)
  203. *> will be ignored. The NN(j) values must be at least zero.
  204. *> \endverbatim
  205. *>
  206. *> \param[in] NTYPES
  207. *> \verbatim
  208. *> NTYPES is INTEGER
  209. *> The number of elements in DOTYPE. If it is zero, CDRVBD
  210. *> does nothing. It must be at least zero. If it is MAXTYP+1
  211. *> and NSIZES is 1, then an additional type, MAXTYP+1 is
  212. *> defined, which is to use whatever matrices are in A and B.
  213. *> This is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
  214. *> DOTYPE(MAXTYP+1) is .TRUE. .
  215. *> \endverbatim
  216. *>
  217. *> \param[in] DOTYPE
  218. *> \verbatim
  219. *> DOTYPE is LOGICAL array, dimension (NTYPES)
  220. *> If DOTYPE(j) is .TRUE., then for each size (m,n), a matrix
  221. *> of type j will be generated. If NTYPES is smaller than the
  222. *> maximum number of types defined (PARAMETER MAXTYP), then
  223. *> types NTYPES+1 through MAXTYP will not be generated. If
  224. *> NTYPES is larger than MAXTYP, DOTYPE(MAXTYP+1) through
  225. *> DOTYPE(NTYPES) will be ignored.
  226. *> \endverbatim
  227. *>
  228. *> \param[in,out] ISEED
  229. *> \verbatim
  230. *> ISEED is INTEGER array, dimension (4)
  231. *> On entry ISEED specifies the seed of the random number
  232. *> generator. The array elements should be between 0 and 4095;
  233. *> if not they will be reduced mod 4096. Also, ISEED(4) must
  234. *> be odd. The random number generator uses a linear
  235. *> congruential sequence limited to small integers, and so
  236. *> should produce machine independent random numbers. The
  237. *> values of ISEED are changed on exit, and can be used in the
  238. *> next call to CDRVBD to continue the same random number
  239. *> sequence.
  240. *> \endverbatim
  241. *>
  242. *> \param[in] THRESH
  243. *> \verbatim
  244. *> THRESH is REAL
  245. *> A test will count as "failed" if the "error", computed as
  246. *> described above, exceeds THRESH. Note that the error
  247. *> is scaled to be O(1), so THRESH should be a reasonably
  248. *> small multiple of 1, e.g., 10 or 100. In particular,
  249. *> it should not depend on the precision (single vs. double)
  250. *> or the size of the matrix. It must be at least zero.
  251. *> \endverbatim
  252. *>
  253. *> \param[out] A
  254. *> \verbatim
  255. *> A is COMPLEX array, dimension (LDA,max(NN))
  256. *> Used to hold the matrix whose singular values are to be
  257. *> computed. On exit, A contains the last matrix actually
  258. *> used.
  259. *> \endverbatim
  260. *>
  261. *> \param[in] LDA
  262. *> \verbatim
  263. *> LDA is INTEGER
  264. *> The leading dimension of A. It must be at
  265. *> least 1 and at least max( MM ).
  266. *> \endverbatim
  267. *>
  268. *> \param[out] U
  269. *> \verbatim
  270. *> U is COMPLEX array, dimension (LDU,max(MM))
  271. *> Used to hold the computed matrix of right singular vectors.
  272. *> On exit, U contains the last such vectors actually computed.
  273. *> \endverbatim
  274. *>
  275. *> \param[in] LDU
  276. *> \verbatim
  277. *> LDU is INTEGER
  278. *> The leading dimension of U. It must be at
  279. *> least 1 and at least max( MM ).
  280. *> \endverbatim
  281. *>
  282. *> \param[out] VT
  283. *> \verbatim
  284. *> VT is COMPLEX array, dimension (LDVT,max(NN))
  285. *> Used to hold the computed matrix of left singular vectors.
  286. *> On exit, VT contains the last such vectors actually computed.
  287. *> \endverbatim
  288. *>
  289. *> \param[in] LDVT
  290. *> \verbatim
  291. *> LDVT is INTEGER
  292. *> The leading dimension of VT. It must be at
  293. *> least 1 and at least max( NN ).
  294. *> \endverbatim
  295. *>
  296. *> \param[out] ASAV
  297. *> \verbatim
  298. *> ASAV is COMPLEX array, dimension (LDA,max(NN))
  299. *> Used to hold a different copy of the matrix whose singular
  300. *> values are to be computed. On exit, A contains the last
  301. *> matrix actually used.
  302. *> \endverbatim
  303. *>
  304. *> \param[out] USAV
  305. *> \verbatim
  306. *> USAV is COMPLEX array, dimension (LDU,max(MM))
  307. *> Used to hold a different copy of the computed matrix of
  308. *> right singular vectors. On exit, USAV contains the last such
  309. *> vectors actually computed.
  310. *> \endverbatim
  311. *>
  312. *> \param[out] VTSAV
  313. *> \verbatim
  314. *> VTSAV is COMPLEX array, dimension (LDVT,max(NN))
  315. *> Used to hold a different copy of the computed matrix of
  316. *> left singular vectors. On exit, VTSAV contains the last such
  317. *> vectors actually computed.
  318. *> \endverbatim
  319. *>
  320. *> \param[out] S
  321. *> \verbatim
  322. *> S is REAL array, dimension (max(min(MM,NN)))
  323. *> Contains the computed singular values.
  324. *> \endverbatim
  325. *>
  326. *> \param[out] SSAV
  327. *> \verbatim
  328. *> SSAV is REAL array, dimension (max(min(MM,NN)))
  329. *> Contains another copy of the computed singular values.
  330. *> \endverbatim
  331. *>
  332. *> \param[out] E
  333. *> \verbatim
  334. *> E is REAL array, dimension (max(min(MM,NN)))
  335. *> Workspace for CGESVD.
  336. *> \endverbatim
  337. *>
  338. *> \param[out] WORK
  339. *> \verbatim
  340. *> WORK is COMPLEX array, dimension (LWORK)
  341. *> \endverbatim
  342. *>
  343. *> \param[in] LWORK
  344. *> \verbatim
  345. *> LWORK is INTEGER
  346. *> The number of entries in WORK. This must be at least
  347. *> MAX(3*MIN(M,N)+MAX(M,N)**2,5*MIN(M,N),3*MAX(M,N)) for all
  348. *> pairs (M,N)=(MM(j),NN(j))
  349. *> \endverbatim
  350. *>
  351. *> \param[out] RWORK
  352. *> \verbatim
  353. *> RWORK is REAL array,
  354. *> dimension ( 5*max(max(MM,NN)) )
  355. *> \endverbatim
  356. *>
  357. *> \param[out] IWORK
  358. *> \verbatim
  359. *> IWORK is INTEGER array, dimension at least 8*min(M,N)
  360. *> \endverbatim
  361. *>
  362. *> \param[in] NOUNIT
  363. *> \verbatim
  364. *> NOUNIT is INTEGER
  365. *> The FORTRAN unit number for printing out error messages
  366. *> (e.g., if a routine returns IINFO not equal to 0.)
  367. *> \endverbatim
  368. *>
  369. *> \param[out] INFO
  370. *> \verbatim
  371. *> INFO is INTEGER
  372. *> If 0, then everything ran OK.
  373. *> -1: NSIZES < 0
  374. *> -2: Some MM(j) < 0
  375. *> -3: Some NN(j) < 0
  376. *> -4: NTYPES < 0
  377. *> -7: THRESH < 0
  378. *> -10: LDA < 1 or LDA < MMAX, where MMAX is max( MM(j) ).
  379. *> -12: LDU < 1 or LDU < MMAX.
  380. *> -14: LDVT < 1 or LDVT < NMAX, where NMAX is max( NN(j) ).
  381. *> -29: LWORK too small.
  382. *> If CLATMS, or CGESVD returns an error code, the
  383. *> absolute value of it is returned.
  384. *> \endverbatim
  385. *
  386. * Authors:
  387. * ========
  388. *
  389. *> \author Univ. of Tennessee
  390. *> \author Univ. of California Berkeley
  391. *> \author Univ. of Colorado Denver
  392. *> \author NAG Ltd.
  393. *
  394. *> \date June 2016
  395. *
  396. *> \ingroup complex_eig
  397. *
  398. * =====================================================================
  399. SUBROUTINE CDRVBD( NSIZES, MM, NN, NTYPES, DOTYPE, ISEED, THRESH,
  400. $ A, LDA, U, LDU, VT, LDVT, ASAV, USAV, VTSAV, S,
  401. $ SSAV, E, WORK, LWORK, RWORK, IWORK, NOUNIT,
  402. $ INFO )
  403. *
  404. * -- LAPACK test routine (version 3.7.0) --
  405. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  406. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  407. * June 2016
  408. *
  409. IMPLICIT NONE
  410. *
  411. * .. Scalar Arguments ..
  412. INTEGER INFO, LDA, LDU, LDVT, LWORK, NOUNIT, NSIZES,
  413. $ NTYPES
  414. REAL THRESH
  415. * ..
  416. * .. Array Arguments ..
  417. LOGICAL DOTYPE( * )
  418. INTEGER ISEED( 4 ), IWORK( * ), MM( * ), NN( * )
  419. REAL E( * ), RWORK( * ), S( * ), SSAV( * )
  420. COMPLEX A( LDA, * ), ASAV( LDA, * ), U( LDU, * ),
  421. $ USAV( LDU, * ), VT( LDVT, * ),
  422. $ VTSAV( LDVT, * ), WORK( * )
  423. * ..
  424. *
  425. * =====================================================================
  426. *
  427. * .. Parameters ..
  428. REAL ZERO, ONE, TWO, HALF
  429. PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0,
  430. $ HALF = 0.5E0 )
  431. COMPLEX CZERO, CONE
  432. PARAMETER ( CZERO = ( 0.0E+0, 0.0E+0 ),
  433. $ CONE = ( 1.0E+0, 0.0E+0 ) )
  434. INTEGER MAXTYP
  435. PARAMETER ( MAXTYP = 5 )
  436. * ..
  437. * .. Local Scalars ..
  438. LOGICAL BADMM, BADNN
  439. CHARACTER JOBQ, JOBU, JOBVT, RANGE
  440. INTEGER I, IINFO, IJQ, IJU, IJVT, IL, IU, ITEMP,
  441. $ IWSPC, IWTMP, J, JSIZE, JTYPE, LSWORK, M,
  442. $ MINWRK, MMAX, MNMAX, MNMIN, MTYPES, N,
  443. $ NERRS, NFAIL, NMAX, NS, NSI, NSV, NTEST,
  444. $ NTESTF, NTESTT, LRWORK
  445. REAL ANORM, DIF, DIV, OVFL, RTUNFL, ULP, ULPINV,
  446. $ UNFL, VL, VU
  447. * ..
  448. * .. Local Scalars for CGESVDQ ..
  449. INTEGER LIWORK, NUMRANK
  450. * ..
  451. * .. Local Arrays ..
  452. CHARACTER CJOB( 4 ), CJOBR( 3 ), CJOBV( 2 )
  453. INTEGER IOLDSD( 4 ), ISEED2( 4 )
  454. REAL RESULT( 39 )
  455. * ..
  456. * .. External Functions ..
  457. REAL SLAMCH, SLARND
  458. EXTERNAL SLAMCH, SLARND
  459. * ..
  460. * .. External Subroutines ..
  461. EXTERNAL ALASVM, XERBLA, CBDT01, CBDT05, CGESDD,
  462. $ CGESVD, CGESVDQ, CGESVJ, CGEJSV, CGESVDX,
  463. $ CLACPY, CLASET, CLATMS, CUNT01, CUNT03
  464. * ..
  465. * .. Intrinsic Functions ..
  466. INTRINSIC ABS, REAL, MAX, MIN
  467. * ..
  468. * .. Scalars in Common ..
  469. CHARACTER*32 SRNAMT
  470. * ..
  471. * .. Common blocks ..
  472. COMMON / SRNAMC / SRNAMT
  473. * ..
  474. * .. Data statements ..
  475. DATA CJOB / 'N', 'O', 'S', 'A' /
  476. DATA CJOBR / 'A', 'V', 'I' /
  477. DATA CJOBV / 'N', 'V' /
  478. * ..
  479. * .. Executable Statements ..
  480. *
  481. * Check for errors
  482. *
  483. INFO = 0
  484. *
  485. * Important constants
  486. *
  487. NERRS = 0
  488. NTESTT = 0
  489. NTESTF = 0
  490. BADMM = .FALSE.
  491. BADNN = .FALSE.
  492. MMAX = 1
  493. NMAX = 1
  494. MNMAX = 1
  495. MINWRK = 1
  496. DO 10 J = 1, NSIZES
  497. MMAX = MAX( MMAX, MM( J ) )
  498. IF( MM( J ).LT.0 )
  499. $ BADMM = .TRUE.
  500. NMAX = MAX( NMAX, NN( J ) )
  501. IF( NN( J ).LT.0 )
  502. $ BADNN = .TRUE.
  503. MNMAX = MAX( MNMAX, MIN( MM( J ), NN( J ) ) )
  504. MINWRK = MAX( MINWRK, MAX( 3*MIN( MM( J ),
  505. $ NN( J ) )+MAX( MM( J ), NN( J ) )**2, 5*MIN( MM( J ),
  506. $ NN( J ) ), 3*MAX( MM( J ), NN( J ) ) ) )
  507. 10 CONTINUE
  508. *
  509. * Check for errors
  510. *
  511. IF( NSIZES.LT.0 ) THEN
  512. INFO = -1
  513. ELSE IF( BADMM ) THEN
  514. INFO = -2
  515. ELSE IF( BADNN ) THEN
  516. INFO = -3
  517. ELSE IF( NTYPES.LT.0 ) THEN
  518. INFO = -4
  519. ELSE IF( LDA.LT.MAX( 1, MMAX ) ) THEN
  520. INFO = -10
  521. ELSE IF( LDU.LT.MAX( 1, MMAX ) ) THEN
  522. INFO = -12
  523. ELSE IF( LDVT.LT.MAX( 1, NMAX ) ) THEN
  524. INFO = -14
  525. ELSE IF( MINWRK.GT.LWORK ) THEN
  526. INFO = -21
  527. END IF
  528. *
  529. IF( INFO.NE.0 ) THEN
  530. CALL XERBLA( 'CDRVBD', -INFO )
  531. RETURN
  532. END IF
  533. *
  534. * Quick return if nothing to do
  535. *
  536. IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 )
  537. $ RETURN
  538. *
  539. * More Important constants
  540. *
  541. UNFL = SLAMCH( 'S' )
  542. OVFL = ONE / UNFL
  543. ULP = SLAMCH( 'E' )
  544. ULPINV = ONE / ULP
  545. RTUNFL = SQRT( UNFL )
  546. *
  547. * Loop over sizes, types
  548. *
  549. NERRS = 0
  550. *
  551. DO 310 JSIZE = 1, NSIZES
  552. M = MM( JSIZE )
  553. N = NN( JSIZE )
  554. MNMIN = MIN( M, N )
  555. *
  556. IF( NSIZES.NE.1 ) THEN
  557. MTYPES = MIN( MAXTYP, NTYPES )
  558. ELSE
  559. MTYPES = MIN( MAXTYP+1, NTYPES )
  560. END IF
  561. *
  562. DO 300 JTYPE = 1, MTYPES
  563. IF( .NOT.DOTYPE( JTYPE ) )
  564. $ GO TO 300
  565. NTEST = 0
  566. *
  567. DO 20 J = 1, 4
  568. IOLDSD( J ) = ISEED( J )
  569. 20 CONTINUE
  570. *
  571. * Compute "A"
  572. *
  573. IF( MTYPES.GT.MAXTYP )
  574. $ GO TO 50
  575. *
  576. IF( JTYPE.EQ.1 ) THEN
  577. *
  578. * Zero matrix
  579. *
  580. CALL CLASET( 'Full', M, N, CZERO, CZERO, A, LDA )
  581. DO 30 I = 1, MIN( M, N )
  582. S( I ) = ZERO
  583. 30 CONTINUE
  584. *
  585. ELSE IF( JTYPE.EQ.2 ) THEN
  586. *
  587. * Identity matrix
  588. *
  589. CALL CLASET( 'Full', M, N, CZERO, CONE, A, LDA )
  590. DO 40 I = 1, MIN( M, N )
  591. S( I ) = ONE
  592. 40 CONTINUE
  593. *
  594. ELSE
  595. *
  596. * (Scaled) random matrix
  597. *
  598. IF( JTYPE.EQ.3 )
  599. $ ANORM = ONE
  600. IF( JTYPE.EQ.4 )
  601. $ ANORM = UNFL / ULP
  602. IF( JTYPE.EQ.5 )
  603. $ ANORM = OVFL*ULP
  604. CALL CLATMS( M, N, 'U', ISEED, 'N', S, 4, REAL( MNMIN ),
  605. $ ANORM, M-1, N-1, 'N', A, LDA, WORK, IINFO )
  606. IF( IINFO.NE.0 ) THEN
  607. WRITE( NOUNIT, FMT = 9996 )'Generator', IINFO, M, N,
  608. $ JTYPE, IOLDSD
  609. INFO = ABS( IINFO )
  610. RETURN
  611. END IF
  612. END IF
  613. *
  614. 50 CONTINUE
  615. CALL CLACPY( 'F', M, N, A, LDA, ASAV, LDA )
  616. *
  617. * Do for minimal and adequate (for blocking) workspace
  618. *
  619. DO 290 IWSPC = 1, 4
  620. *
  621. * Test for CGESVD
  622. *
  623. IWTMP = 2*MIN( M, N )+MAX( M, N )
  624. LSWORK = IWTMP + ( IWSPC-1 )*( LWORK-IWTMP ) / 3
  625. LSWORK = MIN( LSWORK, LWORK )
  626. LSWORK = MAX( LSWORK, 1 )
  627. IF( IWSPC.EQ.4 )
  628. $ LSWORK = LWORK
  629. *
  630. DO 60 J = 1, 35
  631. RESULT( J ) = -ONE
  632. 60 CONTINUE
  633. *
  634. * Factorize A
  635. *
  636. IF( IWSPC.GT.1 )
  637. $ CALL CLACPY( 'F', M, N, ASAV, LDA, A, LDA )
  638. SRNAMT = 'CGESVD'
  639. CALL CGESVD( 'A', 'A', M, N, A, LDA, SSAV, USAV, LDU,
  640. $ VTSAV, LDVT, WORK, LSWORK, RWORK, IINFO )
  641. IF( IINFO.NE.0 ) THEN
  642. WRITE( NOUNIT, FMT = 9995 )'GESVD', IINFO, M, N,
  643. $ JTYPE, LSWORK, IOLDSD
  644. INFO = ABS( IINFO )
  645. RETURN
  646. END IF
  647. *
  648. * Do tests 1--4
  649. *
  650. CALL CBDT01( M, N, 0, ASAV, LDA, USAV, LDU, SSAV, E,
  651. $ VTSAV, LDVT, WORK, RWORK, RESULT( 1 ) )
  652. IF( M.NE.0 .AND. N.NE.0 ) THEN
  653. CALL CUNT01( 'Columns', MNMIN, M, USAV, LDU, WORK,
  654. $ LWORK, RWORK, RESULT( 2 ) )
  655. CALL CUNT01( 'Rows', MNMIN, N, VTSAV, LDVT, WORK,
  656. $ LWORK, RWORK, RESULT( 3 ) )
  657. END IF
  658. RESULT( 4 ) = 0
  659. DO 70 I = 1, MNMIN - 1
  660. IF( SSAV( I ).LT.SSAV( I+1 ) )
  661. $ RESULT( 4 ) = ULPINV
  662. IF( SSAV( I ).LT.ZERO )
  663. $ RESULT( 4 ) = ULPINV
  664. 70 CONTINUE
  665. IF( MNMIN.GE.1 ) THEN
  666. IF( SSAV( MNMIN ).LT.ZERO )
  667. $ RESULT( 4 ) = ULPINV
  668. END IF
  669. *
  670. * Do partial SVDs, comparing to SSAV, USAV, and VTSAV
  671. *
  672. RESULT( 5 ) = ZERO
  673. RESULT( 6 ) = ZERO
  674. RESULT( 7 ) = ZERO
  675. DO 100 IJU = 0, 3
  676. DO 90 IJVT = 0, 3
  677. IF( ( IJU.EQ.3 .AND. IJVT.EQ.3 ) .OR.
  678. $ ( IJU.EQ.1 .AND. IJVT.EQ.1 ) )GO TO 90
  679. JOBU = CJOB( IJU+1 )
  680. JOBVT = CJOB( IJVT+1 )
  681. CALL CLACPY( 'F', M, N, ASAV, LDA, A, LDA )
  682. SRNAMT = 'CGESVD'
  683. CALL CGESVD( JOBU, JOBVT, M, N, A, LDA, S, U, LDU,
  684. $ VT, LDVT, WORK, LSWORK, RWORK, IINFO )
  685. *
  686. * Compare U
  687. *
  688. DIF = ZERO
  689. IF( M.GT.0 .AND. N.GT.0 ) THEN
  690. IF( IJU.EQ.1 ) THEN
  691. CALL CUNT03( 'C', M, MNMIN, M, MNMIN, USAV,
  692. $ LDU, A, LDA, WORK, LWORK, RWORK,
  693. $ DIF, IINFO )
  694. ELSE IF( IJU.EQ.2 ) THEN
  695. CALL CUNT03( 'C', M, MNMIN, M, MNMIN, USAV,
  696. $ LDU, U, LDU, WORK, LWORK, RWORK,
  697. $ DIF, IINFO )
  698. ELSE IF( IJU.EQ.3 ) THEN
  699. CALL CUNT03( 'C', M, M, M, MNMIN, USAV, LDU,
  700. $ U, LDU, WORK, LWORK, RWORK, DIF,
  701. $ IINFO )
  702. END IF
  703. END IF
  704. RESULT( 5 ) = MAX( RESULT( 5 ), DIF )
  705. *
  706. * Compare VT
  707. *
  708. DIF = ZERO
  709. IF( M.GT.0 .AND. N.GT.0 ) THEN
  710. IF( IJVT.EQ.1 ) THEN
  711. CALL CUNT03( 'R', N, MNMIN, N, MNMIN, VTSAV,
  712. $ LDVT, A, LDA, WORK, LWORK,
  713. $ RWORK, DIF, IINFO )
  714. ELSE IF( IJVT.EQ.2 ) THEN
  715. CALL CUNT03( 'R', N, MNMIN, N, MNMIN, VTSAV,
  716. $ LDVT, VT, LDVT, WORK, LWORK,
  717. $ RWORK, DIF, IINFO )
  718. ELSE IF( IJVT.EQ.3 ) THEN
  719. CALL CUNT03( 'R', N, N, N, MNMIN, VTSAV,
  720. $ LDVT, VT, LDVT, WORK, LWORK,
  721. $ RWORK, DIF, IINFO )
  722. END IF
  723. END IF
  724. RESULT( 6 ) = MAX( RESULT( 6 ), DIF )
  725. *
  726. * Compare S
  727. *
  728. DIF = ZERO
  729. DIV = MAX( REAL( MNMIN )*ULP*S( 1 ),
  730. $ SLAMCH( 'Safe minimum' ) )
  731. DO 80 I = 1, MNMIN - 1
  732. IF( SSAV( I ).LT.SSAV( I+1 ) )
  733. $ DIF = ULPINV
  734. IF( SSAV( I ).LT.ZERO )
  735. $ DIF = ULPINV
  736. DIF = MAX( DIF, ABS( SSAV( I )-S( I ) ) / DIV )
  737. 80 CONTINUE
  738. RESULT( 7 ) = MAX( RESULT( 7 ), DIF )
  739. 90 CONTINUE
  740. 100 CONTINUE
  741. *
  742. * Test for CGESDD
  743. *
  744. IWTMP = 2*MNMIN*MNMIN + 2*MNMIN + MAX( M, N )
  745. LSWORK = IWTMP + ( IWSPC-1 )*( LWORK-IWTMP ) / 3
  746. LSWORK = MIN( LSWORK, LWORK )
  747. LSWORK = MAX( LSWORK, 1 )
  748. IF( IWSPC.EQ.4 )
  749. $ LSWORK = LWORK
  750. *
  751. * Factorize A
  752. *
  753. CALL CLACPY( 'F', M, N, ASAV, LDA, A, LDA )
  754. SRNAMT = 'CGESDD'
  755. CALL CGESDD( 'A', M, N, A, LDA, SSAV, USAV, LDU, VTSAV,
  756. $ LDVT, WORK, LSWORK, RWORK, IWORK, IINFO )
  757. IF( IINFO.NE.0 ) THEN
  758. WRITE( NOUNIT, FMT = 9995 )'GESDD', IINFO, M, N,
  759. $ JTYPE, LSWORK, IOLDSD
  760. INFO = ABS( IINFO )
  761. RETURN
  762. END IF
  763. *
  764. * Do tests 1--4
  765. *
  766. CALL CBDT01( M, N, 0, ASAV, LDA, USAV, LDU, SSAV, E,
  767. $ VTSAV, LDVT, WORK, RWORK, RESULT( 8 ) )
  768. IF( M.NE.0 .AND. N.NE.0 ) THEN
  769. CALL CUNT01( 'Columns', MNMIN, M, USAV, LDU, WORK,
  770. $ LWORK, RWORK, RESULT( 9 ) )
  771. CALL CUNT01( 'Rows', MNMIN, N, VTSAV, LDVT, WORK,
  772. $ LWORK, RWORK, RESULT( 10 ) )
  773. END IF
  774. RESULT( 11 ) = 0
  775. DO 110 I = 1, MNMIN - 1
  776. IF( SSAV( I ).LT.SSAV( I+1 ) )
  777. $ RESULT( 11 ) = ULPINV
  778. IF( SSAV( I ).LT.ZERO )
  779. $ RESULT( 11 ) = ULPINV
  780. 110 CONTINUE
  781. IF( MNMIN.GE.1 ) THEN
  782. IF( SSAV( MNMIN ).LT.ZERO )
  783. $ RESULT( 11 ) = ULPINV
  784. END IF
  785. *
  786. * Do partial SVDs, comparing to SSAV, USAV, and VTSAV
  787. *
  788. RESULT( 12 ) = ZERO
  789. RESULT( 13 ) = ZERO
  790. RESULT( 14 ) = ZERO
  791. DO 130 IJQ = 0, 2
  792. JOBQ = CJOB( IJQ+1 )
  793. CALL CLACPY( 'F', M, N, ASAV, LDA, A, LDA )
  794. SRNAMT = 'CGESDD'
  795. CALL CGESDD( JOBQ, M, N, A, LDA, S, U, LDU, VT, LDVT,
  796. $ WORK, LSWORK, RWORK, IWORK, IINFO )
  797. *
  798. * Compare U
  799. *
  800. DIF = ZERO
  801. IF( M.GT.0 .AND. N.GT.0 ) THEN
  802. IF( IJQ.EQ.1 ) THEN
  803. IF( M.GE.N ) THEN
  804. CALL CUNT03( 'C', M, MNMIN, M, MNMIN, USAV,
  805. $ LDU, A, LDA, WORK, LWORK, RWORK,
  806. $ DIF, IINFO )
  807. ELSE
  808. CALL CUNT03( 'C', M, MNMIN, M, MNMIN, USAV,
  809. $ LDU, U, LDU, WORK, LWORK, RWORK,
  810. $ DIF, IINFO )
  811. END IF
  812. ELSE IF( IJQ.EQ.2 ) THEN
  813. CALL CUNT03( 'C', M, MNMIN, M, MNMIN, USAV, LDU,
  814. $ U, LDU, WORK, LWORK, RWORK, DIF,
  815. $ IINFO )
  816. END IF
  817. END IF
  818. RESULT( 12 ) = MAX( RESULT( 12 ), DIF )
  819. *
  820. * Compare VT
  821. *
  822. DIF = ZERO
  823. IF( M.GT.0 .AND. N.GT.0 ) THEN
  824. IF( IJQ.EQ.1 ) THEN
  825. IF( M.GE.N ) THEN
  826. CALL CUNT03( 'R', N, MNMIN, N, MNMIN, VTSAV,
  827. $ LDVT, VT, LDVT, WORK, LWORK,
  828. $ RWORK, DIF, IINFO )
  829. ELSE
  830. CALL CUNT03( 'R', N, MNMIN, N, MNMIN, VTSAV,
  831. $ LDVT, A, LDA, WORK, LWORK,
  832. $ RWORK, DIF, IINFO )
  833. END IF
  834. ELSE IF( IJQ.EQ.2 ) THEN
  835. CALL CUNT03( 'R', N, MNMIN, N, MNMIN, VTSAV,
  836. $ LDVT, VT, LDVT, WORK, LWORK, RWORK,
  837. $ DIF, IINFO )
  838. END IF
  839. END IF
  840. RESULT( 13 ) = MAX( RESULT( 13 ), DIF )
  841. *
  842. * Compare S
  843. *
  844. DIF = ZERO
  845. DIV = MAX( REAL( MNMIN )*ULP*S( 1 ),
  846. $ SLAMCH( 'Safe minimum' ) )
  847. DO 120 I = 1, MNMIN - 1
  848. IF( SSAV( I ).LT.SSAV( I+1 ) )
  849. $ DIF = ULPINV
  850. IF( SSAV( I ).LT.ZERO )
  851. $ DIF = ULPINV
  852. DIF = MAX( DIF, ABS( SSAV( I )-S( I ) ) / DIV )
  853. 120 CONTINUE
  854. RESULT( 14 ) = MAX( RESULT( 14 ), DIF )
  855. 130 CONTINUE
  856. *
  857. * Test CGESVDQ
  858. * Note: CGESVDQ only works for M >= N
  859. *
  860. RESULT( 36 ) = ZERO
  861. RESULT( 37 ) = ZERO
  862. RESULT( 38 ) = ZERO
  863. RESULT( 39 ) = ZERO
  864. *
  865. IF( M.GE.N ) THEN
  866. IWTMP = 2*MNMIN*MNMIN + 2*MNMIN + MAX( M, N )
  867. LSWORK = IWTMP + ( IWSPC-1 )*( LWORK-IWTMP ) / 3
  868. LSWORK = MIN( LSWORK, LWORK )
  869. LSWORK = MAX( LSWORK, 1 )
  870. IF( IWSPC.EQ.4 )
  871. $ LSWORK = LWORK
  872. *
  873. CALL CLACPY( 'F', M, N, ASAV, LDA, A, LDA )
  874. SRNAMT = 'CGESVDQ'
  875. *
  876. LRWORK = MAX(2, M, 5*N)
  877. LIWORK = MAX( N, 1 )
  878. CALL CGESVDQ( 'H', 'N', 'N', 'A', 'A',
  879. $ M, N, A, LDA, SSAV, USAV, LDU,
  880. $ VTSAV, LDVT, NUMRANK, IWORK, LIWORK,
  881. $ WORK, LWORK, RWORK, LRWORK, IINFO )
  882. *
  883. IF( IINFO.NE.0 ) THEN
  884. WRITE( NOUNIT, FMT = 9995 )'CGESVDQ', IINFO, M, N,
  885. $ JTYPE, LSWORK, IOLDSD
  886. INFO = ABS( IINFO )
  887. RETURN
  888. END IF
  889. *
  890. * Do tests 36--39
  891. *
  892. CALL CBDT01( M, N, 0, ASAV, LDA, USAV, LDU, SSAV, E,
  893. $ VTSAV, LDVT, WORK, RWORK, RESULT( 36 ) )
  894. IF( M.NE.0 .AND. N.NE.0 ) THEN
  895. CALL CUNT01( 'Columns', M, M, USAV, LDU, WORK,
  896. $ LWORK, RWORK, RESULT( 37 ) )
  897. CALL CUNT01( 'Rows', N, N, VTSAV, LDVT, WORK,
  898. $ LWORK, RWORK, RESULT( 38 ) )
  899. END IF
  900. RESULT( 39 ) = ZERO
  901. DO 199 I = 1, MNMIN - 1
  902. IF( SSAV( I ).LT.SSAV( I+1 ) )
  903. $ RESULT( 39 ) = ULPINV
  904. IF( SSAV( I ).LT.ZERO )
  905. $ RESULT( 39 ) = ULPINV
  906. 199 CONTINUE
  907. IF( MNMIN.GE.1 ) THEN
  908. IF( SSAV( MNMIN ).LT.ZERO )
  909. $ RESULT( 39 ) = ULPINV
  910. END IF
  911. END IF
  912. *
  913. * Test CGESVJ
  914. * Note: CGESVJ only works for M >= N
  915. *
  916. RESULT( 15 ) = ZERO
  917. RESULT( 16 ) = ZERO
  918. RESULT( 17 ) = ZERO
  919. RESULT( 18 ) = ZERO
  920. *
  921. IF( M.GE.N ) THEN
  922. IWTMP = 2*MNMIN*MNMIN + 2*MNMIN + MAX( M, N )
  923. LSWORK = IWTMP + ( IWSPC-1 )*( LWORK-IWTMP ) / 3
  924. LSWORK = MIN( LSWORK, LWORK )
  925. LSWORK = MAX( LSWORK, 1 )
  926. LRWORK = MAX(6,N)
  927. IF( IWSPC.EQ.4 )
  928. $ LSWORK = LWORK
  929. *
  930. CALL CLACPY( 'F', M, N, ASAV, LDA, USAV, LDA )
  931. SRNAMT = 'CGESVJ'
  932. CALL CGESVJ( 'G', 'U', 'V', M, N, USAV, LDA, SSAV,
  933. & 0, A, LDVT, WORK, LWORK, RWORK,
  934. & LRWORK, IINFO )
  935. *
  936. * CGESVJ returns V not VH
  937. *
  938. DO J=1,N
  939. DO I=1,N
  940. VTSAV(J,I) = CONJG (A(I,J))
  941. END DO
  942. END DO
  943. *
  944. IF( IINFO.NE.0 ) THEN
  945. WRITE( NOUNIT, FMT = 9995 )'GESVJ', IINFO, M, N,
  946. $ JTYPE, LSWORK, IOLDSD
  947. INFO = ABS( IINFO )
  948. RETURN
  949. END IF
  950. *
  951. * Do tests 15--18
  952. *
  953. CALL CBDT01( M, N, 0, ASAV, LDA, USAV, LDU, SSAV, E,
  954. $ VTSAV, LDVT, WORK, RWORK, RESULT( 15 ) )
  955. IF( M.NE.0 .AND. N.NE.0 ) THEN
  956. CALL CUNT01( 'Columns', M, M, USAV, LDU, WORK,
  957. $ LWORK, RWORK, RESULT( 16 ) )
  958. CALL CUNT01( 'Rows', N, N, VTSAV, LDVT, WORK,
  959. $ LWORK, RWORK, RESULT( 17 ) )
  960. END IF
  961. RESULT( 18 ) = ZERO
  962. DO 131 I = 1, MNMIN - 1
  963. IF( SSAV( I ).LT.SSAV( I+1 ) )
  964. $ RESULT( 18 ) = ULPINV
  965. IF( SSAV( I ).LT.ZERO )
  966. $ RESULT( 18 ) = ULPINV
  967. 131 CONTINUE
  968. IF( MNMIN.GE.1 ) THEN
  969. IF( SSAV( MNMIN ).LT.ZERO )
  970. $ RESULT( 18 ) = ULPINV
  971. END IF
  972. END IF
  973. *
  974. * Test CGEJSV
  975. * Note: CGEJSV only works for M >= N
  976. *
  977. RESULT( 19 ) = ZERO
  978. RESULT( 20 ) = ZERO
  979. RESULT( 21 ) = ZERO
  980. RESULT( 22 ) = ZERO
  981. IF( M.GE.N ) THEN
  982. IWTMP = 2*MNMIN*MNMIN + 2*MNMIN + MAX( M, N )
  983. LSWORK = IWTMP + ( IWSPC-1 )*( LWORK-IWTMP ) / 3
  984. LSWORK = MIN( LSWORK, LWORK )
  985. LSWORK = MAX( LSWORK, 1 )
  986. IF( IWSPC.EQ.4 )
  987. $ LSWORK = LWORK
  988. LRWORK = MAX( 7, N + 2*M)
  989. *
  990. CALL CLACPY( 'F', M, N, ASAV, LDA, VTSAV, LDA )
  991. SRNAMT = 'CGEJSV'
  992. CALL CGEJSV( 'G', 'U', 'V', 'R', 'N', 'N',
  993. & M, N, VTSAV, LDA, SSAV, USAV, LDU, A, LDVT,
  994. & WORK, LWORK, RWORK,
  995. & LRWORK, IWORK, IINFO )
  996. *
  997. * CGEJSV returns V not VH
  998. *
  999. DO 133 J=1,N
  1000. DO 132 I=1,N
  1001. VTSAV(J,I) = CONJG (A(I,J))
  1002. 132 END DO
  1003. 133 END DO
  1004. *
  1005. IF( IINFO.NE.0 ) THEN
  1006. WRITE( NOUNIT, FMT = 9995 )'GEJSV', IINFO, M, N,
  1007. $ JTYPE, LSWORK, IOLDSD
  1008. INFO = ABS( IINFO )
  1009. RETURN
  1010. END IF
  1011. *
  1012. * Do tests 19--22
  1013. *
  1014. CALL CBDT01( M, N, 0, ASAV, LDA, USAV, LDU, SSAV, E,
  1015. $ VTSAV, LDVT, WORK, RWORK, RESULT( 19 ) )
  1016. IF( M.NE.0 .AND. N.NE.0 ) THEN
  1017. CALL CUNT01( 'Columns', M, M, USAV, LDU, WORK,
  1018. $ LWORK, RWORK, RESULT( 20 ) )
  1019. CALL CUNT01( 'Rows', N, N, VTSAV, LDVT, WORK,
  1020. $ LWORK, RWORK, RESULT( 21 ) )
  1021. END IF
  1022. RESULT( 22 ) = ZERO
  1023. DO 134 I = 1, MNMIN - 1
  1024. IF( SSAV( I ).LT.SSAV( I+1 ) )
  1025. $ RESULT( 22 ) = ULPINV
  1026. IF( SSAV( I ).LT.ZERO )
  1027. $ RESULT( 22 ) = ULPINV
  1028. 134 CONTINUE
  1029. IF( MNMIN.GE.1 ) THEN
  1030. IF( SSAV( MNMIN ).LT.ZERO )
  1031. $ RESULT( 22 ) = ULPINV
  1032. END IF
  1033. END IF
  1034. *
  1035. * Test CGESVDX
  1036. *
  1037. * Factorize A
  1038. *
  1039. CALL CLACPY( 'F', M, N, ASAV, LDA, A, LDA )
  1040. SRNAMT = 'CGESVDX'
  1041. CALL CGESVDX( 'V', 'V', 'A', M, N, A, LDA,
  1042. $ VL, VU, IL, IU, NS, SSAV, USAV, LDU,
  1043. $ VTSAV, LDVT, WORK, LWORK, RWORK,
  1044. $ IWORK, IINFO )
  1045. IF( IINFO.NE.0 ) THEN
  1046. WRITE( NOUNIT, FMT = 9995 )'GESVDX', IINFO, M, N,
  1047. $ JTYPE, LSWORK, IOLDSD
  1048. INFO = ABS( IINFO )
  1049. RETURN
  1050. END IF
  1051. *
  1052. * Do tests 1--4
  1053. *
  1054. RESULT( 23 ) = ZERO
  1055. RESULT( 24 ) = ZERO
  1056. RESULT( 25 ) = ZERO
  1057. CALL CBDT01( M, N, 0, ASAV, LDA, USAV, LDU, SSAV, E,
  1058. $ VTSAV, LDVT, WORK, RWORK, RESULT( 23 ) )
  1059. IF( M.NE.0 .AND. N.NE.0 ) THEN
  1060. CALL CUNT01( 'Columns', MNMIN, M, USAV, LDU, WORK,
  1061. $ LWORK, RWORK, RESULT( 24 ) )
  1062. CALL CUNT01( 'Rows', MNMIN, N, VTSAV, LDVT, WORK,
  1063. $ LWORK, RWORK, RESULT( 25 ) )
  1064. END IF
  1065. RESULT( 26 ) = ZERO
  1066. DO 140 I = 1, MNMIN - 1
  1067. IF( SSAV( I ).LT.SSAV( I+1 ) )
  1068. $ RESULT( 26 ) = ULPINV
  1069. IF( SSAV( I ).LT.ZERO )
  1070. $ RESULT( 26 ) = ULPINV
  1071. 140 CONTINUE
  1072. IF( MNMIN.GE.1 ) THEN
  1073. IF( SSAV( MNMIN ).LT.ZERO )
  1074. $ RESULT( 26 ) = ULPINV
  1075. END IF
  1076. *
  1077. * Do partial SVDs, comparing to SSAV, USAV, and VTSAV
  1078. *
  1079. RESULT( 27 ) = ZERO
  1080. RESULT( 28 ) = ZERO
  1081. RESULT( 29 ) = ZERO
  1082. DO 170 IJU = 0, 1
  1083. DO 160 IJVT = 0, 1
  1084. IF( ( IJU.EQ.0 .AND. IJVT.EQ.0 ) .OR.
  1085. $ ( IJU.EQ.1 .AND. IJVT.EQ.1 ) ) GO TO 160
  1086. JOBU = CJOBV( IJU+1 )
  1087. JOBVT = CJOBV( IJVT+1 )
  1088. RANGE = CJOBR( 1 )
  1089. CALL CLACPY( 'F', M, N, ASAV, LDA, A, LDA )
  1090. SRNAMT = 'CGESVDX'
  1091. CALL CGESVDX( JOBU, JOBVT, 'A', M, N, A, LDA,
  1092. $ VL, VU, IL, IU, NS, SSAV, U, LDU,
  1093. $ VT, LDVT, WORK, LWORK, RWORK,
  1094. $ IWORK, IINFO )
  1095. *
  1096. * Compare U
  1097. *
  1098. DIF = ZERO
  1099. IF( M.GT.0 .AND. N.GT.0 ) THEN
  1100. IF( IJU.EQ.1 ) THEN
  1101. CALL CUNT03( 'C', M, MNMIN, M, MNMIN, USAV,
  1102. $ LDU, U, LDU, WORK, LWORK, RWORK,
  1103. $ DIF, IINFO )
  1104. END IF
  1105. END IF
  1106. RESULT( 27 ) = MAX( RESULT( 27 ), DIF )
  1107. *
  1108. * Compare VT
  1109. *
  1110. DIF = ZERO
  1111. IF( M.GT.0 .AND. N.GT.0 ) THEN
  1112. IF( IJVT.EQ.1 ) THEN
  1113. CALL CUNT03( 'R', N, MNMIN, N, MNMIN, VTSAV,
  1114. $ LDVT, VT, LDVT, WORK, LWORK,
  1115. $ RWORK, DIF, IINFO )
  1116. END IF
  1117. END IF
  1118. RESULT( 28 ) = MAX( RESULT( 28 ), DIF )
  1119. *
  1120. * Compare S
  1121. *
  1122. DIF = ZERO
  1123. DIV = MAX( REAL( MNMIN )*ULP*S( 1 ),
  1124. $ SLAMCH( 'Safe minimum' ) )
  1125. DO 150 I = 1, MNMIN - 1
  1126. IF( SSAV( I ).LT.SSAV( I+1 ) )
  1127. $ DIF = ULPINV
  1128. IF( SSAV( I ).LT.ZERO )
  1129. $ DIF = ULPINV
  1130. DIF = MAX( DIF, ABS( SSAV( I )-S( I ) ) / DIV )
  1131. 150 CONTINUE
  1132. RESULT( 29) = MAX( RESULT( 29 ), DIF )
  1133. 160 CONTINUE
  1134. 170 CONTINUE
  1135. *
  1136. * Do tests 8--10
  1137. *
  1138. DO 180 I = 1, 4
  1139. ISEED2( I ) = ISEED( I )
  1140. 180 CONTINUE
  1141. IF( MNMIN.LE.1 ) THEN
  1142. IL = 1
  1143. IU = MAX( 1, MNMIN )
  1144. ELSE
  1145. IL = 1 + INT( ( MNMIN-1 )*SLARND( 1, ISEED2 ) )
  1146. IU = 1 + INT( ( MNMIN-1 )*SLARND( 1, ISEED2 ) )
  1147. IF( IU.LT.IL ) THEN
  1148. ITEMP = IU
  1149. IU = IL
  1150. IL = ITEMP
  1151. END IF
  1152. END IF
  1153. CALL CLACPY( 'F', M, N, ASAV, LDA, A, LDA )
  1154. SRNAMT = 'CGESVDX'
  1155. CALL CGESVDX( 'V', 'V', 'I', M, N, A, LDA,
  1156. $ VL, VU, IL, IU, NSI, S, U, LDU,
  1157. $ VT, LDVT, WORK, LWORK, RWORK,
  1158. $ IWORK, IINFO )
  1159. IF( IINFO.NE.0 ) THEN
  1160. WRITE( NOUNIT, FMT = 9995 )'GESVDX', IINFO, M, N,
  1161. $ JTYPE, LSWORK, IOLDSD
  1162. INFO = ABS( IINFO )
  1163. RETURN
  1164. END IF
  1165. *
  1166. RESULT( 30 ) = ZERO
  1167. RESULT( 31 ) = ZERO
  1168. RESULT( 32 ) = ZERO
  1169. CALL CBDT05( M, N, ASAV, LDA, S, NSI, U, LDU,
  1170. $ VT, LDVT, WORK, RESULT( 30 ) )
  1171. IF( M.NE.0 .AND. N.NE.0 ) THEN
  1172. CALL CUNT01( 'Columns', M, NSI, U, LDU, WORK,
  1173. $ LWORK, RWORK, RESULT( 31 ) )
  1174. CALL CUNT01( 'Rows', NSI, N, VT, LDVT, WORK,
  1175. $ LWORK, RWORK, RESULT( 32 ) )
  1176. END IF
  1177. *
  1178. * Do tests 11--13
  1179. *
  1180. IF( MNMIN.GT.0 .AND. NSI.GT.1 ) THEN
  1181. IF( IL.NE.1 ) THEN
  1182. VU = SSAV( IL ) +
  1183. $ MAX( HALF*ABS( SSAV( IL )-SSAV( IL-1 ) ),
  1184. $ ULP*ANORM, TWO*RTUNFL )
  1185. ELSE
  1186. VU = SSAV( 1 ) +
  1187. $ MAX( HALF*ABS( SSAV( NS )-SSAV( 1 ) ),
  1188. $ ULP*ANORM, TWO*RTUNFL )
  1189. END IF
  1190. IF( IU.NE.NS ) THEN
  1191. VL = SSAV( IU ) - MAX( ULP*ANORM, TWO*RTUNFL,
  1192. $ HALF*ABS( SSAV( IU+1 )-SSAV( IU ) ) )
  1193. ELSE
  1194. VL = SSAV( NS ) - MAX( ULP*ANORM, TWO*RTUNFL,
  1195. $ HALF*ABS( SSAV( NS )-SSAV( 1 ) ) )
  1196. END IF
  1197. VL = MAX( VL,ZERO )
  1198. VU = MAX( VU,ZERO )
  1199. IF( VL.GE.VU ) VU = MAX( VU*2, VU+VL+HALF )
  1200. ELSE
  1201. VL = ZERO
  1202. VU = ONE
  1203. END IF
  1204. CALL CLACPY( 'F', M, N, ASAV, LDA, A, LDA )
  1205. SRNAMT = 'CGESVDX'
  1206. CALL CGESVDX( 'V', 'V', 'V', M, N, A, LDA,
  1207. $ VL, VU, IL, IU, NSV, S, U, LDU,
  1208. $ VT, LDVT, WORK, LWORK, RWORK,
  1209. $ IWORK, IINFO )
  1210. IF( IINFO.NE.0 ) THEN
  1211. WRITE( NOUNIT, FMT = 9995 )'GESVDX', IINFO, M, N,
  1212. $ JTYPE, LSWORK, IOLDSD
  1213. INFO = ABS( IINFO )
  1214. RETURN
  1215. END IF
  1216. *
  1217. RESULT( 33 ) = ZERO
  1218. RESULT( 34 ) = ZERO
  1219. RESULT( 35 ) = ZERO
  1220. CALL CBDT05( M, N, ASAV, LDA, S, NSV, U, LDU,
  1221. $ VT, LDVT, WORK, RESULT( 33 ) )
  1222. IF( M.NE.0 .AND. N.NE.0 ) THEN
  1223. CALL CUNT01( 'Columns', M, NSV, U, LDU, WORK,
  1224. $ LWORK, RWORK, RESULT( 34 ) )
  1225. CALL CUNT01( 'Rows', NSV, N, VT, LDVT, WORK,
  1226. $ LWORK, RWORK, RESULT( 35 ) )
  1227. END IF
  1228. *
  1229. * End of Loop -- Check for RESULT(j) > THRESH
  1230. *
  1231. NTEST = 0
  1232. NFAIL = 0
  1233. DO 190 J = 1, 39
  1234. IF( RESULT( J ).GE.ZERO )
  1235. $ NTEST = NTEST + 1
  1236. IF( RESULT( J ).GE.THRESH )
  1237. $ NFAIL = NFAIL + 1
  1238. 190 CONTINUE
  1239. *
  1240. IF( NFAIL.GT.0 )
  1241. $ NTESTF = NTESTF + 1
  1242. IF( NTESTF.EQ.1 ) THEN
  1243. WRITE( NOUNIT, FMT = 9999 )
  1244. WRITE( NOUNIT, FMT = 9998 )THRESH
  1245. NTESTF = 2
  1246. END IF
  1247. *
  1248. DO 200 J = 1, 39
  1249. IF( RESULT( J ).GE.THRESH ) THEN
  1250. WRITE( NOUNIT, FMT = 9997 )M, N, JTYPE, IWSPC,
  1251. $ IOLDSD, J, RESULT( J )
  1252. END IF
  1253. 200 CONTINUE
  1254. *
  1255. NERRS = NERRS + NFAIL
  1256. NTESTT = NTESTT + NTEST
  1257. *
  1258. 290 CONTINUE
  1259. *
  1260. 300 CONTINUE
  1261. 310 CONTINUE
  1262. *
  1263. * Summary
  1264. *
  1265. CALL ALASVM( 'CBD', NOUNIT, NERRS, NTESTT, 0 )
  1266. *
  1267. 9999 FORMAT( ' SVD -- Complex Singular Value Decomposition Driver ',
  1268. $ / ' Matrix types (see CDRVBD for details):',
  1269. $ / / ' 1 = Zero matrix', / ' 2 = Identity matrix',
  1270. $ / ' 3 = Evenly spaced singular values near 1',
  1271. $ / ' 4 = Evenly spaced singular values near underflow',
  1272. $ / ' 5 = Evenly spaced singular values near overflow',
  1273. $ / / ' Tests performed: ( A is dense, U and V are unitary,',
  1274. $ / 19X, ' S is an array, and Upartial, VTpartial, and',
  1275. $ / 19X, ' Spartial are partially computed U, VT and S),', / )
  1276. 9998 FORMAT( ' Tests performed with Test Threshold = ', F8.2,
  1277. $ / ' CGESVD: ', /
  1278. $ ' 1 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
  1279. $ / ' 2 = | I - U**T U | / ( M ulp ) ',
  1280. $ / ' 3 = | I - VT VT**T | / ( N ulp ) ',
  1281. $ / ' 4 = 0 if S contains min(M,N) nonnegative values in',
  1282. $ ' decreasing order, else 1/ulp',
  1283. $ / ' 5 = | U - Upartial | / ( M ulp )',
  1284. $ / ' 6 = | VT - VTpartial | / ( N ulp )',
  1285. $ / ' 7 = | S - Spartial | / ( min(M,N) ulp |S| )',
  1286. $ / ' CGESDD: ', /
  1287. $ ' 8 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
  1288. $ / ' 9 = | I - U**T U | / ( M ulp ) ',
  1289. $ / '10 = | I - VT VT**T | / ( N ulp ) ',
  1290. $ / '11 = 0 if S contains min(M,N) nonnegative values in',
  1291. $ ' decreasing order, else 1/ulp',
  1292. $ / '12 = | U - Upartial | / ( M ulp )',
  1293. $ / '13 = | VT - VTpartial | / ( N ulp )',
  1294. $ / '14 = | S - Spartial | / ( min(M,N) ulp |S| )',
  1295. $ / ' CGESVJ: ', /
  1296. $ / '15 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
  1297. $ / '16 = | I - U**T U | / ( M ulp ) ',
  1298. $ / '17 = | I - VT VT**T | / ( N ulp ) ',
  1299. $ / '18 = 0 if S contains min(M,N) nonnegative values in',
  1300. $ ' decreasing order, else 1/ulp',
  1301. $ / ' CGESJV: ', /
  1302. $ / '19 = | A - U diag(S) VT | / ( |A| max(M,N) ulp )',
  1303. $ / '20 = | I - U**T U | / ( M ulp ) ',
  1304. $ / '21 = | I - VT VT**T | / ( N ulp ) ',
  1305. $ / '22 = 0 if S contains min(M,N) nonnegative values in',
  1306. $ ' decreasing order, else 1/ulp',
  1307. $ / ' CGESVDX(V,V,A): ', /
  1308. $ '23 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
  1309. $ / '24 = | I - U**T U | / ( M ulp ) ',
  1310. $ / '25 = | I - VT VT**T | / ( N ulp ) ',
  1311. $ / '26 = 0 if S contains min(M,N) nonnegative values in',
  1312. $ ' decreasing order, else 1/ulp',
  1313. $ / '27 = | U - Upartial | / ( M ulp )',
  1314. $ / '28 = | VT - VTpartial | / ( N ulp )',
  1315. $ / '29 = | S - Spartial | / ( min(M,N) ulp |S| )',
  1316. $ / ' CGESVDX(V,V,I): ',
  1317. $ / '30 = | U**T A VT**T - diag(S) | / ( |A| max(M,N) ulp )',
  1318. $ / '31 = | I - U**T U | / ( M ulp ) ',
  1319. $ / '32 = | I - VT VT**T | / ( N ulp ) ',
  1320. $ / ' CGESVDX(V,V,V) ',
  1321. $ / '33 = | U**T A VT**T - diag(S) | / ( |A| max(M,N) ulp )',
  1322. $ / '34 = | I - U**T U | / ( M ulp ) ',
  1323. $ / '35 = | I - VT VT**T | / ( N ulp ) ',
  1324. $ ' CGESVDQ(H,N,N,A,A',
  1325. $ / '36 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
  1326. $ / '37 = | I - U**T U | / ( M ulp ) ',
  1327. $ / '38 = | I - VT VT**T | / ( N ulp ) ',
  1328. $ / '39 = 0 if S contains min(M,N) nonnegative values in',
  1329. $ ' decreasing order, else 1/ulp',
  1330. $ / / )
  1331. 9997 FORMAT( ' M=', I5, ', N=', I5, ', type ', I1, ', IWS=', I1,
  1332. $ ', seed=', 4( I4, ',' ), ' test(', I2, ')=', G11.4 )
  1333. 9996 FORMAT( ' CDRVBD: ', A, ' returned INFO=', I6, '.', / 9X, 'M=',
  1334. $ I6, ', N=', I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ),
  1335. $ I5, ')' )
  1336. 9995 FORMAT( ' CDRVBD: ', A, ' returned INFO=', I6, '.', / 9X, 'M=',
  1337. $ I6, ', N=', I6, ', JTYPE=', I6, ', LSWORK=', I6, / 9X,
  1338. $ 'ISEED=(', 3( I5, ',' ), I5, ')' )
  1339. *
  1340. RETURN
  1341. *
  1342. * End of CDRVBD
  1343. *
  1344. END