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.

zdrvls.f 25 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668
  1. *> \brief \b ZDRVLS
  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 ZDRVLS( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL, NNB,
  12. * NBVAL, NXVAL, THRESH, TSTERR, A, COPYA, B,
  13. * COPYB, C, S, COPYS, WORK, RWORK, IWORK, NOUT )
  14. *
  15. * .. Scalar Arguments ..
  16. * LOGICAL TSTERR
  17. * INTEGER NM, NN, NNB, NNS, NOUT
  18. * DOUBLE PRECISION THRESH
  19. * ..
  20. * .. Array Arguments ..
  21. * LOGICAL DOTYPE( * )
  22. * INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
  23. * $ NVAL( * ), NXVAL( * )
  24. * DOUBLE PRECISION COPYS( * ), RWORK( * ), S( * )
  25. * COMPLEX*16 A( * ), B( * ), C( * ), COPYA( * ), COPYB( * ),
  26. * $ WORK( * )
  27. * ..
  28. *
  29. *
  30. *> \par Purpose:
  31. * =============
  32. *>
  33. *> \verbatim
  34. *>
  35. *> ZDRVLS tests the least squares driver routines ZGELS, CGELSS, ZGELSY
  36. *> and CGELSD.
  37. *> \endverbatim
  38. *
  39. * Arguments:
  40. * ==========
  41. *
  42. *> \param[in] DOTYPE
  43. *> \verbatim
  44. *> DOTYPE is LOGICAL array, dimension (NTYPES)
  45. *> The matrix types to be used for testing. Matrices of type j
  46. *> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
  47. *> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
  48. *> The matrix of type j is generated as follows:
  49. *> j=1: A = U*D*V where U and V are random unitary matrices
  50. *> and D has random entries (> 0.1) taken from a uniform
  51. *> distribution (0,1). A is full rank.
  52. *> j=2: The same of 1, but A is scaled up.
  53. *> j=3: The same of 1, but A is scaled down.
  54. *> j=4: A = U*D*V where U and V are random unitary matrices
  55. *> and D has 3*min(M,N)/4 random entries (> 0.1) taken
  56. *> from a uniform distribution (0,1) and the remaining
  57. *> entries set to 0. A is rank-deficient.
  58. *> j=5: The same of 4, but A is scaled up.
  59. *> j=6: The same of 5, but A is scaled down.
  60. *> \endverbatim
  61. *>
  62. *> \param[in] NM
  63. *> \verbatim
  64. *> NM is INTEGER
  65. *> The number of values of M contained in the vector MVAL.
  66. *> \endverbatim
  67. *>
  68. *> \param[in] MVAL
  69. *> \verbatim
  70. *> MVAL is INTEGER array, dimension (NM)
  71. *> The values of the matrix row dimension M.
  72. *> \endverbatim
  73. *>
  74. *> \param[in] NN
  75. *> \verbatim
  76. *> NN is INTEGER
  77. *> The number of values of N contained in the vector NVAL.
  78. *> \endverbatim
  79. *>
  80. *> \param[in] NVAL
  81. *> \verbatim
  82. *> NVAL is INTEGER array, dimension (NN)
  83. *> The values of the matrix column dimension N.
  84. *> \endverbatim
  85. *>
  86. *> \param[in] NNB
  87. *> \verbatim
  88. *> NNB is INTEGER
  89. *> The number of values of NB and NX contained in the
  90. *> vectors NBVAL and NXVAL. The blocking parameters are used
  91. *> in pairs (NB,NX).
  92. *> \endverbatim
  93. *>
  94. *> \param[in] NBVAL
  95. *> \verbatim
  96. *> NBVAL is INTEGER array, dimension (NNB)
  97. *> The values of the blocksize NB.
  98. *> \endverbatim
  99. *>
  100. *> \param[in] NXVAL
  101. *> \verbatim
  102. *> NXVAL is INTEGER array, dimension (NNB)
  103. *> The values of the crossover point NX.
  104. *> \endverbatim
  105. *>
  106. *> \param[in] NNS
  107. *> \verbatim
  108. *> NNS is INTEGER
  109. *> The number of values of NRHS contained in the vector NSVAL.
  110. *> \endverbatim
  111. *>
  112. *> \param[in] NSVAL
  113. *> \verbatim
  114. *> NSVAL is INTEGER array, dimension (NNS)
  115. *> The values of the number of right hand sides NRHS.
  116. *> \endverbatim
  117. *>
  118. *> \param[in] THRESH
  119. *> \verbatim
  120. *> THRESH is DOUBLE PRECISION
  121. *> The threshold value for the test ratios. A result is
  122. *> included in the output file if RESULT >= THRESH. To have
  123. *> every test ratio printed, use THRESH = 0.
  124. *> \endverbatim
  125. *>
  126. *> \param[in] TSTERR
  127. *> \verbatim
  128. *> TSTERR is LOGICAL
  129. *> Flag that indicates whether error exits are to be tested.
  130. *> \endverbatim
  131. *>
  132. *> \param[out] A
  133. *> \verbatim
  134. *> A is COMPLEX*16 array, dimension (MMAX*NMAX)
  135. *> where MMAX is the maximum value of M in MVAL and NMAX is the
  136. *> maximum value of N in NVAL.
  137. *> \endverbatim
  138. *>
  139. *> \param[out] COPYA
  140. *> \verbatim
  141. *> COPYA is COMPLEX*16 array, dimension (MMAX*NMAX)
  142. *> \endverbatim
  143. *>
  144. *> \param[out] B
  145. *> \verbatim
  146. *> B is COMPLEX*16 array, dimension (MMAX*NSMAX)
  147. *> where MMAX is the maximum value of M in MVAL and NSMAX is the
  148. *> maximum value of NRHS in NSVAL.
  149. *> \endverbatim
  150. *>
  151. *> \param[out] COPYB
  152. *> \verbatim
  153. *> COPYB is COMPLEX*16 array, dimension (MMAX*NSMAX)
  154. *> \endverbatim
  155. *>
  156. *> \param[out] C
  157. *> \verbatim
  158. *> C is COMPLEX*16 array, dimension (MMAX*NSMAX)
  159. *> \endverbatim
  160. *>
  161. *> \param[out] S
  162. *> \verbatim
  163. *> S is DOUBLE PRECISION array, dimension
  164. *> (min(MMAX,NMAX))
  165. *> \endverbatim
  166. *>
  167. *> \param[out] COPYS
  168. *> \verbatim
  169. *> COPYS is DOUBLE PRECISION array, dimension
  170. *> (min(MMAX,NMAX))
  171. *> \endverbatim
  172. *>
  173. *> \param[out] WORK
  174. *> \verbatim
  175. *> WORK is COMPLEX*16 array, dimension
  176. *> (MMAX*NMAX + 4*NMAX + MMAX).
  177. *> \endverbatim
  178. *>
  179. *> \param[out] RWORK
  180. *> \verbatim
  181. *> RWORK is DOUBLE PRECISION array, dimension (5*NMAX-1)
  182. *> \endverbatim
  183. *>
  184. *> \param[out] IWORK
  185. *> \verbatim
  186. *> IWORK is INTEGER array, dimension (15*NMAX)
  187. *> \endverbatim
  188. *>
  189. *> \param[in] NOUT
  190. *> \verbatim
  191. *> NOUT is INTEGER
  192. *> The unit number for output.
  193. *> \endverbatim
  194. *
  195. * Authors:
  196. * ========
  197. *
  198. *> \author Univ. of Tennessee
  199. *> \author Univ. of California Berkeley
  200. *> \author Univ. of Colorado Denver
  201. *> \author NAG Ltd.
  202. *
  203. *> \date November 2015
  204. *
  205. *> \ingroup complex16_lin
  206. *
  207. * =====================================================================
  208. SUBROUTINE ZDRVLS( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL, NNB,
  209. $ NBVAL, NXVAL, THRESH, TSTERR, A, COPYA, B,
  210. $ COPYB, C, S, COPYS, WORK, RWORK, IWORK, NOUT )
  211. *
  212. * -- LAPACK test routine (version 3.6.0) --
  213. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  214. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  215. * November 2015
  216. *
  217. * .. Scalar Arguments ..
  218. LOGICAL TSTERR
  219. INTEGER NM, NN, NNB, NNS, NOUT
  220. DOUBLE PRECISION THRESH
  221. * ..
  222. * .. Array Arguments ..
  223. LOGICAL DOTYPE( * )
  224. INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
  225. $ NVAL( * ), NXVAL( * )
  226. DOUBLE PRECISION COPYS( * ), RWORK( * ), S( * )
  227. COMPLEX*16 A( * ), B( * ), C( * ), COPYA( * ), COPYB( * ),
  228. $ WORK( * )
  229. * ..
  230. *
  231. * =====================================================================
  232. *
  233. * .. Parameters ..
  234. INTEGER NTESTS
  235. PARAMETER ( NTESTS = 14 )
  236. INTEGER SMLSIZ
  237. PARAMETER ( SMLSIZ = 25 )
  238. DOUBLE PRECISION ONE, ZERO
  239. PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
  240. COMPLEX*16 CONE, CZERO
  241. PARAMETER ( CONE = ( 1.0D+0, 0.0D+0 ),
  242. $ CZERO = ( 0.0D+0, 0.0D+0 ) )
  243. * ..
  244. * .. Local Scalars ..
  245. CHARACTER TRANS
  246. CHARACTER*3 PATH
  247. INTEGER CRANK, I, IM, IN, INB, INFO, INS, IRANK,
  248. $ ISCALE, ITRAN, ITYPE, J, K, LDA, LDB, LDWORK,
  249. $ LWLSY, LWORK, M, MNMIN, N, NB, NCOLS, NERRS,
  250. $ NFAIL, NRHS, NROWS, NRUN, RANK
  251. DOUBLE PRECISION EPS, NORMA, NORMB, RCOND
  252. * ..
  253. * .. Local Arrays ..
  254. INTEGER ISEED( 4 ), ISEEDY( 4 )
  255. DOUBLE PRECISION RESULT( NTESTS )
  256. * ..
  257. * .. External Functions ..
  258. DOUBLE PRECISION DASUM, DLAMCH, ZQRT12, ZQRT14, ZQRT17
  259. EXTERNAL DASUM, DLAMCH, ZQRT12, ZQRT14, ZQRT17
  260. * ..
  261. * .. External Subroutines ..
  262. EXTERNAL ALAERH, ALAHD, ALASVM, DAXPY, DLASRT, XLAENV,
  263. $ ZDSCAL, ZERRLS, ZGELS, ZGELSD, ZGELSS,
  264. $ ZGELSY, ZGEMM, ZLACPY, ZLARNV, ZQRT13, ZQRT15,
  265. $ ZQRT16
  266. * ..
  267. * .. Intrinsic Functions ..
  268. INTRINSIC DBLE, MAX, MIN, SQRT
  269. * ..
  270. * .. Scalars in Common ..
  271. LOGICAL LERR, OK
  272. CHARACTER*32 SRNAMT
  273. INTEGER INFOT, IOUNIT
  274. * ..
  275. * .. Common blocks ..
  276. COMMON / INFOC / INFOT, IOUNIT, OK, LERR
  277. COMMON / SRNAMC / SRNAMT
  278. * ..
  279. * .. Data statements ..
  280. DATA ISEEDY / 1988, 1989, 1990, 1991 /
  281. * ..
  282. * .. Executable Statements ..
  283. *
  284. * Initialize constants and the random number seed.
  285. *
  286. PATH( 1: 1 ) = 'Zomplex precision'
  287. PATH( 2: 3 ) = 'LS'
  288. NRUN = 0
  289. NFAIL = 0
  290. NERRS = 0
  291. DO 10 I = 1, 4
  292. ISEED( I ) = ISEEDY( I )
  293. 10 CONTINUE
  294. EPS = DLAMCH( 'Epsilon' )
  295. *
  296. * Threshold for rank estimation
  297. *
  298. RCOND = SQRT( EPS ) - ( SQRT( EPS )-EPS ) / 2
  299. *
  300. * Test the error exits
  301. *
  302. CALL XLAENV( 9, SMLSIZ )
  303. IF( TSTERR )
  304. $ CALL ZERRLS( PATH, NOUT )
  305. *
  306. * Print the header if NM = 0 or NN = 0 and THRESH = 0.
  307. *
  308. IF( ( NM.EQ.0 .OR. NN.EQ.0 ) .AND. THRESH.EQ.ZERO )
  309. $ CALL ALAHD( NOUT, PATH )
  310. INFOT = 0
  311. *
  312. DO 140 IM = 1, NM
  313. M = MVAL( IM )
  314. LDA = MAX( 1, M )
  315. *
  316. DO 130 IN = 1, NN
  317. N = NVAL( IN )
  318. MNMIN = MIN( M, N )
  319. LDB = MAX( 1, M, N )
  320. *
  321. DO 120 INS = 1, NNS
  322. NRHS = NSVAL( INS )
  323. LWORK = MAX( 1, ( M+NRHS )*( N+2 ), ( N+NRHS )*( M+2 ),
  324. $ M*N+4*MNMIN+MAX( M, N ), 2*N+M )
  325. *
  326. DO 110 IRANK = 1, 2
  327. DO 100 ISCALE = 1, 3
  328. ITYPE = ( IRANK-1 )*3 + ISCALE
  329. IF( .NOT.DOTYPE( ITYPE ) )
  330. $ GO TO 100
  331. *
  332. IF( IRANK.EQ.1 ) THEN
  333. *
  334. * Test ZGELS
  335. *
  336. * Generate a matrix of scaling type ISCALE
  337. *
  338. CALL ZQRT13( ISCALE, M, N, COPYA, LDA, NORMA,
  339. $ ISEED )
  340. DO 40 INB = 1, NNB
  341. NB = NBVAL( INB )
  342. CALL XLAENV( 1, NB )
  343. CALL XLAENV( 3, NXVAL( INB ) )
  344. *
  345. DO 30 ITRAN = 1, 2
  346. IF( ITRAN.EQ.1 ) THEN
  347. TRANS = 'N'
  348. NROWS = M
  349. NCOLS = N
  350. ELSE
  351. TRANS = 'C'
  352. NROWS = N
  353. NCOLS = M
  354. END IF
  355. LDWORK = MAX( 1, NCOLS )
  356. *
  357. * Set up a consistent rhs
  358. *
  359. IF( NCOLS.GT.0 ) THEN
  360. CALL ZLARNV( 2, ISEED, NCOLS*NRHS,
  361. $ WORK )
  362. CALL ZDSCAL( NCOLS*NRHS,
  363. $ ONE / DBLE( NCOLS ), WORK,
  364. $ 1 )
  365. END IF
  366. CALL ZGEMM( TRANS, 'No transpose', NROWS,
  367. $ NRHS, NCOLS, CONE, COPYA, LDA,
  368. $ WORK, LDWORK, CZERO, B, LDB )
  369. CALL ZLACPY( 'Full', NROWS, NRHS, B, LDB,
  370. $ COPYB, LDB )
  371. *
  372. * Solve LS or overdetermined system
  373. *
  374. IF( M.GT.0 .AND. N.GT.0 ) THEN
  375. CALL ZLACPY( 'Full', M, N, COPYA, LDA,
  376. $ A, LDA )
  377. CALL ZLACPY( 'Full', NROWS, NRHS,
  378. $ COPYB, LDB, B, LDB )
  379. END IF
  380. SRNAMT = 'ZGELS '
  381. CALL ZGELS( TRANS, M, N, NRHS, A, LDA, B,
  382. $ LDB, WORK, LWORK, INFO )
  383. *
  384. IF( INFO.NE.0 )
  385. $ CALL ALAERH( PATH, 'ZGELS ', INFO, 0,
  386. $ TRANS, M, N, NRHS, -1, NB,
  387. $ ITYPE, NFAIL, NERRS,
  388. $ NOUT )
  389. *
  390. * Check correctness of results
  391. *
  392. LDWORK = MAX( 1, NROWS )
  393. IF( NROWS.GT.0 .AND. NRHS.GT.0 )
  394. $ CALL ZLACPY( 'Full', NROWS, NRHS,
  395. $ COPYB, LDB, C, LDB )
  396. CALL ZQRT16( TRANS, M, N, NRHS, COPYA,
  397. $ LDA, B, LDB, C, LDB, RWORK,
  398. $ RESULT( 1 ) )
  399. *
  400. IF( ( ITRAN.EQ.1 .AND. M.GE.N ) .OR.
  401. $ ( ITRAN.EQ.2 .AND. M.LT.N ) ) THEN
  402. *
  403. * Solving LS system
  404. *
  405. RESULT( 2 ) = ZQRT17( TRANS, 1, M, N,
  406. $ NRHS, COPYA, LDA, B, LDB,
  407. $ COPYB, LDB, C, WORK,
  408. $ LWORK )
  409. ELSE
  410. *
  411. * Solving overdetermined system
  412. *
  413. RESULT( 2 ) = ZQRT14( TRANS, M, N,
  414. $ NRHS, COPYA, LDA, B, LDB,
  415. $ WORK, LWORK )
  416. END IF
  417. *
  418. * Print information about the tests that
  419. * did not pass the threshold.
  420. *
  421. DO 20 K = 1, 2
  422. IF( RESULT( K ).GE.THRESH ) THEN
  423. IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
  424. $ CALL ALAHD( NOUT, PATH )
  425. WRITE( NOUT, FMT = 9999 )TRANS, M,
  426. $ N, NRHS, NB, ITYPE, K,
  427. $ RESULT( K )
  428. NFAIL = NFAIL + 1
  429. END IF
  430. 20 CONTINUE
  431. NRUN = NRUN + 2
  432. 30 CONTINUE
  433. 40 CONTINUE
  434. END IF
  435. *
  436. * Generate a matrix of scaling type ISCALE and rank
  437. * type IRANK.
  438. *
  439. CALL ZQRT15( ISCALE, IRANK, M, N, NRHS, COPYA, LDA,
  440. $ COPYB, LDB, COPYS, RANK, NORMA, NORMB,
  441. $ ISEED, WORK, LWORK )
  442. *
  443. * workspace used: MAX(M+MIN(M,N),NRHS*MIN(M,N),2*N+M)
  444. *
  445. LDWORK = MAX( 1, M )
  446. *
  447. * Loop for testing different block sizes.
  448. *
  449. DO 90 INB = 1, NNB
  450. NB = NBVAL( INB )
  451. CALL XLAENV( 1, NB )
  452. CALL XLAENV( 3, NXVAL( INB ) )
  453. *
  454. * Test ZGELSY
  455. *
  456. * ZGELSY: Compute the minimum-norm solution
  457. * X to min( norm( A * X - B ) )
  458. * using the rank-revealing orthogonal
  459. * factorization.
  460. *
  461. CALL ZLACPY( 'Full', M, N, COPYA, LDA, A, LDA )
  462. CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, B,
  463. $ LDB )
  464. *
  465. * Initialize vector IWORK.
  466. *
  467. DO 70 J = 1, N
  468. IWORK( J ) = 0
  469. 70 CONTINUE
  470. *
  471. * Set LWLSY to the adequate value.
  472. *
  473. LWLSY = MNMIN + MAX( 2*MNMIN, NB*( N+1 ),
  474. $ MNMIN+NB*NRHS )
  475. LWLSY = MAX( 1, LWLSY )
  476. *
  477. SRNAMT = 'ZGELSY'
  478. CALL ZGELSY( M, N, NRHS, A, LDA, B, LDB, IWORK,
  479. $ RCOND, CRANK, WORK, LWLSY, RWORK,
  480. $ INFO )
  481. IF( INFO.NE.0 )
  482. $ CALL ALAERH( PATH, 'ZGELSY', INFO, 0, ' ', M,
  483. $ N, NRHS, -1, NB, ITYPE, NFAIL,
  484. $ NERRS, NOUT )
  485. *
  486. * workspace used: 2*MNMIN+NB*NB+NB*MAX(N,NRHS)
  487. *
  488. * Test 3: Compute relative error in svd
  489. * workspace: M*N + 4*MIN(M,N) + MAX(M,N)
  490. *
  491. RESULT( 3 ) = ZQRT12( CRANK, CRANK, A, LDA,
  492. $ COPYS, WORK, LWORK, RWORK )
  493. *
  494. * Test 4: Compute error in solution
  495. * workspace: M*NRHS + M
  496. *
  497. CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, WORK,
  498. $ LDWORK )
  499. CALL ZQRT16( 'No transpose', M, N, NRHS, COPYA,
  500. $ LDA, B, LDB, WORK, LDWORK, RWORK,
  501. $ RESULT( 4 ) )
  502. *
  503. * Test 5: Check norm of r'*A
  504. * workspace: NRHS*(M+N)
  505. *
  506. RESULT( 5 ) = ZERO
  507. IF( M.GT.CRANK )
  508. $ RESULT( 5 ) = ZQRT17( 'No transpose', 1, M,
  509. $ N, NRHS, COPYA, LDA, B, LDB,
  510. $ COPYB, LDB, C, WORK, LWORK )
  511. *
  512. * Test 6: Check if x is in the rowspace of A
  513. * workspace: (M+NRHS)*(N+2)
  514. *
  515. RESULT( 6 ) = ZERO
  516. *
  517. IF( N.GT.CRANK )
  518. $ RESULT( 6 ) = ZQRT14( 'No transpose', M, N,
  519. $ NRHS, COPYA, LDA, B, LDB,
  520. $ WORK, LWORK )
  521. *
  522. * Test ZGELSS
  523. *
  524. * ZGELSS: Compute the minimum-norm solution
  525. * X to min( norm( A * X - B ) )
  526. * using the SVD.
  527. *
  528. CALL ZLACPY( 'Full', M, N, COPYA, LDA, A, LDA )
  529. CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, B,
  530. $ LDB )
  531. SRNAMT = 'ZGELSS'
  532. CALL ZGELSS( M, N, NRHS, A, LDA, B, LDB, S,
  533. $ RCOND, CRANK, WORK, LWORK, RWORK,
  534. $ INFO )
  535. *
  536. IF( INFO.NE.0 )
  537. $ CALL ALAERH( PATH, 'ZGELSS', INFO, 0, ' ', M,
  538. $ N, NRHS, -1, NB, ITYPE, NFAIL,
  539. $ NERRS, NOUT )
  540. *
  541. * workspace used: 3*min(m,n) +
  542. * max(2*min(m,n),nrhs,max(m,n))
  543. *
  544. * Test 7: Compute relative error in svd
  545. *
  546. IF( RANK.GT.0 ) THEN
  547. CALL DAXPY( MNMIN, -ONE, COPYS, 1, S, 1 )
  548. RESULT( 7 ) = DASUM( MNMIN, S, 1 ) /
  549. $ DASUM( MNMIN, COPYS, 1 ) /
  550. $ ( EPS*DBLE( MNMIN ) )
  551. ELSE
  552. RESULT( 7 ) = ZERO
  553. END IF
  554. *
  555. * Test 8: Compute error in solution
  556. *
  557. CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, WORK,
  558. $ LDWORK )
  559. CALL ZQRT16( 'No transpose', M, N, NRHS, COPYA,
  560. $ LDA, B, LDB, WORK, LDWORK, RWORK,
  561. $ RESULT( 8 ) )
  562. *
  563. * Test 9: Check norm of r'*A
  564. *
  565. RESULT( 9 ) = ZERO
  566. IF( M.GT.CRANK )
  567. $ RESULT( 9 ) = ZQRT17( 'No transpose', 1, M,
  568. $ N, NRHS, COPYA, LDA, B, LDB,
  569. $ COPYB, LDB, C, WORK, LWORK )
  570. *
  571. * Test 10: Check if x is in the rowspace of A
  572. *
  573. RESULT( 10 ) = ZERO
  574. IF( N.GT.CRANK )
  575. $ RESULT( 10 ) = ZQRT14( 'No transpose', M, N,
  576. $ NRHS, COPYA, LDA, B, LDB,
  577. $ WORK, LWORK )
  578. *
  579. * Test ZGELSD
  580. *
  581. * ZGELSD: Compute the minimum-norm solution X
  582. * to min( norm( A * X - B ) ) using a
  583. * divide and conquer SVD.
  584. *
  585. CALL XLAENV( 9, 25 )
  586. *
  587. CALL ZLACPY( 'Full', M, N, COPYA, LDA, A, LDA )
  588. CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, B,
  589. $ LDB )
  590. *
  591. SRNAMT = 'ZGELSD'
  592. CALL ZGELSD( M, N, NRHS, A, LDA, B, LDB, S,
  593. $ RCOND, CRANK, WORK, LWORK, RWORK,
  594. $ IWORK, INFO )
  595. IF( INFO.NE.0 )
  596. $ CALL ALAERH( PATH, 'ZGELSD', INFO, 0, ' ', M,
  597. $ N, NRHS, -1, NB, ITYPE, NFAIL,
  598. $ NERRS, NOUT )
  599. *
  600. * Test 11: Compute relative error in svd
  601. *
  602. IF( RANK.GT.0 ) THEN
  603. CALL DAXPY( MNMIN, -ONE, COPYS, 1, S, 1 )
  604. RESULT( 11 ) = DASUM( MNMIN, S, 1 ) /
  605. $ DASUM( MNMIN, COPYS, 1 ) /
  606. $ ( EPS*DBLE( MNMIN ) )
  607. ELSE
  608. RESULT( 11 ) = ZERO
  609. END IF
  610. *
  611. * Test 12: Compute error in solution
  612. *
  613. CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, WORK,
  614. $ LDWORK )
  615. CALL ZQRT16( 'No transpose', M, N, NRHS, COPYA,
  616. $ LDA, B, LDB, WORK, LDWORK, RWORK,
  617. $ RESULT( 12 ) )
  618. *
  619. * Test 13: Check norm of r'*A
  620. *
  621. RESULT( 13 ) = ZERO
  622. IF( M.GT.CRANK )
  623. $ RESULT( 13 ) = ZQRT17( 'No transpose', 1, M,
  624. $ N, NRHS, COPYA, LDA, B, LDB,
  625. $ COPYB, LDB, C, WORK, LWORK )
  626. *
  627. * Test 14: Check if x is in the rowspace of A
  628. *
  629. RESULT( 14 ) = ZERO
  630. IF( N.GT.CRANK )
  631. $ RESULT( 14 ) = ZQRT14( 'No transpose', M, N,
  632. $ NRHS, COPYA, LDA, B, LDB,
  633. $ WORK, LWORK )
  634. *
  635. * Print information about the tests that did not
  636. * pass the threshold.
  637. *
  638. DO 80 K = 3, NTESTS
  639. IF( RESULT( K ).GE.THRESH ) THEN
  640. IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
  641. $ CALL ALAHD( NOUT, PATH )
  642. WRITE( NOUT, FMT = 9998 )M, N, NRHS, NB,
  643. $ ITYPE, K, RESULT( K )
  644. NFAIL = NFAIL + 1
  645. END IF
  646. 80 CONTINUE
  647. NRUN = NRUN + 12
  648. *
  649. 90 CONTINUE
  650. 100 CONTINUE
  651. 110 CONTINUE
  652. 120 CONTINUE
  653. 130 CONTINUE
  654. 140 CONTINUE
  655. *
  656. * Print a summary of the results.
  657. *
  658. CALL ALASVM( PATH, NOUT, NFAIL, NRUN, NERRS )
  659. *
  660. 9999 FORMAT( ' TRANS=''', A1, ''', M=', I5, ', N=', I5, ', NRHS=', I4,
  661. $ ', NB=', I4, ', type', I2, ', test(', I2, ')=', G12.5 )
  662. 9998 FORMAT( ' M=', I5, ', N=', I5, ', NRHS=', I4, ', NB=', I4,
  663. $ ', type', I2, ', test(', I2, ')=', G12.5 )
  664. RETURN
  665. *
  666. * End of ZDRVLS
  667. *
  668. END