C PROGRAM NAME: RODUN08.FOR
C ADAPTED FOR THE FREE GNU G77 FORTRAN COMPILER IN DEC 2008
INTEGER TOTPRO,COUNT,DPSEQ
COMMON /RIT/ IPSET(20,50)
c INDEX NO OF PROJ IN FIXED PROJECT SEQ
C COMMON AREA FOR PROGRAMS RODUN,FINPRI
1,SCALE(20)
c !SCALE FACTOR FOR EACH OF THE FIXED PR SETS
2,AVAER(20)
c !AVERAGE ABSOLUTE ERROR IN HEUR ALG # I
3,AVER(20)
c !AVERAGE RELATIVE ERROR IN HEUR ALG # I
4,FAILR(20)
c !FAILURE RATE IN HEUR ALG # I
CHARACTER*1 ZNAME(60),PIND, IAFLAG, IGFLAG
COMMON /DM/ T(100),D(100),NDEDAT,TAR,DMULT,DKON
COMMON /RI/ DPCTOT,COUNT
c !DISCOUNTED TOTAL COST IN DP SEQ,COUNTER OF SPS
1,DPSEQ(20),PPP,PP
C !DYNPRO SEQUENCE,LOGICAL PRINT FLAGS
1,MSEQ(20)
c !SCRATCH SEQUENCE
1,IDOM
c !dominANCE INDICATOR(IDOM=1:DOMINANCE)
1,ERTO(20)
c !ACCUMULATED ERROR IN SEQ ALG
1,ERRMAX(20)
c !MAX ERROR IN SEQ ALG
1,NERR(20)
c !SP# WHERE MAX ERROR OCCURRED
1,KNT(20)
c !ERROR COUNTER FOR EACH SP
1,CT(20)
c !TOTAL DISCOUNTED COST FOR EACH ALG
1,ISEQ(20,20)
c !INDEX NO OF THE PROJECT IN J-TH PLACE
C !INSEQUENCE FOR ALG # I (I,J)
1,ER(20)
c !PERC ERROR IN EACH ALG
1,NOALG
c !NO OF HEURISTIC ALGORITHMS
1,IBETT(20,20)
c !COUNTER FOR SPS WHEN ALG #I IS SUPERIOR TO ALG#J (I,J)
1,IBEST(20)
c !COUNTER FOR SPS WHERE ALG #I (EXCEPT DP) IS BEST
COMMON /RAN2/LL,MM
c !RANDOM NO SEEDS FOR RAND SEQUENCING
COMMON /RAN/ II,L,M
C !RAND NO SEEDS FOR RAND GEN & SEL OF PROJ
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
LOGICAL PRID,RASEL,P,PP,PPP,PKI,PCI,PSI
IAFLAG ='N'
OPEN(UNIT=9,FILE='DEM.DAT')
OPEN(UNIT=7,FILE='RES.DAT')
OPEN(UNIT=11,FILE='ADEBUG.DAT')
P=.FALSE.
PP=P
PPP=P
PKI=P
PCI=P
PSI=P
NOALG=7
C ! NO OF HEOURISTIC SEQUENCING ALGORITHMS
IGFLAG='N'
WRITE(11,990)
C CALL GETTIM (IHR,IMIN1,ISEC,IS100) OUTDATED CALL FOR CPU TIME
CALL CPU_Time(SEC1)
C WRITE(*,1) OLD VERION WITH SCREEN INPUT
C1 FORMAT(10X,'THE WATER RESOURCES SEQUENCING PROBLEM; '//
C 110X,'ENTER INTEREST RATE IN % ')
990 FORMAT(80('-'))
READ(9,2) XINT
C READ INTEREST RATE FROM DEM.DAT
C WRITE(*,990)
2 FORMAT(F3.0)
write(11,*) xint
C WRITE(*,376) ' DETAILED WRITEOUT TO A FILE DESIRED ? '
C WRITE(*,376) ' R(I)-SEQ DET.FOR ALL SEQ. PROB :A '
C WRITE(*,376) ' MINIMUM OUTPUT FOR ALL SEQ PROB :B '
C WRITE(*,376) ' OUTPUT ONLY WHEN R(I) IS NONOPT :C '
C WRITE(*,376) ' COMPACT DETAILS FOR ALL SEQ PR. :D '
C WRITE(*,376) ' K(I)SEQ DETAILS FOR ALL SEQ PR. :E '
C WRITE(*,376) ' C/X SEQ DETAILS FOR ALL SEQ PR. :F '
C WRITE(*,371) ' S(I)SEQ DETAILS FOR ALL SEQ PR. :G '
371 FORMAT(1X,A50,4X)
READ(9,372) PIND
READ(9,372)
READ(9,372)
READ(9,372)
READ(9,372)
READ(9,372)
READ(9,372)
READ(9,372)
C DUMMY READING OF EXPLANATIONS LINES IN DEM.DAT
372 FORMAT(A1)
write(11,*) pind
376 FORMAT (1X,A50)
IF (PIND.EQ.'A')P=.TRUE.
IF (PIND.EQ.'B')PP=.TRUE.
IF (PIND.EQ.'C')PPP=.TRUE.
IF (PIND.EQ.'E') PKI=.TRUE.
IF (PIND.EQ.'F') PCI=.TRUE.
IF (PIND.EQ.'G') PSI=.TRUE.
PRID=P
C WRITE(*,3)
C3 FORMAT(2X,'ENTER FIRST AND LAST INDEX FOR REQUIRED DEMAND ',
C 1'FUNCTION')
READ(9,4) NBEGD,NENDD
4 FORMAT(2I5)
write(11,*) nbegd,nendd
C WRITE(*,5)
OPEN(UNIT=10,FILE='BOUT.DAT')
5 FORMAT(2X,'DO YOU WANT ACTUAL COST DATA?(Y OR N)')
C READ(9,6) IAFLAG
C YOU ALWAYS WANT ACTUAL COST DATA
IAFLAG='Y'
6 FORMAT(A1)
write(11,*) IAFLAG
IF(IAFLAG.NE.'Y')GO TO 9
OPEN(UNIT=8,FILE='PROCOST.DAT')
READ(8,201)TOTPRO,NOFIX,LFIX,SCA
201 FORMAT(3I3,F10.0)
WRITE(10,202)TOTPRO,NOFIX,LFIX,SCA
C WRITE(*,202) TOTPRO,NOFIX,LFIX,SCA
C DONT WRITE THIS ON THE SCREEN
202 FORMAT(1H1,20X,'TOTAL NO OF PROJECTS IN FILE PROCOST='
1,I3/20X,' NO. OF'
1,' FIXED PROJECT SETS =',I3/20X
1,'LENGTH OF EACH PROJECT SET= ',I3
1,/20X,'SCALE FACTOR FOR ADJUSTING THE SIZE OF ALL PROJECTS='
1,G10.2)
IF(NOFIX.EQ.0)GO TO 228
DO 203 I=1,NOFIX
READ(8,204)SCALE(I),(IPSET(I,J),J=1,LFIX)
IF(P) WRITE(10,207)I,SCALE(I),(IPSET(I,J),J=1,LFIX)
203 CONTINUE
228 IF(P) WRITE(10,208)
IF(P) WRITE(10,910)
204 FORMAT(F10.0,20I3)
DO 205 I=1,TOTPRO
READ(8,206) CC(I),XX(I),(ZNAME(J),J=1,40)
IF(P) WRITE(10,209) I,CC(I),XX(I),(ZNAME(J),J=1,40)
205 CONTINUE
206 FORMAT(2F10.0,40A1)
CLOSE(8)
207 FORMAT(3X,' IN FIXED PROJECT SEQUENCE #',I3,' IS SCALE FACTOR:',
1G10.3,' AND THE PROJECTS:',130(I3,','))
208 FORMAT(//'*****************LIST
1OF ALL ACTUAL PROJECTS*********')
WRITE(10,909)
909 FORMAT(80('*'))
910 FORMAT(1X,T5,'PROJECT NO',T16,'COST',T31,'CATACITY',T50
1 ,'NAME OF PROJECT')
209 FORMAT (1X,T5,I10,T16,F10.0,T31,F10.0,T50,60A1)
WRITE(10,909)
C WRITE(*,7)
7 FORMAT(20X,'ENTER NO OF RANDOMLY SELECTED DATA SETS'
1/25X,' FROM ACTUAL DATA ')
READ(9,46) NCCSET
46 FORMAT(I5)
C WRITE(*,8)
8 FORMAT(20X,'ENTER NO OF PROJECTS IN EACH SET ')
READ(9,46) NPRO
C WRITE(*,22)
22 FORMAT(20X,'ENTER INDEX NO`S OF FIRST AND LAST PROJECT TO BE'/
120X,' INCLUDED IN THE RANGE OF RANDOMLY SELECTED PROJECTS ')
READ(9,4) NBEGP,NENDP
9 CONTINUE
C WRITE(*,10)
write(11,*) NBEGP,NENDP
10 FORMAT(3X,'DO YOU WANT GENERATED COST DATA?(Y OR N) ')
C READ(*,6) IGFLAG
IGFLAG='Y'
IF((IGFLAG.NE.'Y').AND.(IAFLAG.NE.'Y'))STOP 'NO PROCESSING'
IF(IGFLAG.NE.'Y')GO TO 47
C WRITE(*,23)
23 FORMAT(20X,'ENTER # OF RANDOMLY GENERATED DATA SETS ')
READ(9,46) NCGSET
C WRITE(*,24)
24 FORMAT(20X,'ENTER # OF PROJECTS IN EACH SET ')
READ(9,46) NPROG
47 IF (NBEGD.LE.1)GO TO 48
DO 14 I=1,NBEGD-1
CALL DEMAND(PRID)
C THE SUBROUTINE DEAMND READS IN A NEW DEMAND FUNCTION
C FROM THE FILE DEM.DAT
14 CONTINUE
C******************************************************************************
48 DO 15 I=NBEGD,NENDD
RASEL=.FALSE.
L=245791
M=573
LL=167351
MM=12538
DO 36 IK=1,20
ER(IK)=0
IBEST(IK)=0
ERTO(IK)=0
KNT(IK)=0
NERR(IK)=1
ERRMAX(IK)=0
DO 601 JN=1,20
IBETT(IK,JN)=0
601 CONTINUE
36 CONTINUE
KSKR=0
COUNT =0
CALL DEMAND(PRID)
IF(IAFLAG.NE.'Y')GO TO 117
DO 16 II=1,NOFIX+NCCSET
IF(II-NOFIX)301,301,302
301 RASEL=.FALSE.
DO 17 III=1,LFIX
C(III)=CC(IPSET(II,III))
X(III)=XX(IPSET(II,III))*SCALE(II)
D WRITE(10,801) III,C(III),X(III)
D801 FORMAT( ' III=',I5,' C(III)=',F10.0,' X(III)=',F10.0)
17 CONTINUE
CCCCCCCCCCCCCCCCCCCC HERE WE HAVE THE NUCLEUS FOR ACTUAL DATACCCCCCCCCCC
KK=LFIX
303 COUNT=COUNT+1
IF(PP)WRITE(10,990)
IF((.NOT.RASEL).AND.PP)WRITE(10,31) COUNT,I,II
IF(RASEL.AND.PP)WRITE(10,80) COUNT,I,II-NOFIX
31 FORMAT(' SEQUENCING PROBLEM #',I3,' DEMAND FUNCTION #',I3,
1' FIXED PROJECT SEQUENCE #',I3)
80 FORMAT(' SEQUENCING PROBLEM #',I3,' DEMAND FUNCTION #',I3,
1' RANDOMLY SELECTED SEQUENCE #',I3)
write(11,*) 'DP byrjun'
CALL DOMIN(IDOM)
C------------------------ DYNAMIC PROGRAMING SEQUENCING-----------------
CALL XDPROD(DPCTOT,DPSEQ)
C--------------------------K(I) SEQUENCING------------------------------
CALL KISEQ(CT(3),MSEQ,3,PKI)
CALL COMPR(DPSEQ,MSEQ,KOMP)
DO 86 IK=1,KK
ISEQ(3,IK)=MSEQ(IK)
86 CONTINUE
CALL SEQNOT(3,KOMP)
C--------------------------TMR SEQUENCING--------------------------
CALL KISEQ(CT(2),MSEQ,2,P)
DO 87 IK=1,KK
ISEQ(2,IK)=MSEQ(IK)
87 CONTINUE
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(2,KOMP)
KSKR=KOMP
C-------------------------- C/X SEQUENGING-------------------------
CALL KISEQ(CT(1),MSEQ,1,PCI)
DO 41 IK=1,KK
ISEQ(1,IK)=MSEQ(IK)
41 CONTINUE
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(1,KOMP)
C-------------------------- S(I) SEQUENCING------------------------
CALL KISEQ(CT(4),MSEQ,4,PSI)
DO 42 IK=1,KK
ISEQ(4,IK)=MSEQ(IK)
42 CONTINUE
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(4,KOMP)
C---------------------------RANDOM SEQUENCING----------------------
CALL KISEQ(CT(5),MSEQ,5,.FALSE.)
DO 33 IK=1,KK
ISEQ(5,IK)=MSEQ(IK)
33 CONTINUE
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(5,KOMP)
C---------------------------MAX(K,R) SEQUENCING-------------------
CALL MAXALG(3,2,6)
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(6,KOMP)
C---------------------------MAX(S,R) SEQUENCING--------------------
CALL MAXALG(4,2,7)
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(7,KOMP)
C---------------------------PERFORMANCE EVALUATION-------------------
CALL PERFOR
IF(PPP.AND.(KSKR.EQ.1))CALL SKRIFT
C !WRITE WHEN R(I) NONOPT
IF(PIND.EQ.'D')CALL SKRIFT
C !COMPACT WRITEOUT
GO TO 43
302 CALL RANPAC(PRID)
KK=NPRO
DO 76 I5=1,KK
X(I5)=X(I5)*SCA
76 CONTINUE
C !SCALE ALL PROJECTS
RASEL=.TRUE.
GO TO 303
43 CONTINUE
16 CONTINUE
C SEQUENCING STARTS FOR GENERATED DATA
117 IF(IGFLAG.NE.'Y')GO TO 15
WRITE(10,990)
WRITE(10,63)
DO 61 II=1,NCGSET
CALL PROGEN(0,PRID)
KK=NPROG
DO 75 I5=1,KK
C !SCALE ALL PROJECTS
X(I5)=X(I5)*SCA
75 CONTINUE
COUNT=COUNT+1
IF(PP)WRITE(10,990)
IF(PP)WRITE(10,62) COUNT,I,II
63 FORMAT(10X,'*************RANDOMLY GENERATED DATA*************')
62 FORMAT(' SEQUENCING PROBLEM #:',I3,
1 ' DEMAND FUNCTION INDEX NO:',
1 I3,' RANDOMLY GENERATED SEQUENCE NO:',I3)
CALL DOMIN(IDOM)
C------------------------- DYNAMIC PROGRAMING SEQUENCING---------------------
CALL XDPROD(DPCTOT,DPSEQ)
C---------------------------K(I) SEQUENCING---------------------------------
CALL KISEQ(CT(3),MSEQ,3,PKI)
DO 88 IK=1,KK
ISEQ(3,IK)=MSEQ(IK)
88 CONTINUE
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(3,KOMP)
C--------------------------TMR SEQUENCING-----------------------------------
CALL KISEQ(CT(2),MSEQ,2,P)
DO 89 IK=1,KK
ISEQ(2,IK)=MSEQ(IK)
89 CONTINUE
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(2,KOMP)
KSKR=KOMP
C-------------------------- C/X SEQUENGING-------------------------
CALL KISEQ(CT(1),MSEQ,1,PCI)
DO 81 IK=1,KK
ISEQ(1,IK)=MSEQ(IK)
81 CONTINUE
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(1,KOMP)
C-------------------------- S(I) SEQUENCING------------------------
CALL KISEQ(CT(4),MSEQ,4,PSI)
DO 82 IK=1,KK
ISEQ(4,IK)=MSEQ(IK)
82 CONTINUE
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(4,KOMP)
C---------------------------RANDOM SEQUENCING----------------------
CALL KISEQ(CT(5),MSEQ,5,.FALSE.)
DO 83 IK=1,KK
ISEQ(5,IK)=MSEQ(IK)
83 CONTINUE
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(5,KOMP)
C---------------------------MAX(K,R) SEQUENCING-------------------
CALL MAXALG(3,2,6)
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(6,KOMP)
C---------------------------MAX(S,R) SEQUENCING--------------------
CALL MAXALG(4,2,7)
CALL COMPR(DPSEQ,MSEQ,KOMP)
CALL SEQNOT(7,KOMP)
C----------------------------PERFORMANCE EVALUATION--------------------
CALL PERFOR
C----------------------------PRINTOUT OF RESULTS ---------------------
IF(PPP.AND.(KSKR.EQ.1))CALL SKRIFT
C !WRITE WHEN R(I) NONOPT
IF(PIND.EQ.'D')CALL SKRIFT
C !COMPACT WRITEOUT
61 CONTINUE
C---------------------------PREPARE AND PRINT DATA-------------------------
XCOUNT=COUNT
DO 761 IK=1,NOALG
AVAER(IK)=ERTO(IK)/XCOUNT
XKNT=KNT(IK)
IF(KNT(IK).NE.0)AVER(IK)=ERTO(IK)/XKNT
IF(KNT(IK).EQ.0)AVER(IK)=0
FAILR(IK)=XKNT/XCOUNT
761 CONTINUE
CALL FINPRI(I)
C------------------NEXT DEMAND FUNCTION-----------------------------------
15 CONTINUE
C CALL GETTIM (IHR,IMIN2,ISEC,IS100) old version
CALL CPU_Time(SEC2)
DELTA=(SEC2-SEC1)/60
C DELTA=IMIN2-IMIN1 old version
WRITE(10,261) DELTA
261 FORMAT(5X,'RUNNING TIME=',F15.2,' MINUTES')
CLOSE(6)
CLOSE(7)
CLOSE(9)
STOP'FINISHED OK!'
END
SUBROUTINE MAXALG(I1,I2,IOUT)
C----------------------------------------------------------------------
C THE SUBROUTINE MAXALG COMBINES ALGORITHM #I1 AND ALGO-
C RITHM #I2 AND CREATES A NEW ALGO #IOUT. FOR EACH SP THE
C BETTER ALGO OF THE 2 (#I1 & #I2) IS SELQCTED.
C----------------------------------------------------------------------
COMMON /RI/ DPCTOT,COUNT
c !DISCOUNTED TOTAL COST IN DP SEQ,COUNTER OF SPS
1,DPSEQ(20),PPP,PP
C !DYNPRO SEQUENCE,LOGICAL PRINT FLAGS
1,MSEQ(20)
c !SCRATCH SEQUENCE
1,IDOM
c !DOMINANCE INDICATOR(IDOM=1:DOMINANCE)
1,ERTO(20)
c !ACCUMULATED ERROR IN SEQ ALG
1,ERRMAX(20)
c !MAX ERROR IN SEQ ALG
1,NERR(20)
c !SP# WHERE MAX ERROR OCCURRED
1,KNT(20)
c !ERROR COUNTER FOR EACH SP
1,CT(20)
c !TOTAL DISCOUNTED COST FOR EACH ALG
1,ISEQ(20,20)
c !INDEX NO OF THE PROJECT IN J-TH PLACE
C !INSEQUENCE FOR ALG # I (I,J)
1,ER(20)
c !PERC ERROR IN EACH ALG
1,NOALG
c !NO OF HEURISTIC ALGORITHMS
1,IBETT(20,20)
c !COUNTER FOR SPS WHEN ALG #I IS SUPERIOR TO ALG#J (I,J)
1,IBEST(20)
c !COUNTER FOR SPS WHERE ALG #I (EXCEPT DP) IS BEST
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
INTEGER DPSEQ,COUNT
I=I2
IF(ER(I1).LT.ER(I2))I=I1
CT(IOUT)=CT(I)
DO 10 IK=1,KK
ISEQ(IOUT,IK)=ISEQ(I,IK)
MSEQ(IK)=ISEQ(I,IK)
10 CONTINUE
RETURN
END
SUBROUTINE SKRIFT
C-------------------------------------------------------------------------
C THIS SUB WRITES COMPACT DETAILS ABOUT ALL SP
C-------------------------------------------------------------------------
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
INTEGER DPSEQ,COUNT
COMMON /RI/ DPCTOT,COUNT
c !DISCOUNTED TOTAL COST IN DP SEQ,COUNTER OF SPS
1,DPSEQ(20),PPP,PP
C !DYNPRO SEQUENCE,LOGICAL PRINT FLAGS
1,MSEQ(20)
c !SCRATCH SEQUENCE
1,IDOM
c !DOMINANCE INDICATOR(IDOM=1:DOMINANCE)
1,ERTO(20)
c !ACCUMULATED ERROR IN SEQ ALG
1,ERRMAX(20)
c !MAX ERROR IN SEQ ALG
1,NERR(20)
c !SP# WHERE MAX ERROR OCCURRED
1,KNT(20)
c !ERROR COUNTER FOR EACH SP
1,CT(20)
c !TOTAL DISCOUNTED COST FOR EACH ALG
1,ISEQ(20,20)
c !INDEX NO OF THE PROJECT IN J-TH PLACE
C !INSEQUENCE FOR ALG # I (I,J)
1,ER(20)
c !PERC ERROR IN EACH ALG
1,NOALG
c !NO OF HEURISTIC ALGORITHMS
1,IBETT(20,20)
c !COUNTER FOR SPS WHEN ALG #I IS SUPERIOR TO ALG#J (I,J)
1,IBEST(20)
c !COUNTER FOR SPS WHERE ALG #I (EXCEPT DP) IS BEST
LOGICAL PPP,PP
WRITE(10,10) COUNT,(C(I),I=1,KK)
10 FORMAT(1X,'SP#: ',I4,' C(I): ',20F10.0)
WRITE(10,11)(X(I),I=1,KK)
11 FORMAT(10X,' X(I): ',20F10.0)
WRITE(10,12) DPCTOT,(DPSEQ(I),I=1,KK)
12 FORMAT(10X,' DPCTOT= ',F10.0,17X,' DPSEQ= ',20(I3,','))
13 FORMAT(11X,'R(I)SEQ;CT= ',F10.0,' ER= ',F5.2,
1' % RISEQ= ',20(I3,','))
14 FORMAT(11X,'K(I)SEQ;CT= ',F10.0,' ER= ',F5.2,
1' % KISEQ= ',20(I3,','))
15 FORMAT(11X,'C/X SEQ;CT= ',F10.0,' ER= ',F6.2,
1' % CXSEQ= ',20(I3,','))
16 FORMAT(11X,'S(I)SEQ;CT= ',F10.0,' ER= ',F6.2,
1' % SISEQ= ',20(I3,','))
17 FORMAT(11X,'RAN SEQ;CT= ',F10.0,' ER= ',F6.2,
1' % RANSEQ= ',20(I3,','))
18 FORMAT(11X,'MAX(K,R)CT= ',F10.0,' ER= ',F6.2,
1' % SEQ= ',20(I3,','))
19 FORMAT(11X,'MAX(S,R)CT= ',F10.0,' ER= ',F6.2,
1' % SEQ= ',20(I3,','))
WRITE(10,15) CT(1),ER(1),(ISEQ(1,I),I=1,KK)
WRITE(10,13)CT(2),ER(2),(ISEQ(2,I),I=1,KK)
WRITE(10,14)CT(3),ER(3),(ISEQ(3,I),I=1,KK)
WRITE(10,16)CT(4),ER(4),(ISEQ(4,I),I=1,KK)
WRITE(10,17)CT(5),ER(5),(ISEQ(5,I),I=1,KK)
WRITE(10,18)CT(6),ER(6),(ISEQ(6,I),I=1,KK)
WRITE(10,19)CT(7),ER(7),(ISEQ(7,I),I=1,KK)
50 RETURN
END
SUBROUTINE SEQNOT(JALG,KOMP)
C--------------------------------------------------------------------------
C THIS SUBROUTINE PERFORMS PROCESSING OF THE ERROR OCCURRING
C IN HEURISTIC ALGORITHM # JALG.THE VARIABLE KOMP IS AN INPUT
C TO THE ROUTINE. KOMP=1 IF THE ALGORITHM IS NOT OPTIMAL,BUT
C ZERO OTHERWISE.
C-------------------------------------------------------------------------
COMMON /RI/ DPCTOT,COUNT
c !DISCOUNTED TOTAL COST IN DP SEQ,COUNTER OF SPS
1,DPSEQ(20),PPP,PP
C !DYNPRO SEQUENCE,LOGICAL PRINT FLAGS
1,MSEQ(20)
c !SCRATCH SEQUENCE
1,IDOM
c !DOMINANCE INDICATOR(IDOM=1:DOMINANCE)
1,ERTO(20)
c !ACCUMULATED ERROR IN SEQ ALG
1,ERRMAX(20)
c !MAX ERROR IN SEQ ALG
1,NERR(20)
c !SP# WHERE MAX ERROR OCCURRED
1,KNT(20)
c !ERROR COUNTER FOR EACH SP
1,CT(20)
c !TOTAL DISCOUNTED COST FOR EACH ALG
1,ISEQ(20,20)
c !INDEX NO OF THE PROJECT IN J-TH PLACE
C !INSEQUENCE FOR ALG # I (I,J)
1,ER(20)
c !PERC ERROR IN EACH ALG
1,NOALG
c !NO OF HEURISTIC ALGORITHMS
1,IBETT(20,20)
c !COUNTER FOR SPS WHEN ALG #I IS SUPERIOR TO ALG#J (I,J)
1,IBEST(20)
c !COUNTER FOR SPS WHERE ALG #I (EXCEPT DP) IS BEST
INTEGER DPSEQ,COUNT
LOGICAL PPP,PP
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
IF(KOMP.EQ.0)GO TO 30
ER(JALG)=(CT(JALG)-DPCTOT)/DPCTOT*100
KNT(JALG)=KNT(JALG)+1
ERTO(JALG)=ERTO(JALG)+ER(JALG)
IF(ER(JALG).LE.ERRMAX(JALG))GO TO 20
NERR(JALG)=COUNT
ERRMAX(JALG)=ER(JALG)
GO TO 20
30 ER(JALG)=0
20 RETURN
END
SUBROUTINE PERFOR
C----------------------------------------------------------------------------
C THE SUBROUTINE PERFOR COMPARES THE PERFORMANCE OF ALL HEURISTIC
C ALGORITHMS AGAINST EACH OTHER.
C IT INCREMENTS THE COUNTERS IN THE MATRIX IBETT(I,J) FOR AN SP WHEN
C ALGORITHM #I IS SUPERIOR TO ALGORITHM #J.
C IT INCREMENTS THE VECTOR IBEST, WHEN ALGORITHM # I GIVES THE BEST
C SOLUTION OUT OF THE HEURISTIC ALGORITHMS ONLY.
C-------------------------------------------------------------------------------
COMMON /RI/ DPCTOT,COUNT
c !DISCOUNTED TOTAL COST IN DP SEQ,COUNTER OF SPS
1,DPSEQ(20),PPP,PP
C !DYNPRO SEQUENCE,LOGICAL PRINT FLAGS
1,MSEQ(20)
c !SCRATCH SEQUENCE
1,IDOM
c !DOMINANCE INDICATOR(IDOM=1:DOMINANCE)
1,ERTO(20)
c !ACCUMULATED ERROR IN SEQ ALG
1,ERRMAX(20)
c !MAX ERROR IN SEQ ALG
1,NERR(20)
c !SP# WHERE MAX ERROR OCCURRED
1,KNT(20)
c !ERROR COUNTER FOR EACH SP
1,CT(20)
c !TOTAL DISCOUNTED COST FOR EACH ALG
1,ISEQ(20,20)
c !INDEX NO OF THE PROJECT IN J-TH PLACE
C !INSEQUENCE FOR ALG # I (I,J)
1,ER(20)
c !PERC ERROR IN EACH ALG
1,NOALG
c !NO OF HEURISTIC ALGORITHMS
1,IBETT(20,20)
c !COUNTER FOR SPS WHEN ALG #I IS SUPERIOR TO ALG#J (I,J)
1,IBEST(20)
c !COUNTER FOR SPS WHERE ALG #I (EXCEPT DP) IS BEST
INTEGER DPSEQ,COUNT
LOGICAL PPP,PP
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
BESS=100
DO 10 I=1,NOALG
DO 20 J=1,NOALG
IF(ER(I).LT.ER(J))IBETT(I,J)=IBETT(I,J)+1
20 CONTINUE
10 CONTINUE
DO 30 I=1,NOALG
IF (ER(I).LE.BESS) BESS=ER(I)
C !FIND MINIMUM ERROR
30 CONTINUE
DO 40 I=1,NOALG
IF(ER(I).LE.BESS) IBEST(I)=IBEST(I)+1
C !INCREMENT COUNTER FOR ALL ALG
40 CONTINUE
C ! WITH ERROR EQ TO MIN ERROR
RETURN
END
SUBROUTINE FINPRI(ID)
C--------------------------------------------------------------------------
C THE SUBROUTINE FINPRI WRITES FINAL RESULTS ONTO FILE 'RES.DAT'
C FOR DEMAND FUNCTION # ID
C---------------------------------------------------------------------------
COMMON /RIT/ IPSET(20,50)
C !INDEX NO OF PROJ IN FIXED PROJECT SEQ
C COMMON AREA FOR PROGRAMS RODUN,FINPRI
1,SCALE(20)
C !SCALE FACTOR FOR EACH OF THE FIXED PROJECT SETS
2,AVAER(20)
C !AVERAGE ABSOLUTE ERROR IN HEUR ALG # I
3,AVER(20)
C !AVERAGE RELATIVE ERROR IN HEUR ALG # I
4,FAILR(20)
C !FAILURE RATE IN HEUR ALG # I
COMMON /RI/ DPCTOT,COUNT
c !DISCOUNTED TOTAL COST IN DP SEQ,COUNTER OF SPS
1,DPSEQ(20),PPP,PP
C !DYNPRO SEQUENCE,LOGICAL PRINT FLAGS
1,MSEQ(20)
c !SCRATCH SEQUENCE
1,IDOM
c !DOMINANCE INDICATOR(IDOM=1:DOMINANCE)
1,ERTO(20)
c !ACCUMULATED ERROR IN SEQ ALG
1,ERRMAX(20)
c !MAX ERROR IN SEQ ALG
1,NERR(20)
c !SP# WHERE MAX ERROR OCCURRED
1,KNT(20)
c !ERROR COUNTER FOR EACH SP
1,CT(20)
c !TOTAL DISCOUNTED COST FOR EACH ALG
1,ISEQ(20,20)
c !INDEX NO OF THE PROJECT IN J-TH PLACE
C !INSEQUENCE FOR ALG # I (I,J)
1,ER(20)
c !PERC ERROR IN EACH ALG
1,NOALG
c !NO OF HEURISTIC ALGORITHMS
1,IBETT(20,20)
c !COUNTER FOR SPS WHEN ALG #I IS SUPERIOR TO ALG#J (I,J)
1,IBEST(20)
c !COUNTER FOR SPS WHERE ALG #I (EXCEPT DP) IS BEST
INTEGER DPSEQ,COUNT
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
LOGICAL PPP,PP
WRITE(7,3)
3 FORMAT(1H1)
WRITE(7,909)
909 FORMAT(132('*'))
WRITE(7,2)
2 FORMAT(1X)
WRITE(7,1)ID,COUNT,KK,XINT
1 FORMAT(9X,'RESULTS FOR DEMAND FUNCTION # ',
1 I6,' NO OF SEQ PROB.:',I5
1 ,' PROBLEM SIZE:',I3,' INTEREST RATE:',F5.1,' %')
WRITE(7,2)
WRITE(7,4)
4 FORMAT(2X,'ALGORITHM',6X,' C(I)/X(I)'
1,'TMR K(I) S('
1,'I) RAND SEQ MAX(K,R) MAX(S,R)')
WRITE(7,5)(I,I=1,20)
5 FORMAT(2X,'ALGORITHM NO',11X,20(I10))
WRITE(7,2)
WRITE(7,6)(KNT(I),I=1,NOALG)
6 FORMAT(1X,'NO OF SEQ FAILURES',6X,20I10)
WRITE(7,7)(FAILR(I),I=1,NOALG)
7 FORMAT(1X,'FAILURE RATE %',10X,20G10.3)
WRITE(7,2)
8 FORMAT(1X,'AVE REL FAILURE %',10X,20G10.3)
9 FORMAT(1X,'AVE ABS FAILURE %',10X,20G10.3)
10 FORMAT(1X,' MAX FAILURE %',10X,20G10.3)
WRITE(7,8) (AVER(I),I=1,NOALG)
WRITE(7,9) (AVAER(I),I=1,NOALG)
WRITE(7,10)(ERRMAX(I),I=1,NOALG)
WRITE(7,2)
WRITE(7,12)(IBEST(I),I=1,NOALG)
12 FORMAT(1X,'NO OF SPS WHEN BEST (EXC DP)',20(I6,4X))
WRITE(7,13)
13 FORMAT(1X,'NO OF SPS WHEN:')
21 FORMAT(1X,'C(I)/X(I) BETTER (',I2,')',3X,20I10)
23 FORMAT(1X,' K(I) BETTER (',I2,')',3X,20I10)
24 FORMAT(1X,' S(I) BETTER (',I2,')',3X,20I10)
22 FORMAT(1X,' TMR BETTER (',I2,')',3X,20I10)
25 FORMAT(1X,' RAN SEQ BETTER (',I2,')',3X,20I10)
26 FORMAT(1X,' MAX(K,R) BETTER (',I2,')',3X,20I10)
27 FORMAT(1X,' MAX(S,R) BETTER (',I2,')',3X,20I10)
WRITE(7,21)1,(IBETT(1,I),I=1,NOALG)
WRITE(7,22)2,(IBETT(2,I),I=1,NOALG)
WRITE(7,23)3,(IBETT(3,I),I=1,NOALG)
WRITE(7,24)4,(IBETT(4,I),I=1,NOALG)
WRITE(7,25)5,(IBETT(5,I),I=1,NOALG)
WRITE(7,26)6,(IBETT(6,I),I=1,NOALG)
WRITE(7,27)7,(IBETT(7,I),I=1,NOALG)
RETURN
END
SUBROUTINE DOMIN(IDOM)
C----------------------------------------------------------------------------
C THE SUBROUTINE DOMIN CHECKS FOR DOMINANCE AMONG THE KK
C PROJECTS IN THE SP.
C IF DOMINANCE :IDOM=1
C IF NO DOMINANCE: IDOM=0
C----------------------------------------------------------------------------
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
IDOM=0
DO 10 I=1,KK
DO 20 J=1,KK
IF((C(I).LT.C(J)).AND.(X(I).GT.X(J))) IDOM=1
C write(*,*) 'domin loop', i,j
20 CONTINUE
10 CONTINUE
C write(*,*) ' £t £r domin', KK
RETURN
END
SUBROUTINE COMPR(OPTSEQ,TSEQ,KOMP)
C----------------------------------------------------------------------
C THIS SUBROUTINE COMPARES THE SEQUENCE TSEQ AGAINST
C THE OPTIMAL SEQUENCE OPTSEQ.
C IF OPTSEQ=TSEQ :KOMP=0
C IF OPTSEQ NOT= TSEQ:KOMP=1
C------------------------------------------------------------------------
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
INTEGER TSEQ(20),OPTSEQ(20)
KOMP=0
DO 10 I=1,KK
IF(TSEQ(I).NE.OPTSEQ(I))KOMP=1
10 CONTINUE
RETURN
END
C FUNCTION XRAN(ISEED)
C--------------------------------
C RANDOMNUMBER GENERATOR
C (In RODUN08 this generator XRAN(ISEED) has in G77
c FOR NUMBERS BETWEEN 0 AND 1
C (been replaced by G77 intrinsic function RAND(0)
C-----------------------------
C ISEED=IRAND(ISEED)
C IRN =ISEED/100
C XRAN = FLOAT(IRN)/1000000
C RETURN
C END
C INTEGER FUNCTION IRAND(ISEED)
C INTEGER ISEED1, ISEED2, M1,M2, MULT1,MULT2, I
C MULT1= 3141
C MULT2 = 5821
C M1 = 100000000
C M2 = 10000
C ISEED1 = ISEED / M2
C ISEED2 = MOD( ISEED, M2 )
C I = MOD( ( ( MOD( (ISEED2*MULT1 + ISEED1*MULT2), M2)
C 1 * M2 ) + (ISEED2*MULT2) ), M1)
C IRAND = MOD(I+1,M1)
C RETURN
C END
C
SUBROUTINE XDPROD(CTOT,DPSEQ)
C------------------------------------------------------------------
C THIS SUBROUTINE SEQUENCES KK PROJECTS INTO AN OPTIMAL
C SEQUENCE USING DYNAMIC PROGRAMING.
C OPTIMAL COST: CTOT
C OPTIMAL SEQUENCE DPSEQ
C-------------------------------------------------------------------
COMMON /DM/ T(100),D(100),NDEDAT,TAR,DMULT,DKON
COMMON /DP/ CO(3003,2)
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
DIMENSION CI(3003,2)
C COST AND CAPACITY FOR I-STAGE
INTEGER DPSEQ(20)
INTEGER*2 STI(3003,20)
C PARTIAL SEQUENCE FOR I STAGE
INTEGER*2 STO(3003,20)
C --------------- FOR O-STAGE
INTEGER*2 NSTO(20)
C SCRATCH PARTIAL SEQUENCE
CHARACTER*1 CH
DO 10 I=1,3003
DO 20 J=1,KK
STI(I,J)=0
20 CONTINUE
10 CONTINUE
DO 30 I=1,KK
CI(I,1)=C(I)
C ESTABLISH THE FIRST COST,
CI(I,2)=X(I)
c AND CAPACITY VECTORS.
STI(I,1)=I
c ESTABLISH THE FIRST STAGE
30 CONTINUE
c VECTOR
DO 100 IST=1,KK-1
c GO THRU THE STAGES
DO 13 I=1,3003
DO 12 J=1,KK
STO(I,J)=0
CO(I,1)=0
CO(I,2)=0
12 CONTINUE
13 CONTINUE
CALL COMB(KK,IST,NCOM)
CALL COMB(KK,IST+1,NCO)
c CALCULATE # OF COMBINATIONS
DO 40 I=1,NCOM
C GO DOWN THE STAGE VECTOR
DO 50 IP=1,KK
c CHECK ALL PROJECTS
CH='N'
c LOWER 'ALREADY IN PSEQ-FLAG'
DO 60 IX=1,IST
IF (STI(I,IX).EQ.IP)CH='Y'
c THE FLAG IS RAISED
60 CONTINUE
IF(CH.EQ.'Y') GO TO 50
c PROJECT HAS ALREADY BEEN SELECTED
DO 70 IX=1,IST
c MOVE PARTIAL SEQUENCE
NSTO(IX)=STI(I,IX)
c TO NSTO AND ADD PROJECT # IP
70 CONTINUE
CALL INVDMD(CI(I,2),TCH)
CR=CI(I,1)+C(IP)/(1.+XINT/100)**TCH
XR=CI(I,2)+X(IP)
NSTO(IST+1)=IP
D WRITE(11,320) NSTO,I,IST,NCO,NCOM,IP,CR,XR
D320 FORMAT(' NSTO=',20I2,' I=',I2,' IST=',I2,' NCO=',I5,
D 1' NCOM=',I5,' IP=',I2/' CR=',F10.1,' XR=',F10.1)
CALL NEWST(STO,NSTO,IST,NCO,CR,XR)
C------- MATCH WITH NEXT STAGE VECTOR AND ADJUST THE VECTOR
50 CONTINUE
c NEW UNSELECTED PROJECT
40 CONTINUE
c STAGE VECTOR FINISHED
DO 80 I=1,NCO
IF(STO(I,1).EQ.0)GO TO 11
c CHECK # OF ELEMENTS IN STAGE
80 CONTINUE
c VECTOR STO
I=NCO
11 IF(I.NE.NCO) STOP'XDPROD 11'
DO 200 I=1,NCO
DO 201 J=1,IST+1
c MOVE STO TO STI AS WELL AS COST
STI(I,J)=STO(I,J)
c AND CAPACITY VECTORS
201 CONTINUE
CI(I,1)=CO(I,1)
CI(I,2)=CO(I,2)
200 CONTINUE
100 CONTINUE
c NEXT STAGE
DO 300 I=1,KK
DPSEQ(I)=STO(1,I)
300 CONTINUE
CTOT=CO(1,1)
c SET UP THE RESULTS
RETURN
END
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
SUBROUTINE NEWST(STO,NSTO,IST,NC,CR,XR)
C---------------------------------------------------------
C
C THIS SUBROUTINE TRIES TO MATCH THE PARTIAL SEQUENCE
C NSTO TO ANY PARTIAL SEQUENCE(ROW) IN THE OUT-STAGE VECTOR
C STO.
C IF A MATCH IS FOUND, A CHECK IS PERFORMED TO COMPARE THE
C TOTAL DISCOUNTED COST OF NSTO AGAINST THE APPROPRIATE ROW IN
C STO. SHOULD NSTO HAVE A LOWER COST, IT REPLACES THAT ROW
C IF NO MATCH IS FOUND, NSTO IS ADDED AS A NEW ROW TO THE
C STAGE VECTOR STO
C-----------------------------------------------------------------------
C CR= TOTAL DISCOUNTED COST OF NSTO
C XR= ACCUMULATED CAPACITY OF NSTO
C IST= STAGE # OF INPUT STAGE
C NC= # OF ROWS IN STO
C
COMMON/DM/ T(100),D(100),NDEDAT,TAR,DMULT,DKON
COMMON/DP/ CO(3003,2)
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
INTEGER*2 STO(3003,20),NSTO(20)
DO 10 I=1,NC
c GO DOWN STO
IF(STO(I,1).EQ.0)GO TO 100
DO 20 J=1,IST+1
NCH=0
DO 30 K=1,IST+1
c CHECK FOR A COMPLETE MATCH
IF(STO(I,J).EQ.NSTO(K)) NCH=1
30 CONTINUE
IF(NCH.EQ.0) GO TO 10
20 CONTINUE
IF(CR.GT.CO(I,1))GO TO 70
c COMPARE DISCOUNTED COST
CO(I,1)=CR
CO(I,2)=XR
DO 50 L=1,IST+1
STO(I,L)=NSTO(L)
50 CONTINUE
D WRITE(11,310) I,(STO(I,JA),JA=1,IST+1)
D310 FORMAT (' C STO(',I1,')= ',20I2)
RETURN
10 CONTINUE
8888 STOP' NEWST 8888'
100 CO(I,1)=CR
CO(I,2)=XR
DO 60 J=1,IST+1
STO(I,J)=NSTO(J)
60 CONTINUE
D WRITE(11,320) I,(STO(I,JA),JA=1,IST+1)
D320 FORMAT (' A STO(',I1,')= ',20I2)
RETURN
70 CONTINUE
D WRITE(11,330) I,(STO(I,JA),JA=1,IST+1)
D330 FORMAT (' B STO(',I1,')= ',20I2)
RETURN
END
SUBROUTINE SINDX(XA,L,SIND)
C--------------------------------------------------------------
C 0OOTHIS SUBROUTINE CALCULATES THE INDEX SIND, THAT
c RESULTS
C FROM ADDING PROJECT # L TO THE DEMAND FUNCTION
c WHEN
C ACCUMULATED CAPACITY OF PREVIOUSLY CONSTRUCTED
C PROJECTS
C EQUALS XA.THE INDEX SIND IS PROPORTIONAL TO THE RATIO
C OF
C CAPITAL COST TO PRESENT VALUE OF THE ACTUAL OUTPUT
C FROM THE PROJECT
C---------------------------------------------------------
COMMON /DM/ T(100),D(100),NDEDAT,TAR,DMULT,DKON
COMMON /RI/ DPCTOT,COUNT
c !DISCOUNTED TOTAL COST IN DP SEQ,COUNTER OF SPS
1,DPSEQ(20),PPP,PP
C !DYNPRO SEQUENCE,LOGICAL PRINT FLAGS
1,MSEQ(20)
c !SCRATCH SEQUENCE
1,IDOM
c !dominANCE INDICATOR(IDOM=1:DOMINANCE)
1,ERTO(20)
c !ACCUMULATED ERROR IN SEQ ALG
1,ERRMAX(20)
c !MAX ERROR IN SEQ ALG
1,NERR(20)
c !SP# WHERE MAX ERROR OCCURRED
1,KNT(20)
c !ERROR COUNTER FOR EACH SP
1,CT(20)
c !TOTAL DISCOUNTED COST FOR EACH ALG
1,ISEQ(20,20)
c !INDEX NO OF THE PROJECT IN J-TH PLACE
C !INSEQUENCE FOR ALG # I (I,J)
1,ER(20)
c !PERC ERROR IN EACH ALG
1,NOALG
c !NO OF HEURISTIC ALGORITHMS
1,IBETT(20,20)
c !COUNTER FOR SPS WHEN ALG #I IS SUPERIOR TO ALG#J (I,J)
1,IBEST(20)
c !COUNTER FOR SPS WHERE ALG #I (EXCEPT DP) IS BEST
COMMON /RAN2/ LL,MM
C 'RANDOM NO SEEDS FOR RAND '
C SEQUENCING
COMMON /RAN/ II,L1,M
C 'RAND NO SEEDS FOR RAND GEN and SEQ OF '
C PROJ
INTEGER DPSEQ,SEQ,RICNT,COUNT
COMMON C(20),X(20),CC(100),XX(100),KK,
1N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
REAL KICTOT
G=0
CALL INVDMD (XA,TA)
CALL INVDMD (XA+X(L),TOA)
C 'CALC. TIMING FOR PROJ L'
DO 10 I=1,NDEDAT
IF(XA.LT.D(I)) GO TO 11
C 'IN WHAT LINEAR SEGM OF D(T) IS PROJ L started up'
10 CONTINUE
GO TO 12
C ' STARTED LATER THAN POINT no NDEDAT'
11 IR = I
C ' FIRST POINT IN D(I) AFTER STARTUP'
DO 20 I=1,NDEDAT
IF (XA+X(L).LT.D(I)) GO TO 13
C 'IN WHAT LINEaR SEGM IS P no L EXHAUSTED'
20 CONTINUE
G=ENU(T(NDEDAT),TOA,D(NDEDAT),XA+X(L))
IS=NDEDAT
GO TO 14
13 IS=I-1
IF(IS.LT.IR) GO TO 12
G=ENU(T(IS),TOA,D(IS),XA+X(L))
14 G=G+ENU(TA,T(IR),XA,D(IR))
IF(IS.EQ.IR) GO TO 31
DO 30 I=IR,IS-1
G=G+ENU(T(I),T(I+1),D(I),D(I+1))
30 CONTINUE
31 GO TO 32
12 G=ENU(TA,TOA,XA,XA+X(L))
32 SIND=C(L)/G
RETURN
END
FUNCTION ENU(T1,T2,D1,D2)
C-------------------------------------------------------
C THIS FUNCTION CALCULATES THE DISCOUNTED 'AREA'
C BELOW A LINEAR S EGMENT OF THE DEMAND CURVE.
C-------------------------------------------------------
COMMON /DM/ T(100),D(100),NDEDAT,TAR,DMULT,DKON
COMMON C(20),X(20),CC(100),
1XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
IF (T2.EQ.T1) GO TO 10
ENU=(D2-D1)/(T2-T1)*((1+XINT/100)**(-T1)-(1+XINT/100)
1**(-T2))
10 IF(T2.EQ.T1)ENU = 0
RETURN
END
SUBROUTINE COMB(N,K,NCO)
C--------------------------------------------------------------------------
C THIS SUBROUTINE CALCULATES THR NUMBER OF COMBINATIONS
C EQUAL TO N!/(K!*(N-K)!)=NCO
C---------------------------------------------------------------------------
IF(K.EQ.N)GO TO 21
NK=N/2
K1=K
IF(K.GT.NK)K1=N-K
FNC=1
FNCO=1
DO 10 I=1,K1
FNC=FNC*I
FNCO=FNCO*(N-I+1)
10 CONTINUE
NCO=FNCO/FNC+0.01
GO TO 22
21 NCO=1
22 RETURN
END
C
SUBROUTINE PROGEN(ISEL,PRI)
C-----------------------------------------------------------------------
C THIS SUBROUTINE GENERATES FROM A RANDOM DISTRIBUTION,
C COST AND CAPACITY DATA FOR WATER RESOURCES PROJECTS.
C IF ISEL=0 NORMAL DISTRIB. & C=A*X**B
C----------------------------------------------------------------------
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
COMMON/RAN/II,L,M
LOGICAL PRI
REAL NORM
IF (ISEL) 2,3,4
3 XSD=1950
XMED=2050
A=89
B=.82
DO 10 I=1,NPROG
5 X(I)=NORM(XMED,XSD)
IF(X(I).LE.0)GO TO 5
CMED=A*X(I)**B
C(I)=NORM(CMED,CMED*.5)
IF(PRI) WRITE(10,1) I,C(I),I,X(I)
1 FORMAT(' PROGEN: C(',I2,')=',F10.2,' X(',I2,')=',F10.2)
10 CONTINUE
2 CONTINUE
4 RETURN
END
C
SUBROUTINE RANPAC(PRI)
LOGICAL PRI
C-----------------------------------------------------------------
C THIS SUBROUTINE RANDOMLY SELSCTS A SUBSET OF "NPRO"
C PROJECTS OUT OF A SUPERSET OF"NENDP-NBEGP" ACTUAL PROJECTS
C---------------------------------------------------------------------
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
COMMON/RAN/II,L,M
INTEGER*1 IE(100)
IF(PRI)WRITE(11,2) II-NOFIX
DO 20 I=NBEGP,NENDP
20 IE(I)=0
DO 10 I=1,NPRO
IGEN=0
40 J=RAND(0)*(NENDP-NBEGP+1)+NBEGP
IGEN=IGEN+1
IF(IGEN.GT.15700)STOP'PROJECT SELECTION NOT SUCCESSFUL'
IF(IE(J).EQ.1)GO TO 40
C(I)=CC(J)
X(I)=XX(J)
IE(J)=1
IF(PRI)WRITE(11,1)J
1 FORMAT(' PROJECT NO',I3,' SELECTED')
10 CONTINUE
2 FORMAT(' RANDOMLY SELECTED SET #',I3)
write(11,*) 'OUT OF ransel NEW NEW'
RETURN
END
C
SUBROUTINE KISEQ(CTOT,NSEQ,IND,P)
C------------------------------------------------------------------------------
C THIS PROGRAM SEQUENCES KK PROJECTS ACCORDING TO A INDEX
C OR BY RANDOM SELECTION IN THE SUBROUTINE RANSEL.
C THE INDEX IS CHOSEN AS FOLLOWS:
C IF IND=1 C(I)/X(I) IS THE INDEX
C IND=2 R(I) IS THE INDEX(TMR ALGORITHM)
C IND=3 K(I) IS THE INDEX
C IND=4 S(I) IS THE INDEX
C IND=5 RANDOM SELECTION IS PERFORMED
C TOTAL DISCOUNTED COST OF THE MINIMUM SEQUENCE IS RETURNED
C IN THE VARIABLE "CTOT"
C THE OPTIMAL SEQUENCE IS RETURNED IN THE ARRAY "NSEQ".
C IF THE LOGICAL VARIABLE P=.TRUE. DETAILED PRINTING IS
C PERFORMED.
C ------------------------------------------------------------------------------
DIMENSION TIMING(20),NSEQ(20)
COMMON /RAN2/ LL,MM
COMMON/DM/ T(100),D(100),NDEDAT,TAR,DMULT,DKON
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
REAL K(20,4)
c !INDEX FOR EACH PROJ AND EACH ALG.
LOGICAL P
IF((IND.EQ.1).AND.P) WRITE(10,24)' C(I)/X(I) SEQUENCING'
IF((IND.EQ.2).AND.P) WRITE(10,24)' R(I) SEQUENCING'
IF((IND.EQ.3).AND.P) WRITE(10,24)' K(I) SEQUENCING'
IF((IND.EQ.4).AND.P) WRITE(10,24)' S(I) SEQUENCING'
IF((IND.EQ.5).AND.P) WRITE(10,24)' RANDOM SEQUENCING '
24 FORMAT(15X,5('-'),A25,5('-'))
IF(P) WRITE(10,5)
5 FORMAT(1X,'PROJ.NO.,I |COST,C(I)|CAPACITY'
1 ,' | K(I) | KF(I) |'
1 ,' R(I) | S(I) |')
CTOT=0.
DO 30 JJ=1,KK
N(JJ)=0
30 CONTINUE
XACCUM=0.
DO 20 J=1,KK
XKMIN=1E10
IF(IND.EQ.5) CALL RANSEL(J,NMIN)
IF(IND.EQ.5) GO TO 15
DO 10 I=1,KK
IF(N(I).EQ.1) GO TO 10
K(I,1)=XINT*C(I)/(X(I)*100)
CALL INVDMD(XACCUM,T1)
CALL INVDMD(XACCUM+X(I),T2)
DELT=T2-T1
K(I,2)=C(I)/(1-(1+XINT/100)**(-DELT))
K(I,3)=K(I,2)*DELT/X(I)
IF(IND.EQ.4) CALL SINDX(XACCUM,I,SIND)
IF(IND.EQ.4) K(I,4)=SIND
IF(IND.NE.4) K(I,4)=1E12
IF(P) WRITE(10,6)I,C(I),X(I),K(I,3),K(I,1),K(I,2),K(I,4)
6 FORMAT (2X,I10,20F10.2)
IF(K(I,IND).GE.XKMIN)GO TO 10
XKMIN=K(I,IND)
NMIN=I
10 CONTINUE
15 N(NMIN)=1
CALL INVDMD(XACCUM,TIMING(J))
XACCUM=XACCUM +X(NMIN)
NSEQ(J)=NMIN
CTOT=CTOT+C(NMIN)/(1+XINT/100)**TIMING(J)
20 CONTINUE
IF(P) WRITE(10,21)(NSEQ(J),J=1,KK)
21 FORMAT(' PROJECT SEQUENCE:',20I6)
IF(P) WRITE(10,23)(TIMING(J),J=1,KK)
23 FORMAT(' TIMING OF PROJECTS',20F6.2)
IF(P) WRITE(10,22)CTOT
22 FORMAT(' DISCOUNTED TOTAL COST:',F15.2)
RETURN
END
SUBROUTINE RANSEL(JPRO,NSEL)
C----------------------------------------------------------------
C THIS SUBROUTINE RANDOMLY SELECTS A PROJECT TO THE JPRO-TH
C PLACE IN THE SEQUENCE OF KK PROJECTS. THE INDEX NO OF THE
C SELECTED PROJECT IS RETURNED AS NSEL.
C-----------------------------------------------------------------
COMMON /RAN2/ LL,MM
COMMON C(20),X(20),CC(100),XX(100),KK,N(20),NPRO,XINT,NPROG
1,NENDP,NBEGP,NOFIX
1 IRAN=INT(RAND(0)*(KK-JPRO+1))
IF(IRAN.EQ.(KK-JPRO+1)) GO TO 1
IRAN=IRAN+1
ICO=1
DO 10 I=1,KK
c ! MATCH IRAN TO UNSELECTED PROJ
IF(N(I).EQ.1)GO TO 10
c !THIS PROJ ALREADY SELECTED
IF(ICO .EQ. IRAN) GO TO 11
c ! MATCH FOUND
ICO=ICO+1
10 CONTINUE
STOP ' ERROR IN RANSEL'
11 NSEL=I
RETURN
END
C
REAL FUNCTION NORM(X,S)
COMMON/RAN/II,L,M
C-------------------------------------------------------------------------
C THIS ROUTINE GENERATES
C NORMALLY DISTRIBUTED RANDOM NUMBERS.
C------------------------------------------------------------------------
ASNI =RAND(0)
NORM=2.*S*ASNI + X-S
C--------------TEMPORARY SQUARE DISTRIBUTION---------------------
RETURN
END
SUBROUTINE INVDMD(XACC,TACC)
COMMON/DM/T(100),D(100),NDEDAT,TAR,DMULT,DKON
C
C THIS ROUTINE CALCULATES THE INVERSE DEMAND FUNCTION
C
DO 1 I=1,NDEDAT-1
IF((XACC.GE.D(I)).AND.(XACC.LT.D(I+1))) GO TO 20
1 CONTINUE
I=NDEDAT-1
20 TACC=T(I)+ (XACC-D(I))/(D(I+1)-D(I))*(T(I+1)-T(I))
D WRITE(10,2)D(I+1),D(I),XACC,TACC,I
D2 FORMAT(' INVDMD: ',4F12.2,I5)
RETURN
END
C
SUBROUTINE DEMAND(PRI)
C
C------------------------------------------------------------------
C THIS SUBROUTINE READS IN A NEW DEMAND FUNCTION
C FROM THE FILE DEM.DAT
C
LOGICAL PRI
COMMON/DM/ T(100),D(100),NDEDAT,TAR,DMULT,DKON
READ(9,1,END=20) TAR,NDEDAT,DMULT,DKON
IF(PRI)WRITE(10,4)TAR,NDEDAT,DMULT,DKON
DO 10 I = 1,NDEDAT
READ(9,2,END=30)T(I),D(I)
D(I)=D(I)*DMULT+DKON
IF (PRI) WRITE(11,3) I,T(I),D(I)
10 CONTINUE
GO TO 20
1 FORMAT(F5.0,I5,2F5.0)
2 FORMAT(2F10.0)
3 FORMAT(1X,'DEM. POINT # (I)=',I3,' T(I)=',F10.0,'D(I)=',F10.0)
4 FORMAT(20X,' NEW DEMAND FUNCTION '/ ' TIME HORIZON=',F10.0
1/' NO OF DEMAND POINTS',
2I3,' MULTIPLICATION FACTOR=',F10.2,' ADDING CONSTANT=',F10.2)
30 PAUSE ' ERROR IN FILE DEM.DAT'
20 CONTINUE
RETURN
END