csyrk.f 9.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299
  1. *DECK CSYRK
  2. SUBROUTINE CSYRK (UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
  3. C***BEGIN PROLOGUE CSYRK
  4. C***PURPOSE Perform symmetric rank k update of a complex symmetric
  5. C matrix.
  6. C***LIBRARY SLATEC (BLAS)
  7. C***CATEGORY D1B6
  8. C***TYPE COMPLEX (SSYRK-S, DSYRK-D, CSYRK-C)
  9. C***KEYWORDS LEVEL 3 BLAS, LINEAR ALGEBRA
  10. C***AUTHOR Dongarra, J., (ANL)
  11. C Duff, I., (AERE)
  12. C Du Croz, J., (NAG)
  13. C Hammarling, S. (NAG)
  14. C***DESCRIPTION
  15. C
  16. C CSYRK performs one of the symmetric rank k operations
  17. C
  18. C C := alpha*A*A' + beta*C,
  19. C
  20. C or
  21. C
  22. C C := alpha*A'*A + beta*C,
  23. C
  24. C where alpha and beta are scalars, C is an n by n symmetric matrix
  25. C and A is an n by k matrix in the first case and a k by n matrix
  26. C in the second case.
  27. C
  28. C Parameters
  29. C ==========
  30. C
  31. C UPLO - CHARACTER*1.
  32. C On entry, UPLO specifies whether the upper or lower
  33. C triangular part of the array C is to be referenced as
  34. C follows:
  35. C
  36. C UPLO = 'U' or 'u' Only the upper triangular part of C
  37. C is to be referenced.
  38. C
  39. C UPLO = 'L' or 'l' Only the lower triangular part of C
  40. C is to be referenced.
  41. C
  42. C Unchanged on exit.
  43. C
  44. C TRANS - CHARACTER*1.
  45. C On entry, TRANS specifies the operation to be performed as
  46. C follows:
  47. C
  48. C TRANS = 'N' or 'n' C := alpha*A*A' + beta*C.
  49. C
  50. C TRANS = 'T' or 't' C := alpha*A'*A + beta*C.
  51. C
  52. C Unchanged on exit.
  53. C
  54. C N - INTEGER.
  55. C On entry, N specifies the order of the matrix C. N must be
  56. C at least zero.
  57. C Unchanged on exit.
  58. C
  59. C K - INTEGER.
  60. C On entry with TRANS = 'N' or 'n', K specifies the number
  61. C of columns of the matrix A, and on entry with
  62. C TRANS = 'T' or 't', K specifies the number of rows of the
  63. C matrix A. K must be at least zero.
  64. C Unchanged on exit.
  65. C
  66. C ALPHA - COMPLEX .
  67. C On entry, ALPHA specifies the scalar alpha.
  68. C Unchanged on exit.
  69. C
  70. C A - COMPLEX array of DIMENSION ( LDA, ka ), where ka is
  71. C k when TRANS = 'N' or 'n', and is n otherwise.
  72. C Before entry with TRANS = 'N' or 'n', the leading n by k
  73. C part of the array A must contain the matrix A, otherwise
  74. C the leading k by n part of the array A must contain the
  75. C matrix A.
  76. C Unchanged on exit.
  77. C
  78. C LDA - INTEGER.
  79. C On entry, LDA specifies the first dimension of A as declared
  80. C in the calling (sub) program. When TRANS = 'N' or 'n'
  81. C then LDA must be at least max( 1, n ), otherwise LDA must
  82. C be at least max( 1, k ).
  83. C Unchanged on exit.
  84. C
  85. C BETA - COMPLEX .
  86. C On entry, BETA specifies the scalar beta.
  87. C Unchanged on exit.
  88. C
  89. C C - COMPLEX array of DIMENSION ( LDC, n ).
  90. C Before entry with UPLO = 'U' or 'u', the leading n by n
  91. C upper triangular part of the array C must contain the upper
  92. C triangular part of the symmetric matrix and the strictly
  93. C lower triangular part of C is not referenced. On exit, the
  94. C upper triangular part of the array C is overwritten by the
  95. C upper triangular part of the updated matrix.
  96. C Before entry with UPLO = 'L' or 'l', the leading n by n
  97. C lower triangular part of the array C must contain the lower
  98. C triangular part of the symmetric matrix and the strictly
  99. C upper triangular part of C is not referenced. On exit, the
  100. C lower triangular part of the array C is overwritten by the
  101. C lower triangular part of the updated matrix.
  102. C
  103. C LDC - INTEGER.
  104. C On entry, LDC specifies the first dimension of C as declared
  105. C in the calling (sub) program. LDC must be at least
  106. C max( 1, n ).
  107. C Unchanged on exit.
  108. C
  109. C***REFERENCES Dongarra, J., Du Croz, J., Duff, I., and Hammarling, S.
  110. C A set of level 3 basic linear algebra subprograms.
  111. C ACM TOMS, Vol. 16, No. 1, pp. 1-17, March 1990.
  112. C***ROUTINES CALLED LSAME, XERBLA
  113. C***REVISION HISTORY (YYMMDD)
  114. C 890208 DATE WRITTEN
  115. C 910605 Modified to meet SLATEC prologue standards. Only comment
  116. C lines were modified. (BKS)
  117. C***END PROLOGUE CSYRK
  118. C .. Scalar Arguments ..
  119. CHARACTER*1 UPLO, TRANS
  120. INTEGER N, K, LDA, LDC
  121. COMPLEX ALPHA, BETA
  122. C .. Array Arguments ..
  123. COMPLEX A( LDA, * ), C( LDC, * )
  124. C .. External Functions ..
  125. LOGICAL LSAME
  126. EXTERNAL LSAME
  127. C .. External Subroutines ..
  128. EXTERNAL XERBLA
  129. C .. Intrinsic Functions ..
  130. INTRINSIC MAX
  131. C .. Local Scalars ..
  132. LOGICAL UPPER
  133. INTEGER I, INFO, J, L, NROWA
  134. COMPLEX TEMP
  135. C .. Parameters ..
  136. COMPLEX ONE
  137. PARAMETER ( ONE = ( 1.0E+0, 0.0E+0 ) )
  138. COMPLEX ZERO
  139. PARAMETER ( ZERO = ( 0.0E+0, 0.0E+0 ) )
  140. C***FIRST EXECUTABLE STATEMENT CSYRK
  141. C
  142. C Test the input parameters.
  143. C
  144. IF( LSAME( TRANS, 'N' ) )THEN
  145. NROWA = N
  146. ELSE
  147. NROWA = K
  148. END IF
  149. UPPER = LSAME( UPLO, 'U' )
  150. C
  151. INFO = 0
  152. IF( ( .NOT.UPPER ).AND.
  153. $ ( .NOT.LSAME( UPLO , 'L' ) ) )THEN
  154. INFO = 1
  155. ELSE IF( ( .NOT.LSAME( TRANS, 'N' ) ).AND.
  156. $ ( .NOT.LSAME( TRANS, 'T' ) ) )THEN
  157. INFO = 2
  158. ELSE IF( N .LT.0 )THEN
  159. INFO = 3
  160. ELSE IF( K .LT.0 )THEN
  161. INFO = 4
  162. ELSE IF( LDA.LT.MAX( 1, NROWA ) )THEN
  163. INFO = 7
  164. ELSE IF( LDC.LT.MAX( 1, N ) )THEN
  165. INFO = 10
  166. END IF
  167. IF( INFO.NE.0 )THEN
  168. CALL XERBLA( 'CSYRK ', INFO )
  169. RETURN
  170. END IF
  171. C
  172. C Quick return if possible.
  173. C
  174. IF( ( N.EQ.0 ).OR.
  175. $ ( ( ( ALPHA.EQ.ZERO ).OR.( K.EQ.0 ) ).AND.( BETA.EQ.ONE ) ) )
  176. $ RETURN
  177. C
  178. C And when alpha.eq.zero.
  179. C
  180. IF( ALPHA.EQ.ZERO )THEN
  181. IF( UPPER )THEN
  182. IF( BETA.EQ.ZERO )THEN
  183. DO 20, J = 1, N
  184. DO 10, I = 1, J
  185. C( I, J ) = ZERO
  186. 10 CONTINUE
  187. 20 CONTINUE
  188. ELSE
  189. DO 40, J = 1, N
  190. DO 30, I = 1, J
  191. C( I, J ) = BETA*C( I, J )
  192. 30 CONTINUE
  193. 40 CONTINUE
  194. END IF
  195. ELSE
  196. IF( BETA.EQ.ZERO )THEN
  197. DO 60, J = 1, N
  198. DO 50, I = J, N
  199. C( I, J ) = ZERO
  200. 50 CONTINUE
  201. 60 CONTINUE
  202. ELSE
  203. DO 80, J = 1, N
  204. DO 70, I = J, N
  205. C( I, J ) = BETA*C( I, J )
  206. 70 CONTINUE
  207. 80 CONTINUE
  208. END IF
  209. END IF
  210. RETURN
  211. END IF
  212. C
  213. C Start the operations.
  214. C
  215. IF( LSAME( TRANS, 'N' ) )THEN
  216. C
  217. C Form C := alpha*A*A' + beta*C.
  218. C
  219. IF( UPPER )THEN
  220. DO 130, J = 1, N
  221. IF( BETA.EQ.ZERO )THEN
  222. DO 90, I = 1, J
  223. C( I, J ) = ZERO
  224. 90 CONTINUE
  225. ELSE IF( BETA.NE.ONE )THEN
  226. DO 100, I = 1, J
  227. C( I, J ) = BETA*C( I, J )
  228. 100 CONTINUE
  229. END IF
  230. DO 120, L = 1, K
  231. IF( A( J, L ).NE.ZERO )THEN
  232. TEMP = ALPHA*A( J, L )
  233. DO 110, I = 1, J
  234. C( I, J ) = C( I, J ) + TEMP*A( I, L )
  235. 110 CONTINUE
  236. END IF
  237. 120 CONTINUE
  238. 130 CONTINUE
  239. ELSE
  240. DO 180, J = 1, N
  241. IF( BETA.EQ.ZERO )THEN
  242. DO 140, I = J, N
  243. C( I, J ) = ZERO
  244. 140 CONTINUE
  245. ELSE IF( BETA.NE.ONE )THEN
  246. DO 150, I = J, N
  247. C( I, J ) = BETA*C( I, J )
  248. 150 CONTINUE
  249. END IF
  250. DO 170, L = 1, K
  251. IF( A( J, L ).NE.ZERO )THEN
  252. TEMP = ALPHA*A( J, L )
  253. DO 160, I = J, N
  254. C( I, J ) = C( I, J ) + TEMP*A( I, L )
  255. 160 CONTINUE
  256. END IF
  257. 170 CONTINUE
  258. 180 CONTINUE
  259. END IF
  260. ELSE
  261. C
  262. C Form C := alpha*A'*A + beta*C.
  263. C
  264. IF( UPPER )THEN
  265. DO 210, J = 1, N
  266. DO 200, I = 1, J
  267. TEMP = ZERO
  268. DO 190, L = 1, K
  269. TEMP = TEMP + A( L, I )*A( L, J )
  270. 190 CONTINUE
  271. IF( BETA.EQ.ZERO )THEN
  272. C( I, J ) = ALPHA*TEMP
  273. ELSE
  274. C( I, J ) = ALPHA*TEMP + BETA*C( I, J )
  275. END IF
  276. 200 CONTINUE
  277. 210 CONTINUE
  278. ELSE
  279. DO 240, J = 1, N
  280. DO 230, I = J, N
  281. TEMP = ZERO
  282. DO 220, L = 1, K
  283. TEMP = TEMP + A( L, I )*A( L, J )
  284. 220 CONTINUE
  285. IF( BETA.EQ.ZERO )THEN
  286. C( I, J ) = ALPHA*TEMP
  287. ELSE
  288. C( I, J ) = ALPHA*TEMP + BETA*C( I, J )
  289. END IF
  290. 230 CONTINUE
  291. 240 CONTINUE
  292. END IF
  293. END IF
  294. C
  295. RETURN
  296. C
  297. C End of CSYRK .
  298. C
  299. END