$offsymlist offsymxref offuellist offuelxref
option limrow=0, limcol=0;

sets
 I / 1 * 3 /,
 J / 1 * 4 /;

table A(I,J)
	1	2	3	4
1	1		1	
2	1	1	1	1
3		1			;

parameters
c(J) /
1	4,
2	1,
3	3,
4	1
/,
b(I) /
1	2,
2	4,
3	2
/;

free variable z;

positive variables
x(J);

equations
obj,
f(I);

obj ..	z =e= sum(J, c(J)*x(J) );

f(I) ..	sum(J, A(I,J)*x(J) ) =l= b(I);

model testlp / obj, f /;

solve testlp using lp maximizing z;

* here, we have all constraints active

* however, since the dual f.m("3") is 0, we know b("3") can be
* increased without bound, giving a "range" of [0 .. inf] for this solution

b("3") = b("3") + 8;

solve testlp using lp maximizing z;

* this demonstrates that the solution need not change from that above.

* in a LP of this form, we can decrease the c(j)'s for which x(j) = 0

c(J)$(x.l(J) eq x.lo(J)) = c(J) - 5;

solve testlp using lp maximizing z;

* now, note that since the duals (reduced costs)
* for c("3") and c("4") are negative,
* we know we can increase these c's by the negative of their reduced
* cost and decrease them without bound,
* while keeping the same solution

* this will introduce some degeneracy!!
c(J) = c(J) - x.m(J);

solve testlp using lp maximizing z;


************************************************************
* The complementarity slackness conditions for the LP can
* be written down and solved for explicitly using an LCP
************************************************************


positive variables
y(I);

equations
delx(J),
dely(I);

delx(J) ..	sum(I, y(I)*A(I,J) ) =g= c(J);
dely(I) ..	b(I) =g= sum(J, A(I,J)*x(J) );

model testlcp / delx.x, dely.y /;

* avoid repeating work already done
y.l(I) = f.m(I);
solve testlcp using mcp;