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.

sdrvbd.f 43 kB

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