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.

daxpy.f 3.5 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152
  1. *> \brief \b DAXPY
  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 DAXPY(N,DA,DX,INCX,DY,INCY)
  12. *
  13. * .. Scalar Arguments ..
  14. * DOUBLE PRECISION DA
  15. * INTEGER INCX,INCY,N
  16. * ..
  17. * .. Array Arguments ..
  18. * DOUBLE PRECISION DX(*),DY(*)
  19. * ..
  20. *
  21. *
  22. *> \par Purpose:
  23. * =============
  24. *>
  25. *> \verbatim
  26. *>
  27. *> DAXPY constant times a vector plus a vector.
  28. *> uses unrolled loops for increments equal to one.
  29. *> \endverbatim
  30. *
  31. * Arguments:
  32. * ==========
  33. *
  34. *> \param[in] N
  35. *> \verbatim
  36. *> N is INTEGER
  37. *> number of elements in input vector(s)
  38. *> \endverbatim
  39. *>
  40. *> \param[in] DA
  41. *> \verbatim
  42. *> DA is DOUBLE PRECISION
  43. *> On entry, DA specifies the scalar alpha.
  44. *> \endverbatim
  45. *>
  46. *> \param[in] DX
  47. *> \verbatim
  48. *> DX is DOUBLE PRECISION array, dimension ( 1 + ( N - 1 )*abs( INCX ) )
  49. *> \endverbatim
  50. *>
  51. *> \param[in] INCX
  52. *> \verbatim
  53. *> INCX is INTEGER
  54. *> storage spacing between elements of DX
  55. *> \endverbatim
  56. *>
  57. *> \param[in,out] DY
  58. *> \verbatim
  59. *> DY is DOUBLE PRECISION array, dimension ( 1 + ( N - 1 )*abs( INCY ) )
  60. *> \endverbatim
  61. *>
  62. *> \param[in] INCY
  63. *> \verbatim
  64. *> INCY is INTEGER
  65. *> storage spacing between elements of DY
  66. *> \endverbatim
  67. *
  68. * Authors:
  69. * ========
  70. *
  71. *> \author Univ. of Tennessee
  72. *> \author Univ. of California Berkeley
  73. *> \author Univ. of Colorado Denver
  74. *> \author NAG Ltd.
  75. *
  76. *> \date November 2017
  77. *
  78. *> \ingroup double_blas_level1
  79. *
  80. *> \par Further Details:
  81. * =====================
  82. *>
  83. *> \verbatim
  84. *>
  85. *> jack dongarra, linpack, 3/11/78.
  86. *> modified 12/3/93, array(1) declarations changed to array(*)
  87. *> \endverbatim
  88. *>
  89. * =====================================================================
  90. SUBROUTINE DAXPY(N,DA,DX,INCX,DY,INCY)
  91. *
  92. * -- Reference BLAS level1 routine (version 3.8.0) --
  93. * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
  94. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  95. * November 2017
  96. *
  97. * .. Scalar Arguments ..
  98. DOUBLE PRECISION DA
  99. INTEGER INCX,INCY,N
  100. * ..
  101. * .. Array Arguments ..
  102. DOUBLE PRECISION DX(*),DY(*)
  103. * ..
  104. *
  105. * =====================================================================
  106. *
  107. * .. Local Scalars ..
  108. INTEGER I,IX,IY,M,MP1
  109. * ..
  110. * .. Intrinsic Functions ..
  111. INTRINSIC MOD
  112. * ..
  113. IF (N.LE.0) RETURN
  114. IF (DA.EQ.0.0d0) RETURN
  115. IF (INCX.EQ.1 .AND. INCY.EQ.1) THEN
  116. *
  117. * code for both increments equal to 1
  118. *
  119. *
  120. * clean-up loop
  121. *
  122. M = MOD(N,4)
  123. IF (M.NE.0) THEN
  124. DO I = 1,M
  125. DY(I) = DY(I) + DA*DX(I)
  126. END DO
  127. END IF
  128. IF (N.LT.4) RETURN
  129. MP1 = M + 1
  130. DO I = MP1,N,4
  131. DY(I) = DY(I) + DA*DX(I)
  132. DY(I+1) = DY(I+1) + DA*DX(I+1)
  133. DY(I+2) = DY(I+2) + DA*DX(I+2)
  134. DY(I+3) = DY(I+3) + DA*DX(I+3)
  135. END DO
  136. ELSE
  137. *
  138. * code for unequal increments or equal increments
  139. * not equal to 1
  140. *
  141. IX = 1
  142. IY = 1
  143. IF (INCX.LT.0) IX = (-N+1)*INCX + 1
  144. IF (INCY.LT.0) IY = (-N+1)*INCY + 1
  145. DO I = 1,N
  146. DY(IY) = DY(IY) + DA*DX(IX)
  147. IX = IX + INCX
  148. IY = IY + INCY
  149. END DO
  150. END IF
  151. RETURN
  152. END