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.

ddrvac.f 14 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455
  1. *> \brief \b DDRVAC
  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 DDRVAC( DOTYPE, NM, MVAL, NNS, NSVAL, THRESH, NMAX,
  12. * A, AFAC, B, X, WORK,
  13. * RWORK, SWORK, NOUT )
  14. *
  15. * .. Scalar Arguments ..
  16. * INTEGER NMAX, NM, NNS, NOUT
  17. * DOUBLE PRECISION THRESH
  18. * ..
  19. * .. Array Arguments ..
  20. * LOGICAL DOTYPE( * )
  21. * INTEGER MVAL( * ), NSVAL( * )
  22. * REAL SWORK(*)
  23. * DOUBLE PRECISION A( * ), AFAC( * ), B( * ),
  24. * $ RWORK( * ), WORK( * ), X( * )
  25. * ..
  26. *
  27. *
  28. *> \par Purpose:
  29. * =============
  30. *>
  31. *> \verbatim
  32. *>
  33. *> DDRVAC tests DSPOSV.
  34. *> \endverbatim
  35. *
  36. * Arguments:
  37. * ==========
  38. *
  39. *> \param[in] DOTYPE
  40. *> \verbatim
  41. *> DOTYPE is LOGICAL array, dimension (NTYPES)
  42. *> The matrix types to be used for testing. Matrices of type j
  43. *> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
  44. *> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
  45. *> \endverbatim
  46. *>
  47. *> \param[in] NM
  48. *> \verbatim
  49. *> NM is INTEGER
  50. *> The number of values of N contained in the vector MVAL.
  51. *> \endverbatim
  52. *>
  53. *> \param[in] MVAL
  54. *> \verbatim
  55. *> MVAL is INTEGER array, dimension (NM)
  56. *> The values of the matrix dimension N.
  57. *> \endverbatim
  58. *>
  59. *> \param[in] NNS
  60. *> \verbatim
  61. *> NNS is INTEGER
  62. *> The number of values of NRHS contained in the vector NSVAL.
  63. *> \endverbatim
  64. *>
  65. *> \param[in] NSVAL
  66. *> \verbatim
  67. *> NSVAL is INTEGER array, dimension (NNS)
  68. *> The values of the number of right hand sides NRHS.
  69. *> \endverbatim
  70. *>
  71. *> \param[in] THRESH
  72. *> \verbatim
  73. *> THRESH is DOUBLE PRECISION
  74. *> The threshold value for the test ratios. A result is
  75. *> included in the output file if RESULT >= THRESH. To have
  76. *> every test ratio printed, use THRESH = 0.
  77. *> \endverbatim
  78. *>
  79. *> \param[in] NMAX
  80. *> \verbatim
  81. *> NMAX is INTEGER
  82. *> The maximum value permitted for N, used in dimensioning the
  83. *> work arrays.
  84. *> \endverbatim
  85. *>
  86. *> \param[out] A
  87. *> \verbatim
  88. *> A is DOUBLE PRECISION array, dimension (NMAX*NMAX)
  89. *> \endverbatim
  90. *>
  91. *> \param[out] AFAC
  92. *> \verbatim
  93. *> AFAC is DOUBLE PRECISION array, dimension (NMAX*NMAX)
  94. *> \endverbatim
  95. *>
  96. *> \param[out] B
  97. *> \verbatim
  98. *> B is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
  99. *> \endverbatim
  100. *>
  101. *> \param[out] X
  102. *> \verbatim
  103. *> X is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
  104. *> \endverbatim
  105. *>
  106. *> \param[out] WORK
  107. *> \verbatim
  108. *> WORK is DOUBLE PRECISION array, dimension
  109. *> (NMAX*max(3,NSMAX))
  110. *> \endverbatim
  111. *>
  112. *> \param[out] RWORK
  113. *> \verbatim
  114. *> RWORK is DOUBLE PRECISION array, dimension
  115. *> (max(2*NMAX,2*NSMAX+NWORK))
  116. *> \endverbatim
  117. *>
  118. *> \param[out] SWORK
  119. *> \verbatim
  120. *> SWORK is REAL array, dimension
  121. *> (NMAX*(NSMAX+NMAX))
  122. *> \endverbatim
  123. *>
  124. *> \param[in] NOUT
  125. *> \verbatim
  126. *> NOUT is INTEGER
  127. *> The unit number for output.
  128. *> \endverbatim
  129. *
  130. * Authors:
  131. * ========
  132. *
  133. *> \author Univ. of Tennessee
  134. *> \author Univ. of California Berkeley
  135. *> \author Univ. of Colorado Denver
  136. *> \author NAG Ltd.
  137. *
  138. *> \ingroup double_lin
  139. *
  140. * =====================================================================
  141. SUBROUTINE DDRVAC( DOTYPE, NM, MVAL, NNS, NSVAL, THRESH, NMAX,
  142. $ A, AFAC, B, X, WORK,
  143. $ RWORK, SWORK, NOUT )
  144. *
  145. * -- LAPACK test routine --
  146. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  147. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  148. *
  149. * .. Scalar Arguments ..
  150. INTEGER NMAX, NM, NNS, NOUT
  151. DOUBLE PRECISION THRESH
  152. * ..
  153. * .. Array Arguments ..
  154. LOGICAL DOTYPE( * )
  155. INTEGER MVAL( * ), NSVAL( * )
  156. REAL SWORK(*)
  157. DOUBLE PRECISION A( * ), AFAC( * ), B( * ),
  158. $ RWORK( * ), WORK( * ), X( * )
  159. * ..
  160. *
  161. * =====================================================================
  162. *
  163. * .. Parameters ..
  164. DOUBLE PRECISION ZERO
  165. PARAMETER ( ZERO = 0.0D+0 )
  166. INTEGER NTYPES
  167. PARAMETER ( NTYPES = 9 )
  168. INTEGER NTESTS
  169. PARAMETER ( NTESTS = 1 )
  170. * ..
  171. * .. Local Scalars ..
  172. LOGICAL ZEROT
  173. CHARACTER DIST, TYPE, UPLO, XTYPE
  174. CHARACTER*3 PATH
  175. INTEGER I, IM, IMAT, INFO, IOFF, IRHS, IUPLO,
  176. $ IZERO, KL, KU, LDA, MODE, N,
  177. $ NERRS, NFAIL, NIMAT, NRHS, NRUN
  178. DOUBLE PRECISION ANORM, CNDNUM
  179. * ..
  180. * .. Local Arrays ..
  181. CHARACTER UPLOS( 2 )
  182. INTEGER ISEED( 4 ), ISEEDY( 4 )
  183. DOUBLE PRECISION RESULT( NTESTS )
  184. * ..
  185. * .. Local Variables ..
  186. INTEGER ITER, KASE
  187. * ..
  188. * .. External Functions ..
  189. LOGICAL LSAME
  190. EXTERNAL LSAME
  191. * ..
  192. * .. External Subroutines ..
  193. EXTERNAL ALAERH, DLACPY,
  194. $ DLARHS, DLASET, DLATB4, DLATMS,
  195. $ DPOT06, DSPOSV
  196. * ..
  197. * .. Intrinsic Functions ..
  198. INTRINSIC DBLE, MAX, SQRT
  199. * ..
  200. * .. Scalars in Common ..
  201. LOGICAL LERR, OK
  202. CHARACTER*32 SRNAMT
  203. INTEGER INFOT, NUNIT
  204. * ..
  205. * .. Common blocks ..
  206. COMMON / INFOC / INFOT, NUNIT, OK, LERR
  207. COMMON / SRNAMC / SRNAMT
  208. * ..
  209. * .. Data statements ..
  210. DATA ISEEDY / 1988, 1989, 1990, 1991 /
  211. DATA UPLOS / 'U', 'L' /
  212. * ..
  213. * .. Executable Statements ..
  214. *
  215. * Initialize constants and the random number seed.
  216. *
  217. KASE = 0
  218. PATH( 1: 1 ) = 'Double precision'
  219. PATH( 2: 3 ) = 'PO'
  220. NRUN = 0
  221. NFAIL = 0
  222. NERRS = 0
  223. DO 10 I = 1, 4
  224. ISEED( I ) = ISEEDY( I )
  225. 10 CONTINUE
  226. *
  227. INFOT = 0
  228. *
  229. * Do for each value of N in MVAL
  230. *
  231. DO 120 IM = 1, NM
  232. N = MVAL( IM )
  233. LDA = MAX( N, 1 )
  234. NIMAT = NTYPES
  235. IF( N.LE.0 )
  236. $ NIMAT = 1
  237. *
  238. DO 110 IMAT = 1, NIMAT
  239. *
  240. * Do the tests only if DOTYPE( IMAT ) is true.
  241. *
  242. IF( .NOT.DOTYPE( IMAT ) )
  243. $ GO TO 110
  244. *
  245. * Skip types 3, 4, or 5 if the matrix size is too small.
  246. *
  247. ZEROT = IMAT.GE.3 .AND. IMAT.LE.5
  248. IF( ZEROT .AND. N.LT.IMAT-2 )
  249. $ GO TO 110
  250. *
  251. * Do first for UPLO = 'U', then for UPLO = 'L'
  252. *
  253. DO 100 IUPLO = 1, 2
  254. UPLO = UPLOS( IUPLO )
  255. *
  256. * Set up parameters with DLATB4 and generate a test matrix
  257. * with DLATMS.
  258. *
  259. CALL DLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, MODE,
  260. $ CNDNUM, DIST )
  261. *
  262. SRNAMT = 'DLATMS'
  263. CALL DLATMS( N, N, DIST, ISEED, TYPE, RWORK, MODE,
  264. $ CNDNUM, ANORM, KL, KU, UPLO, A, LDA, WORK,
  265. $ INFO )
  266. *
  267. * Check error code from DLATMS.
  268. *
  269. IF( INFO.NE.0 ) THEN
  270. CALL ALAERH( PATH, 'DLATMS', INFO, 0, UPLO, N, N, -1,
  271. $ -1, -1, IMAT, NFAIL, NERRS, NOUT )
  272. GO TO 100
  273. END IF
  274. *
  275. * For types 3-5, zero one row and column of the matrix to
  276. * test that INFO is returned correctly.
  277. *
  278. IF( ZEROT ) THEN
  279. IF( IMAT.EQ.3 ) THEN
  280. IZERO = 1
  281. ELSE IF( IMAT.EQ.4 ) THEN
  282. IZERO = N
  283. ELSE
  284. IZERO = N / 2 + 1
  285. END IF
  286. IOFF = ( IZERO-1 )*LDA
  287. *
  288. * Set row and column IZERO of A to 0.
  289. *
  290. IF( IUPLO.EQ.1 ) THEN
  291. DO 20 I = 1, IZERO - 1
  292. A( IOFF+I ) = ZERO
  293. 20 CONTINUE
  294. IOFF = IOFF + IZERO
  295. DO 30 I = IZERO, N
  296. A( IOFF ) = ZERO
  297. IOFF = IOFF + LDA
  298. 30 CONTINUE
  299. ELSE
  300. IOFF = IZERO
  301. DO 40 I = 1, IZERO - 1
  302. A( IOFF ) = ZERO
  303. IOFF = IOFF + LDA
  304. 40 CONTINUE
  305. IOFF = IOFF - IZERO
  306. DO 50 I = IZERO, N
  307. A( IOFF+I ) = ZERO
  308. 50 CONTINUE
  309. END IF
  310. ELSE
  311. IZERO = 0
  312. END IF
  313. *
  314. DO 60 IRHS = 1, NNS
  315. NRHS = NSVAL( IRHS )
  316. XTYPE = 'N'
  317. *
  318. * Form an exact solution and set the right hand side.
  319. *
  320. SRNAMT = 'DLARHS'
  321. CALL DLARHS( PATH, XTYPE, UPLO, ' ', N, N, KL, KU,
  322. $ NRHS, A, LDA, X, LDA, B, LDA,
  323. $ ISEED, INFO )
  324. *
  325. * Compute the L*L' or U'*U factorization of the
  326. * matrix and solve the system.
  327. *
  328. SRNAMT = 'DSPOSV '
  329. KASE = KASE + 1
  330. *
  331. CALL DLACPY( 'All', N, N, A, LDA, AFAC, LDA)
  332. *
  333. CALL DSPOSV( UPLO, N, NRHS, AFAC, LDA, B, LDA, X, LDA,
  334. $ WORK, SWORK, ITER, INFO )
  335. IF (ITER.LT.0) THEN
  336. CALL DLACPY( 'All', N, N, A, LDA, AFAC, LDA )
  337. ENDIF
  338. *
  339. * Check error code from DSPOSV .
  340. *
  341. IF( INFO.NE.IZERO ) THEN
  342. *
  343. IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
  344. $ CALL ALAHD( NOUT, PATH )
  345. NERRS = NERRS + 1
  346. *
  347. IF( INFO.NE.IZERO .AND. IZERO.NE.0 ) THEN
  348. WRITE( NOUT, FMT = 9988 )'DSPOSV',INFO,IZERO,N,
  349. $ IMAT
  350. ELSE
  351. WRITE( NOUT, FMT = 9975 )'DSPOSV',INFO,N,IMAT
  352. END IF
  353. END IF
  354. *
  355. * Skip the remaining test if the matrix is singular.
  356. *
  357. IF( INFO.NE.0 )
  358. $ GO TO 110
  359. *
  360. * Check the quality of the solution
  361. *
  362. CALL DLACPY( 'All', N, NRHS, B, LDA, WORK, LDA )
  363. *
  364. CALL DPOT06( UPLO, N, NRHS, A, LDA, X, LDA, WORK,
  365. $ LDA, RWORK, RESULT( 1 ) )
  366. *
  367. * Check if the test passes the testing.
  368. * Print information about the tests that did not
  369. * pass the testing.
  370. *
  371. * If iterative refinement has been used and claimed to
  372. * be successful (ITER>0), we want
  373. * NORM1(B - A*X)/(NORM1(A)*NORM1(X)*EPS*SRQT(N)) < 1
  374. *
  375. * If double precision has been used (ITER<0), we want
  376. * NORM1(B - A*X)/(NORM1(A)*NORM1(X)*EPS) < THRES
  377. * (Cf. the linear solver testing routines)
  378. *
  379. IF ((THRESH.LE.0.0E+00)
  380. $ .OR.((ITER.GE.0).AND.(N.GT.0)
  381. $ .AND.(RESULT(1).GE.SQRT(DBLE(N))))
  382. $ .OR.((ITER.LT.0).AND.(RESULT(1).GE.THRESH))) THEN
  383. *
  384. IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) THEN
  385. WRITE( NOUT, FMT = 8999 )'DPO'
  386. WRITE( NOUT, FMT = '( '' Matrix types:'' )' )
  387. WRITE( NOUT, FMT = 8979 )
  388. WRITE( NOUT, FMT = '( '' Test ratios:'' )' )
  389. WRITE( NOUT, FMT = 8960 )1
  390. WRITE( NOUT, FMT = '( '' Messages:'' )' )
  391. END IF
  392. *
  393. WRITE( NOUT, FMT = 9998 )UPLO, N, NRHS, IMAT, 1,
  394. $ RESULT( 1 )
  395. *
  396. NFAIL = NFAIL + 1
  397. *
  398. END IF
  399. *
  400. NRUN = NRUN + 1
  401. *
  402. 60 CONTINUE
  403. 100 CONTINUE
  404. 110 CONTINUE
  405. 120 CONTINUE
  406. *
  407. * Print a summary of the results.
  408. *
  409. IF( NFAIL.GT.0 ) THEN
  410. WRITE( NOUT, FMT = 9996 )'DSPOSV', NFAIL, NRUN
  411. ELSE
  412. WRITE( NOUT, FMT = 9995 )'DSPOSV', NRUN
  413. END IF
  414. IF( NERRS.GT.0 ) THEN
  415. WRITE( NOUT, FMT = 9994 )NERRS
  416. END IF
  417. *
  418. 9998 FORMAT( ' UPLO=''', A1, ''', N =', I5, ', NRHS=', I3, ', type ',
  419. $ I2, ', test(', I2, ') =', G12.5 )
  420. 9996 FORMAT( 1X, A6, ': ', I6, ' out of ', I6,
  421. $ ' tests failed to pass the threshold' )
  422. 9995 FORMAT( /1X, 'All tests for ', A6,
  423. $ ' routines passed the threshold ( ', I6, ' tests run)' )
  424. 9994 FORMAT( 6X, I6, ' error messages recorded' )
  425. *
  426. * SUBNAM, INFO, INFOE, N, IMAT
  427. *
  428. 9988 FORMAT( ' *** ', A6, ' returned with INFO =', I5, ' instead of ',
  429. $ I5, / ' ==> N =', I5, ', type ',
  430. $ I2 )
  431. *
  432. * SUBNAM, INFO, N, IMAT
  433. *
  434. 9975 FORMAT( ' *** Error code from ', A6, '=', I5, ' for M=', I5,
  435. $ ', type ', I2 )
  436. 8999 FORMAT( / 1X, A3, ': positive definite dense matrices' )
  437. 8979 FORMAT( 4X, '1. Diagonal', 24X, '7. Last n/2 columns zero', / 4X,
  438. $ '2. Upper triangular', 16X,
  439. $ '8. Random, CNDNUM = sqrt(0.1/EPS)', / 4X,
  440. $ '3. Lower triangular', 16X, '9. Random, CNDNUM = 0.1/EPS',
  441. $ / 4X, '4. Random, CNDNUM = 2', 13X,
  442. $ '10. Scaled near underflow', / 4X, '5. First column zero',
  443. $ 14X, '11. Scaled near overflow', / 4X,
  444. $ '6. Last column zero' )
  445. 8960 FORMAT( 3X, I2, ': norm_1( B - A * X ) / ',
  446. $ '( norm_1(A) * norm_1(X) * EPS * SQRT(N) ) > 1 if ITERREF',
  447. $ / 4x, 'or norm_1( B - A * X ) / ',
  448. $ '( norm_1(A) * norm_1(X) * EPS ) > THRES if DPOTRF' )
  449. RETURN
  450. *
  451. * End of DDRVAC
  452. *
  453. END