ctpmv.f 11 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345
  1. *DECK CTPMV
  2. SUBROUTINE CTPMV (UPLO, TRANS, DIAG, N, AP, X, INCX)
  3. C***BEGIN PROLOGUE CTPMV
  4. C***PURPOSE Perform one of the matrix-vector operations.
  5. C***LIBRARY SLATEC (BLAS)
  6. C***CATEGORY D1B4
  7. C***TYPE COMPLEX (STPMV-S, DTPMV-D, CTPMV-C)
  8. C***KEYWORDS LEVEL 2 BLAS, LINEAR ALGEBRA
  9. C***AUTHOR Dongarra, J. J., (ANL)
  10. C Du Croz, J., (NAG)
  11. C Hammarling, S., (NAG)
  12. C Hanson, R. J., (SNLA)
  13. C***DESCRIPTION
  14. C
  15. C CTPMV performs one of the matrix-vector operations
  16. C
  17. C x := A*x, or x := A'*x, or x := conjg( A')*x,
  18. C
  19. C where x is an n element vector and A is an n by n unit, or non-unit,
  20. C upper or lower triangular matrix, supplied in packed form.
  21. C
  22. C Parameters
  23. C ==========
  24. C
  25. C UPLO - CHARACTER*1.
  26. C On entry, UPLO specifies whether the matrix is an upper or
  27. C lower triangular matrix as follows:
  28. C
  29. C UPLO = 'U' or 'u' A is an upper triangular matrix.
  30. C
  31. C UPLO = 'L' or 'l' A is a lower triangular matrix.
  32. C
  33. C Unchanged on exit.
  34. C
  35. C TRANS - CHARACTER*1.
  36. C On entry, TRANS specifies the operation to be performed as
  37. C follows:
  38. C
  39. C TRANS = 'N' or 'n' x := A*x.
  40. C
  41. C TRANS = 'T' or 't' x := A'*x.
  42. C
  43. C TRANS = 'C' or 'c' x := conjg( A' )*x.
  44. C
  45. C Unchanged on exit.
  46. C
  47. C DIAG - CHARACTER*1.
  48. C On entry, DIAG specifies whether or not A is unit
  49. C triangular as follows:
  50. C
  51. C DIAG = 'U' or 'u' A is assumed to be unit triangular.
  52. C
  53. C DIAG = 'N' or 'n' A is not assumed to be unit
  54. C triangular.
  55. C
  56. C Unchanged on exit.
  57. C
  58. C N - INTEGER.
  59. C On entry, N specifies the order of the matrix A.
  60. C N must be at least zero.
  61. C Unchanged on exit.
  62. C
  63. C AP - COMPLEX array of DIMENSION at least
  64. C ( ( n*( n + 1 ) )/2 ).
  65. C Before entry with UPLO = 'U' or 'u', the array AP must
  66. C contain the upper triangular matrix packed sequentially,
  67. C column by column, so that AP( 1 ) contains a( 1, 1 ),
  68. C AP( 2 ) and AP( 3 ) contain a( 1, 2 ) and a( 2, 2 )
  69. C respectively, and so on.
  70. C Before entry with UPLO = 'L' or 'l', the array AP must
  71. C contain the lower triangular matrix packed sequentially,
  72. C column by column, so that AP( 1 ) contains a( 1, 1 ),
  73. C AP( 2 ) and AP( 3 ) contain a( 2, 1 ) and a( 3, 1 )
  74. C respectively, and so on.
  75. C Note that when DIAG = 'U' or 'u', the diagonal elements of
  76. C A are not referenced, but are assumed to be unity.
  77. C Unchanged on exit.
  78. C
  79. C X - COMPLEX array of dimension at least
  80. C ( 1 + ( n - 1 )*abs( INCX ) ).
  81. C Before entry, the incremented array X must contain the n
  82. C element vector x. On exit, X is overwritten with the
  83. C transformed vector x.
  84. C
  85. C INCX - INTEGER.
  86. C On entry, INCX specifies the increment for the elements of
  87. C X. INCX must not be zero.
  88. C Unchanged on exit.
  89. C
  90. C***REFERENCES Dongarra, J. J., Du Croz, J., Hammarling, S., and
  91. C Hanson, R. J. An extended set of Fortran basic linear
  92. C algebra subprograms. ACM TOMS, Vol. 14, No. 1,
  93. C pp. 1-17, March 1988.
  94. C***ROUTINES CALLED LSAME, XERBLA
  95. C***REVISION HISTORY (YYMMDD)
  96. C 861022 DATE WRITTEN
  97. C 910605 Modified to meet SLATEC prologue standards. Only comment
  98. C lines were modified. (BKS)
  99. C***END PROLOGUE CTPMV
  100. C .. Scalar Arguments ..
  101. INTEGER INCX, N
  102. CHARACTER*1 DIAG, TRANS, UPLO
  103. C .. Array Arguments ..
  104. COMPLEX AP( * ), X( * )
  105. C .. Parameters ..
  106. COMPLEX ZERO
  107. PARAMETER ( ZERO = ( 0.0E+0, 0.0E+0 ) )
  108. C .. Local Scalars ..
  109. COMPLEX TEMP
  110. INTEGER I, INFO, IX, J, JX, K, KK, KX
  111. LOGICAL NOCONJ, NOUNIT
  112. C .. External Functions ..
  113. LOGICAL LSAME
  114. EXTERNAL LSAME
  115. C .. External Subroutines ..
  116. EXTERNAL XERBLA
  117. C .. Intrinsic Functions ..
  118. INTRINSIC CONJG
  119. C***FIRST EXECUTABLE STATEMENT CTPMV
  120. C
  121. C Test the input parameters.
  122. C
  123. INFO = 0
  124. IF ( .NOT.LSAME( UPLO , 'U' ).AND.
  125. $ .NOT.LSAME( UPLO , 'L' ) )THEN
  126. INFO = 1
  127. ELSE IF( .NOT.LSAME( TRANS, 'N' ).AND.
  128. $ .NOT.LSAME( TRANS, 'T' ).AND.
  129. $ .NOT.LSAME( TRANS, 'C' ) )THEN
  130. INFO = 2
  131. ELSE IF( .NOT.LSAME( DIAG , 'U' ).AND.
  132. $ .NOT.LSAME( DIAG , 'N' ) )THEN
  133. INFO = 3
  134. ELSE IF( N.LT.0 )THEN
  135. INFO = 4
  136. ELSE IF( INCX.EQ.0 )THEN
  137. INFO = 7
  138. END IF
  139. IF( INFO.NE.0 )THEN
  140. CALL XERBLA( 'CTPMV ', INFO )
  141. RETURN
  142. END IF
  143. C
  144. C Quick return if possible.
  145. C
  146. IF( N.EQ.0 )
  147. $ RETURN
  148. C
  149. NOCONJ = LSAME( TRANS, 'T' )
  150. NOUNIT = LSAME( DIAG , 'N' )
  151. C
  152. C Set up the start point in X if the increment is not unity. This
  153. C will be ( N - 1 )*INCX too small for descending loops.
  154. C
  155. IF( INCX.LE.0 )THEN
  156. KX = 1 - ( N - 1 )*INCX
  157. ELSE IF( INCX.NE.1 )THEN
  158. KX = 1
  159. END IF
  160. C
  161. C Start the operations. In this version the elements of AP are
  162. C accessed sequentially with one pass through AP.
  163. C
  164. IF( LSAME( TRANS, 'N' ) )THEN
  165. C
  166. C Form x:= A*x.
  167. C
  168. IF( LSAME( UPLO, 'U' ) )THEN
  169. KK = 1
  170. IF( INCX.EQ.1 )THEN
  171. DO 20, J = 1, N
  172. IF( X( J ).NE.ZERO )THEN
  173. TEMP = X( J )
  174. K = KK
  175. DO 10, I = 1, J - 1
  176. X( I ) = X( I ) + TEMP*AP( K )
  177. K = K + 1
  178. 10 CONTINUE
  179. IF( NOUNIT )
  180. $ X( J ) = X( J )*AP( KK + J - 1 )
  181. END IF
  182. KK = KK + J
  183. 20 CONTINUE
  184. ELSE
  185. JX = KX
  186. DO 40, J = 1, N
  187. IF( X( JX ).NE.ZERO )THEN
  188. TEMP = X( JX )
  189. IX = KX
  190. DO 30, K = KK, KK + J - 2
  191. X( IX ) = X( IX ) + TEMP*AP( K )
  192. IX = IX + INCX
  193. 30 CONTINUE
  194. IF( NOUNIT )
  195. $ X( JX ) = X( JX )*AP( KK + J - 1 )
  196. END IF
  197. JX = JX + INCX
  198. KK = KK + J
  199. 40 CONTINUE
  200. END IF
  201. ELSE
  202. KK = ( N*( N + 1 ) )/2
  203. IF( INCX.EQ.1 )THEN
  204. DO 60, J = N, 1, -1
  205. IF( X( J ).NE.ZERO )THEN
  206. TEMP = X( J )
  207. K = KK
  208. DO 50, I = N, J + 1, -1
  209. X( I ) = X( I ) + TEMP*AP( K )
  210. K = K - 1
  211. 50 CONTINUE
  212. IF( NOUNIT )
  213. $ X( J ) = X( J )*AP( KK - N + J )
  214. END IF
  215. KK = KK - ( N - J + 1 )
  216. 60 CONTINUE
  217. ELSE
  218. KX = KX + ( N - 1 )*INCX
  219. JX = KX
  220. DO 80, J = N, 1, -1
  221. IF( X( JX ).NE.ZERO )THEN
  222. TEMP = X( JX )
  223. IX = KX
  224. DO 70, K = KK, KK - ( N - ( J + 1 ) ), -1
  225. X( IX ) = X( IX ) + TEMP*AP( K )
  226. IX = IX - INCX
  227. 70 CONTINUE
  228. IF( NOUNIT )
  229. $ X( JX ) = X( JX )*AP( KK - N + J )
  230. END IF
  231. JX = JX - INCX
  232. KK = KK - ( N - J + 1 )
  233. 80 CONTINUE
  234. END IF
  235. END IF
  236. ELSE
  237. C
  238. C Form x := A'*x or x := conjg( A' )*x.
  239. C
  240. IF( LSAME( UPLO, 'U' ) )THEN
  241. KK = ( N*( N + 1 ) )/2
  242. IF( INCX.EQ.1 )THEN
  243. DO 110, J = N, 1, -1
  244. TEMP = X( J )
  245. K = KK - 1
  246. IF( NOCONJ )THEN
  247. IF( NOUNIT )
  248. $ TEMP = TEMP*AP( KK )
  249. DO 90, I = J - 1, 1, -1
  250. TEMP = TEMP + AP( K )*X( I )
  251. K = K - 1
  252. 90 CONTINUE
  253. ELSE
  254. IF( NOUNIT )
  255. $ TEMP = TEMP*CONJG( AP( KK ) )
  256. DO 100, I = J - 1, 1, -1
  257. TEMP = TEMP + CONJG( AP( K ) )*X( I )
  258. K = K - 1
  259. 100 CONTINUE
  260. END IF
  261. X( J ) = TEMP
  262. KK = KK - J
  263. 110 CONTINUE
  264. ELSE
  265. JX = KX + ( N - 1 )*INCX
  266. DO 140, J = N, 1, -1
  267. TEMP = X( JX )
  268. IX = JX
  269. IF( NOCONJ )THEN
  270. IF( NOUNIT )
  271. $ TEMP = TEMP*AP( KK )
  272. DO 120, K = KK - 1, KK - J + 1, -1
  273. IX = IX - INCX
  274. TEMP = TEMP + AP( K )*X( IX )
  275. 120 CONTINUE
  276. ELSE
  277. IF( NOUNIT )
  278. $ TEMP = TEMP*CONJG( AP( KK ) )
  279. DO 130, K = KK - 1, KK - J + 1, -1
  280. IX = IX - INCX
  281. TEMP = TEMP + CONJG( AP( K ) )*X( IX )
  282. 130 CONTINUE
  283. END IF
  284. X( JX ) = TEMP
  285. JX = JX - INCX
  286. KK = KK - J
  287. 140 CONTINUE
  288. END IF
  289. ELSE
  290. KK = 1
  291. IF( INCX.EQ.1 )THEN
  292. DO 170, J = 1, N
  293. TEMP = X( J )
  294. K = KK + 1
  295. IF( NOCONJ )THEN
  296. IF( NOUNIT )
  297. $ TEMP = TEMP*AP( KK )
  298. DO 150, I = J + 1, N
  299. TEMP = TEMP + AP( K )*X( I )
  300. K = K + 1
  301. 150 CONTINUE
  302. ELSE
  303. IF( NOUNIT )
  304. $ TEMP = TEMP*CONJG( AP( KK ) )
  305. DO 160, I = J + 1, N
  306. TEMP = TEMP + CONJG( AP( K ) )*X( I )
  307. K = K + 1
  308. 160 CONTINUE
  309. END IF
  310. X( J ) = TEMP
  311. KK = KK + ( N - J + 1 )
  312. 170 CONTINUE
  313. ELSE
  314. JX = KX
  315. DO 200, J = 1, N
  316. TEMP = X( JX )
  317. IX = JX
  318. IF( NOCONJ )THEN
  319. IF( NOUNIT )
  320. $ TEMP = TEMP*AP( KK )
  321. DO 180, K = KK + 1, KK + N - J
  322. IX = IX + INCX
  323. TEMP = TEMP + AP( K )*X( IX )
  324. 180 CONTINUE
  325. ELSE
  326. IF( NOUNIT )
  327. $ TEMP = TEMP*CONJG( AP( KK ) )
  328. DO 190, K = KK + 1, KK + N - J
  329. IX = IX + INCX
  330. TEMP = TEMP + CONJG( AP( K ) )*X( IX )
  331. 190 CONTINUE
  332. END IF
  333. X( JX ) = TEMP
  334. JX = JX + INCX
  335. KK = KK + ( N - J + 1 )
  336. 200 CONTINUE
  337. END IF
  338. END IF
  339. END IF
  340. C
  341. RETURN
  342. C
  343. C End of CTPMV .
  344. C
  345. END