lanczos_iteration_d Subroutine

subroutine lanczos_iteration_d(MatVec,iter,vin,vout,alfa,beta)
  interface
     subroutine MatVec(Nloc,vin,vout)
       integer                 :: Nloc
       real(8),dimension(Nloc) :: vin
       real(8),dimension(Nloc) :: vout
     end subroutine MatVec
  end interface
  integer                                    :: iter
  real(8),dimension(:),intent(inout)         :: vin
  real(8),dimension(size(vin)),intent(inout) :: vout
  real(8),dimension(size(vin))               :: tmp
  real(8),intent(inout)                      :: alfa,beta
  integer                                    :: nloc
  real(8)                                    :: norm
  !
  nloc=size(vin)
  !
  if(iter==1)then
     norm=sqrt(dot_product(vin,vin))
     if(norm==0.d0)stop "LANCZOS_ITERATION_D: norm =0!!"
     vin    = vin/norm
  else
     tmp = vin
     vin = vout/beta
     vout= -beta*tmp
  endif
  call MatVec(nloc,vin,tmp)
  vout = vout + tmp
  alfa = dot_product(vin,vout)
  vout = vout - alfa*vin
  beta = sqrt(dot_product(vout,vout))
end subroutine lanczos_iteration_d