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.

cher2f.f 8.4 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249
  1. SUBROUTINE CHER2F ( UPLO, N, ALPHA, X, INCX, Y, INCY, A, LDA )
  2. * .. Scalar Arguments ..
  3. COMPLEX ALPHA
  4. INTEGER INCX, INCY, LDA, N
  5. CHARACTER*1 UPLO
  6. * .. Array Arguments ..
  7. COMPLEX A( LDA, * ), X( * ), Y( * )
  8. * ..
  9. *
  10. * Purpose
  11. * =======
  12. *
  13. * CHER2 performs the hermitian rank 2 operation
  14. *
  15. * A := alpha*x*conjg( y' ) + conjg( alpha )*y*conjg( x' ) + A,
  16. *
  17. * where alpha is a scalar, x and y are n element vectors and A is an n
  18. * by n hermitian matrix.
  19. *
  20. * Parameters
  21. * ==========
  22. *
  23. * UPLO - CHARACTER*1.
  24. * On entry, UPLO specifies whether the upper or lower
  25. * triangular part of the array A is to be referenced as
  26. * follows:
  27. *
  28. * UPLO = 'U' or 'u' Only the upper triangular part of A
  29. * is to be referenced.
  30. *
  31. * UPLO = 'L' or 'l' Only the lower triangular part of A
  32. * is to be referenced.
  33. *
  34. * Unchanged on exit.
  35. *
  36. * N - INTEGER.
  37. * On entry, N specifies the order of the matrix A.
  38. * N must be at least zero.
  39. * Unchanged on exit.
  40. *
  41. * ALPHA - COMPLEX .
  42. * On entry, ALPHA specifies the scalar alpha.
  43. * Unchanged on exit.
  44. *
  45. * X - COMPLEX array of dimension at least
  46. * ( 1 + ( n - 1 )*abs( INCX ) ).
  47. * Before entry, the incremented array X must contain the n
  48. * element vector x.
  49. * Unchanged on exit.
  50. *
  51. * INCX - INTEGER.
  52. * On entry, INCX specifies the increment for the elements of
  53. * X. INCX must not be zero.
  54. * Unchanged on exit.
  55. *
  56. * Y - COMPLEX array of dimension at least
  57. * ( 1 + ( n - 1 )*abs( INCY ) ).
  58. * Before entry, the incremented array Y must contain the n
  59. * element vector y.
  60. * Unchanged on exit.
  61. *
  62. * INCY - INTEGER.
  63. * On entry, INCY specifies the increment for the elements of
  64. * Y. INCY must not be zero.
  65. * Unchanged on exit.
  66. *
  67. * A - COMPLEX array of DIMENSION ( LDA, n ).
  68. * Before entry with UPLO = 'U' or 'u', the leading n by n
  69. * upper triangular part of the array A must contain the upper
  70. * triangular part of the hermitian matrix and the strictly
  71. * lower triangular part of A is not referenced. On exit, the
  72. * upper triangular part of the array A is overwritten by the
  73. * upper triangular part of the updated matrix.
  74. * Before entry with UPLO = 'L' or 'l', the leading n by n
  75. * lower triangular part of the array A must contain the lower
  76. * triangular part of the hermitian matrix and the strictly
  77. * upper triangular part of A is not referenced. On exit, the
  78. * lower triangular part of the array A is overwritten by the
  79. * lower triangular part of the updated matrix.
  80. * Note that the imaginary parts of the diagonal elements need
  81. * not be set, they are assumed to be zero, and on exit they
  82. * are set to zero.
  83. *
  84. * LDA - INTEGER.
  85. * On entry, LDA specifies the first dimension of A as declared
  86. * in the calling (sub) program. LDA must be at least
  87. * max( 1, n ).
  88. * Unchanged on exit.
  89. *
  90. *
  91. * Level 2 Blas routine.
  92. *
  93. * -- Written on 22-October-1986.
  94. * Jack Dongarra, Argonne National Lab.
  95. * Jeremy Du Croz, Nag Central Office.
  96. * Sven Hammarling, Nag Central Office.
  97. * Richard Hanson, Sandia National Labs.
  98. *
  99. *
  100. * .. Parameters ..
  101. COMPLEX ZERO
  102. PARAMETER ( ZERO = ( 0.0E+0, 0.0E+0 ) )
  103. * .. Local Scalars ..
  104. COMPLEX TEMP1, TEMP2
  105. INTEGER I, INFO, IX, IY, J, JX, JY, KX, KY
  106. * .. External Functions ..
  107. LOGICAL LSAME
  108. EXTERNAL LSAME
  109. * .. External Subroutines ..
  110. EXTERNAL XERBLA
  111. * .. Intrinsic Functions ..
  112. INTRINSIC CONJG, MAX, REAL
  113. * ..
  114. * .. Executable Statements ..
  115. *
  116. * Test the input parameters.
  117. *
  118. INFO = 0
  119. IF ( .NOT.LSAME( UPLO, 'U' ).AND.
  120. $ .NOT.LSAME( UPLO, 'L' ) )THEN
  121. INFO = 1
  122. ELSE IF( N.LT.0 )THEN
  123. INFO = 2
  124. ELSE IF( INCX.EQ.0 )THEN
  125. INFO = 5
  126. ELSE IF( INCY.EQ.0 )THEN
  127. INFO = 7
  128. ELSE IF( LDA.LT.MAX( 1, N ) )THEN
  129. INFO = 9
  130. END IF
  131. IF( INFO.NE.0 )THEN
  132. CALL XERBLA( 'CHER2 ', INFO )
  133. RETURN
  134. END IF
  135. *
  136. * Quick return if possible.
  137. *
  138. IF( ( N.EQ.0 ).OR.( ALPHA.EQ.ZERO ) )
  139. $ RETURN
  140. *
  141. * Set up the start points in X and Y if the increments are not both
  142. * unity.
  143. *
  144. IF( ( INCX.NE.1 ).OR.( INCY.NE.1 ) )THEN
  145. IF( INCX.GT.0 )THEN
  146. KX = 1
  147. ELSE
  148. KX = 1 - ( N - 1 )*INCX
  149. END IF
  150. IF( INCY.GT.0 )THEN
  151. KY = 1
  152. ELSE
  153. KY = 1 - ( N - 1 )*INCY
  154. END IF
  155. JX = KX
  156. JY = KY
  157. END IF
  158. *
  159. * Start the operations. In this version the elements of A are
  160. * accessed sequentially with one pass through the triangular part
  161. * of A.
  162. *
  163. IF( LSAME( UPLO, 'U' ) )THEN
  164. *
  165. * Form A when A is stored in the upper triangle.
  166. *
  167. IF( ( INCX.EQ.1 ).AND.( INCY.EQ.1 ) )THEN
  168. DO 20, J = 1, N
  169. IF( ( X( J ).NE.ZERO ).OR.( Y( J ).NE.ZERO ) )THEN
  170. TEMP1 = ALPHA*CONJG( Y( J ) )
  171. TEMP2 = CONJG( ALPHA*X( J ) )
  172. DO 10, I = 1, J - 1
  173. A( I, J ) = A( I, J ) + X( I )*TEMP1 + Y( I )*TEMP2
  174. 10 CONTINUE
  175. A( J, J ) = REAL( A( J, J ) ) +
  176. $ REAL( X( J )*TEMP1 + Y( J )*TEMP2 )
  177. ELSE
  178. A( J, J ) = REAL( A( J, J ) )
  179. END IF
  180. 20 CONTINUE
  181. ELSE
  182. DO 40, J = 1, N
  183. IF( ( X( JX ).NE.ZERO ).OR.( Y( JY ).NE.ZERO ) )THEN
  184. TEMP1 = ALPHA*CONJG( Y( JY ) )
  185. TEMP2 = CONJG( ALPHA*X( JX ) )
  186. IX = KX
  187. IY = KY
  188. DO 30, I = 1, J - 1
  189. A( I, J ) = A( I, J ) + X( IX )*TEMP1
  190. $ + Y( IY )*TEMP2
  191. IX = IX + INCX
  192. IY = IY + INCY
  193. 30 CONTINUE
  194. A( J, J ) = REAL( A( J, J ) ) +
  195. $ REAL( X( JX )*TEMP1 + Y( JY )*TEMP2 )
  196. ELSE
  197. A( J, J ) = REAL( A( J, J ) )
  198. END IF
  199. JX = JX + INCX
  200. JY = JY + INCY
  201. 40 CONTINUE
  202. END IF
  203. ELSE
  204. *
  205. * Form A when A is stored in the lower triangle.
  206. *
  207. IF( ( INCX.EQ.1 ).AND.( INCY.EQ.1 ) )THEN
  208. DO 60, J = 1, N
  209. IF( ( X( J ).NE.ZERO ).OR.( Y( J ).NE.ZERO ) )THEN
  210. TEMP1 = ALPHA*CONJG( Y( J ) )
  211. TEMP2 = CONJG( ALPHA*X( J ) )
  212. A( J, J ) = REAL( A( J, J ) ) +
  213. $ REAL( X( J )*TEMP1 + Y( J )*TEMP2 )
  214. DO 50, I = J + 1, N
  215. A( I, J ) = A( I, J ) + X( I )*TEMP1 + Y( I )*TEMP2
  216. 50 CONTINUE
  217. ELSE
  218. A( J, J ) = REAL( A( J, J ) )
  219. END IF
  220. 60 CONTINUE
  221. ELSE
  222. DO 80, J = 1, N
  223. IF( ( X( JX ).NE.ZERO ).OR.( Y( JY ).NE.ZERO ) )THEN
  224. TEMP1 = ALPHA*CONJG( Y( JY ) )
  225. TEMP2 = CONJG( ALPHA*X( JX ) )
  226. A( J, J ) = REAL( A( J, J ) ) +
  227. $ REAL( X( JX )*TEMP1 + Y( JY )*TEMP2 )
  228. IX = JX
  229. IY = JY
  230. DO 70, I = J + 1, N
  231. IX = IX + INCX
  232. IY = IY + INCY
  233. A( I, J ) = A( I, J ) + X( IX )*TEMP1
  234. $ + Y( IY )*TEMP2
  235. 70 CONTINUE
  236. ELSE
  237. A( J, J ) = REAL( A( J, J ) )
  238. END IF
  239. JX = JX + INCX
  240. JY = JY + INCY
  241. 80 CONTINUE
  242. END IF
  243. END IF
  244. *
  245. RETURN
  246. *
  247. * End of CHER2 .
  248. *
  249. END