*==>freebert.gms $offsymxref offsymlist offuellist offuelxref $TITLE Bertsekas problem as a 15-variable MCP $ontext This problem is based on a traffic assigment problem, found in: D.P. Bertsekas and E.M. Gafni, "Projection Methods for Variational Inequalities with Application to the Traffic Assignment Problem", Mathematical Programming Study 17(1982), 139-159. Just so you know, x(i) is the flow (from node(i) to node(i++3)) on the (longer) outside loop, z(i) the flow on the inside loop, and y the dual variables. The source problem here is a VI, where the constraints consist of non-negative bounds and x(I) + z(I) = demand(I). Without much fiddling, we can write the problem as a MCP. We can write the problem in (at least) 2 different ways. In the first, y is free, and the demand constraints are satisfied exactly. In the second, we change the demand constraint to an inequality, where extra flow is allowed. This doesn't occur in solutions, and we can write the problem as an NCP. $offtext sets run / run1 * run7 /, I nodes in the network / 1 * 5 /, types types of arcs / ihy, ion, ioff, iby, ohy, oon, ooff, oby / ; parameters includerun (run) / run1 01 run2 01 run3 01 run4 01 run5 01 run6 01 run7 01 /, demand(I), c (I, types); scalar gamma; $include 'bert_gaf.dat' * demand and gamma data included in bertsekas and gafni files parameters ipty(I,run), iptx(I,run), iptz(I,run); table iptx(I,run) run1 run2 run3 run4 run5 run6 1 2 3 4 5 ; table iptz(I,run) run1 run2 run3 run4 run5 run6 1 2 3 4 5 ; iptz (I,"run1") = demand (I) - iptx (I,"run1"); iptz (I,"run2") = demand (I) - iptx (I,"run2"); iptz (I,"run3") = demand (I) - iptx (I,"run3"); ipty (I,"run2") = 4000; ipty (I,"run5") = 4000; ipty (I,"run3") = 1; ipty (I,"run6") = 1; iptx (I,"run7") = 0.5; ipty (I,"run7") = 0.5; iptz (I,"run7") = 0.5; positive variable x(I) outside flow, z(I) inside flow; free variable y(I) dual vars; equations delay_o(I), delay_i(I), dem_cons(I); delay_o(I) .. c(I,"oon")*(1 + x(I) + power(x(I),2)) + gamma * (1 + x(I++3) + x(I++4) + power(x(I++3)+x(I++4),2)) + 10 * c(I,"ohy")*(1 + x(I) + x(I++3) + x(I++4) + power(x(I) + x(I++3) + x(I++4),2)) + 2 * gamma * (1 + x(I++3) + power(x(I++3),2)) + c(I++1,"oby") * (1 + x(I) + x(I++4) + power(x(I) + x(I++4),2)) + 10 * c(I++1,"ohy") * (1 + x(I) + x(I++1) + x(I++4) + power(x(I) + x(I++1) + x(I++4), 2)) + 2 * gamma * (1 + x(I++4) + power(x(I++4),2)) + c(I++2,"oby") * (1 + x(I) + x(I++1) + power(x(I)+x(I++1),2)) + 10 * c(I++2,"ohy") * (1 + x(I) + x(I++1) + x(I++2) + power(x(I) + x(I++1) + x(I++2), 2)) + 2 * gamma * (1 + x(I) + power(x(I),2)) + c(I++3,"ooff") * (1 + x(I) + power(x(I),2)) =g= y(I); delay_i(I) .. c(I,"ion")*(1 + z(I) + power(z(I),2)) + gamma * (1 + z(I++1) + power(z(I++1),2)) + 10 * c(I,"ihy")*(1 + z(I) + z(I++1) + power(z(I) + z(I++1),2)) + 2 * gamma * (1 + z(I++1) + power(z(I++1),2)) + c(I--1,"iby")*(1 + z(I) + power(z(I),2)) + 10 * c(I--1,"ihy") * (1 + z(I) + z(I++4) + power(z(I)+z(I++4), 2)) + 2 * gamma * (1 + z(I) + power(z(I),2)) + c(I--2,"ioff") * (1 + z(I) + power(z(I),2)) =g= y(I); dem_cons (I) .. x(I) + z(I) =e= demand(I); MODEL freebert / delay_o.x, delay_i.z, dem_cons.y /; option limrow = 0; option limcol = 0; option iterlim = 1000; option reslim = 120; loop (run, if (includerun(run), x.l(I) = iptx(I,run); z.l(I) = iptz(I,run); y.l(I) = ipty(I,run); solve freebert using mcp; ));