daxpyf.f

subroutinedaxpyf(n,da,dx,incx,dy,incy)
c
c constant times a vector plus a vector.
c uses unrolled loops for increments equal to one.
c jack dongarra, linpack, 3/11/78.
c modified 12/3/93, array(1) declarations changed to array(*)
c
double precision dx(*),dy(*),da
integer i,incx,incy,ix,iy,m,mp1,n
c
if(n.le.0)return
if (da .eq.0.0d0) return
if(incx.eq.1.and.incy.eq.1)go to20
c
c code for unequal increments or equal increments
c not equal to 1
c
 ix =1
 iy =1
if(incx.lt.0)ix = (-n+1)*incx +1
if(incy.lt.0)iy = (-n+1)*incy +1
do10 i =1,n
 dy(iy) = dy(iy) + da*dx(ix)
 ix = ix + incx
 iy = iy + incy
10continue
return
c
c code for both increments equal to 1
c
c
c clean-up loop
c
20 m =mod(n,4)
if( m .eq.0 ) go to40
do30 i =1,m
 dy(i) = dy(i) + da*dx(i)
30continue
if( n .lt.4 ) return
40 mp1 = m +1
do50 i = mp1,n,4
 dy(i) = dy(i) + da*dx(i)
 dy(i +1) = dy(i +1) + da*dx(i +1)
 dy(i +2) = dy(i +2) + da*dx(i +2)
 dy(i +3) = dy(i +3) + da*dx(i +3)
50continue
return
end