Juice Problem ILP Formulation
Juice Problem ILP (GurobiPy)
Author: Ben Rosenberg
Imports
We begin by importing some relevant libraries. We import GurobiPy
to formulate and solve the ILP, and we import time
to time the entire process of supplying the constraints solving the ILP. Lastly, we'll use pandas
to work with the input data to the ILP.
import gurobipy as gp
from gurobipy import GRB
import time
import pandas as pd
Problem definition
The "Juice Problem" is pretty vague. This is basically just a simple example of an ILP with constraints based on nutritional values of various ingredients. We'll define an input here:
Ingredient | Calories | Fat | Cholesterol | Sodium | Sugar |
---|---|---|---|---|---|
A | 110 | 4 | 30 | 340 | 15 |
B | 150 | 8 | 40 | 120 | 35 |
C | 250 | 2 | 60 | 450 | 40 |
D | 270 | 10 | 20 | 240 | 60 |
E | 350 | 15 | 70 | 760 | 45 |
And let's make some constraints. Our juice can have:
- No fewer than 2000 units of calories
- No more than 70 units of fat
- No more than 750 units of cholesterol
- No more than 3000 units of sodium
Furthermore, let's say that the ratio of cholesterol to sodium must be at most 1:10.
Lastly, let's make the objective to minimize the amount of sugar.
Input entry
We'll use pandas
to store our input in a more manageable form (a DataFrame). This allows this notebook to also serve as a minimal example of how one might use the pandas
library.
# define ingredients using a dictionary
data = {
'Ingredients' : ['A', 'B', 'C', 'D', 'E'],
'Calories' : [110, 150, 250, 270, 350],
'Fat' : [4, 8, 2, 10, 15],
'Cholesterol' : [30, 40, 60, 20, 70],
'Sodium' : [340, 120, 450, 240, 760],
'Sugar' : [15, 35, 40, 60, 45]
}
df = pd.DataFrame.from_dict(data).set_index('Ingredients')
df.head()
print(df['Calories']['A'])
110
start_time = time.time()
model = model = gp.Model()
# set output level to max
model.Params.TuneOutput = 3
# add variables for each ingredient
I = model.addVars(df.index, vtype=GRB.CONTINUOUS, name='I')
# add constraint on calories
model.addConstr(gp.quicksum(I[ingredient] * df['Calories'][ingredient] for ingredient in df.index) >= 2000)
# add constraint on fat
model.addConstr(gp.quicksum(I[ingredient] * df['Fat'][ingredient] for ingredient in df.index) <= 70)
# add constraint on cholesterol
model.addConstr(gp.quicksum(I[ingredient] * df['Cholesterol'][ingredient] for ingredient in df.index) <= 750)
# add constraint on sodium
model.addConstr(gp.quicksum(I[ingredient] * df['Sodium'][ingredient] for ingredient in df.index) <= 3000)
# add constraint on cholesterol:sodium ratio
model.addConstr(gp.quicksum(I[ingredient] * df['Cholesterol'][ingredient] for ingredient in df.index) * 10
<= gp.quicksum(I[ingredient] * df['Sodium'][ingredient] for ingredient in df.index))
# set objective
model.setObjective(gp.quicksum(I[ingredient] * df['Sugar'][ingredient] for ingredient in df.index))
model.optimize()
end_time = time.time()
diff = time.gmtime(end_time - start_time)
print('\n[Total time used: {} minutes, {} seconds]'.format(diff.tm_min, diff.tm_sec))
try:
print(f'\nObjective value found: {model.objVal}')
except AttributeError as e:
print(f'\nCould not find an objective value. \nTraceback:\n\t{e}')
Set parameter Username Academic license - for non-commercial use only - expires 2023-10-13 Set parameter TuneOutput to value 3 Gurobi Optimizer version 9.5.2 build v9.5.2rc0 (win64) Thread count: 4 physical cores, 8 logical processors, using up to 8 threads Optimize a model with 5 rows, 5 columns and 25 nonzeros Model fingerprint: 0xcf744ec6 Coefficient statistics: Matrix range [2e+00, 8e+02] Objective range [2e+01, 6e+01] Bounds range [0e+00, 0e+00] RHS range [7e+01, 3e+03] Presolve time: 0.01s Presolved: 5 rows, 5 columns, 25 nonzeros Iteration Objective Primal Inf. Dual Inf. Time 0 0.0000000e+00 1.250000e+02 0.000000e+00 0s 3 3.5644712e+02 0.000000e+00 0.000000e+00 0s Solved in 3 iterations and 0.02 seconds (0.00 work units) Optimal objective 3.564471197e+02 [Total time used: 0 minutes, 0 seconds] Objective value found: 356.44711968958313
Now that we have our objective, let's print out our solution values:
for ingredient in df.index:
print('{} units of {}'.format(
I[ingredient].x, ingredient
))
0.0 units of A 0.0 units of B 1.1600835737737538 units of C 3.5668092727091834 units of D 2.1341160083573776 units of E
Integer Program
Say that we have the additional constraint that we can only use integer amounts of each ingredient. This means we need to use integer variables (m.IntVar(...)
) rather than continuous ones (m.NumVar(...)
):
start_time = time.time()
model = model = gp.Model()
# set output level to max
model.Params.TuneOutput = 3
# add variables for each ingredient
I = model.addVars(df.index, vtype=GRB.INTEGER, name='I') # <-- change
# add constraint on calories
model.addConstr(gp.quicksum(I[ingredient] * df['Calories'][ingredient] for ingredient in df.index) >= 2000)
# add constraint on fat
model.addConstr(gp.quicksum(I[ingredient] * df['Fat'][ingredient] for ingredient in df.index) <= 70)
# add constraint on cholesterol
model.addConstr(gp.quicksum(I[ingredient] * df['Cholesterol'][ingredient] for ingredient in df.index) <= 750)
# add constraint on sodium
model.addConstr(gp.quicksum(I[ingredient] * df['Sodium'][ingredient] for ingredient in df.index) <= 3000)
# add constraint on cholesterol:sodium ratio
model.addConstr(gp.quicksum(I[ingredient] * df['Cholesterol'][ingredient] for ingredient in df.index) * 10
<= gp.quicksum(I[ingredient] * df['Sodium'][ingredient] for ingredient in df.index))
# set objective
model.setObjective(gp.quicksum(I[ingredient] * df['Sugar'][ingredient] for ingredient in df.index))
model.optimize()
end_time = time.time()
diff = time.gmtime(end_time - start_time)
print('\n[Total time used: {} minutes, {} seconds]'.format(diff.tm_min, diff.tm_sec))
try:
print(f'\nObjective value found: {model.objVal}')
except AttributeError as e:
print(f'\nCould not find an objective value. \nTraceback:\n\t{e}')
Set parameter TuneOutput to value 3 Gurobi Optimizer version 9.5.2 build v9.5.2rc0 (win64) Thread count: 4 physical cores, 8 logical processors, using up to 8 threads Optimize a model with 5 rows, 5 columns and 25 nonzeros Model fingerprint: 0x9f3dd579 Variable types: 0 continuous, 5 integer (0 binary) Coefficient statistics: Matrix range [2e+00, 8e+02] Objective range [2e+01, 6e+01] Bounds range [0e+00, 0e+00] RHS range [7e+01, 3e+03] Presolve removed 1 rows and 0 columns Presolve time: 0.01s Presolved: 4 rows, 5 columns, 20 nonzeros Variable types: 0 continuous, 5 integer (0 binary) Found heuristic solution: objective 430.0000000 Root relaxation: objective 3.564471e+02, 3 iterations, 0.00 seconds (0.00 work units) Nodes | Current Node | Objective Bounds | Work Expl Unexpl | Obj Depth IntInf | Incumbent BestBd Gap | It/Node Time 0 0 356.44712 0 3 430.00000 356.44712 17.1% - 0s 0 0 392.28571 0 3 430.00000 392.28571 8.77% - 0s Cutting planes: Gomory: 1 StrongCG: 2 Explored 1 nodes (8 simplex iterations) in 0.04 seconds (0.00 work units) Thread count was 8 (of 8 available processors) Solution count 1: 430 Optimal solution found (tolerance 1.00e-04) Best objective 4.300000000000e+02, best bound 4.300000000000e+02, gap 0.0000% [Total time used: 0 minutes, 0 seconds] Objective value found: 430.0
And here are the solution values:
for ingredient in df.index:
print('{} units of {}'.format(
I[ingredient].x, ingredient
))
2.0 units of A -0.0 units of B 1.0 units of C 6.0 units of D -0.0 units of E
Note how adding the constraint that the amounts of each ingredient must be integer results in a higher (worse) objective value. In general, adding constraints to a problem will result in an equal or worse result than the previously obtained one.
Juice Problem ILP
Author: Ben Rosenberg
Imports
We begin by importing some relevant libraries. We import ortools
to formulate and solve the ILP, and we import time
to time the entire process of supplying the constraints solving the ILP. Lastly, we'll use pandas
to work with the input data to the ILP.
from ortools.linear_solver import pywraplp as OR
import time
import pandas as pd
Problem definition
The "Juice Problem" is pretty vague. This is basically just a simple example of an ILP with constraints based on nutritional values of various ingredients. We'll define an input here:
Ingredient | Calories | Fat | Cholesterol | Sodium | Sugar |
---|---|---|---|---|---|
A | 110 | 4 | 30 | 340 | 15 |
B | 150 | 8 | 40 | 120 | 35 |
C | 250 | 2 | 60 | 450 | 40 |
D | 270 | 10 | 20 | 240 | 60 |
E | 350 | 15 | 70 | 760 | 45 |
And let's make some constraints. Our juice can have:
- No fewer than 2000 units of calories
- No more than 70 units of fat
- No more than 750 units of cholesterol
- No more than 3000 units of sodium
Furthermore, let's say that the ratio of cholesterol to sodium must be at most 1:10.
Lastly, let's make the objective to minimize the amount of sugar.
Input entry
We'll use pandas
to store our input in a more manageable form (a DataFrame). This allows this notebook to also serve as a minimal example of how one might use the pandas
library.
# define ingredients using a dictionary
data = {
'Ingredients' : ['A', 'B', 'C', 'D', 'E'],
'Calories' : [110, 150, 250, 270, 350],
'Fat' : [4, 8, 2, 10, 15],
'Cholesterol' : [30, 40, 60, 20, 70],
'Sodium' : [340, 120, 450, 240, 760],
'Sugar' : [15, 35, 40, 60, 45]
}
df = pd.DataFrame.from_dict(data).set_index('Ingredients')
df.head()
print(df['Calories']['A'])
110
start_time = time.time()
m = OR.Solver('Juice Problem', OR.Solver.CBC_MIXED_INTEGER_PROGRAMMING)
# add variables for each ingredient
I = {}
for ingredient in df.index:
I[ingredient] = m.NumVar(0, m.infinity(), ingredient)
# add constraint on calories
m.Add(sum(I[ingredient] * df['Calories'][ingredient] for ingredient in df.index) >= 2000)
# add constraint on fat
m.Add(sum(I[ingredient] * df['Fat'][ingredient] for ingredient in df.index) <= 70)
# add constraint on cholesterol
m.Add(sum(I[ingredient] * df['Cholesterol'][ingredient] for ingredient in df.index) <= 750)
# add constraint on sodium
m.Add(sum(I[ingredient] * df['Sodium'][ingredient] for ingredient in df.index) <= 3000)
# add constraint on cholesterol:sodium ratio
m.Add(sum(I[ingredient] * df['Cholesterol'][ingredient] for ingredient in df.index) * 10
<= sum(I[ingredient] * df['Sodium'][ingredient] for ingredient in df.index))
# set objective
m.Minimize(sum(I[ingredient] * df['Sugar'][ingredient] for ingredient in df.index))
m.Solve()
end_time = time.time()
diff = time.gmtime(end_time - start_time)
print('\n[Total time used: {} minutes, {} seconds]'.format(diff.tm_min, diff.tm_sec))
print('Objective:', m.Objective().Value())
[Total time used: 0 minutes, 0 seconds] Objective: 356.44711968958313
Now that we have our objective, let's print out our solution values:
for ingredient in df.index:
print('{} units of {}'.format(
I[ingredient].solution_value(), ingredient
))
0.0 units of A 0.0 units of B 1.1600835737737498 units of C 3.5668092727091842 units of D 2.1341160083573794 units of E
Integer Program
Say that we have the additional constraint that we can only use integer amounts of each ingredient. This means we need to use integer variables (m.IntVar(...)
) rather than continuous ones (m.NumVar(...)
):
start_time = time.time()
m = OR.Solver('Juice Problem', OR.Solver.CBC_MIXED_INTEGER_PROGRAMMING)
# add variables for each ingredient
I = {}
for ingredient in df.index:
I[ingredient] = m.IntVar(0, m.infinity(), ingredient) # <-- change
# add constraint on calories
m.Add(sum(I[ingredient] * df['Calories'][ingredient] for ingredient in df.index) >= 2000)
# add constraint on fat
m.Add(sum(I[ingredient] * df['Fat'][ingredient] for ingredient in df.index) <= 70)
# add constraint on cholesterol
m.Add(sum(I[ingredient] * df['Cholesterol'][ingredient] for ingredient in df.index) <= 750)
# add constraint on sodium
m.Add(sum(I[ingredient] * df['Sodium'][ingredient] for ingredient in df.index) <= 3000)
# add constraint on cholesterol:sodium ratio
m.Add(sum(I[ingredient] * df['Cholesterol'][ingredient] for ingredient in df.index) * 10
<= sum(I[ingredient] * df['Sodium'][ingredient] for ingredient in df.index))
# set objective
m.Minimize(sum(I[ingredient] * df['Sugar'][ingredient] for ingredient in df.index))
m.Solve()
end_time = time.time()
diff = time.gmtime(end_time - start_time)
print('\n[Total time used: {} minutes, {} seconds]'.format(diff.tm_min, diff.tm_sec))
print('Objective:', m.Objective().Value())
[Total time used: 0 minutes, 0 seconds] Objective: 430.0
And here are the solution values:
for ingredient in df.index:
print('{} units of {}'.format(
I[ingredient].solution_value(), ingredient
))
2.0 units of A 0.0 units of B 1.0 units of C 6.0 units of D 0.0 units of E
Note how adding the constraint that the amounts of each ingredient must be integer results in a higher (worse) objective value. In general, adding constraints to a problem will result in an equal or worse result than the previously obtained one.