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.

dsyrf.f 6.0 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197
  1. SUBROUTINE DSYRF ( UPLO, N, ALPHA, X, INCX, A, LDA )
  2. * .. Scalar Arguments ..
  3. DOUBLE PRECISION ALPHA
  4. INTEGER INCX, LDA, N
  5. CHARACTER*1 UPLO
  6. * .. Array Arguments ..
  7. DOUBLE PRECISION A( LDA, * ), X( * )
  8. * ..
  9. *
  10. * Purpose
  11. * =======
  12. *
  13. * DSYR performs the symmetric rank 1 operation
  14. *
  15. * A := alpha*x*x' + A,
  16. *
  17. * where alpha is a real scalar, x is an n element vector and A is an
  18. * n by n symmetric 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 - DOUBLE PRECISION.
  42. * On entry, ALPHA specifies the scalar alpha.
  43. * Unchanged on exit.
  44. *
  45. * X - DOUBLE PRECISION 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. * A - DOUBLE PRECISION array of DIMENSION ( LDA, n ).
  57. * Before entry with UPLO = 'U' or 'u', the leading n by n
  58. * upper triangular part of the array A must contain the upper
  59. * triangular part of the symmetric matrix and the strictly
  60. * lower triangular part of A is not referenced. On exit, the
  61. * upper triangular part of the array A is overwritten by the
  62. * upper triangular part of the updated matrix.
  63. * Before entry with UPLO = 'L' or 'l', the leading n by n
  64. * lower triangular part of the array A must contain the lower
  65. * triangular part of the symmetric matrix and the strictly
  66. * upper triangular part of A is not referenced. On exit, the
  67. * lower triangular part of the array A is overwritten by the
  68. * lower triangular part of the updated matrix.
  69. *
  70. * LDA - INTEGER.
  71. * On entry, LDA specifies the first dimension of A as declared
  72. * in the calling (sub) program. LDA must be at least
  73. * max( 1, n ).
  74. * Unchanged on exit.
  75. *
  76. *
  77. * Level 2 Blas routine.
  78. *
  79. * -- Written on 22-October-1986.
  80. * Jack Dongarra, Argonne National Lab.
  81. * Jeremy Du Croz, Nag Central Office.
  82. * Sven Hammarling, Nag Central Office.
  83. * Richard Hanson, Sandia National Labs.
  84. *
  85. *
  86. * .. Parameters ..
  87. DOUBLE PRECISION ZERO
  88. PARAMETER ( ZERO = 0.0D+0 )
  89. * .. Local Scalars ..
  90. DOUBLE PRECISION TEMP
  91. INTEGER I, INFO, IX, J, JX, KX
  92. * .. External Functions ..
  93. LOGICAL LSAME
  94. EXTERNAL LSAME
  95. * .. External Subroutines ..
  96. EXTERNAL XERBLA
  97. * .. Intrinsic Functions ..
  98. INTRINSIC MAX
  99. * ..
  100. * .. Executable Statements ..
  101. *
  102. * Test the input parameters.
  103. *
  104. INFO = 0
  105. IF ( .NOT.LSAME( UPLO, 'U' ).AND.
  106. $ .NOT.LSAME( UPLO, 'L' ) )THEN
  107. INFO = 1
  108. ELSE IF( N.LT.0 )THEN
  109. INFO = 2
  110. ELSE IF( INCX.EQ.0 )THEN
  111. INFO = 5
  112. ELSE IF( LDA.LT.MAX( 1, N ) )THEN
  113. INFO = 7
  114. END IF
  115. IF( INFO.NE.0 )THEN
  116. CALL XERBLA( 'DSYR ', INFO )
  117. RETURN
  118. END IF
  119. *
  120. * Quick return if possible.
  121. *
  122. IF( ( N.EQ.0 ).OR.( ALPHA.EQ.ZERO ) )
  123. $ RETURN
  124. *
  125. * Set the start point in X if the increment is not unity.
  126. *
  127. IF( INCX.LE.0 )THEN
  128. KX = 1 - ( N - 1 )*INCX
  129. ELSE IF( INCX.NE.1 )THEN
  130. KX = 1
  131. END IF
  132. *
  133. * Start the operations. In this version the elements of A are
  134. * accessed sequentially with one pass through the triangular part
  135. * of A.
  136. *
  137. IF( LSAME( UPLO, 'U' ) )THEN
  138. *
  139. * Form A when A is stored in upper triangle.
  140. *
  141. IF( INCX.EQ.1 )THEN
  142. DO 20, J = 1, N
  143. IF( X( J ).NE.ZERO )THEN
  144. TEMP = ALPHA*X( J )
  145. DO 10, I = 1, J
  146. A( I, J ) = A( I, J ) + X( I )*TEMP
  147. 10 CONTINUE
  148. END IF
  149. 20 CONTINUE
  150. ELSE
  151. JX = KX
  152. DO 40, J = 1, N
  153. IF( X( JX ).NE.ZERO )THEN
  154. TEMP = ALPHA*X( JX )
  155. IX = KX
  156. DO 30, I = 1, J
  157. A( I, J ) = A( I, J ) + X( IX )*TEMP
  158. IX = IX + INCX
  159. 30 CONTINUE
  160. END IF
  161. JX = JX + INCX
  162. 40 CONTINUE
  163. END IF
  164. ELSE
  165. *
  166. * Form A when A is stored in lower triangle.
  167. *
  168. IF( INCX.EQ.1 )THEN
  169. DO 60, J = 1, N
  170. IF( X( J ).NE.ZERO )THEN
  171. TEMP = ALPHA*X( J )
  172. DO 50, I = J, N
  173. A( I, J ) = A( I, J ) + X( I )*TEMP
  174. 50 CONTINUE
  175. END IF
  176. 60 CONTINUE
  177. ELSE
  178. JX = KX
  179. DO 80, J = 1, N
  180. IF( X( JX ).NE.ZERO )THEN
  181. TEMP = ALPHA*X( JX )
  182. IX = JX
  183. DO 70, I = J, N
  184. A( I, J ) = A( I, J ) + X( IX )*TEMP
  185. IX = IX + INCX
  186. 70 CONTINUE
  187. END IF
  188. JX = JX + INCX
  189. 80 CONTINUE
  190. END IF
  191. END IF
  192. *
  193. RETURN
  194. *
  195. * End of DSYR .
  196. *
  197. END