daxpy Subroutine

public subroutine daxpy(n, da, dx, incx, dy, incy)

DAXPY constant times a vector plus a vector. uses unrolled loops for increments equal to one.

Arguments

Type IntentOptional Attributes Name
integer(kind=ip) :: n
real(kind=wp) :: da
real(kind=wp) :: dx(*)
integer(kind=ip) :: incx
real(kind=wp) :: dy(*)
integer(kind=ip) :: incy

Called by

proc~~daxpy~~CalledByGraph proc~daxpy bspline_blas_module::daxpy proc~dfcmn bspline_defc_module::dfcmn proc~dfcmn->proc~daxpy proc~dlsei bspline_defc_module::dlsei proc~dfcmn->proc~dlsei proc~dbndac bspline_defc_module::dbndac proc~dfcmn->proc~dbndac proc~dh12 bspline_defc_module::dh12 proc~dh12->proc~daxpy proc~dlsei->proc~daxpy proc~dlsei->proc~dh12 proc~dlsi bspline_defc_module::dlsi proc~dlsei->proc~dlsi proc~dlsi->proc~daxpy proc~dlsi->proc~dh12 proc~dhfti bspline_defc_module::dhfti proc~dlsi->proc~dhfti proc~dlpdp bspline_defc_module::dlpdp proc~dlsi->proc~dlpdp proc~dwnlsm bspline_defc_module::dwnlsm proc~dwnlsm->proc~daxpy proc~dwnlsm->proc~dh12 proc~dwnlit bspline_defc_module::dwnlit proc~dwnlsm->proc~dwnlit proc~dbndac->proc~dh12 proc~dfc bspline_defc_module::dfc proc~dfc->proc~dfcmn proc~dhfti->proc~dh12 proc~dwnlit->proc~dh12 proc~dwnnls bspline_defc_module::dwnnls proc~dwnnls->proc~dwnlsm proc~defcmn bspline_defc_module::defcmn proc~defcmn->proc~dbndac proc~dlpdp->proc~dwnnls proc~defc bspline_defc_module::defc proc~defc->proc~defcmn

Source Code

    subroutine daxpy(n, da, dx, incx, dy, incy)
         !! DAXPY constant times a vector plus a vector.
         !! uses unrolled loops for increments equal to one.

      real(wp) :: da
      integer(ip) :: incx, incy, n
      real(wp) :: dx(*), dy(*)

      integer(ip) :: i, ix, iy, m, mp1

      if (n <= 0_ip) return
      if (da == 0.0_wp) return
      if (incx == 1_ip .and. incy == 1_ip) then
         ! code for both increments equal to 1
         ! clean-up loop
         m = mod(n, 4_ip)
         if (m /= 0_ip) then
            do i = 1_ip, m
               dy(i) = dy(i) + da*dx(i)
            end do
         end if
         if (n < 4_ip) return
         mp1 = m + 1_ip
         do i = mp1, n, 4_ip
            dy(i) = dy(i) + da*dx(i)
            dy(i + 1_ip) = dy(i + 1_ip) + da*dx(i + 1_ip)
            dy(i + 2_ip) = dy(i + 2_ip) + da*dx(i + 2_ip)
            dy(i + 3_ip) = dy(i + 3_ip) + da*dx(i + 3_ip)
         end do
      else
         ! code for unequal increments or equal increments
         ! not equal to 1
         ix = 1_ip
         iy = 1_ip
         if (incx < 0_ip) ix = (-n + 1_ip)*incx + 1_ip
         if (incy < 0_ip) iy = (-n + 1_ip)*incy + 1_ip
         do i = 1_ip, n
            dy(iy) = dy(iy) + da*dx(ix)
            ix = ix + incx
            iy = iy + incy
         end do
      end if

   end subroutine daxpy