WS 2.5: Profit vs Planet¶

No description has been provided for this image No description has been provided for this image

CEGM1000 MUDE: Week 2.5, Optimization. For: December 11, 2024

Overview¶

A civil engineering company wants to decide on the projects that they should do. Their objective is to minimize the environmental impact of their projects while making enough profit to keep the company running.

They have a portfolio of 6 possible projects to invest in, where A, B , and C are new infrastructure projects (so-called type 1), and D, E, F are refurbishment projects (so-called type 2).

The environmental impact of each project is given by $I_i$ where $i \in [1,(...),6]$ is the index of the project. $I_i=[90,45,78,123,48,60]$

The profit of each project is given by $P_i$ where $i\in [1,(...),6]$ is the index of the project: $P_i=[120,65,99,110,33,99]$

The company wants to do 3 out of the 6 projects, therefore please formulate the mathematical program that allows solving the problem, also knowing that the projects of type 2 must be at least as many as the ones of type 1 and that the profit of all projects together must be greater or equal than $250$ ($\beta$)

You are not allowed to use ChatGPT for this task otherwise you won’t learn ;)

Task 1: Writting the mathematical formulation

Write down every formulation and constrain that is relevant to solve this optimization problem.

Solution $$Min(Z) \sum_{i=1}^{6} x_i I_i $$

subject to

$$ \sum_{i=1}^{6} x_i = 3 $$$$ \sum_{i=4}^{6} x_i \ge \sum_{i=1}^{3} x_i $$$$ \sum_{i=1}^{6} x_i P_i \ge β $$$$x\in(1,0), i\in(1,....,6)$$

Extra Challenge (Task 6) $$ \sum_{i=1}^{6} x_i I_i \le x_1 γ + (1-x_1) M $$

where M is a sufficiently big number such as 9999

Task 2: Setting up the software

We'll continue using Gurobi this week, which you've set up in last week's PA. We'll use some other special packages as well. Therefore, be sure to use the special conda environment created for this week.

In [ ]:
import gurobipy as gp

Task 3: Setting up the problem

Define any variables you might need to setup your model.

In [ ]:
# YOUR_CODE_HERE

# SOLUTION
# Project data
I = [90, 45, 78, 123, 48, 60]  # Environmental impact
P = [120, 65, 99, 110, 33, 99]  # Profit

# Minimum required profit
beta = 250
M = 100000

# Number of projects and types
num_projects = len(I)
num_type1_projects = 3
num_type2_projects = num_projects - num_type1_projects

Create model with Gurobi¶

Remember that examples of using Gurobi to create and optimize a model are provided in the online textbook, and generally consist of the following steps (the first instantiates a class and the rest are executed as methods of the class):

  1. Define the model (instantiate the class)
  2. Define variables
  3. Define objective function
  4. Add constraints
  5. Optimize the model

Remember, you can always ask for help to understand a function of gurobi

help(gurobipy.model.addVars)

Task 4: Create the Gurobi model

Create the Gurobi model, set your decision variables, your function and your constrains. Take a look at the book for an example implementation in Python if you don't know where to start.

Solution.

Each line of code does the following:

  1. Create a Gurobi model
  2. Add binary variables
  3. Objective function: Minimize environmental impact
  4. Constraint: Select exactly 3 projects
  5. Constraint: Number of type 2 projects must be at least as many as type 1 projects selected
  6. Constraint: Minimum profit requirement
  7. Optimize the model

In [ ]:
# YOUR_CODE_HERE

# SOLUTION
model = gp.Model("Project_Selection")
x = model.addVars(num_projects, vtype=gp.GRB.BINARY, name="x")
model.setObjective(sum(I[i] * x[i] for i in range(num_projects)),
                   gp.GRB.MINIMIZE)
model.addConstr(x.sum() == 3, "Select_Projects")
model.addConstr((sum(x[i] for i in range(num_type2_projects, num_projects))
                 - sum(x[i] for i in range(num_type1_projects)) >= 0),
                 "Type_Constraint")
model.addConstr(sum(P[i] * x[i] for i in range(num_projects)) >= beta,
                "Minimum_Profit")
model.optimize()
Gurobi Optimizer version 9.5.1 build v9.5.1rc2 (win64)
Thread count: 4 physical cores, 8 logical processors, using up to 8 threads
Optimize a model with 3 rows, 6 columns and 18 nonzeros
Model fingerprint: 0xa9bc0e77
Variable types: 0 continuous, 6 integer (6 binary)
Coefficient statistics:
  Matrix range     [1e+00, 1e+02]
  Objective range  [5e+01, 1e+02]
  Bounds range     [1e+00, 1e+00]
  RHS range        [3e+00, 3e+02]
Presolve time: 0.00s
Presolved: 3 rows, 6 columns, 17 nonzeros
Variable types: 0 continuous, 6 integer (6 binary)
Found heuristic solution: objective 228.0000000
Found heuristic solution: objective 198.0000000

Root relaxation: objective 1.855909e+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 infeasible    0       198.00000  198.00000  0.00%     -    0s

Explored 1 nodes (3 simplex iterations) in 0.01 seconds (0.00 work units)
Thread count was 8 (of 8 available processors)

Solution count 2: 198 228 

Optimal solution found (tolerance 1.00e-04)
Best objective 1.980000000000e+02, best bound 1.980000000000e+02, gap 0.0000%

Task 5: Display your results

Display the model in a good way to interpret and print the solution of the optimization.

In [ ]:
# YOUR_CODE_HERE

# SOLUTION
print("Model structure:")        
# see the model that you have built in a nice why to interpret
model.display()  

if model.status == gp.GRB.OPTIMAL:
    print("Optimal Solution:")
    for i in range(num_projects):
        if x[i].x > 0.9:
            print(f"Project {i+1}: Selected")
else:
    print("No optimal solution found.")

    
print("Optimal Objective function Value", model.objVal)    
Model structure:
Minimize
  <gurobi.LinExpr: 90.0 x[0] + 45.0 x[1] + 78.0 x[2] + 123.0 x[3] + 48.0 x[4] + 60.0 x[5]>
Subject To
  Select_Projects: <gurobi.LinExpr: x[0] + x[1] + x[2] + x[3] + x[4] + x[5]> = 3
Type_Constraint: <gurobi.LinExpr: -1.0 x[0] + -1.0 x[1] + -1.0 x[2] + x[3] + x[4] +
 x[5]> >= 0
Minimum_Profit: <gurobi.LinExpr: 120.0 x[0] + 65.0 x[1] + 99.0 x[2] + 110.0 x[3] + 33.0
 x[4] + 99.0 x[5]> >= 250
Binaries
  ['x[0]', 'x[1]', 'x[2]', 'x[3]', 'x[4]', 'x[5]']
Optimal Solution:
Project 1: Selected
Project 5: Selected
Project 6: Selected
Optimal Objective function Value 198.0

Task 6: Additional constraint

Solve the model with an additional constraint: if project 1 is done then the impact of all projects together should be lower than $\gamma$ with $\gamma=130$.

In the first cell you should add the constraint, then in a second cell optimize the model.

In [ ]:
# YOUR_CODE_HERE

# SOLUTION
gamma = 130
model.addConstr((sum(I[i] * x[i] for i in range(num_projects))
                 <= gamma * x[0]+ M * (1 - x[0])),
                 "Impact_Constraint") 
Out[ ]:
<gurobi.Constr *Awaiting Model Update*>
In [ ]:
# YOUR_CODE_HERE

# SOLUTION
model.optimize()

print("Model structure:")        
# see the model that you have built in a nice why to interpret
model.display()  

# Display the solution
if model.status == gp.GRB.OPTIMAL:
    print("Optimal Solution:")
    for i in range(num_projects):
        if x[i].x > 0.9:
            print(f"Project {i+1}: Selected")
else:
    print("No optimal solution found.")

    
print("Optimal Objective function Value", model.objVal)   
Gurobi Optimizer version 9.5.1 build v9.5.1rc2 (win64)
Thread count: 4 physical cores, 8 logical processors, using up to 8 threads
Optimize a model with 4 rows, 6 columns and 24 nonzeros
Model fingerprint: 0x37d09666
Variable types: 0 continuous, 6 integer (6 binary)
Coefficient statistics:
  Matrix range     [1e+00, 1e+05]
  Objective range  [5e+01, 1e+02]
  Bounds range     [1e+00, 1e+00]
  RHS range        [3e+00, 1e+05]

MIP start from previous solve did not produce a new incumbent solution
MIP start from previous solve violates constraint Impact_Constraint by 68.000000000

Presolve removed 4 rows and 6 columns
Presolve time: 0.00s
Presolve: All rows and columns removed

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

Solution count 1: 228 

Optimal solution found (tolerance 1.00e-04)
Best objective 2.280000000000e+02, best bound 2.280000000000e+02, gap 0.0000%
Model structure:
Minimize
  <gurobi.LinExpr: 90.0 x[0] + 45.0 x[1] + 78.0 x[2] + 123.0 x[3] + 48.0 x[4] + 60.0 x[5]>
Subject To
  Select_Projects: <gurobi.LinExpr: x[0] + x[1] + x[2] + x[3] + x[4] + x[5]> = 3
Type_Constraint: <gurobi.LinExpr: -1.0 x[0] + -1.0 x[1] + -1.0 x[2] + x[3] + x[4] +
 x[5]> >= 0
Minimum_Profit: <gurobi.LinExpr: 120.0 x[0] + 65.0 x[1] + 99.0 x[2] + 110.0 x[3] + 33.0
 x[4] + 99.0 x[5]> >= 250
Impact_Constraint: <gurobi.LinExpr: 99960.0 x[0] + 45.0 x[1] + 78.0 x[2] + 123.0 x[3] +
 48.0 x[4] + 60.0 x[5]> <= 100000
Binaries
  ['x[0]', 'x[1]', 'x[2]', 'x[3]', 'x[4]', 'x[5]']
Optimal Solution:
Project 2: Selected
Project 4: Selected
Project 6: Selected
Optimal Objective function Value 228.0

End of notebook.

Creative Commons License TU Delft MUDE

© Copyright 2024 MUDE TU Delft. This work is licensed under a CC BY 4.0 License.