Kidney Donor Exchange (dcmpe10)
/***************************************************************/
/* */
/* S A S S A M P L E L I B R A R Y */
/* */
/* NAME: dcmpe10 */
/* TITLE: Kidney Donor Exchange (dcmpe10) */
/* PRODUCT: OR */
/* SYSTEM: ALL */
/* KEYS: OR */
/* PROCS: OPTMODEL, OPTNET */
/* DATA: */
/* */
/* SUPPORT: UPDATE: */
/* REF: */
/* MISC: Example 10 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;
do i = 0 to &n - 1;
do j = 0 to &n - 1;
if i eq j then continue;
else if ranuni(1) < &p then do;
weight = ranuni(2);
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 node_packing {i in NODES}:
sum {m in MATCHINGS} UseNode[i,m] + Slack[i] = 1;
/* at most one donee 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 donee */
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;
/* decompose by matching (aggregate formulation) */
for {i in NODES, m in MATCHINGS} donate[i,m].block = m;
for {j in NODES, m in MATCHINGS} receive[j,m].block = m;
for {m in MATCHINGS} cardinality[m].block = m;
solve with milp / presolver=basic decomp=();
/* 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;
proc optnet
direction = directed
data_links = Solution;
data_links_var
from = i
to = j;
cycle
mode = all_cycles
out = Cycles;
run;
proc print data=Cycles noobs label;
by cycle;
run;