********************************************************************************
*                                                                              * 
*  2DHF version 1-2003                                                         *
*  Copyright (C) 1996  Jacek Kobus, Leif Laaksonen, Dage Sundholm              *
*                                                                              *
*  This software may be used and distributed according to the terms            *
*  of the GNU General Public License, see README and COPYING.                  *
*                                                                              *
********************************************************************************
c ### preprelexch ###
c
c     Determines the number of exchange potentials associated with a given
c     orbital (iorb) and stores in array iexchorb orbitals needed to
c     evaluate these potentials 

      subroutine preprelexch(iorb) 
      implicit integer*4 (i-n)
      implicit real*8 (a-h,o-z)
c_mpi
      include 'mpif.h'
      include 'commons_mpi.inc'
c_mpi
      include 'commons8.inc'
      parameter(iprint=0)

      call MPI_COMM_RANK( MPI_COMM_WORLD, myid, IERROR )

      if     (exlexp.eq.0.or.exlexp.eq.1) then
	  ifirst=1
      elseif (exlexp.eq.2) then
c        recalculate only those exchange potentials which depend on orbitals 
c        modified so far (note the reverse order of relaxation in scf)
	  ifirst=iorb
      endif

c     zeroise working arrays
      nexchorb(iorb)=0
      do iorb1=1,norb
         iexchorb(iorb,iorb1)=0
         iexchadd(iorb,iorb1)=0
         iexchlen(iorb,iorb1)=0
      enddo

      do iorb1=ifirst,norb
         iswtch=0
         
         if (iorb.eq.iorb1.and.mgx(6,iorb).eq.0 ) goto 10
         if ((iorb.eq.iorb1).and.(ilc(iorb*(iorb+1)/2).lt.1)) goto 10

c        orbitals in increasing order
         if (iorb.lt.iorb1) then
	    in1=iorb
	    in2=iorb1	  
         else	    
	    in1=iorb1
	    in2=iorb
         endif	    

         ipc=in1+in2*(in2-1)/2
         iwexch(ipc)=-1
         ibexp=i3b(ipc)
         ngexp=i3si(ipc)

         idel=iabs(mgx(6,in1)-mgx(6,in2))
         if (in1.eq.in2) idel=2*mgx(6,in1)
         idel1=mgx(6,in1)+mgx(6,in2)

         ipex=idel-2*(idel/2)
         if (ipex.eq.0) then	  
	    isym= 1
         else
	    isym=-1
         endif
 1000    continue

c        nexchorb(iorb) is used to index array iexchorb(iorb,i) which
c        contains orbitals forming exchange potentials with the orbital
c        iorb. Eventually nexchorb will be equal to the number of
c        exchange potentials for the orbital iorb.

         nexchorb(iorb)=nexchorb(iorb)+1
         if (iswtch.eq.0) then 
            iexchorb(iorb,nexchorb(iorb))=iorb1
            iexchadd(iorb,iorb1)=ibexp
            iexchlen(iorb,iorb1)=ngexp
         elseif (iswtch.eq.1) then
            iexchorb(iorb,nexchorb(iorb))=(-1)*iorb1
         endif

         if (iswtch.eq.1) goto 10
         if (ilc(ipc).lt.2) goto 10
         iwexch(ipc)=-2
         idel=idel1
         ibexp=ibexp+ngexp
         ipc=ipc+norb*(norb+1)/2
         iswtch=1
         goto 1000
 10      continue
      enddo
      if (myid.eq.0.and.iprint.eq.1) then 
         write(*,*) 'preprelexch'
         write(*,'("iorb =",i10,/,20i5)') 
     &        iorb,(iexchorb(iorb,i),i=1,nexchorb(iorb)) 
      endif
      return
      end