outlets/iso.f

      PROGRAM OUTLET
      IMPLICIT NONE
      INTEGER NMMAX,NMAX
      PARAMETER (NMMAX=50000,NMAX=100)
      DOUBLE PRECISION X(NMAX)
      INTEGER I, J, N, M, NM, IT
      INTEGER IPERM(NMAX)
      INTEGER L1(NMAX), L2(NMAX)
      INTEGER PATTERNS ( NMMAX, NMAX )
      INTEGER CONF(NMAX,6), FREQ(NMMAX), IPAT
      DOUBLE PRECISION DFREQ(NMMAX)
      INTEGER PERM ( NMAX, 6), P1(NMAX), P2(NMAX), P(NMAX)
      INTEGER PDISTANCE
      INTEGER NITER
      CHARACTER (LEN=10) STR
      CHARACTER (LEN=20) ARGV
      INTEGER NARG

      CALL GETARG(1, ARGV)
      READ (ARGV,*) N
      CALL GETARG(2, ARGV)
      READ (ARGV,*) M
      CALL GETARG(3, ARGV)
      READ (ARGV,*) NITER

    
      NM = N * M


      DFREQ = 0.D0
      FREQ  = 0
      CALL GRAPHS (N,M,PATTERNS,IPAT,FREQ,NITER)
      PRINT*, N, M, IPAT
      DO 10, I = 1, IPAT
      DFREQ ( I ) = DBLE( FREQ(I) ) / DBLE(NITER) * 100.D0 
      WRITE(STR,'(I5)') NM
10    WRITE(6, "(I5,F8.2," // ADJUSTL(STR) // "I5)" )  
     &                 I, DFREQ( I ), (PATTERNS(I,J), J=1,NM)

 
      END



      SUBROUTINE GRAPHS ( N, M, PATTERNS, IPAT, FREQ, ITMAX )
      IMPLICIT NONE
      INTEGER NMMAX,NMAX
      PARAMETER (NMMAX=50000,NMAX=100)
      INTEGER N, M, NM, IPAT, ITMAX
      INTEGER PATTERNS ( NMMAX, NMAX )
      INTEGER FREQ (NMAX)
      INTEGER I, J, K, IT
      INTEGER PERM ( NMAX, 6), P1(NMAX), P2(NMAX), P(NMAX)
      LOGICAL ISOMORPHICQ
      LOGICAL Q

      NM = N * M

      CALL PERMUTATIONS ( N, 100, PERM )
C      DO 1, I = 1, NM
C1     PRINT * , ( PERM ( I , J ), J = 1 , N )

      CALL GENER ( N, M, P1 )
C     PRINT *,  (P1(I), I=1,NM)

      DO 2, J = 1, NM
2     PATTERNS(1,J) = P1(J)

      FREQ    = 0

      IPAT    = 1
      FREQ(1) = 1


C MAIN LOOP
      IT = 1
10    CONTINUE

      Q = .FALSE.

      CALL GENER ( N, M, P2 )

C Check isomorphy against all existing patterns
      K = 0
20    CONTINUE
      K = K + 1

      DO 12, J = 1, NM
 12   P1(J) = PATTERNS(K,J)


      IF (ISOMORPHICQ(N,M,P1,P2,PERM)) THEN
         Q = .TRUE.
         FREQ ( K ) = FREQ ( K ) + 1
         IT = IT + 1
C        PRINT*, (P2(J), J=1,NM)
         IF (IT.LT.ITMAX) GOTO 10
      ELSE
         IF (K.LT.IPAT) GOTO 20
      ENDIF

C If it doesn't match any existing pattern add it
      IF (.NOT. Q) THEN
        IPAT = IPAT + 1
        FREQ(IPAT) = 1
        DO 14, J=1, NM
 14     PATTERNS(IPAT,J) = P2(J)
C       WRITE(*,'(1H+,A,I5,A)') 'PROGRESS: ', IPAT, '%' 
C       PRINT*, (P2(J), J=1,NM)
      ENDIF
      IT = IT + 1
C     IF (ISOMORPHICQ(N,M,P1,P2,PERM)) THEN
C       FREQ(K) = FREQ(K)+1
C       IT = IT + 1
C       IF (IT.LT.ITMAX) GOTO 10
C       GOTO 60
C     ENDIF
C50   CONTINUE

      IF (IT.LT.ITMAX) GOTO 10

C55   CONTINUE

C     IF ( Q ) THEN 
C       IPAT = IPAT + 1
C       FREQ(IPAT) = 1
C       DO 14, J=1, NM
C14     PATTERNS(IPAT,J) = P2(J)
C       PRINT*, (P2(J), J=1,NM)
C     ELSE


C     IT = IT + 1

C     IF (IT.LT.ITMAX) GOTO 10 

C60   CONTINUE

      RETURN

      END



      FUNCTION ISOMORPHICQ (N, M, P1, P2, PERM)
      IMPLICIT NONE
      LOGICAL ISOMORPHICQ, Q
      INTEGER N, M
      INTEGER P1( N*M ), P2( N*M ), P3( N*M ), P4( N*M ), P11( N*M )
      INTEGER PERM ( 100, N * M )
      INTEGER I, J, K,  NM, FACT, DIFF1, DIFF11

      Q = .FALSE.

      NM = N * M

      FACT = 1
      DO 1, I = 1, N
 1    FACT = I * FACT

      P11 (1) = P1 (1)
      DO 10, I = 2, NM
 10   P11(I) = P1(NM+2-I) 


      DO 2, I = 1, FACT

      DO 3, J = 1, NM    
 3    P3 ( J ) = PERM ( I , P2 ( J ) )
C     PRINT *, "P3 -> " , (P3(K), K=1,NM) 

      DO 2, J = 0, NM - 1
      CALL PADLIST ( P3, P4, NM, J )
      DIFF1  = 0
      DIFF11 = 0
      DO 4, K = 1, NM
      DIFF1  = DIFF1  + ABS ( P4 ( K ) - P1  ( K ) )
 4    DIFF11 = DIFF11 + ABS ( P4 ( K ) - P11 ( K ) )
C     PRINT *,  (P4(K), K=1,NM), DIFF
      
      IF ((DIFF1.EQ.0).OR.(DIFF11.EQ.0)) THEN
        Q = .TRUE.
        GOTO 100
      ENDIF 

 2    CONTINUE 

 100  CONTINUE
      ISOMORPHICQ = Q

      RETURN
      END


      SUBROUTINE GENER ( N, M, P )
      IMPLICIT NONE
      INTEGER NMAX
      PARAMETER (NMAX=100)
      INTEGER N, M, NM
      INTEGER P ( * )
      DOUBLE PRECISION X ( NMAX )
      INTEGER IPERM (NMAX) , IER
      INTEGER I

      NM = N * M 

      DO 1, I=1, NM
 1    X(I) = RAND() 

      CALL DPSORT(X, NM, IPERM, 1, IER) 

      DO 2, I = 1, NM
 2    P ( I ) = ( IPERM ( I ) - 1 ) / M + 1


      RETURN
      END



      SUBROUTINE PERMUTATIONS ( M, MMAX, PERM )
      IMPLICIT NONE
      INTEGER M, MMAX, MFAC
      INTEGER PERM ( MMAX, M )
      INTEGER I, J, A ( M )
      LOGICAL NEXTP
      EXTERNAL NEXTP

      MFAC = 1
      DO 1, I = 1, M
      A ( I ) = I
      PERM ( 1, I ) = A ( I )
      MFAC = MFAC * I
  1   CONTINUE

      DO 2, I = 2, MFAC
      IF(NEXTP(M, A)) THEN
        DO 3, J = 1, M
  3     PERM ( I, J ) = A ( J )
      ENDIF
  2   CONTINUE

      RETURN
      END


C     SUBROUTINE PATTERNS ( CONF, NMAX, N, M, IT, DIST )
C     IMPLICIT NONE
C     INTEGER NMAX, N, M,  IT
C     INTEGER CONF ( NMAX, N ), DIST ( NMAX )
C     INTEGER A(M)
C     DOUBLE PRECISION X ( N )
C     INTEGER PERMUT ( NMAX, N )
C     INTEGER I, J

C     DO 1, I=1, M
C     A(I) = I
C 1   CONTINUE

C     DO 2, I = 1, M


C     END 





      SUBROUTINE PADLIST ( L1, L2, N, INC )
      IMPLICIT NONE
      INTEGER N, INC, I
      INTEGER L1 ( N ), L2 ( N )


      IF (INC.LE.-1) GOTO 50

      IF (ABS(INC).GE.N) THEN
        PRINT*, "Incorrect value of INC, try smaller value"
        GOTO 100
      ENDIF

      DO 10, I=1, N
      IF(I+INC.LE.N) THEN
        L2(I+INC) = L1(I)
      ELSE
        L2(I+INC-N) = L1(I)
      ENDIF 
 10   CONTINUE

      GOTO 100 

 50   CONTINUE

      DO 11, I=1, N
      IF(I+INC.GE.1) THEN
        L2(I+INC) = L1(I)
      ELSE
        L2(I+INC+N) = L1(I)
      ENDIF 
 11   CONTINUE

100   CONTINUE

      RETURN
      END





      FUNCTION NEXTP(N,A)
      INTEGER N,A,I,J,K,T
      LOGICAL NEXTP
      DIMENSION A(N)
      I=N-1
   10 IF(A(I).LT.A(I+1)) GO TO 20
      I=I-1
      IF(I.EQ.0) GO TO 20
      GO TO 10
   20 J=I+1
      K=N
   30 T=A(J)
      A(J)=A(K)
      A(K)=T
      J=J+1
      K=K-1
      IF(J.LT.K) GO TO 30
      J=I
      IF(J.NE.0) GO TO 40
      NEXTP=.FALSE.
      RETURN
   40 J=J+1
      IF(A(J).LT.A(I)) GO TO 40
      T=A(I)
      A(I)=A(J)
      A(J)=T
      NEXTP=.TRUE.
      END


      FUNCTION PDISTANCE (N, M, P1, P2, PERM)
      IMPLICIT NONE
      INTEGER N, M
      INTEGER P1( N*M ), P2( N*M ), P3( N*M ), P4( N*M ), P11( N*M )
      INTEGER PERM ( 100, N * M )
      INTEGER I, J, K,  NM, FACT, DIFF1, DIFF11
      INTEGER Q

      Q = 1000

      NM = N * M

      FACT = 1
      DO 1, I = 1, N
 1    FACT = I * FACT

      P11 (1) = P1 (1)
      DO 10, I = 2, NM
 10   P11(I) = P1(NM+2-I) 


      DO 2, I = 1, FACT

      DO 3, J = 1, NM    
 3    P3 ( J ) = PERM ( I , P2 ( J ) )
C     PRINT *, "P3 -> " , (P3(K), K=1,NM) 

      DO 2, J = 0, NM - 1
      CALL PADLIST ( P3, P4, NM, J )
      DIFF1  = 0
      DIFF11 = 0
      DO 4, K = 1, NM
      DIFF1  = DIFF1  + ABS ( P4 ( K ) - P1  ( K ) )
 4    DIFF11 = DIFF11 + ABS ( P4 ( K ) - P11 ( K ) )
C     PRINT *,  (P4(K), K=1,NM), DIFF
      
      DIFFMIN = MIN ( DIFF1, DIFF11 )

 2    CONTINUE 

      PDISTANCE = DIFFMIN

      RETURN
      END