504 lines
9.7 KiB
Rexx
504 lines
9.7 KiB
Rexx
/* REXX ***************************************************************
|
|
* Solve the Transportation Problem using Vogel's Approximation
|
|
Default Input
|
|
2 3 # of sources / # of demands
|
|
25 35 sources
|
|
20 30 10 demands
|
|
3 5 7 cost matrix <
|
|
3 2 5
|
|
* 20201210 support no input file -courtesy GS
|
|
* Note: correctness of input is not checked
|
|
* 20210102 restored Vogel's Approximation and added Optimization
|
|
* 20210103 eliminated debug code
|
|
**********************************************************************/
|
|
Signal On Halt
|
|
Signal On Novalue
|
|
Signal On Syntax
|
|
|
|
Parse Arg fid
|
|
If fid='' Then
|
|
fid='input1.txt'
|
|
Call init
|
|
m.=0
|
|
Do Forever
|
|
dmax.=0
|
|
dmax=0
|
|
Do r=1 To rr
|
|
dr.r=''
|
|
Do c=1 To cc
|
|
If cost.r.c<>'*' Then
|
|
dr.r=dr.r cost.r.c
|
|
End
|
|
dr.r=words(dr.r) dr.r
|
|
dr.r=diff(dr.r)
|
|
If dr.r>dmax Then Do; dmax=dr.r; dmax.0='R'; dmax.1=r; dmax.2=dr.r; End
|
|
End
|
|
Do c=1 To cc
|
|
dc.c=''
|
|
Do r=1 To rr
|
|
If cost.r.c<>'*' Then
|
|
dc.c=dc.c cost.r.c
|
|
End
|
|
dc.c=words(dc.c) dc.c
|
|
dc.c=diff(dc.c)
|
|
If dc.c>dmax Then Do; dmax=dc.c; dmax.0='C'; dmax.1=c; dmax.2=dc.c; End
|
|
End
|
|
cmin=999
|
|
Select
|
|
When dmax.0='R' Then Do
|
|
r=dmax.1
|
|
Do c=1 To cc
|
|
If cost.r.c<>'*' &,
|
|
cost.r.c<cmin Then Do
|
|
cmin=cost.r.c
|
|
cx=c
|
|
End
|
|
End
|
|
Call allocate r cx
|
|
End
|
|
When dmax.0='C' Then Do
|
|
c=dmax.1
|
|
Do r=1 To rr
|
|
If cost.r.c<>'*' &,
|
|
cost.r.c<cmin Then Do
|
|
cmin=cost.r.c
|
|
rx=r
|
|
End
|
|
End
|
|
Call allocate rx c
|
|
End
|
|
Otherwise
|
|
Leave
|
|
End
|
|
End
|
|
|
|
Do r=1 To rr
|
|
Do c=1 To cc
|
|
If cost.r.c<>'*' Then Do
|
|
Call allocate r c
|
|
cost.r.c='*'
|
|
End
|
|
End
|
|
End
|
|
|
|
Call show_alloc 'Vogel''s Approximation'
|
|
|
|
Do r=1 To rr
|
|
Do c=1 To cc
|
|
cost.r.c=word(matrix.r.c,3) /* restore cost.*.* */
|
|
End
|
|
End
|
|
|
|
Call steppingstone
|
|
Exit
|
|
|
|
/**********************************************************************
|
|
* Subroutines for Vogel's Approximation
|
|
**********************************************************************/
|
|
|
|
init:
|
|
If lines(fid)=0 Then Do
|
|
Say 'Input file not specified or not found. Using default input instead.'
|
|
fid='Default input'
|
|
in.1=sourceline(4)
|
|
Parse Var in.1 numSources .
|
|
Do i=2 To numSources+3
|
|
in.i=sourceline(i+3)
|
|
End
|
|
End
|
|
Else Do
|
|
Do i=1 By 1 while lines(fid)>0
|
|
in.i=linein(fid)
|
|
End
|
|
End
|
|
Parse Var in.1 numSources numDestinations . 1 rr cc .
|
|
source.=0
|
|
demand.=0
|
|
source_sum=0
|
|
Do i=1 To numSources
|
|
Parse Var in.2 source.i in.2
|
|
ss.i=source.i
|
|
source_in.i=source.i
|
|
source_sum=source_sum+source.i
|
|
End
|
|
l=linein(fid)
|
|
demand_sum=0
|
|
Do i=1 To numDestinations
|
|
Parse Var in.3 demand.i in.3
|
|
dd.i=demand.i
|
|
demand_in.i=demand.i
|
|
demand_sum=demand_sum+demand.i
|
|
End
|
|
Do i=1 To numSources
|
|
j=i+3
|
|
l=in.j
|
|
Do j=1 To numDestinations
|
|
Parse Var l cost.i.j l
|
|
End
|
|
End
|
|
Do i=1 To numSources
|
|
ol=format(source.i,3)
|
|
Do j=1 To numDestinations
|
|
ol=ol format(cost.i.j,4)
|
|
End
|
|
End
|
|
ol=' '
|
|
Do j=1 To numDestinations
|
|
ol=ol format(demand.j,4)
|
|
End
|
|
|
|
Select
|
|
When source_sum=demand_sum Then Nop /* balanced */
|
|
When source_sum>demand_sum Then Do /* unbalanced - add dummy demand */
|
|
Say 'This is an unbalanced case (sources exceed demands). We add a dummy consumer.'
|
|
cc=cc+1
|
|
demand.cc=source_sum-demand_sum
|
|
demand_in.cc=demand.cc
|
|
dd.cc=demand.cc
|
|
Do r=1 To rr
|
|
cost.r.cc=0
|
|
End
|
|
End
|
|
Otherwise /* demand_sum>source_sum */ Do /* unbalanced - add dummy source */
|
|
Say 'This is an unbalanced case (demands exceed sources). We add a dummy source.'
|
|
rr=rr+1
|
|
source.rr=demand_sum-source_sum
|
|
source_in.rr=source.rr
|
|
ss.rr=source.rr
|
|
Do c=1 To cc
|
|
cost.rr.c=0
|
|
End
|
|
End
|
|
End
|
|
|
|
Say 'Sources / Demands / Cost'
|
|
ol=' '
|
|
Do c=1 To cc
|
|
ol=ol format(demand.c,3)
|
|
End
|
|
Say ol
|
|
|
|
Do r=1 To rr
|
|
ol=format(source.r,4)
|
|
Do c=1 To cc
|
|
ol=ol format(cost.r.c,3)
|
|
matrix.r.c=r c cost.r.c 0
|
|
End
|
|
Say ol
|
|
End
|
|
Return
|
|
|
|
allocate: Procedure Expose m. source. demand. cost. rr cc matrix.
|
|
Parse Arg r c
|
|
sh=min(source.r,demand.c)
|
|
source.r=source.r-sh
|
|
demand.c=demand.c-sh
|
|
m.r.c=sh
|
|
matrix.r.c=subword(matrix.r.c,1,3) sh
|
|
If source.r=0 Then Do
|
|
Do c=1 To cc
|
|
cost.r.c='*'
|
|
End
|
|
End
|
|
If demand.c=0 Then Do
|
|
Do r=1 To rr
|
|
cost.r.c='*'
|
|
End
|
|
End
|
|
Return
|
|
|
|
diff: Procedure
|
|
Parse Value arg(1) With n list
|
|
If n<2 Then Return 0
|
|
list=wordsort(list)
|
|
Return word(list,2)-word(list,1)
|
|
|
|
wordsort: Procedure
|
|
/**********************************************************************
|
|
* Sort the list of words supplied as argument. Return the sorted list
|
|
**********************************************************************/
|
|
Parse Arg wl
|
|
wa.=''
|
|
wa.0=0
|
|
Do While wl<>''
|
|
Parse Var wl w wl
|
|
Do i=1 To wa.0
|
|
If wa.i>w Then Leave
|
|
End
|
|
If i<=wa.0 Then Do
|
|
Do j=wa.0 To i By -1
|
|
ii=j+1
|
|
wa.ii=wa.j
|
|
End
|
|
End
|
|
wa.i=w
|
|
wa.0=wa.0+1
|
|
End
|
|
swl=''
|
|
Do i=1 To wa.0
|
|
swl=swl wa.i
|
|
End
|
|
/* Say swl */
|
|
Return strip(swl)
|
|
|
|
show_alloc: Procedure Expose matrix. rr cc demand_in. source_in.
|
|
Parse Arg header
|
|
If header='' Then
|
|
Return
|
|
Say ''
|
|
Say header
|
|
total=0
|
|
ol=' '
|
|
Do c=1 to cc
|
|
ol=ol format(demand_in.c,3)
|
|
End
|
|
Say ol
|
|
as=''
|
|
Do r=1 to rr
|
|
ol=format(source_in.r,4)
|
|
a=word(matrix.r.1,4)
|
|
If a=0.0000000001 Then a=0
|
|
If a>0 Then
|
|
ol=ol format(a,3)
|
|
Else
|
|
ol=ol ' - '
|
|
total=total+word(matrix.r.1,4)*word(matrix.r.1,3)
|
|
Do c=2 To cc
|
|
a=word(matrix.r.c,4)
|
|
If a=0.0000000001 Then a=0
|
|
If a>0 Then
|
|
ol=ol format(a,3)
|
|
Else
|
|
ol=ol ' - '
|
|
total=total+word(matrix.r.c,4)*word(matrix.r.c,3)
|
|
as=as a
|
|
End
|
|
Say ol
|
|
End
|
|
Say 'Total costs:' format(total,4,1)
|
|
Return
|
|
|
|
|
|
/**********************************************************************
|
|
* Subroutines for Optimization
|
|
**********************************************************************/
|
|
|
|
steppingstone: Procedure Expose matrix. cost. rr cc matrix. demand_in.,
|
|
source_in. ms fid move cnt.
|
|
maxReduction=0
|
|
move=''
|
|
Call fixDegenerateCase
|
|
Do r=1 To rr
|
|
Do c=1 To cc
|
|
Parse Var matrix.r.c r c cost qrc
|
|
If qrc=0 Then Do
|
|
path=getclosedpath(r,c)
|
|
If pelems(path)<4 Then Do
|
|
Iterate
|
|
End
|
|
reduction = 0
|
|
lowestQuantity = 1e10
|
|
leavingCandidate = ''
|
|
plus=1
|
|
pathx=path
|
|
Do While pathx<>''
|
|
Parse Var pathx s '|' pathx
|
|
If plus Then
|
|
reduction=reduction+word(s,3)
|
|
Else Do
|
|
reduction=reduction-word(s,3)
|
|
If word(s,4)<lowestQuantity Then Do
|
|
leavingCandidate = s
|
|
lowestQuantity = word(s,4)
|
|
End
|
|
End
|
|
plus=\plus
|
|
End
|
|
If reduction < maxreduction Then Do
|
|
move=path
|
|
leaving=leavingCandidate
|
|
maxReduction = reduction
|
|
End
|
|
End
|
|
End
|
|
End
|
|
if move<>'' Then Do
|
|
quant=word(leaving,4)
|
|
If quant=0 Then Do
|
|
Call show_alloc 'Optimum'
|
|
Exit
|
|
End
|
|
plus=1
|
|
Do While move<>''
|
|
Parse Var move m '|' move
|
|
Parse Var m r c cpu qrc
|
|
Parse Var matrix.r.c vr vc vcost vquant
|
|
If plus Then
|
|
nquant=vquant+quant
|
|
Else
|
|
nquant=vquant-quant
|
|
matrix.r.c = vr vc vcost nquant
|
|
plus=\plus
|
|
End
|
|
move=''
|
|
Call steppingStone
|
|
End
|
|
Else
|
|
Call show_alloc 'Optimal Solution' fid
|
|
Return
|
|
|
|
getclosedpath: Procedure Expose matrix. cost. rr cc matrix.
|
|
Parse Arg rd,cd
|
|
path=rd cd cost.rd.cd word(matrix.rd.cd,4)
|
|
do r=1 To rr
|
|
Do c=1 To cc
|
|
If word(matrix.r.c,4)>0 Then Do
|
|
path=path'|'r c cost.r.c word(matrix.r.c,4)
|
|
End
|
|
End
|
|
End
|
|
path=magic(path)
|
|
Return stones(path)
|
|
|
|
magic: Procedure
|
|
Parse Arg list
|
|
Do Forever
|
|
list_1=remove_1(list)
|
|
If list_1=list Then Leave
|
|
list=list_1
|
|
End
|
|
Return list_1
|
|
|
|
remove_1: Procedure
|
|
Parse Arg list
|
|
cntr.=0
|
|
cntc.=0
|
|
Do i=1 By 1 While list<>''
|
|
parse Var list e.i '|' list
|
|
Parse Var e.i r c .
|
|
cntr.r=cntr.r+1
|
|
cntc.c=cntc.c+1
|
|
End
|
|
n=i-1
|
|
keep.=1
|
|
Do i=1 To n
|
|
Parse Var e.i r c .
|
|
If cntr.r<2 |,
|
|
cntc.c<2 Then Do
|
|
keep.i=0
|
|
End
|
|
End
|
|
list=e.1
|
|
Do i=2 To n
|
|
If keep.i Then
|
|
list=list'|'e.i
|
|
End
|
|
Return list
|
|
|
|
stones: Procedure
|
|
Parse Arg lst
|
|
tstc=lst
|
|
Do i=1 By 1 While tstc<>''
|
|
Parse Var tstc o.i '|' tstc
|
|
end
|
|
stones=lst
|
|
o.0=i-1
|
|
prev=o.1
|
|
Do i=1 To o.0
|
|
st.i=prev
|
|
k=i//2
|
|
nbrs=getNeighbors(prev,lst)
|
|
Parse Var nbrs n.1 '|' n.2
|
|
If k=0 Then
|
|
prev=n.2
|
|
Else
|
|
prev=n.1
|
|
End
|
|
stones=st.1
|
|
Do i=2 To o.0
|
|
stones=stones'|'st.i
|
|
End
|
|
Return stones
|
|
|
|
getNeighbors: Procedure Expose o.
|
|
parse Arg s, lst
|
|
Do i=1 To 4
|
|
Parse Var lst o.i '|' lst
|
|
End
|
|
nbrs.=''
|
|
sr=word(s,1)
|
|
sc=word(s,2)
|
|
Do i=1 To o.0
|
|
If o.i<>s Then Do
|
|
or=word(o.i,1)
|
|
oc=word(o.i,2)
|
|
If or=sr & nbrs.0='' Then
|
|
nbrs.0 = o.i
|
|
else if oc=sc & nbrs.1='' Then
|
|
nbrs.1 = o.i
|
|
If nbrs.0<>'' & nbrs.1<>'' Then
|
|
Leave
|
|
End
|
|
End
|
|
return nbrs.0'|'nbrs.1
|
|
|
|
m1: Procedure
|
|
Parse Arg z
|
|
Return z-1
|
|
|
|
pelems: Procedure
|
|
Call Trace 'O'
|
|
Parse Arg p
|
|
n=0
|
|
Do While p<>''
|
|
Parse Var p x '|' p
|
|
If x<>'' Then n=n+1
|
|
End
|
|
Return n
|
|
|
|
fixDegenerateCase: Procedure Expose matrix. rr cc ms
|
|
Call matrixtolist
|
|
If (rr+cc-1)<>ms Then Do
|
|
Do r=1 To rr
|
|
Do c=1 To cc
|
|
If word(matrix.r.c,4)=0 Then Do
|
|
matrix.r.c=subword(matrix.r.c,1,3) 1.e-10
|
|
Return
|
|
End
|
|
End
|
|
End
|
|
End
|
|
Return
|
|
|
|
matrixtolist: Procedure Expose matrix. rr cc ms
|
|
ms=0
|
|
list=''
|
|
Do r=1 To rr
|
|
Do c=1 To cc
|
|
If word(matrix.r.c,4)>0 Then Do
|
|
list=list'|'matrix.r.c
|
|
ms=ms+1
|
|
End
|
|
End
|
|
End
|
|
Return strip(list,,'|')
|
|
|
|
Novalue:
|
|
Say 'Novalue raised in line' sigl
|
|
Say sourceline(sigl)
|
|
Say 'Variable' condition('D')
|
|
Signal lookaround
|
|
|
|
Syntax:
|
|
Say 'Syntax raised in line' sigl
|
|
Say sourceline(sigl)
|
|
Say 'rc='rc '('errortext(rc)')'
|
|
|
|
halt:
|
|
lookaround:
|
|
If fore() Then Do
|
|
Say 'You can look around now.'
|
|
Trace ?R
|
|
Nop
|
|
End
|
|
Exit 12
|