PROGRAM MAIN
      IMPLICIT INTEGER (A-Z)
      PARAMETER (NEXMX=10)

      DIMENSION ILEFT(2*NEXMX),IRIGHT(2*NEXMX)

      WRITE(6,*) ' give NOCC '
      READ(5,*) NOCC
    
100   CONTINUE
      WRITE(6,*) ' give the left determinant (excitation '
     $     ,'degree, indices, -1 stops the program)' 
      READ(5,*,IOSTAT=KK) NEXL,(ILEFT(I),I=1,2*NEXL)
      IF (KK.NE.0) GO TO 200
      if (nexl.eq.-1) stop

      WRITE(6,*) ' give the right determinant (excitation '
     $     ,'degree, indices)' 
      READ(5,*) NEXR,(IRIGHT(I),I=1,2*NEXR)

      CALL XMATEL(ILEFT,IRIGHT,NEXL,NEXR,NOCC)
      WRITE(6,*) 
      WRITE(6,*)  ' - - - - - - - - - - - - - - - - - -'
      WRITE(6,*) 
      WRITE(6,*) 
      GO TO 100
 200  CONTINUE

      END

      SUBROUTINE XMATEL(ILEFT,IRIGHT,NEXL,NEXR,NOCC)
C maximally doubly excited
      IMPLICIT INTEGER (A-Z)
      PARAMETER (NEXMX=10)
      CHARACTER*8 PNAME
      CHARACTER*1 CSIGN(-1:1)
      DIMENSION ILEFT(2*NEXMX),IRIGHT(2*NEXMX)
      DIMENSION ISPINL(2*NEXMX),ISPINR(2*NEXMX)
      DIMENSION ITYPEL(2*NEXMX),ITYPER(2*NEXMX)
      DIMENSION ILEFT2(2*NEXMX),IRIGHT2(2*NEXMX)
      DIMENSION IHOLES(2*NEXMX)
      DATA CSIGN /'-',' ','+'/
      
      IPRINT=5
            
      DO I=1,2*MAX(NEXL,NEXR)
       ILEFT2(I)=0
       IRIGHT2(I)=0
      END DO

      WRITE(6,*) ' left string :  ',(ILEFT(I),I=1,2*NEXL)
      WRITE(6,*) ' right string : ',(IRIGHT(I),I=1,2*NEXR)
C
C how many holes in the LHS ?
      ILHA=0
      ILHB=0
      ILPA=0
      ILPB=0
      DO I=1,2*NEXL
       IF (ILEFT(I).LT.0) THEN
        IF (ABS(ILEFT(I)).LE.NOCC) THEN
         ILHB=ILHB+1
         ITYPEL(I)=1
         ISPINL(I)=1
        ELSE 
         ILPB=ILPB+1
         ITYPEL(I)=2
         ISPINL(I)=1
        END IF 
       ELSE IF (ILEFT(I).GT.0) THEN
        IF (ILEFT(I).LE.NOCC) THEN
         ILHA=ILHA+1 
         ITYPEL(I)=1
         ISPINL(I)=2
        ELSE
         ILPA=ILPA+1 
         ITYPEL(I)=2
         ISPINL(I)=2
        END IF 
       END IF
      END DO 
      
      IRHA=0
      IRHB=0
      IRPA=0
      IRPB=0
      DO I=1,2*NEXR
       IF (IRIGHT(I).LT.0) THEN
        IF (ABS(IRIGHT(I)).LE.NOCC) THEN
         IRHB=IRHB+1
         ITYPER(I)=1
         ISPINR(I)=1
        ELSE
         IRPB=IRPB+1
         ITYPER(I)=2
         ISPINR(I)=1
        END IF 
       ELSE IF (IRIGHT(I).GT.0) THEN
        IF (IRIGHT(I).LE.NOCC) THEN
         IRHA=IRHA+1 
         ITYPER(I)=1
         ISPINR(I)=2
        ELSE
         IRPA=IRPA+1 
         ITYPER(I)=2
         ISPINR(I)=2
        END IF 
       END IF 
      END DO 
      NLEFT=ILHA+ILHB+ILPA+ILPB
      NRIGHT=IRHA+IRHB+IRPA+IRPB
      
      WRITE(6,*) ' left  determinant : holes alpha, beta ',ILHA,ILHB
      WRITE(6,*) ' left  determinant : part. alpha, beta ',ILPA,ILPB
      WRITE(6,*) ' right determinant : holes alpha, beta ',IRHA,IRHB
      WRITE(6,*) ' right determinant : part. alpha, beta ',IRPA,IRPB

      IF (ILHA.NE.ILPA) STOP ' error in holes/particles alpha LHS'
      IF (ILHB.NE.ILPB) STOP ' error in holes/particles beta  LHS'
      IF (IRHA.NE.IRPA) STOP ' error in holes/particles alpha RHS'
      IF (IRHB.NE.IRPB) STOP ' error in holes/particles beta  RHS'

C transform the strings to occupations
      INDX=0
      DO I=1,2*NEXL
       IF (ILEFT(I).NE.0.AND.ABS(ILEFT(I)).LE.NOCC) THEN
        INDX=INDX+1
        ILEFT2(INDX)=ILEFT(I)
       END IF
      END DO
      DO I=1,2*NEXR
       IF (IRIGHT(I).NE.0.AND.ABS(IRIGHT(I)).LE.NOCC) THEN
        INDX=INDX+1
        ILEFT2(INDX)=IRIGHT(I)
       END IF
      END DO

      WRITE(6,*) ' ILEFT2 ',(ILEFT2(JJ),JJ=1,NEXL+NEXR)
C
C eliminate same holes
C
      DO I=1,NEXL+NEXR
       DO J=I+1,NEXL+NEXR
        IF (ILEFT2(I).EQ.ILEFT2(J).AND.ILEFT2(I).NE.0) THEN
         ILEFT2(J)=0 
        END IF
       END DO
      END DO

      ISIGN=1
C
C we may have NEXL+NEXR different holes
C 1 2 3 4
C
C they may be reduced to 3 or 2, and during this action the RHS holes
C may have changed their order
C
      NORBI=NEXL+NEXR
      DO I=1,NEXL+NEXR
       IF (ILEFT2(I).EQ.0) THEN
        LFOUND=0
        DO J=I+1,NEXL+NEXR
         IF (LFOUND.EQ.0) THEN
          IF (ILEFT2(J).NE.0) THEN
           ILEFT2(I)=ILEFT2(J)
           ILEFT2(J)=0
           LFOUND=1
          END IF
         END IF
        END DO
       END IF
      END DO

      NORBI=0
      DO I=1,NEXL+NEXR
       IF (ILEFT2(I).NE.0) NORBI=NORBI+1
      END DO

C check order of the RHS holes
      IF (IRHA.GE.2.OR.IRHB.GE.2) THEN
       DO I=1,NORBI
        IF (ILEFT2(I).EQ.IRIGHT(1)) IRP1=I
        IF (ILEFT2(I).EQ.IRIGHT(2)) IRP2=I
       END DO
       IF (IRP1.GT.IRP2) ISIGN=-ISIGN
      END IF

      WRITE(6,9902) NORBI,(ILEFT2(J),J=1,NORBI)
 9902 FORMAT('   NORBI = ',I4,' STRING ',20I4)
C save the string of the holes
      DO I=1,NORBI
       IHOLES(I)=ILEFT2(I)
      END DO
      NHOLES=NORBI

      DO I=1,NORBI
       IRIGHT2(I)=ILEFT2(I)
      END DO

C now we have to place the particles

      DO I=1,NORBI
C if this hole appears in the left string, we look for an adequate
c particle 
       LFOUND=0
       DO J=1,2*NEXL
        IF (ILEFT2(I).EQ.ILEFT(J)) THEN
         LFOUND=1
         LSPIN=ISPINL(J)
        END IF
       END DO
C and replace the hole
       IF (LFOUND.EQ.1) THEN
        LSUB=0
        DO J=1,2*NEXL
         IF (LSUB.EQ.0) THEN
C look for a particle with the same spin
          IF ((LSPIN.EQ.ISPINL(J)).AND.(ITYPEL(J).EQ.2)) THEN
           ILEFT2(I)=ILEFT(J)
           ITYPEL(J)=0
           LSUB=1
          END IF
         END IF
        END DO
       END IF
      END DO
      
      DO I=1,NORBI
C if this hole appears in the right string, we look for an adequate
c particle 
       LFOUND=0
       DO J=1,2*NEXR
        IF (IRIGHT2(I).EQ.IRIGHT(J)) THEN
         LFOUND=1
         LSPIN=ISPINR(J)
        END IF
       END DO
C and replace the hole
       IF (LFOUND.EQ.1) THEN
        LSUB=0
        DO J=1,2*NEXR
         IF (LSUB.EQ.0) THEN
C look for a particle with the same spin
          IF ((LSPIN.EQ.ISPINR(J)).AND.(ITYPER(J).EQ.2)) THEN
           IRIGHT2(I)=IRIGHT(J)
           ITYPER(J)=0
           LSUB=1
          END IF
         END IF
        END DO
       END IF
      END DO

      WRITE(6,9903) 'left ',(ILEFT2(J),J=1,NORBI)
      WRITE(6,9903) 'right',(IRIGHT2(J),J=1,NORBI)
      WRITE(6,9804) 'holes ',(IHOLES(J),J=1,NORBI)
 9903 FORMAT('   ',A5,' determinant as occupation : ',10I4) 
 9804 FORMAT('   ',A5,' in the determinants       : ',10I4)

C
C now we have to determine how many of the indices are equal
C
      NEQUAL=0
      DO I=1,NORBI
       LFOUND=0
       DO J=1,NORBI
        IF (ILEFT2(I).EQ.IRIGHT2(J).AND.LFOUND.EQ.0) THEN
         NEQUAL=NEQUAL+1
         LFOUND=1
         IF (I.NE.J) THEN
C exchange I and J in the right string
          ISIGN=-ISIGN
          IDUM=IRIGHT2(J)
          IRIGHT2(J)=IRIGHT2(I)
          IRIGHT2(I)=IDUM          
         END IF
        END IF
       END DO
      END DO
C
C get the different particles to the left-hand side
      DO I=1,NORBI
       LFOUND=0
       DO J=I,NORBI
C are these already different ?
        IF (ILEFT2(J).NE.IRIGHT2(J).AND.LFOUND.EQ.0) THEN
         IDUM=ILEFT2(I)
         ILEFT2(I)=ILEFT2(J)
         ILEFT2(J)=IDUM
         IDUM=IRIGHT2(I)
         IRIGHT2(I)=IRIGHT2(J)
         IRIGHT2(J)=IDUM
         LFOUND=1 
        END IF
       END DO
      END DO

      WRITE(6,*) '-left ',(ILEFT2(J),J=1,NORBI)
      WRITE(6,*) '-right',(IRIGHT2(J),J=1,NORBI)
 
      WRITE(6,*) 'left ',(ILEFT2(J),J=1,NORBI)
      WRITE(6,*) 'right',(IRIGHT2(J),J=1,NORBI)
      
      NDIFF=NORBI-NEQUAL

      WRITE(6,*) ' implied orbitals : ',NORBI
      WRITE(6,*) ' equal indices    : ',NEQUAL
      WRITE(6,*) ' different indices: ',NDIFF
      WRITE(6,*) 
      WRITE(6,*) ' sign of the permutation: ',ISIGN
      
      WRITE(6,9903) 'left ',(ILEFT2(J),J=1,NDIFF)
      WRITE(6,9903) 'right',(IRIGHT2(J),J=1,NDIFF)

      IF (NDIFF.LE.2) THEN
       
       DO I=1,NORBI
        IF (ILEFT2(I).LT.0) THEN
         ISPINL(I)=1
        ELSE
         ISPINL(I)=2
        END IF
        IF (IRIGHT2(I).LT.0) THEN
         ISPINR(I)=1
        ELSE
         ISPINR(I)=2
        END IF
       END DO
C
       IF (NDIFF.EQ.0) THEN
        WRITE(6,*) ' NDIFF = 0 : diagonal element '
        WRITE(6,9924) CSIGN(ISIGN)
 9924   FORMAT('  calculating matrix element : ',A1,'E_HF ')
        DO I=1,NORBI
         I1=ABS(IHOLES(I))
         J1=ABS(ILEFT2(I))
         WRITE(6,9921) CSIGN(-ISIGN),I1,I1,CSIGN(ISIGN),J1,J1
 9921    FORMAT(' diagonal element : ',A1,'F (',I3,',',I3,' ) ',A1
     $        ,' F (',I3,',',I3,' )')
        END DO

C
C hole-hole, hole-particle and particle-particle interactions
C hole-hole and part-part enter with +
C hole-part enter with -
C same-spin combinations are J-K
C different-spin combinations are J only
C
        DO I=1,NORBI
         DO J=1,NORBI

C hole-part and part-hole

          I1=ABS(IHOLES(I))
          J1=ABS(ILEFT2(J))
          
          JSIGN=1
          IF (I1.LE.NOCC.AND.J1.GT.NOCC) JSIGN=-1 
          IF (I1.GT.NOCC.AND.J1.LE.NOCC) JSIGN=-1
          IF (ISPINL(I).EQ.ISPINL(J)) THEN
           WRITE(6,9922) CSIGN(ISIGN*JSIGN), I1,I1,J1,J1,CSIGN(
     $          -ISIGN*JSIGN),I1,J1,I1,J1
          ELSE
           WRITE(6,9923) CSIGN(ISIGN*JSIGN),I1,I1,J1,J1
          END IF 

         END DO
        END DO

        DO I=1,NORBI-1
         DO J=I+1,NORBI
C hole-hole and part-part
          I1=ABS(IHOLES(I))
          J1=ABS(IHOLES(J))
          
          JSIGN=1
          IF (I1.LE.NOCC.AND.J1.GT.NOCC) JSIGN=-1 
          IF (I1.GT.NOCC.AND.J1.LE.NOCC) JSIGN=-1
          IF (ISPINL(I).EQ.ISPINL(J)) THEN
           WRITE(6,9922) CSIGN(ISIGN*JSIGN), I1,I1,J1,J1,CSIGN(
     $          -ISIGN*JSIGN),I1,J1,I1,J1
          ELSE
           WRITE(6,9923) CSIGN(ISIGN*JSIGN),I1,I1,J1,J1
          END IF 


          I1=ABS(ILEFT2(I))
          J1=ABS(ILEFT2(J))
          
          JSIGN=1
          IF (I1.LE.NOCC.AND.J1.GT.NOCC) JSIGN=-1 
          IF (I1.GT.NOCC.AND.J1.LE.NOCC) JSIGN=-1
          IF (ISPINL(I).EQ.ISPINL(J)) THEN
           WRITE(6,9922) CSIGN(ISIGN*JSIGN), I1,I1,J1,J1,CSIGN(
     $          -ISIGN*JSIGN),I1,J1,I1,J1
          ELSE
           WRITE(6,9923) CSIGN(ISIGN*JSIGN),I1,I1,J1,J1
          END IF 

         END DO
        END DO
 9922   FORMAT(' diagonal element : ',A1,' (',2I4,'|',2I4,' ) ',A1,' ('
     $       ,2I4,'|',2I4,' ) ')
 9923   FORMAT(' diagonal element : ',A1,' (',2I4,'|',2I4,' )')
       ELSE IF (NDIFF.EQ.1) THEN
        I=MIN(ABS(ILEFT2(1)),ABS(IRIGHT2(1)))
        J=MAX(ABS(ILEFT2(1)),ABS(IRIGHT2(1)))
        WRITE(6,9912) ILEFT2(1),IRIGHT2(1),CSIGN(ISIGN),I,J 
C subtract all hole interactions
        DO IH=1,NHOLES
         IHOLE=IHOLES(IH)
C        WRITE(6,*) ' IHOLE = ',IHOLE,ILEFT2(1)
         IF (SIGN(1,IHOLE).EQ.SIGN(1,ILEFT2(1))) THEN
          IHOLE=ABS(IHOLE)
          II1=ABS(ILEFT2(1))
          JJ1=ABS(IRIGHT2(1))
          KK1=IHOLE
          LL1=IHOLE
          II2=ABS(ILEFT2(1))
          JJ2=IHOLE
          KK2=ABS(IRIGHT2(1))
          LL2=IHOLE
          CALL OCAN(II1,JJ1,KK1,LL1)
          CALL OCAN(II2,JJ2,KK2,LL2)
          WRITE(6,9913) ILEFT2(1),IRIGHT2(1),CSIGN(
     $         -ISIGN),II1,JJ1,KK1,LL1,CSIGN(ISIGN),II2,JJ2,KK2,LL2
         ELSE
          IHOLE=ABS(IHOLE)
          II1=ABS(ILEFT2(1))
          JJ1=ABS(IRIGHT2(1))
          KK1=IHOLE
          LL1=IHOLE
          CALL OCAN(II1,JJ1,KK1,LL1)
          WRITE(6,9914) ILEFT2(1),IRIGHT2(1),CSIGN(
     $         -ISIGN),II1,JJ1,KK1,LL1
         END IF 
        END DO  
C add all particle interactions
        DO IEQUAL=1,NEQUAL
         IPART=ILEFT2(1+IEQUAL)
C        WRITE(6,*) ' IPART = ',IPART
         IF (SIGN(1,IPART).EQ.SIGN(1,ILEFT2(1))) THEN
          IPART=ABS(IPART)
          II1=ABS(ILEFT2(1))
          JJ1=ABS(IRIGHT2(1))
          KK1=IPART
          LL1=IPART
          II2=ABS(ILEFT2(1))
          JJ2=IPART
          KK2=ABS(IRIGHT2(1))
          LL2=IPART
          CALL OCAN(II1,JJ1,KK1,LL1)
          CALL OCAN(II2,JJ2,KK2,LL2)
          WRITE(6,9913) ILEFT2(1),IRIGHT2(1)
     $         ,CSIGN(ISIGN),II1,JJ1,KK1,LL1,CSIGN(-ISIGN),II2,JJ2,KK2 
     $         ,LL2
         ELSE
          IPART=ABS(IPART)
          II1=ABS(ILEFT2(1))
          JJ1=ABS(IRIGHT2(1))
          KK1=IPART
          LL1=IPART
          CALL OCAN(II1,JJ1,KK1,LL1)
          WRITE(6,9914) ILEFT2(1),IRIGHT2(1)
     $         ,CSIGN(ISIGN),II1,JJ1,KK1,LL1
         END IF
        END DO
C we have the interaction with the holes 
       ELSE IF (NDIFF.EQ.2) THEN
        I1=ILEFT2(1)
        I2=ILEFT2(2)
        J1=IRIGHT2(1)
        J2=IRIGHT2(2)
        DO I=1,2
         ILEFT2(I)=ABS(ILEFT2(I))
         IRIGHT2(I)=ABS(IRIGHT2(I))
        END DO
C same spins
        IF (ISPINL(1).EQ.ISPINL(2)) THEN
         
         II1=ILEFT2(1)
         JJ1=IRIGHT2(1)
         KK1=ILEFT2(2)
         LL1=IRIGHT2(2)
         
         CALL OCAN(II1,JJ1,KK1,LL1)
         
         II2=ILEFT2(1)
         JJ2=IRIGHT2(2)
         KK2=ILEFT2(2)
         LL2=IRIGHT2(1)
         CALL OCAN(II2,JJ2,KK2,LL2)
         WRITE(6,9915) I1,I2,J1,J2,CSIGN(ISIGN),II1,JJ1,KK1,LL1,CSIGN(
     $        -ISIGN),II2,JJ2,KK2,LL2  
        ELSE
         II1=ILEFT2(1)
         JJ1=IRIGHT2(1)
         KK1=ILEFT2(2)
         LL1=IRIGHT2(2)
         CALL OCAN(II1,JJ1,KK1,LL1)
         WRITE(6,9916)  I1,I2,J1,J2,CSIGN(ISIGN),II1,JJ1,KK1,LL1
        END IF
       END IF
      ELSE
       WRITE(6,9917)
      END IF
      
 9904 FORMAT(' calculating diagonal element <',4I4,'||',4I4
     $     ,'> ')
 9912 FORMAT(' calculating matrix element <',I4,'||',I4
     $     ,'> : interaction ',A1,'F(',I3,',',I3,')',F20.12)
 9913 FORMAT(' calculating matrix element <',I4,'||',I4
     $     ,'> : interaction ',A1,'(',2I4,'|',2I4,') ',A1,' (',2I4,'|'
     $     ,2I4,') ',F20.12) 
 9914 FORMAT(' calculating matrix element <',I4,'||',I4
     $     ,'> : interaction ',A1,'(',2I4,'|',2I4,') ',F20.12) 
 9915 FORMAT(' calculating matrix element <',2I4,'||',2I4
     $     ,'> : interaction ',A1,'(',2I4,'|',2I4,') ',A1,' (',2I4,'|'
     $     ,2I4,') ',F20.12) 
 9916 FORMAT(' calculating matrix element <',2I4,'||',2I4
     $     ,'> : interaction ',A1,'(',2I4,'|',2I4,') ',F20.12) 
 9917 FORMAT(' calculating matrix element : no interaction ') 

      RETURN
      END

      SUBROUTINE OCAN(I,J,K,L)
      IMPLICIT INTEGER (A-Z)

      I1=I
      J1=J
      K1=K
      L1=L

      IF (I1.GT.J1) THEN
       IDUM=I1
       I1=J1
       J1=IDUM
      END IF
      IF (K1.GT.L1) THEN
       IDUM=K1
       K1=L1
       L1=IDUM
      END IF

      IF (I1.GT.K1) THEN
       IDUM=I1
       I1=K1
       K1=IDUM
       IDUM=J1
       J1=L1
       L1=IDUM
      ELSE IF (I1.EQ.K1) THEN
       IF (J1.GT.L1) THEN
        IDUM=J1
        J1=L1
        L1=IDUM
       END IF
      END IF

      I=I1
      J=J1
      K=K1
      L=L1

      RETURN
      END