pulp constraint problem

927 views
Skip to first unread message

Christopher Cardinal

unread,
Aug 2, 2010, 7:47:13 PM8/2/10
to pulp-or...@googlegroups.com
Hi. I'm a newbie to PuLP and OR, in general (geography student). I am trying to implement a binary linear integer programming model to optimize a collection of parcels, each of which has a conservation value(CV), cost(C) and an environmental benefit (EB).

CV = a real number representing the conservation value of a parcel, e.g. 3.4
C = the cost of a parcel, e.g. $45,000
EB = the env. benefits of the parcel, a real number, e.g. 6.5
B = the annual budget of the land conservation org.

In the first case, I simply Maximize CV wrt B:

X = LpVariable.dict('X',UID,0,1,LpInteger)
prob = LpProblem('Parcel Optimization -- Max CV', LpMaximize)
prob += lpSum([(X[uidcv[0]] * uidcv[1]) for uidcv in UIDCV])
prob += lpSum([X[uidc[0]] * uidc[1] for uidc in UIDC]) <= dLPLimits['BUDGET']

UID is a list of parcel IDs: [56,499,500 ...]
UIDCV is a list of tuples, parcel ID and CV: [(56, 0.273553404245), (499, 0.974708856126), (500, 0.139934542203) ...]
UIDC is also a list of tuples, parcel ID and C: [(56, 354392.468123), (499, 336757.595564), (500, 318850.964382) ...]
This works for many different values of B, always returning a Max CV value spending as much of B as possible. I've confirmed results in Open Office Solver.

Example result for a $5,000,000 budget:
Total CV from Problem Objective = 308.980117485
Sum of cost = 4,998,813.48679
Result of problem = 1

My problem concerns adding another constraint. I want to set a minimum value for EB:

X = LpVariable.dict('X',UID,0,1,LpInteger)
prob = LpProblem('Parcel Optimization -- With EB Limits', LpMaximize)
prob += lpSum([(X[uidcv[0]] * uidcv[1]) for uidcv in UIDCV])
prob += lpSum([X[uidc[0]] * uidc[1] for uidc in UIDC]) <= dLPLimits['BUDGET']
prob += lpSum([X[uideb[0]] * uideb[1] for uideb in UIDEB]) >= dLPLimits['EBMIN']

UIDEB is a list of tuples, parcel ID and EB: [(56, 7.44573307698), (499, 26.5301833495), (500, 5.56020769487) ...]

The results returned, however, seem to ignore the budget constraint:
Total CV from Problem Objective = 664.317056443
Sum of cost = 65,035,395.1858 (too big!)
Result of problem = 1

If I change the sense of the inequality of the EB constraint, the results look good again:
Total CV from Problem Objective = 309.105734267
Sum of cost = 4,998,873.95198
Result of problem = 1

Am I missing something here?

Thanks,
C


Stuart Mitchell

unread,
Aug 3, 2010, 12:51:34 AM8/3/10
to pulp-or...@googlegroups.com
On Tue, 2010-08-03 at 11:47 +1200, Christopher Cardinal wrote:
> Hi. I'm a newbie to PuLP and OR, in general (geography student). I am
> trying to implement a binary linear integer programming model to
> optimize a collection of parcels, each of which has a conservation
> value(CV), cost(C) and an environmental benefit (EB).

Cool the more people that use it the better


> UID is a list of parcel IDs: [56,499,500 ...]
> UIDCV is a list of tuples, parcel ID and CV: [(56, 0.273553404245),
> (499, 0.974708856126), (500, 0.139934542203) ...]
> UIDC is also a list of tuples, parcel ID and C: [(56, 354392.468123),
> (499, 336757.595564), (500, 318850.964382) ...]

Please try (so that my eyes don't hurt)
>>> UIDCV = dict(UIDCV) # turn it into a dictionary
>>> UIDC = dict(UIDC)

Then the following is

> In the first case, I simply Maximize CV wrt B:
>
X = LpVariable.dict('X',UID,0,1,LpInteger)
prob = LpProblem('Parcel Optimization -- Max CV', LpMaximize)

prob += lpSum([(X[id] * uidcv[id]) for id in UID])
prob += lpSum([X[id] * uidc[id] for id in UID]) <= dLPLimits['BUDGET']


> This works for many different values of B, always returning a Max CV
> value spending as much of B as possible. I've confirmed results in
> Open Office Solver.
>
> Example result for a $5,000,000 budget:
> Total CV from Problem Objective = 308.980117485
> Sum of cost = 4,998,813.48679
> Result of problem = 1
>
> My problem concerns adding another constraint. I want to set a minimum
> value for EB:
>
> X = LpVariable.dict('X',UID,0,1,LpInteger)
> prob = LpProblem('Parcel Optimization -- With EB Limits', LpMaximize)
> prob += lpSum([(X[uidcv[0]] * uidcv[1]) for uidcv in UIDCV])
> prob += lpSum([X[uidc[0]] * uidc[1] for uidc in UIDC]) <=
> dLPLimits['BUDGET']
> prob += lpSum([X[uideb[0]] * uideb[1] for uideb in UIDEB]) >=
> dLPLimits['EBMIN']
>
> UIDEB is a list of tuples, parcel ID and EB: [(56, 7.44573307698),
> (499, 26.5301833495), (500, 5.56020769487) ...]

Please check the status of your problem after the solve
>>> print "Status:", LpStatus[prob.status]
I have a feeling that your problem may be infeasible i.e. there is not a
solution that satisfies all of the constraints.
Common practice in the OR world is to return a status though I always
think that it would be more user friendly to raise an exception.


Vinaka
Stu

--
___________________________________
Dr Stuart Mitchell
Research Fellow
Light Metals Research Centre (LMRC)
University of Auckland
Private Bag 92019
Auckland
New Zealand

Ph (wk) +64 9 3737599 ext 84867
(ddi) +64 9 9234867
(fax) +64 9 3737925
(mb) +64 21 441331
___________________________________


Christopher Cardinal

unread,
Aug 3, 2010, 1:41:25 PM8/3/10
to pulp-or...@googlegroups.com
Pardon the extra email, but I wanted to add something I found. If I hold budget constant at $5million, and gradually increase the EB minimum, i.e. move towards infeasibility, I see some odd results. First, is it normal that running the LP process repeatedly would produce different results?  I see different values for the objective function when the same process is run repeatedly. Other than that, results looked good until EB = 9000.

EB: 1000
Problem Objective = 308.988173812  --- (value returned by value(prob.objective))
Sum of CV = 308.988173812 --- (my calc of the sum of CV values)
Sum of cost = 4999645.63167 --- (my calc of sum of costs)
Sum of EB = 6788.08587788 --- (my calc of sum of EB)
Result of problem = Optimal

Problem Objective = 309.040460187
Sum of CV = 309.040460187
Sum of cost = 4999220.35968
Sum of EB = 6791.26154025
Result of problem = Optimal

Problem Objective = 309.105734267
Sum of CV = 309.105734267
Sum of cost = 4998873.95198
Sum of EB = 6699.56096171
Result of problem = Optimal

EB: 5000
Problem Objective = 309.040460187
Sum of CV = 309.040460187
Sum of cost = 4999220.35968
Sum of EB = 6791.26154025
Result of problem = Optimal

Problem Objective = 309.105734267
Sum of CV = 309.105734267
Sum of cost = 4998873.95198
Sum of EB = 6699.56096171
Result of problem = Optimal

EB: 8000
Problem Objective = 299.351831582
Sum of CV = 299.351831582
Sum of cost = 4998428.36565
Sum of EB = 8020.24120641
Result of problem = Optimal

Problem Objective = 299.546541364
Sum of CV = 299.546541364
Sum of cost = 4999121.10279
Sum of EB = 8013.31549288
Result of problem = Optimal

EB: 9000
Problem Objective = 100.43116751
Sum of CV = 293.128137386
Sum of cost = 5557999.96889
Sum of EB = 9395.10977153
Result of problem = Optimal

With EB = 9,000, the solution should be infeasible: the cost is greater than the budget. I have checked the value of the costs a few times to be sure I haven't made a calculation error. However, when I scan the varValues for each solution variable, two values are odd:
-416.285769049
453.505096974

In the previous iterations, all solution variables are zero or one, aside from a couple of variables with a value of 0.999999999998.

C


On Tue, Aug 3, 2010 at 11:57 AM, Christopher Cardinal <hobo...@gmail.com> wrote:
Thanks for the suggestions. I've converted everything to a dict -- much clearer. I have been expecting an infeasibility result. However, in every case, the solution status is Optimal.

If I set a minimum constraint for EB (70,000) and maximum for Budget ($1million):

    X = LpVariable.dict('X',UID,0,1,LpInteger)
    prob = LpProblem('Parcel Optimization -- Min EB', LpMaximize)
    prob += lpSum([(X[uid] * UIDCV[uid]) for uid in UID])
    prob += lpSum([X[uid] * UIDC[uid]  for uid in UID]) <= dLPLimits['BUDGET']
    prob += lpSum([X[uid] * UIDEB[uid]  for uid in UID]) >= dLPLimits['EBMIN']

>>> prob.solve()

1

>>> print "Status:", LpStatus[prob.status]
Status: Optimal

and a weird (to me), negative value for the objective function:
>>> value(prob.objective)
-125322.07153962553

But 26 parcels are chosen, based on varValue:
    i = 0
    for v in prob.variables():
        if v.varValue > 0.0:
            i += 1
    print i
    26

- 22 solution variables equal 1.0, and the other four equal 0.99999999999999989
The amount of $$ spent is well over the budget.

If I change the sense of the EB constraint:
    prob += lpSum([X[uid] * UIDEB[uid]  for uid in UID]) <= landcons.dLPLimits['EBMIN']

>>> prob.solve()

1

>>> print "Status:", LpStatus[prob.status]
Status: Optimal

>>> value(prob.objective)
105.572096431254

and the process spends $999,654 of a $1million budget.

It seems like the negative value for the value of the objective function indicates infeasibility. If I increase the budget and decrease minimum EB, $20million and 10,000, respectively:
>>> prob.solve()

1

>>> print "Status:", LpStatus[prob.status]
Status: Optimal
>>> value(prob.objective)
789.10768481136404

and $19,998,804 is spent to purchase the parcels.This is an expected result.

Any more thoughts? Thanks a lot!

C




--
You received this message because you are subscribed to the Google Groups "pulp-or-discuss" group.
To post to this group, send email to pulp-or...@googlegroups.com.
To unsubscribe from this group, send email to pulp-or-discu...@googlegroups.com.
For more options, visit this group at http://groups.google.com/group/pulp-or-discuss?hl=en.



Christopher Cardinal

unread,
Aug 3, 2010, 11:57:23 AM8/3/10
to pulp-or...@googlegroups.com
Thanks for the suggestions. I've converted everything to a dict -- much clearer. I have been expecting an infeasibility result. However, in every case, the solution status is Optimal.

If I set a minimum constraint for EB (70,000) and maximum for Budget ($1million):
    X = LpVariable.dict('X',UID,0,1,LpInteger)
    prob = LpProblem('Parcel Optimization -- Min EB', LpMaximize)
    prob += lpSum([(X[uid] * UIDCV[uid]) for uid in UID])
    prob += lpSum([X[uid] * UIDC[uid]  for uid in UID]) <= dLPLimits['BUDGET']
    prob += lpSum([X[uid] * UIDEB[uid]  for uid in UID]) >= dLPLimits['EBMIN']

>>> prob.solve()
1

>>> print "Status:", LpStatus[prob.status]
Status: Optimal

and a weird (to me), negative value for the objective function:
>>> value(prob.objective)
-125322.07153962553

But 26 parcels are chosen, based on varValue:
    i = 0
    for v in prob.variables():
        if v.varValue > 0.0:
            i += 1
    print i
    26

- 22 solution variables equal 1.0, and the other four equal 0.99999999999999989
The amount of $$ spent is well over the budget.

If I change the sense of the EB constraint:
    prob += lpSum([X[uid] * UIDEB[uid]  for uid in UID]) <= landcons.dLPLimits['EBMIN']

>>> prob.solve()
1
>>> print "Status:", LpStatus[prob.status]
Status: Optimal

>>> value(prob.objective)
105.572096431254

and the process spends $999,654 of a $1million budget.

It seems like the negative value for the value of the objective function indicates infeasibility. If I increase the budget and decrease minimum EB, $20million and 10,000, respectively:
>>> prob.solve()
1
>>> print "Status:", LpStatus[prob.status]
Status: Optimal
>>> value(prob.objective)
789.10768481136404

and $19,998,804 is spent to purchase the parcels.This is an expected result.

Any more thoughts? Thanks a lot!

C

On Tue, Aug 3, 2010 at 12:51 AM, Stuart Mitchell <s.mit...@auckland.ac.nz> wrote:

Stuart Mitchell

unread,
Aug 3, 2010, 5:30:14 PM8/3/10
to pulp-or...@googlegroups.com
I think I have your answer,

There is an error in the CoinMP.dll solver that does not report integer
infeasible solutions correctly.

If you look at the output (with msg=True) you will probably see output
that shows this.

To fix the problem you can use another solver glpk or cbc command line,
check the values of your X variables (strange according to your comment
about their values you should be fine).
Alternatively, use the LpProblem.writeLP() method to actually write out
the lp file to check that your problem is correct.
Or, send me your complete code for a quick look (off list if it is
confidential) and I'll see what I can do.

Stu

> To unsubscribe from this group, send email to pulp-or-discuss
> +unsub...@googlegroups.com.


> For more options, visit this group at
> http://groups.google.com/group/pulp-or-discuss?hl=en.

Christopher Cardinal

unread,
Aug 7, 2010, 8:11:26 PM8/7/10
to pulp-or...@googlegroups.com
Thanks. I haven't had a chance to play with other solvers, but hopefully in the near future. I'd like to ask another question. The question is whether or not linear programming is suitable for solving the following problem.


UID is a list of parcel IDs
UIDCV is a dict linking UID to a conservation value of the parcel
UIDC is a dict linking UID to its cost

I can maximize CV for a given BUDGET with
   X = LpVariable.dict('X',UID,0,1,LpInteger)
   prob = LpProblem('Parcel Optimization -- Max CV',LpMaximize)
   prob += lpSum([(X[uid] * UIDCV[uid]) for uid in UID])
   prob += lpSum([X[uid] * UIDC[uid]  for uid in UID]) <= dLPLimits['BUDGET']

Now, suppose I want the CV of each parcel to increase if it's neighbor is purchased as well. This would be determined by the program, if possible.

UIDSIZE is a dict linking each parcel to its acreage

So, for instance, if parcel 1 has neighbors 2, 3, and 4, then the CV of parcel 1 =
X_1 * UIDCV[1] + X_2 * UIDSIZE[2] + X_3 * UIDSIZE[3] + X_4 * UIDSIZE[4]
etc., etc.

If you catch my drift, does this seem solvable by BLP? As far as formulating the problem, that's another story.

C

Stuart Mitchell

unread,
Aug 8, 2010, 5:52:47 PM8/8/10
to pulp-or...@googlegroups.com
On Sun, 2010-08-08 at 12:11 +1200, Christopher Cardinal wrote:

> Now, suppose I want the CV of each parcel to increase if it's neighbor
> is purchased as well. This would be determined by the program, if
> possible.

Ok I think you will find some papers on this subject try a google
scholar search for
conservation value integer programming


>
> UIDSIZE is a dict linking each parcel to its acreage
>
> So, for instance, if parcel 1 has neighbors 2, 3, and 4, then the CV
> of parcel 1 =
> X_1 * UIDCV[1] + X_2 * UIDSIZE[2] + X_3 * UIDSIZE[3] + X_4 *
> UIDSIZE[4]
> etc., etc.

Hmm not really able to model it as a linear program that way as it is
implicitly
X_1 * UIDCV[1] + X_1*X_2 * UIDSIZE[2] + X_1*X_3 * UIDSIZE[3] + X_1*X_4 *
UIDSIZE[4]

but I would recommend as a simple exercise (and if the groupings are not
that big) to model the problem as a group of larger perhaps overlapping
groups of parcels. Then you can use a set packing formulation. You are
free to use whatever method you like for pricing the groups.

Things will get a lot more complicated if you want to automatically
generate the groups and guarantee an optimal solution. But for the
moment see if that approximation works.

Reply all
Reply to author
Forward
0 new messages