anchor.diff:    !      write(6,*)'               reduk: end'
     !      write(6,*)'------------------------------------------------------'
     return
  +contains
  +subroutine r3frac(eps,v,iv)
  +! Copyright (C) 2002-2005 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
  +! This file is distributed under the terms of the GNU Lesser General Public
  +! License.
  +implicit none
  +! arguments
  +real(8), intent(in) :: eps
  +real(8), intent(inout) :: v(3)
  +integer, intent(out) :: iv(3)
  +! local variables
  +integer i
  +do i=1,3
  +  iv(i)=int(v(i))
  +  v(i)=v(i)-dble(iv(i))
  +  if (v(i).lt.0.d0) then
  +    v(i)=v(i)+1.d0
  +    iv(i)=iv(i)-1
  +  end if
  +  if (1.d0-v(i).lt.eps) then
  +    v(i)=0.d0
  +    iv(i)=iv(i)+1
  +  end if
  +  if (v(i).lt.eps) then
  +    v(i)=0.d0
  +  end if
  +end do
  +return
  +end subroutine
   end subroutine reduk
   !EOC
  
positive.code:   implicit none
    !<sag>
    logical, allocatable :: done(:)
    !</sag>
  
  
  
negative.code:                  !              write(24,'(i6,3i4)')divsh,div
                   !              write(24,'(3i4,"  ",3i4)')onek,kk
                   do i=1,nsymt
                      ! from internal coordinates to lattice coordinates
                      kpr(:)=dble(kk(:))/dble(divsh*div(:))
                      ! rotate k-point
                      kprt(:)=iio(:,1,i)*kpr(1)+iio(:,2,i)*kpr(2)+&
                           iio(:,3,i)*kpr(3)
                      ! map to reciprocal unit cell
                      call r3frac(1.d-6,kprt,iv)
                      ! from lattice coordinates to internal coordinates
                      kprt(:)=kprt(:)*divsh*div(:)
                      nkp(:)=nint(kprt(:))
                      ! determine fractional part of rotated k-point
                      call r3frac(1.d-6,kprt,iv)
                      ! if fractional part is present discard k-point
                      if (any(kprt>1.d-6)) nkp(:)=-1
  
                      !               write(24,'(i6,3i4,"  ",3i4)')i,ktp,nkp
                      nkpid=idkp(nkp)
                      !                write(24,'(4x,2i4)')nirkp,nkpid
                      !---                call jset(kpav,nkpid,1)
                      if (.not.all(nkp.eq.-1)) then
                         done(nkpid)=.true. !+++
                         redkp(nkpid)=kpid
                         starr(i)=nkpid
                      end if
  write(1234,'(a,2i9,3i6,i9)') '  kpid,sym,nkp,nkpid',kpid,i,nkp,nkpid
                   enddo
                   members=0
                   do i=1,nsymt
  