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.

slattp.f 24 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809
  1. *> \brief \b SLATTP
  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 SLATTP( IMAT, UPLO, TRANS, DIAG, ISEED, N, A, B, WORK,
  12. * INFO )
  13. *
  14. * .. Scalar Arguments ..
  15. * CHARACTER DIAG, TRANS, UPLO
  16. * INTEGER IMAT, INFO, N
  17. * ..
  18. * .. Array Arguments ..
  19. * INTEGER ISEED( 4 )
  20. * REAL A( * ), B( * ), WORK( * )
  21. * ..
  22. *
  23. *
  24. *> \par Purpose:
  25. * =============
  26. *>
  27. *> \verbatim
  28. *>
  29. *> SLATTP generates a triangular test matrix in packed storage.
  30. *> IMAT and UPLO uniquely specify the properties of the test
  31. *> matrix, which is returned in the array AP.
  32. *> \endverbatim
  33. *
  34. * Arguments:
  35. * ==========
  36. *
  37. *> \param[in] IMAT
  38. *> \verbatim
  39. *> IMAT is INTEGER
  40. *> An integer key describing which matrix to generate for this
  41. *> path.
  42. *> \endverbatim
  43. *>
  44. *> \param[in] UPLO
  45. *> \verbatim
  46. *> UPLO is CHARACTER*1
  47. *> Specifies whether the matrix A will be upper or lower
  48. *> triangular.
  49. *> = 'U': Upper triangular
  50. *> = 'L': Lower triangular
  51. *> \endverbatim
  52. *>
  53. *> \param[in] TRANS
  54. *> \verbatim
  55. *> TRANS is CHARACTER*1
  56. *> Specifies whether the matrix or its transpose will be used.
  57. *> = 'N': No transpose
  58. *> = 'T': Transpose
  59. *> = 'C': Conjugate transpose (= Transpose)
  60. *> \endverbatim
  61. *>
  62. *> \param[out] DIAG
  63. *> \verbatim
  64. *> DIAG is CHARACTER*1
  65. *> Specifies whether or not the matrix A is unit triangular.
  66. *> = 'N': Non-unit triangular
  67. *> = 'U': Unit triangular
  68. *> \endverbatim
  69. *>
  70. *> \param[in,out] ISEED
  71. *> \verbatim
  72. *> ISEED is INTEGER array, dimension (4)
  73. *> The seed vector for the random number generator (used in
  74. *> SLATMS). Modified on exit.
  75. *> \endverbatim
  76. *>
  77. *> \param[in] N
  78. *> \verbatim
  79. *> N is INTEGER
  80. *> The order of the matrix to be generated.
  81. *> \endverbatim
  82. *>
  83. *> \param[out] A
  84. *> \verbatim
  85. *> A is REAL array, dimension (N*(N+1)/2)
  86. *> The upper or lower triangular matrix A, packed columnwise in
  87. *> a linear array. The j-th column of A is stored in the array
  88. *> AP as follows:
  89. *> if UPLO = 'U', AP((j-1)*j/2 + i) = A(i,j) for 1<=i<=j;
  90. *> if UPLO = 'L',
  91. *> AP((j-1)*(n-j) + j*(j+1)/2 + i-j) = A(i,j) for j<=i<=n.
  92. *> \endverbatim
  93. *>
  94. *> \param[out] B
  95. *> \verbatim
  96. *> B is REAL array, dimension (N)
  97. *> The right hand side vector, if IMAT > 10.
  98. *> \endverbatim
  99. *>
  100. *> \param[out] WORK
  101. *> \verbatim
  102. *> WORK is REAL array, dimension (3*N)
  103. *> \endverbatim
  104. *>
  105. *> \param[out] INFO
  106. *> \verbatim
  107. *> INFO is INTEGER
  108. *> = 0: successful exit
  109. *> < 0: if INFO = -k, the k-th argument had an illegal value
  110. *> \endverbatim
  111. *
  112. * Authors:
  113. * ========
  114. *
  115. *> \author Univ. of Tennessee
  116. *> \author Univ. of California Berkeley
  117. *> \author Univ. of Colorado Denver
  118. *> \author NAG Ltd.
  119. *
  120. *> \ingroup single_lin
  121. *
  122. * =====================================================================
  123. SUBROUTINE SLATTP( IMAT, UPLO, TRANS, DIAG, ISEED, N, A, B, WORK,
  124. $ INFO )
  125. *
  126. * -- LAPACK test routine --
  127. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  128. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  129. *
  130. * .. Scalar Arguments ..
  131. CHARACTER DIAG, TRANS, UPLO
  132. INTEGER IMAT, INFO, N
  133. * ..
  134. * .. Array Arguments ..
  135. INTEGER ISEED( 4 )
  136. REAL A( * ), B( * ), WORK( * )
  137. * ..
  138. *
  139. * =====================================================================
  140. *
  141. * .. Parameters ..
  142. REAL ONE, TWO, ZERO
  143. PARAMETER ( ONE = 1.0E+0, TWO = 2.0E+0, ZERO = 0.0E+0 )
  144. * ..
  145. * .. Local Scalars ..
  146. LOGICAL UPPER
  147. CHARACTER DIST, PACKIT, TYPE
  148. CHARACTER*3 PATH
  149. INTEGER I, IY, J, JC, JCNEXT, JCOUNT, JJ, JL, JR, JX,
  150. $ KL, KU, MODE
  151. REAL ANORM, BIGNUM, BNORM, BSCAL, C, CNDNUM, PLUS1,
  152. $ PLUS2, RA, RB, REXP, S, SFAC, SMLNUM, STAR1,
  153. $ STEMP, T, TEXP, TLEFT, TSCAL, ULP, UNFL, X, Y,
  154. $ Z
  155. * ..
  156. * .. External Functions ..
  157. LOGICAL LSAME
  158. INTEGER ISAMAX
  159. REAL SLAMCH, SLARND
  160. EXTERNAL LSAME, ISAMAX, SLAMCH, SLARND
  161. * ..
  162. * .. External Subroutines ..
  163. EXTERNAL SLABAD, SLARNV, SLATB4, SLATMS, SROT, SROTG,
  164. $ SSCAL
  165. * ..
  166. * .. Intrinsic Functions ..
  167. INTRINSIC ABS, MAX, REAL, SIGN, SQRT
  168. * ..
  169. * .. Executable Statements ..
  170. *
  171. PATH( 1: 1 ) = 'Single precision'
  172. PATH( 2: 3 ) = 'TP'
  173. UNFL = SLAMCH( 'Safe minimum' )
  174. ULP = SLAMCH( 'Epsilon' )*SLAMCH( 'Base' )
  175. SMLNUM = UNFL
  176. BIGNUM = ( ONE-ULP ) / SMLNUM
  177. CALL SLABAD( SMLNUM, BIGNUM )
  178. IF( ( IMAT.GE.7 .AND. IMAT.LE.10 ) .OR. IMAT.EQ.18 ) THEN
  179. DIAG = 'U'
  180. ELSE
  181. DIAG = 'N'
  182. END IF
  183. INFO = 0
  184. *
  185. * Quick return if N.LE.0.
  186. *
  187. IF( N.LE.0 )
  188. $ RETURN
  189. *
  190. * Call SLATB4 to set parameters for SLATMS.
  191. *
  192. UPPER = LSAME( UPLO, 'U' )
  193. IF( UPPER ) THEN
  194. CALL SLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, MODE,
  195. $ CNDNUM, DIST )
  196. PACKIT = 'C'
  197. ELSE
  198. CALL SLATB4( PATH, -IMAT, N, N, TYPE, KL, KU, ANORM, MODE,
  199. $ CNDNUM, DIST )
  200. PACKIT = 'R'
  201. END IF
  202. *
  203. * IMAT <= 6: Non-unit triangular matrix
  204. *
  205. IF( IMAT.LE.6 ) THEN
  206. CALL SLATMS( N, N, DIST, ISEED, TYPE, B, MODE, CNDNUM, ANORM,
  207. $ KL, KU, PACKIT, A, N, WORK, INFO )
  208. *
  209. * IMAT > 6: Unit triangular matrix
  210. * The diagonal is deliberately set to something other than 1.
  211. *
  212. * IMAT = 7: Matrix is the identity
  213. *
  214. ELSE IF( IMAT.EQ.7 ) THEN
  215. IF( UPPER ) THEN
  216. JC = 1
  217. DO 20 J = 1, N
  218. DO 10 I = 1, J - 1
  219. A( JC+I-1 ) = ZERO
  220. 10 CONTINUE
  221. A( JC+J-1 ) = J
  222. JC = JC + J
  223. 20 CONTINUE
  224. ELSE
  225. JC = 1
  226. DO 40 J = 1, N
  227. A( JC ) = J
  228. DO 30 I = J + 1, N
  229. A( JC+I-J ) = ZERO
  230. 30 CONTINUE
  231. JC = JC + N - J + 1
  232. 40 CONTINUE
  233. END IF
  234. *
  235. * IMAT > 7: Non-trivial unit triangular matrix
  236. *
  237. * Generate a unit triangular matrix T with condition CNDNUM by
  238. * forming a triangular matrix with known singular values and
  239. * filling in the zero entries with Givens rotations.
  240. *
  241. ELSE IF( IMAT.LE.10 ) THEN
  242. IF( UPPER ) THEN
  243. JC = 0
  244. DO 60 J = 1, N
  245. DO 50 I = 1, J - 1
  246. A( JC+I ) = ZERO
  247. 50 CONTINUE
  248. A( JC+J ) = J
  249. JC = JC + J
  250. 60 CONTINUE
  251. ELSE
  252. JC = 1
  253. DO 80 J = 1, N
  254. A( JC ) = J
  255. DO 70 I = J + 1, N
  256. A( JC+I-J ) = ZERO
  257. 70 CONTINUE
  258. JC = JC + N - J + 1
  259. 80 CONTINUE
  260. END IF
  261. *
  262. * Since the trace of a unit triangular matrix is 1, the product
  263. * of its singular values must be 1. Let s = sqrt(CNDNUM),
  264. * x = sqrt(s) - 1/sqrt(s), y = sqrt(2/(n-2))*x, and z = x**2.
  265. * The following triangular matrix has singular values s, 1, 1,
  266. * ..., 1, 1/s:
  267. *
  268. * 1 y y y ... y y z
  269. * 1 0 0 ... 0 0 y
  270. * 1 0 ... 0 0 y
  271. * . ... . . .
  272. * . . . .
  273. * 1 0 y
  274. * 1 y
  275. * 1
  276. *
  277. * To fill in the zeros, we first multiply by a matrix with small
  278. * condition number of the form
  279. *
  280. * 1 0 0 0 0 ...
  281. * 1 + * 0 0 ...
  282. * 1 + 0 0 0
  283. * 1 + * 0 0
  284. * 1 + 0 0
  285. * ...
  286. * 1 + 0
  287. * 1 0
  288. * 1
  289. *
  290. * Each element marked with a '*' is formed by taking the product
  291. * of the adjacent elements marked with '+'. The '*'s can be
  292. * chosen freely, and the '+'s are chosen so that the inverse of
  293. * T will have elements of the same magnitude as T. If the *'s in
  294. * both T and inv(T) have small magnitude, T is well conditioned.
  295. * The two offdiagonals of T are stored in WORK.
  296. *
  297. * The product of these two matrices has the form
  298. *
  299. * 1 y y y y y . y y z
  300. * 1 + * 0 0 . 0 0 y
  301. * 1 + 0 0 . 0 0 y
  302. * 1 + * . . . .
  303. * 1 + . . . .
  304. * . . . . .
  305. * . . . .
  306. * 1 + y
  307. * 1 y
  308. * 1
  309. *
  310. * Now we multiply by Givens rotations, using the fact that
  311. *
  312. * [ c s ] [ 1 w ] [ -c -s ] = [ 1 -w ]
  313. * [ -s c ] [ 0 1 ] [ s -c ] [ 0 1 ]
  314. * and
  315. * [ -c -s ] [ 1 0 ] [ c s ] = [ 1 0 ]
  316. * [ s -c ] [ w 1 ] [ -s c ] [ -w 1 ]
  317. *
  318. * where c = w / sqrt(w**2+4) and s = 2 / sqrt(w**2+4).
  319. *
  320. STAR1 = 0.25
  321. SFAC = 0.5
  322. PLUS1 = SFAC
  323. DO 90 J = 1, N, 2
  324. PLUS2 = STAR1 / PLUS1
  325. WORK( J ) = PLUS1
  326. WORK( N+J ) = STAR1
  327. IF( J+1.LE.N ) THEN
  328. WORK( J+1 ) = PLUS2
  329. WORK( N+J+1 ) = ZERO
  330. PLUS1 = STAR1 / PLUS2
  331. REXP = SLARND( 2, ISEED )
  332. STAR1 = STAR1*( SFAC**REXP )
  333. IF( REXP.LT.ZERO ) THEN
  334. STAR1 = -SFAC**( ONE-REXP )
  335. ELSE
  336. STAR1 = SFAC**( ONE+REXP )
  337. END IF
  338. END IF
  339. 90 CONTINUE
  340. *
  341. X = SQRT( CNDNUM ) - ONE / SQRT( CNDNUM )
  342. IF( N.GT.2 ) THEN
  343. Y = SQRT( TWO / REAL( N-2 ) )*X
  344. ELSE
  345. Y = ZERO
  346. END IF
  347. Z = X*X
  348. *
  349. IF( UPPER ) THEN
  350. *
  351. * Set the upper triangle of A with a unit triangular matrix
  352. * of known condition number.
  353. *
  354. JC = 1
  355. DO 100 J = 2, N
  356. A( JC+1 ) = Y
  357. IF( J.GT.2 )
  358. $ A( JC+J-1 ) = WORK( J-2 )
  359. IF( J.GT.3 )
  360. $ A( JC+J-2 ) = WORK( N+J-3 )
  361. JC = JC + J
  362. 100 CONTINUE
  363. JC = JC - N
  364. A( JC+1 ) = Z
  365. DO 110 J = 2, N - 1
  366. A( JC+J ) = Y
  367. 110 CONTINUE
  368. ELSE
  369. *
  370. * Set the lower triangle of A with a unit triangular matrix
  371. * of known condition number.
  372. *
  373. DO 120 I = 2, N - 1
  374. A( I ) = Y
  375. 120 CONTINUE
  376. A( N ) = Z
  377. JC = N + 1
  378. DO 130 J = 2, N - 1
  379. A( JC+1 ) = WORK( J-1 )
  380. IF( J.LT.N-1 )
  381. $ A( JC+2 ) = WORK( N+J-1 )
  382. A( JC+N-J ) = Y
  383. JC = JC + N - J + 1
  384. 130 CONTINUE
  385. END IF
  386. *
  387. * Fill in the zeros using Givens rotations
  388. *
  389. IF( UPPER ) THEN
  390. JC = 1
  391. DO 150 J = 1, N - 1
  392. JCNEXT = JC + J
  393. RA = A( JCNEXT+J-1 )
  394. RB = TWO
  395. CALL SROTG( RA, RB, C, S )
  396. *
  397. * Multiply by [ c s; -s c] on the left.
  398. *
  399. IF( N.GT.J+1 ) THEN
  400. JX = JCNEXT + J
  401. DO 140 I = J + 2, N
  402. STEMP = C*A( JX+J ) + S*A( JX+J+1 )
  403. A( JX+J+1 ) = -S*A( JX+J ) + C*A( JX+J+1 )
  404. A( JX+J ) = STEMP
  405. JX = JX + I
  406. 140 CONTINUE
  407. END IF
  408. *
  409. * Multiply by [-c -s; s -c] on the right.
  410. *
  411. IF( J.GT.1 )
  412. $ CALL SROT( J-1, A( JCNEXT ), 1, A( JC ), 1, -C, -S )
  413. *
  414. * Negate A(J,J+1).
  415. *
  416. A( JCNEXT+J-1 ) = -A( JCNEXT+J-1 )
  417. JC = JCNEXT
  418. 150 CONTINUE
  419. ELSE
  420. JC = 1
  421. DO 170 J = 1, N - 1
  422. JCNEXT = JC + N - J + 1
  423. RA = A( JC+1 )
  424. RB = TWO
  425. CALL SROTG( RA, RB, C, S )
  426. *
  427. * Multiply by [ c -s; s c] on the right.
  428. *
  429. IF( N.GT.J+1 )
  430. $ CALL SROT( N-J-1, A( JCNEXT+1 ), 1, A( JC+2 ), 1, C,
  431. $ -S )
  432. *
  433. * Multiply by [-c s; -s -c] on the left.
  434. *
  435. IF( J.GT.1 ) THEN
  436. JX = 1
  437. DO 160 I = 1, J - 1
  438. STEMP = -C*A( JX+J-I ) + S*A( JX+J-I+1 )
  439. A( JX+J-I+1 ) = -S*A( JX+J-I ) - C*A( JX+J-I+1 )
  440. A( JX+J-I ) = STEMP
  441. JX = JX + N - I + 1
  442. 160 CONTINUE
  443. END IF
  444. *
  445. * Negate A(J+1,J).
  446. *
  447. A( JC+1 ) = -A( JC+1 )
  448. JC = JCNEXT
  449. 170 CONTINUE
  450. END IF
  451. *
  452. * IMAT > 10: Pathological test cases. These triangular matrices
  453. * are badly scaled or badly conditioned, so when used in solving a
  454. * triangular system they may cause overflow in the solution vector.
  455. *
  456. ELSE IF( IMAT.EQ.11 ) THEN
  457. *
  458. * Type 11: Generate a triangular matrix with elements between
  459. * -1 and 1. Give the diagonal norm 2 to make it well-conditioned.
  460. * Make the right hand side large so that it requires scaling.
  461. *
  462. IF( UPPER ) THEN
  463. JC = 1
  464. DO 180 J = 1, N
  465. CALL SLARNV( 2, ISEED, J, A( JC ) )
  466. A( JC+J-1 ) = SIGN( TWO, A( JC+J-1 ) )
  467. JC = JC + J
  468. 180 CONTINUE
  469. ELSE
  470. JC = 1
  471. DO 190 J = 1, N
  472. CALL SLARNV( 2, ISEED, N-J+1, A( JC ) )
  473. A( JC ) = SIGN( TWO, A( JC ) )
  474. JC = JC + N - J + 1
  475. 190 CONTINUE
  476. END IF
  477. *
  478. * Set the right hand side so that the largest value is BIGNUM.
  479. *
  480. CALL SLARNV( 2, ISEED, N, B )
  481. IY = ISAMAX( N, B, 1 )
  482. BNORM = ABS( B( IY ) )
  483. BSCAL = BIGNUM / MAX( ONE, BNORM )
  484. CALL SSCAL( N, BSCAL, B, 1 )
  485. *
  486. ELSE IF( IMAT.EQ.12 ) THEN
  487. *
  488. * Type 12: Make the first diagonal element in the solve small to
  489. * cause immediate overflow when dividing by T(j,j).
  490. * In type 12, the offdiagonal elements are small (CNORM(j) < 1).
  491. *
  492. CALL SLARNV( 2, ISEED, N, B )
  493. TSCAL = ONE / MAX( ONE, REAL( N-1 ) )
  494. IF( UPPER ) THEN
  495. JC = 1
  496. DO 200 J = 1, N
  497. CALL SLARNV( 2, ISEED, J-1, A( JC ) )
  498. CALL SSCAL( J-1, TSCAL, A( JC ), 1 )
  499. A( JC+J-1 ) = SIGN( ONE, SLARND( 2, ISEED ) )
  500. JC = JC + J
  501. 200 CONTINUE
  502. A( N*( N+1 ) / 2 ) = SMLNUM
  503. ELSE
  504. JC = 1
  505. DO 210 J = 1, N
  506. CALL SLARNV( 2, ISEED, N-J, A( JC+1 ) )
  507. CALL SSCAL( N-J, TSCAL, A( JC+1 ), 1 )
  508. A( JC ) = SIGN( ONE, SLARND( 2, ISEED ) )
  509. JC = JC + N - J + 1
  510. 210 CONTINUE
  511. A( 1 ) = SMLNUM
  512. END IF
  513. *
  514. ELSE IF( IMAT.EQ.13 ) THEN
  515. *
  516. * Type 13: Make the first diagonal element in the solve small to
  517. * cause immediate overflow when dividing by T(j,j).
  518. * In type 13, the offdiagonal elements are O(1) (CNORM(j) > 1).
  519. *
  520. CALL SLARNV( 2, ISEED, N, B )
  521. IF( UPPER ) THEN
  522. JC = 1
  523. DO 220 J = 1, N
  524. CALL SLARNV( 2, ISEED, J-1, A( JC ) )
  525. A( JC+J-1 ) = SIGN( ONE, SLARND( 2, ISEED ) )
  526. JC = JC + J
  527. 220 CONTINUE
  528. A( N*( N+1 ) / 2 ) = SMLNUM
  529. ELSE
  530. JC = 1
  531. DO 230 J = 1, N
  532. CALL SLARNV( 2, ISEED, N-J, A( JC+1 ) )
  533. A( JC ) = SIGN( ONE, SLARND( 2, ISEED ) )
  534. JC = JC + N - J + 1
  535. 230 CONTINUE
  536. A( 1 ) = SMLNUM
  537. END IF
  538. *
  539. ELSE IF( IMAT.EQ.14 ) THEN
  540. *
  541. * Type 14: T is diagonal with small numbers on the diagonal to
  542. * make the growth factor underflow, but a small right hand side
  543. * chosen so that the solution does not overflow.
  544. *
  545. IF( UPPER ) THEN
  546. JCOUNT = 1
  547. JC = ( N-1 )*N / 2 + 1
  548. DO 250 J = N, 1, -1
  549. DO 240 I = 1, J - 1
  550. A( JC+I-1 ) = ZERO
  551. 240 CONTINUE
  552. IF( JCOUNT.LE.2 ) THEN
  553. A( JC+J-1 ) = SMLNUM
  554. ELSE
  555. A( JC+J-1 ) = ONE
  556. END IF
  557. JCOUNT = JCOUNT + 1
  558. IF( JCOUNT.GT.4 )
  559. $ JCOUNT = 1
  560. JC = JC - J + 1
  561. 250 CONTINUE
  562. ELSE
  563. JCOUNT = 1
  564. JC = 1
  565. DO 270 J = 1, N
  566. DO 260 I = J + 1, N
  567. A( JC+I-J ) = ZERO
  568. 260 CONTINUE
  569. IF( JCOUNT.LE.2 ) THEN
  570. A( JC ) = SMLNUM
  571. ELSE
  572. A( JC ) = ONE
  573. END IF
  574. JCOUNT = JCOUNT + 1
  575. IF( JCOUNT.GT.4 )
  576. $ JCOUNT = 1
  577. JC = JC + N - J + 1
  578. 270 CONTINUE
  579. END IF
  580. *
  581. * Set the right hand side alternately zero and small.
  582. *
  583. IF( UPPER ) THEN
  584. B( 1 ) = ZERO
  585. DO 280 I = N, 2, -2
  586. B( I ) = ZERO
  587. B( I-1 ) = SMLNUM
  588. 280 CONTINUE
  589. ELSE
  590. B( N ) = ZERO
  591. DO 290 I = 1, N - 1, 2
  592. B( I ) = ZERO
  593. B( I+1 ) = SMLNUM
  594. 290 CONTINUE
  595. END IF
  596. *
  597. ELSE IF( IMAT.EQ.15 ) THEN
  598. *
  599. * Type 15: Make the diagonal elements small to cause gradual
  600. * overflow when dividing by T(j,j). To control the amount of
  601. * scaling needed, the matrix is bidiagonal.
  602. *
  603. TEXP = ONE / MAX( ONE, REAL( N-1 ) )
  604. TSCAL = SMLNUM**TEXP
  605. CALL SLARNV( 2, ISEED, N, B )
  606. IF( UPPER ) THEN
  607. JC = 1
  608. DO 310 J = 1, N
  609. DO 300 I = 1, J - 2
  610. A( JC+I-1 ) = ZERO
  611. 300 CONTINUE
  612. IF( J.GT.1 )
  613. $ A( JC+J-2 ) = -ONE
  614. A( JC+J-1 ) = TSCAL
  615. JC = JC + J
  616. 310 CONTINUE
  617. B( N ) = ONE
  618. ELSE
  619. JC = 1
  620. DO 330 J = 1, N
  621. DO 320 I = J + 2, N
  622. A( JC+I-J ) = ZERO
  623. 320 CONTINUE
  624. IF( J.LT.N )
  625. $ A( JC+1 ) = -ONE
  626. A( JC ) = TSCAL
  627. JC = JC + N - J + 1
  628. 330 CONTINUE
  629. B( 1 ) = ONE
  630. END IF
  631. *
  632. ELSE IF( IMAT.EQ.16 ) THEN
  633. *
  634. * Type 16: One zero diagonal element.
  635. *
  636. IY = N / 2 + 1
  637. IF( UPPER ) THEN
  638. JC = 1
  639. DO 340 J = 1, N
  640. CALL SLARNV( 2, ISEED, J, A( JC ) )
  641. IF( J.NE.IY ) THEN
  642. A( JC+J-1 ) = SIGN( TWO, A( JC+J-1 ) )
  643. ELSE
  644. A( JC+J-1 ) = ZERO
  645. END IF
  646. JC = JC + J
  647. 340 CONTINUE
  648. ELSE
  649. JC = 1
  650. DO 350 J = 1, N
  651. CALL SLARNV( 2, ISEED, N-J+1, A( JC ) )
  652. IF( J.NE.IY ) THEN
  653. A( JC ) = SIGN( TWO, A( JC ) )
  654. ELSE
  655. A( JC ) = ZERO
  656. END IF
  657. JC = JC + N - J + 1
  658. 350 CONTINUE
  659. END IF
  660. CALL SLARNV( 2, ISEED, N, B )
  661. CALL SSCAL( N, TWO, B, 1 )
  662. *
  663. ELSE IF( IMAT.EQ.17 ) THEN
  664. *
  665. * Type 17: Make the offdiagonal elements large to cause overflow
  666. * when adding a column of T. In the non-transposed case, the
  667. * matrix is constructed to cause overflow when adding a column in
  668. * every other step.
  669. *
  670. TSCAL = UNFL / ULP
  671. TSCAL = ( ONE-ULP ) / TSCAL
  672. DO 360 J = 1, N*( N+1 ) / 2
  673. A( J ) = ZERO
  674. 360 CONTINUE
  675. TEXP = ONE
  676. IF( UPPER ) THEN
  677. JC = ( N-1 )*N / 2 + 1
  678. DO 370 J = N, 2, -2
  679. A( JC ) = -TSCAL / REAL( N+1 )
  680. A( JC+J-1 ) = ONE
  681. B( J ) = TEXP*( ONE-ULP )
  682. JC = JC - J + 1
  683. A( JC ) = -( TSCAL / REAL( N+1 ) ) / REAL( N+2 )
  684. A( JC+J-2 ) = ONE
  685. B( J-1 ) = TEXP*REAL( N*N+N-1 )
  686. TEXP = TEXP*TWO
  687. JC = JC - J + 2
  688. 370 CONTINUE
  689. B( 1 ) = ( REAL( N+1 ) / REAL( N+2 ) )*TSCAL
  690. ELSE
  691. JC = 1
  692. DO 380 J = 1, N - 1, 2
  693. A( JC+N-J ) = -TSCAL / REAL( N+1 )
  694. A( JC ) = ONE
  695. B( J ) = TEXP*( ONE-ULP )
  696. JC = JC + N - J + 1
  697. A( JC+N-J-1 ) = -( TSCAL / REAL( N+1 ) ) / REAL( N+2 )
  698. A( JC ) = ONE
  699. B( J+1 ) = TEXP*REAL( N*N+N-1 )
  700. TEXP = TEXP*TWO
  701. JC = JC + N - J
  702. 380 CONTINUE
  703. B( N ) = ( REAL( N+1 ) / REAL( N+2 ) )*TSCAL
  704. END IF
  705. *
  706. ELSE IF( IMAT.EQ.18 ) THEN
  707. *
  708. * Type 18: Generate a unit triangular matrix with elements
  709. * between -1 and 1, and make the right hand side large so that it
  710. * requires scaling.
  711. *
  712. IF( UPPER ) THEN
  713. JC = 1
  714. DO 390 J = 1, N
  715. CALL SLARNV( 2, ISEED, J-1, A( JC ) )
  716. A( JC+J-1 ) = ZERO
  717. JC = JC + J
  718. 390 CONTINUE
  719. ELSE
  720. JC = 1
  721. DO 400 J = 1, N
  722. IF( J.LT.N )
  723. $ CALL SLARNV( 2, ISEED, N-J, A( JC+1 ) )
  724. A( JC ) = ZERO
  725. JC = JC + N - J + 1
  726. 400 CONTINUE
  727. END IF
  728. *
  729. * Set the right hand side so that the largest value is BIGNUM.
  730. *
  731. CALL SLARNV( 2, ISEED, N, B )
  732. IY = ISAMAX( N, B, 1 )
  733. BNORM = ABS( B( IY ) )
  734. BSCAL = BIGNUM / MAX( ONE, BNORM )
  735. CALL SSCAL( N, BSCAL, B, 1 )
  736. *
  737. ELSE IF( IMAT.EQ.19 ) THEN
  738. *
  739. * Type 19: Generate a triangular matrix with elements between
  740. * BIGNUM/(n-1) and BIGNUM so that at least one of the column
  741. * norms will exceed BIGNUM.
  742. *
  743. TLEFT = BIGNUM / MAX( ONE, REAL( N-1 ) )
  744. TSCAL = BIGNUM*( REAL( N-1 ) / MAX( ONE, REAL( N ) ) )
  745. IF( UPPER ) THEN
  746. JC = 1
  747. DO 420 J = 1, N
  748. CALL SLARNV( 2, ISEED, J, A( JC ) )
  749. DO 410 I = 1, J
  750. A( JC+I-1 ) = SIGN( TLEFT, A( JC+I-1 ) ) +
  751. $ TSCAL*A( JC+I-1 )
  752. 410 CONTINUE
  753. JC = JC + J
  754. 420 CONTINUE
  755. ELSE
  756. JC = 1
  757. DO 440 J = 1, N
  758. CALL SLARNV( 2, ISEED, N-J+1, A( JC ) )
  759. DO 430 I = J, N
  760. A( JC+I-J ) = SIGN( TLEFT, A( JC+I-J ) ) +
  761. $ TSCAL*A( JC+I-J )
  762. 430 CONTINUE
  763. JC = JC + N - J + 1
  764. 440 CONTINUE
  765. END IF
  766. CALL SLARNV( 2, ISEED, N, B )
  767. CALL SSCAL( N, TWO, B, 1 )
  768. END IF
  769. *
  770. * Flip the matrix across its counter-diagonal if the transpose will
  771. * be used.
  772. *
  773. IF( .NOT.LSAME( TRANS, 'N' ) ) THEN
  774. IF( UPPER ) THEN
  775. JJ = 1
  776. JR = N*( N+1 ) / 2
  777. DO 460 J = 1, N / 2
  778. JL = JJ
  779. DO 450 I = J, N - J
  780. T = A( JR-I+J )
  781. A( JR-I+J ) = A( JL )
  782. A( JL ) = T
  783. JL = JL + I
  784. 450 CONTINUE
  785. JJ = JJ + J + 1
  786. JR = JR - ( N-J+1 )
  787. 460 CONTINUE
  788. ELSE
  789. JL = 1
  790. JJ = N*( N+1 ) / 2
  791. DO 480 J = 1, N / 2
  792. JR = JJ
  793. DO 470 I = J, N - J
  794. T = A( JL+I-J )
  795. A( JL+I-J ) = A( JR )
  796. A( JR ) = T
  797. JR = JR - I
  798. 470 CONTINUE
  799. JL = JL + N - J + 1
  800. JJ = JJ - J - 1
  801. 480 CONTINUE
  802. END IF
  803. END IF
  804. *
  805. RETURN
  806. *
  807. * End of SLATTP
  808. *
  809. END