program ggraph c --------------------------------------------------------------- c generator: max-flow min-cost directed grid graph c c Author: Mauricio G.C. Resende c AT&T Bell Laboratories c mgcr@research.att.com c c Date: 23 Sep 91 c 01 Feb 92 - modified to generate COST in [1,MAXCOST] c CAP in [1,MAXCAP] c c standard input: h, w, MAXCAP, MAXCOST, SEED c c graph has h * w + 2 nodes: source s, sink t and a grid of c h by w nodes. c c arcs go from s to h nodes of first grid column c from h nodes of last grid column to t c from (i,j) to (i+1,j) and (j+1,i) except c for last row that go from (i,h) to (i+1,h) and last c column that go from (w,j) to t c c arc capacities for grid arcs uniformly distributed [0,MAXCAP] c for s - (1,j) = INFINITY c for (w,j) - t = INFINITY c c arc costs for grid arcs uniformly distributed [0,MAXCOST] c for s - (1,j) = 0 c for (w,j) - t = 0 c c Supply at node s = MAXFLO(s-t) c Demand at node t = MAXFLO(s-t) c Supply and demands for all other nodes = 0 c c Output: DIMACS min format c c Uses subroutines: dnfwd and dnsub to compute max-flow. c Reference: D. Goldbarb, M.D. Grigoriadis, c "A computational comparison of the dinic c and network simplex methods for maximum flow", c Annals of Operations Research 7, 1988. c c rand - portable pseudo-random number generator. c Reference: L. Schrage, "A More Portable FORTRAN c Random Number Generator", ACM Transactions on c Mathematical Software, June, 1979. c c --------------------------------------------------------------- integer maxh,maxw,maxarc parameter (maxh=1000,maxw=1000) parameter (maxarc=524273) integer caps(maxh),costs(maxh),capt(maxh),costt(maxh) integer h,w,maxf,maxc,n1,n2,cap,unif,i,j integer nnodes,narcs,snode,tnode,seed,seed0 integer cost(524273) integer v10(maxarc),v20(maxarc),cap0(maxarc),cost0(maxarc) C C************START DNSUB STANDARD USER XFACE DECLARATION BLOCK********* C C USER MUST INCLUDE THIS BLOCK INHIS/HER CALLING PROGRAM. C C THESE DECLARATIONS ACCOUNT FOR THE ENTIRE ARRAY STORAGE USED BY C THE DNSUB SUBROUTINES, INCLUDING WHAT IS NECESSARY TO STORETHE C INPUT DATA. IT AMOUNTS TO 6*NODES +5*ARCS WORDS, BUT IT MAY BE C STATED MORE PRECISELY AS FOLLOWS: C C 1) FOR SOLVING PROBLEMS WITH UP TO 32,765 NODES AND 2,147,483,647 C ARCS, AND WITH FLOW VALUES OF UP TO2,147,483,647, IT SUFFICES C TO DECLARE ARRAYS AS SHOWN UNDER 'STD ARRANGEMENT' IN THE C TABLE BELOW. THIS ARRANGEMENT IS THE ONE USED IN THE RELEASE C VERSION OF THE DNSUB SUBROUTINES. C C 2) IN ORDER TO HANDLE PROBLEMS WITHMORE THAN 32,765 NODES, ALL C INTEGER*2 DECLARATIONS BELOW MUST BE CHANGED TO INTEGER*4 AND C THE DNSUB SUBROUTINES AND USER CALLING PROGRAMS MUST BE COMPILED C AND LINKED.(SEE ALTERNATE 1 BELOW.) C C 3) FOR SOLVING ONLY PROBLEMS WITH AT MOST 32,765 NODES, 32,767 C ARCS, AND WITH FLOWS NOT EXCEEDING 32,767, ALL INTEGER*4 C DECLARATIONS BELOW MAY BE CHANGED TO INTEGER*2, AND THE DNSUB C SUBROUTINES AND USER CALLING PROGRAMS MUSTBE COMPILED AND C LINKED. (SEE ALTERNATE 2 BELOW.) IT IS ALSO NECESSARY TO C ACTIVATE THE STATEMENT DNIBIG=32767 in subroutine dnsub, and to C include 'implicit INTEGER*2 (i-n)' statements in each program C unit. THIS ARRANGEMENT, WHICH USESTHE LEAST AMOUNT OF MEMORY C AND IS THE FASTEST when executed on16-BIT PROCESSORS, IS ONLY C USEFUL FOR A SET OF VERY SPECIAL APPLICATIONS. ITSgeneral C use is not recommended. OTHER ALTERNATES ARE ALSO POSSIBLE. C C ----------------------------------------------------------------- C C STD ARRANGEMENT ALTERNATE 1 ALTERNATE 2 C (not recommended) C --------------- --------------- ---------------- C ARRAY LENGTH TYPE EQV'D TYPE TYPE EQV'D TYPE TYPE EQV'D TYPE C ------ ----- ---- ------- ---- ----- ---- ---- ----- ----- C DNFAPT N 4 4 2 C DNFADJ A 2 4 2 C DNCAP A 4 4 2 C DNBAPT N 4 4 2 C DNBADJ A 2 DNFROM 2 4 DNFROM 4 2 DNFROM 2 C DNFLOW A 4 DNTO 2 4 DNTO 4 2 DNTO 2 C DNBTOF A 4 DNGCAP 4 4 DNGCAP 4 2 DNGCAP 2 C DNPTRF N 4 4 2 C DNPTRB N 4 4 2 C DNLIST N 2 4 2 C DNFLAB N 4 DNDIST 2 4 DNDIST 4 2 DNDIST 2 C ----------------------------------------------------------------- C C TOTAL BYTES: 22*N + 16*A 24*N + 20*A 12*N + 10*A C (N=NODES, A=ARCS) C ----------------------------------------------------------------- C C DIMENSION THESE TO AT LEAST MAX NUMBER OF NODES + 2: INTEGER DNLIST(65538),DNDIST(65538) INTEGER DNFAPT(65538),DNPTRF(65538),DNBAPT(65538), * DNFLAB(65538),DNPTRB(65538) C DIMENSION THESE TO AT LEAST MAX NUMBER OF ARCS + 1: INTEGER DNFADJ(524273),DNFROM(524273),DNTO(524273),DNBADJ(524273) INTEGER DNCAP(524273),DNGCAP(524273),DNFLOW(524273),DNBTOF(524273) C COMMON /DN01/DNFAPT/DN02/DNFADJ/DN03/DNCAP/DN04/DNFROM COMMON /DN05/DNFLOW/DN06/DNGCAP/DN08/DNPTRF/DN09/DNLIST COMMON /DN10/DNBAPT/DN11/DNFLAB/DN12/DNPTRB C INTEGER DNNODE,NODP12,DNSRCE,DNSINK,I2,II2,K2,J2,IHEAD2, * ITAIL2,QHEAD2,QTAIL2,KP2,KSAT2,DNNOP2,LAYMAX,K2DST,K2DP1,MAXD2 INTEGER DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE,DNIBIG COMMON /DN00/DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE, * DNELTM,DNIBIG,DNNODE,DNSRCE,DNSINK,DNNOP2 EQUIVALENCE (DNFLOW(1),DNTO(1)),(DNFLAB(1),DNDIST(1)), * (DNFROM(1),DNBADJ(1)),(DNGCAP(1),DNBTOF(1)) C C************END DNSUB STANDARD USER XFACE DECLARATION BLOCK*********** read (5,*) h,w,maxf,maxc,seed seed0=seed nnodes=h*w+2 narcs=2*(h-1)*w+w+h snode=h*w+1 tnode=h*w+2 c for dinic interface c ------------------------------------------------------ nodes = nnodes isou=snode isnk=tnode c ------------------------------------------------------ numarc=0 do 1010 i=1,h-1 numarc=numarc+1 dnfrom(numarc)=(i-1)*w+1 dnto(numarc)=dnfrom(numarc)+1 dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) caps(i)=dngcap(numarc) numarc=numarc+1 dnfrom(numarc)=(i-1)*w+1 dnto(numarc)=dnfrom(numarc)+w dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) caps(i)=caps(i)+dngcap(numarc) costs(i)=unif(maxc,seed) 1010 continue numarc=numarc+1 dnfrom(numarc)=(h-1)*w+1 dnto(numarc)=dnfrom(numarc)+1 dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) caps(h)=dngcap(numarc) costs(h)=unif(maxc,seed) do 1020 i=1,h-1 do 1030 j=2,w-2 numarc=numarc+1 dnfrom(numarc)=(i-1)*w+j dnto(numarc)=dnfrom(numarc)+1 dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) numarc=numarc+1 dnfrom(numarc)=(i-1)*w+j dnto(numarc)=dnfrom(numarc)+w dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) 1030 continue 1020 continue do 1040 j=2,w-2 numarc=numarc+1 dnfrom(numarc)=(h-1)*w+j dnto(numarc)=dnfrom(numarc)+1 dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) 1040 continue do 1050 i=1,h-1 numarc=numarc+1 dnfrom(numarc)=(i-1)*w+w-1 dnto(numarc)=dnfrom(numarc)+1 dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) capt(i)=dngcap(numarc) numarc=numarc+1 dnfrom(numarc)=(i-1)*w+w-1 dnto(numarc)=dnfrom(numarc)+w dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) costt(i)=unif(maxc,seed) 1050 continue numarc=numarc+1 dnfrom(numarc)=(h-1)*w+w-1 dnto(numarc)=dnfrom(numarc)+1 dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) capt(h)=dngcap(numarc) costt(h)=unif(maxc,seed) do 1060 i=1,h-1 numarc=numarc+1 dnfrom(numarc)=(i-1)*w+w dnto(numarc)=dnfrom(numarc)+w dngcap(numarc)=unif(maxf,seed) cost(numarc)=unif(maxc,seed) capt(i+1)=capt(i+1)+dngcap(numarc) 1060 continue do 1070 i=1,h numarc=numarc+1 dnfrom(numarc)=h*w+1 dnto(numarc)=(i-1)*w+1 dngcap(numarc)=caps(i) cost(numarc)=costs(i) 1070 continue do 1080 i=1,h numarc=numarc+1 dnfrom(numarc)=i*w dnto(numarc)=h*w+2 dngcap(numarc)=capt(i) cost(numarc)=costt(i) 1080 continue do 90101 i=1,numarc v10(i)=dnfrom(i) v20(i)=dnto(i) cap0(i)=dngcap(i) cost0(i)=cost(i) 90101 continue call dnfwd(nodes,narcs,iretn) call dnsub(nodes,narcs,isou,isnk,maxflo,numaug,nphase, + nncut,nacut,etime0,iretn) write(6,100) h,w 100 format("c max-flow min-cost on a",i3," by ",i3," grid") write(6,105) maxf,maxc,seed0 105 format("c cap: [0,",i8,"] cost: [0,",i8,"] seed:", i10) write(6,101) nodes,narcs 101 format("p min ",2i15) write(6,103) isou,maxflo 103 format("n ",2i15) write(6,103) isnk,-maxflo do 9010 i=1,numarc write(6,102) v10(i),v20(i),cap0(i),cost0(i) 102 format("a ",2i10," 0 ",2i10) 9010 continue stop end integer function unif(max,seed) integer max,seed real rand unif=1+rand(seed)*(max-1) return end real function rand(ix) c ================================================================= c Portable pseudo-random number generator. c Reference: l. schrage, "A More Portable FORTRAN c Random Number Generator", ACM Transactions on c Mathematical Software, June, 1979. c ================================================================= integer*4 a,p,ix,b15,b16,xhi,xalo,leftlo,fhi,k data a/16807/,b15/32768/,b16/65536/,p/2147483647/ xhi=ix/b16 xalo=(ix-xhi*b16)*a leftlo=xalo/b16 fhi=xhi*a+leftlo k=fhi/b15 ix=(((xalo-leftlo*b16)-p)+(fhi-k*b15)*b16)+k if (ix.lt.0) ix=ix+p rand=float(ix)*4.656612875e-10 return end SUBROUTINE DNSUB(NODES,NARCS,ISOURC,ISINK,MAXFLO, * NUMAUG,NUMSTG,NNCUT,NACUT,ELTIM,IRETN) C C************START DNSUB STANDARD USER XFACE DECLARATION BLOCK********* C C USER MUST INCLUDE THIS BLOCK INHIS/HER CALLING PROGRAM. C C THESE DECLARATIONS ACCOUNT FOR THE ENTIRE ARRAY STORAGE USED BY C THE DNSUB SUBROUTINES, INCLUDING WHAT IS NECESSARY TO STORETHE C INPUT DATA. IT AMOUNTS TO 6*NODES +5*ARCS WORDS, BUT IT MAY BE C STATED MORE PRECISELY AS FOLLOWS: C C 1) FOR SOLVING PROBLEMS WITH UP TO 32,765 NODES AND 2,147,483,647 C ARCS, AND WITH FLOW VALUES OF UP TO2,147,483,647, IT SUFFICES C TO DECLARE ARRAYS AS SHOWN UNDER 'STD ARRANGEMENT' IN THE C TABLE BELOW. THIS ARRANGEMENT IS THE ONE USED IN THE RELEASE C VERSION OF THE DNSUB SUBROUTINES. C C 2) IN ORDER TO HANDLE PROBLEMS WITHMORE THAN 32,765 NODES, ALL C INTEGER*2 DECLARATIONS BELOW MUST BE CHANGED TO INTEGER*4 AND C THE DNSUB SUBROUTINES AND USER CALLING PROGRAMS MUST BE COMPILED C AND LINKED.(SEE ALTERNATE 1 BELOW.) C C 3) FOR SOLVING ONLY PROBLEMS WITH AT MOST 32,765 NODES, 32,767 C ARCS, AND WITH FLOWS NOT EXCEEDING 32,767, ALL INTEGER*4 C DECLARATIONS BELOW MAY BE CHANGED TO INTEGER*2, AND THE DNSUB C SUBROUTINES AND USER CALLING PROGRAMS MUSTBE COMPILED AND C LINKED. (SEE ALTERNATE 2 BELOW.) IT IS ALSO NECESSARY TO C ACTIVATE THE STATEMENT DNIBIG=32767 in subroutine dnsub, and to C include 'implicit INTEGER*2 (i-n)' statements in each program C unit. THIS ARRANGEMENT, WHICH USESTHE LEAST AMOUNT OF MEMORY C AND IS THE FASTEST when executed on16-BIT PROCESSORS, IS ONLY C USEFUL FOR A SET OF VERY SPECIAL APPLICATIONS. ITSgeneral C use is not recommended. OTHER ALTERNATES ARE ALSO POSSIBLE. C C ----------------------------------------------------------------- C C STD ARRANGEMENT ALTERNATE 1 ALTERNATE 2 C (not recommended) C --------------- --------------- ---------------- C ARRAY LENGTH TYPE EQV'D TYPE TYPE EQV'D TYPE TYPE EQV'D TYPE C ------ ----- ---- ------- ---- ----- ---- ---- ----- ----- C DNFAPT N 4 4 2 C DNFADJ A 2 4 2 C DNCAP A 4 4 2 C DNBAPT N 4 4 2 C DNBADJ A 2 DNFROM 2 4 DNFROM 4 2 DNFROM 2 C DNFLOW A 4 DNTO 2 4 DNTO 4 2 DNTO 2 C DNBTOF A 4 DNGCAP 4 4 DNGCAP 4 2 DNGCAP 2 C DNPTRF N 4 4 2 C DNPTRB N 4 4 2 C DNLIST N 2 4 2 C DNFLAB N 4 DNDIST 2 4 DNDIST 4 2 DNDIST 2 C ----------------------------------------------------------------- C C TOTAL BYTES: 22*N + 16*A 24*N + 20*A 12*N + 10*A C (N=NODES, A=ARCS) C ----------------------------------------------------------------- C C DIMENSION THESE TO AT LEAST MAX NUMBER OF NODES + 2: INTEGER DNLIST(65538),DNDIST(65538) INTEGER DNFAPT(65538),DNPTRF(65538),DNBAPT(65538), * DNFLAB(65538),DNPTRB(65538) C DIMENSION THESE TO AT LEAST MAX NUMBER OF ARCS + 1: INTEGER DNFADJ(524273),DNFROM(524273),DNTO(524273),DNBADJ(524273) INTEGER DNCAP(524273),DNGCAP(524273),DNFLOW(524273),DNBTOF(524273) C COMMON /DN01/DNFAPT/DN02/DNFADJ/DN03/DNCAP/DN04/DNFROM COMMON /DN05/DNFLOW/DN06/DNGCAP/DN08/DNPTRF/DN09/DNLIST COMMON /DN10/DNBAPT/DN11/DNFLAB/DN12/DNPTRB C INTEGER DNNODE,NODP12,DNSRCE,DNSINK,I2,II2,K2,J2,IHEAD2, * ITAIL2,QHEAD2,QTAIL2,KP2,KSAT2,DNNOP2,LAYMAX,K2DST,K2DP1,MAXD2 INTEGER DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE,DNIBIG COMMON /DN00/DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE, * DNELTM,DNIBIG,DNNODE,DNSRCE,DNSINK,DNNOP2 EQUIVALENCE (DNFLOW(1),DNTO(1)),(DNFLAB(1),DNDIST(1)), * (DNFROM(1),DNBADJ(1)),(DNGCAP(1),DNBTOF(1)) C C************END DNSUB STANDARD USER XFACE DECLARATION BLOCK*********** C*********************************************************************** C C CALLING CONDITIONS: C USER CALLABLE WITH INPUT DATA AS FOLLOWS. C INPUT: C SCALARS (IN CALLING SEQUENCE; ALL INTEGER*4): C NODES: NUMBER OF NODES (INCLUDING SOURCE AND SINK) C NARCS: NUMBER OF ARCS C ISOURC: NODE NUMBER FOR SOURCE C ISINK: NODE NUMBER FOR SINK C C ARRAYS (IN COMMON): C DNFAPT: (NODES+1)-LONG INTEGER*4 POINTER ARRAY FOR FORWARD ADJA- C CENCY LISTS (I.E. THE FORWARD ADJACENCY LIST OF A NODE I C IS THE SET OF ARCS C {(I,DNFADJ(J)) :J=DNFAPT(I),...,DNFAPT(I+1)-1 } C NOTE: MUST HAVE DNFAPT(NODES+1) = NARCS+1. C C DNFADJ: NARCS-LONG INTEGER*2 ARRAY GIVING THE LIST OF NODES IN C THE FORWARD ADJACENCY LIST DNFADJ(J) FOR EACH NODE J, C AS DESCRIBED ABOVE. C C DNCAP: NARCS-LONG INTEGER*4 ARRAY GIVING THE ARC CAPACITIES, IN C THE ORDER PRESCRIBED BY DNFAPT(.) AND DNFADJ(.). ALL ARC C CAPACITIES MUST BE GIVEN AS POSITIVE INTEGERS. C C NOTE: FOR TRANSFORMING UNORDERED ARC LISTS TO THE ABOVE INPUT C DATA STRUCTURE,USE SUBROUTINE DNFWD BEFORE CALLING DNSUB C C OUTPUT: C C SCALARS (IN CALLING SEQUENCE): C MAXFLO: VALUE OF MAXIMUM FLOW (INTEGER*4). C NUMAUG: NUMBER OF FLOW AUGMENTATIONS (INTEGER*4). C NUMSTG: NUMBER OF STAGES (LAYERED NETWORKS CREATED; INTEGER*4). C NNCUT: NUMBER OF NODES ON source SIDE OF FINAL CUT (INTEGER*4) C NACUT: NUMBER OF SATURATED ARCS IN THE FINAL CUT (INTEGER*4) C ELTIM: EXECUTION TIME, IN SECONDS (SEE NOTE 3 BELOW) (REAL*4). C IRETN: NONZERO IF THERE ARE ERRORS IN INPUT DATA. C C ARRAYS (IN COMMON:): C DNFAPT: THE ORIGINAL DNFAPT(.) WITH SOME OF ITS ELEMENTS C NEGATED TO MARK THE NODES THAT ARE ON THE SINK SIDE OF C THE MINIMUM CUT, I.E. DNFAPT(I)<0 IMPLIES THAT NODE I IS C ON THE SINK SIDE OF CUT, ELSE, IT IS ON THE source SIDE. C C DNFADJ: THE ORIGINAL DNFADJ(.) WITH SOME OF ITS ELEMENTS NEGATED C TO MARK THOSE SATURATED ARCS THAT ARE IN THE MIN CUT C FOUND BY THE ALGORITHM, I.E. DNFADJ(I)<0 IMPLIES THAT C THE I-TH ARC IN FORWARD ADJACENY ORDER IS IN THE CUT; C THE FOLLOWING CODE SEGMENT PRINTS THESE ARCS: C C DO 2 K=1,NODES C IBEG=IABS(DNFAPT(K)) C IEND=IABS(DNFAPT(K+1))-1 C DO 1 I=IBEG,IEND C IF(DNFADJ(I).LT.0)PRINT THE ARC C 1 CONTINUE C 2 CONTINUE C C (SEE NOTE 1 BELOW). C C DNCAP: THE ORIGINAL CAPACITIES, UNALTERED. C C DNFLOW: NARCS-LONG INTEGER*4 ARRAY WHICH GIVES FLOWS ON ARCS, IN C FORWARD ADJACENCY ORDER. C C DNBAPT: (NODES+1)-LONG INTEGER*4 POINTER ARRAY FOR BACKWARD C ADJACENCY LISTS (I.E. THE FORWARD ADJACENCY LIST OF C NODE I IS THE SET OF ARCS C {(DNBADJ(J),J): J=DNBAPT(I),...,DNBAPT(I+1)-1 }. C NOTE: WE MUST HAVE DNBAPT(NODES+1) = NARCS+1. C C DNBADJ: NARCS-LONG INTEGER*2 ARRAY GIVING THE BACKWARD ADJACENCY C LIST OF EACH NODE J, ONE AFTER THE OTHER, AS DESCRIBED C ABOVE (ALSO SEE NOTE 2 BELOW). C C DNBTOF: NARCS-LONG INTEGER*4 ARRAY GIVING THE BACKWARD ADJACENCY C TO FORWARD ADJACENCY MAPPING, I.E. THE J-TH ARC C IN THE BACKWARD ADJACENCY ORDER IS THE DNBTOF(J)-TH C ARC IN THE FORWARD ADJACENCY ORDER. C C*********************************************************************** C C INTERNAL SCRATCH ARRAYS (IN COMMON): C C DNLIST: (NODES+1)-LONG INTEGER*2 ARRAY USED AS FOLLOWS: C 1) AS THE QUEUE IN BFS SEARCH FOR CONSTRUCTING THE C LAYERED GRAPH (SUBR. DNBFS). C 2) AS THE "PARENT" ARRAY IN DFS SEARCH OF LAYERED C GRAPH (SUBR: DNDFS, DNPUSH). C 3) AS SCRATCH IN CONSTRUCTING FWD ADJACENCIES C (SUBR DNFWD), AND BWD ADJACENCIES (SUBR. DNSUB). C C DNFLAB: NODES-LONG INTEGER*4 ARRAY USED IN DFS SEARCH OF C LAYERED GRAPH TO STORE THE FLOW LABEL FOR EACH NODE C (SEE ARRAY DNDIST() ). C C DNDIST: NODES-LONG INTEGER*2 ARRAY USED IN CONSTRUCTING THE C LAYERED GRAPH BY BFS, STORING THE DISTANCE OF EACH C NODE FROM THE SOURCE; DNDIST(ISOURC)=0 (EQUIVALENCED C TO THE INTEGER*4 ARRAY DNFLAB() ). C C DNPTRF: NODES-LONG INTEGER*4 STATUS AND POINTER ARRAY, I.E. C DNPTRF(K)=0 IF THERE ARE NO OPEN ARCS LEAVING NODE K C IN THE REPRESENTATION OF THE LAYERED C NETWORK. C >0 THE ARC INDEX, IN 'FA' ORDER, THAT C CORRESPONDS TO THE LAST ARC IN THE C FWD ADJACENCY OF K WHICH IS SCANNED C BACKWARD. C C DNPTRB: NODES-LONG INTEGER*4 STATUS AND POINTER ARRAY, I.E. C DNPTRB(K)=0 IF THERE ARE NO OPEN ARCS ENTERING NODE K C IN THE REPRESENTATION OF THE LAYERED C NETWORK. C >0 THE ARC INDEX (IN 'BA' ORDER) OF THE FIRST C ARC IN THE BWD ADJACENCY OF K TO BE C SCANNED IN LIFO ORDER. C C*********************************************************************** C C NOTES:1) DNFADJ(J) IS INTERNALLY NEGATED TO MARK THE J-TH ARC (IN C FORWARD ADJACENCY ORDER) AS 'CLOSED'. BEFORE RETURNING C TO THE CALLING PGM, THESE MARKINGS ARE DISCARDED, AND SOME C ELEMENTS OF DNFADJ(.) ARE NEGATED TO REFLECT THE OUTPUT C SPECIFICATION DESCRIBED ABOVE. C C 2) DNBADJ(J) IS INTERNALLY NEGATED TO MARK THE J-TH ARC C (IN BACKWARD ADJACENCY ORDER) AS 'CLOSED'. C C 3) THE USER MUST PROVIDE A TIMING SUBROUTINE THAT RETURNS C THE CURRENT TIME IN VARIABLE ' T ', IN SECONDS. SUCH C A SUBROUTINE IS INSTALLATION-DEPENDENT AND IS NOT GIVEN HERE. C IF THIS CANNOT BE DONE, PROVIDE THE FOLLOWING SUBROUTINE AND C FORGET THE EXECUTION TIMING INFORMATION: C C SUBROUTINE GETIME(T) C REAL T C T=0.0 C RETURN C END C C*********************************************************************** C C SUBROUTINES CALLED: DNBFS, DNDFS, DNCUT, GETIME,DNCLEA C C*********************************************************************** C C NAMES COLLECTIVELY RESERVED BY ALL DNSUB SUBROUTINES: C C DN00, DN01, DN02, DN03, DN04, DN05, DN06, DN08, DN09, DN10, DN11, C DN12, DNARC, DNAUG, DNBADJ, DNBAPT, DNBFS, DNBTOF, DNCAP, DNCLEA, C DNCUT, DNDFS, DNDIST, DNFADJ, DNPUSH, DNFAPT, DNFLAB, DNFLOW, C DNFROM, DNFVA, DNFWD, DNGCAP, DNIBIG, DNLAUG, DNLFVA, DNLIST, C DNNODE, DNNOP2, DNOUT, DNPTRB, DNPTRF, DNSINK, DNSRCE, DNSTGE, C DNSUB, DNELTM, DNTO. C C********************************************************************** C INITIALIZATION: C C CALL GETIME(TIMBEG) TIMBEG=etime(it) IRETN=0 DNIBIG=2147483647 C DNIBIG=32767 DNARC=NARCS DNNODE=NODES DNNOP2=DNNODE+2 DNSRCE=ISOURC DNSINK=ISINK DO 200 I4=1,DNARC 200 DNFLOW(I4)=0 C*********************************************************************** C C USING THE FORWARD ADJACENCY ARRAYS DNFAPT(.) AND DNFADJ(.) (IN C COMMON), CREATE THE BACKWARD ADJACENCY ARRAYS DNBAPT(.), DNBADJ(.) C AND THE BACKWARD TO FORWARD ADJACENCY MAPPING DNBTOF(.). DNLIST(.) IS C USED HERE AS A SCRATCH ARRAY. C NODP12=DNNODE+1 DO 210 I2=1,NODP12 DNLIST(I2)=0 210 DNBAPT(I2)=0 C TEMPORARILY STORE IN DNBAPT(.) NUMBER OF ARCS INTO EACH NODE: DO 230 I2=1,DNNODE IBEG4=DNFAPT(I2) IEND4=DNFAPT(I2+1)-1 DO 220 I4=IBEG4,IEND4 JJDN=DNFADJ(I4) 220 DNBAPT(JJDN)=DNBAPT(JJDN)+1 230 CONTINUE C CONSTRUCT DNBAPT(.): IHPI4=DNBAPT(1) DNBAPT(1)=1 DO 240 I2=1,DNNODE IHPSV4=IHPI4+DNBAPT(I2) IHPI4=DNBAPT(I2+1) 240 DNBAPT(I2+1)=IHPSV4 C CONSTRUCT DNBADJ(.) AND DNBTOF(.): DO 260 I2=1,DNNODE IBEG4=DNFAPT(I2) IEND4=DNFAPT(I2+1)-1 DO 250 I4=IBEG4,IEND4 IHEAD2=DNFADJ(I4) IHPUT4=DNBAPT(IHEAD2)+DNLIST(IHEAD2) IF(IHPUT4.LE.0)GOTO1900 DNBADJ(IHPUT4)=I2 DNBTOF(IHPUT4)=I4 250 DNLIST(IHEAD2)=DNLIST(IHEAD2)+1 260 CONTINUE C C*********************************************************************** C DNFVA=0 DNAUG=0 DNSTGE=1 C C STAGE LOOP: 300 CONTINUE C FORM NEW LAYERED NETWORK W.R.T. CURRENT FEASIBLE FLOW. CALL DNBFS(LAYMAX) IF(LAYMAX.GT.0)GOTO1500 C COMPUTE A MAXIMAL FLOW IN THIS LAYERED NET: CALL DNDFS C UPDATE FLOW VALUE: DNFVA=DNFVA+DNLFVA DNAUG=DNAUG+DNLAUG DNSTGE=DNSTGE+1 C ERASE MARKINGS THAT DEFINE LAYERED NETWORK: CALL DNCLEA GOTO300 C*********************************************************************** C C CURRENT FLOW IS MAXIMUM: 1500 MAXFLO=DNFVA NUMAUG=DNAUG NUMSTG=DNSTGE CALL DNCLEA C MARK EXECUTION TIME: C CALL GETIME(TIMEND) TIMEND=etime(it) DNELTM=TIMEND-TIMBEG ELTIM=DNELTM C C MARK THE MIN CUT-SET IN DNFAPT(.) AND DNFADJ(.): CALL DNCUT(NNCUT,NACUT) C RETURN C 1900 IRETN=1 RETURN END SUBROUTINE DNBFS(LAYMAX) C C C DIMENSION THESE TO AT LEAST MAX NUMBER OF NODES + 2: INTEGER DNLIST(65538),DNDIST(65538) INTEGER DNFAPT(65538),DNPTRF(65538),DNBAPT(65538), * DNFLAB(65538),DNPTRB(65538) C DIMENSION THESE TO AT LEAST MAX NUMBER OF ARCS + 1: INTEGER DNFADJ(524273),DNFROM(524273),DNTO(524273),DNBADJ(524273) INTEGER DNCAP(524273),DNGCAP(524273),DNFLOW(524273),DNBTOF(524273) C COMMON /DN01/DNFAPT/DN02/DNFADJ/DN03/DNCAP/DN04/DNFROM COMMON /DN05/DNFLOW/DN06/DNGCAP/DN08/DNPTRF/DN09/DNLIST COMMON /DN10/DNBAPT/DN11/DNFLAB/DN12/DNPTRB C INTEGER DNNODE,NODP12,DNSRCE,DNSINK,I2,II2,K2,J2,IHEAD2, * ITAIL2,QHEAD2,QTAIL2,KP2,KSAT2,DNNOP2,LAYMAX,K2DST,K2DP1,MAXD2 INTEGER DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE,DNIBIG COMMON /DN00/DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE, * DNELTM,DNIBIG,DNNODE,DNSRCE,DNSINK,DNNOP2 EQUIVALENCE (DNFLOW(1),DNTO(1)),(DNFLAB(1),DNDIST(1)), * (DNFROM(1),DNBADJ(1)),(DNGCAP(1),DNBTOF(1)) C C******************************************************************* C C CALLING CONDITIONS: C INTERNAL CALL FROM DNSUB ONLY. C INPUT ARRAYS: C DNFAPT, DNBAPT, DNFADJ, DNBADJ, DNBTOF, DNCAP, DNFLOW C SCRATCH ARRAYS: C DNLIST C OUTPUT: C IN CALLING SEQUENCE (INTEGER*2): C LAYMAX: 0 (LAYERED NET COMPLETE WITH SINK IN LAST LAYER.) C 1 (LAST LAYER IS EMPTY; CURRENT FLOW IN ORIGINAL C NETWORK IS MAXIMUM.) C IN COMMON (ARRAYS): C DNPTRF, DNPTRB, DNDIST C*********************************************************************** C C THIS SUBROUTINE CONSTRUCTS THE LAYERED NETWORK OF THE RESIDUAL C NETWORK W.R.T. THE CURRENT FLOW. THE RESIDUAL NETWORK IS C INFERRED, AND THE LAYERED NETWORK IS ONLY RECORDED BY TEMPORARY C MARKERS ON THE ORIGINAL GRAPH. THE SUBROUTINE IMPLEMENTS A C BREADTH-FIRST SEARCH ON THE RESIDUAL NETWORK, I.E. NODES ARE C PROCESSED AS IN A QUEUE. SUCH A SEARCH PARTITIONS THE NODES INTO C A SET OF 'LAYERS', ACCORDING TO THEIR CARDNALITY DISTANCE FROM C THE source. THE DNDIST(.) ARRAY IS USED TO RECORD THESE DISTANCES, C AND THUS DEFINES THE LAYERS (NOTE: DNDIST(ISOURC)=0). C C IN SELECTING ARCS OF THE RESIDUAL NETWORK THAT ARE DIRECTED FROM C ONE LAYER TO THE NEXT, WE SHALL MARK THE ARCS IN THE ORIGINAL C NETWORK AS 'OPEN' OR 'CLOSED'. A 'CLOSED' ARC (V,W) IN THE C LAYERED NETWORK IS EITHER AN ARC I=(V,W) THAT IS SATURATED IN THE C ORIGINAL NETWORK, OR AN ARC I=(W,V) OF THE ORIGINAL NETWORK WITH C ZERO FLOW. IN THE FORMER CASE, THE ABSENCE OF THIS ARC FROM THE C LAYERED NETWORK IS MARKED BY NEGATING DNFADJ(I), AND IN THE LATTER C CASE, BY NEGATING DNBADJ(I). AN 'OPEN' ARC I=(V,W) FOR THE LAYERED C NETWORK EITHER HAS POSITIVE RESIDUAL CAPACITY, OR ARC (W,V) HAS C POSITIVE FLOW, IN THE ORIGINAL NETWORK. IF A NODE V HAS NO OPEN C ARCS IN THE LAYERED NETWORK THAT ARE IN THE FORWARD ADJACENCY LIST C OF V IN THE ORIGINAL NETWORK, THEN WE SET DNPTRF(V)=0; AND IF IT HAS C NO OPEN ARCS THAT ARE IN THE BACKWARD ADJACENCY LIST OF V IN THE C ORIGINAL NETWORK, WE SET DNPTRB(V)=0. C C IN THIS BFS SEARCH WE ONLY USE OPEN ARCS THAT LEAD US TO NEW C NODES, OR TO NODES THAT HAVE ALREADY BEEN PLACED IN THE NEXT C LAYER. INITIALLY, THE QUEUE CONTAINS ONLY THE source. AS THE C VERTICES ARE POPPED FROM THE QUEUE AND SCANNED, NEW NODES ARE C INJECTED INTO THE QUEUE. EVENTUALLY, EITHER THE SINK IS REACHED, C OR SOME VERTICES ARE NOT REACHABLE WITH OPEN ARCS FROM THE C CURRENT LAYER. IN THE FORMER CASE THE LAYERED NETWORK FOR THE C CURRENT STAGE IS COMPLETE, AND THUS A FLOW AUGMENTATION IS C POSSIBLE. IN THE LATTER CASE, THE CURRENT FLOW ON THE ORIGINAL C NETWORK IS MAXIMUM, AND THE RUN TERMINATES. C C THE QUEUE IS MAINTAINED IN THE LIST DNLIST(.) WITH TWO POINTERS, C AS SHOWN: C ------------------------------------------- C ARRAY | NODES TO BE SCANNED | C DNLIST(.): ... | | C ------------------------------------------- C ^ ^ C | | C QHEAD2 QTAIL2 C******************************************************************* C C INITIALIZE: DO 10 I2=1,DNNODE DNPTRF(I2)=0 DNPTRB(I2)=0 10 DNDIST(I2)=DNNOP2 DNDIST(DNSRCE)=0 C PUT source INTO QUEUE: QHEAD2=0 QTAIL2=1 DNLIST(1)=DNSRCE MAXD2=DNNODE-1 C C SCAN EACH NODE IN QUEUE: C C--------------- C 100 IF(QHEAD2.EQ.QTAIL2)GOTO2000 QHEAD2=QHEAD2+1 C POP NODE IN FRONT OF QUEUE: K2=DNLIST(QHEAD2) K2DST=DNDIST(K2) IF(K2DST.GE.MAXD2)GOTO100 K2DP1=K2DST+1 C C SCAN NODE K2, I.E. SEARCH OVER FWD ADJACENCY OF K2 AND FOR ARCS C (K2,J2) SUCH THAT J2 IS UNSCANNED AND J2 IS NOT IN THE QUEUE AND C ARC (K2,J2) HAS POSITIVE RESIDUAL CAPACITY IN ORIGINAL NETWORK. C 'J2 UNSCANNED' IS CHECKED BY THE CONDITION: 'DNDIST(J2)>=DNDIST(K2)'. C 'J2 NOT IN QUEUE' IS CHECKED BY THE CONDITION C 'DNDIST(J2)=DNDIST(K2)+1 FOR A SCANNED J2'. C IBEG4=DNFAPT(K2) IEND4=DNFAPT(K2+1)-1 IF(IBEG4.GT.IEND4)GOTO1000 DO 400 I4=IBEG4,IEND4 J2=DNFADJ(I4) IF(DNDIST(J2).LE.K2DST)GOTO300 IF(DNCAP(I4).LE.DNFLOW(I4))GOTO300 C (IMPLICITLY) MARK ARC (K2,J2) AS OPEN IN LAYERED NETWORK; ALSO, C SAVE ITS INDEX IN DNPTRF(.): DNPTRF(K2)=I4 IF(J2.NE.DNSINK)GOTO200 C J2=SINK: MAXD2=K2DP1 DNDIST(J2)=MAXD2 GOTO400 C APPEND NODE J2 TO THE QUEUE, IF NOT ALREADY THERE: 200 IF(DNDIST(J2).EQ.K2DP1)GOTO400 DNDIST(J2)=K2DP1 QTAIL2=QTAIL2+1 DNLIST(QTAIL2)=J2 GOTO400 C (EXPLICITLY) MARK ARC (K2,J2) AS CLOSED IN LAYERED NETWORK: 300 DNFADJ(I4)=-DNFADJ(I4) 400 CONTINUE C C--------------- C C IF SINK WAS REACHED, DON'T NEED TO SCAN BACKWARD ARCS INTO ITS LAYER: 1000 IF(DNDIST(DNSINK).NE.DNNOP2)GOTO100 C C ELSE, CONTINUE SCAN OVER FORWARD ARCSOUT OF K2 IN THE RESIDUAL C NETWORK (WHICH ARE BACKWARD ARCS OF THE ORIGINAL NETWORK WITH C POSITIVE FLOW). THIS CODE SEGMENT IS ANALOGOUS TO THE ONE ABOVE: C IBEG4=DNBAPT(K2) IEND4=DNBAPT(K2+1)-1 IF(IBEG4.GT.IEND4)GOTO100 DO 1400 I4=IBEG4,IEND4 J2=DNBADJ(I4) IF(DNDIST(J2).LE.K2DST)GOTO1300 JJDN=DNBTOF(I4) IF(DNFLOW(JJDN).LE.0)GOTO1300 DNPTRB(K2)=I4 IF(DNDIST(J2).EQ.K2DP1)GOTO1400 DNDIST(J2)=K2DP1 QTAIL2=QTAIL2+1 DNLIST(QTAIL2)=J2 GOTO1400 1300 DNBADJ(I4)=-DNBADJ(I4) 1400 CONTINUE GOTO100 C C--------------- C 2000 LAYER2=LAYER2+1 LAYP12=LAYER2+1 C HERE, ALL NODES IN QUEUE HAVEBEEN PROCESSED. C IF THE SINK WAS REACHED, THE LAYERED NETWORK IS COMPLETE: LAYMAX=0 IF(DNDIST(DNSINK).NE.DNNOP2)RETURN C ELSE, THE CURRENT FLOW IS MAXIMUM: LAYMAX=1 RETURN C END SUBROUTINE DNDFS C C C DIMENSION THESE TO AT LEAST MAX NUMBER OF NODES + 2: INTEGER DNLIST(65538),DNDIST(65538) INTEGER DNFAPT(65538),DNPTRF(65538),DNBAPT(65538), * DNFLAB(65538),DNPTRB(65538) C DIMENSION THESE TO AT LEAST MAX NUMBER OF ARCS + 1: INTEGER DNFADJ(524273),DNFROM(524273),DNTO(524273),DNBADJ(524273) INTEGER DNCAP(524273),DNGCAP(524273),DNFLOW(524273),DNBTOF(524273) C COMMON /DN01/DNFAPT/DN02/DNFADJ/DN03/DNCAP/DN04/DNFROM COMMON /DN05/DNFLOW/DN06/DNGCAP/DN08/DNPTRF/DN09/DNLIST COMMON /DN10/DNBAPT/DN11/DNFLAB/DN12/DNPTRB C INTEGER DNNODE,NODP12,DNSRCE,DNSINK,I2,II2,K2,J2,IHEAD2, * ITAIL2,QHEAD2,QTAIL2,KP2,KSAT2,DNNOP2,LAYMAX,K2DST,K2DP1,MAXD2 INTEGER DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE,DNIBIG COMMON /DN00/DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE, * DNELTM,DNIBIG,DNNODE,DNSRCE,DNSINK,DNNOP2 EQUIVALENCE (DNFLOW(1),DNTO(1)),(DNFLAB(1),DNDIST(1)), * (DNFROM(1),DNBADJ(1)),(DNGCAP(1),DNBTOF(1)) C C*********************************************************************** C C THIS SUBROUTINE FINDS A MAXIMAL FLOWIN THE LAYERED NETWORK. IT C USES A DEPTH-FIRST SEARCH STARTING FROM THE SOURCE AND RECORDS C THE DFS TREE BY THE PREDECESSOR ARRAY DNLIST(.) AND FLOW LABEL C ARRAY DNFLAB(.). LAYERED NETWORK ARCS THAT ARE FORWARD IN THE C ORIGINAL NETWORK ARE SEARCHED BEFOREBACKWARD ARCS. THE C FORWARD AND BACKWARD ADJACENCIES OF A NODE K2 ARE SCANNED INLIFO C ORDER, STARTING FROM THE LAST ARC THAT WAS VISITED DURING A C PREVIOUS SCAN OF K2. THUS, WITH THE AID OF THE TWO POINTER LISTS C DNPTRF(.) AND DNPTRB(.), THE ADJACENCY LISTSOF THE ORIGINAL C NETWORK ARE SCANNED ONLY ONCE BY THIS SUBROUTINE. DNPTRF(K2)=0 C IMPLIES THAT NODE K2HAS NO OUTGOING UNSCANNED ARCS; DNPTRB(K2)=0 C IMPLIES THAT NODE K2HAS NO INCOMING UNSCANNED ARCS; A 'CLOSED' C NODE IS DETECTED BY THE CONDITION DNPTRF(K2)=DNPTRB(K2)=0. C C*********************************************************************** C C CALLING CONDITIONS: C INTERNAL CALL FROM DNSUB ONLY. C INPUT ARRAYS: C DNPTRB,DNPTR, DNBTOF, DNCAP, DNFLOW, DNBADJ, DNFADJ C OUTPUT ARRAYS: C DNPTRF,DNPTRB, DNLIST, DNFLAB. C*********************************************************************** C C SUBROUTINES CALLED:DNPUSH C C*********************************************************************** C K2=DNSRCE DNLFVA=0 DNLAUG=0 DNFLAB(K2)=DNIBIG C C----------- C C SCAN NODE K2: C 100 IF(DNPTRF(K2).EQ.0)GOTO1000 C FIND AN OPEN ARC FROM K2 TO SOME NODEJ2: SCAN THE C FORWARD ADJACENCY LIST OF K2, STARTING FROM DNPTRF(K2) AND C PROCEEDING BACKWARD ON DNFADJ() TOWARD DNFAPT(K2): I4=DNPTRF(K2) IPT4=DNFAPT(K2) 300 J2=DNFADJ(I4) IF(J2.GE.0)GOTO600 400 I4=I4-1 IF(I4.GE.IPT4)GOTO300 C NO OPEN ARCS OUT OF K2; MARK K2 AS CLOSED: DNPTRF(K2)=0 GOTO1000 C AN OPEN ARC FOUND; MOVE POINTER DNPTRF() TO THIS ARC; C EXTEND DFS TREE TO NODE J2: 600 DNPTRF(K2)=I4 DNFLAB(J2)=DNFLAB(K2) IF(DNCAP(I4)-DNFLOW(I4).LT.DNFLAB(J2)) * DNFLAB(J2)=DNCAP(I4)-DNFLOW(I4) DNLIST(J2)=K2 K2=J2 C IF NODE K2 (FORMERLY J2) IS THE SINK, AUGMENTFLOW; GET NEW K2: IF(K2.EQ.DNSINK)CALL DNPUSH(K2) GOTO100 C C----------- C 1000 IF(DNPTRB(K2).EQ.0)GOTO2000 C SCAN BACKWARD ARCS INTO K2 TOFIND AN OPEN ARC: I4=DNPTRB(K2) IPT4=DNBAPT(K2) 1300 J2=DNBADJ(I4) IF(J2.GE.0)GOTO1600 I4=I4-1 IF(I4.GE.IPT4)GOTO1300 C NO OPEN ARC INTO K2: DNPTRB(K2)=0 GOTO2000 C AN OPEN ARC FOUND; EXTEND DFS TREE TONODE J2: 1600 DNPTRB(K2)=I4 DNFLAB(J2)=DNFLAB(K2) JJDN=DNBTOF(I4) IFL4=DNFLOW(JJDN) IF(IFL4.LT.DNFLAB(J2))DNFLAB(J2)=IFL4 DNLIST(J2)=-K2 K2=J2 GOTO100 C C----------- C C K2 IS 'CLOSED' NODE; IF IT ISTHE source, THEN WE HAVE A MAXIMAL C FLOW IN LAYERED NETWORK: 2000 IF(K2.EQ.DNSRCE)RETURN C ELSE, BACK UPONE NODE FROM K2 IN DFS TREE: KP2=DNLIST(K2) K2=KP2 IF(K2.LT.0)K2=-K2 IF(KP2.GE.0)GOTO2200 I4=DNPTRB(K2) DNBADJ(I4)=-DNBADJ(I4) GOTO100 2200 I4=DNPTRF(K2) DNFADJ(I4)=-DNFADJ(I4) GOTO100 C C----------- END SUBROUTINE DNPUSH(K2) C C C DIMENSION THESE TO AT LEAST MAX NUMBER OF NODES + 2: INTEGER DNLIST(65538),DNDIST(65538) INTEGER DNFAPT(65538),DNPTRF(65538),DNBAPT(65538), * DNFLAB(65538),DNPTRB(65538) C DIMENSION THESE TO AT LEAST MAX NUMBER OF ARCS + 1: INTEGER DNFADJ(524273),DNFROM(524273),DNTO(524273),DNBADJ(524273) INTEGER DNCAP(524273),DNGCAP(524273),DNFLOW(524273),DNBTOF(524273) C COMMON /DN01/DNFAPT/DN02/DNFADJ/DN03/DNCAP/DN04/DNFROM COMMON /DN05/DNFLOW/DN06/DNGCAP/DN08/DNPTRF/DN09/DNLIST COMMON /DN10/DNBAPT/DN11/DNFLAB/DN12/DNPTRB C INTEGER DNNODE,NODP12,DNSRCE,DNSINK,I2,II2,K2,J2,IHEAD2, * ITAIL2,QHEAD2,QTAIL2,KP2,KSAT2,DNNOP2,LAYMAX,K2DST,K2DP1,MAXD2 INTEGER DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE,DNIBIG COMMON /DN00/DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE, * DNELTM,DNIBIG,DNNODE,DNSRCE,DNSINK,DNNOP2 EQUIVALENCE (DNFLOW(1),DNTO(1)),(DNFLAB(1),DNDIST(1)), * (DNFROM(1),DNBADJ(1)),(DNGCAP(1),DNBTOF(1)) C C*********************************************************************** C C AUGMENT FLOW ALONG FLOW AUGMENTING PATH DEFINED SUBROUTINE DNDFS, I.E. C USING THE PREDECESSORARRAY DNLIST(.), START FROM THE SINK AND C TRAVERSE TO THE SOURCE ARC FLOWS AND FLOW LABELS BY AN AMOUNT EQUAL C TO DNFLAB(DNSINK). MARK SATURATED FORWARD (FROM S TO T) ARCS, AND C BACKWARD ARCS HAVING ZERO FLOW, AS CLOSED. C C*********************************************************************** C C CALLING CONDITIONS: C INTERNAL CALL FROM DNDFS ONLY. C INPUT ARRAYS: C DNFLAB,DNLIST,DNPTRB,DNPTRF,DNBTOF,DNCAP,DNFLOW,DNBADJ,DNFADJ C OUTPUT ARRAYS: C DNFLAB,DNFLOW,DNBADJ,DNFADJ C*********************************************************************** C KSAT2=0 J2=DNSINK INCRE4=DNFLAB(J2) DNLFVA=DNLFVA+INCRE4 DNLAUG=DNLAUG+1 C 100 K2=J2 J2=DNLIST(K2) DNFLAB(K2)=DNFLAB(K2)-INCRE4 IF(J2.GT.0)GOTO200 C C DECREASE FLOW ON BACKWARD ARC(K2,J2): J2=-J2 I4=DNPTRB(J2) II4=DNBTOF(I4) DNFLOW(II4)=DNFLOW(II4)-INCRE4 IF(DNFLOW(II4).NE.0)GOTO100 C FLOW IS ZERO; MARK ARC AS 'CLOSED' INLAYERED NETWORK: DNBADJ(I4)=-DNBADJ(I4) KSAT2=J2 GOTO100 C C INCREASE FLOW ON FORWARD ARC (J2,K2): 200 I4=DNPTRF(J2) DNFLOW(I4)=DNFLOW(I4)+INCRE4 IF(DNCAP(I4).NE.DNFLOW(I4))GOTO300 C ARC IS NOW SATURATED; MARK ISAS 'CLOSED' IN LAYERED NET DNFADJ(I4)=-DNFADJ(I4) KSAT2=J2 300 IF(J2.NE.DNSRCE)GOTO100 C C RETURN TO RESUME SEARCH FROM K2 (NODE CLOSESTTO source OF C THE ARC CLOSED LAST) IF(KSAT2.EQ.0)STOP K2=KSAT2 RETURN END SUBROUTINE DNFWD(NODES,NARCS,IRETN) C C C DIMENSION THESE TO AT LEAST MAX NUMBER OF NODES + 2: INTEGER DNLIST(65538),DNDIST(65538) INTEGER DNFAPT(65538),DNPTRF(65538),DNBAPT(65538), * DNFLAB(65538),DNPTRB(65538) C DIMENSION THESE TO AT LEAST MAX NUMBER OF ARCS + 1: INTEGER DNFADJ(524273),DNFROM(524273),DNTO(524273),DNBADJ(524273) INTEGER DNCAP(524273),DNGCAP(524273),DNFLOW(524273),DNBTOF(524273) C COMMON /DN01/DNFAPT/DN02/DNFADJ/DN03/DNCAP/DN04/DNFROM COMMON /DN05/DNFLOW/DN06/DNGCAP/DN08/DNPTRF/DN09/DNLIST COMMON /DN10/DNBAPT/DN11/DNFLAB/DN12/DNPTRB C INTEGER DNNODE,NODP12,DNSRCE,DNSINK,I2,II2,K2,J2,IHEAD2, * ITAIL2,QHEAD2,QTAIL2,KP2,KSAT2,DNNOP2,LAYMAX,K2DST,K2DP1,MAXD2 INTEGER DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE,DNIBIG COMMON /DN00/DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE, * DNELTM,DNIBIG,DNNODE,DNSRCE,DNSINK,DNNOP2 EQUIVALENCE (DNFLOW(1),DNTO(1)),(DNFLAB(1),DNDIST(1)), * (DNFROM(1),DNBADJ(1)),(DNGCAP(1),DNBTOF(1)) C C*********************************************************************** C C CALLING CONDITIONS: C USER CALLABLE; BEFORE CALLING DNSUB. C INPUT: C IN CALLING SEQUENCE (ALL INTEGER*4): C NODES=NUMBER OF NODES, INCLUDING source AND SINK C NARCS=NUMBER OF ARCS. C C IN COMMON, PROVIDE THREE ARC LISTS WITH ARCS IN ANY ORDER: C DNFROM=THE LIST OF TAILS. (INTEGER*2) C DNTO=THE LIST OF HEADS. (INTEGER*2) C DNGCAP=THE LIST OF CAPACITIES. (INTEGER*4). C C SCRATCH ARRAY (IN COMMON): C DNLIST C OUTPUT: C IN CALLING SEQUENCE (INTEGER*4): C IRETN=0 NO ERRORS; =1 ERROR DETECTED, CHECK INPUT. C IN COMMON, ARRAYS: C DNFAPT,DNFADJ,DNCAP (I.E. THE FWD ADJACENCY INPUT DATA STRUCTURE C THAT IS REQUIRED AS INPUT TO SUBROUTINE DNSUB). C C*********************************************************************** C C USING THE THREE ARC LISTS DNFROM(.), DNTO(.), DNGCAP(.), THIS C SUBROUTINE CONSTRUCTS THE FORWARD ADJACENCY ARRAYS DNFAPT(.), C DNFADJ(.), AND DNCAP(.), AS THE INPUT REQUIRED BY SUBROUTINE DNSUB. C C*********************************************************************** C C INITIALIZE: IRETN=0 NODP1=NODES+1 DO 100 I2=1,NODP1 DNFAPT(I2)=0 100 DNLIST(I2)=0 C TEMPORARILY STORE IN DNFAPT(.) NUMBER OF ARCSOUT OF EACH NODE: DO 200 I4=1,NARCS JJDN=DNFROM(I4) 200 DNFAPT(JJDN)=DNFAPT(JJDN)+1 C CONSTRUCT DNFAPT(.): ITPI4=DNFAPT(1) DNFAPT(1)=1 DO 300 I2=1,NODES ITPSV4=ITPI4+DNFAPT(I2) ITPI4=DNFAPT(I2+1) 300 DNFAPT(I2+1)=ITPSV4 C CONSTRUCT DNFADJ(.) AND DNCAP(.): DO 400 I4=1,NARCS ITAIL2=DNFROM(I4) ITPUT4=DNFAPT(ITAIL2)+DNLIST(ITAIL2) IF(ITPUT4.LE.0)GOTO500 DNFADJ(ITPUT4)=DNTO(I4) DNLIST(ITAIL2)=DNLIST(ITAIL2)+1 DNCAP(ITPUT4)=DNGCAP(I4) 400 CONTINUE RETURN 500 IRETN=1 RETURN END SUBROUTINE DNCLEA C C C DIMENSION THESE TO AT LEAST MAX NUMBER OF NODES + 2: INTEGER DNLIST(65538),DNDIST(65538) INTEGER DNFAPT(65538),DNPTRF(65538),DNBAPT(65538), * DNFLAB(65538),DNPTRB(65538) C DIMENSION THESE TO AT LEAST MAX NUMBER OF ARCS + 1: INTEGER DNFADJ(524273),DNFROM(524273),DNTO(524273),DNBADJ(524273) INTEGER DNCAP(524273),DNGCAP(524273),DNFLOW(524273),DNBTOF(524273) C COMMON /DN01/DNFAPT/DN02/DNFADJ/DN03/DNCAP/DN04/DNFROM COMMON /DN05/DNFLOW/DN06/DNGCAP/DN08/DNPTRF/DN09/DNLIST COMMON /DN10/DNBAPT/DN11/DNFLAB/DN12/DNPTRB C INTEGER DNNODE,NODP12,DNSRCE,DNSINK,I2,II2,K2,J2,IHEAD2, * ITAIL2,QHEAD2,QTAIL2,KP2,KSAT2,DNNOP2,LAYMAX,K2DST,K2DP1,MAXD2 INTEGER DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE,DNIBIG COMMON /DN00/DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE, * DNELTM,DNIBIG,DNNODE,DNSRCE,DNSINK,DNNOP2 EQUIVALENCE (DNFLOW(1),DNTO(1)),(DNFLAB(1),DNDIST(1)), * (DNFROM(1),DNBADJ(1)),(DNGCAP(1),DNBTOF(1)) C C******************************************************************* C CALLING CONDITIONS: C INTERNAL CALL FROM DNSUB ONLY. C INPUT (IN COMMON): C SCALARS: DNARC C ARRAYS: DNFADJ, DNBADJ C OUTPUT: C ARRAYS (IN COMMON): DNFADJ,DNBADJ C******************************************************************** C C CLEAR MARKINGS IN ARRAYS DNFADJ() ANDDNBADJ().: DO 100 I4=1,DNARC IF(DNBADJ(I4).LT.0)DNBADJ(I4)=-DNBADJ(I4) IF(DNFADJ(I4).LT.0)DNFADJ(I4)=-DNFADJ(I4) 100 CONTINUE RETURN END SUBROUTINE DNCUT(NNCUT,NACUT) C C C DIMENSION THESE TO AT LEAST MAX NUMBER OF NODES + 2: INTEGER DNLIST(65538),DNDIST(65538) INTEGER DNFAPT(65538),DNPTRF(65538),DNBAPT(65538), * DNFLAB(65538),DNPTRB(65538) C DIMENSION THESE TO AT LEAST MAX NUMBER OF ARCS + 1: INTEGER DNFADJ(524273),DNFROM(524273),DNTO(524273),DNBADJ(524273) INTEGER DNCAP(524273),DNGCAP(524273),DNFLOW(524273),DNBTOF(524273) C COMMON /DN01/DNFAPT/DN02/DNFADJ/DN03/DNCAP/DN04/DNFROM COMMON /DN05/DNFLOW/DN06/DNGCAP/DN08/DNPTRF/DN09/DNLIST COMMON /DN10/DNBAPT/DN11/DNFLAB/DN12/DNPTRB C INTEGER DNNODE,NODP12,DNSRCE,DNSINK,I2,II2,K2,J2,IHEAD2, * ITAIL2,QHEAD2,QTAIL2,KP2,KSAT2,DNNOP2,LAYMAX,K2DST,K2DP1,MAXD2 INTEGER DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE,DNIBIG COMMON /DN00/DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE, * DNELTM,DNIBIG,DNNODE,DNSRCE,DNSINK,DNNOP2 EQUIVALENCE (DNFLOW(1),DNTO(1)),(DNFLAB(1),DNDIST(1)), * (DNFROM(1),DNBADJ(1)),(DNGCAP(1),DNBTOF(1)) C C*********************************************************************** C CALLING CONDITIONS: C ONLY AFTER A SUCCESSFULRETURN FROM A CALL TO DNSUB. C INPUT ARRAYS: C DNFAPT,DNFADJ, DNCAP,DNFLOW, DNDIST C OUTPUT ARRAYS: C DNFAPT,DNFADJ C PASSED IN CALLING SEQUENCE: C NNCUT:NUMBER OF NODES ON source SIDE OF FINAL CUT (INTEGER*4) C NACUT:NUMBER OF ARCS IN THE FINAL CUT (INTEGER*4) C*********************************************************************** C C NEGATE DNFAPT(.) FOR NODES ON SINK SIDE OF CUT: NNCUT=0 DO 1700I2=1,DNNODE IF(DNDIST(I2).NE.DNNOP2)GOTO1650 DNFAPT(I2)=-IABS(DNFAPT(I2)) GOTO1700 1650 NNCUT=NNCUT+1 1700 CONTINUE C NEGATE DNFADJ(.) FOR THOSE SATURATED ARCS IN MIN CUT: NACUT=0 DO 1900I2=1,DNNODE IF(DNFAPT(I2).LT.0)GOTO1900 C NODE I2 IS ON source SIDE: IBEG4=IABS(DNFAPT(I2)) IEND4=IABS(DNFAPT(I2+1))-1 DO 1800 I4=IBEG4,IEND4 JJDN=DNFADJ(I4) IF(DNFAPT(JJDN).GT.0)GOTO1800 C NODE DNFADJ(I4) ON SINK SIDE: IF(DNFLOW(I4).NE.DNCAP(I4))GOTO1800 II2=DNFADJ(I4) DNFADJ(I4)=-IABS(II2) NACUT=NACUT+1 1800 CONTINUE 1900 CONTINUE RETURN END SUBROUTINE DNOUT(IUNREP) C CHARACTER*1 T , JJS, JJT , ISIDE, JSIDE CHARACTER*3 BL3, BSAT, TAG1 CHARACTER*4 BL4, BCUT, TAG2 C C C DIMENSION THESE TO AT LEAST MAX NUMBER OF NODES + 2: INTEGER DNLIST(65538),DNDIST(65538) INTEGER DNFAPT(65538),DNPTRF(65538),DNBAPT(65538), * DNFLAB(65538),DNPTRB(65538) C DIMENSION THESE TO AT LEAST MAX NUMBER OF ARCS + 1: INTEGER DNFADJ(524273),DNFROM(524273),DNTO(524273),DNBADJ(524273) INTEGER DNCAP(524273),DNGCAP(524273),DNFLOW(524273),DNBTOF(524273) C COMMON /DN01/DNFAPT/DN02/DNFADJ/DN03/DNCAP/DN04/DNFROM COMMON /DN05/DNFLOW/DN06/DNGCAP/DN08/DNPTRF/DN09/DNLIST COMMON /DN10/DNBAPT/DN11/DNFLAB/DN12/DNPTRB C INTEGER DNNODE,NODP12,DNSRCE,DNSINK,I2,II2,K2,J2,IHEAD2, * ITAIL2,QHEAD2,QTAIL2,KP2,KSAT2,DNNOP2,LAYMAX,K2DST,K2DP1,MAXD2 INTEGER DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE,DNIBIG COMMON /DN00/DNARC,DNFVA,DNAUG,DNLFVA,DNLAUG,DNSTGE, * DNELTM,DNIBIG,DNNODE,DNSRCE,DNSINK,DNNOP2 EQUIVALENCE (DNFLOW(1),DNTO(1)),(DNFLAB(1),DNDIST(1)), * (DNFROM(1),DNBADJ(1)),(DNGCAP(1),DNBTOF(1)) DATA BL3/3H /, BL4/4H /, BSAT/3HSAT/, BCUT/4H-CUT/, 1 JJT/1HT/ DATA JJS/1HS/ C C*********************************************************************** C C PREPARES OUTPUT REPORT. MAY ONLY BE CALLED AFTER A CALL TO DNSUB. C C*********************************************************************** C C CALLING CONDITIONS: C USER CALLABLE; ONLY AFTER CALLING DNSUB. C INPUT ARRAYS: C DNFAPT,DNFADJ, DNCAP, DNFLOW C SCRATCH ARRAYS: C DNLIST C OUTPUT ARRAYS: NONE C*********************************************************************** C C PRINT SUMMARYSECTION: WRITE(IUNREP,1000) * DNNODE,DNARC,DNSRCE,DNSINK,DNFVA,DNELTM,DNAUG,DNSTGE C PRINT CUT-SETNODES: K2=0 DO 200 I2=1,DNNODE IF(DNFAPT(I2).LT.0)GOTO200 K2=K2+1 DNLIST(K2)=I2 200 CONTINUE c WRITE(IUNREP,1003)K2 c WRITE(IUNREP,1004)(DNLIST(I2),I2=1,K2) c WRITE(IUNREP,1001) C PRINT ARC LIST: NTOTC4=0 NSAINC=0 NSA=0 NPOSFL=0 DO 400 I2=1,DNNODE IBEG4=IABS(DNFAPT(I2)) IEND4=IABS(DNFAPT(I2+1))-1 DO 300 I4=IBEG4,IEND4 MFLO4=DNFLOW(I4) IF(MFLO4.EQ.0)GOTO300 NPOSFL=NPOSFL+1 TAG1=BL3 IF(MFLO4.NE.DNCAP(I4))GOTO270 TAG1=BSAT NSA=NSA+1 270 TAG2=BL4 C IF(DNFADJ(I4).LT.0)TAG2=BCUT IF(DNFADJ(I4).GE.0)GOTO280 NTOTC4=NTOTC4+MFLO4 NSAINC=NSAINC+1 280 ISIDE=JJS JSIDE=ISIDE IF(DNFAPT(I2).LT.0)ISIDE=JJT IFAPLU=DNFADJ(I4) IFAPLU=IABS(IFAPLU) IF(DNFAPT(IFAPLU).LT.0)JSIDE=JJT IF(ISIDE.NE.JSIDE)TAG2=BCUT C c WRITE(IUNREP,1002)I2,ISIDE,IFAPLU,JSIDE, c * MFLO4,DNCAP(I4),TAG1,TAG2 300 CONTINUE 400 CONTINUE C PRINT MAXFLO VALUE AND COMPUTED MINCUT c WRITE(IUNREP,1006)NPOSFL,NSA,NSAINC c IF(DNFVA.NE.NTOTC4)WRITE(IUNREP,1007)DNFVA,NTOTC4 c WRITE(IUNREP,1008) RETURN C C FORMATS: 1000 FORMAT(/13H DINIC MAXFLO, */13H PROBLEM HAS ,I6,7H NODES,,I6,24H ARCS, SOURCE/SINK PAIR,,2I6 */16H MAX FLOW VALUE=,I16,21H, PC/AT ELAPSED TIME ,F8.2,6H SECS., */21H TOTAL AUGMENTATIONS=,I6,4H IN ,I6,7H STAGES) 1001 FORMAT(/25H ARCS WITH POSITIVE FLOW , *//23H TAIL HEAD, * 40H FLOW CAPACITY STATUS/) 1002 FORMAT(I10,1X,A1,I10,1X,A1,2I15,4X,A3,A4) 1003 FORMAT(/10H THERE ARE,I8,29H NODES IN SOURCE SIDE OF CUT,, * 10H THEY ARE /) 1004 FORMAT(10I6) 1006 FORMAT(/7H NOTES , */10X,13H 1) THERE ARE,I6,26H ARCS WITH POSITIVE FLOW, ,I5,3H OF */10X,22H THEM SATURATED (,I5,26H OF WHICH ARE IN MINCUT). , */10X,52H 2) S / T REFER TO SOURCE / SINK SIDES OF CUT.) 1007 FORMAT( * 10X,41H 3) MAX FLOW VALUE FROM SOURCE TO SINK IS,I14,1H, */10X,41H BUT, CUT CAPACITY WAS COMPUTED AS ,I14) 1008 FORMAT(/15H END OF REPORT.) END