PROGRAM D10R13 real etime ! Declare the type of etime() real elapsed(2) ! For receiving user and system time real total ! For receiving total time PARAMETER (NCITY=900) DIMENSION X(NCITY),Y(NCITY),IORDER(NCITY) C CREATE POINTS OF SALE IDUM=-111 DO 11 I=1,NCITY X(I)=RAN3(IDUM) Y(I)=RAN3(IDUM) IORDER(I)=I 11 CONTINUE CALL ANNEAL(X,Y,IORDER,NCITY) OPEN(UNIT=6,FILE='ANNEAL.OUT',STATUS='NEW') WRITE(6,*) '*** SYSTEM FROZEN ***' WRITE(6,*) 'FINAL PATH:' WRITE(6,'(1X,T3,A,T13,A,T23,A)') 'CITY','X','Y' DO 12 I=1,NCITY II=IORDER(I) WRITE(6,'(1X,I4,2F10.4)') II,X(II),Y(II) 12 CONTINUE total = etime(elapsed) WRITE(6,*) 'End: total=', total, ' user=', elapsed(1), ' system=', elapsed(2) END SUBROUTINE ANNEAL(X,Y,IORDER,NCITY) C REV. 12/13/85 DIMENSION X(NCITY),Y(NCITY),IORDER(NCITY),N(6) LOGICAL ANS ALEN(X1,X2,Y1,Y2)=SQRT((X2-X1)**2+(Y2-Y1)**2) NOVER=100*NCITY NLIMIT=10*NCITY TFACTR=0.9 PATH=0.0 T=0.5 DO 11 I=1,NCITY-1 I1=IORDER(I) I2=IORDER(I+1) PATH=PATH+ALEN(X(I1),X(I2),Y(I1),Y(I2)) 11 CONTINUE I1=IORDER(NCITY) I2=IORDER(1) PATH=PATH+ALEN(X(I1),X(I2),Y(I1),Y(I2)) IDUM=-1 ISEED=111 DO 15 J=1,100 NSUCC=0 DO 13 K=1,NOVER 12 N(1)=1+INT(NCITY*RAN3(IDUM)) N(2)=1+INT((NCITY-1)*RAN3(IDUM)) IF (N(2).GE.N(1)) N(2)=N(2)+1 IDEC=IRBIT1(ISEED) NN=1+MOD((N(1)-N(2)+NCITY-1),NCITY) IF (NN.LT.3) GOTO 12 IF (IDEC.EQ.0) THEN N(3)=N(2)+INT(ABS(NN-2)*RAN3(IDUM))+1 N(3)=1+MOD(N(3)-1,NCITY) CALL TRNCST(X,Y,IORDER,NCITY,N,DE) CALL METROP(DE,T,ANS) IF (ANS) THEN NSUCC=NSUCC+1 PATH=PATH+DE CALL TRNSPT(IORDER,NCITY,N) ENDIF ELSE CALL REVCST(X,Y,IORDER,NCITY,N,DE) CALL METROP(DE,T,ANS) IF (ANS) THEN NSUCC=NSUCC+1 PATH=PATH+DE CALL REVERS(IORDER,NCITY,N) ENDIF ENDIF IF (NSUCC.GE.NLIMIT) GOTO 14 13 CONTINUE 14 WRITE(6,*) WRITE(6,'(1X,A,F10.6,A,F10.6)') 'T =',T, * ' PATH LENGTH =',PATH WRITE(6,'(1X,A,I6)') 'SUCCESSFUL MOVES: ',NSUCC T=T*TFACTR IF (NSUCC.EQ.0) RETURN 15 CONTINUE RETURN END SUBROUTINE REVCST(X,Y,IORDER,NCITY,N,DE) C REV. 12/13/85 DIMENSION X(NCITY),Y(NCITY),IORDER(NCITY),N(6),XX(4),YY(4) ALEN(X1,X2,Y1,Y2)=SQRT((X2-X1)**2+(Y2-Y1)**2) N(3)=1+MOD((N(1)+NCITY-2),NCITY) N(4)=1+MOD(N(2),NCITY) DO 11 J=1,4 II=IORDER(N(J)) XX(J)=X(II) YY(J)=Y(II) 11 CONTINUE DE=-ALEN(XX(1),XX(3),YY(1),YY(3))-ALEN(XX(2),XX(4),YY(2),YY(4)) * +ALEN(XX(1),XX(4),YY(1),YY(4))+ALEN(XX(2),XX(3),YY(2),YY(3)) RETURN END SUBROUTINE REVERS(IORDER,NCITY,N) DIMENSION IORDER(NCITY),N(6) NN=(1+MOD(N(2)-N(1)+NCITY,NCITY))/2 DO 11 J=1,NN K=1+MOD((N(1)+J-2),NCITY) L=1+MOD((N(2)-J+NCITY),NCITY) ITMP=IORDER(K) IORDER(K)=IORDER(L) IORDER(L)=ITMP 11 CONTINUE RETURN END SUBROUTINE TRNCST(X,Y,IORDER,NCITY,N,DE) DIMENSION X(NCITY),Y(NCITY),IORDER(NCITY),N(6),XX(6),YY(6) ALEN(X1,X2,Y1,Y2)=SQRT((X2-X1)**2+(Y2-Y1)**2) N(4)=1+MOD(N(3),NCITY) N(5)=1+MOD((N(1)+NCITY-2),NCITY) N(6)=1+MOD(N(2),NCITY) DO 11 J=1,6 II=IORDER(N(J)) XX(J)=X(II) YY(J)=Y(II) 11 CONTINUE DE=-ALEN(XX(2),XX(6),YY(2),YY(6))-ALEN(XX(1),XX(5),YY(1),YY(5)) * -ALEN(XX(3),XX(4),YY(3),YY(4))+ALEN(XX(1),XX(3),YY(1),YY(3)) * +ALEN(XX(2),XX(4),YY(2),YY(4))+ALEN(XX(5),XX(6),YY(5),YY(6)) RETURN END SUBROUTINE TRNSPT(IORDER,NCITY,N) PARAMETER(MAXCIT=1000) DIMENSION IORDER(NCITY),JORDER(MAXCIT),N(6) M1=1+MOD((N(2)-N(1)+NCITY),NCITY) M2=1+MOD((N(5)-N(4)+NCITY),NCITY) M3=1+MOD((N(3)-N(6)+NCITY),NCITY) NN=1 DO 11 J=1,M1 JJ=1+MOD((J+N(1)-2),NCITY) JORDER(NN)=IORDER(JJ) NN=NN+1 11 CONTINUE IF (M2.GT.0) THEN DO 12 J=1,M2 JJ=1+MOD((J+N(4)-2),NCITY) JORDER(NN)=IORDER(JJ) NN=NN+1 12 CONTINUE ENDIF IF (M3.GT.0) THEN DO 13 J=1,M3 JJ=1+MOD((J+N(6)-2),NCITY) JORDER(NN)=IORDER(JJ) NN=NN+1 13 CONTINUE ENDIF DO 14 J=1,NCITY IORDER(J)=JORDER(J) 14 CONTINUE RETURN END SUBROUTINE METROP(DE,T,ANS) PARAMETER(JDUM=1) LOGICAL ANS ANS=(DE.LT.0.0).OR.(RAN3(JDUM).LT.EXP(-DE/T)) RETURN END FUNCTION RAN3(IDUM) C IMPLICIT REAL*4(M) C PARAMETER (MBIG=4000000.,MSEED=1618033.,MZ=0.,FAC=2.5E-7) PARAMETER (MBIG=1000000000,MSEED=161803398,MZ=0,FAC=1.E-9) DIMENSION MA(55) DATA IFF /0/ IF(IDUM.LT.0.OR.IFF.EQ.0)THEN IFF=1 MJ=MSEED-IABS(IDUM) MJ=MOD(MJ,MBIG) MA(55)=MJ MK=1 DO 11 I=1,54 II=MOD(21*I,55) MA(II)=MK MK=MJ-MK IF(MK.LT.MZ)MK=MK+MBIG 11 CONTINUE DO 13 K=1,4 DO 12 I=1,55 MA(I)=MA(I)-MA(1+MOD(I+30,55)) IF(MA(I).LT.MZ)MA(I)=MA(I)+MBIG 12 CONTINUE 13 CONTINUE INEXT=0 INEXTP=31 IDUM=1 ENDIF INEXT=INEXT+1 IF(INEXT.EQ.56)INEXT=1 INEXTP=INEXTP+1 IF(INEXTP.EQ.56)INEXTP=1 MJ=MA(INEXT)-MA(INEXTP) IF(MJ.LT.MZ)MJ=MJ+MBIG MA(INEXT)=MJ RAN3=MJ*FAC RETURN END FUNCTION IRBIT1(ISEED) LOGICAL NEWBIT PARAMETER (IB1=1,IB2=2,IB5=16,IB18=131072) NEWBIT=IAND(ISEED,IB18).NE.0 IF(IAND(ISEED,IB5).NE.0)NEWBIT=.NOT.NEWBIT IF(IAND(ISEED,IB2).NE.0)NEWBIT=.NOT.NEWBIT IF(IAND(ISEED,IB1).NE.0)NEWBIT=.NOT.NEWBIT IRBIT1=0 ISEED=IAND(ISHFT(ISEED,1),NOT(IB1)) IF(NEWBIT)THEN IRBIT1=1 ISEED=IOR(ISEED,IB1) ENDIF RETURN END