stdlib_saxpy Subroutine

public pure subroutine stdlib_saxpy(n, sa, sx, incx, sy, incy)

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

Arguments

Type IntentOptional Attributes Name
integer(kind=ilp), intent(in) :: n
real(kind=sp), intent(in) :: sa
real(kind=sp), intent(in) :: sx(*)
integer(kind=ilp), intent(in) :: incx
real(kind=sp), intent(inout) :: sy(*)
integer(kind=ilp), intent(in) :: incy

Source Code

     pure subroutine stdlib_saxpy(n,sa,sx,incx,sy,incy)
     !! SAXPY constant times a vector plus a vector.
     !! uses unrolled loops for increments equal to one.
        ! -- reference blas level1 routine --
        ! -- reference blas is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           real(sp), intent(in) :: sa
           integer(ilp), intent(in) :: incx, incy, n
           ! Array Arguments 
           real(sp), intent(in) :: sx(*)
           real(sp), intent(inout) :: sy(*)
        ! =====================================================================
           ! Local Scalars 
           integer(ilp) :: i, ix, iy, m, mp1
           ! Intrinsic Functions 
           intrinsic :: mod
           if (n<=0) return
           if (sa==0.0_sp) return
           if (incx==1 .and. incy==1) then
              ! code for both increments equal to 1
              ! clean-up loop
              m = mod(n,4)
              if (m/=0) then
                 do i = 1,m
                    sy(i) = sy(i) + sa*sx(i)
                 end do
              end if
              if (n<4) return
              mp1 = m + 1
              do i = mp1,n,4
                 sy(i) = sy(i) + sa*sx(i)
                 sy(i+1) = sy(i+1) + sa*sx(i+1)
                 sy(i+2) = sy(i+2) + sa*sx(i+2)
                 sy(i+3) = sy(i+3) + sa*sx(i+3)
              end do
           else
              ! code for unequal increments or equal increments
                ! not equal to 1
              ix = 1
              iy = 1
              if (incx<0) ix = (-n+1)*incx + 1
              if (incy<0) iy = (-n+1)*incy + 1
              do i = 1,n
               sy(iy) = sy(iy) + sa*sx(ix)
               ix = ix + incx
               iy = iy + incy
              end do
           end if
           return
     end subroutine stdlib_saxpy