Skip to main content

Model is infeasible or unbounded

Answered
William Praat
  • Gurobi-versary
  • First Comment
  • First Question

Dear Gurobi community. For an educational assignment I made a mathematical model and wrote it in python and am trying to use Gurobi as a solver.

Although I think I have implemented the decision variables, constraints, and objective function correctly, Gurobi seems to disagree with that as it states that the "Model is infeasible or unbounded".

I've tried many things and try to figure out the issue by reading the documentation but with no success so far. Note that I am a beginner, excuse me for the foolish mistakes that might be in my code.

T = 40 #10 years, 4 quarters per year
M = [[0, 0],[1, 4],[2, 8],[3,12],[4,16],[5,20],[6,24],[7, 28],[8, 32],[9, 36]] #10 technology levels: each year a new system is introduced

StartDemand = 0.5 #initial demand
GrowthRate = 0.1 #Growth per period (exponential)
randomness = 0.25 #0-1, percentage of randomness to demand growth


S = 2000 #setup costs (static for now)
A = 10000 # Acquesitioin costs (static for now)
O = 200 # operational costs (static for now)
L = 500 # salvage profit (static for now)

def generatePurchaseDecisionVars(T, M):
"""
In this simple version of the model the age of the machines is not determined
So only the purchasing and salvaging decision variables have to be made

A decision var is made for each tech level and period if the tech level is available in period t
"""
decisionVars = []

for m in M:
t = 0
while t < T:
if t >= m[1]:
tempTuple = (m[0],t)
decisionVars.append(tempTuple)
t = t + 1
return decisionVars

def generateSalvageDecisionVars(T,M):
"""
similar to generatePurchaseDecisionVars but I kept it in a different function as this might change in the future
"""
decisionVars = []

for m in M:
t = 0
while t < T:
if t >= m[1]:
tempTuple = (m[0],t)
decisionVars.append(tempTuple)
t = t + 1
return decisionVars

def generateOperatingSystemsCountVars(T, M):
"""
similar to generatePurchaseDecisionVars but I kept it in a different function as this might change in the future
"""
decisionVars = []
for m in M:
t = 0
techlevel = []
while t < T:
techlevel.append(0)
t = t + 1
decisionVars.append(techlevel)
return decisionVars

def generateDeteriorationOfSystems(T, M):
Umt = []
for m in M:
techlevel = []
for t in range(T):
deter = 0
if t >= m[1]:
# per period deteriroration increases by 0.5%, each tech level m deteriration is increasing less fast
deter = (0.995+(1-0.995)*0.1*m[0])**(t - m[1])
techlevel.append(deter)
Umt.append(techlevel)
return Umt

def getDemandGrowth(T, startDemand, growthRate, randomness):
#Generating the exponential growth
DeamandRange = range(0, T)
demandGrowth = [StartDemand*exp(GrowthRate*t) for t in DeamandRange]

randomizedDemand = []
#print(demandGrowth)
for d in demandGrowth:
randomRange = d*randomness
randomizedDemand.append(d + N.random.uniform(-randomRange, randomRange))

return randomizedDemand

purchaseDecisionVars = generatePurchaseDecisionVars(T, M)
salvageDecisionVars = generateSalvageDecisionVars(T, M)

model = Model("HydrogenSystemOptimization")

"""
for each t and m a X (purchase system m in period t) and Y (Salvage system in period t) decision variables have been made
All are integer values and can be equal or greater than 0
"""
X = model.addVars(purchaseDecisionVars, lb=0, vtype=GRB.INTEGER)
Y = model.addVars(salvageDecisionVars, lb=0, vtype=GRB.INTEGER)

"""
for each period t a K variable is created. K(t) should be 1 if a system is purchased in period t and 0 if this is not the case.
"""
K = []
for t in range(T):
K.append(0)
model.update()
"""
### Equation 1
if in period t 1 or more systems are purchased K(t) is forced to be 1, otherwise the minimize objective will force
it to be 0.
"""
for t in range(T):
model.addConstr((K[t]*10000) >= (quicksum(X[m[0],t] for m in M if t >= m[1])))

"""
### Equation 2 + 3
Number of operational machines F(m,t) is calculated here which is later used in the next constraint and in the objective function. non-negativity is added here as well.
"""
operatingSystemsVars = generateOperatingSystemsCountVars(T, M)
for t in range(T):
for m in M:
if t >= m[1]:
model.addConstr(operatingSystemsVars[m[0]][t] == quicksum(X[m[0],i] - Y[m[0],i] for i in range(t) if i >= m[1]), name="operating")

D = getDemandGrowth(T, StartDemand, GrowthRate, randomness)
U = generateDeteriorationOfSystems(T, M)


"""
### Equation 4
The total capacity in period t must be greater or equal to the demand in that period.
"""
for t in range(T):
model.addConstr((quicksum(U[m[0]][t] * operatingSystemsVars[m[0]-1][t] for m in M if t >= m[1])) >= D[t])

"""
### Equation 5, objective funciton
Minimize: Set-up costs*K(t) + Acquesition cost * number of machines purchased(m,t) + #of operating machines * operating costs - Salvage profits * systems salvaged(m,t)
"""
model.setObjective((quicksum((quicksum(S * K[t] + A * X[m[0],t] + operatingSystemsVars[m[0]-1][t] * O - L * Y[m[0],t] for t in range(T) if t >= m[1])) for m in M)), GRB.MINIMIZE)

#model.optimize()

 

 

Gurobi creates the following output:

Optimize a model with 300 rows, 440 columns and 6160 nonzeros
Variable types: 0 continuous, 440 integer (0 binary)
Coefficient statistics:
Matrix range [1e+00, 1e+00]
Objective range [5e+02, 1e+04]
Bounds range [0e+00, 0e+00]
RHS range [5e-01, 3e+01]
Presolve removed 10 rows and 0 columns
Presolve time: 0.00s

Explored 0 nodes (0 simplex iterations) in 0.01 seconds
Thread count was 1 (of 8 available processors)

Solution count 0

Model is infeasible or unbounded
Best objective -, best bound -, gap -

IIS computed: 1 constraints, 0 bounds
IIS runtime: 0.01 seconds

0

Comments

5 comments

Date Votes
  • Official comment
    Simranjit Kaur
    • Gurobi Staff
    Edited
    This post is more than three years old. Some information may not be up to date. For current information, please check the Gurobi Documentation or Knowledge Base. If you need more help, please create a new post in the community forum. Or why not try our AI Gurobot?.
  • Eli Towle
    • Gurobi Staff

    Hi William,

    You can write a model to the human-readable LP (.lp) format with the Model.write() method. It looks like you are computing an IIS; you can also write an ILP (.ilp) file to determine which subsystem of constraints is infeasible (in your case, it is a single constraint). Both of these are more helpful if you name your variables and constraints.

    I suspect the issue is with the generateOperatingSystemsCountVars() function that you define. It returns a list of lists, each element of which is zero. Thus, when you add your Equation 4 constraints using operatingSystemsVars, the LHSs all evaluate to zero:

    # The following are all equivalent (Equation 4)
    model.addConstr((quicksum(U[m[0]][t] * operatingSystemsVars[m[0]-1][t] for m in M if t >= m[1])) >= D[t])
    model.addConstr(quicksum(U[m[0]][t] * 0 for m in M if t >= m[1]) >= D[t])
    model.addConstr(0 >= D[t])

    If the RHS is positive, this constraint is infeasible. Also, note that the Equation 4 constraints do not include any of the X or Y model variables.

    Eli

    1
  • William Praat
    Edited
    • Gurobi-versary
    • First Comment
    • First Question

    Eli Towle Thank you for your reply.

    Thank you for explaining what is going wrong. However, does this also mean that it is not possible to let Gurobi change the value of for example K(t) or operatingSystemsVars. I guess it is technically possible by adding them as vars to the model but to me this does not seem logical as they are not decision variables in the original mathematical model.

     

    William

    Update:

    I've made a very basic version of the model and did a test with it and it does find a solution. However, this is not how I expected it to go. I did the test with the following code and with T = 10 for now:

    purchaseDecisionVars = generatePurchaseDecisionVars(T, M)
    salvageDecisionVars = generateSalvageDecisionVars(T, M)

    model = Model("HydrogenSystemOptimization")


    X = model.addVars(purchaseDecisionVars, lb=0, vtype=GRB.INTEGER, name="systems purchases: ")
    Y = model.addVars(salvageDecisionVars, lb=0, vtype=GRB.INTEGER, name="systems salvaged: ")

    for t in range(T):
    for m in M:
    print("T= ", t, " m[0]: ", m[1])
    model.addConstr(quicksum(X[m[0],i] - Y[m[0],i] for i in range(t) if i >= m[1]) >= 0, name="non-zero")

    model.setObjective((quicksum((quicksum(A * X[m[0],t] - L * Y[m[0],t] for t in range(T) if t >= m[1])) for m in M)), GRB.MINIMIZE)

     

     

    This is the result:

    systems purchases: [0,0] -0
    systems purchases: [0,1] -0
    systems purchases: [0,2] -0
    systems purchases: [0,3] -0
    systems purchases: [0,4] -0
    systems purchases: [0,5] -0
    systems purchases: [0,6] -0
    systems purchases: [0,7] -0
    systems purchases: [0,8] -0
    systems purchases: [0,9] -0
    systems purchases: [1,4] -0
    systems purchases: [1,5] -0
    systems purchases: [1,6] -0
    systems purchases: [1,7] -0
    systems purchases: [1,8] -0
    systems purchases: [1,9] -0
    systems purchases: [2,8] -0
    systems purchases: [2,9] -0
    systems salvaged: [0,0] -0
    systems salvaged: [0,1] -0
    systems salvaged: [0,2] -0
    systems salvaged: [0,3] -0
    systems salvaged: [0,4] -0
    systems salvaged: [0,5] -0
    systems salvaged: [0,6] -0
    systems salvaged: [0,7] -0
    systems salvaged: [0,8] -0
    systems salvaged: [0,9] 1e+30
    systems salvaged: [1,4] -0
    systems salvaged: [1,5] -0
    systems salvaged: [1,6] -0
    systems salvaged: [1,7] -0
    systems salvaged: [1,8] -0
    systems salvaged: [1,9] 1e+30
    systems salvaged: [2,8] -0
    systems salvaged: [2,9] 1e+30
    Obj: total costs: -1.5e+33

     

    So for some reason the constraint I added does not work for the last t but I don't understand why. Also interesting, if I change the non-zero constraint from 0 to 1 (simulating some demand) it cannot find a solution.

    0
  • Eli Towle
    • Gurobi Staff
    Edited

    Hi William,

    Gurobi searches for an optimal assignment of values to the model variables. If the elements of operatingSystemsVars aren't model variables, then they are fixed constants/coefficients that will not be changed by Gurobi.

    Have you written out an LP file to inspect the model visually? With this code, a lot of your "non-zero" constraints are \( 0 \geq 0 \). Changing the RHS to 1 will result in a constraint that is never satisfied, and thus an infeasible model. Note that every constraint should include a variable. If this isn't the case, you end up with a constraint that is always false (so the model is automatically infeasible), or a constraint that is always true (meaning the constraint doesn't affect the feasible region).

    For t=9, you add constraints for every i in range(t). However, range(9) is the numbers 0 through 8. Thus, "systems salvaged" variables with a second index of t=9 do not appear in any constraints. Perhaps you should be using range(t+1) in the quicksum() function of the "non-zero" constraints?

    Eli

     

    1
  • William Praat
    • Gurobi-versary
    • First Comment
    • First Question

    Dear Eli,

    Thank you for this clarification. Your description gave me a much better understanding of how Gurobi works and I was able to get my model working. There were some more issues I had to resolve but I got my model working as intended and it is generating results that make sense.

    Especially the LP file helped me to find the problems, thank you for pointing that out.

    Again, thank you for the support!

    William

    0

Post is closed for comments.

Didn't find what you were looking for?

New post