Kidney Donor Exchange (dcmpe11)
/***************************************************************/
/* */
/* S A S S A M P L E L I B R A R Y */
/* */
/* NAME: dcmpe11 */
/* TITLE: Kidney Donor Exchange (dcmpe11) */
/* PRODUCT: OR */
/* SYSTEM: ALL */
/* KEYS: OR */
/* PROCS: OPTMODEL, OPTNETWORK */
/* DATA: */
/* */
/* SUPPORT: UPDATE: */
/* REF: */
/* MISC: Example 11 from the Decomposition Algorithm */
/* chapter of Mathematical Programming. */
/* */
/***************************************************************/
/* create random graph on n nodes with arc probability p
and uniform(0,1) weight */
%let n = 100;
%let p = 0.02;
data ArcData;
call streaminit(1);
do i = 0 to &n - 1;
do j = 0 to &n - 1;
if i eq j then continue;
else if rand('UNIFORM') < &p then do;
weight = rand('UNIFORM');
output;
end;
end;
end;
run;
%let max_length = 10;
proc optmodel;
set <num,num> ARCS;
num weight {ARCS};
read data ArcData into ARCS=[i j] weight;
print weight;
set NODES = union {<i,j> in ARCS} {i,j};
set MATCHINGS = 1..card(NODES)/2;
/* UseNode[i,m] = 1 if node i is used in matching m, 0 otherwise */
var UseNode {NODES, MATCHINGS} binary;
/* UseArc[i,j,m] = 1 if arc (i,j) is used in matching m, 0 otherwise */
var UseArc {ARCS, MATCHINGS} binary;
/* maximize total weight of arcs used */
max TotalWeight
= sum {<i,j> in ARCS, m in MATCHINGS} weight[i,j] * UseArc[i,j,m];
/* each node appears in at most one matching */
/* rewrite as set partitioning (so decomp uses identical blocks)
sum{} x <= 1 => sum{} x + s = 1, s >= 0 with no associated cost */
var Slack {NODES} binary;
con Packing {i in NODES}:
sum {m in MATCHINGS} UseNode[i,m] + Slack[i] = 1;
/* at most one recipient for each donor */
con Donate {i in NODES, m in MATCHINGS}:
sum {<(i),j> in ARCS} UseArc[i,j,m] = UseNode[i,m];
/* at most one donor for each recipient */
con Receive {j in NODES, m in MATCHINGS}:
sum {<i,(j)> in ARCS} UseArc[i,j,m] = UseNode[j,m];
/* exclude long matchings */
con Cardinality {m in MATCHINGS}:
sum {<i,j> in ARCS} UseArc[i,j,m] <= &max_length;
/* automatically decompose using METHOD=SET */
solve with milp / presolver=basic decomp=(method=set);
/* save solution to a data set */
create data Solution from
[m i j]={m in MATCHINGS, <i,j> in ARCS: UseArc[i,j,m].sol > 0.5}
weight[i,j];
quit;
data Solution;
set Solution;
run;
proc optnet
direction = directed
links = Solution;
links_var
from = i
to = j;
cycle
mode = all_cycles
out = Cycles;
run;
data Cycles; set Cycles; run;
proc sort data=Cycles; by cycle order; run;
proc print data=Cycles noobs label;
by cycle;
run;