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.

clatm4.f 12 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404
  1. *> \brief \b CLATM4
  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 CLATM4( ITYPE, N, NZ1, NZ2, RSIGN, AMAGN, RCOND,
  12. * TRIANG, IDIST, ISEED, A, LDA )
  13. *
  14. * .. Scalar Arguments ..
  15. * LOGICAL RSIGN
  16. * INTEGER IDIST, ITYPE, LDA, N, NZ1, NZ2
  17. * REAL AMAGN, RCOND, TRIANG
  18. * ..
  19. * .. Array Arguments ..
  20. * INTEGER ISEED( 4 )
  21. * COMPLEX A( LDA, * )
  22. * ..
  23. *
  24. *
  25. *> \par Purpose:
  26. * =============
  27. *>
  28. *> \verbatim
  29. *>
  30. *> CLATM4 generates basic square matrices, which may later be
  31. *> multiplied by others in order to produce test matrices. It is
  32. *> intended mainly to be used to test the generalized eigenvalue
  33. *> routines.
  34. *>
  35. *> It first generates the diagonal and (possibly) subdiagonal,
  36. *> according to the value of ITYPE, NZ1, NZ2, RSIGN, AMAGN, and RCOND.
  37. *> It then fills in the upper triangle with random numbers, if TRIANG is
  38. *> non-zero.
  39. *> \endverbatim
  40. *
  41. * Arguments:
  42. * ==========
  43. *
  44. *> \param[in] ITYPE
  45. *> \verbatim
  46. *> ITYPE is INTEGER
  47. *> The "type" of matrix on the diagonal and sub-diagonal.
  48. *> If ITYPE < 0, then type abs(ITYPE) is generated and then
  49. *> swapped end for end (A(I,J) := A'(N-J,N-I).) See also
  50. *> the description of AMAGN and RSIGN.
  51. *>
  52. *> Special types:
  53. *> = 0: the zero matrix.
  54. *> = 1: the identity.
  55. *> = 2: a transposed Jordan block.
  56. *> = 3: If N is odd, then a k+1 x k+1 transposed Jordan block
  57. *> followed by a k x k identity block, where k=(N-1)/2.
  58. *> If N is even, then k=(N-2)/2, and a zero diagonal entry
  59. *> is tacked onto the end.
  60. *>
  61. *> Diagonal types. The diagonal consists of NZ1 zeros, then
  62. *> k=N-NZ1-NZ2 nonzeros. The subdiagonal is zero. ITYPE
  63. *> specifies the nonzero diagonal entries as follows:
  64. *> = 4: 1, ..., k
  65. *> = 5: 1, RCOND, ..., RCOND
  66. *> = 6: 1, ..., 1, RCOND
  67. *> = 7: 1, a, a^2, ..., a^(k-1)=RCOND
  68. *> = 8: 1, 1-d, 1-2*d, ..., 1-(k-1)*d=RCOND
  69. *> = 9: random numbers chosen from (RCOND,1)
  70. *> = 10: random numbers with distribution IDIST (see CLARND.)
  71. *> \endverbatim
  72. *>
  73. *> \param[in] N
  74. *> \verbatim
  75. *> N is INTEGER
  76. *> The order of the matrix.
  77. *> \endverbatim
  78. *>
  79. *> \param[in] NZ1
  80. *> \verbatim
  81. *> NZ1 is INTEGER
  82. *> If abs(ITYPE) > 3, then the first NZ1 diagonal entries will
  83. *> be zero.
  84. *> \endverbatim
  85. *>
  86. *> \param[in] NZ2
  87. *> \verbatim
  88. *> NZ2 is INTEGER
  89. *> If abs(ITYPE) > 3, then the last NZ2 diagonal entries will
  90. *> be zero.
  91. *> \endverbatim
  92. *>
  93. *> \param[in] RSIGN
  94. *> \verbatim
  95. *> RSIGN is LOGICAL
  96. *> = .TRUE.: The diagonal and subdiagonal entries will be
  97. *> multiplied by random numbers of magnitude 1.
  98. *> = .FALSE.: The diagonal and subdiagonal entries will be
  99. *> left as they are (usually non-negative real.)
  100. *> \endverbatim
  101. *>
  102. *> \param[in] AMAGN
  103. *> \verbatim
  104. *> AMAGN is REAL
  105. *> The diagonal and subdiagonal entries will be multiplied by
  106. *> AMAGN.
  107. *> \endverbatim
  108. *>
  109. *> \param[in] RCOND
  110. *> \verbatim
  111. *> RCOND is REAL
  112. *> If abs(ITYPE) > 4, then the smallest diagonal entry will be
  113. *> RCOND. RCOND must be between 0 and 1.
  114. *> \endverbatim
  115. *>
  116. *> \param[in] TRIANG
  117. *> \verbatim
  118. *> TRIANG is REAL
  119. *> The entries above the diagonal will be random numbers with
  120. *> magnitude bounded by TRIANG (i.e., random numbers multiplied
  121. *> by TRIANG.)
  122. *> \endverbatim
  123. *>
  124. *> \param[in] IDIST
  125. *> \verbatim
  126. *> IDIST is INTEGER
  127. *> On entry, DIST specifies the type of distribution to be used
  128. *> to generate a random matrix .
  129. *> = 1: real and imaginary parts each UNIFORM( 0, 1 )
  130. *> = 2: real and imaginary parts each UNIFORM( -1, 1 )
  131. *> = 3: real and imaginary parts each NORMAL( 0, 1 )
  132. *> = 4: complex number uniform in DISK( 0, 1 )
  133. *> \endverbatim
  134. *>
  135. *> \param[in,out] ISEED
  136. *> \verbatim
  137. *> ISEED is INTEGER array, dimension (4)
  138. *> On entry ISEED specifies the seed of the random number
  139. *> generator. The values of ISEED are changed on exit, and can
  140. *> be used in the next call to CLATM4 to continue the same
  141. *> random number sequence.
  142. *> Note: ISEED(4) should be odd, for the random number generator
  143. *> used at present.
  144. *> \endverbatim
  145. *>
  146. *> \param[out] A
  147. *> \verbatim
  148. *> A is COMPLEX array, dimension (LDA, N)
  149. *> Array to be computed.
  150. *> \endverbatim
  151. *>
  152. *> \param[in] LDA
  153. *> \verbatim
  154. *> LDA is INTEGER
  155. *> Leading dimension of A. Must be at least 1 and at least N.
  156. *> \endverbatim
  157. *
  158. * Authors:
  159. * ========
  160. *
  161. *> \author Univ. of Tennessee
  162. *> \author Univ. of California Berkeley
  163. *> \author Univ. of Colorado Denver
  164. *> \author NAG Ltd.
  165. *
  166. *> \date November 2011
  167. *
  168. *> \ingroup complex_eig
  169. *
  170. * =====================================================================
  171. SUBROUTINE CLATM4( ITYPE, N, NZ1, NZ2, RSIGN, AMAGN, RCOND,
  172. $ TRIANG, IDIST, ISEED, A, LDA )
  173. *
  174. * -- LAPACK test routine (version 3.4.0) --
  175. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  176. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  177. * November 2011
  178. *
  179. * .. Scalar Arguments ..
  180. LOGICAL RSIGN
  181. INTEGER IDIST, ITYPE, LDA, N, NZ1, NZ2
  182. REAL AMAGN, RCOND, TRIANG
  183. * ..
  184. * .. Array Arguments ..
  185. INTEGER ISEED( 4 )
  186. COMPLEX A( LDA, * )
  187. * ..
  188. *
  189. * =====================================================================
  190. *
  191. * .. Parameters ..
  192. REAL ZERO, ONE
  193. PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 )
  194. COMPLEX CZERO, CONE
  195. PARAMETER ( CZERO = ( 0.0E+0, 0.0E+0 ),
  196. $ CONE = ( 1.0E+0, 0.0E+0 ) )
  197. * ..
  198. * .. Local Scalars ..
  199. INTEGER I, ISDB, ISDE, JC, JD, JR, K, KBEG, KEND, KLEN
  200. REAL ALPHA
  201. COMPLEX CTEMP
  202. * ..
  203. * .. External Functions ..
  204. REAL SLARAN
  205. COMPLEX CLARND
  206. EXTERNAL SLARAN, CLARND
  207. * ..
  208. * .. External Subroutines ..
  209. EXTERNAL CLASET
  210. * ..
  211. * .. Intrinsic Functions ..
  212. INTRINSIC ABS, CMPLX, EXP, LOG, MAX, MIN, MOD, REAL
  213. * ..
  214. * .. Executable Statements ..
  215. *
  216. IF( N.LE.0 )
  217. $ RETURN
  218. CALL CLASET( 'Full', N, N, CZERO, CZERO, A, LDA )
  219. *
  220. * Insure a correct ISEED
  221. *
  222. IF( MOD( ISEED( 4 ), 2 ).NE.1 )
  223. $ ISEED( 4 ) = ISEED( 4 ) + 1
  224. *
  225. * Compute diagonal and subdiagonal according to ITYPE, NZ1, NZ2,
  226. * and RCOND
  227. *
  228. IF( ITYPE.NE.0 ) THEN
  229. IF( ABS( ITYPE ).GE.4 ) THEN
  230. KBEG = MAX( 1, MIN( N, NZ1+1 ) )
  231. KEND = MAX( KBEG, MIN( N, N-NZ2 ) )
  232. KLEN = KEND + 1 - KBEG
  233. ELSE
  234. KBEG = 1
  235. KEND = N
  236. KLEN = N
  237. END IF
  238. ISDB = 1
  239. ISDE = 0
  240. GO TO ( 10, 30, 50, 80, 100, 120, 140, 160,
  241. $ 180, 200 )ABS( ITYPE )
  242. *
  243. * abs(ITYPE) = 1: Identity
  244. *
  245. 10 CONTINUE
  246. DO 20 JD = 1, N
  247. A( JD, JD ) = CONE
  248. 20 CONTINUE
  249. GO TO 220
  250. *
  251. * abs(ITYPE) = 2: Transposed Jordan block
  252. *
  253. 30 CONTINUE
  254. DO 40 JD = 1, N - 1
  255. A( JD+1, JD ) = CONE
  256. 40 CONTINUE
  257. ISDB = 1
  258. ISDE = N - 1
  259. GO TO 220
  260. *
  261. * abs(ITYPE) = 3: Transposed Jordan block, followed by the
  262. * identity.
  263. *
  264. 50 CONTINUE
  265. K = ( N-1 ) / 2
  266. DO 60 JD = 1, K
  267. A( JD+1, JD ) = CONE
  268. 60 CONTINUE
  269. ISDB = 1
  270. ISDE = K
  271. DO 70 JD = K + 2, 2*K + 1
  272. A( JD, JD ) = CONE
  273. 70 CONTINUE
  274. GO TO 220
  275. *
  276. * abs(ITYPE) = 4: 1,...,k
  277. *
  278. 80 CONTINUE
  279. DO 90 JD = KBEG, KEND
  280. A( JD, JD ) = CMPLX( JD-NZ1 )
  281. 90 CONTINUE
  282. GO TO 220
  283. *
  284. * abs(ITYPE) = 5: One large D value:
  285. *
  286. 100 CONTINUE
  287. DO 110 JD = KBEG + 1, KEND
  288. A( JD, JD ) = CMPLX( RCOND )
  289. 110 CONTINUE
  290. A( KBEG, KBEG ) = CONE
  291. GO TO 220
  292. *
  293. * abs(ITYPE) = 6: One small D value:
  294. *
  295. 120 CONTINUE
  296. DO 130 JD = KBEG, KEND - 1
  297. A( JD, JD ) = CONE
  298. 130 CONTINUE
  299. A( KEND, KEND ) = CMPLX( RCOND )
  300. GO TO 220
  301. *
  302. * abs(ITYPE) = 7: Exponentially distributed D values:
  303. *
  304. 140 CONTINUE
  305. A( KBEG, KBEG ) = CONE
  306. IF( KLEN.GT.1 ) THEN
  307. ALPHA = RCOND**( ONE / REAL( KLEN-1 ) )
  308. DO 150 I = 2, KLEN
  309. A( NZ1+I, NZ1+I ) = CMPLX( ALPHA**REAL( I-1 ) )
  310. 150 CONTINUE
  311. END IF
  312. GO TO 220
  313. *
  314. * abs(ITYPE) = 8: Arithmetically distributed D values:
  315. *
  316. 160 CONTINUE
  317. A( KBEG, KBEG ) = CONE
  318. IF( KLEN.GT.1 ) THEN
  319. ALPHA = ( ONE-RCOND ) / REAL( KLEN-1 )
  320. DO 170 I = 2, KLEN
  321. A( NZ1+I, NZ1+I ) = CMPLX( REAL( KLEN-I )*ALPHA+RCOND )
  322. 170 CONTINUE
  323. END IF
  324. GO TO 220
  325. *
  326. * abs(ITYPE) = 9: Randomly distributed D values on ( RCOND, 1):
  327. *
  328. 180 CONTINUE
  329. ALPHA = LOG( RCOND )
  330. DO 190 JD = KBEG, KEND
  331. A( JD, JD ) = EXP( ALPHA*SLARAN( ISEED ) )
  332. 190 CONTINUE
  333. GO TO 220
  334. *
  335. * abs(ITYPE) = 10: Randomly distributed D values from DIST
  336. *
  337. 200 CONTINUE
  338. DO 210 JD = KBEG, KEND
  339. A( JD, JD ) = CLARND( IDIST, ISEED )
  340. 210 CONTINUE
  341. *
  342. 220 CONTINUE
  343. *
  344. * Scale by AMAGN
  345. *
  346. DO 230 JD = KBEG, KEND
  347. A( JD, JD ) = AMAGN*REAL( A( JD, JD ) )
  348. 230 CONTINUE
  349. DO 240 JD = ISDB, ISDE
  350. A( JD+1, JD ) = AMAGN*REAL( A( JD+1, JD ) )
  351. 240 CONTINUE
  352. *
  353. * If RSIGN = .TRUE., assign random signs to diagonal and
  354. * subdiagonal
  355. *
  356. IF( RSIGN ) THEN
  357. DO 250 JD = KBEG, KEND
  358. IF( REAL( A( JD, JD ) ).NE.ZERO ) THEN
  359. CTEMP = CLARND( 3, ISEED )
  360. CTEMP = CTEMP / ABS( CTEMP )
  361. A( JD, JD ) = CTEMP*REAL( A( JD, JD ) )
  362. END IF
  363. 250 CONTINUE
  364. DO 260 JD = ISDB, ISDE
  365. IF( REAL( A( JD+1, JD ) ).NE.ZERO ) THEN
  366. CTEMP = CLARND( 3, ISEED )
  367. CTEMP = CTEMP / ABS( CTEMP )
  368. A( JD+1, JD ) = CTEMP*REAL( A( JD+1, JD ) )
  369. END IF
  370. 260 CONTINUE
  371. END IF
  372. *
  373. * Reverse if ITYPE < 0
  374. *
  375. IF( ITYPE.LT.0 ) THEN
  376. DO 270 JD = KBEG, ( KBEG+KEND-1 ) / 2
  377. CTEMP = A( JD, JD )
  378. A( JD, JD ) = A( KBEG+KEND-JD, KBEG+KEND-JD )
  379. A( KBEG+KEND-JD, KBEG+KEND-JD ) = CTEMP
  380. 270 CONTINUE
  381. DO 280 JD = 1, ( N-1 ) / 2
  382. CTEMP = A( JD+1, JD )
  383. A( JD+1, JD ) = A( N+1-JD, N-JD )
  384. A( N+1-JD, N-JD ) = CTEMP
  385. 280 CONTINUE
  386. END IF
  387. *
  388. END IF
  389. *
  390. * Fill in upper triangle
  391. *
  392. IF( TRIANG.NE.ZERO ) THEN
  393. DO 300 JC = 2, N
  394. DO 290 JR = 1, JC - 1
  395. A( JR, JC ) = TRIANG*CLARND( IDIST, ISEED )
  396. 290 CONTINUE
  397. 300 CONTINUE
  398. END IF
  399. *
  400. RETURN
  401. *
  402. * End of CLATM4
  403. *
  404. END