ddotf.f

double precision function ddotf(n,dx,incx,dy,incy)
c
c forms the dot product of two vectors.
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(*),dtemp
integer i,incx,incy,ix,iy,m,mp1,n
c
 ddotf =0.0d0
 dtemp =0.0d0
if(n.le.0)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
 dtemp = dtemp + dx(ix)*dy(iy)
 ix = ix + incx
 iy = iy + incy
10continue
 ddotf = dtemp
return
c
c code for both increments equal to 1
c
c
c clean-up loop
c
20 m =mod(n,5)
if( m .eq.0 ) go to40
do30 i =1,m
 dtemp = dtemp + dx(i)*dy(i)
30continue
if( n .lt.5 ) go to60
40 mp1 = m +1
do50 i = mp1,n,5
 dtemp = dtemp + dx(i)*dy(i) + dx(i +1)*dy(i +1) +
* dx(i +2)*dy(i +2) + dx(i +3)*dy(i +3) + dx(i +4)*dy(i +4)
50continue
60 ddotf = dtemp
return
 end