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.

schkqp3rk.f 29 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796
  1. *> \brief \b SCHKQP3RK
  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 SCHKQP3RK( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL,
  12. * $ NNB, NBVAL, NXVAL, THRESH, A, COPYA,
  13. * $ B, COPYB, S, TAU,
  14. * $ WORK, IWORK, NOUT )
  15. * IMPLICIT NONE
  16. *
  17. * -- LAPACK test routine --
  18. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  19. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  20. *
  21. * .. Scalar Arguments ..
  22. * INTEGER NM, NN, NNS, NNB, NOUT
  23. * REAL THRESH
  24. * ..
  25. * .. Array Arguments ..
  26. * LOGICAL DOTYPE( * )
  27. * INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
  28. * $ NVAL( * ), NXVAL( * )
  29. * REAL A( * ), COPYA( * ), B( * ), COPYB( * ),
  30. * $ S( * ), TAU( * ), WORK( * )
  31. * ..
  32. *
  33. *
  34. *> \par Purpose:
  35. * =============
  36. *>
  37. *> \verbatim
  38. *>
  39. *> SCHKQP3RK tests SGEQP3RK.
  40. *> \endverbatim
  41. *
  42. * Arguments:
  43. * ==========
  44. *
  45. *> \param[in] DOTYPE
  46. *> \verbatim
  47. *> DOTYPE is LOGICAL array, dimension (NTYPES)
  48. *> The matrix types to be used for testing. Matrices of type j
  49. *> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
  50. *> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
  51. *> \endverbatim
  52. *>
  53. *> \param[in] NM
  54. *> \verbatim
  55. *> NM is INTEGER
  56. *> The number of values of M contained in the vector MVAL.
  57. *> \endverbatim
  58. *>
  59. *> \param[in] MVAL
  60. *> \verbatim
  61. *> MVAL is INTEGER array, dimension (NM)
  62. *> The values of the matrix row dimension M.
  63. *> \endverbatim
  64. *>
  65. *> \param[in] NN
  66. *> \verbatim
  67. *> NN is INTEGER
  68. *> The number of values of N contained in the vector NVAL.
  69. *> \endverbatim
  70. *>
  71. *> \param[in] NVAL
  72. *> \verbatim
  73. *> NVAL is INTEGER array, dimension (NN)
  74. *> The values of the matrix column dimension N.
  75. *> \endverbatim
  76. *>
  77. *> \param[in] NNS
  78. *> \verbatim
  79. *> NNS is INTEGER
  80. *> The number of values of NRHS contained in the vector NSVAL.
  81. *> \endverbatim
  82. *>
  83. *> \param[in] NSVAL
  84. *> \verbatim
  85. *> NSVAL is INTEGER array, dimension (NNS)
  86. *> The values of the number of right hand sides NRHS.
  87. *> \endverbatim
  88. *>
  89. *> \param[in] NNB
  90. *> \verbatim
  91. *> NNB is INTEGER
  92. *> The number of values of NB and NX contained in the
  93. *> vectors NBVAL and NXVAL. The blocking parameters are used
  94. *> in pairs (NB,NX).
  95. *> \endverbatim
  96. *>
  97. *> \param[in] NBVAL
  98. *> \verbatim
  99. *> NBVAL is INTEGER array, dimension (NNB)
  100. *> The values of the blocksize NB.
  101. *> \endverbatim
  102. *>
  103. *> \param[in] NXVAL
  104. *> \verbatim
  105. *> NXVAL is INTEGER array, dimension (NNB)
  106. *> The values of the crossover point NX.
  107. *> \endverbatim
  108. *>
  109. *> \param[in] THRESH
  110. *> \verbatim
  111. *> THRESH is REAL
  112. *> The threshold value for the test ratios. A result is
  113. *> included in the output file if RESULT >= THRESH. To have
  114. *> every test ratio printed, use THRESH = 0.
  115. *> \endverbatim
  116. *>
  117. *> \param[out] A
  118. *> \verbatim
  119. *> A is REAL array, dimension (MMAX*NMAX)
  120. *> where MMAX is the maximum value of M in MVAL and NMAX is the
  121. *> maximum value of N in NVAL.
  122. *> \endverbatim
  123. *>
  124. *> \param[out] COPYA
  125. *> \verbatim
  126. *> COPYA is REAL array, dimension (MMAX*NMAX)
  127. *> \endverbatim
  128. *>
  129. *> \param[out] B
  130. *> \verbatim
  131. *> B is REAL array, dimension (MMAX*NSMAX)
  132. *> where MMAX is the maximum value of M in MVAL and NSMAX is the
  133. *> maximum value of NRHS in NSVAL.
  134. *> \endverbatim
  135. *>
  136. *> \param[out] COPYB
  137. *> \verbatim
  138. *> COPYB is REAL array, dimension (MMAX*NSMAX)
  139. *> \endverbatim
  140. *>
  141. *> \param[out] S
  142. *> \verbatim
  143. *> S is REAL array, dimension
  144. *> (min(MMAX,NMAX))
  145. *> \endverbatim
  146. *>
  147. *> \param[out] TAU
  148. *> \verbatim
  149. *> TAU is REAL array, dimension (MMAX)
  150. *> \endverbatim
  151. *>
  152. *> \param[out] WORK
  153. *> \verbatim
  154. *> WORK is REAL array, dimension
  155. *> (MMAX*NMAX + 4*NMAX + MMAX)
  156. *> \endverbatim
  157. *>
  158. *> \param[out] IWORK
  159. *> \verbatim
  160. *> IWORK is INTEGER array, dimension (2*NMAX)
  161. *> \endverbatim
  162. *>
  163. *> \param[in] NOUT
  164. *> \verbatim
  165. *> NOUT is INTEGER
  166. *> The unit number for output.
  167. *> \endverbatim
  168. *
  169. * Authors:
  170. * ========
  171. *
  172. *> \author Univ. of Tennessee
  173. *> \author Univ. of California Berkeley
  174. *> \author Univ. of Colorado Denver
  175. *> \author NAG Ltd.
  176. *
  177. *> \ingroup single_lin
  178. *
  179. * =====================================================================
  180. SUBROUTINE SCHKQP3RK( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL,
  181. $ NNB, NBVAL, NXVAL, THRESH, A, COPYA,
  182. $ B, COPYB, S, TAU,
  183. $ WORK, IWORK, NOUT )
  184. IMPLICIT NONE
  185. *
  186. * -- LAPACK test routine --
  187. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  188. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  189. *
  190. * .. Scalar Arguments ..
  191. INTEGER NM, NN, NNB, NNS, NOUT
  192. REAL THRESH
  193. * ..
  194. * .. Array Arguments ..
  195. LOGICAL DOTYPE( * )
  196. INTEGER IWORK( * ), NBVAL( * ), MVAL( * ), NVAL( * ),
  197. $ NSVAL( * ), NXVAL( * )
  198. REAL A( * ), COPYA( * ), B( * ), COPYB( * ),
  199. $ S( * ), TAU( * ), WORK( * )
  200. * ..
  201. *
  202. * =====================================================================
  203. *
  204. * .. Parameters ..
  205. INTEGER NTYPES
  206. PARAMETER ( NTYPES = 19 )
  207. INTEGER NTESTS
  208. PARAMETER ( NTESTS = 5 )
  209. REAL ONE, ZERO, BIGNUM
  210. PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0,
  211. $ BIGNUM = 1.0E+38 )
  212. * ..
  213. * .. Local Scalars ..
  214. CHARACTER DIST, TYPE
  215. CHARACTER*3 PATH
  216. INTEGER I, IHIGH, ILOW, IM, IMAT, IN, INC_ZERO,
  217. $ INB, IND_OFFSET_GEN,
  218. $ IND_IN, IND_OUT, INS, INFO,
  219. $ ISTEP, J, J_INC, J_FIRST_NZ, JB_ZERO,
  220. $ KFACT, KL, KMAX, KU, LDA, LW, LWORK,
  221. $ LWORK_MQR, M, MINMN, MINMNB_GEN, MODE, N,
  222. $ NB, NB_ZERO, NERRS, NFAIL, NB_GEN, NRHS,
  223. $ NRUN, NX, T
  224. REAL ANORM, CNDNUM, EPS, ABSTOL, RELTOL,
  225. $ DTEMP, MAXC2NRMK, RELMAXC2NRMK
  226. * ..
  227. * .. Local Arrays ..
  228. INTEGER ISEED( 4 ), ISEEDY( 4 )
  229. REAL RESULT( NTESTS ), RDUMMY( 1 )
  230. * ..
  231. * .. External Functions ..
  232. REAL SLAMCH, SQPT01, SQRT11, SQRT12, SLANGE
  233. EXTERNAL SLAMCH, SQPT01, SQRT11, SQRT12, SLANGE
  234. * ..
  235. * .. External Subroutines ..
  236. EXTERNAL ALAERH, ALAHD, ALASUM, SAXPY, SGEQP3RK,
  237. $ SLACPY, SLAORD, SLASET, SLATB4, SLATMS,
  238. $ SORMQR, SSWAP, ICOPY, XLAENV
  239. * ..
  240. * .. Intrinsic Functions ..
  241. INTRINSIC ABS, MAX, MIN, MOD, REAL
  242. * ..
  243. * .. Scalars in Common ..
  244. LOGICAL LERR, OK
  245. CHARACTER*32 SRNAMT
  246. INTEGER INFOT, IOUNIT
  247. * ..
  248. * .. Common blocks ..
  249. COMMON / INFOC / INFOT, IOUNIT, OK, LERR
  250. COMMON / SRNAMC / SRNAMT
  251. * ..
  252. * .. Data statements ..
  253. DATA ISEEDY / 1988, 1989, 1990, 1991 /
  254. * ..
  255. * .. Executable Statements ..
  256. *
  257. * Initialize constants and the random number seed.
  258. *
  259. PATH( 1: 1 ) = 'Single precision'
  260. PATH( 2: 3 ) = 'QK'
  261. NRUN = 0
  262. NFAIL = 0
  263. NERRS = 0
  264. DO I = 1, 4
  265. ISEED( I ) = ISEEDY( I )
  266. END DO
  267. EPS = SLAMCH( 'Epsilon' )
  268. INFOT = 0
  269. *
  270. DO IM = 1, NM
  271. *
  272. * Do for each value of M in MVAL.
  273. *
  274. M = MVAL( IM )
  275. LDA = MAX( 1, M )
  276. *
  277. DO IN = 1, NN
  278. *
  279. * Do for each value of N in NVAL.
  280. *
  281. N = NVAL( IN )
  282. MINMN = MIN( M, N )
  283. LWORK = MAX( 1, M*MAX( M, N )+4*MINMN+MAX( M, N ),
  284. $ M*N + 2*MINMN + 4*N )
  285. *
  286. DO INS = 1, NNS
  287. NRHS = NSVAL( INS )
  288. *
  289. * Set up parameters with SLATB4 and generate
  290. * M-by-NRHS B matrix with SLATMS.
  291. * IMAT = 14:
  292. * Random matrix, CNDNUM = 2, NORM = ONE,
  293. * MODE = 3 (geometric distribution of singular values).
  294. *
  295. CALL SLATB4( PATH, 14, M, NRHS, TYPE, KL, KU, ANORM,
  296. $ MODE, CNDNUM, DIST )
  297. *
  298. SRNAMT = 'SLATMS'
  299. CALL SLATMS( M, NRHS, DIST, ISEED, TYPE, S, MODE,
  300. $ CNDNUM, ANORM, KL, KU, 'No packing',
  301. $ COPYB, LDA, WORK, INFO )
  302. *
  303. * Check error code from SLATMS.
  304. *
  305. IF( INFO.NE.0 ) THEN
  306. CALL ALAERH( PATH, 'SLATMS', INFO, 0, ' ', M,
  307. $ NRHS, -1, -1, -1, 6, NFAIL, NERRS,
  308. $ NOUT )
  309. CYCLE
  310. END IF
  311. *
  312. DO IMAT = 1, NTYPES
  313. *
  314. * Do the tests only if DOTYPE( IMAT ) is true.
  315. *
  316. IF( .NOT.DOTYPE( IMAT ) )
  317. $ CYCLE
  318. *
  319. * The type of distribution used to generate the random
  320. * eigen-/singular values:
  321. * ( 'S' for symmetric distribution ) => UNIFORM( -1, 1 )
  322. *
  323. * Do for each type of NON-SYMMETRIC matrix: CNDNUM NORM MODE
  324. * 1. Zero matrix
  325. * 2. Random, Diagonal, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  326. * 3. Random, Upper triangular, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  327. * 4. Random, Lower triangular, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  328. * 5. Random, First column is zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  329. * 6. Random, Last MINMN column is zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  330. * 7. Random, Last N column is zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  331. * 8. Random, Middle column in MINMN is zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  332. * 9. Random, First half of MINMN columns are zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  333. * 10. Random, Last columns are zero starting from MINMN/2+1, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  334. * 11. Random, Half MINMN columns in the middle are zero starting
  335. * from MINMN/2-(MINMN/2)/2+1, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  336. * 12. Random, Odd columns are ZERO, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  337. * 13. Random, Even columns are ZERO, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  338. * 14. Random, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
  339. * 15. Random, CNDNUM = sqrt(0.1/EPS) CNDNUM = BADC1 = sqrt(0.1/EPS) ONE 3 ( geometric distribution of singular values )
  340. * 16. Random, CNDNUM = 0.1/EPS CNDNUM = BADC2 = 0.1/EPS ONE 3 ( geometric distribution of singular values )
  341. * 17. Random, CNDNUM = 0.1/EPS, CNDNUM = BADC2 = 0.1/EPS ONE 2 ( one small singular value, S(N)=1/CNDNUM )
  342. * one small singular value S(N)=1/CNDNUM
  343. * 18. Random, CNDNUM = 2, scaled near underflow CNDNUM = 2 SMALL = SAFMIN
  344. * 19. Random, CNDNUM = 2, scaled near overflow CNDNUM = 2 LARGE = 1.0/( 0.25 * ( SAFMIN / EPS ) ) 3 ( geometric distribution of singular values )
  345. *
  346. IF( IMAT.EQ.1 ) THEN
  347. *
  348. * Matrix 1: Zero matrix
  349. *
  350. CALL SLASET( 'Full', M, N, ZERO, ZERO, COPYA, LDA )
  351. DO I = 1, MINMN
  352. S( I ) = ZERO
  353. END DO
  354. *
  355. ELSE IF( (IMAT.GE.2 .AND. IMAT.LE.4 )
  356. $ .OR. (IMAT.GE.14 .AND. IMAT.LE.19 ) ) THEN
  357. *
  358. * Matrices 2-5.
  359. *
  360. * Set up parameters with SLATB4 and generate a test
  361. * matrix with SLATMS.
  362. *
  363. CALL SLATB4( PATH, IMAT, M, N, TYPE, KL, KU, ANORM,
  364. $ MODE, CNDNUM, DIST )
  365. *
  366. SRNAMT = 'SLATMS'
  367. CALL SLATMS( M, N, DIST, ISEED, TYPE, S, MODE,
  368. $ CNDNUM, ANORM, KL, KU, 'No packing',
  369. $ COPYA, LDA, WORK, INFO )
  370. *
  371. * Check error code from SLATMS.
  372. *
  373. IF( INFO.NE.0 ) THEN
  374. CALL ALAERH( PATH, 'SLATMS', INFO, 0, ' ', M, N,
  375. $ -1, -1, -1, IMAT, NFAIL, NERRS,
  376. $ NOUT )
  377. CYCLE
  378. END IF
  379. *
  380. CALL SLAORD( 'Decreasing', MINMN, S, 1 )
  381. *
  382. ELSE IF( MINMN.GE.2
  383. $ .AND. IMAT.GE.5 .AND. IMAT.LE.13 ) THEN
  384. *
  385. * Rectangular matrices 5-13 that contain zero columns,
  386. * only for matrices MINMN >=2.
  387. *
  388. * JB_ZERO is the column index of ZERO block.
  389. * NB_ZERO is the column block size of ZERO block.
  390. * NB_GEN is the column blcok size of the
  391. * generated block.
  392. * J_INC in the non_zero column index increment
  393. * for matrix 12 and 13.
  394. * J_FIRS_NZ is the index of the first non-zero
  395. * column.
  396. *
  397. IF( IMAT.EQ.5 ) THEN
  398. *
  399. * First column is zero.
  400. *
  401. JB_ZERO = 1
  402. NB_ZERO = 1
  403. NB_GEN = N - NB_ZERO
  404. *
  405. ELSE IF( IMAT.EQ.6 ) THEN
  406. *
  407. * Last column MINMN is zero.
  408. *
  409. JB_ZERO = MINMN
  410. NB_ZERO = 1
  411. NB_GEN = N - NB_ZERO
  412. *
  413. ELSE IF( IMAT.EQ.7 ) THEN
  414. *
  415. * Last column N is zero.
  416. *
  417. JB_ZERO = N
  418. NB_ZERO = 1
  419. NB_GEN = N - NB_ZERO
  420. *
  421. ELSE IF( IMAT.EQ.8 ) THEN
  422. *
  423. * Middle column in MINMN is zero.
  424. *
  425. JB_ZERO = MINMN / 2 + 1
  426. NB_ZERO = 1
  427. NB_GEN = N - NB_ZERO
  428. *
  429. ELSE IF( IMAT.EQ.9 ) THEN
  430. *
  431. * First half of MINMN columns is zero.
  432. *
  433. JB_ZERO = 1
  434. NB_ZERO = MINMN / 2
  435. NB_GEN = N - NB_ZERO
  436. *
  437. ELSE IF( IMAT.EQ.10 ) THEN
  438. *
  439. * Last columns are zero columns,
  440. * starting from (MINMN / 2 + 1) column.
  441. *
  442. JB_ZERO = MINMN / 2 + 1
  443. NB_ZERO = N - JB_ZERO + 1
  444. NB_GEN = N - NB_ZERO
  445. *
  446. ELSE IF( IMAT.EQ.11 ) THEN
  447. *
  448. * Half of the columns in the middle of MINMN
  449. * columns is zero, starting from
  450. * MINMN/2 - (MINMN/2)/2 + 1 column.
  451. *
  452. JB_ZERO = MINMN / 2 - (MINMN / 2) / 2 + 1
  453. NB_ZERO = MINMN / 2
  454. NB_GEN = N - NB_ZERO
  455. *
  456. ELSE IF( IMAT.EQ.12 ) THEN
  457. *
  458. * Odd-numbered columns are zero,
  459. *
  460. NB_GEN = N / 2
  461. NB_ZERO = N - NB_GEN
  462. J_INC = 2
  463. J_FIRST_NZ = 2
  464. *
  465. ELSE IF( IMAT.EQ.13 ) THEN
  466. *
  467. * Even-numbered columns are zero.
  468. *
  469. NB_ZERO = N / 2
  470. NB_GEN = N - NB_ZERO
  471. J_INC = 2
  472. J_FIRST_NZ = 1
  473. *
  474. END IF
  475. *
  476. *
  477. * 1) Set the first NB_ZERO columns in COPYA(1:M,1:N)
  478. * to zero.
  479. *
  480. CALL SLASET( 'Full', M, NB_ZERO, ZERO, ZERO,
  481. $ COPYA, LDA )
  482. *
  483. * 2) Generate an M-by-(N-NB_ZERO) matrix with the
  484. * chosen singular value distribution
  485. * in COPYA(1:M,NB_ZERO+1:N).
  486. *
  487. CALL SLATB4( PATH, IMAT, M, NB_GEN, TYPE, KL, KU,
  488. $ ANORM, MODE, CNDNUM, DIST )
  489. *
  490. SRNAMT = 'SLATMS'
  491. *
  492. IND_OFFSET_GEN = NB_ZERO * LDA
  493. *
  494. CALL SLATMS( M, NB_GEN, DIST, ISEED, TYPE, S, MODE,
  495. $ CNDNUM, ANORM, KL, KU, 'No packing',
  496. $ COPYA( IND_OFFSET_GEN + 1 ), LDA,
  497. $ WORK, INFO )
  498. *
  499. * Check error code from SLATMS.
  500. *
  501. IF( INFO.NE.0 ) THEN
  502. CALL ALAERH( PATH, 'SLATMS', INFO, 0, ' ', M,
  503. $ NB_GEN, -1, -1, -1, IMAT, NFAIL,
  504. $ NERRS, NOUT )
  505. CYCLE
  506. END IF
  507. *
  508. * 3) Swap the gererated colums from the right side
  509. * NB_GEN-size block in COPYA into correct column
  510. * positions.
  511. *
  512. IF( IMAT.EQ.6
  513. $ .OR. IMAT.EQ.7
  514. $ .OR. IMAT.EQ.8
  515. $ .OR. IMAT.EQ.10
  516. $ .OR. IMAT.EQ.11 ) THEN
  517. *
  518. * Move by swapping the generated columns
  519. * from the right NB_GEN-size block from
  520. * (NB_ZERO+1:NB_ZERO+JB_ZERO)
  521. * into columns (1:JB_ZERO-1).
  522. *
  523. DO J = 1, JB_ZERO-1, 1
  524. CALL SSWAP( M,
  525. $ COPYA( ( NB_ZERO+J-1)*LDA+1), 1,
  526. $ COPYA( (J-1)*LDA + 1 ), 1 )
  527. END DO
  528. *
  529. ELSE IF( IMAT.EQ.12 .OR. IMAT.EQ.13 ) THEN
  530. *
  531. * ( IMAT = 12, Odd-numbered ZERO columns. )
  532. * Swap the generated columns from the right
  533. * NB_GEN-size block into the even zero colums in the
  534. * left NB_ZERO-size block.
  535. *
  536. * ( IMAT = 13, Even-numbered ZERO columns. )
  537. * Swap the generated columns from the right
  538. * NB_GEN-size block into the odd zero colums in the
  539. * left NB_ZERO-size block.
  540. *
  541. DO J = 1, NB_GEN, 1
  542. IND_OUT = ( NB_ZERO+J-1 )*LDA + 1
  543. IND_IN = ( J_INC*(J-1)+(J_FIRST_NZ-1) )*LDA
  544. $ + 1
  545. CALL SSWAP( M,
  546. $ COPYA( IND_OUT ), 1,
  547. $ COPYA( IND_IN), 1 )
  548. END DO
  549. *
  550. END IF
  551. *
  552. * 5) Order the singular values generated by
  553. * DLAMTS in decreasing order and add trailing zeros
  554. * that correspond to zero columns.
  555. * The total number of singular values is MINMN.
  556. *
  557. MINMNB_GEN = MIN( M, NB_GEN )
  558. *
  559. DO I = MINMNB_GEN+1, MINMN
  560. S( I ) = ZERO
  561. END DO
  562. *
  563. ELSE
  564. *
  565. * IF(MINMN.LT.2) skip this size for this matrix type.
  566. *
  567. CYCLE
  568. END IF
  569. *
  570. * Initialize a copy array for a pivot array for SGEQP3RK.
  571. *
  572. DO I = 1, N
  573. IWORK( I ) = 0
  574. END DO
  575. *
  576. DO INB = 1, NNB
  577. *
  578. * Do for each pair of values (NB,NX) in NBVAL and NXVAL.
  579. *
  580. NB = NBVAL( INB )
  581. CALL XLAENV( 1, NB )
  582. NX = NXVAL( INB )
  583. CALL XLAENV( 3, NX )
  584. *
  585. * We do MIN(M,N)+1 because we need a test for KMAX > N,
  586. * when KMAX is larger than MIN(M,N), KMAX should be
  587. * KMAX = MIN(M,N)
  588. *
  589. DO KMAX = 0, MIN(M,N)+1
  590. *
  591. * Get a working copy of COPYA into A( 1:M,1:N ).
  592. * Get a working copy of COPYB into A( 1:M, (N+1):NRHS ).
  593. * Get a working copy of COPYB into into B( 1:M, 1:NRHS ).
  594. * Get a working copy of IWORK(1:N) awith zeroes into
  595. * which is going to be used as pivot array IWORK( N+1:2N ).
  596. * NOTE: IWORK(2N+1:3N) is going to be used as a WORK array
  597. * for the routine.
  598. *
  599. CALL SLACPY( 'All', M, N, COPYA, LDA, A, LDA )
  600. CALL SLACPY( 'All', M, NRHS, COPYB, LDA,
  601. $ A( LDA*N + 1 ), LDA )
  602. CALL SLACPY( 'All', M, NRHS, COPYB, LDA,
  603. $ B, LDA )
  604. CALL ICOPY( N, IWORK( 1 ), 1, IWORK( N+1 ), 1 )
  605. DO I = 1, NTESTS
  606. RESULT( I ) = ZERO
  607. END DO
  608. *
  609. ABSTOL = -1.0
  610. RELTOL = -1.0
  611. *
  612. * Compute the QR factorization with pivoting of A
  613. *
  614. LW = MAX( 1, MAX( 2*N + NB*( N+NRHS+1 ),
  615. $ 3*N + NRHS - 1 ) )
  616. *
  617. * Compute SGEQP3RK factorization of A.
  618. *
  619. SRNAMT = 'SGEQP3RK'
  620. CALL SGEQP3RK( M, N, NRHS, KMAX, ABSTOL, RELTOL,
  621. $ A, LDA, KFACT, MAXC2NRMK,
  622. $ RELMAXC2NRMK, IWORK( N+1 ), TAU,
  623. $ WORK, LW, IWORK( 2*N+1 ), INFO )
  624. *
  625. * Check error code from SGEQP3RK.
  626. *
  627. IF( INFO.LT.0 )
  628. $ CALL ALAERH( PATH, 'SGEQP3RK', INFO, 0, ' ',
  629. $ M, N, NX, -1, NB, IMAT,
  630. $ NFAIL, NERRS, NOUT )
  631. *
  632. * Compute test 1:
  633. *
  634. * This test in only for the full rank factorization of
  635. * the matrix A.
  636. *
  637. * Array S(1:min(M,N)) contains svd(A) the sigular values
  638. * of the original matrix A in decreasing absolute value
  639. * order. The test computes svd(R), the vector sigular
  640. * values of the upper trapezoid of A(1:M,1:N) that
  641. * contains the factor R, in decreasing order. The test
  642. * returns the ratio:
  643. *
  644. * 2-norm(svd(R) - svd(A)) / ( max(M,N) * 2-norm(svd(A)) * EPS )
  645. *
  646. IF( KFACT.EQ.MINMN ) THEN
  647. *
  648. RESULT( 1 ) = SQRT12( M, N, A, LDA, S, WORK,
  649. $ LWORK )
  650. *
  651. NRUN = NRUN + 1
  652. *
  653. * End test 1
  654. *
  655. END IF
  656. *
  657. * Compute test 2:
  658. *
  659. * The test returns the ratio:
  660. *
  661. * 1-norm( A*P - Q*R ) / ( max(M,N) * 1-norm(A) * EPS )
  662. *
  663. RESULT( 2 ) = SQPT01( M, N, KFACT, COPYA, A, LDA, TAU,
  664. $ IWORK( N+1 ), WORK, LWORK )
  665. *
  666. * Compute test 3:
  667. *
  668. * The test returns the ratio:
  669. *
  670. * 1-norm( Q**T * Q - I ) / ( M * EPS )
  671. *
  672. RESULT( 3 ) = SQRT11( M, KFACT, A, LDA, TAU, WORK,
  673. $ LWORK )
  674. *
  675. NRUN = NRUN + 2
  676. *
  677. * Compute test 4:
  678. *
  679. * This test is only for the factorizations with the
  680. * rank greater than 2.
  681. * The elements on the diagonal of R should be non-
  682. * increasing.
  683. *
  684. * The test returns the ratio:
  685. *
  686. * Returns 1.0D+100 if abs(R(K+1,K+1)) > abs(R(K,K)),
  687. * K=1:KFACT-1
  688. *
  689. IF( MIN(KFACT, MINMN).GE.2 ) THEN
  690. *
  691. DO J = 1, KFACT-1, 1
  692. DTEMP = (( ABS( A( (J-1)*LDA+J ) ) -
  693. $ ABS( A( (J)*LDA+J+1 ) ) ) /
  694. $ ABS( A(1) ) )
  695. *
  696. IF( DTEMP.LT.ZERO ) THEN
  697. RESULT( 4 ) = BIGNUM
  698. END IF
  699. *
  700. END DO
  701. *
  702. NRUN = NRUN + 1
  703. *
  704. * End test 4.
  705. *
  706. END IF
  707. *
  708. * Compute test 5:
  709. *
  710. * This test in only for matrix A with min(M,N) > 0.
  711. *
  712. * The test returns the ratio:
  713. *
  714. * 1-norm(Q**T * B - Q**T * B ) /
  715. * ( M * EPS )
  716. *
  717. * (1) Compute B:=Q**T * B in the matrix B.
  718. *
  719. IF( MINMN.GT.0 ) THEN
  720. *
  721. LWORK_MQR = MAX(1, NRHS)
  722. CALL SORMQR( 'Left', 'Transpose',
  723. $ M, NRHS, KFACT, A, LDA, TAU, B, LDA,
  724. $ WORK, LWORK_MQR, INFO )
  725. *
  726. DO I = 1, NRHS
  727. *
  728. * Compare N+J-th column of A and J-column of B.
  729. *
  730. CALL SAXPY( M, -ONE, A( ( N+I-1 )*LDA+1 ), 1,
  731. $ B( ( I-1 )*LDA+1 ), 1 )
  732. END DO
  733. *
  734. RESULT( 5 ) = ABS(
  735. $ SLANGE( 'One-norm', M, NRHS, B, LDA, RDUMMY ) /
  736. $ ( REAL( M )*SLAMCH( 'Epsilon' ) ) )
  737. *
  738. NRUN = NRUN + 1
  739. *
  740. * End compute test 5.
  741. *
  742. END IF
  743. *
  744. * Print information about the tests that did not pass
  745. * the threshold.
  746. *
  747. DO T = 1, NTESTS
  748. IF( RESULT( T ).GE.THRESH ) THEN
  749. IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
  750. $ CALL ALAHD( NOUT, PATH )
  751. WRITE( NOUT, FMT = 9999 ) 'SGEQP3RK', M, N,
  752. $ NRHS, KMAX, ABSTOL, RELTOL,
  753. $ NB, NX, IMAT, T, RESULT( T )
  754. NFAIL = NFAIL + 1
  755. END IF
  756. END DO
  757. *
  758. * END DO KMAX = 1, MIN(M,N)+1
  759. *
  760. END DO
  761. *
  762. * END DO for INB = 1, NNB
  763. *
  764. END DO
  765. *
  766. * END DO for IMAT = 1, NTYPES
  767. *
  768. END DO
  769. *
  770. * END DO for INS = 1, NNS
  771. *
  772. END DO
  773. *
  774. * END DO for IN = 1, NN
  775. *
  776. END DO
  777. *
  778. * END DO for IM = 1, NM
  779. *
  780. END DO
  781. *
  782. * Print a summary of the results.
  783. *
  784. CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS )
  785. *
  786. 9999 FORMAT( 1X, A, ' M =', I5, ', N =', I5, ', NRHS =', I5,
  787. $ ', KMAX =', I5, ', ABSTOL =', G12.5,
  788. $ ', RELTOL =', G12.5, ', NB =', I4, ', NX =', I4,
  789. $ ', type ', I2, ', test ', I2, ', ratio =', G12.5 )
  790. *
  791. * End of SCHKQP3RK
  792. *
  793. END