I'm trying to set up MIP problem using LpSolve and I'm at a loss on how to set it up.
We are trying to allocate resources to various investments, to maximize revenue, subject to a budget constraint. Each business is located in a particular state, and We also have to spread our allocations across 5 states (this is the one that's causing me issues). My current linear equations are this:
Objective function: maximize following where Xi represents the number of dollars allocated to ith business, and Ri represents revenue generated by the ith business
X1 * R1 +....+ X35 * R35
Subject to the constraint:
X1 * C1 +....+ X35 * C35 <= 1,000,000
Where Ci represents the costs of the ith business. Setting this objective function and constraint matrix is very easy. However each business is located in some state, and what I need to do is setup a constraint that we must allocate to business in at least 5 states.
What I initially considered was a set of constraints where each state was represented as a row: so if this row represents the state of California for instance:
X1 * 1 +....+ X2 * 1 + X3 *0 + X35 * 1 <= Some Constraint
Then each business in California is assigned a value of 1, then this would give me the sum of $s (from the Xi assigned to each business) allocated to California. I could do this for all states, but that doesn't quite get me what I want as I thought about it mroe. I could constrain each state using this methodology, but what I want is something that lets me assign Unique_States >5
But this doesn't quite get me what I want. What I really want is something that gives me a sum of total number of states "used".
Is this even possible?
I've spent a while trying to rack my brain on how to do this but I can't quite figure it out.
There's one other questions like this on SO, but I couldn't follow any of the answers: LpSolve R conditional constraint
This question is awfully similar albeit a totally different context, but I couldn't follow it and I don't have enough reputation to ask any questions.
If anyone could give any sort of guidance
Here is a skeleton of your model in pyomo. This is basically a "Knapsack Problem" with a twist (the min number of states reqt.).
import pyomo.environ as pyo
from collections import defaultdict
# some data... this could be from dictionary (as below) or file read or ....
# name state C R
investments = { 'Burger Joint': ('CA', 1.2, 3.1),
'car wash': ('CA', 1.1, 2.6),
'gem store': ('NV', 1.0, 2.5),
'warehouse': ('UT', 1.2, 3.6)}
budget = 10
min_investment = 1
min_states_to_use = 2
states = {v[0] for v in investments.values()}
investments_by_state = defaultdict(list)
for k, v in investments.items():
investments_by_state[v[0]].append(k)
# set up the model
m = pyo.ConcreteModel()
# SETS
m.I = pyo.Set(initialize=investments.keys())
m.S = pyo.Set(initialize=list(states))
# the below is an "indexed set" which is indexed by STATE and contains the
# subset of I within that state
m.IS = pyo.Set(m.S, within=m.I, initialize=investments_by_state)
# PARAMETERS
m.cost = pyo.Param(m.I, initialize={k:v[1] for k, v in investments.items()})
m.rev = pyo.Param(m.I, initialize={k:v[2] for k, v in investments.items()})
# VARIABLES
m.X = pyo.Var(m.I, domain=pyo.NonNegativeReals) # amount invested in i
m.Y = pyo.Var(m.S, domain=pyo.Binary) # 1 if state y has an investment
# OBJECTIVE
m.obj = pyo.Objective(expr=sum(m.X[i]*m.rev[i] for i in m.I), sense=pyo.maximize)
# CONSTRAINTS
# stay in budget...
m.C1 = pyo.Constraint(expr=sum(m.X[i]*m.cost[i] for i in m.I) <= budget)
# connect indicator Y to X...
def state_invest(m, s):
return m.Y[s] * min_investment <= sum(m.X[i] for i in m.IS[s])
m.C2 = pyo.Constraint(m.S, rule=state_invest)
# use needed number of states
m.C3 = pyo.Constraint(expr=sum(m.Y[s] for s in m.S) >= min_states_to_use)
m.pprint() # <-- show the whole model
# instantiate a solver and solve...
solver = pyo.SolverFactory('cbc')
result = solver.solve(m)
print(result)
m.X.display() # <-- "display" substitutes in the variable values
m.Y.display()
3 Set Declarations
I : Size=1, Index=None, Ordered=Insertion
Key : Dimen : Domain : Size : Members
None : 1 : Any : 4 : {'Burger Joint', 'car wash', 'gem store', 'warehouse'}
IS : Size=3, Index=S, Ordered=Insertion
Key : Dimen : Domain : Size : Members
CA : 1 : I : 2 : {'Burger Joint', 'car wash'}
NV : 1 : I : 1 : {'gem store',}
UT : 1 : I : 1 : {'warehouse',}
S : Size=1, Index=None, Ordered=Insertion
Key : Dimen : Domain : Size : Members
None : 1 : Any : 3 : {'CA', 'NV', 'UT'}
2 Param Declarations
cost : Size=4, Index=I, Domain=Any, Default=None, Mutable=False
Key : Value
Burger Joint : 1.2
car wash : 1.1
gem store : 1.0
warehouse : 1.2
rev : Size=4, Index=I, Domain=Any, Default=None, Mutable=False
Key : Value
Burger Joint : 3.1
car wash : 2.6
gem store : 2.5
warehouse : 3.6
2 Var Declarations
X : Size=4, Index=I
Key : Lower : Value : Upper : Fixed : Stale : Domain
Burger Joint : 0 : None : None : False : True : NonNegativeReals
car wash : 0 : None : None : False : True : NonNegativeReals
gem store : 0 : None : None : False : True : NonNegativeReals
warehouse : 0 : None : None : False : True : NonNegativeReals
Y : Size=3, Index=S
Key : Lower : Value : Upper : Fixed : Stale : Domain
CA : 0 : None : 1 : False : True : Binary
NV : 0 : None : 1 : False : True : Binary
UT : 0 : None : 1 : False : True : Binary
1 Objective Declarations
obj : Size=1, Index=None, Active=True
Key : Active : Sense : Expression
None : True : maximize : 3.1*X[Burger Joint] + 2.6*X[car wash] + 2.5*X[gem store] + 3.6*X[warehouse]
3 Constraint Declarations
C1 : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : -Inf : 1.2*X[Burger Joint] + 1.1*X[car wash] + X[gem store] + 1.2*X[warehouse] : 10.0 : True
C2 : Size=3, Index=S, Active=True
Key : Lower : Body : Upper : Active
CA : -Inf : Y[CA] - (X[Burger Joint] + X[car wash]) : 0.0 : True
NV : -Inf : Y[NV] - X[gem store] : 0.0 : True
UT : -Inf : Y[UT] - X[warehouse] : 0.0 : True
C3 : Size=1, Index=None, Active=True
Key : Lower : Body : Upper : Active
None : 2.0 : Y[CA] + Y[NV] + Y[UT] : +Inf : True
11 Declarations: I S IS cost rev X Y obj C1 C2 C3
Problem:
- Name: unknown
Lower bound: -29.5
Upper bound: -29.5
Number of objectives: 1
Number of constraints: 5
Number of variables: 7
Number of binary variables: 3
Number of integer variables: 3
Number of nonzeros: 4
Sense: maximize
Solver:
- Status: ok
User time: -1.0
System time: 0.0
Wallclock time: 0.0
Termination condition: optimal
Termination message: Model was solved to optimality (subject to tolerances), and an optimal solution is available.
Statistics:
Branch and bound:
Number of bounded subproblems: 0
Number of created subproblems: 0
Black box:
Number of iterations: 0
Error rc: 0
Time: 0.014362096786499023
Solution:
- number of solutions: 0
number of solutions displayed: 0
X : Size=4, Index=I
Key : Lower : Value : Upper : Fixed : Stale : Domain
Burger Joint : 0 : 1.0 : None : False : False : NonNegativeReals
car wash : 0 : 0.0 : None : False : False : NonNegativeReals
gem store : 0 : 0.0 : None : False : False : NonNegativeReals
warehouse : 0 : 7.3333333 : None : False : False : NonNegativeReals
Y : Size=3, Index=S
Key : Lower : Value : Upper : Fixed : Stale : Domain
CA : 0 : 1.0 : 1 : False : False : Binary
NV : 0 : 0.0 : 1 : False : False : Binary
UT : 0 : 1.0 : 1 : False : False : Binary
[Finished in 696ms]